Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

[SEDONA-299] Add ST_FrechetDistance #869

Merged
merged 6 commits into from
Jun 26, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
Original file line number Diff line number Diff line change
Expand Up @@ -955,6 +955,10 @@ public static Geometry geometricMedian(Geometry geometry) throws Exception {
return geometricMedian(geometry, DEFAULT_TOLERANCE, DEFAULT_MAX_ITER, false);
}

public static double frechetDistance(Geometry g1, Geometry g2) {
return GeomUtils.getFrechetDistance(g1, g2);
}

public static Geometry boundingDiagonal(Geometry geometry) {
if (geometry.isEmpty()) {
return GEOMETRY_FACTORY.createLineString();
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -26,6 +26,7 @@
import org.locationtech.jts.io.WKTWriter;
import org.locationtech.jts.operation.polygonize.Polygonizer;
import org.locationtech.jts.operation.union.UnaryUnionOp;
import org.locationtech.jts.algorithm.distance.DiscreteFrechetDistance;
import org.locationtech.jts.algorithm.distance.DiscreteHausdorffDistance;

import java.nio.ByteOrder;
Expand Down Expand Up @@ -424,7 +425,6 @@ public static Geometry[] getSubGeometries(Geometry geom) {
return geometries;
}


public static Geometry get3DGeom(Geometry geometry, double zValue) {
Coordinate[] coordinates = geometry.getCoordinates();
if (coordinates.length == 0) return geometry;
Expand Down Expand Up @@ -481,6 +481,11 @@ public static void affineGeom(Geometry geometry, Double a, Double b, Double d, D
geometry.geometryChanged();
}

public static double getFrechetDistance(Geometry g1, Geometry g2) {
if (g1.isEmpty() || g2.isEmpty()) return 0.0;
return DiscreteFrechetDistance.distance(g1, g2);
}

public static Double getHausdorffDistance(Geometry g1, Geometry g2, double densityFrac) throws Exception {
if (g1.isEmpty() || g2.isEmpty()) return 0.0;
DiscreteHausdorffDistance hausdorffDistanceObj = new DiscreteHausdorffDistance(g1, g2);
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -918,10 +918,10 @@ public void translateHybridGeomCollectionNoDeltaZ() {
public void translateHybridGeomCollectionDeltaZ() {
Polygon polygon = GEOMETRY_FACTORY.createPolygon(coordArray(1, 0, 1, 1, 2, 1, 2, 0, 1, 0));
Polygon polygon3D = GEOMETRY_FACTORY.createPolygon(coordArray3d(1, 0, 1, 2, 0, 2, 2, 1, 2, 1, 0, 1));
MultiPolygon multiPolygon = GEOMETRY_FACTORY.createMultiPolygon(new Polygon[] {polygon3D, polygon});
MultiPolygon multiPolygon = GEOMETRY_FACTORY.createMultiPolygon(new Polygon[]{polygon3D, polygon});
Point point3D = GEOMETRY_FACTORY.createPoint(new Coordinate(1, 1, 1));
LineString emptyLineString = GEOMETRY_FACTORY.createLineString();
Geometry geomCollection = GEOMETRY_FACTORY.createGeometryCollection(new Geometry[] {GEOMETRY_FACTORY.createGeometryCollection(new Geometry[] {multiPolygon, point3D, emptyLineString})});
Geometry geomCollection = GEOMETRY_FACTORY.createGeometryCollection(new Geometry[]{GEOMETRY_FACTORY.createGeometryCollection(new Geometry[]{multiPolygon, point3D, emptyLineString})});
Polygon expectedPolygon = GEOMETRY_FACTORY.createPolygon(coordArray(2, 3, 2, 4, 3, 4, 3, 3, 2, 3));
Polygon expectedPolygon3D = GEOMETRY_FACTORY.createPolygon(coordArray3d(2, 3, 6, 3, 3, 7, 3, 4, 7, 2, 3, 6));
Point expectedPoint3D = GEOMETRY_FACTORY.createPoint(new Coordinate(2, 4, 6));
Expand All @@ -934,6 +934,45 @@ public void translateHybridGeomCollectionDeltaZ() {
assertEquals(emptyLineString.toText(), actualGeometry.getGeometryN(0).getGeometryN(2).toText());
}

@Test
public void testFrechetGeom2D() {
LineString lineString1 = GEOMETRY_FACTORY.createLineString(coordArray(0, 0, 100, 0));
LineString lineString2 = GEOMETRY_FACTORY.createLineString(coordArray(0, 0, 50, 50, 100, 0));
double expected = 70.7106781186548;
double actual = Functions.frechetDistance(lineString1, lineString2);
assertEquals(expected, actual, 1e-9);
}

@Test
public void testFrechetGeom3D() {
LineString lineString = GEOMETRY_FACTORY.createLineString(coordArray3d(1, 0, 1, 2, 2, 2, 3, 3, 3));
Polygon polygon = GEOMETRY_FACTORY.createPolygon(coordArray3d(1, 0, 0, 1, 1, 1, 2, 1, 1, 1, 0, 0));
double expected = 3.605551275463989;
double actual = Functions.frechetDistance(lineString, polygon);
assertEquals(expected, actual, 1e-9);
}

@Test
public void testFrechetGeomCollection() {
Geometry point = GEOMETRY_FACTORY.createPoint(new Coordinate(1, 2));
Geometry lineString1 = GEOMETRY_FACTORY.createLineString(coordArray(2, 2, 3, 3, 4, 4));
Geometry lineString2 = GEOMETRY_FACTORY.createLineString(coordArray(-1, -1, -4, -4, -10, -10));
Geometry geometryCollection = GEOMETRY_FACTORY.createGeometryCollection(new Geometry[] {point, lineString1, lineString2});
Polygon polygon = GEOMETRY_FACTORY.createPolygon(coordArray(1, 0, 1, 1, 2, 1, 2, 0, 1, 0));
double expected = 14.866068747318506;
double actual = Functions.frechetDistance(polygon, geometryCollection);
assertEquals(expected, actual, 1e-9);
}

@Test
public void testFrechetGeomEmpty() {
Polygon p1 = GEOMETRY_FACTORY.createPolygon(coordArray(1, 0, 1, 1, 2, 1, 2, 0, 1, 0));
LineString emptyPoint = GEOMETRY_FACTORY.createLineString();
double expected = 0.0;
double actual = Functions.frechetDistance(p1, emptyPoint);
assertEquals(expected, actual, 1e-9);
}

@Test
public void boundingDiagonalGeom2D() {
Polygon polygon = GEOMETRY_FACTORY.createPolygon(coordArray(1, 0, 1, 1, 2, 1, 2, 2, 2, 0, 1, 0));
Expand Down
20 changes: 20 additions & 0 deletions docs/api/flink/Function.md
Original file line number Diff line number Diff line change
Expand Up @@ -544,6 +544,26 @@ Input: `LINESTRING EMPTY`

Output: `LINESTRING EMPTY`

## ST_FrechetDistance

Introduction: Computes and returns discrete [Frechet Distance](https://en.wikipedia.org/wiki/Fr%C3%A9chet_distance) between the given two geometrie,
based on [Computing Discrete Frechet Distance](http://www.kr.tuwien.ac.at/staff/eiter/et-archive/cdtr9464.pdf)

If any of the geometries is empty, returns 0.0

Format: `ST_FrechetDistance(g1: geomtry, g2: geometry)`

Since: `1.5.0`

Example:
```sql
SELECT ST_FrechetDistance(g1, g2)
```

Input: `g1: POINT (0 1), g2: LINESTRING (0 0, 1 0, 2 0, 3 0, 4 0, 5 0)`

Output: `5.0990195135927845`

## ST_GeoHash

Introduction: Returns GeoHash of the geometry with given precision
Expand Down
20 changes: 20 additions & 0 deletions docs/api/sql/Function.md
Original file line number Diff line number Diff line change
Expand Up @@ -730,6 +730,26 @@ Input: `LINESTRING EMPTY`

Output: `LINESTRING EMPTY`

## ST_FrechetDistance

Introduction: Computes and returns discrete [Frechet Distance](https://en.wikipedia.org/wiki/Fr%C3%A9chet_distance) between the given two geometrie,
based on [Computing Discrete Frechet Distance](http://www.kr.tuwien.ac.at/staff/eiter/et-archive/cdtr9464.pdf)

If any of the geometries is empty, returns 0.0

Format: `ST_FrechetDistance(g1: geomtry, g2: geometry)`

Since: `1.5.0`

Example:
```sql
SELECT ST_FrechetDistance(g1, g2)
```

Input: `g1: POINT (0 1), g2: LINESTRING (0 0, 1 0, 2 0, 3 0, 4 0, 5 0)`

Output: `5.0990195135927845`

## ST_GeoHash

Introduction: Returns GeoHash of the geometry with given precision
Expand Down
1 change: 1 addition & 0 deletions flink/src/main/java/org/apache/sedona/flink/Catalog.java
Original file line number Diff line number Diff line change
Expand Up @@ -100,6 +100,7 @@ public static UserDefinedFunction[] getFuncs() {
new Functions.ST_Force3D(),
new Functions.ST_NRings(),
new Functions.ST_Translate(),
new Functions.ST_FrechetDistance(),
new Functions.ST_Affine(),
new Functions.ST_BoundingDiagonal(),
new Functions.ST_HausdorffDistance(),
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -581,7 +581,16 @@ public Geometry eval(@DataTypeHint(value = "RAW", bridgedTo = org.locationtech.j
Geometry geometry = (Geometry) o;
return org.apache.sedona.common.Functions.geometricMedian(geometry, tolerance, maxIter, failIfNotConverged);
}
}

public static class ST_FrechetDistance extends ScalarFunction {
@DataTypeHint("Double")
public Double eval(@DataTypeHint(value = "RAW", bridgedTo = org.locationtech.jts.geom.Geometry.class) Object g1,
@DataTypeHint(value = "RAW", bridgedTo = org.locationtech.jts.geom.Geometry.class) Object g2) {
Geometry geom1 = (Geometry) g1;
Geometry geom2 = (Geometry) g2;
return org.apache.sedona.common.Functions.frechetDistance(geom1, geom2);
}
}

public static class ST_NumPoints extends ScalarFunction {
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -747,6 +747,14 @@ public void testTranslate() {
}

@Test
public void testFrechet() {
Table polyTable = tableEnv.sqlQuery("SELECT ST_GeomFromWKT('POINT (1 2)') AS g1, ST_GeomFromWKT('POINT (10 10)') as g2");
polyTable = polyTable.select(call(Functions.ST_FrechetDistance.class.getSimpleName(), $("g1"), $("g2")));
Double expected = 12.041594578792296;
Double actual = (Double) first(polyTable).getField(0);
assertEquals(expected, actual);
}

public void testBoundingDiagonal() {
Table polyTable = tableEnv.sqlQuery("SELECT ST_BoundingDiagonal(ST_GeomFromWKT('POLYGON ((1 0, 1 1, 2 1, 2 0, 1 0))'))" +" AS " + polygonColNames[0]);
polyTable = polyTable.select(call(Functions.ST_AsText.class.getSimpleName(), $(polygonColNames[0])));
Expand Down
13 changes: 13 additions & 0 deletions python/sedona/sql/st_functions.py
Original file line number Diff line number Diff line change
Expand Up @@ -113,6 +113,7 @@
"ST_Force3D",
"ST_NRings",
"ST_Translate",
"ST_FrechetDistance",
"ST_Affine",
"ST_BoundingDiagonal"
]
Expand Down Expand Up @@ -1290,6 +1291,18 @@ def ST_Translate(geometry: ColumnOrName, deltaX: Union[ColumnOrName, float], del
args = (geometry, deltaX, deltaY, deltaZ)
return _call_st_function("ST_Translate", args)

def ST_FrechetDistance(g1: ColumnOrName, g2: ColumnOrName) -> Column:
"""
Computes discrete frechet distance between the two geometries.
If any of the geometry is empty, ST_FrechetDistance returns 0
:param g1:
:param g2:
:return: Computed Discrete Frechet Distance between g1 and g2
"""

args = (g1, g2)
return _call_st_function("ST_FrechetDistance", args)

@validate_argument_types
def ST_Affine(geometry: ColumnOrName, a: Union[ColumnOrName, float], b: Union[ColumnOrName, float], d: Union[ColumnOrName, float],
e: Union[ColumnOrName, float], xOff: Union[ColumnOrName, float], yOff: Union[ColumnOrName, float], c: Optional[Union[ColumnOrName, float]] = None, f: Optional[Union[ColumnOrName, float]] = None,
Expand Down
1 change: 1 addition & 0 deletions python/tests/sql/test_dataframe_api.py
Original file line number Diff line number Diff line change
Expand Up @@ -89,6 +89,7 @@
(stf.ST_FlipCoordinates, ("point",), "point_geom", "", "POINT (1 0)"),
(stf.ST_Force_2D, ("point",), "point_geom", "", "POINT (0 1)"),
(stf.ST_Force3D, ("point", 1.0), "point_geom", "", "POINT Z (0 1 1)"),
(stf.ST_FrechetDistance, ("point", "line",), "point_and_line", "", 5.0990195135927845),
(stf.ST_GeometricMedian, ("multipoint",), "multipoint_geom", "", "POINT (22.500002656424286 21.250001168173426)"),
(stf.ST_GeometryN, ("geom", 0), "multipoint", "", "POINT (0 0)"),
(stf.ST_GeometryType, ("point",), "point_geom", "", "ST_Point"),
Expand Down
7 changes: 7 additions & 0 deletions python/tests/sql/test_function.py
Original file line number Diff line number Diff line change
Expand Up @@ -1129,6 +1129,13 @@ def test_translate(self):
actual = actual_df.selectExpr("ST_AsText(geom)").take(1)[0][0]
assert expected == actual

def test_frechetDistance(self):
expected = 5.0990195135927845
actual_df = self.spark.sql("SELECT ST_FrechetDistance(ST_GeomFromText('LINESTRING (0 0, 1 0, 2 0, 3 0, 4 0, "
"5 0)'), ST_GeomFromText('POINT (0 1)'))")
actual = actual_df.take(1)[0][0]
assert expected == actual

def test_affine(self):
expected = "POLYGON Z((2 3 1, 4 5 1, 7 8 2, 2 3 1))"
actual_df = self.spark.sql("SELECT ST_Affine(ST_GeomFromText('POLYGON ((1 0 1, 1 1 1, 2 2 2, 1 0 1))'), 1, 2, "
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -152,6 +152,7 @@ object Catalog {
function[ST_Force3D](0.0),
function[ST_NRings](),
function[ST_Translate](0.0),
function[ST_FrechetDistance](),
function[ST_Affine](null, null, null, null, null, null),
function[ST_BoundingDiagonal](),
function[ST_HausdorffDistance](-1),
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -1005,6 +1005,13 @@ case class ST_Translate(inputExpressions: Seq[Expression])
}
}

case class ST_FrechetDistance(inputExpressions: Seq[Expression])
extends InferredExpression(Functions.frechetDistance _) with FoldableExpression {
protected def withNewChildrenInternal(newChildren: IndexedSeq[Expression]) = {
copy(inputExpressions = newChildren)
}
}

case class ST_Affine(inputExpressions: Seq[Expression])
extends InferredExpression(InferrableFunction.allowSixRightNull(Functions.affine _)) with FoldableExpression {
protected def withNewChildrenInternal(newChildren: IndexedSeq[Expression]) = {
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -329,6 +329,9 @@ object st_functions extends DataFrameAPI {

def ST_Translate(geometry: String, deltaX: Double, deltaY: Double): Column = wrapExpression[ST_Translate](geometry, deltaX, deltaY, 0.0)

def ST_FrechetDistance(g1: Column, g2: Column): Column = wrapExpression[ST_FrechetDistance](g1, g2)

def ST_FrechetDistance(g1: String, g2: String): Column = wrapExpression[ST_FrechetDistance](g1, g2)
def ST_Affine(geometry: Column, a: Column, b: Column, d: Column, e: Column, xOff: Column, yOff: Column, c: Column, f: Column, g: Column, h: Column, i: Column, zOff: Column): Column =
wrapExpression[ST_Affine](geometry, a, b, d, e, xOff, yOff, c, f, g, h, i, zOff)

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -1003,6 +1003,14 @@ class dataFrameAPITestScala extends TestBaseScala {
assert(expectedDefaultValue == actualDefaultValue)
}

it("Passed ST_FrechetDistance") {
val polyDf = sparkSession.sql("SELECT ST_GeomFromWKT('POINT (1 2)') as g1, ST_GeomFromWKT('POINT (100 230)') as g2")
val df = polyDf.select(ST_FrechetDistance("g1", "g2"))
val expected = 248.5658866377283
val actual = df.take(1)(0).get(0).asInstanceOf[Double]
assertEquals(expected, actual, 1e-9)
}

it("Passed ST_Affine") {
val polyDf = sparkSession.sql("SELECT ST_GeomFromWKT('POLYGON ((2 3 1, 4 5 1, 7 8 2, 2 3 1))') AS geom")
val df = polyDf.select(ST_Affine("geom", 1, 2, 3, 4, 1, 2, 3, 4, 1, 4, 2, 1));
Expand All @@ -1029,7 +1037,7 @@ class dataFrameAPITestScala extends TestBaseScala {

it("Passed ST_HausdorffDistance") {
val polyDf = sparkSession.sql("SELECT ST_GeomFromWKT('POLYGON ((1 2, 2 1, 2 0, 4 1, 1 2))') AS g1, " +
"ST_GeomFromWKT('MULTILINESTRING ((1 1, 2 1, 4 4, 5 5), (10 10, 11 11, 12 12, 14 14), (-11 -20, -11 -21, -15 -19))') AS g2")
"ST_GeomFromWKT('MULTILINESTRING ((1 1, 2 1, 4 4, 5 5), (10 10, 11 11, 12 12, 14 14), (-11 -20, -11 -21, -15 -19))') AS g2")
val df = polyDf.select(ST_HausdorffDistance("g1", "g2", 0.05))
val dfDefaultValue = polyDf.select(ST_HausdorffDistance("g1", "g2"))
val expected = 25.495097567963924
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -2005,6 +2005,22 @@ class functionTestScala extends TestBaseScala with Matchers with GeometrySample
}
}

it ("should pass ST_FrechetDistance") {
val geomTestCases = Map(
("'POINT (1 2)'", "'POINT (10 10)'") -> 12.041594578792296d,
("'LINESTRING (0 0, 100 0)'", "'LINESTRING (0 0, 50 50, 100 0)'") -> 70.7106781186548d
)
for (((geom), expectedResult) <- geomTestCases) {
val g1 = geom._1
val g2 = geom._2
val df = sparkSession.sql(s"SELECT ST_FrechetDistance(ST_GeomFromWKT($g1), ST_GeomFromWKT($g2))")
val actual = df.take(1)(0).get(0).asInstanceOf[Double]
val expected = expectedResult
assertEquals(expected, actual, 1e-9)

}
}

it ("should pass ST_Affine") {
val geomTestCases = Map (
("'POLYGON ((1 0 1, 1 1 1, 2 2 2, 1 0 1))'")-> ("'POLYGON Z((5 8 16, 7 9 20, 13 16 37, 5 8 16))'", "'POLYGON Z((2 3 1, 4 5 1, 7 8 2, 2 3 1))'"),
Expand Down
Loading