From 971ad408ceb8de1c3a5bc55278032fb2301079d5 Mon Sep 17 00:00:00 2001 From: Oliver Tan Date: Mon, 27 Jul 2020 18:03:40 -0700 Subject: [PATCH] builtins: implement ST_Intersection and ST_Buffer for geography types PostGIS uses _ST_BestSRID to approximate a geography as a geometry to perform ST_Buffer and ST_Intersection. We're obliged to do the same - comments are inline. Release note (sql change): Implements the ST_Buffer and ST_Intersection functions for Geography types. Release note (sql change): Implement ST_Intersection for string types. --- docs/generated/sql/functions.md | 44 +++ pkg/geo/geogfn/best_projection.go | 139 +++++++++ pkg/geo/geogfn/best_projection_test.go | 94 ++++++ .../logictest/testdata/logic_test/geospatial | 284 ++++++++++-------- pkg/sql/sem/builtins/geo_builtins.go | 205 ++++++++++++- 5 files changed, 641 insertions(+), 125 deletions(-) create mode 100644 pkg/geo/geogfn/best_projection.go create mode 100644 pkg/geo/geogfn/best_projection_test.go diff --git a/docs/generated/sql/functions.md b/docs/generated/sql/functions.md index 822f00221a52..9d75b591252f 100644 --- a/docs/generated/sql/functions.md +++ b/docs/generated/sql/functions.md @@ -930,6 +930,39 @@ has no relationship with the commit order of concurrent transactions.

st_azimuth(geometry_a: geometry, geometry_b: geometry) → float

Returns the azimuth in radians of the segment defined by the given point geometries, or NULL if the two points are coincident.

The azimuth is angle is referenced from north, and is positive clockwise: North = 0; East = π/2; South = π; West = 3π/2.

 +st_buffer(geography: geography, distance: float) → geography

Returns a Geometry that represents all points whose distance is less than or equal to the given distance +from the given Geometry.

+

This function utilizes the GEOS module.

+

This operation is done by transforming the object into a Geometry. This occurs by translating +the Geography objects into Geometry objects before applying an LAEA, UTM or Web Mercator +based projection based on the bounding boxes of the given Geography objects. When the result is +calculated, the result is transformed back into a Geography with SRID 4326.

+ +st_buffer(geography: geography, distance: float, buffer_style_params: string) → geography

Returns a Geometry that represents all points whose distance is less than or equal to the given distance from the +given Geometry.

+

This variant takes in a space separate parameter string, which will augment the buffer styles. Valid parameters are:

+ +

This function utilizes the GEOS module.

+

This operation is done by transforming the object into a Geometry. This occurs by translating +the Geography objects into Geometry objects before applying an LAEA, UTM or Web Mercator +based projection based on the bounding boxes of the given Geography objects. When the result is +calculated, the result is transformed back into a Geography with SRID 4326.

+ +st_buffer(geography: geography, distance: float, quad_segs: int) → geography

Returns a Geometry that represents all points whose distance is less than or equal to the given distance from the +given Geometry.

+

This variant approximates the circle into quad_seg segments per line (the default is 8).

+

This function utilizes the GEOS module.

+

This operation is done by transforming the object into a Geometry. This occurs by translating +the Geography objects into Geometry objects before applying an LAEA, UTM or Web Mercator +based projection based on the bounding boxes of the given Geography objects. When the result is +calculated, the result is transformed back into a Geography with SRID 4326.

+ st_buffer(geometry: geometry, distance: decimal) → geometry

Returns a Geometry that represents all points whose distance is less than or equal to the given distance from the given Geometry.

This function utilizes the GEOS module.

@@ -1171,9 +1204,20 @@ Bottom Left.

 st_interiorringn(geometry: geometry, n: int) → geometry

Returns the n-th (1-indexed) interior ring of a Polygon as a LineString. Returns NULL if the shape is not a Polygon, or the ring does not exist.

 +st_intersection(geography_a: geography, geography_b: geography) → geography

Returns the point intersections of the given geographies.

+

This operation is done by transforming the object into a Geometry. This occurs by translating +the Geography objects into Geometry objects before applying an LAEA, UTM or Web Mercator +based projection based on the bounding boxes of the given Geography objects. When the result is +calculated, the result is transformed back into a Geography with SRID 4326.

+

This function utilizes the GEOS module.

+ st_intersection(geometry_a: geometry, geometry_b: geometry) → geometry

Returns the point intersections of the given geometries.

This function utilizes the GEOS module.

 +st_intersection(geometry_a_str: string, geometry_b_str: string) → geometry

Returns the point intersections of the given geometries.

+

This function utilizes the GEOS module.

+

This variant will cast all geometry_str arguments into Geometry types.

+ st_intersects(geography_a: geography, geography_b: geography) → bool

Returns true if geography_a shares any portion of space with geography_b.

The calculations performed are have a precision of 1cm.

This function utilizes the S2 library for spherical calculations.

diff --git a/pkg/geo/geogfn/best_projection.go b/pkg/geo/geogfn/best_projection.go new file mode 100644 index 000000000000..91e6c0718220 --- /dev/null +++ b/pkg/geo/geogfn/best_projection.go @@ -0,0 +1,139 @@ +// Copyright 2020 The Cockroach Authors. +// +// Use of this software is governed by the Business Source License +// included in the file licenses/BSL.txt. +// +// As of the Change Date specified in that file, in accordance with +// the Business Source License, use of this software will be governed +// by the Apache License, Version 2.0, included in the file +// licenses/APL.txt. + +package geogfn + +import ( + "fmt" + "math" + + "github.com/cockroachdb/cockroach/pkg/geo/geopb" + "github.com/cockroachdb/cockroach/pkg/geo/geoprojbase" + "github.com/cockroachdb/errors" + "github.com/golang/geo/s1" + "github.com/golang/geo/s2" +) + +// BestGeomProjection translates roughly to the ST_BestSRID function in PostGIS. +// It attempts to find the best projection for a bounding box into an accurate +// geometry-type projection. +// +// The algorithm is described by ST_Buffer/ST_Intersection documentation (paraphrased): +// It first determines the best SRID that fits the bounding box of the 2 geography objects (ST_Intersection only). +// It favors a north/south pole projection, then UTM, then LAEA for smaller zones, otherwise falling back +// to web mercator. +// If geography objects are within one half zone UTM but not the same UTM it will pick one of those. +// After the calculation is complete, it will fall back to WGS84 Geography. +func BestGeomProjection(boundingRect s2.Rect) (geoprojbase.Proj4Text, error) { + center := boundingRect.Center() + + latWidth := s1.Angle(boundingRect.Lat.Length()) + lngWidth := s1.Angle(boundingRect.Lng.Length()) + + // Check if these fit either the North Pole or South Pole areas. + // If the center has latitude greater than 70 (an arbitrary polar threshold), and it is + // within the polar ranges, return that. + if center.Lat.Degrees() > 70 && boundingRect.Lo().Lat.Degrees() > 45 { + // See: https://epsg.io/3574. + return getGeomProjection(3574) + } + // Same for south pole. + if center.Lat.Degrees() < -70 && boundingRect.Hi().Lat.Degrees() < -45 { + // See: https://epsg.io/3409 + return getGeomProjection(3409) + } + + // Each UTM zone is 6 degrees wide and distortion is low for geometries that fit within the zone. We use + // UTM if the width is lower than the UTM zone width, even though the geometry may span 2 zones -- using + // the geometry center to pick the UTM zone should result in most of the geometry being in the picked zone. + if lngWidth.Degrees() < 6 { + // Determine the offset of the projection. + // Offset longitude -180 to 0 and divide by 6 to get the zone. + // Note that we treat 180 degree longtitudes as offset 59. + // TODO(#geo): do we care about https://en.wikipedia.org/wiki/Universal_Transverse_Mercator_coordinate_system#Exceptions? + // PostGIS's _ST_BestSRID function doesn't seem to care: . + sridOffset := geopb.SRID(math.Min(math.Floor((center.Lng.Degrees()+180)/6), 59)) + if center.Lat.Degrees() >= 0 { + // Start at the north UTM SRID. + return getGeomProjection(32601 + sridOffset) + } + // Start at the south UTM SRID. + // This should make no difference in end result compared to using the north UTMs, + // but for completeness we do it. + return getGeomProjection(32701 + sridOffset) + } + + // Attempt to fit into LAEA areas if the width is less than 25 degrees (we can go up to 30 + // but want to leave some room for precision issues). + // + // LAEA areas are separated into 3 latitude zones between 0 and 90 and 3 latitude zones + // between -90 and 0. Within each latitude zones, they have different longitude bands: + // * The bands closest to the equator have 12x30 degree longitude zones. + // * The bands in the temperate area 8x45 degree longitude zones. + // * The bands near the poles have 4x90 degree longitude zones. + // + // For each of these bands, we custom define a LAEA area with the center of the LAEA area + // as the lat/lon offset. + // + // See also: https://en.wikipedia.org/wiki/Lambert_azimuthal_equal-area_projection. + if latWidth.Degrees() < 25 { + // Convert lat to a known 30 degree zone.. + // -3 represents [-90, -60), -2 represents [-60, -30) ... and 2 represents [60, 90]. + // (note: 90 is inclusive at the end). + // Treat a 90 degree latitude as band 2. + latZone := math.Min(math.Floor(center.Lat.Degrees()/30), 2) + latZoneCenterDegrees := (latZone * 30) + 15 + // Equator bands - 30 degree zones. + if (latZone == 0 || latZone == -1) && lngWidth.Degrees() <= 30 { + lngZone := math.Floor(center.Lng.Degrees() / 30) + return geoprojbase.MakeProj4Text( + fmt.Sprintf( + "+proj=laea +ellps=WGS84 +datum=WGS84 +lat_0=%g +lon_0=%g +units=m +no_defs", + latZoneCenterDegrees, + (lngZone*30)+15, + ), + ), nil + } + // Temperate bands - 45 degree zones. + if (latZone == -2 || latZone == 1) && lngWidth.Degrees() <= 45 { + lngZone := math.Floor(center.Lng.Degrees() / 45) + return geoprojbase.MakeProj4Text( + fmt.Sprintf( + "+proj=laea +ellps=WGS84 +datum=WGS84 +lat_0=%g +lon_0=%g +units=m +no_defs", + latZoneCenterDegrees, + (lngZone*45)+22.5, + ), + ), nil + } + // Polar bands -- 90 degree zones. + if (latZone == -3 || latZone == 2) && lngWidth.Degrees() <= 90 { + lngZone := math.Floor(center.Lng.Degrees() / 90) + return geoprojbase.MakeProj4Text( + fmt.Sprintf( + "+proj=laea +ellps=WGS84 +datum=WGS84 +lat_0=%g +lon_0=%g +units=m +no_defs", + latZoneCenterDegrees, + (lngZone*90)+45, + ), + ), nil + } + } + + // Default to Web Mercator. + return getGeomProjection(3857) +} + +// getGeomProjection returns the Proj4Text associated with an SRID. +func getGeomProjection(srid geopb.SRID) (geoprojbase.Proj4Text, error) { + proj, ok := geoprojbase.Projection(srid) + if !ok { + return geoprojbase.Proj4Text{}, errors.Newf("unexpected SRID %d", srid) + } + return proj.Proj4Text, nil +} diff --git a/pkg/geo/geogfn/best_projection_test.go b/pkg/geo/geogfn/best_projection_test.go new file mode 100644 index 000000000000..f339690dbe06 --- /dev/null +++ b/pkg/geo/geogfn/best_projection_test.go @@ -0,0 +1,94 @@ +// Copyright 2020 The Cockroach Authors. +// +// Use of this software is governed by the Business Source License +// included in the file licenses/BSL.txt. +// +// As of the Change Date specified in that file, in accordance with +// the Business Source License, use of this software will be governed +// by the Apache License, Version 2.0, included in the file +// licenses/APL.txt. + +package geogfn + +import ( + "testing" + + "github.com/cockroachdb/cockroach/pkg/geo/geoprojbase" + "github.com/golang/geo/s2" + "github.com/stretchr/testify/require" +) + +func TestBestGeomProjection(t *testing.T) { + testCases := []struct { + desc string + rect s2.Rect + expected geoprojbase.Proj4Text + }{ + { + "north pole", + s2.RectFromLatLng(s2.LatLngFromDegrees(75, 75)), + geoprojbase.Projections[3574].Proj4Text, + }, + { + "south pole", + s2.RectFromLatLng(s2.LatLngFromDegrees(-75, -75)), + geoprojbase.Projections[3409].Proj4Text}, + { + "utm 15 on top hemisphere", + s2.RectFromLatLng(s2.LatLngFromDegrees(15, 93)), + geoprojbase.Projections[32646].Proj4Text}, + { + "utm -16 on bottom hemisphere", + s2.RectFromLatLng(s2.LatLngFromDegrees(-15, -111)), + geoprojbase.Projections[32712].Proj4Text, + }, + { + "LAEA at equator bands (north half)", + s2.RectFromCenterSize(s2.LatLngFromDegrees(12, 22), s2.LatLngFromDegrees(10, 10)), + geoprojbase.MakeProj4Text("+proj=laea +ellps=WGS84 +datum=WGS84 +lat_0=15 +lon_0=15 +units=m +no_defs"), + }, + { + "LAEA at equator bands (south half)", + s2.RectFromCenterSize(s2.LatLngFromDegrees(-13, 123), s2.LatLngFromDegrees(10, 10)), + geoprojbase.MakeProj4Text("+proj=laea +ellps=WGS84 +datum=WGS84 +lat_0=-15 +lon_0=135 +units=m +no_defs"), + }, + { + "LAEA at temperate band (north half)", + s2.RectFromCenterSize(s2.LatLngFromDegrees(33, 87), s2.LatLngFromDegrees(10, 10)), + geoprojbase.MakeProj4Text("+proj=laea +ellps=WGS84 +datum=WGS84 +lat_0=45 +lon_0=67.5 +units=m +no_defs"), + }, + { + "LAEA at temperate band (south half)", + s2.RectFromCenterSize(s2.LatLngFromDegrees(-53, -120.5), s2.LatLngFromDegrees(10, 10)), + geoprojbase.MakeProj4Text("+proj=laea +ellps=WGS84 +datum=WGS84 +lat_0=-45 +lon_0=-112.5 +units=m +no_defs"), + }, + { + "LAEA at polar band (north half)", + s2.RectFromCenterSize(s2.LatLngFromDegrees(63, 87), s2.LatLngFromDegrees(10, 10)), + geoprojbase.MakeProj4Text("+proj=laea +ellps=WGS84 +datum=WGS84 +lat_0=75 +lon_0=45 +units=m +no_defs"), + }, + { + "LAEA at polar band (south half)", + s2.RectFromCenterSize(s2.LatLngFromDegrees(-66, -120.5), s2.LatLngFromDegrees(10, 10)), + geoprojbase.MakeProj4Text("+proj=laea +ellps=WGS84 +datum=WGS84 +lat_0=-75 +lon_0=-135 +units=m +no_defs"), + }, + { + "UTM which should be 32V, but we return 31V as we do not handle UTM exceptions", + s2.RectFromLatLng(s2.LatLngFromDegrees(59.4136, 5.26)), + geoprojbase.Projections[32631].Proj4Text, // Should be 32632 + }, + { + "wide example", + s2.RectFromCenterSize(s2.LatLngFromDegrees(0, 0), s2.LatLngFromDegrees(50, 50)), + geoprojbase.Projections[3857].Proj4Text, + }, + } + + for _, tc := range testCases { + t.Run(tc.desc, func(t *testing.T) { + result, err := BestGeomProjection(tc.rect) + require.NoError(t, err) + require.Equal(t, tc.expected.String(), result.String()) + }) + } +} diff --git a/pkg/sql/logictest/testdata/logic_test/geospatial b/pkg/sql/logictest/testdata/logic_test/geospatial index 919c961c3a4b..ec9a6b54cede 100644 --- a/pkg/sql/logictest/testdata/logic_test/geospatial +++ b/pkg/sql/logictest/testdata/logic_test/geospatial @@ -782,136 +782,137 @@ Square overlapping left and right square Square (left) Square overlapping left and right square Square (right) POLYGON ((0 0, -0.1 0, -0.1 1, 0 1, 1 1, 1 0, 0 0)) Square overlapping left and right square Square overlapping left and right square POLYGON ((1 0, -0.1 0, -0.1 1, 1 1, 1 0)) -query TTT +query TTTB SELECT a.dsc, b.dsc, - ST_AsEWKT(ST_Intersection(a.geom, b.geom)) + ST_AsEWKT(ST_Intersection(a.geom, b.geom)), + ST_Intersection(a.geom, b.geom) = ST_Intersection(a.geom::string, b.geom::string) FROM geom_operators_test a JOIN geom_operators_test b ON (1=1) ORDER BY a.dsc, b.dsc ---- -Empty GeometryCollection Empty GeometryCollection GEOMETRYCOLLECTION EMPTY -Empty GeometryCollection Empty LineString GEOMETRYCOLLECTION EMPTY -Empty GeometryCollection Empty Point GEOMETRYCOLLECTION EMPTY -Empty GeometryCollection Faraway point GEOMETRYCOLLECTION EMPTY -Empty GeometryCollection Line going through left and right square GEOMETRYCOLLECTION EMPTY -Empty GeometryCollection NULL NULL -Empty GeometryCollection Point middle of Left Square GEOMETRYCOLLECTION EMPTY -Empty GeometryCollection Point middle of Right Square GEOMETRYCOLLECTION EMPTY -Empty GeometryCollection Square (left) GEOMETRYCOLLECTION EMPTY -Empty GeometryCollection Square (right) GEOMETRYCOLLECTION EMPTY -Empty GeometryCollection Square overlapping left and right square GEOMETRYCOLLECTION EMPTY -Empty LineString Empty GeometryCollection GEOMETRYCOLLECTION EMPTY -Empty LineString Empty LineString GEOMETRYCOLLECTION EMPTY -Empty LineString Empty Point GEOMETRYCOLLECTION EMPTY -Empty LineString Faraway point GEOMETRYCOLLECTION EMPTY -Empty LineString Line going through left and right square GEOMETRYCOLLECTION EMPTY -Empty LineString NULL NULL -Empty LineString Point middle of Left Square GEOMETRYCOLLECTION EMPTY -Empty LineString Point middle of Right Square GEOMETRYCOLLECTION EMPTY -Empty LineString Square (left) GEOMETRYCOLLECTION EMPTY -Empty LineString Square (right) GEOMETRYCOLLECTION EMPTY -Empty LineString Square overlapping left and right square GEOMETRYCOLLECTION EMPTY -Empty Point Empty GeometryCollection GEOMETRYCOLLECTION EMPTY -Empty Point Empty LineString GEOMETRYCOLLECTION EMPTY -Empty Point Empty Point GEOMETRYCOLLECTION EMPTY -Empty Point Faraway point GEOMETRYCOLLECTION EMPTY -Empty Point Line going through left and right square GEOMETRYCOLLECTION EMPTY -Empty Point NULL NULL -Empty Point Point middle of Left Square GEOMETRYCOLLECTION EMPTY -Empty Point Point middle of Right Square GEOMETRYCOLLECTION EMPTY -Empty Point Square (left) GEOMETRYCOLLECTION EMPTY -Empty Point Square (right) GEOMETRYCOLLECTION EMPTY -Empty Point Square overlapping left and right square GEOMETRYCOLLECTION EMPTY -Faraway point Empty GeometryCollection GEOMETRYCOLLECTION EMPTY -Faraway point Empty LineString GEOMETRYCOLLECTION EMPTY -Faraway point Empty Point GEOMETRYCOLLECTION EMPTY -Faraway point Faraway point POINT (5 5) -Faraway point Line going through left and right square POINT EMPTY -Faraway point NULL NULL -Faraway point Point middle of Left Square POINT EMPTY -Faraway point Point middle of Right Square POINT EMPTY -Faraway point Square (left) POINT EMPTY -Faraway point Square (right) POINT EMPTY -Faraway point Square overlapping left and right square POINT EMPTY -Line going through left and right square Empty GeometryCollection GEOMETRYCOLLECTION EMPTY -Line going through left and right square Empty LineString GEOMETRYCOLLECTION EMPTY -Line going through left and right square Empty Point GEOMETRYCOLLECTION EMPTY -Line going through left and right square Faraway point POINT EMPTY -Line going through left and right square Line going through left and right square LINESTRING (-0.5 0.5, 0.5 0.5) -Line going through left and right square NULL NULL -Line going through left and right square Point middle of Left Square POINT (-0.5 0.5) -Line going through left and right square Point middle of Right Square POINT (0.5 0.5) -Line going through left and right square Square (left) LINESTRING (-0.5 0.5, 0 0.5) -Line going through left and right square Square (right) LINESTRING (0 0.5, 0.5 0.5) -Line going through left and right square Square overlapping left and right square LINESTRING (-0.1 0.5, 0.5 0.5) -NULL Empty GeometryCollection NULL -NULL Empty LineString NULL -NULL Empty Point NULL -NULL Faraway point NULL -NULL Line going through left and right square NULL -NULL NULL NULL -NULL Point middle of Left Square NULL -NULL Point middle of Right Square NULL -NULL Square (left) NULL -NULL Square (right) NULL -NULL Square overlapping left and right square NULL -Point middle of Left Square Empty GeometryCollection GEOMETRYCOLLECTION EMPTY -Point middle of Left Square Empty LineString GEOMETRYCOLLECTION EMPTY -Point middle of Left Square Empty Point GEOMETRYCOLLECTION EMPTY -Point middle of Left Square Faraway point POINT EMPTY -Point middle of Left Square Line going through left and right square POINT (-0.5 0.5) -Point middle of Left Square NULL NULL -Point middle of Left Square Point middle of Left Square POINT (-0.5 0.5) -Point middle of Left Square Point middle of Right Square POINT EMPTY -Point middle of Left Square Square (left) POINT (-0.5 0.5) -Point middle of Left Square Square (right) POINT EMPTY -Point middle of Left Square Square overlapping left and right square POINT EMPTY -Point middle of Right Square Empty GeometryCollection GEOMETRYCOLLECTION EMPTY -Point middle of Right Square Empty LineString GEOMETRYCOLLECTION EMPTY -Point middle of Right Square Empty Point GEOMETRYCOLLECTION EMPTY -Point middle of Right Square Faraway point POINT EMPTY -Point middle of Right Square Line going through left and right square POINT (0.5 0.5) -Point middle of Right Square NULL NULL -Point middle of Right Square Point middle of Left Square POINT EMPTY -Point middle of Right Square Point middle of Right Square POINT (0.5 0.5) -Point middle of Right Square Square (left) POINT EMPTY -Point middle of Right Square Square (right) POINT (0.5 0.5) -Point middle of Right Square Square overlapping left and right square POINT (0.5 0.5) -Square (left) Empty GeometryCollection GEOMETRYCOLLECTION EMPTY -Square (left) Empty LineString GEOMETRYCOLLECTION EMPTY -Square (left) Empty Point GEOMETRYCOLLECTION EMPTY -Square (left) Faraway point POINT EMPTY -Square (left) Line going through left and right square LINESTRING (-0.5 0.5, 0 0.5) -Square (left) NULL NULL -Square (left) Point middle of Left Square POINT (-0.5 0.5) -Square (left) Point middle of Right Square POINT EMPTY -Square (left) Square (left) POLYGON ((0 0, -1 0, -1 1, 0 1, 0 0)) -Square (left) Square (right) LINESTRING (0 0, 0 1) -Square (left) Square overlapping left and right square POLYGON ((0 0, -0.1 0, -0.1 1, 0 1, 0 0)) -Square (right) Empty GeometryCollection GEOMETRYCOLLECTION EMPTY -Square (right) Empty LineString GEOMETRYCOLLECTION EMPTY -Square (right) Empty Point GEOMETRYCOLLECTION EMPTY -Square (right) Faraway point POINT EMPTY -Square (right) Line going through left and right square LINESTRING (0 0.5, 0.5 0.5) -Square (right) NULL NULL -Square (right) Point middle of Left Square POINT EMPTY -Square (right) Point middle of Right Square POINT (0.5 0.5) -Square (right) Square (left) LINESTRING (0 1, 0 0) -Square (right) Square (right) POLYGON ((1 0, 0 0, 0 1, 1 1, 1 0)) -Square (right) Square overlapping left and right square POLYGON ((1 0, 0 0, 0 1, 1 1, 1 0)) -Square overlapping left and right square Empty GeometryCollection GEOMETRYCOLLECTION EMPTY -Square overlapping left and right square Empty LineString GEOMETRYCOLLECTION EMPTY -Square overlapping left and right square Empty Point GEOMETRYCOLLECTION EMPTY -Square overlapping left and right square Faraway point POINT EMPTY -Square overlapping left and right square Line going through left and right square LINESTRING (-0.1 0.5, 0.5 0.5) -Square overlapping left and right square NULL NULL -Square overlapping left and right square Point middle of Left Square POINT EMPTY -Square overlapping left and right square Point middle of Right Square POINT (0.5 0.5) -Square overlapping left and right square Square (left) POLYGON ((0 0, -0.1 0, -0.1 1, 0 1, 0 0)) -Square overlapping left and right square Square (right) POLYGON ((1 0, 0 0, 0 1, 1 1, 1 0)) -Square overlapping left and right square Square overlapping left and right square POLYGON ((1 0, -0.1 0, -0.1 1, 1 1, 1 0)) +Empty GeometryCollection Empty GeometryCollection GEOMETRYCOLLECTION EMPTY true +Empty GeometryCollection Empty LineString GEOMETRYCOLLECTION EMPTY true +Empty GeometryCollection Empty Point GEOMETRYCOLLECTION EMPTY true +Empty GeometryCollection Faraway point GEOMETRYCOLLECTION EMPTY true +Empty GeometryCollection Line going through left and right square GEOMETRYCOLLECTION EMPTY true +Empty GeometryCollection NULL NULL NULL +Empty GeometryCollection Point middle of Left Square GEOMETRYCOLLECTION EMPTY true +Empty GeometryCollection Point middle of Right Square GEOMETRYCOLLECTION EMPTY true +Empty GeometryCollection Square (left) GEOMETRYCOLLECTION EMPTY true +Empty GeometryCollection Square (right) GEOMETRYCOLLECTION EMPTY true +Empty GeometryCollection Square overlapping left and right square GEOMETRYCOLLECTION EMPTY true +Empty LineString Empty GeometryCollection GEOMETRYCOLLECTION EMPTY true +Empty LineString Empty LineString GEOMETRYCOLLECTION EMPTY true +Empty LineString Empty Point GEOMETRYCOLLECTION EMPTY true +Empty LineString Faraway point GEOMETRYCOLLECTION EMPTY true +Empty LineString Line going through left and right square GEOMETRYCOLLECTION EMPTY true +Empty LineString NULL NULL NULL +Empty LineString Point middle of Left Square GEOMETRYCOLLECTION EMPTY true +Empty LineString Point middle of Right Square GEOMETRYCOLLECTION EMPTY true +Empty LineString Square (left) GEOMETRYCOLLECTION EMPTY true +Empty LineString Square (right) GEOMETRYCOLLECTION EMPTY true +Empty LineString Square overlapping left and right square GEOMETRYCOLLECTION EMPTY true +Empty Point Empty GeometryCollection GEOMETRYCOLLECTION EMPTY true +Empty Point Empty LineString GEOMETRYCOLLECTION EMPTY true +Empty Point Empty Point GEOMETRYCOLLECTION EMPTY true +Empty Point Faraway point GEOMETRYCOLLECTION EMPTY true +Empty Point Line going through left and right square GEOMETRYCOLLECTION EMPTY true +Empty Point NULL NULL NULL +Empty Point Point middle of Left Square GEOMETRYCOLLECTION EMPTY true +Empty Point Point middle of Right Square GEOMETRYCOLLECTION EMPTY true +Empty Point Square (left) GEOMETRYCOLLECTION EMPTY true +Empty Point Square (right) GEOMETRYCOLLECTION EMPTY true +Empty Point Square overlapping left and right square GEOMETRYCOLLECTION EMPTY true +Faraway point Empty GeometryCollection GEOMETRYCOLLECTION EMPTY true +Faraway point Empty LineString GEOMETRYCOLLECTION EMPTY true +Faraway point Empty Point GEOMETRYCOLLECTION EMPTY true +Faraway point Faraway point POINT (5 5) true +Faraway point Line going through left and right square POINT EMPTY true +Faraway point NULL NULL NULL +Faraway point Point middle of Left Square POINT EMPTY true +Faraway point Point middle of Right Square POINT EMPTY true +Faraway point Square (left) POINT EMPTY true +Faraway point Square (right) POINT EMPTY true +Faraway point Square overlapping left and right square POINT EMPTY true +Line going through left and right square Empty GeometryCollection GEOMETRYCOLLECTION EMPTY true +Line going through left and right square Empty LineString GEOMETRYCOLLECTION EMPTY true +Line going through left and right square Empty Point GEOMETRYCOLLECTION EMPTY true +Line going through left and right square Faraway point POINT EMPTY true +Line going through left and right square Line going through left and right square LINESTRING (-0.5 0.5, 0.5 0.5) true +Line going through left and right square NULL NULL NULL +Line going through left and right square Point middle of Left Square POINT (-0.5 0.5) true +Line going through left and right square Point middle of Right Square POINT (0.5 0.5) true +Line going through left and right square Square (left) LINESTRING (-0.5 0.5, 0 0.5) true +Line going through left and right square Square (right) LINESTRING (0 0.5, 0.5 0.5) true +Line going through left and right square Square overlapping left and right square LINESTRING (-0.1 0.5, 0.5 0.5) true +NULL Empty GeometryCollection NULL NULL +NULL Empty LineString NULL NULL +NULL Empty Point NULL NULL +NULL Faraway point NULL NULL +NULL Line going through left and right square NULL NULL +NULL NULL NULL NULL +NULL Point middle of Left Square NULL NULL +NULL Point middle of Right Square NULL NULL +NULL Square (left) NULL NULL +NULL Square (right) NULL NULL +NULL Square overlapping left and right square NULL NULL +Point middle of Left Square Empty GeometryCollection GEOMETRYCOLLECTION EMPTY true +Point middle of Left Square Empty LineString GEOMETRYCOLLECTION EMPTY true +Point middle of Left Square Empty Point GEOMETRYCOLLECTION EMPTY true +Point middle of Left Square Faraway point POINT EMPTY true +Point middle of Left Square Line going through left and right square POINT (-0.5 0.5) true +Point middle of Left Square NULL NULL NULL +Point middle of Left Square Point middle of Left Square POINT (-0.5 0.5) true +Point middle of Left Square Point middle of Right Square POINT EMPTY true +Point middle of Left Square Square (left) POINT (-0.5 0.5) true +Point middle of Left Square Square (right) POINT EMPTY true +Point middle of Left Square Square overlapping left and right square POINT EMPTY true +Point middle of Right Square Empty GeometryCollection GEOMETRYCOLLECTION EMPTY true +Point middle of Right Square Empty LineString GEOMETRYCOLLECTION EMPTY true +Point middle of Right Square Empty Point GEOMETRYCOLLECTION EMPTY true +Point middle of Right Square Faraway point POINT EMPTY true +Point middle of Right Square Line going through left and right square POINT (0.5 0.5) true +Point middle of Right Square NULL NULL NULL +Point middle of Right Square Point middle of Left Square POINT EMPTY true +Point middle of Right Square Point middle of Right Square POINT (0.5 0.5) true +Point middle of Right Square Square (left) POINT EMPTY true +Point middle of Right Square Square (right) POINT (0.5 0.5) true +Point middle of Right Square Square overlapping left and right square POINT (0.5 0.5) true +Square (left) Empty GeometryCollection GEOMETRYCOLLECTION EMPTY true +Square (left) Empty LineString GEOMETRYCOLLECTION EMPTY true +Square (left) Empty Point GEOMETRYCOLLECTION EMPTY true +Square (left) Faraway point POINT EMPTY true +Square (left) Line going through left and right square LINESTRING (-0.5 0.5, 0 0.5) true +Square (left) NULL NULL NULL +Square (left) Point middle of Left Square POINT (-0.5 0.5) true +Square (left) Point middle of Right Square POINT EMPTY true +Square (left) Square (left) POLYGON ((0 0, -1 0, -1 1, 0 1, 0 0)) true +Square (left) Square (right) LINESTRING (0 0, 0 1) true +Square (left) Square overlapping left and right square POLYGON ((0 0, -0.1 0, -0.1 1, 0 1, 0 0)) true +Square (right) Empty GeometryCollection GEOMETRYCOLLECTION EMPTY true +Square (right) Empty LineString GEOMETRYCOLLECTION EMPTY true +Square (right) Empty Point GEOMETRYCOLLECTION EMPTY true +Square (right) Faraway point POINT EMPTY true +Square (right) Line going through left and right square LINESTRING (0 0.5, 0.5 0.5) true +Square (right) NULL NULL NULL +Square (right) Point middle of Left Square POINT EMPTY true +Square (right) Point middle of Right Square POINT (0.5 0.5) true +Square (right) Square (left) LINESTRING (0 1, 0 0) true +Square (right) Square (right) POLYGON ((1 0, 0 0, 0 1, 1 1, 1 0)) true +Square (right) Square overlapping left and right square POLYGON ((1 0, 0 0, 0 1, 1 1, 1 0)) true +Square overlapping left and right square Empty GeometryCollection GEOMETRYCOLLECTION EMPTY true +Square overlapping left and right square Empty LineString GEOMETRYCOLLECTION EMPTY true +Square overlapping left and right square Empty Point GEOMETRYCOLLECTION EMPTY true +Square overlapping left and right square Faraway point POINT EMPTY true +Square overlapping left and right square Line going through left and right square LINESTRING (-0.1 0.5, 0.5 0.5) true +Square overlapping left and right square NULL NULL NULL +Square overlapping left and right square Point middle of Left Square POINT EMPTY true +Square overlapping left and right square Point middle of Right Square POINT (0.5 0.5) true +Square overlapping left and right square Square (left) POLYGON ((0 0, -0.1 0, -0.1 1, 0 1, 0 0)) true +Square overlapping left and right square Square (right) POLYGON ((1 0, 0 0, 0 1, 1 1, 1 0)) true +Square overlapping left and right square Square overlapping left and right square POLYGON ((1 0, -0.1 0, -0.1 1, 1 1, 1 0)) true query TTTT SELECT @@ -3162,6 +3163,43 @@ Square (left) POLYGON ((-1 0, -0.5 0, 0 0, 0 0.5, 0 Square (right) POLYGON ((0 0, 0.5 0, 1 0, 1 0.5, 1 1, 0.5 1.000038070652873, 0 1, 0 0.5, 0 0)) POLYGON ((0 0, 0.25 0, 0.5 0, 0.75 0, 1 0, 1 0.25, 1 0.5, 1 0.75, 1 1, 0.750000001449807 1.000028552944327, 0.5 1.000038070652873, 0.249999998550193 1.000028552944327, 0 1, 0 0.75, 0 0.5, 0 0.25, 0 0)) Square overlapping left and right square POLYGON ((-0.1 0, 0.45 0, 1 0, 1 0.5, 1 1, 0.45 1.000046065796834, -0.1 1, -0.1 0.5, -0.1 0)) POLYGON ((-0.1 0, 0.175 0, 0.45 0, 0.725 0, 1 0, 1 0.25, 1 0.5, 1 0.75, 1 1, 0.725000001929716 1.000034549281259, 0.45 1.000046065796834, 0.174999998070284 1.000034549281259, -0.1 1, -0.1 0.75, -0.1 0.5, -0.1 0.25, -0.1 0)) +# ST_Buffer +query TTTT +SELECT + a.dsc, + ST_AsEWKT(ST_Buffer(a.geog, 10), 5), + ST_AsEWKT(ST_Buffer(a.geog, 10, 2), 5), + ST_AsEWKT(ST_Buffer(a.geog, 10, 'quad_segs=4 endcap=flat'), 5) +FROM geog_operators_test a +ORDER BY a.dsc +---- +Empty GeometryCollection SRID=4326;POLYGON EMPTY SRID=4326;POLYGON EMPTY SRID=4326;POLYGON EMPTY +Empty LineString SRID=4326;POLYGON EMPTY SRID=4326;POLYGON EMPTY SRID=4326;POLYGON EMPTY +Empty Point SRID=4326;POLYGON EMPTY SRID=4326;POLYGON EMPTY SRID=4326;POLYGON EMPTY +Faraway point SRID=4326;POLYGON ((5.00009 5, 5.00009 4.99998, 5.00008 4.99997, 5.00007 4.99995, 5.00006 4.99994, 5.00005 4.99992, 5.00003 4.99992, 5.00002 4.99991, 5 4.99991, 4.99998 4.99991, 4.99997 4.99992, 4.99995 4.99992, 4.99994 4.99994, 4.99992 4.99995, 4.99992 4.99997, 4.99991 4.99998, 4.99991 5, 4.99991 5.00002, 4.99992 5.00003, 4.99993 5.00005, 4.99994 5.00006, 4.99995 5.00008, 4.99997 5.00008, 4.99998 5.00009, 5 5.00009, 5.00002 5.00009, 5.00003 5.00008, 5.00005 5.00008, 5.00006 5.00006, 5.00008 5.00005, 5.00008 5.00003, 5.00009 5.00002, 5.00009 5)) SRID=4326;POLYGON ((5.00009 5, 5.00006 4.99994, 5 4.99991, 4.99994 4.99994, 4.99991 5, 4.99994 5.00006, 5 5.00009, 5.00006 5.00006, 5.00009 5)) SRID=4326;POLYGON EMPTY +Line going through left and right square SRID=4326;POLYGON ((0.5 0.50009, 0.50002 0.50009, 0.50003 0.50008, 0.50005 0.50008, 0.50006 0.50006, 0.50007 0.50005, 0.50008 0.50003, 0.50009 0.50002, 0.50009 0.5, 0.50009 0.49998, 0.50008 0.49997, 0.50007 0.49995, 0.50006 0.49994, 0.50005 0.49992, 0.50003 0.49992, 0.50002 0.49991, 0.5 0.49991, -0.5 0.49991, -0.50002 0.49991, -0.50003 0.49992, -0.50005 0.49992, -0.50006 0.49994, -0.50007 0.49995, -0.50008 0.49997, -0.50009 0.49998, -0.50009 0.5, -0.50009 0.50002, -0.50008 0.50003, -0.50007 0.50005, -0.50006 0.50006, -0.50005 0.50008, -0.50003 0.50008, -0.50002 0.50009, -0.5 0.50009, 0.5 0.50009)) SRID=4326;POLYGON ((0.5 0.50009, 0.50006 0.50006, 0.50009 0.5, 0.50006 0.49994, 0.5 0.49991, -0.5 0.49991, -0.50006 0.49994, -0.50009 0.5, -0.50006 0.50006, -0.5 0.50009, 0.5 0.50009)) SRID=4326;POLYGON ((0.5 0.50009, 0.5 0.49991, -0.5 0.49991, -0.5 0.50009, 0.5 0.50009)) +NULL NULL NULL NULL +Point middle of Left Square SRID=4326;POLYGON ((-0.49991 0.5, -0.49991 0.49998, -0.49992 0.49997, -0.49993 0.49995, -0.49994 0.49994, -0.49995 0.49992, -0.49997 0.49992, -0.49998 0.49991, -0.5 0.49991, -0.50002 0.49991, -0.50003 0.49992, -0.50005 0.49992, -0.50006 0.49994, -0.50007 0.49995, -0.50008 0.49997, -0.50009 0.49998, -0.50009 0.5, -0.50009 0.50002, -0.50008 0.50003, -0.50007 0.50005, -0.50006 0.50006, -0.50005 0.50008, -0.50003 0.50008, -0.50002 0.50009, -0.5 0.50009, -0.49998 0.50009, -0.49997 0.50008, -0.49995 0.50008, -0.49994 0.50006, -0.49993 0.50005, -0.49992 0.50003, -0.49991 0.50002, -0.49991 0.5)) SRID=4326;POLYGON ((-0.49991 0.5, -0.49994 0.49994, -0.5 0.49991, -0.50006 0.49994, -0.50009 0.5, -0.50006 0.50006, -0.5 0.50009, -0.49994 0.50006, -0.49991 0.5)) SRID=4326;POLYGON EMPTY +Point middle of Right Square SRID=4326;POLYGON ((0.50009 0.5, 0.50009 0.49998, 0.50008 0.49997, 0.50007 0.49995, 0.50006 0.49994, 0.50005 0.49992, 0.50003 0.49992, 0.50002 0.49991, 0.5 0.49991, 0.49998 0.49991, 0.49997 0.49992, 0.49995 0.49992, 0.49994 0.49994, 0.49993 0.49995, 0.49992 0.49997, 0.49991 0.49998, 0.49991 0.5, 0.49991 0.50002, 0.49992 0.50003, 0.49993 0.50005, 0.49994 0.50006, 0.49995 0.50008, 0.49997 0.50008, 0.49998 0.50009, 0.5 0.50009, 0.50002 0.50009, 0.50003 0.50008, 0.50005 0.50008, 0.50006 0.50006, 0.50007 0.50005, 0.50008 0.50003, 0.50009 0.50002, 0.50009 0.5)) SRID=4326;POLYGON ((0.50009 0.5, 0.50006 0.49994, 0.5 0.49991, 0.49994 0.49994, 0.49991 0.5, 0.49994 0.50006, 0.5 0.50009, 0.50006 0.50006, 0.50009 0.5)) SRID=4326;POLYGON EMPTY +Square (left) SRID=4326;POLYGON ((-1.00009 -0, -1.00009 1, -1.00009 1.00002, -1.00008 1.00003, -1.00007 1.00005, -1.00006 1.00006, -1.00005 1.00008, -1.00003 1.00008, -1.00002 1.00009, -1 1.00009, 0 1.00009, 0.00002 1.00009, 0.00003 1.00008, 0.00005 1.00008, 0.00006 1.00006, 0.00007 1.00005, 0.00008 1.00003, 0.00009 1.00002, 0.00009 1, 0.00009 0, 0.00009 -0.00002, 0.00008 -0.00003, 0.00007 -0.00005, 0.00006 -0.00006, 0.00005 -0.00008, 0.00003 -0.00008, 0.00002 -0.00009, 0 -0.00009, -1 -0.00009, -1.00002 -0.00009, -1.00003 -0.00008, -1.00005 -0.00008, -1.00006 -0.00006, -1.00007 -0.00005, -1.00008 -0.00003, -1.00009 -0.00002, -1.00009 -0)) SRID=4326;POLYGON ((-1.00009 -0, -1.00009 1, -1.00006 1.00006, -1 1.00009, 0 1.00009, 0.00006 1.00006, 0.00009 1, 0.00009 0, 0.00006 -0.00006, 0 -0.00009, -1 -0.00009, -1.00006 -0.00006, -1.00009 -0)) SRID=4326;POLYGON ((-1.00009 -0, -1.00009 1, -1.00008 1.00003, -1.00006 1.00006, -1.00003 1.00008, -1 1.00009, 0 1.00009, 0.00003 1.00008, 0.00006 1.00006, 0.00008 1.00003, 0.00009 1, 0.00009 0, 0.00008 -0.00003, 0.00006 -0.00006, 0.00003 -0.00008, 0 -0.00009, -1 -0.00009, -1.00003 -0.00008, -1.00006 -0.00006, -1.00008 -0.00003, -1.00009 -0)) +Square (right) SRID=4326;POLYGON ((-0.00009 0, -0.00009 1, -0.00009 1.00002, -0.00008 1.00003, -0.00007 1.00005, -0.00006 1.00006, -0.00005 1.00008, -0.00003 1.00008, -0.00002 1.00009, -0 1.00009, 1 1.00009, 1.00002 1.00009, 1.00003 1.00008, 1.00005 1.00008, 1.00006 1.00006, 1.00007 1.00005, 1.00008 1.00003, 1.00009 1.00002, 1.00009 1, 1.00009 -0, 1.00009 -0.00002, 1.00008 -0.00003, 1.00007 -0.00005, 1.00006 -0.00006, 1.00005 -0.00008, 1.00003 -0.00008, 1.00002 -0.00009, 1 -0.00009, -0 -0.00009, -0.00002 -0.00009, -0.00003 -0.00008, -0.00005 -0.00008, -0.00006 -0.00006, -0.00007 -0.00005, -0.00008 -0.00003, -0.00009 -0.00002, -0.00009 0)) SRID=4326;POLYGON ((-0.00009 0, -0.00009 1, -0.00006 1.00006, -0 1.00009, 1 1.00009, 1.00006 1.00006, 1.00009 1, 1.00009 -0, 1.00006 -0.00006, 1 -0.00009, -0 -0.00009, -0.00006 -0.00006, -0.00009 0)) SRID=4326;POLYGON ((-0.00009 0, -0.00009 1, -0.00008 1.00003, -0.00006 1.00006, -0.00003 1.00008, -0 1.00009, 1 1.00009, 1.00003 1.00008, 1.00006 1.00006, 1.00008 1.00003, 1.00009 1, 1.00009 -0, 1.00008 -0.00003, 1.00006 -0.00006, 1.00003 -0.00008, 1 -0.00009, -0 -0.00009, -0.00003 -0.00008, -0.00006 -0.00006, -0.00008 -0.00003, -0.00009 0)) +Square overlapping left and right square SRID=4326;POLYGON ((-0.10009 0, -0.10009 1, -0.10009 1.00002, -0.10008 1.00003, -0.10007 1.00005, -0.10006 1.00006, -0.10005 1.00008, -0.10003 1.00008, -0.10002 1.00009, -0.1 1.00009, 1 1.00009, 1.00002 1.00009, 1.00003 1.00008, 1.00005 1.00008, 1.00006 1.00006, 1.00007 1.00005, 1.00008 1.00003, 1.00009 1.00002, 1.00009 1, 1.00009 -0, 1.00009 -0.00002, 1.00008 -0.00003, 1.00007 -0.00005, 1.00006 -0.00006, 1.00005 -0.00008, 1.00003 -0.00008, 1.00002 -0.00009, 1 -0.00009, -0.1 -0.00009, -0.10002 -0.00009, -0.10003 -0.00008, -0.10005 -0.00008, -0.10006 -0.00006, -0.10007 -0.00005, -0.10008 -0.00003, -0.10009 -0.00002, -0.10009 0)) SRID=4326;POLYGON ((-0.10009 0, -0.10009 1, -0.10006 1.00006, -0.1 1.00009, 1 1.00009, 1.00006 1.00006, 1.00009 1, 1.00009 -0, 1.00006 -0.00006, 1 -0.00009, -0.1 -0.00009, -0.10006 -0.00006, -0.10009 0)) SRID=4326;POLYGON ((-0.10009 0, -0.10009 1, -0.10008 1.00003, -0.10006 1.00006, -0.10003 1.00008, -0.1 1.00009, 1 1.00009, 1.00003 1.00008, 1.00006 1.00006, 1.00008 1.00003, 1.00009 1, 1.00009 -0, 1.00008 -0.00003, 1.00006 -0.00006, 1.00003 -0.00008, 1 -0.00009, -0.1 -0.00009, -0.10003 -0.00008, -0.10006 -0.00006, -0.10008 -0.00003, -0.10009 0)) + +# ST_Intersection for geography +query TT +SELECT + dsc, + st_asewkt(st_intersection(a, b), 5) +FROM ( VALUES + ('empty', 'POINT EMPTY'::geography, 'POINT EMPTY'::geography), + ('non intersecting', 'POINT (1.5 1.5)'::geography, 'POINT (1.6 1.6)'::geography), + ('intersecting', 'LINESTRING (0 0, 1 1)'::geography, 'LINESTRING (0 1, 1 0)'::geography) +) t(dsc, a, b) +---- +empty SRID=4326;GEOMETRYCOLLECTION EMPTY +non intersecting SRID=4326;POINT EMPTY +intersecting SRID=4326;POINT (0.50019 0.50006) + # ST_Centroid for Geography query TTT SELECT diff --git a/pkg/sql/sem/builtins/geo_builtins.go b/pkg/sql/sem/builtins/geo_builtins.go index 25274c40a73e..706d7b408fc9 100644 --- a/pkg/sql/sem/builtins/geo_builtins.go +++ b/pkg/sql/sem/builtins/geo_builtins.go @@ -292,6 +292,56 @@ This variant approximates the circle into quad_seg segments per line (the defaul libraryUsage: usesGEOS, } +var usingBestGeomProjectionWarning = ` + +This operation is done by transforming the object into a Geometry. This occurs by translating +the Geography objects into Geometry objects before applying an LAEA, UTM or Web Mercator +based projection based on the bounding boxes of the given Geography objects. When the result is +calculated, the result is transformed back into a Geography with SRID 4326.` + +// performGeographyOperationUsingBestGeomProjection performs an operation on a +// Geography by transforming it to a relevant Geometry SRID and applying the closure, +// before retransforming it back into a geopb.DefaultGeographySRID geometry. +func performGeographyOperationUsingBestGeomProjection( + g *geo.Geography, f func(*geo.Geometry) (*geo.Geometry, error), +) (*geo.Geography, error) { + proj, err := geogfn.BestGeomProjection(g.BoundingRect()) + if err != nil { + return nil, err + } + + inLatLonGeom, err := g.AsGeometry() + if err != nil { + return nil, err + } + + inProjectedGeom, err := geotransform.Transform( + inLatLonGeom, + geoprojbase.Projections[g.SRID()].Proj4Text, + proj, + g.SRID(), + ) + if err != nil { + return nil, err + } + + outProjectedGeom, err := f(inProjectedGeom) + if err != nil { + return nil, err + } + + outGeom, err := geotransform.Transform( + outProjectedGeom, + proj, + geoprojbase.Projections[geopb.DefaultGeographySRID].Proj4Text, + geopb.DefaultGeographySRID, + ) + if err != nil { + return nil, err + } + return outGeom.AsGeography() +} + // fitMaxDecimalDigitsToBounds ensures maxDecimalDigits falls within the bounds that // is permitted by strconv.FormatFloat. func fitMaxDecimalDigitsToBounds(maxDecimalDigits int) int { @@ -2721,6 +2771,68 @@ For flags=1, validity considers self-intersecting rings forming holes as valid a }, tree.VolatilityImmutable, ), + geographyOverload2( + func(ctx *tree.EvalContext, a *tree.DGeography, b *tree.DGeography) (tree.Datum, error) { + proj, err := geogfn.BestGeomProjection(a.Geography.BoundingRect().Union(b.Geography.BoundingRect())) + if err != nil { + return nil, err + } + + aInGeom, err := a.Geography.AsGeometry() + if err != nil { + return nil, err + } + bInGeom, err := b.Geography.AsGeometry() + if err != nil { + return nil, err + } + + aInProjected, err := geotransform.Transform( + aInGeom, + geoprojbase.Projections[a.Geography.SRID()].Proj4Text, + proj, + a.Geography.SRID(), + ) + if err != nil { + return nil, err + } + bInProjected, err := geotransform.Transform( + bInGeom, + geoprojbase.Projections[b.Geography.SRID()].Proj4Text, + proj, + b.Geography.SRID(), + ) + if err != nil { + return nil, err + } + + projectedIntersection, err := geomfn.Intersection(aInProjected, bInProjected) + if err != nil { + return nil, err + } + + outGeom, err := geotransform.Transform( + projectedIntersection, + proj, + geoprojbase.Projections[geopb.DefaultGeographySRID].Proj4Text, + geopb.DefaultGeographySRID, + ) + if err != nil { + return nil, err + } + ret, err := outGeom.AsGeography() + if err != nil { + return nil, err + } + return tree.NewDGeography(ret), nil + }, + types.Geography, + infoBuilder{ + info: "Returns the point intersections of the given geographies." + usingBestGeomProjectionWarning, + libraryUsage: usesGEOS, + }, + tree.VolatilityImmutable, + ), ), "st_union": makeBuiltin( defProps(), @@ -3079,6 +3191,96 @@ The calculations are done on a sphere.`, Info: stBufferWithParamsInfoBuilder.String(), Volatility: tree.VolatilityImmutable, }, + tree.Overload{ + Types: tree.ArgTypes{ + {"geography", types.Geography}, + {"distance", types.Float}, + }, + ReturnType: tree.FixedReturnType(types.Geography), + Fn: func(_ *tree.EvalContext, args tree.Datums) (tree.Datum, error) { + g := args[0].(*tree.DGeography) + distance := *args[1].(*tree.DFloat) + + ret, err := performGeographyOperationUsingBestGeomProjection( + g.Geography, + func(g *geo.Geometry) (*geo.Geometry, error) { + return geomfn.Buffer(g, geomfn.MakeDefaultBufferParams(), float64(distance)) + }, + ) + if err != nil { + return nil, err + } + + return tree.NewDGeography(ret), nil + }, + Info: stBufferInfoBuilder.String() + usingBestGeomProjectionWarning, + Volatility: tree.VolatilityImmutable, + }, + tree.Overload{ + Types: tree.ArgTypes{ + {"geography", types.Geography}, + {"distance", types.Float}, + {"quad_segs", types.Int}, + }, + ReturnType: tree.FixedReturnType(types.Geography), + Fn: func(_ *tree.EvalContext, args tree.Datums) (tree.Datum, error) { + g := args[0].(*tree.DGeography) + distance := *args[1].(*tree.DFloat) + quadSegs := *args[2].(*tree.DInt) + + ret, err := performGeographyOperationUsingBestGeomProjection( + g.Geography, + func(g *geo.Geometry) (*geo.Geometry, error) { + return geomfn.Buffer( + g, + geomfn.MakeDefaultBufferParams().WithQuadrantSegments(int(quadSegs)), + float64(distance), + ) + }, + ) + if err != nil { + return nil, err + } + return tree.NewDGeography(ret), nil + }, + Info: stBufferWithQuadSegInfoBuilder.String() + usingBestGeomProjectionWarning, + Volatility: tree.VolatilityImmutable, + }, + tree.Overload{ + Types: tree.ArgTypes{ + {"geography", types.Geography}, + {"distance", types.Float}, + {"buffer_style_params", types.String}, + }, + ReturnType: tree.FixedReturnType(types.Geography), + Fn: func(_ *tree.EvalContext, args tree.Datums) (tree.Datum, error) { + g := args[0].(*tree.DGeography) + distance := *args[1].(*tree.DFloat) + paramsString := *args[2].(*tree.DString) + + params, modifiedDistance, err := geomfn.ParseBufferParams(string(paramsString), float64(distance)) + if err != nil { + return nil, err + } + + ret, err := performGeographyOperationUsingBestGeomProjection( + g.Geography, + func(g *geo.Geometry) (*geo.Geometry, error) { + return geomfn.Buffer( + g, + params, + modifiedDistance, + ) + }, + ) + if err != nil { + return nil, err + } + return tree.NewDGeography(ret), nil + }, + Info: stBufferWithParamsInfoBuilder.String() + usingBestGeomProjectionWarning, + Volatility: tree.VolatilityImmutable, + }, ), "st_envelope": makeBuiltin( defProps(), @@ -3628,8 +3830,7 @@ func initGeoBuiltins() { "st_distance", "st_dwithin", "st_intersects", - // TODO(#48398): uncomment - // "st_intersection", + "st_intersection", "st_length", } { builtin, exists := geoBuiltins[builtinName]