github.com/cockroachdb/cockroachdb-parser@v0.23.3-0.20240213214944-911057d40c9a/pkg/geo/encode.go (about)

     1  // Copyright 2020 The Cockroach Authors.
     2  //
     3  // Use of this software is governed by the Business Source License
     4  // included in the file licenses/BSL.txt.
     5  //
     6  // As of the Change Date specified in that file, in accordance with
     7  // the Business Source License, use of this software will be governed
     8  // by the Apache License, Version 2.0, included in the file
     9  // licenses/APL.txt.
    10  
    11  package geo
    12  
    13  import (
    14  	"bytes"
    15  	"encoding/binary"
    16  	"fmt"
    17  	"strings"
    18  
    19  	"github.com/cockroachdb/cockroachdb-parser/pkg/geo/geopb"
    20  	"github.com/cockroachdb/cockroachdb-parser/pkg/geo/geoprojbase"
    21  	"github.com/cockroachdb/cockroachdb-parser/pkg/sql/pgwire/pgcode"
    22  	"github.com/cockroachdb/cockroachdb-parser/pkg/sql/pgwire/pgerror"
    23  	"github.com/golang/geo/s1"
    24  	"github.com/pierrre/geohash"
    25  	"github.com/twpayne/go-geom"
    26  	"github.com/twpayne/go-geom/encoding/ewkb"
    27  	"github.com/twpayne/go-geom/encoding/geojson"
    28  	"github.com/twpayne/go-geom/encoding/kml"
    29  	"github.com/twpayne/go-geom/encoding/wkb"
    30  	"github.com/twpayne/go-geom/encoding/wkbcommon"
    31  	"github.com/twpayne/go-geom/encoding/wkbhex"
    32  	"github.com/twpayne/go-geom/encoding/wkt"
    33  )
    34  
    35  // DefaultGeoJSONDecimalDigits is the default number of digits coordinates in GeoJSON.
    36  const DefaultGeoJSONDecimalDigits = 9
    37  
    38  // SpatialObjectToWKT transforms a given SpatialObject to WKT.
    39  func SpatialObjectToWKT(so geopb.SpatialObject, maxDecimalDigits int) (geopb.WKT, error) {
    40  	t, err := ewkb.Unmarshal([]byte(so.EWKB))
    41  	if err != nil {
    42  		return "", err
    43  	}
    44  	ret, err := wkt.Marshal(t, wkt.EncodeOptionWithMaxDecimalDigits(maxDecimalDigits))
    45  	return geopb.WKT(ret), err
    46  }
    47  
    48  // SpatialObjectToEWKT transforms a given SpatialObject to EWKT.
    49  func SpatialObjectToEWKT(so geopb.SpatialObject, maxDecimalDigits int) (geopb.EWKT, error) {
    50  	t, err := ewkb.Unmarshal([]byte(so.EWKB))
    51  	if err != nil {
    52  		return "", err
    53  	}
    54  	ret, err := wkt.Marshal(t, wkt.EncodeOptionWithMaxDecimalDigits(maxDecimalDigits))
    55  	if err != nil {
    56  		return "", err
    57  	}
    58  	if t.SRID() != 0 {
    59  		ret = fmt.Sprintf("SRID=%d;%s", t.SRID(), ret)
    60  	}
    61  	return geopb.EWKT(ret), err
    62  }
    63  
    64  // SpatialObjectToWKB transforms a given SpatialObject to WKB.
    65  func SpatialObjectToWKB(so geopb.SpatialObject, byteOrder binary.ByteOrder) (geopb.WKB, error) {
    66  	t, err := ewkb.Unmarshal([]byte(so.EWKB))
    67  	if err != nil {
    68  		return nil, err
    69  	}
    70  	ret, err := wkb.Marshal(t, byteOrder, wkbcommon.WKBOptionEmptyPointHandling(wkbcommon.EmptyPointHandlingNaN))
    71  	return geopb.WKB(ret), err
    72  }
    73  
    74  // SpatialObjectToGeoJSONFlag maps to the ST_AsGeoJSON flags for PostGIS.
    75  type SpatialObjectToGeoJSONFlag int
    76  
    77  // These should be kept with ST_AsGeoJSON in PostGIS.
    78  // 0: means no option
    79  // 1: GeoJSON BBOX
    80  // 2: GeoJSON Short CRS (e.g EPSG:4326)
    81  // 4: GeoJSON Long CRS (e.g urn:ogc:def:crs:EPSG::4326)
    82  // 8: GeoJSON Short CRS if not EPSG:4326 (default)
    83  const (
    84  	SpatialObjectToGeoJSONFlagIncludeBBox SpatialObjectToGeoJSONFlag = 1 << (iota)
    85  	SpatialObjectToGeoJSONFlagShortCRS
    86  	SpatialObjectToGeoJSONFlagLongCRS
    87  	SpatialObjectToGeoJSONFlagShortCRSIfNot4326
    88  
    89  	SpatialObjectToGeoJSONFlagZero = 0
    90  )
    91  
    92  // geomToGeoJSONCRS converts a geom to its CRS GeoJSON form.
    93  func geomToGeoJSONCRS(t geom.T, long bool) (*geojson.CRS, error) {
    94  	projection, err := geoprojbase.Projection(geopb.SRID(t.SRID()))
    95  	if err != nil {
    96  		return nil, err
    97  	}
    98  	var prop string
    99  	if long {
   100  		prop = fmt.Sprintf("urn:ogc:def:crs:%s::%d", projection.AuthName, projection.AuthSRID)
   101  	} else {
   102  		prop = fmt.Sprintf("%s:%d", projection.AuthName, projection.AuthSRID)
   103  	}
   104  	crs := &geojson.CRS{
   105  		Type: "name",
   106  		Properties: map[string]interface{}{
   107  			"name": prop,
   108  		},
   109  	}
   110  	return crs, nil
   111  }
   112  
   113  // SpatialObjectToGeoJSON transforms a given SpatialObject to GeoJSON.
   114  func SpatialObjectToGeoJSON(
   115  	so geopb.SpatialObject, maxDecimalDigits int, flag SpatialObjectToGeoJSONFlag,
   116  ) ([]byte, error) {
   117  	t, err := ewkb.Unmarshal([]byte(so.EWKB))
   118  	if err != nil {
   119  		return nil, err
   120  	}
   121  	options := []geojson.EncodeGeometryOption{
   122  		geojson.EncodeGeometryWithMaxDecimalDigits(maxDecimalDigits),
   123  	}
   124  	if flag&SpatialObjectToGeoJSONFlagIncludeBBox != 0 {
   125  		// Do not encoding empty bounding boxes.
   126  		if so.BoundingBox != nil {
   127  			options = append(
   128  				options,
   129  				geojson.EncodeGeometryWithBBox(),
   130  			)
   131  		}
   132  	}
   133  	// Take CRS flag in order of precedence.
   134  	if t.SRID() != 0 {
   135  		if flag&SpatialObjectToGeoJSONFlagLongCRS != 0 {
   136  			crs, err := geomToGeoJSONCRS(t, true /* long */)
   137  			if err != nil {
   138  				return nil, err
   139  			}
   140  			options = append(options, geojson.EncodeGeometryWithCRS(crs))
   141  		} else if flag&SpatialObjectToGeoJSONFlagShortCRS != 0 {
   142  			crs, err := geomToGeoJSONCRS(t, false /* long */)
   143  			if err != nil {
   144  				return nil, err
   145  			}
   146  			options = append(options, geojson.EncodeGeometryWithCRS(crs))
   147  		} else if flag&SpatialObjectToGeoJSONFlagShortCRSIfNot4326 != 0 {
   148  			if t.SRID() != 4326 {
   149  				crs, err := geomToGeoJSONCRS(t, false /* long */)
   150  				if err != nil {
   151  					return nil, err
   152  				}
   153  				options = append(options, geojson.EncodeGeometryWithCRS(crs))
   154  			}
   155  		}
   156  	}
   157  
   158  	return geojson.Marshal(t, options...)
   159  }
   160  
   161  // SpatialObjectToWKBHex transforms a given SpatialObject to WKBHex.
   162  func SpatialObjectToWKBHex(so geopb.SpatialObject) (string, error) {
   163  	t, err := ewkb.Unmarshal([]byte(so.EWKB))
   164  	if err != nil {
   165  		return "", err
   166  	}
   167  	ret, err := wkbhex.Encode(t, DefaultEWKBEncodingFormat, wkbcommon.WKBOptionEmptyPointHandling(wkbcommon.EmptyPointHandlingNaN))
   168  	return strings.ToUpper(ret), err
   169  }
   170  
   171  // SpatialObjectToKML transforms a given SpatialObject to KML.
   172  func SpatialObjectToKML(so geopb.SpatialObject) (string, error) {
   173  	t, err := ewkb.Unmarshal([]byte(so.EWKB))
   174  	if err != nil {
   175  		return "", err
   176  	}
   177  	kmlElement, err := kml.Encode(t)
   178  	if err != nil {
   179  		return "", err
   180  	}
   181  	var buf bytes.Buffer
   182  	if err := kmlElement.Write(&buf); err != nil {
   183  		return "", err
   184  	}
   185  	return buf.String(), nil
   186  }
   187  
   188  // GeoHashAutoPrecision means to calculate the precision of SpatialObjectToGeoHash
   189  // based on input, up to 32 characters.
   190  const GeoHashAutoPrecision = 0
   191  
   192  // GeoHashMaxPrecision is the maximum precision for GeoHashes.
   193  // 20 is picked as doubles have 51 decimals of precision, and each base32 position
   194  // can contain 5 bits of data. As we have two points, we use floor((2 * 51) / 5) = 20.
   195  const GeoHashMaxPrecision = 20
   196  
   197  // SpatialObjectToGeoHash transforms a given SpatialObject to a GeoHash.
   198  func SpatialObjectToGeoHash(so geopb.SpatialObject, p int) (string, error) {
   199  	if so.BoundingBox == nil {
   200  		return "", nil
   201  	}
   202  	bbox := so.BoundingBox
   203  	if so.Type == geopb.SpatialObjectType_GeographyType {
   204  		// Convert bounding box back to degrees.
   205  		bbox = &geopb.BoundingBox{
   206  			LoX: s1.Angle(bbox.LoX).Degrees(),
   207  			HiX: s1.Angle(bbox.HiX).Degrees(),
   208  			LoY: s1.Angle(bbox.LoY).Degrees(),
   209  			HiY: s1.Angle(bbox.HiY).Degrees(),
   210  		}
   211  	}
   212  	if bbox.LoX < -180 || bbox.HiX > 180 || bbox.LoY < -90 || bbox.HiY > 90 {
   213  		return "", pgerror.Newf(
   214  			pgcode.InvalidParameterValue,
   215  			"object has bounds greater than the bounds of lat/lng, got (%f %f, %f %f)",
   216  			bbox.LoX, bbox.LoY,
   217  			bbox.HiX, bbox.HiY,
   218  		)
   219  	}
   220  
   221  	// Get precision using the bounding box if required.
   222  	if p <= GeoHashAutoPrecision {
   223  		p = getPrecisionForBBox(bbox)
   224  	}
   225  
   226  	// Support up to 20, which is the same as PostGIS.
   227  	if p > GeoHashMaxPrecision {
   228  		p = GeoHashMaxPrecision
   229  	}
   230  
   231  	bbCenterLng := bbox.LoX + (bbox.HiX-bbox.LoX)/2.0
   232  	bbCenterLat := bbox.LoY + (bbox.HiY-bbox.LoY)/2.0
   233  
   234  	return geohash.Encode(bbCenterLat, bbCenterLng, p), nil
   235  }
   236  
   237  // getPrecisionForBBox is a function imitating PostGIS's ability to go from
   238  // a world bounding box and truncating a GeoHash to fit the given bounding box.
   239  // The algorithm halves the world bounding box until it intersects with the
   240  // feature bounding box to get a precision that will encompass the entire
   241  // bounding box.
   242  func getPrecisionForBBox(bbox *geopb.BoundingBox) int {
   243  	bitPrecision := 0
   244  
   245  	// This is a point, for points we use the full bitPrecision.
   246  	if bbox.LoX == bbox.HiX && bbox.LoY == bbox.HiY {
   247  		return GeoHashMaxPrecision
   248  	}
   249  
   250  	// Starts from a world bounding box:
   251  	lonMin := -180.0
   252  	lonMax := 180.0
   253  	latMin := -90.0
   254  	latMax := 90.0
   255  
   256  	// Each iteration shrinks the world bounding box by half in the dimension that
   257  	// does not fit, making adjustments each iteration until it intersects with
   258  	// the object bbox.
   259  	for {
   260  		lonWidth := lonMax - lonMin
   261  		latWidth := latMax - latMin
   262  		latMaxDelta, lonMaxDelta, latMinDelta, lonMinDelta := 0.0, 0.0, 0.0, 0.0
   263  
   264  		// Look at whether the longitudes of the bbox are to the left or
   265  		// the right of the world bbox longitudes, shrinks it and makes adjustments
   266  		// for the next iteration.
   267  		if bbox.LoX > lonMin+lonWidth/2.0 {
   268  			lonMinDelta = lonWidth / 2.0
   269  		} else if bbox.HiX < lonMax-lonWidth/2.0 {
   270  			lonMaxDelta = lonWidth / -2.0
   271  		}
   272  		// Look at whether the latitudes of the bbox are to the left or
   273  		// the right of the world bbox latitudes, shrinks it and makes adjustments
   274  		// for the next iteration.
   275  		if bbox.LoY > latMin+latWidth/2.0 {
   276  			latMinDelta = latWidth / 2.0
   277  		} else if bbox.HiY < latMax-latWidth/2.0 {
   278  			latMaxDelta = latWidth / -2.0
   279  		}
   280  
   281  		// Every change we make that splits the box up adds precision.
   282  		// If we detect no change, we've intersected a box and so must exit.
   283  		precisionDelta := 0
   284  		if lonMinDelta != 0.0 || lonMaxDelta != 0.0 {
   285  			lonMin += lonMinDelta
   286  			lonMax += lonMaxDelta
   287  			precisionDelta++
   288  		} else {
   289  			break
   290  		}
   291  		if latMinDelta != 0.0 || latMaxDelta != 0.0 {
   292  			latMin += latMinDelta
   293  			latMax += latMaxDelta
   294  			precisionDelta++
   295  		} else {
   296  			break
   297  		}
   298  		bitPrecision += precisionDelta
   299  	}
   300  	// Each character can represent 5 bits of bitPrecision.
   301  	// As such, divide by 5 to get GeoHash precision.
   302  	return bitPrecision / 5
   303  }
   304  
   305  // StringToByteOrder returns the byte order of string.
   306  func StringToByteOrder(s string) binary.ByteOrder {
   307  	switch strings.ToLower(s) {
   308  	case "ndr":
   309  		return binary.LittleEndian
   310  	case "xdr":
   311  		return binary.BigEndian
   312  	default:
   313  		return DefaultEWKBEncodingFormat
   314  	}
   315  }