diff --git a/jmeos-core/src/main/java/utils/spatial/Haversine.java b/jmeos-core/src/main/java/utils/spatial/Haversine.java new file mode 100644 index 00000000..2ade89d5 --- /dev/null +++ b/jmeos-core/src/main/java/utils/spatial/Haversine.java @@ -0,0 +1,54 @@ +package utils.spatial; + +import functions.functions; +import jnr.ffi.Pointer; + +/** + * Geodesic ("haversine") distance between two WGS84 (longitude, latitude) + * points, computed via MEOS {@code geog_distance} on ephemeral geography + * point objects. This wrapper is the JMEOS-backed replacement for + * pure-Java haversine formulas used by stream-side consumers (MobilityFlink + * / MobilityKafka per-event spatial-predicate call sites). + * + *
Semantics match MEOS-on-PostgreSQL {@code ST_Distance(geography, + * geography)}: WGS84 spheroidal (not spherical) geodesic distance in + * meters. This is the canonical reference shared with every other MEOS + * spatial operator in the pipeline; the wrapper exists so consumers do + * not maintain a parallel pure-Java haversine that semantically drifts. + * + *
Cost: two {@code geog_in} allocations per call, one + * {@code geog_distance} invocation. The JNR-FFI pointers go to the JVM + * for GC; the underlying GSERIALIZED bytes follow MEOS's existing + * memory-management contract (see {@code functions.functions.geog_in}'s + * pattern; no explicit free path is required by current JMEOS callers). + */ +public final class Haversine { + + private Haversine() { + // Utility class — not instantiable. + } + + /** + * Geodesic distance in meters between two WGS84 points. + * + *
The coordinates are interpreted as (longitude, latitude) in + * degrees, in SRID 4326 (WGS84). + * + * @param lon1 longitude of point 1, degrees in {@code (-180, 180]} + * @param lat1 latitude of point 1, degrees in {@code [-90, 90]} + * @param lon2 longitude of point 2, degrees in {@code (-180, 180]} + * @param lat2 latitude of point 2, degrees in {@code [-90, 90]} + * @return geodesic distance in meters, non-negative + */ + public static double distance(double lon1, double lat1, + double lon2, double lat2) { + // Build EWKT strings with Java's locale-independent Double.toString + // (avoids "%f" which depends on Locale.getDefault() for the + // decimal separator and would corrupt WKT in fr-FR / de-DE etc.). + String wkt1 = "SRID=4326;POINT(" + lon1 + " " + lat1 + ")"; + String wkt2 = "SRID=4326;POINT(" + lon2 + " " + lat2 + ")"; + Pointer g1 = functions.geog_in(wkt1, -1); + Pointer g2 = functions.geog_in(wkt2, -1); + return functions.geog_distance(g1, g2); + } +} diff --git a/jmeos-core/src/main/java/utils/spatial/PointToSegment.java b/jmeos-core/src/main/java/utils/spatial/PointToSegment.java new file mode 100644 index 00000000..7f99e128 --- /dev/null +++ b/jmeos-core/src/main/java/utils/spatial/PointToSegment.java @@ -0,0 +1,62 @@ +package utils.spatial; + +import functions.functions; +import jnr.ffi.Pointer; + +/** + * Minimum geodesic distance between a WGS84 point and a 2-vertex line + * segment (s1 → s2), computed via MEOS {@code geog_distance} on + * ephemeral geography point + linestring objects. This wrapper is the + * JMEOS-backed replacement for the planar equirectangular + * point-to-segment fallback used by stream-side consumers + * (MobilityFlink {@code SegmentDistance.java}, MobilityKafka). + * + *
Semantics match MEOS-on-PostgreSQL {@code ST_Distance(POINT, + * LINESTRING)} on the geography flavour: WGS84 spheroidal geodesic + * distance in meters, including the segment-perpendicular case + * (closest point falls between the endpoints). This is the canonical + * reference shared with every other MEOS spatial operator in the + * pipeline; the wrapper exists so consumers do not maintain a parallel + * planar fallback that semantically drifts. + * + *
Cost: two {@code geog_in} allocations per call (point + 2-vertex + * linestring), one {@code geog_distance} invocation. The JNR-FFI + * pointers go to the JVM for GC; the underlying GSERIALIZED bytes + * follow MEOS's existing memory-management contract. + */ +public final class PointToSegment { + + private PointToSegment() { + // Utility class — not instantiable. + } + + /** + * Minimum geodesic distance in meters from {@code p} to the line + * segment {@code s1 → s2}. + * + *
All coordinates are interpreted as (longitude, latitude) in + * degrees, SRID 4326 (WGS84). The segment is degenerate if + * {@code s1 == s2}; in that case this reduces to point-to-point + * geodesic distance. + * + * @param pLon longitude of the query point, degrees + * @param pLat latitude of the query point, degrees + * @param s1Lon longitude of segment endpoint 1, degrees + * @param s1Lat latitude of segment endpoint 1, degrees + * @param s2Lon longitude of segment endpoint 2, degrees + * @param s2Lat latitude of segment endpoint 2, degrees + * @return minimum geodesic distance in meters, non-negative + */ + public static double distance(double pLon, double pLat, + double s1Lon, double s1Lat, + double s2Lon, double s2Lat) { + String ptWkt = + "SRID=4326;POINT(" + pLon + " " + pLat + ")"; + String segWkt = + "SRID=4326;LINESTRING(" + s1Lon + " " + s1Lat + + "," + s2Lon + " " + s2Lat + ")"; + Pointer pt = functions.geog_in(ptWkt, -1); + Pointer seg = functions.geog_in(segWkt, -1); + return functions.geog_distance(pt, seg); + } +} diff --git a/jmeos-core/src/test/java/utils/spatial/HaversineTest.java b/jmeos-core/src/test/java/utils/spatial/HaversineTest.java new file mode 100644 index 00000000..1c3a9110 --- /dev/null +++ b/jmeos-core/src/test/java/utils/spatial/HaversineTest.java @@ -0,0 +1,74 @@ +package utils.spatial; + +import functions.error_handler; +import functions.error_handler_fn; +import functions.functions; +import org.junit.jupiter.api.BeforeAll; +import org.junit.jupiter.api.Test; +import org.junit.jupiter.api.extension.ExtendWith; +import utils.TestLogger; + +import static org.junit.jupiter.api.Assertions.assertEquals; +import static org.junit.jupiter.api.Assertions.assertTrue; + +/** + * Unit tests for {@link utils.spatial.Haversine}. Numeric expected + * values are MEOS-on-PostgreSQL ground truth (WGS84 spheroidal). + */ +@ExtendWith(TestLogger.class) +public class HaversineTest { + + private static final double METER_TOLERANCE = 1.0; + + private static final error_handler_fn errorHandler = new error_handler(); + + @BeforeAll + static void initMeos() { + functions.meos_initialize_timezone("UTC"); + functions.meos_initialize_error_handler(errorHandler); + } + + /** Same point twice = zero distance. */ + @Test + void zeroDistanceForIdenticalPoints() { + double d = Haversine.distance(4.35, 50.85, 4.35, 50.85); + assertEquals(0.0, d, 1e-9); + } + + /** ~5.5 km Brussels-area segment along a meridian. */ + @Test + void shortMeridianSegment() { + // (4.35, 50.85) to (4.35, 50.90) — 0.05° latitude ≈ 5.56 km on WGS84 + double d = Haversine.distance(4.35, 50.85, 4.35, 50.90); + assertEquals(5562.0, d, METER_TOLERANCE); + } + + /** ~264 km Brussels-to-Paris diagonal. */ + @Test + void brusselsToParis() { + double d = Haversine.distance(4.35, 50.85, 2.35, 48.86); + // MEOS-on-PostgreSQL ground truth (WGS84 spheroidal): + // SELECT ST_Distance( + // 'SRID=4326;POINT(4.35 50.85)'::geography, + // 'SRID=4326;POINT(2.35 48.86)'::geography + // ) -> ~263538 m + assertEquals(263538.0, d, 500.0); + } + + /** Symmetric: d(a, b) == d(b, a). */ + @Test + void symmetric() { + double d1 = Haversine.distance(4.35, 50.85, 2.35, 48.86); + double d2 = Haversine.distance(2.35, 48.86, 4.35, 50.85); + assertEquals(d1, d2, 1e-9); + } + + /** Non-negative for any input. */ + @Test + void nonNegative() { + double d = Haversine.distance(-122.4, 37.8, 139.7, 35.7); // SF -> Tokyo + assertTrue(d > 0.0); + // Sanity: SF-Tokyo great-circle distance is ~8270 km + assertTrue(d > 8.0e6 && d < 9.0e6, "expected ~8.27e6 m, got " + d); + } +} diff --git a/jmeos-core/src/test/java/utils/spatial/PointToSegmentTest.java b/jmeos-core/src/test/java/utils/spatial/PointToSegmentTest.java new file mode 100644 index 00000000..c3454306 --- /dev/null +++ b/jmeos-core/src/test/java/utils/spatial/PointToSegmentTest.java @@ -0,0 +1,88 @@ +package utils.spatial; + +import functions.error_handler; +import functions.error_handler_fn; +import functions.functions; +import org.junit.jupiter.api.BeforeAll; +import org.junit.jupiter.api.Test; +import org.junit.jupiter.api.extension.ExtendWith; +import utils.TestLogger; + +import static org.junit.jupiter.api.Assertions.assertEquals; +import static org.junit.jupiter.api.Assertions.assertTrue; + +/** + * Unit tests for {@link utils.spatial.PointToSegment}. Numeric expected + * values are MEOS-on-PostgreSQL ground truth (WGS84 spheroidal). + */ +@ExtendWith(TestLogger.class) +public class PointToSegmentTest { + + private static final double METER_TOLERANCE = 1.0; + + private static final error_handler_fn errorHandler = new error_handler(); + + @BeforeAll + static void initMeos() { + functions.meos_initialize_timezone("UTC"); + functions.meos_initialize_error_handler(errorHandler); + } + + /** Point on the segment endpoint = zero distance. */ + @Test + void zeroAtEndpoint() { + double d = PointToSegment.distance( + 4.35, 50.85, + 4.35, 50.85, 4.40, 50.90); + assertEquals(0.0, d, METER_TOLERANCE); + } + + /** Point on the segment interior (midpoint) = zero distance. */ + @Test + void zeroAtInterior() { + // Midpoint of the 5.5 km meridian segment is on the segment exactly. + double d = PointToSegment.distance( + 4.35, 50.875, + 4.35, 50.85, 4.35, 50.90); + assertEquals(0.0, d, METER_TOLERANCE); + } + + /** + * Point off to the side of the segment: perpendicular distance. + * + *
A 5.5 km meridian segment along longitude 4.35 between latitudes + * 50.85 and 50.90; a query point at (4.40, 50.875) — same latitude + * as the midpoint but 0.05° east, which is ~3.5 km at this latitude. + */ + @Test + void perpendicularDistance() { + double d = PointToSegment.distance( + 4.40, 50.875, + 4.35, 50.85, 4.35, 50.90); + // MEOS-on-PostgreSQL ground truth: WGS84 spheroidal distance from + // (4.40, 50.875) to the meridian segment ≈ 3520 m + assertTrue(d > 3400.0 && d < 3700.0, "expected ~3520 m, got " + d); + } + + /** Beyond the segment: distance to the nearest endpoint. */ + @Test + void beyondEndpointFallsBackToEndpoint() { + // Query point far north of the segment (north of latitude 50.90) + // → nearest is endpoint (4.35, 50.90). + double dToSegment = PointToSegment.distance( + 4.35, 51.00, + 4.35, 50.85, 4.35, 50.90); + double dToEndpoint = Haversine.distance(4.35, 51.00, 4.35, 50.90); + assertEquals(dToEndpoint, dToSegment, METER_TOLERANCE); + } + + /** Degenerate segment (s1 == s2) reduces to point-to-point distance. */ + @Test + void degenerateSegmentReducesToHaversine() { + double dSeg = PointToSegment.distance( + 4.35, 50.85, + 4.40, 50.90, 4.40, 50.90); + double dPt = Haversine.distance(4.35, 50.85, 4.40, 50.90); + assertEquals(dPt, dSeg, METER_TOLERANCE); + } +}