[ANN] Meshes.jl - Computational Geometry in Julia

Package Announcements
, , , ,
41

Check the Geospatial Data Science with Julia book:

https://juliaearth.github.io/geospatial-data-science-with-julia

Use a GeoTable.

3 Likes
42

Thanks for the tip ! I just went through the book, and I feel like I’m missing some things.

  1. It feels like everything is based on the data values at geometry points, but no data can be linked directly with the elements, is that so ?
  2. How do you dispatch plotting recipes based on a given 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 ?
43

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.

1 Like
44

Meshes.jl v0.43

@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.

Breaking changes

  • The coordinates function is no longer present. Use the to function instead to return the Vec from the origin to a given Point.
  • Geometry and domain types no longer have a type parameter T . They still have a Dim parameter as in Geometry{Dim}.
  • All points and vectors have unitful coordinates. End-users are not supposed to handle these coordinates directly, but we provide a function coords(point) for that. Most times, the to function is what is needed.

Roadmap

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.

16 Likes
45

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?

1 Like
46

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.

1 Like
47

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.

2 Likes
48

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.

5 Likes
49

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?

1 Like
50

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
1 Like
51

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
52

Check CoordRefSystems.jl for other coordinate types with angle units.

Point(Polar(…)) like in the example above.

53

Meshes.jl v0.45

New transforms:

New discretization/tesselation:

Breaking changes

(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
10 Likes
54

Quick updates:

  • Addition of Wedge geometry as a 3-Polytope

  • Discretization of Cylinder into Hexahedron and Wedge geometries

image

  • New 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:

image

  • Angles with explicit rad units:
julia> ∠(Vec(1, 0), Vec(0, 1))
1.5707963267948966 rad
  • Addition of various CRS for georeferencing data in Brazil, including the Aratu, SAD69, SAD96 and SIRGAS2000 datums:
9 Likes
55

Meshes.jl v0.47

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.

15 Likes
56

For those interested in geographic applications:

5 Likes