Check the Geospatial Data Science with Julia book:
https://juliaearth.github.io/geospatial-data-science-with-julia
Use a GeoTable.
Check the Geospatial Data Science with Julia book:
https://juliaearth.github.io/geospatial-data-science-with-julia
Use a GeoTable.
Thanks for the tip ! I just went through the book, and I feel like I’m missing some things.
georef
. In my current practice, I have a dedicated Mesh
-like type that stores all the info and then I use a custom recipe based on various Makie
basic recipes to do the job. What would be the Meshes
to achieve that ? Dispatch on the column names of the geotable
?The geometry column can contain any geometry, not just points.
The recipes are dispatched on Domain types which are collections of geometries.
Please let us know if something is not clear. Better ask on Zulip to avoid spamming everyone in this thread.
@eliascarv and I are super excited to share a new release of the project, after months of hard work to formalize coordinate reference systems (CRS) with Julia types. This release comes with many breaking changes, which are needed for real-world usage of meshes in industrial applications.
The main innovation that puts us years ahead of other mesh/GIS software out there is the ability to abstract CRS information in type parameters, without increasing memory requirements.
Every point now has a well-defined set of coordinates with Unitful.jl units thanks to CoordRefSystems.jl, and we can carry this important information forward in algorithms without runtime nor storage penalties:
julia> Point(1,2,3)
Point with Cartesian{NoDatum} coordinates
├─ x: 1.0 m
├─ y: 2.0 m
└─ z: 3.0 m
julia> using CoordRefSystems
julia> Point(Polar(0,0))
Point with Polar{NoDatum} coordinates
├─ ρ: 0.0 m
└─ ϕ: 0.0 rad
All geometries and domain types also have this CRS information, which is useful to dispatch specialized algorithms on non-Cartesian coordinates.
We are effectively solving the “implicit” approach to CRS handling where users are responsible for magically finding the CRS of their data.
coordinates
function is no longer present. Use the to
function instead to return the Vec
from the origin to a given Point
.T
. They still have a Dim
parameter as in Geometry{Dim}
.coords(point)
for that. Most times, the to
function is what is needed.The roadmap for the following weeks consists of updating the rest of the GeoStats.jl stack to accommodate this new release with CRS support, and fix miscellaneous issues that are listed here.
This is great to see! What was the motivation for the naming of to
function? It seems a bit strange to me. I would have never guessed it returns the vector (Vec
) from the origin. Why not just allow people to use Vec
?
This function takes any point and returns the geometric vector point - origin
, i.e. our Vec
. The motivation in the name comes from
origin to point
.
This is a very specific operation that is not very useful in curved manifolds for example, where the vector connecting the origin of the CRS to the point does not belong to the surface. We will add more features like this in the future, and associated functions.
The name to
was the shortest name that we found. Given that it is used a lot in implementations, it is convenient to keep it short.
Awesome news! Even though I’m not the target audience for this feature, I’m happy to see the active development of the package.
I’m wondering if there are any plans for Makie.jl to use Meshes.jl under the hood? Or Meshes.jl doesn’t want to be the new GeometryBasics.jl package.
I can answer from the Meshes.jl side that we don’t see any issue with this idea of providing the geometric processing infrastructure for Makie.jl. Both projects would certainly benefit from tons of features, and human resources would be used more wisely.
v0.43 automatically adds units to all the points, which are by default [m] . Is it possible to use a customized unit system, or simply turn off the units?
Hi @henry2004y , the following is true now:
Point(xyz...) = Point(Cartesian(xyz...))
This means that Point(1, 2, 3)
is syntactic sugar to construct Cartesian
coordinates from CoordRefSystems.jl If the coordinates don’t have appropriate units, then [m] will be added as you said. But if the coordinates already have length (L) units, they will pass just fine:
julia> using Meshes
julia> using Unitful
julia> Point(1u"km", 2u"km", 3u"ft")
Point with Cartesian{NoDatum} coordinates
├─ x: 1000.0 m
├─ y: 2000.0 m
└─ z: 0.9144 m
We don’t allow unit-less coordinates anymore because that is the root of all evil in geospatial software. Pretend that the default [m] units always existed implicitly, and that they are just explicit now. Algorithms will be more robust as a consequence, and everything else such as distances, measures, etc. will have a clear physical interpretation:
julia> area(rand(Triangle{2}))
0.07232380768033886 m^2
if the coordinates already have length (L) units, they will pass just fine
Meshes can be useful in other spaces
(for instance in reciprocal space, units may be rad / nm
).
julia> using Meshes
julia> using Unitful
julia> Point(1u"rad/nm", 2u"rad/nm")
ERROR: ArgumentError: invalid units for coordinates, please check the documentation
Check CoordRefSystems.jl for other coordinate types with angle units.
Point(Polar(…)) like in the example above.
Proj
for changing the CRS of geometries, domains and geotablesShadow
for computing the “shadow” of a 3D domain onto XY, XZ or YZLenghtUnit
for changing the lenght unit of geometries and domainsDelaunayTriangulation
of geometriesDelaunayTesselation
of point setsVoronoiTesselation
of point sets(Multi-)polygons with more than 5000 vertices are now discretized with DelaunayTriangulation
by default. DehnTriangulation
is preferred for smaller polygons because it is faster.
FIST
renamed to HeldTriangulation
Dehn1899
renamed to DehnTriangulation
Quick updates:
Addition of Wedge
geometry as a 3-Polytope
Discretization of Cylinder
into Hexahedron
and Wedge
geometries
sideof(point, ring)
and point ∈ poly
with Hao et al. algorithm, which we implemented with multiple-threads. Speedup of ~12x. The algorithm is quite robust, and can handle all sorts of degenerate cases in 2D polygonal areas:rad
units:julia> ∠(Vec(1, 0), Vec(0, 1))
1.5707963267948966 rad
Aratu
, SAD69
, SAD96
and SIRGAS2000
datums:This release introduces support for non-Euclidean manifolds, which means that geometries can be “overlaid” over surfaces such as ellipsoids of revolution representing the Earth.
To give an example, consider the construction of a triangle using LatLon
coordinates for the vertices:
julia> using Meshes
julia> using CoordRefSystems
julia> triangle = Triangle(Point(LatLon(0, 0)), Point(LatLon(0, 90)), Point(LatLon(90, 0)))
Triangle
├─ Point(lat: 0.0°, lon: 0.0°)
├─ Point(lat: 0.0°, lon: 90.0°)
└─ Point(lat: 90.0°, lon: 0.0°)
julia> manifold(triangle)
🌐
The manifold is automatically detected from the CRS of the points at compile-time, leading to a geodesic triangle over the globe.
Since all our geometry and domain types have an explicit manifold, we can now dispatch specialized algorithms from geodesic geometry.
This is the last major refactor that we had planned before v1.0. We will test the design over the rest of the year before making the first major release official.
We invite all potential users and collaborators to give it a try. Now is the time to report issues and suggest breaking changes.
For those interested in geographic applications:
For those interested in differential forms and integrals, @mike.ingold @JoshuaLampert and @kylebeggs are working on AD rules for Meshes.jl geometries: