dms polygon operators

Back to Geometric functions and polygon operators

The dms_ operators are the polygon set operations that GeoDMS implements itself, as a snap rounding sweep, next to the four geometry libraries (geos, boost geometry, boost polygon and cgal). They differ from those in the one respect that matters when the data is not perfect: they accept polygons that are not valid. Rings may cross themselves, touch themselves, overlap each other, touch each other or fail to close, and an optional grid size merges vertices and slivers that are too close to be meant apart.

Since GeoDMS 20.20.0.

the operators

Element by element, on two polygon arrays:

Each takes two polygon data items of the same Value type and Values unit, over the same Domain unit or with one of them a Parameter, and each takes an optional third argument, the grid size.

Cleaning, dissolving and splitting:

  • dms_polygon - read each element under the even-odd rule and write it back as a valid polygon
  • dms_union_polygon - dissolve all polygons into one, optionally grouped by a partition attribute
  • dms_split_polygon - split multi-polygons into their single polygons, each a separate domain entry
  • dms_split_union_polygon - dissolve, then split into single polygons

Overlay and connectivity:

Growing and shrinking:

  • dms_minkowski_sum - grow each geometry by a kernel: the kernel as a polygon argument, or one of six named shapes at a given size
  • dms_minkowski_difference - shrink each geometry by a kernel; the eroding counterpart

All six polygon value types are supported throughout: Spoint, Ipoint, Wpoint, Upoint, Fpoint and Dpoint. That is wider than any of the libraries, where geos is dpoint-only and boost polygon integer-only.

what an operand may look like: the even-odd rule

The Polygon value of one element is a sequence of points, with the rings of a multi polygon strung together by zero-width corridors. The dms_ operators read that sequence as one closed walk: consecutive points are line segments, and the last point is connected back to the first. Which part of the plane is then inside is decided by the even-odd rule:

a point lies inside the polygon when a ray from that point to infinity crosses the boundary an odd number of times.

Three consequences a modeller can use:

  • Ring orientation does not matter on input. A shell given counter-clockwise, a hole given clockwise, or a feature whose rings are only partly flipped, are all read the same way. This is what has_correct_winding and fix_winding_order exist for on the other families.
  • A self-intersecting ring has a defined meaning. A ring in the shape of a bow tie is two triangles, since the middle is crossed twice. The usual definition of a valid polygon, the one the geometry libraries implement, leaves such a ring without an interior, which is why they refuse it or repair it first.
  • An unclosed ring is closed and a corridor walked twice cancels, because a segment traversed in both directions does not change the crossing parity.

The result is always a valid polygon: clockwise shells and counter-clockwise holes, the order that Area is positive for.

the grid, and the tolerance

Every vertex of the input and every computed intersection is snapped to a grid, and that grid is the tolerance. Two vertices in the same cell become one vertex, and a sliver thinner than a cell collapses and is dropped: polygon in, polygon out, so a collapsed line or point is not part of the result.

form the grid use
two arguments, integer coordinates the integer grid itself no vertex is moved and every vertex of the result is an integer point
two arguments, float coordinates a power of two, the coarser of the last-place unit of the largest coordinate and the extent divided by 2^36: the extent of the element pair for the binary operators, and for a dissolve the declared range of the Values unit, or the extent of all the data when the unit declares none far below the precision of geographic data, so the answer is the exact one
three arguments the given grid size in coordinate units the tolerance at which near-coincident vertices merge and slivers disappear

The grid size of the three-argument form is a Parameter, one number for the whole attribute, so that two neighbours snap a shared boundary the same way. It must be positive, a whole number for integer coordinates, and coarse enough for the operands: at most 2^36 cells in either direction over their joint extent, which a finer grid is refused for by name.

A dissolve, dms_union_polygon or dms_split_union_polygon, derives its grid once for the whole call rather than per element pair: from the declared range of the values unit when it has one, and from the extent of all the elements otherwise. Every element and every intermediate then lands on one lattice, so the result depends neither on the order of the elements nor on how the attribute is tiled. A values unit ranged over the whole of the Netherlands in metres gives a cell of about 8 micrometres.

what the result looks like

  • Rings carry no vertex that lies on the straight line between its neighbours, so a shared edge that disappears in a union takes its collinear corners with it.
  • Shells and holes start at their lexicographically first vertex and are written in that order, so the result does not depend on the order of the operands or on where the input rings started: dms_union(A, B) and dms_union(B, A) are the same sequence, and dms_xor(A, A) is empty.
  • An element whose result is empty is an empty sequence, not a zero-area remnant; an element with an undefined operand is undefined.

which family to use

  valid input required coordinates fault tolerance
dms_ no all six point types the even-odd rule, plus a grid size as tolerance
geos_ yes, repaired with MakeValid when needed dpoint a repair pass before the operation
bg_ yes all six point types a repair pass through GEOS
bp_ yes integer only, below 2^25 none
cgal_ yes all six point types exact arithmetic, slow

For data that is already valid, the geos_ operators remain the first choice; they are faster and carry the full OGC operation set. Reach for dms_ when the source is not clean, or when a tolerance is what the model needs.

diagnostics

At minor trace level each call reports how much snapping it did, for instance:

dms_union: 1240 sequences; 87 crossings snapped to the grid, 0 extra noding rounds, 0 undefined results

The crossing count is how many intersections were moved onto the grid, the noding rounds are the passes in which that snapping created a new incidence that had to be resolved, and the last number counts the elements that came out undefined.

what is not in this family

The infix operators *, +, - and ^ keep their boost polygon binding for integer coordinates (bp_intersect and its family) and their geos binding for float coordinates (geos_intersect and its family). To get the fault-tolerant reading, call the operator by name.

There are no dms_ buffer or simplify operators; use the geos ones, which need valid input, or clean first with dms_polygon.

examples

  • Cleaning invalid polygons - cleaning a source, reading a self-intersecting ring, closing slivers with a grid size, cutting a hole
  • The configurations that pin these operators ship with GeoDMS in examples\testcases: oper_dms_overlay.dms (the areas of the four operations, the algebraic identities, agreement with the bp_ and geos_ results), oper_dms_overlay_types.dms (all six point types), oper_dms_overlay_topology.dms (the multi polygon layout: several shells, a hole, an island in a lake, rings that touch) oper_dms_overlay_neg1.dms to _neg5.dms (the four grid refusals and the values unit check), oper_dms_union_counts.dms (the dissolve: overlap is union, a ring wound the wrong way adds, partitions, the grid from a ranged and from an unranged dpoint unit), oper_dms_union_steep.dms (a snapped steep edge, against the pairwise result) and oper_dms_crossings.dms (the crossing sweep on its degenerate cases: several edges through one crossing, crossings on lattice points, vertical edges, collinear edges, a self-crossing ring).

see also