diff --git a/CHANGELOG.md b/CHANGELOG.md index 09bd7ef7..040fd7d7 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -2,6 +2,11 @@ ## Unreleased +- Add `geos.ClipByRect` function that clips a geometry to an axis-aligned + rectangle (defined by a `geom.Envelope`). This wraps the GEOS + `GEOSClipByRect` operation, which is faster than computing a full + `Intersection` with a rectangular polygon. + ## v0.58.0 2026-02-15 diff --git a/geos/entrypoints.go b/geos/entrypoints.go index b635e7ce..5b10bcd1 100644 --- a/geos/entrypoints.go +++ b/geos/entrypoints.go @@ -331,3 +331,16 @@ func UnaryUnion(g geom.Geometry) (geom.Geometry, error) { func ConcaveHull(g geom.Geometry, concavenessRatio float64, allowHoles bool) (geom.Geometry, error) { return rawgeos.ConcaveHull(g, concavenessRatio, allowHoles) } + +// ClipByRect clips a geometry to an axis-aligned rectangle defined by the +// given [geom.Envelope]. If the envelope is empty, then an empty +// [geom.GeometryCollection] is returned. +// +// The validity of the result is not checked. +func ClipByRect(g geom.Geometry, rect geom.Envelope) (geom.Geometry, error) { + lo, hi, ok := rect.MinMaxXYs() + if !ok { + return geom.Geometry{}, nil + } + return rawgeos.ClipByRect(g, lo.X, lo.Y, hi.X, hi.Y) +} diff --git a/geos/entrypoints_test.go b/geos/entrypoints_test.go index bbea0a2c..aaebea45 100644 --- a/geos/entrypoints_test.go +++ b/geos/entrypoints_test.go @@ -1016,3 +1016,86 @@ func TestCoverageIsValid(t *testing.T) { }) } } + +func TestClipByRect(t *testing.T) { + for _, tc := range []struct { + name string + input string + rect geom.Envelope + want string + }{ + { + name: "polygon fully inside rect", + input: "POLYGON((1 1,1 2,2 2,2 1,1 1))", + rect: geom.NewEnvelope(geom.XY{X: 0, Y: 0}, geom.XY{X: 3, Y: 3}), + want: "POLYGON((1 1,1 2,2 2,2 1,1 1))", + }, + { + name: "polygon partially overlapping rect", + input: "POLYGON((0 0,0 4,4 4,4 0,0 0))", + rect: geom.NewEnvelope(geom.XY{X: 1, Y: 1}, geom.XY{X: 3, Y: 3}), + want: "POLYGON((1 1,1 3,3 3,3 1,1 1))", + }, + { + name: "polygon fully outside rect", + input: "POLYGON((0 0,0 1,1 1,1 0,0 0))", + rect: geom.NewEnvelope(geom.XY{X: 5, Y: 5}, geom.XY{X: 6, Y: 6}), + want: "GEOMETRYCOLLECTION EMPTY", + }, + { + name: "linestring clipped by rect", + input: "LINESTRING(0 0,4 4)", + rect: geom.NewEnvelope(geom.XY{X: 1, Y: 1}, geom.XY{X: 3, Y: 3}), + want: "LINESTRING(1 1,3 3)", + }, + { + name: "point inside rect", + input: "POINT(2 2)", + rect: geom.NewEnvelope(geom.XY{X: 1, Y: 1}, geom.XY{X: 3, Y: 3}), + want: "POINT(2 2)", + }, + { + name: "point outside rect", + input: "POINT(0 0)", + rect: geom.NewEnvelope(geom.XY{X: 1, Y: 1}, geom.XY{X: 3, Y: 3}), + want: "GEOMETRYCOLLECTION EMPTY", + }, + { + name: "empty input geometry", + input: "GEOMETRYCOLLECTION EMPTY", + rect: geom.NewEnvelope(geom.XY{X: 0, Y: 0}, geom.XY{X: 1, Y: 1}), + want: "GEOMETRYCOLLECTION EMPTY", + }, + { + name: "u-shaped polygon clipped through both arms produces multipolygon", + input: "POLYGON((0 0,4 0,4 3,3 3,3 1,1 1,1 3,0 3,0 0))", + rect: geom.NewEnvelope(geom.XY{X: 0, Y: 2}, geom.XY{X: 4, Y: 4}), + want: "MULTIPOLYGON(((0 2,0 3,1 3,1 2,0 2)),((3 2,3 3,4 3,4 2,3 2)))", + }, + { + name: "polygon with hole inside rect", + input: "POLYGON((0 0,0 6,6 6,6 0,0 0),(2 2,4 2,4 4,2 4,2 2))", + rect: geom.NewEnvelope(geom.XY{X: 1, Y: 1}, geom.XY{X: 5, Y: 5}), + want: "POLYGON((1 1,1 5,5 5,5 1,1 1),(2 2,4 2,4 4,2 4,2 2))", + }, + { + name: "polygon with hole partially outside rect removes hole", + input: "POLYGON((0 0,0 6,6 6,6 0,0 0),(1 1,3 1,3 3,1 3,1 1))", + rect: geom.NewEnvelope(geom.XY{X: 2, Y: 2}, geom.XY{X: 5, Y: 5}), + want: "POLYGON((2 3,2 5,5 5,5 2,3 2,3 3,2 3))", + }, + { + name: "empty envelope", + input: "POLYGON((0 0,0 1,1 1,1 0,0 0))", + rect: geom.Envelope{}, + want: "GEOMETRYCOLLECTION EMPTY", + }, + } { + t.Run(tc.name, func(t *testing.T) { + got, err := geos.ClipByRect(geomFromWKT(t, tc.input), tc.rect) + skipIfUnsupported(t, err) + expectNoErr(t, err) + expectGeomEq(t, got, geomFromWKT(t, tc.want), geom.IgnoreOrder) + }) + } +} diff --git a/internal/rawgeos/entrypoints.go b/internal/rawgeos/entrypoints.go index e512a69b..94ed640a 100644 --- a/internal/rawgeos/entrypoints.go +++ b/internal/rawgeos/entrypoints.go @@ -44,6 +44,16 @@ GEOSGeometry *GEOSCoverageSimplifyVW_r(GEOSContextHandle_t handle, const GEOSGeo int GEOSCoverageIsValid_r(GEOSContextHandle_t handle, const GEOSGeometry* g, double gapWidth, GEOSGeometry** invalidEdges) { return 2; } #endif +#define CLIP_BY_RECT_MIN_VERSION "3.5.0" +#define CLIP_BY_RECT_MISSING ( \ + GEOS_VERSION_MAJOR < 3 || \ + (GEOS_VERSION_MAJOR == 3 && GEOS_VERSION_MINOR < 5) \ +) +#if CLIP_BY_RECT_MISSING +// This stub implementation always fails: +GEOSGeometry *GEOSClipByRect_r(GEOSContextHandle_t handle, const GEOSGeometry* g, double xmin, double ymin, double xmax, double ymax) { return NULL; } +#endif + #define CONCAVE_HULL_MIN_VERSION "3.11.0" #define CONCAVE_HULL_MISSING ( \ GEOS_VERSION_MAJOR < 3 || \ @@ -441,6 +451,18 @@ func Envelope(g geom.Geometry) (geom.Geometry, error) { return result, wrap(err, "executing GEOSEnvelope_r") } +func ClipByRect(g geom.Geometry, xmin, ymin, xmax, ymax float64) (geom.Geometry, error) { + if C.CLIP_BY_RECT_MISSING != 0 { + return geom.Geometry{}, UnsupportedGEOSVersionError{ + C.CLIP_BY_RECT_MIN_VERSION, "ClipByRect", + } + } + result, err := unaryOpG(g, func(ctx C.GEOSContextHandle_t, g *C.GEOSGeometry) *C.GEOSGeometry { + return C.GEOSClipByRect_r(ctx, g, C.double(xmin), C.double(ymin), C.double(xmax), C.double(ymax)) + }) + return result, wrap(err, "executing GEOSClipByRect_r") +} + func Area(g geom.Geometry) (float64, error) { result, err := unaryOpF(g, func(h C.GEOSContextHandle_t, g *C.GEOSGeometry, d *C.double) C.int { return C.GEOSArea_r(h, g, d)