Skip to content

Commit

Permalink
Merge #49742
Browse files Browse the repository at this point in the history
49742: geo/geomfn: implement ST_LineInterpolatePoint(s) r=otan a=abhishek20123g

Fixes #48971
Fixes #48972

This PR adds following builtin functions
* ST_LineInterpolatePoint{{geometry, float8}}
* ST_LineInterpolatePoints{{geometry, float8, bool}}

which works for LineString only, allows us to determine one or more
interpolated points in the LineString which at an integral multiple
of given fraction of LineString's total length.

Release note (sql change): This PR implement adds the following built-in functions.
* ST_LineInterpolatePoint{{geometry, float8}}
* ST_LineInterpolatePoints{{geometry, float8, bool}}

Co-authored-by: abhishek20123g <abhishek20123g@gmail.com>
  • Loading branch information
craig[bot] and abhishek20123g committed Jun 2, 2020
2 parents d65b224 + cf81e9f commit cbe7fdc
Show file tree
Hide file tree
Showing 8 changed files with 397 additions and 1 deletion.
11 changes: 11 additions & 0 deletions docs/generated/sql/functions.md
Original file line number Diff line number Diff line change
Expand Up @@ -960,6 +960,17 @@ given Geometry.</p>
</span></td></tr>
<tr><td><a name="st_linefromwkb"></a><code>st_linefromwkb(wkb: <a href="bytes.html">bytes</a>, srid: <a href="int.html">int</a>) &rarr; geometry</code></td><td><span class="funcdesc"><p>Returns the Geometry from a WKB representation with an SRID. If the shape underneath is not LineString, NULL is returned.</p>
</span></td></tr>
<tr><td><a name="st_lineinterpolatepoint"></a><code>st_lineinterpolatepoint(geometry: geometry, fraction: <a href="float.html">float</a>) &rarr; geometry</code></td><td><span class="funcdesc"><p>Returns a point along the given LineString which is at given fraction of LineString’s total length.</p>
<p>This function utilizes the GEOS module.</p>
</span></td></tr>
<tr><td><a name="st_lineinterpolatepoints"></a><code>st_lineinterpolatepoints(geometry: geometry, fraction: <a href="float.html">float</a>) &rarr; geometry</code></td><td><span class="funcdesc"><p>Returns one or more points along the LineString which is at an integral multiples of given fraction of LineString’s total length.</p>
<p>Note If the result has zero or one points, it will be returned as a POINT. If it has two or more points, it will be returned as a MULTIPOINT.</p>
<p>This function utilizes the GEOS module.</p>
</span></td></tr>
<tr><td><a name="st_lineinterpolatepoints"></a><code>st_lineinterpolatepoints(geometry: geometry, fraction: <a href="float.html">float</a>, repeat: <a href="bool.html">bool</a>) &rarr; geometry</code></td><td><span class="funcdesc"><p>Returns one or more points along the LineString which is at an integral multiples of given fraction of LineString’s total length. If repeat is false (default true) then it returns first point.</p>
<p>Note If the result has zero or one points, it will be returned as a POINT. If it has two or more points, it will be returned as a MULTIPOINT.</p>
<p>This function utilizes the GEOS module.</p>
</span></td></tr>
<tr><td><a name="st_linestringfromtext"></a><code>st_linestringfromtext(str: <a href="string.html">string</a>, srid: <a href="int.html">int</a>) &rarr; geometry</code></td><td><span class="funcdesc"><p>Returns the Geometry from a WKT or EWKT representation with an SRID. If the shape underneath is not LineString, NULL is returned. If the SRID is present in both the EWKT and the argument, the argument value is used.</p>
</span></td></tr>
<tr><td><a name="st_linestringfromtext"></a><code>st_linestringfromtext(val: <a href="string.html">string</a>) &rarr; geometry</code></td><td><span class="funcdesc"><p>Returns the Geometry from a WKT or EWKT representation. If the shape underneath is not LineString, NULL is returned.</p>
Expand Down
71 changes: 71 additions & 0 deletions pkg/geo/geomfn/linear_reference.go
Original file line number Diff line number Diff line change
@@ -0,0 +1,71 @@
// 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 geomfn

import (
"github.com/cockroachdb/cockroach/pkg/geo"
"github.com/cockroachdb/cockroach/pkg/geo/geos"
"github.com/cockroachdb/errors"
"github.com/twpayne/go-geom"
"github.com/twpayne/go-geom/encoding/ewkb"
)

// LineInterpolatePoints returns one or more points along the given
// LineString which are at an integral multiples of given fraction of
// LineString's total length. When repeat is set to false, it returns
// the first point.
func LineInterpolatePoints(g *geo.Geometry, fraction float64, repeat bool) (*geo.Geometry, error) {
if fraction < 0 || fraction > 1 {
return nil, errors.Newf("fraction %f should be within [0 1] range", fraction)
}
geomRepr, err := g.AsGeomT()
if err != nil {
return nil, err
}
// Empty geometries do not react well in GEOS, so we have to
// convert and check beforehand.
// Remove after #49209 is resolved.
if geomRepr.Empty() {
return geo.NewGeometryFromGeom(geom.NewPointEmpty(geom.XY))
}
switch geomRepr := geomRepr.(type) {
case *geom.LineString:
// In case fraction is greater than 0.5 or equal to 0 or repeat is false,
// then we will have only one interpolated point.
lengthOfLineString := geomRepr.Length()
if repeat && fraction <= 0.5 && fraction != 0 {
numberOfInterpolatedPoints := int(1 / fraction)
interpolatedPoints := geom.NewMultiPoint(geom.XY).SetSRID(geomRepr.SRID())
for pointInserted := 1; pointInserted <= numberOfInterpolatedPoints; pointInserted++ {
pointEWKB, err := geos.InterpolateLine(g.EWKB(), float64(pointInserted)*fraction*lengthOfLineString)
if err != nil {
return nil, err
}
point, err := ewkb.Unmarshal(pointEWKB)
if err != nil {
return nil, err
}
err = interpolatedPoints.Push(point.(*geom.Point))
if err != nil {
return nil, err
}
}
return geo.NewGeometryFromGeom(interpolatedPoints)
}
interpolatedPointEWKB, err := geos.InterpolateLine(g.EWKB(), fraction*lengthOfLineString)
if err != nil {
return nil, err
}
return geo.ParseGeometryFromEWKB(interpolatedPointEWKB)
default:
return nil, errors.Newf("geometry %s should be LineString", g.Shape())
}
}
107 changes: 107 additions & 0 deletions pkg/geo/geomfn/linear_reference_test.go
Original file line number Diff line number Diff line change
@@ -0,0 +1,107 @@
// 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 geomfn

import (
"fmt"
"testing"

"github.com/cockroachdb/cockroach/pkg/geo"
"github.com/stretchr/testify/require"
)

func TestLineInterpolatePoints(t *testing.T) {
var testCasesForLineInterpolate = []struct {
wkb string
errMsg string
fraction float64
expectedWKTForRepeatTrue string
expectedWKTForRepeatFalse string
}{
{
wkb: "LINESTRING (0 0, 1 1)",
errMsg: "fraction -0.200000 should be within [0 1] range",
fraction: -0.2,
},
{
wkb: "LINESTRING (0 0, 1 1, 2 5)",
errMsg: "fraction 1.500000 should be within [0 1] range",
fraction: 1.5,
},
{
wkb: "MULTILINESTRING ((0 0, 1 1, 2 5), (0 0, 1 1))",
errMsg: "geometry MultiLineString should be LineString",
fraction: 0.3,
},
{
wkb: "POINT (0 0)",
errMsg: "geometry Point should be LineString",
fraction: 0.3,
},
{
wkb: "POLYGON((-1.0 0.0, 0.0 0.0, 0.0 1.0, -1.0 1.0, -1.0 0.0))",
errMsg: "geometry Polygon should be LineString",
fraction: 0.3,
},
{
wkb: "LINESTRING (0 0, 1 1, 2 5)",
fraction: 0.51,
expectedWKTForRepeatTrue: "POINT (1.3419313865603413 2.367725546241365)",
expectedWKTForRepeatFalse: "POINT (1.3419313865603413 2.367725546241365)",
},
{
wkb: "LINESTRING (0 0, 1 1, 2 5)",
fraction: 0.5,
expectedWKTForRepeatTrue: "MULTIPOINT (1.3285014148574912 2.3140056594299647, 2 5)",
expectedWKTForRepeatFalse: "POINT (1.3285014148574912 2.3140056594299647)",
},
{
wkb: "LINESTRING (0 0, 1 1, 2 5)",
fraction: 0.2,
expectedWKTForRepeatTrue: "MULTIPOINT (0.78309518948453 0.78309518948453, 1.1942016978289893 1.7768067913159575, 1.462801131885993 2.851204527543972, 1.7314005659429965 3.925602263771986, 2 5)",
expectedWKTForRepeatFalse: "POINT (0.78309518948453 0.78309518948453)",
},
{
wkb: "LINESTRING (0 0, 1 1, 2 5)",
fraction: 0,
expectedWKTForRepeatTrue: "POINT (0 0)",
expectedWKTForRepeatFalse: "POINT (0 0)",
},
{
wkb: "LINESTRING (0 0, 1 1, 2 5)",
fraction: 1,
expectedWKTForRepeatTrue: "POINT (2 5)",
expectedWKTForRepeatFalse: "POINT (2 5)",
},
}
for _, test := range testCasesForLineInterpolate {
for _, repeat := range []bool{false, true} {
t.Run(fmt.Sprintf("%s for fraction %f where repeat is %t", test.wkb, test.fraction, repeat),
func(t *testing.T) {
geometry, err := geo.ParseGeometry(test.wkb)
require.NoError(t, err)
interpolatedPoint, err := LineInterpolatePoints(geometry, test.fraction, repeat)
if test.errMsg == "" {
require.NoError(t, err)
expectedWKTForRepeat := test.expectedWKTForRepeatFalse
if repeat {
expectedWKTForRepeat = test.expectedWKTForRepeatTrue
}
expectedInterpolatedPoint, err := geo.ParseGeometry(expectedWKTForRepeat)
require.NoError(t, err)
require.Equal(t, expectedInterpolatedPoint, interpolatedPoint)
} else {
require.EqualError(t, err, test.errMsg)
}
})
}
}
}
34 changes: 33 additions & 1 deletion pkg/geo/geos/geos.cc
Original file line number Diff line number Diff line change
Expand Up @@ -83,6 +83,8 @@ typedef int (*CR_GEOS_Area_r)(CR_GEOS_Handle, CR_GEOS_Geometry, double*);
typedef int (*CR_GEOS_Length_r)(CR_GEOS_Handle, CR_GEOS_Geometry, double*);
typedef CR_GEOS_Geometry (*CR_GEOS_Centroid_r)(CR_GEOS_Handle, CR_GEOS_Geometry);

typedef CR_GEOS_Geometry (*CR_GEOS_Interpolate_r)(CR_GEOS_Handle, CR_GEOS_Geometry, double);

typedef int (*CR_GEOS_Distance_r)(CR_GEOS_Handle, CR_GEOS_Geometry, CR_GEOS_Geometry, double*);

typedef char (*CR_GEOS_Covers_r)(CR_GEOS_Handle, CR_GEOS_Geometry, CR_GEOS_Geometry);
Expand Down Expand Up @@ -147,6 +149,8 @@ struct CR_GEOS {
CR_GEOS_Length_r GEOSLength_r;
CR_GEOS_Centroid_r GEOSGetCentroid_r;

CR_GEOS_Interpolate_r GEOSInterpolate_r;

CR_GEOS_Distance_r GEOSDistance_r;

CR_GEOS_Covers_r GEOSCovers_r;
Expand Down Expand Up @@ -205,6 +209,7 @@ struct CR_GEOS {
INIT(GEOSArea_r);
INIT(GEOSLength_r);
INIT(GEOSGetCentroid_r);
INIT(GEOSInterpolate_r);
INIT(GEOSDistance_r);
INIT(GEOSCovers_r);
INIT(GEOSCoveredBy_r);
Expand Down Expand Up @@ -455,7 +460,34 @@ CR_GEOS_Status CR_GEOS_Centroid(CR_GEOS* lib, CR_GEOS_Slice a, CR_GEOS_String* c
return toGEOSString(error.data(), error.length());
}

CR_GEOS_Status CR_GEOS_Distance(CR_GEOS* lib, CR_GEOS_Slice a, CR_GEOS_Slice b, double* ret) {
//
// Linear Reference
//

CR_GEOS_Status CR_GEOS_Interpolate(CR_GEOS* lib, CR_GEOS_Slice a, double distance,
CR_GEOS_String* interpolatedPointEWKB) {
std::string error;
auto handle = initHandleWithErrorBuffer(lib, &error);
auto geom = CR_GEOS_GeometryFromSlice(lib, handle, a);
*interpolatedPointEWKB = {.data = NULL, .len = 0};
if (geom != nullptr) {
auto interpolatedPoint = lib->GEOSInterpolate_r(handle, geom, distance);
if (interpolatedPoint != nullptr) {
auto srid = lib->GEOSGetSRID_r(handle, geom);
CR_GEOS_writeGeomToEWKB(lib, handle, interpolatedPoint, interpolatedPointEWKB, srid);
lib->GEOSGeom_destroy_r(handle, interpolatedPoint);
}
lib->GEOSGeom_destroy_r(handle, geom);
}
lib->GEOS_finish_r(handle);
return toGEOSString(error.data(), error.length());
}

//
// Binary operators
//

CR_GEOS_Status CR_GEOS_Distance(CR_GEOS* lib, CR_GEOS_Slice a, CR_GEOS_Slice b, double *ret) {
return CR_GEOS_BinaryOperator(lib, lib->GEOSDistance_r, a, b, ret);
}

Expand Down
19 changes: 19 additions & 0 deletions pkg/geo/geos/geos.go
Original file line number Diff line number Diff line change
Expand Up @@ -320,6 +320,25 @@ func Centroid(ewkb geopb.EWKB) (geopb.EWKB, error) {
return cStringToSafeGoBytes(cEWKB), nil
}

// InterpolateLine returns the point along the given LineString which is at
// a given distance from starting point.
// Note: For distance less than 0 it returns start point similarly for distance
// greater LineString's length.
// InterpolateLine also works with (Multi)LineString. However, the result is
// not appropriate as it combines all the LineString present in (MULTI)LineString,
// considering all the corner points of LineString overlaps each other.
func InterpolateLine(ewkb geopb.EWKB, distance float64) (geopb.EWKB, error) {
g, err := ensureInitInternal()
if err != nil {
return nil, err
}
var cEWKB C.CR_GEOS_String
if err := statusToError(C.CR_GEOS_Interpolate(g, goToCSlice(ewkb), C.double(distance), &cEWKB)); err != nil {
return nil, err
}
return cStringToSafeGoBytes(cEWKB), nil
}

// MinDistance returns the minimum distance between two EWKBs.
func MinDistance(a geopb.EWKB, b geopb.EWKB) (float64, error) {
g, err := ensureInitInternal()
Expand Down
6 changes: 6 additions & 0 deletions pkg/geo/geos/geos.h
Original file line number Diff line number Diff line change
Expand Up @@ -77,6 +77,12 @@ CR_GEOS_Status CR_GEOS_Area(CR_GEOS* lib, CR_GEOS_Slice a, double* ret);
CR_GEOS_Status CR_GEOS_Length(CR_GEOS* lib, CR_GEOS_Slice a, double* ret);
CR_GEOS_Status CR_GEOS_Centroid(CR_GEOS* lib, CR_GEOS_Slice a, CR_GEOS_String* centroidEWKB);

//
// Linear reference.
//
CR_GEOS_Status CR_GEOS_Interpolate(CR_GEOS *lib, CR_GEOS_Slice a, double distance,
CR_GEOS_String *interpolatedPoint);

//
// Binary operators.
//
Expand Down
Loading

0 comments on commit cbe7fdc

Please sign in to comment.