From 3e58d9ebe5b89f3ed635bf2dd61df978d42acae9 Mon Sep 17 00:00:00 2001 From: t-ober <63147366+t-ober@users.noreply.github.com> Date: Wed, 1 Dec 2021 14:49:00 +0100 Subject: [PATCH 01/17] add GeoUtils --- docs/OsmModel.puml | 2 +- .../scala/edu/ie3/util/geo/GeoUtils.scala | 83 +++++++++++++++++++ .../scala/edu/ie3/util/osm/OsmEntities.scala | 24 +++++- .../scala/edu/ie3/util/osm/OsmModel.scala | 5 +- .../scala/edu/ie3/util/osm/TagUtils.scala | 1 + 5 files changed, 109 insertions(+), 6 deletions(-) create mode 100644 src/main/scala/edu/ie3/util/geo/GeoUtils.scala diff --git a/docs/OsmModel.puml b/docs/OsmModel.puml index 6a5973a6..c1d01cdd 100644 --- a/docs/OsmModel.puml +++ b/docs/OsmModel.puml @@ -28,7 +28,7 @@ Class Node { } abstract Class Way { - nodes: List[UUID] + nodes: List[Node] } Class OpenWay { diff --git a/src/main/scala/edu/ie3/util/geo/GeoUtils.scala b/src/main/scala/edu/ie3/util/geo/GeoUtils.scala new file mode 100644 index 00000000..7d35ac5f --- /dev/null +++ b/src/main/scala/edu/ie3/util/geo/GeoUtils.scala @@ -0,0 +1,83 @@ +/* + * © 2021. TU Dortmund University, + * Institute of Energy Systems, Energy Efficiency and Energy Economics, + * Research group Distribution grid planning and operation +*/ +package edu.ie3.util.geo + +import edu.ie3.util.osm.OsmEntities.{ClosedWay, Way} +import edu.ie3.util.quantities.PowerSystemUnits.KILOMETRE +import org.locationtech.jts.geom.{GeometryFactory, Point, PrecisionModel} +import tech.units.indriya.ComparableQuantity +import tech.units.indriya.quantity.Quantities +import tech.units.indriya.unit.Units.{METRE, RADIAN} + +import java.lang.Math.{atan2, cos, sin, sqrt, toRadians} +import javax.measure.quantity.Length +import scala.math.pow + +object GeoUtils { + val DEFAULT_GEOMETRY_FACTORY: GeometryFactory = + new GeometryFactory(new PrecisionModel(), 4326) + + val EARTH_RADIUS: ComparableQuantity[Length] = + Quantities.getQuantity(6378137.0, METRE); + + /** Calculates the great circle disteance between two coordinates + * + * @param coordinateA + * coordinate a + * @param coordinateB + * coordinate b + * @return + * the distance between the coordinates as a quantity + */ + def calcHaversine( + coordinateA: Point, + coordinateB: Point + ): ComparableQuantity[Length] = { + calcHaversine( + coordinateA.getY, + coordinateA.getX, + coordinateB.getY, + coordinateB.getX + ) + } + + /** Calculates the great circle disteance between two coordinates + * + * @param latA + * latitude of coordinate a + * @param longA + * longitude of coordinate a + * @param latB + * latitude of coordinate b + * @param longB + * longitude of coordinate b + * @return + * the distance between the coordinates as a quantity + */ + def calcHaversine( + latA: Double, + longA: Double, + latB: Double, + longB: Double + ): ComparableQuantity[Length] = { + val r = EARTH_RADIUS.to(KILOMETRE) + val dLat = Quantities.getQuantity(toRadians(latB - latA), RADIAN) + val dLon = Quantities.getQuantity(toRadians(longB - longA), RADIAN) + val a = pow(sin(dLat.getValue.doubleValue / 2), 2) + cos( + toRadians(latA) + ) * cos(toRadians(latB)) * pow(sin(dLon.getValue.doubleValue / 2), 2) + val c = 2 * atan2(sqrt(a), sqrt(1 - a)) + r.multiply(c) + } + + def isInsideLandUse(coordinate: Point, landuses: List[ClosedWay]): Boolean = { + landuses.foreach(landuse => + if (landuse.toPolygon.convexHull().contains(coordinate)) return true + ) + false + } + +} diff --git a/src/main/scala/edu/ie3/util/osm/OsmEntities.scala b/src/main/scala/edu/ie3/util/osm/OsmEntities.scala index 704eb670..89e61fa2 100644 --- a/src/main/scala/edu/ie3/util/osm/OsmEntities.scala +++ b/src/main/scala/edu/ie3/util/osm/OsmEntities.scala @@ -5,13 +5,17 @@ */ package edu.ie3.util.osm -import org.locationtech.jts.geom.Point +import edu.ie3.util.geo.GeoUtils.DEFAULT_GEOMETRY_FACTORY +import org.locationtech.jts.geom.impl.CoordinateArraySequence +import org.locationtech.jts.geom.{Coordinate, GeometryFactory, LinearRing, Polygon, PrecisionModel} import java.time.ZonedDateTime import java.util.UUID +import scala.jdk.CollectionConverters._ object OsmEntities { + /** Common trait to all OpenStreetMap entities */ sealed trait OsmEntity { @@ -63,7 +67,7 @@ object OsmEntities { * time stamp of the elements last edit * @param tags * associated tags - * @param coordinates + * @param coordinate * the coordinates of the specific node */ final case class Node( @@ -71,7 +75,7 @@ object OsmEntities { osmId: Int, lastEdited: ZonedDateTime, tags: Map[String, String], - coordinates: Point + coordinate: Coordinate ) extends OsmEntity /** Common trait to all OSM ways @@ -121,7 +125,19 @@ object OsmEntities { lastEdited: ZonedDateTime, tags: Map[String, String], nodes: List[Node] - ) extends Way + ) extends Way { + + def getCoordinates: List[Coordinate] = { + nodes.map(_.coordinate) + } + + def toPolygon: Polygon = { + val arrayCoordinates = new CoordinateArraySequence(getCoordinates.toArray) + val linearRing = new LinearRing(arrayCoordinates, DEFAULT_GEOMETRY_FACTORY) + new Polygon(linearRing, Array[LinearRing](), DEFAULT_GEOMETRY_FACTORY) + } + + } /** An OSM relation. * diff --git a/src/main/scala/edu/ie3/util/osm/OsmModel.scala b/src/main/scala/edu/ie3/util/osm/OsmModel.scala index 8d924b2a..9e03a20f 100644 --- a/src/main/scala/edu/ie3/util/osm/OsmModel.scala +++ b/src/main/scala/edu/ie3/util/osm/OsmModel.scala @@ -48,7 +48,10 @@ object OsmModel { * @return * all ways that represent highways */ - def extractHighways(ways: List[Way]): List[Way] = { + def extractHighways( + ways: List[Way], + highWayValues: Set[String] + ): List[Way] = { ways.filter(way => way.containsKeyValuePair(highway, highWayValues)) } diff --git a/src/main/scala/edu/ie3/util/osm/TagUtils.scala b/src/main/scala/edu/ie3/util/osm/TagUtils.scala index 8dac3e07..18a215b6 100644 --- a/src/main/scala/edu/ie3/util/osm/TagUtils.scala +++ b/src/main/scala/edu/ie3/util/osm/TagUtils.scala @@ -13,6 +13,7 @@ object TagUtils { val landuse = "landuse" } + // Fixme: Rausschmeißen object Values { val highWayValues = Set( "residential", From cf29710df3663bf721768b4e067fd8bf350efac3 Mon Sep 17 00:00:00 2001 From: t-ober <63147366+t-ober@users.noreply.github.com> Date: Wed, 1 Dec 2021 14:51:21 +0100 Subject: [PATCH 02/17] rename deprecated GeoUtils --- .../edu/ie3/util/geo/CoordinateDistance.java | 4 +- ...{GeoUtils.java => DeprecatedGeoUtils.java} | 8 ++-- .../scala/edu/ie3/util/osm/OsmEntities.scala | 12 ++++-- .../util/geo/CoordinateDistanceTest.groovy | 16 ++++---- .../edu/ie3/util/geo/GeoUtilsTest.groovy | 38 +++++++++---------- 5 files changed, 42 insertions(+), 36 deletions(-) rename src/main/java/edu/ie3/util/geo/{GeoUtils.java => DeprecatedGeoUtils.java} (99%) diff --git a/src/main/java/edu/ie3/util/geo/CoordinateDistance.java b/src/main/java/edu/ie3/util/geo/CoordinateDistance.java index 142d8c1d..2be7b554 100644 --- a/src/main/java/edu/ie3/util/geo/CoordinateDistance.java +++ b/src/main/java/edu/ie3/util/geo/CoordinateDistance.java @@ -20,7 +20,7 @@ public class CoordinateDistance implements Comparable { /** * Calculates the distance from the first to the second coordinate using {@link - * GeoUtils#calcHaversine(double, double, double, double)} + * DeprecatedGeoUtils#calcHaversine(double, double, double, double)} * * @param coordinateA The first coordinate * @param coordinateB The second coordinate @@ -30,7 +30,7 @@ public CoordinateDistance( this( coordinateA, coordinateB, - GeoUtils.calcHaversine( + DeprecatedGeoUtils.calcHaversine( coordinateA.getY(), coordinateA.getX(), coordinateB.getY(), coordinateB.getX())); } diff --git a/src/main/java/edu/ie3/util/geo/GeoUtils.java b/src/main/java/edu/ie3/util/geo/DeprecatedGeoUtils.java similarity index 99% rename from src/main/java/edu/ie3/util/geo/GeoUtils.java rename to src/main/java/edu/ie3/util/geo/DeprecatedGeoUtils.java index 39857f30..fcaf11c9 100644 --- a/src/main/java/edu/ie3/util/geo/GeoUtils.java +++ b/src/main/java/edu/ie3/util/geo/DeprecatedGeoUtils.java @@ -33,8 +33,8 @@ import tech.units.indriya.quantity.Quantities; /** Functionality to deal with geographical and geometric information */ -public class GeoUtils { - private static final Logger logger = LoggerFactory.getLogger(GeoUtils.class); +public class DeprecatedGeoUtils { + private static final Logger logger = LoggerFactory.getLogger(DeprecatedGeoUtils.class); public static final ComparableQuantity EARTH_RADIUS = Quantities.getQuantity(6378137.0, METRE); @@ -43,7 +43,7 @@ public class GeoUtils { public static final GeometryFactory DEFAULT_GEOMETRY_FACTORY = new GeometryFactory(new PrecisionModel(), 4326); - protected GeoUtils() { + protected DeprecatedGeoUtils() { throw new IllegalStateException("Utility classes cannot be instantiated"); } @@ -607,7 +607,7 @@ > calcHaversine( /** * Calculates the area, which is surrounded by a closed way by the help of {@link - * GeoUtils#calcArea(Polygon)} + * DeprecatedGeoUtils#calcArea(Polygon)} * * @param w Closed way, that surrounds the area * @return The covered area in {@link edu.ie3.util.quantities.PowerSystemUnits#SQUARE_METRE} diff --git a/src/main/scala/edu/ie3/util/osm/OsmEntities.scala b/src/main/scala/edu/ie3/util/osm/OsmEntities.scala index 89e61fa2..6b6f5510 100644 --- a/src/main/scala/edu/ie3/util/osm/OsmEntities.scala +++ b/src/main/scala/edu/ie3/util/osm/OsmEntities.scala @@ -7,7 +7,13 @@ package edu.ie3.util.osm import edu.ie3.util.geo.GeoUtils.DEFAULT_GEOMETRY_FACTORY import org.locationtech.jts.geom.impl.CoordinateArraySequence -import org.locationtech.jts.geom.{Coordinate, GeometryFactory, LinearRing, Polygon, PrecisionModel} +import org.locationtech.jts.geom.{ + Coordinate, + GeometryFactory, + LinearRing, + Polygon, + PrecisionModel +} import java.time.ZonedDateTime import java.util.UUID @@ -15,7 +21,6 @@ import scala.jdk.CollectionConverters._ object OsmEntities { - /** Common trait to all OpenStreetMap entities */ sealed trait OsmEntity { @@ -133,7 +138,8 @@ object OsmEntities { def toPolygon: Polygon = { val arrayCoordinates = new CoordinateArraySequence(getCoordinates.toArray) - val linearRing = new LinearRing(arrayCoordinates, DEFAULT_GEOMETRY_FACTORY) + val linearRing = + new LinearRing(arrayCoordinates, DEFAULT_GEOMETRY_FACTORY) new Polygon(linearRing, Array[LinearRing](), DEFAULT_GEOMETRY_FACTORY) } diff --git a/src/test/groovy/edu/ie3/util/geo/CoordinateDistanceTest.groovy b/src/test/groovy/edu/ie3/util/geo/CoordinateDistanceTest.groovy index 6f57f31b..e1e8a009 100644 --- a/src/test/groovy/edu/ie3/util/geo/CoordinateDistanceTest.groovy +++ b/src/test/groovy/edu/ie3/util/geo/CoordinateDistanceTest.groovy @@ -10,9 +10,9 @@ import spock.lang.Specification class CoordinateDistanceTest extends Specification { def "The constructor without a distance parameter calculates the distance as expected"() { given: - def pointA = GeoUtils.xyToPoint(49d, 7d) - def pointB = GeoUtils.xyToPoint(50d, 7d) - def expectedDistance = GeoUtils.calcHaversine(pointA.y, pointA.x, pointB.y, pointB.x) + def pointA = DeprecatedGeoUtils.xyToPoint(49d, 7d) + def pointB = DeprecatedGeoUtils.xyToPoint(50d, 7d) + def expectedDistance = DeprecatedGeoUtils.calcHaversine(pointA.y, pointA.x, pointB.y, pointB.x) when: def coordinateDistance = new CoordinateDistance(pointA, pointB) @@ -26,11 +26,11 @@ class CoordinateDistanceTest extends Specification { def "CoordinateDistances are sortable using their distance field"() { given: - def basePoint = GeoUtils.xyToPoint(49d, 7d) - def distA = new CoordinateDistance(basePoint, GeoUtils.xyToPoint(50d, 7d)) - def distB = new CoordinateDistance(basePoint, GeoUtils.xyToPoint(50d, 7.1d)) - def distC = new CoordinateDistance(basePoint, GeoUtils.xyToPoint(49d, 7.1d)) - def distD = new CoordinateDistance(basePoint, GeoUtils.xyToPoint(52d, 9d)) + def basePoint = DeprecatedGeoUtils.xyToPoint(49d, 7d) + def distA = new CoordinateDistance(basePoint, DeprecatedGeoUtils.xyToPoint(50d, 7d)) + def distB = new CoordinateDistance(basePoint, DeprecatedGeoUtils.xyToPoint(50d, 7.1d)) + def distC = new CoordinateDistance(basePoint, DeprecatedGeoUtils.xyToPoint(49d, 7.1d)) + def distD = new CoordinateDistance(basePoint, DeprecatedGeoUtils.xyToPoint(52d, 9d)) def coordinateDistances = [distA, distB, distC, distD] when: diff --git a/src/test/groovy/edu/ie3/util/geo/GeoUtilsTest.groovy b/src/test/groovy/edu/ie3/util/geo/GeoUtilsTest.groovy index e1d312cc..5401cb2c 100644 --- a/src/test/groovy/edu/ie3/util/geo/GeoUtilsTest.groovy +++ b/src/test/groovy/edu/ie3/util/geo/GeoUtilsTest.groovy @@ -42,7 +42,7 @@ class GeoUtilsTest extends Specification { ComparableQuantity expected = Quantities.getQuantity(450.18011568984845, METRE) when: - ComparableQuantity actual = GeoUtils.calcHaversine(start.lat, start.lon, end.lat, end.lon) + ComparableQuantity actual = DeprecatedGeoUtils.calcHaversine(start.lat, start.lon, end.lat, end.lon) then: Math.abs(actual.subtract(expected).to(METRE).value.doubleValue()) < tolerance.value.doubleValue() @@ -50,14 +50,14 @@ class GeoUtilsTest extends Specification { def "Total length of LineString is correctly calculated"() { given: - def lineString = GeoUtils.DEFAULT_GEOMETRY_FACTORY.createLineString([ + def lineString = DeprecatedGeoUtils.DEFAULT_GEOMETRY_FACTORY.createLineString([ new Coordinate(22.69962d, 11.13038d, 0), new Coordinate(20.84247d, 28.14743d, 0), new Coordinate(24.21942d, 12.04265d, 0) ] as Coordinate[]) when: - ComparableQuantity y = GeoUtils.totalLengthOfLineString(lineString) + ComparableQuantity y = DeprecatedGeoUtils.totalLengthOfLineString(lineString) then: QuantityUtil.isEquivalentAbs(y, Quantities.getQuantity(3463.37, PowerSystemUnits.KILOMETRE), 10) @@ -66,12 +66,12 @@ class GeoUtilsTest extends Specification { def "The GridAndGeoUtils should get the CoordinateDistances between a base point and a collection of other points correctly"() { given: - def basePoint = GeoUtils.xyToPoint(49d, 7d) + def basePoint = DeprecatedGeoUtils.xyToPoint(49d, 7d) def points = [ - GeoUtils.xyToPoint(50d, 7d), - GeoUtils.xyToPoint(50d, 7.1d), - GeoUtils.xyToPoint(49d, 7.1d), - GeoUtils.xyToPoint(52d, 9d) + DeprecatedGeoUtils.xyToPoint(50d, 7d), + DeprecatedGeoUtils.xyToPoint(50d, 7.1d), + DeprecatedGeoUtils.xyToPoint(49d, 7.1d), + DeprecatedGeoUtils.xyToPoint(52d, 9d) ] def coordinateDistances = [ new CoordinateDistance(basePoint, points[0]), @@ -80,7 +80,7 @@ class GeoUtilsTest extends Specification { new CoordinateDistance(basePoint, points[3]) ] expect: - GeoUtils.getCoordinateDistances(basePoint, points) == new TreeSet(coordinateDistances) + DeprecatedGeoUtils.getCoordinateDistances(basePoint, points) == new TreeSet(coordinateDistances) } def "The GridAndGeoUtils should build a line string with two exact equal geo coordinates correctly avoiding the known bug in jts geometry"() { @@ -88,7 +88,7 @@ class GeoUtilsTest extends Specification { def line = geoJsonReader.read(lineString) as LineString when: - def safeLineString = GeoUtils.buildSafeLineString(line) + def safeLineString = DeprecatedGeoUtils.buildSafeLineString(line) def actualCoordinates = safeLineString.coordinates then: @@ -125,10 +125,10 @@ class GeoUtilsTest extends Specification { new Coordinate(51.49391, 7.41172), new Coordinate(51.49404, 7.41279) ] as Coordinate[] - def lineString = GeoUtils.DEFAULT_GEOMETRY_FACTORY.createLineString(coordinates) + def lineString = DeprecatedGeoUtils.DEFAULT_GEOMETRY_FACTORY.createLineString(coordinates) when: - def actual = GeoUtils.buildSafeLineString(lineString).coordinates + def actual = DeprecatedGeoUtils.buildSafeLineString(lineString).coordinates then: coordinates.length == actual.length @@ -142,23 +142,23 @@ class GeoUtilsTest extends Specification { def coord = new Coordinate(1, 1, 0) expect: - GeoUtils.buildSafeCoord(coord) == new Coordinate(1.0000000000001, 1.0000000000001, 1.0E-13) + DeprecatedGeoUtils.buildSafeCoord(coord) == new Coordinate(1.0000000000001, 1.0000000000001, 1.0E-13) } def "The GridAndGeoUtils should build a safe instance of a LineString between two provided coordinates correctly"() { expect: - GeoUtils.buildSafeLineStringBetweenCoords(coordA, coordB) == resLineString + DeprecatedGeoUtils.buildSafeLineStringBetweenCoords(coordA, coordB) == resLineString // do not change or remove the following line, it is NOT equal to the line above in this case! - GeoUtils.buildSafeLineStringBetweenCoords(coordA, coordB).equals(resLineString) + DeprecatedGeoUtils.buildSafeLineStringBetweenCoords(coordA, coordB).equals(resLineString) where: coordA | coordB || resLineString - new Coordinate(1, 1, 0) | new Coordinate(1, 1, 0) || GeoUtils.DEFAULT_GEOMETRY_FACTORY.createLineString([ + new Coordinate(1, 1, 0) | new Coordinate(1, 1, 0) || DeprecatedGeoUtils.DEFAULT_GEOMETRY_FACTORY.createLineString([ new Coordinate(1.0000000000001, 1.0000000000001, 1.0E-13), new Coordinate(1, 1, 0) ] as Coordinate[]) - new Coordinate(1, 1, 0) | new Coordinate(2, 2, 0) || GeoUtils.DEFAULT_GEOMETRY_FACTORY.createLineString([ + new Coordinate(1, 1, 0) | new Coordinate(2, 2, 0) || DeprecatedGeoUtils.DEFAULT_GEOMETRY_FACTORY.createLineString([ new Coordinate(1, 1, 0), new Coordinate(2, 2, 0) ] as Coordinate[]) @@ -171,7 +171,7 @@ class GeoUtilsTest extends Specification { Quantity radius = Quantities.getQuantity(50, METRE) when: - Polygon poly = GeoUtils.radiusWithCircleAsPolygon(center, radius) + Polygon poly = DeprecatedGeoUtils.radiusWithCircleAsPolygon(center, radius) List circlePoints = poly.getCoords() then: @@ -183,7 +183,7 @@ class GeoUtilsTest extends Specification { circlePoints.size() == 361 // rounded distance should be 50 meters circlePoints.forEach({ point -> - Double distance = GeoUtils.calcHaversine(center.lat, center.lon, point.lat, point.lon).to(METRE).value.doubleValue() + Double distance = DeprecatedGeoUtils.calcHaversine(center.lat, center.lon, point.lat, point.lon).to(METRE).value.doubleValue() Math.round(distance) == 50 }) From 01b65f0c18bc5107987ed3dd315bbeae772adc92 Mon Sep 17 00:00:00 2001 From: t-ober <63147366+t-ober@users.noreply.github.com> Date: Thu, 2 Dec 2021 15:08:01 +0100 Subject: [PATCH 03/17] add various GeoUtils and introduce RichGeometries --- .../ie3/util/exceptions/GeoException.scala | 11 ++ .../ie3/util/exceptions/OsmException.scala | 13 ++ .../scala/edu/ie3/util/geo/GeoUtils.scala | 138 ++++++++++++++---- .../scala/edu/ie3/util/osm/OsmEntities.scala | 26 ++-- .../scala/edu/ie3/util/osm/OsmUtils.scala | 20 +++ .../edu/ie3/util/osm/RichGeometries.scala | 135 +++++++++++++++++ .../ie3/util/quantities/QuantityUtils.scala | 14 +- 7 files changed, 316 insertions(+), 41 deletions(-) create mode 100644 src/main/scala/edu/ie3/util/exceptions/GeoException.scala create mode 100644 src/main/scala/edu/ie3/util/exceptions/OsmException.scala create mode 100644 src/main/scala/edu/ie3/util/osm/OsmUtils.scala create mode 100644 src/main/scala/edu/ie3/util/osm/RichGeometries.scala diff --git a/src/main/scala/edu/ie3/util/exceptions/GeoException.scala b/src/main/scala/edu/ie3/util/exceptions/GeoException.scala new file mode 100644 index 00000000..8ec57bf6 --- /dev/null +++ b/src/main/scala/edu/ie3/util/exceptions/GeoException.scala @@ -0,0 +1,11 @@ +/* + * © 2021. TU Dortmund University, + * Institute of Energy Systems, Energy Efficiency and Energy Economics, + * Research group Distribution grid planning and operation +*/ +package edu.ie3.util.exceptions + +class GeoException( + private val msg: String, + private val cause: Throwable = None.orNull +) extends Exception(msg, cause) diff --git a/src/main/scala/edu/ie3/util/exceptions/OsmException.scala b/src/main/scala/edu/ie3/util/exceptions/OsmException.scala new file mode 100644 index 00000000..67bf3cc5 --- /dev/null +++ b/src/main/scala/edu/ie3/util/exceptions/OsmException.scala @@ -0,0 +1,13 @@ +/* + * © 2021. TU Dortmund University, + * Institute of Energy Systems, Energy Efficiency and Energy Economics, + * Research group Distribution grid planning and operation +*/ +package edu.ie3.util.exceptions + +/** Base class for grouping OSM related exceptions + */ +class OsmException( + private val msg: String, + private val cause: Throwable = None.orNull +) extends Exception(msg, cause) diff --git a/src/main/scala/edu/ie3/util/geo/GeoUtils.scala b/src/main/scala/edu/ie3/util/geo/GeoUtils.scala index 7d35ac5f..a03b7bd4 100644 --- a/src/main/scala/edu/ie3/util/geo/GeoUtils.scala +++ b/src/main/scala/edu/ie3/util/geo/GeoUtils.scala @@ -5,16 +5,29 @@ */ package edu.ie3.util.geo -import edu.ie3.util.osm.OsmEntities.{ClosedWay, Way} +import edu.ie3.util.exceptions.GeoException import edu.ie3.util.quantities.PowerSystemUnits.KILOMETRE -import org.locationtech.jts.geom.{GeometryFactory, Point, PrecisionModel} +import org.locationtech.jts.algorithm.ConvexHull +import org.locationtech.jts.geom.impl.CoordinateArraySequence +import org.locationtech.jts.geom.{ + Coordinate, + GeometryCollection, + GeometryFactory, + LineString, + LinearRing, + Point, + Polygon, + PrecisionModel +} import tech.units.indriya.ComparableQuantity import tech.units.indriya.quantity.Quantities import tech.units.indriya.unit.Units.{METRE, RADIAN} import java.lang.Math.{atan2, cos, sin, sqrt, toRadians} +import javax.measure.Quantity import javax.measure.quantity.Length import scala.math.pow +import scala.util.{Failure, Success, Try} object GeoUtils { val DEFAULT_GEOMETRY_FACTORY: GeometryFactory = @@ -23,27 +36,6 @@ object GeoUtils { val EARTH_RADIUS: ComparableQuantity[Length] = Quantities.getQuantity(6378137.0, METRE); - /** Calculates the great circle disteance between two coordinates - * - * @param coordinateA - * coordinate a - * @param coordinateB - * coordinate b - * @return - * the distance between the coordinates as a quantity - */ - def calcHaversine( - coordinateA: Point, - coordinateB: Point - ): ComparableQuantity[Length] = { - calcHaversine( - coordinateA.getY, - coordinateA.getX, - coordinateB.getY, - coordinateB.getX - ) - } - /** Calculates the great circle disteance between two coordinates * * @param latA @@ -73,11 +65,101 @@ object GeoUtils { r.multiply(c) } - def isInsideLandUse(coordinate: Point, landuses: List[ClosedWay]): Boolean = { - landuses.foreach(landuse => - if (landuse.toPolygon.convexHull().contains(coordinate)) return true - ) - false + /** Builds a convex hull from a set of coordinates. + * + * @param coordinates + * the coordinates to consider + * @return + * a Try of the resulting polygon + */ + def buildConvexHull(coordinates: Set[Coordinate]): Try[Polygon] = { + new ConvexHull( + coordinates.toArray, + DEFAULT_GEOMETRY_FACTORY + ).getConvexHull match { + case polygon: Polygon => Success(polygon) + case _: LineString => + Failure( + new GeoException( + s"Got a line string as a convex hull. Probably coordinates: $coordinates only contains two different coordinates." + ) + ) + case _: Point => + Failure( + new GeoException( + s"Got a point as a convex hull. Probably coordinates: $coordinates only contains one coordinate." + ) + ) + case _: GeometryCollection => + Failure( + new GeoException( + s"Got a GeometryCollection. Probably $coordinates was empty." + ) + ) + case _ => Failure(new GeoException(s"Got an unexpected return type.")) + } + } + + /** Build a coordinate from latitude and longitude values. + * + * @param lat + * latitude value + * @param long + * longitude value + * @return + * the built [[Coordinate]] + */ + def buildCoordinate(lat: Double, long: Double): Coordinate = { + new Coordinate(long, lat) + } + + /** Builds a polygon from a List of coordinates + * + * @param coordinates + * the coordinates for building the polygon + * @return + * a [[Polygon]] + */ + def buildPolygon(coordinates: List[Coordinate]): Polygon = { + val arrayCoordinates = new CoordinateArraySequence(coordinates.toArray) + val linearRing = + new LinearRing(arrayCoordinates, DEFAULT_GEOMETRY_FACTORY) + new Polygon(linearRing, Array[LinearRing](), DEFAULT_GEOMETRY_FACTORY) + } + + /** Credits to Joe Kington + * (https://stackoverflow.com/questions/4681737/how-to-calculate-the-area-of-a-polygon-on-the-earths-surface-using-python#:~:text=Basically%2C%20you%20just%20multiply%20the,the%20cosine%20of%20the%20latitude.) + * Sinusoidal equal-area projection of latitude and longitude values to a 2d + * surface. This is an approximation gets worse when lines become very long + * + * @param lat + * latitude to project + * @param long + * longitude to project + * @return + * a projected Coordinate with values in metre + */ + def equalAreaProjection(lat: Double, long: Double): Coordinate = { + val latDist = + Math.PI * (EARTH_RADIUS.to(METRE).getValue.doubleValue() / 180.0) + val y = lat * latDist + val x = long * latDist * cos(toRadians(lat)) + new Coordinate(x, y) + } + + /** Draws a circle with a radius of the provided distance around the provided + * center coordinates and returns the result as a drawable polygon (one point + * per degree) + * + * @param center + * coordinate of the circle's center + * @param radius + * radius of the circle + * @return + * a polygon without the center, but with all points of the circle + */ + def buildCircle(center: Coordinate, radius: Quantity[Length]): Polygon = { + ??? } } diff --git a/src/main/scala/edu/ie3/util/osm/OsmEntities.scala b/src/main/scala/edu/ie3/util/osm/OsmEntities.scala index 6b6f5510..08d7d76b 100644 --- a/src/main/scala/edu/ie3/util/osm/OsmEntities.scala +++ b/src/main/scala/edu/ie3/util/osm/OsmEntities.scala @@ -5,18 +5,21 @@ */ package edu.ie3.util.osm -import edu.ie3.util.geo.GeoUtils.DEFAULT_GEOMETRY_FACTORY -import org.locationtech.jts.geom.impl.CoordinateArraySequence -import org.locationtech.jts.geom.{ - Coordinate, - GeometryFactory, - LinearRing, - Polygon, - PrecisionModel +import edu.ie3.util.geo.GeoUtils.{ + DEFAULT_GEOMETRY_FACTORY, + buildPolygon, + calcAreaOnEarth } +import edu.ie3.util.quantities.PowerSystemUnits +import edu.ie3.util.quantities.QuantityUtils.RichQuantityDouble +import org.locationtech.jts.geom.impl.CoordinateArraySequence +import org.locationtech.jts.geom.{Coordinate, LinearRing, Polygon} +import tech.units.indriya.ComparableQuantity +import tech.units.indriya.quantity.Quantities import java.time.ZonedDateTime import java.util.UUID +import javax.measure.quantity.Area import scala.jdk.CollectionConverters._ object OsmEntities { @@ -137,12 +140,11 @@ object OsmEntities { } def toPolygon: Polygon = { - val arrayCoordinates = new CoordinateArraySequence(getCoordinates.toArray) - val linearRing = - new LinearRing(arrayCoordinates, DEFAULT_GEOMETRY_FACTORY) - new Polygon(linearRing, Array[LinearRing](), DEFAULT_GEOMETRY_FACTORY) + buildPolygon(getCoordinates) } + def calculateArea: ComparableQuantity[Area] = + calcAreaOnEarth(toPolygon) } /** An OSM relation. diff --git a/src/main/scala/edu/ie3/util/osm/OsmUtils.scala b/src/main/scala/edu/ie3/util/osm/OsmUtils.scala new file mode 100644 index 00000000..f2bb4f72 --- /dev/null +++ b/src/main/scala/edu/ie3/util/osm/OsmUtils.scala @@ -0,0 +1,20 @@ +/* + * © 2021. TU Dortmund University, + * Institute of Energy Systems, Energy Efficiency and Energy Economics, + * Research group Distribution grid planning and operation +*/ +package edu.ie3.util.osm + +import edu.ie3.util.osm.OsmEntities.ClosedWay +import org.locationtech.jts.geom.Point + +object OsmUtils { + + def isInsideLandUse(coordinate: Point, landuses: List[ClosedWay]): Boolean = { + landuses.foreach(landuse => + if (landuse.toPolygon.convexHull().contains(coordinate)) return true + ) + false + } + +} diff --git a/src/main/scala/edu/ie3/util/osm/RichGeometries.scala b/src/main/scala/edu/ie3/util/osm/RichGeometries.scala new file mode 100644 index 00000000..0c539820 --- /dev/null +++ b/src/main/scala/edu/ie3/util/osm/RichGeometries.scala @@ -0,0 +1,135 @@ +/* + * © 2021. TU Dortmund University, + * Institute of Energy Systems, Energy Efficiency and Energy Economics, + * Research group Distribution grid planning and operation +*/ +package edu.ie3.util.osm + +import edu.ie3.util.exceptions.GeoException +import edu.ie3.util.geo.GeoUtils +import edu.ie3.util.geo.GeoUtils.{ + buildPolygon, + calcHaversine, + equalAreaProjection +} +import edu.ie3.util.quantities.PowerSystemUnits.KILOMETRE +import edu.ie3.util.quantities.QuantityUtils.RichQuantityDouble +import org.locationtech.jts.geom.{Coordinate, LineString, Polygon} +import tech.units.indriya.ComparableQuantity +import tech.units.indriya.quantity.Quantities + +import javax.measure.quantity.{Area, Length} +import scala.math.abs +import scala.util.{Failure, Success, Try} + +object RichGeometries { + + implicit class EarthCoordinate(coordinate: Coordinate) { + + /** Calculates the great circle disteance between two coordinates + * + * @param coordinateB + * coordinate b + * @return + * the distance between the coordinates as a quantity + */ + def haversineDistance( + coordinateB: Coordinate + ): ComparableQuantity[Length] = { + calcHaversine( + coordinate.getY, + coordinate.getX, + coordinateB.getY, + coordinateB.getX + ) + } + + /** Checks if the coordinate lies between two coordinates a and b. + * + * @param a + * coordinate a + * @param b + * coordinate b + * @param epsilon + * permitted deviation + * @return + * whether or not the coordinate lies between + */ + def isBetween( + a: Coordinate, + b: Coordinate, + epsilon: Double = 1e-11 + ): Boolean = { + val distance = a.haversineDistance(b) + val distancePassingMe = a + .haversineDistance(coordinate) + .add(coordinate.haversineDistance(b)) + .getValue + .doubleValue + if (abs(1 - (distancePassingMe / distance)) < epsilon) true + else false + } + } + + implicit class EarthLineString(lineString: LineString) { + + /** Compute length of a [[LineString]] on earth's surface. + * + * @return + * the length in kilometre as a quantity + */ + def haversineLength: ComparableQuantity[Length] = { + val coordinates = lineString.getCoordinates.toVector + val coordinateSize = coordinates.size + coordinates.zipWithIndex + .foldLeft(Quantities.getQuantity(0, KILOMETRE))((acc, current) => { + if (current._2 < coordinateSize - 1) { + val currentCoordinate = current._1 + val nextCoordinate = coordinates(current._2 + 1) + acc.add(currentCoordinate.haversineDistance(nextCoordinate)) + } else acc + }) + } + } + + implicit class EarthPolygon(polygon: Polygon) { + + // Fixme: Is not precisely right because earths curvature - use equal area projection ? + def intersect(polygonB: Polygon): Try[Set[Coordinate]] = { + Try(polygon.intersection(polygonB)) match { + case Failure(exception) => + Failure( + new GeoException( + s"Couldn't calculate intersection of polygons: ${polygon.toString} and ${polygonB.toString}. Reason:", + exception + ) + ) + case Success(geometry) => Success(geometry.getCoordinates.toSet) + } + } + + /** Calculates the area of a polygon on earth's surface. + * + * @return + * a Quantity of area in metre + */ + def calcAreaOnEarth: ComparableQuantity[Area] = { + equalAreaProjection.getArea.asSquareMetre + } + + /** Does an equal area projection of the polygon onto a two-dimensional + * surface to account for earth's curvature when calculating the polygon's + * area. + * + * @return + * the projected polygon + */ + def equalAreaProjection: Polygon = { + val projectedCoordinates = polygon.getCoordinates.map(coordinate => + GeoUtils.equalAreaProjection(coordinate.y, coordinate.getX) + ) + buildPolygon(projectedCoordinates.toList) + } + + } +} diff --git a/src/main/scala/edu/ie3/util/quantities/QuantityUtils.scala b/src/main/scala/edu/ie3/util/quantities/QuantityUtils.scala index 7d9993e0..701a0a19 100644 --- a/src/main/scala/edu/ie3/util/quantities/QuantityUtils.scala +++ b/src/main/scala/edu/ie3/util/quantities/QuantityUtils.scala @@ -24,11 +24,20 @@ import edu.ie3.util.quantities.interfaces.{ } import tech.units.indriya.ComparableQuantity import tech.units.indriya.quantity.Quantities -import tech.units.indriya.unit.Units.{AMPERE, OHM, PERCENT, SIEMENS, VOLT, WATT} +import tech.units.indriya.unit.Units.{ + AMPERE, + OHM, + PERCENT, + SIEMENS, + SQUARE_METRE, + VOLT, + WATT +} import javax.measure.MetricPrefix import javax.measure.quantity.{ Angle, + Area, Dimensionless, ElectricConductance, ElectricCurrent, @@ -57,6 +66,9 @@ object QuantityUtils { /* indriya units */ + def asSquareMetre: ComparableQuantity[Area] = + Quantities.getQuantity(value, SQUARE_METRE) + def asVolt: ComparableQuantity[ElectricPotential] = Quantities.getQuantity(value, VOLT) From 4b262d5cd4cb3aa2b31751937b19e7a1d9dfe395 Mon Sep 17 00:00:00 2001 From: t-ober <63147366+t-ober@users.noreply.github.com> Date: Mon, 6 Dec 2021 11:30:45 +0100 Subject: [PATCH 04/17] refactor GeoUtils --- ...java => DeprecatedCoordinateDistance.java} | 10 +- .../edu/ie3/util/geo/DeprecatedGeoUtils.java | 4 +- .../edu/ie3/util/geo/CoordinateDistance.scala | 96 ++++++++++++ .../scala/edu/ie3/util/geo/GeoUtils.scala | 148 +++++++++++++++++- .../util/{osm => geo}/RichGeometries.scala | 35 +++-- .../scala/edu/ie3/util/osm/OsmEntities.scala | 15 +- ...> DeprecatedCoordinateDistanceTest.groovy} | 14 +- .../edu/ie3/util/geo/GeoUtilsTest.groovy | 8 +- 8 files changed, 289 insertions(+), 41 deletions(-) rename src/main/java/edu/ie3/util/geo/{CoordinateDistance.java => DeprecatedCoordinateDistance.java} (90%) create mode 100644 src/main/scala/edu/ie3/util/geo/CoordinateDistance.scala rename src/main/scala/edu/ie3/util/{osm => geo}/RichGeometries.scala (81%) rename src/test/groovy/edu/ie3/util/geo/{CoordinateDistanceTest.groovy => DeprecatedCoordinateDistanceTest.groovy} (61%) diff --git a/src/main/java/edu/ie3/util/geo/CoordinateDistance.java b/src/main/java/edu/ie3/util/geo/DeprecatedCoordinateDistance.java similarity index 90% rename from src/main/java/edu/ie3/util/geo/CoordinateDistance.java rename to src/main/java/edu/ie3/util/geo/DeprecatedCoordinateDistance.java index 2be7b554..ff671b4c 100644 --- a/src/main/java/edu/ie3/util/geo/CoordinateDistance.java +++ b/src/main/java/edu/ie3/util/geo/DeprecatedCoordinateDistance.java @@ -13,7 +13,7 @@ * Wraps two coordinates with the distance between the first one and the second one, can be compared * by distance to another CoordinateDistance */ -public class CoordinateDistance implements Comparable { +public class DeprecatedCoordinateDistance implements Comparable { private final org.locationtech.jts.geom.Point coordinateA; private final org.locationtech.jts.geom.Point coordinateB; private final ComparableQuantity distance; @@ -25,7 +25,7 @@ public class CoordinateDistance implements Comparable { * @param coordinateA The first coordinate * @param coordinateB The second coordinate */ - public CoordinateDistance( + public DeprecatedCoordinateDistance( org.locationtech.jts.geom.Point coordinateA, org.locationtech.jts.geom.Point coordinateB) { this( coordinateA, @@ -39,7 +39,7 @@ public CoordinateDistance( * @param coordinateB The second coordinate * @param distance The distance from A to B */ - private CoordinateDistance( + private DeprecatedCoordinateDistance( org.locationtech.jts.geom.Point coordinateA, org.locationtech.jts.geom.Point coordinateB, ComparableQuantity distance) { @@ -72,7 +72,7 @@ public ComparableQuantity getDistance() { * number higher than 0 if that has a lower distance */ @Override - public int compareTo(CoordinateDistance that) { + public int compareTo(DeprecatedCoordinateDistance that) { return this.distance.compareTo(that.distance); } @@ -80,7 +80,7 @@ public int compareTo(CoordinateDistance that) { public boolean equals(Object o) { if (this == o) return true; if (o == null || getClass() != o.getClass()) return false; - CoordinateDistance that = (CoordinateDistance) o; + DeprecatedCoordinateDistance that = (DeprecatedCoordinateDistance) o; return coordinateA.equals(that.coordinateA) && coordinateB.equals(that.coordinateB) && distance.equals(that.distance); diff --git a/src/main/java/edu/ie3/util/geo/DeprecatedGeoUtils.java b/src/main/java/edu/ie3/util/geo/DeprecatedGeoUtils.java index fcaf11c9..a136aae2 100644 --- a/src/main/java/edu/ie3/util/geo/DeprecatedGeoUtils.java +++ b/src/main/java/edu/ie3/util/geo/DeprecatedGeoUtils.java @@ -103,11 +103,11 @@ public static ComparableQuantity totalLengthOfLineString(LineString line * @param coordinates the points to calculate the distance from the base point for * @return a sorted set of distances between the base and other coordinates */ - public static SortedSet getCoordinateDistances( + public static SortedSet getCoordinateDistances( org.locationtech.jts.geom.Point baseCoordinate, Collection coordinates) { return coordinates.stream() - .map(coordinate -> new CoordinateDistance(baseCoordinate, coordinate)) + .map(coordinate -> new DeprecatedCoordinateDistance(baseCoordinate, coordinate)) .collect(Collectors.toCollection(TreeSet::new)); } diff --git a/src/main/scala/edu/ie3/util/geo/CoordinateDistance.scala b/src/main/scala/edu/ie3/util/geo/CoordinateDistance.scala new file mode 100644 index 00000000..b1190fd1 --- /dev/null +++ b/src/main/scala/edu/ie3/util/geo/CoordinateDistance.scala @@ -0,0 +1,96 @@ +/* + * © 2021. TU Dortmund University, + * Institute of Energy Systems, Energy Efficiency and Energy Economics, + * Research group Distribution grid planning and operation +*/ +package edu.ie3.util.geo + +import org.locationtech.jts.geom.Point + +import java.util.Objects +import javax.measure.quantity.Length +import tech.units.indriya.ComparableQuantity + +/** Wraps two coordinates with the distance between the first one and the second + * one, can be compared by distance to another CoordinateDistance + * @param coordinateA + * The first coordinate + * @param coordinateB + * The second coordinate + * @param distance + * The distance from A to B + */ +class CoordinateDistance private ( + val coordinateA: Point, + val coordinateB: Point, + val distance: ComparableQuantity[Length] +) extends Comparable[CoordinateDistance] { + + /** Calculates the distance from the first to the second coordinate using + * [[GeoUtils.calcHaversine(double, double, double, double)]] + * + * @param coordinateA + * The first coordinate + * @param coordinateB + * The second coordinate + */ + def this(coordinateA: Point, coordinateB: Point) = { + this( + coordinateA, + coordinateB, + GeoUtils.calcHaversine( + coordinateA.getY, + coordinateA.getX, + coordinateB.getY, + coordinateB.getX + ) + ) + } + + /** @return The first coordinate */ + def getCoordinateA: Point = { + coordinateA + } + + /** @return The second coordinate */ + def getCoordinateB: Point = { + coordinateB + } + + /** @return + * The distance from the first coordinate to the second coordinate in km + */ + def getDistance: ComparableQuantity[Length] = { + distance + } + + /** Compares two coordinate distances on the length of the distance alone, + * thus having a natural ordering that is inconsistent with equals + * + * @param that + * the distance to compare + * @return + * a number lower than 0 if this has a lower distance than that, 0 if they + * are the same, a number higher than 0 if that has a lower distance + */ + override def compareTo(that: CoordinateDistance): Int = { + this.distance.compareTo(that.distance) + } + + override def equals(other: Any): Boolean = other match { + case other: CoordinateDistance if this == other => true + case other: CoordinateDistance => + coordinateA.equals(other.coordinateA) && coordinateB.equals( + other.coordinateB + ) && distance == other.distance + case _ => false + } + + override def hashCode: Int = { + Objects.hash(coordinateA, coordinateB, distance) + } + + override def toString: String = { + "CoordinateDistance{" + "coordinateA=" + coordinateA + ", coordinateB=" + coordinateB + ", distance=" + distance + '}' + } +} diff --git a/src/main/scala/edu/ie3/util/geo/GeoUtils.scala b/src/main/scala/edu/ie3/util/geo/GeoUtils.scala index a03b7bd4..d305a32f 100644 --- a/src/main/scala/edu/ie3/util/geo/GeoUtils.scala +++ b/src/main/scala/edu/ie3/util/geo/GeoUtils.scala @@ -7,6 +7,7 @@ package edu.ie3.util.geo import edu.ie3.util.exceptions.GeoException import edu.ie3.util.quantities.PowerSystemUnits.KILOMETRE +import org.apache.commons.lang3.ArrayUtils import org.locationtech.jts.algorithm.ConvexHull import org.locationtech.jts.geom.impl.CoordinateArraySequence import org.locationtech.jts.geom.{ @@ -24,8 +25,11 @@ import tech.units.indriya.quantity.Quantities import tech.units.indriya.unit.Units.{METRE, RADIAN} import java.lang.Math.{atan2, cos, sin, sqrt, toRadians} +import java.util import javax.measure.Quantity import javax.measure.quantity.Length +import scala.collection.immutable.{SortedSet, TreeSet} +import scala.jdk.CollectionConverters.IterableHasAsJava import scala.math.pow import scala.util.{Failure, Success, Try} @@ -36,7 +40,133 @@ object GeoUtils { val EARTH_RADIUS: ComparableQuantity[Length] = Quantities.getQuantity(6378137.0, METRE); - /** Calculates the great circle disteance between two coordinates + /** Build an instance of [[LineString]] between two points that is safe to be + * compared even if the provided two points consist of exactly the same + * coordinates. This is done by increasing the coordinate of the provided + * Point `p1` by a small amount to make it different from Point `p2`. For + * details on the bug inside [[LineString]] that is addressed here, see + * https://github.com/locationtech/jts/issues/531 + * + * @param p1 + * start point of the linestring + * @param p2 + * end point of the linestring + * @return + * a [[LineString]] between the provided points + */ + def buildSafeLineStringBetweenPoints(p1: Point, p2: Point): LineString = { + val safePoint1 = if (p1 == p2) buildSafePoint(p1) else p1 + DEFAULT_GEOMETRY_FACTORY.createLineString( + safePoint1.getCoordinates ++ p2.getCoordinates + ) + } + + /** Build an instance of [[LineString]] between two coordinates that is safe + * to be compared even if the provided two coordinates are exactly the same + * coordinates. This is done by increasing the coordinate of the provided + * Point `code c1` by a small amount to make it different from Point `code + * c2`. For details on the bug inside [[LineString]] that is addressed here, + * see https://github.com/locationtech/jts/issues/531 + * + * @param c1 + * start coordinate of the linestring + * @param c2 + * end coordinate of the linestring + * @return + * A safely build line string + */ + def buildSafeLineStringBetweenCoords( + c1: Coordinate, + c2: Coordinate + ): LineString = { + val safeCoord1: Coordinate = if (c1 == c2) buildSafeCoord(c1) else c1 + DEFAULT_GEOMETRY_FACTORY.createLineString( + ArrayUtils.addAll(Array[Coordinate](safeCoord1), c2) + ) + } + + /** Adapt the provided point as described in buildSafeCoord ( Coordinate )} + * and return a new, adapted instance of [[Point]] + * + * @param p1 + * the point that should be adapted + * @return + * the adapted point with a slightly changed coordinate + */ + private def buildSafePoint(p1: Point): Point = { + val safeCoord = buildSafeCoord(p1.getCoordinate) + val safeCoordSeq = new CoordinateArraySequence(Array[Coordinate](safeCoord)) + new Point(safeCoordSeq, p1.getFactory) + } + + /** Adapted [[Coordinate]] x, [[Coordinate]] y and [[Coordinate]] z of the + * provided [[Coordinate]] by 1e-13 and return a new, adapted instance of + * [[Coordinate]] + * + * @param coord + * the coordinate that should be adapted + * @return + * the adapted coordinate with slightly changed x,y,z values + */ + private def buildSafeCoord(coord: Coordinate): Coordinate = { + val modVal: Double = 1e-13 + val p1X: Double = coord.getX + modVal + val p1Y: Double = coord.getY + modVal + val p1Z: Double = coord.getZ + modVal + new Coordinate(p1X, p1Y, p1Z) + } + + /** Convert a given [[LineString]] with at least two points into a 'safe to be + * compared' [[LineString]] This is done by removing duplicates in the points + * in the provided linestring as well as a small change of the start + * coordinate if the linestring only consists of two coordinates. For details + * on the bug inside [[LineString]] that is addressed here, see + * https://github.com/locationtech/jts/issues/531 + * + * @param lineString + * the linestring that should be checked and maybe converted to a 'safe to + * be compared' linestring + * @return + * a 'safe to be compared' linestring + */ + def buildSafeLineString(lineString: LineString): LineString = + if (lineString.getCoordinates.length == 2) + buildSafeLineStringBetweenPoints( + lineString.getStartPoint, + lineString.getEndPoint + ) + // rebuild line with unique points + else { + val uniqueCoords: Array[Coordinate] = lineString.getCoordinates.distinct + if (uniqueCoords.length == 1) + buildSafeLineStringBetweenPoints( + lineString.getStartPoint, + lineString.getEndPoint + ) + else DEFAULT_GEOMETRY_FACTORY.createLineString(uniqueCoords) + } + + /** Calculates and orders the coordinate distances from a base coordinate to a + * list of coordinates + * + * @param baseCoordinate + * the base coordinate + * @param coordinates + * the coordinates to calculate the distance for + * @return + * a sorted set of [[CoordinateDistance]] + */ + def calcOrderedCoordinateDistances( + baseCoordinate: Point, + coordinates: List[Point] + ): SortedSet[CoordinateDistance] = { + val coordinateDistances = coordinates.map(coordinate => + new CoordinateDistance(baseCoordinate, coordinate) + ) + TreeSet(coordinateDistances: _*) + } + + /** Calculates the great circle distance between two coordinates * * @param latA * latitude of coordinate a @@ -65,6 +195,7 @@ object GeoUtils { r.multiply(c) } + // Fixme: Might be slightly off as it does not account for earth's curvature /** Builds a convex hull from a set of coordinates. * * @param coordinates @@ -158,8 +289,19 @@ object GeoUtils { * @return * a polygon without the center, but with all points of the circle */ - def buildCircle(center: Coordinate, radius: Quantity[Length]): Polygon = { - ??? + def buildCirclePolygon( + center: Coordinate, + radius: Quantity[Length] + ): Polygon = { + val coordinates = (0 to 359) + .map(deg => { + val angle = deg.toRadians + val x = center.x + radius.getValue.doubleValue() * cos(angle) + val y = center.y + radius.getValue.doubleValue() * sin(angle) + new Coordinate(x, y) + }) + .toList + buildPolygon(coordinates) } } diff --git a/src/main/scala/edu/ie3/util/osm/RichGeometries.scala b/src/main/scala/edu/ie3/util/geo/RichGeometries.scala similarity index 81% rename from src/main/scala/edu/ie3/util/osm/RichGeometries.scala rename to src/main/scala/edu/ie3/util/geo/RichGeometries.scala index 0c539820..eae27655 100644 --- a/src/main/scala/edu/ie3/util/osm/RichGeometries.scala +++ b/src/main/scala/edu/ie3/util/geo/RichGeometries.scala @@ -3,18 +3,16 @@ * Institute of Energy Systems, Energy Efficiency and Energy Economics, * Research group Distribution grid planning and operation */ -package edu.ie3.util.osm +package edu.ie3.util.geo import edu.ie3.util.exceptions.GeoException -import edu.ie3.util.geo.GeoUtils import edu.ie3.util.geo.GeoUtils.{ + DEFAULT_GEOMETRY_FACTORY, buildPolygon, - calcHaversine, - equalAreaProjection + calcHaversine } import edu.ie3.util.quantities.PowerSystemUnits.KILOMETRE -import edu.ie3.util.quantities.QuantityUtils.RichQuantityDouble -import org.locationtech.jts.geom.{Coordinate, LineString, Polygon} +import org.locationtech.jts.geom.{Coordinate, LineString, Point, Polygon} import tech.units.indriya.ComparableQuantity import tech.units.indriya.quantity.Quantities @@ -24,7 +22,7 @@ import scala.util.{Failure, Success, Try} object RichGeometries { - implicit class EarthCoordinate(coordinate: Coordinate) { + implicit class RichCoordinate(coordinate: Coordinate) { /** Calculates the great circle disteance between two coordinates * @@ -69,9 +67,13 @@ object RichGeometries { if (abs(1 - (distancePassingMe / distance)) < epsilon) true else false } + + def toPoint: Point = { + DEFAULT_GEOMETRY_FACTORY.createPoint(coordinate) + } } - implicit class EarthLineString(lineString: LineString) { + implicit class RichLineString(lineString: LineString) { /** Compute length of a [[LineString]] on earth's surface. * @@ -92,10 +94,11 @@ object RichGeometries { } } - implicit class EarthPolygon(polygon: Polygon) { + implicit class RichPolygon(polygon: Polygon) { // Fixme: Is not precisely right because earths curvature - use equal area projection ? def intersect(polygonB: Polygon): Try[Set[Coordinate]] = { + Try(polygon.intersection(polygonB)) match { case Failure(exception) => Failure( @@ -131,5 +134,19 @@ object RichGeometries { buildPolygon(projectedCoordinates.toList) } + // Fixme: Is not precisely right because earths curvatur + + /** Checks whether the polygon contains the coordinate. Uses "covers()" + * insted of "contains()" so borders are included. + * + * @param coordinate + * the coordinate to check + * @return + * whether the polygon contains the coordinate + */ + def containsCoordinate(coordinate: Coordinate): Boolean = { + polygon.covers(coordinate.toPoint) + } + } } diff --git a/src/main/scala/edu/ie3/util/osm/OsmEntities.scala b/src/main/scala/edu/ie3/util/osm/OsmEntities.scala index 08d7d76b..80b6ccf4 100644 --- a/src/main/scala/edu/ie3/util/osm/OsmEntities.scala +++ b/src/main/scala/edu/ie3/util/osm/OsmEntities.scala @@ -5,17 +5,10 @@ */ package edu.ie3.util.osm -import edu.ie3.util.geo.GeoUtils.{ - DEFAULT_GEOMETRY_FACTORY, - buildPolygon, - calcAreaOnEarth -} -import edu.ie3.util.quantities.PowerSystemUnits -import edu.ie3.util.quantities.QuantityUtils.RichQuantityDouble -import org.locationtech.jts.geom.impl.CoordinateArraySequence -import org.locationtech.jts.geom.{Coordinate, LinearRing, Polygon} +import edu.ie3.util.geo.GeoUtils.buildPolygon +import edu.ie3.util.geo.RichGeometries.RichPolygon +import org.locationtech.jts.geom.{Coordinate, Polygon} import tech.units.indriya.ComparableQuantity -import tech.units.indriya.quantity.Quantities import java.time.ZonedDateTime import java.util.UUID @@ -144,7 +137,7 @@ object OsmEntities { } def calculateArea: ComparableQuantity[Area] = - calcAreaOnEarth(toPolygon) + toPolygon.calcAreaOnEarth } /** An OSM relation. diff --git a/src/test/groovy/edu/ie3/util/geo/CoordinateDistanceTest.groovy b/src/test/groovy/edu/ie3/util/geo/DeprecatedCoordinateDistanceTest.groovy similarity index 61% rename from src/test/groovy/edu/ie3/util/geo/CoordinateDistanceTest.groovy rename to src/test/groovy/edu/ie3/util/geo/DeprecatedCoordinateDistanceTest.groovy index e1e8a009..691e988b 100644 --- a/src/test/groovy/edu/ie3/util/geo/CoordinateDistanceTest.groovy +++ b/src/test/groovy/edu/ie3/util/geo/DeprecatedCoordinateDistanceTest.groovy @@ -7,7 +7,7 @@ package edu.ie3.util.geo import spock.lang.Specification -class CoordinateDistanceTest extends Specification { +class DeprecatedCoordinateDistanceTest extends Specification { def "The constructor without a distance parameter calculates the distance as expected"() { given: def pointA = DeprecatedGeoUtils.xyToPoint(49d, 7d) @@ -15,8 +15,8 @@ class CoordinateDistanceTest extends Specification { def expectedDistance = DeprecatedGeoUtils.calcHaversine(pointA.y, pointA.x, pointB.y, pointB.x) when: - def coordinateDistance = new CoordinateDistance(pointA, pointB) - def expectedCoordinateDistance = new CoordinateDistance(coordinateDistance.coordinateA, coordinateDistance.coordinateB, expectedDistance) + def coordinateDistance = new DeprecatedCoordinateDistance(pointA, pointB) + def expectedCoordinateDistance = new DeprecatedCoordinateDistance(coordinateDistance.coordinateA, coordinateDistance.coordinateB, expectedDistance) then: coordinateDistance.distance == expectedDistance @@ -27,10 +27,10 @@ class CoordinateDistanceTest extends Specification { def "CoordinateDistances are sortable using their distance field"() { given: def basePoint = DeprecatedGeoUtils.xyToPoint(49d, 7d) - def distA = new CoordinateDistance(basePoint, DeprecatedGeoUtils.xyToPoint(50d, 7d)) - def distB = new CoordinateDistance(basePoint, DeprecatedGeoUtils.xyToPoint(50d, 7.1d)) - def distC = new CoordinateDistance(basePoint, DeprecatedGeoUtils.xyToPoint(49d, 7.1d)) - def distD = new CoordinateDistance(basePoint, DeprecatedGeoUtils.xyToPoint(52d, 9d)) + def distA = new DeprecatedCoordinateDistance(basePoint, DeprecatedGeoUtils.xyToPoint(50d, 7d)) + def distB = new DeprecatedCoordinateDistance(basePoint, DeprecatedGeoUtils.xyToPoint(50d, 7.1d)) + def distC = new DeprecatedCoordinateDistance(basePoint, DeprecatedGeoUtils.xyToPoint(49d, 7.1d)) + def distD = new DeprecatedCoordinateDistance(basePoint, DeprecatedGeoUtils.xyToPoint(52d, 9d)) def coordinateDistances = [distA, distB, distC, distD] when: diff --git a/src/test/groovy/edu/ie3/util/geo/GeoUtilsTest.groovy b/src/test/groovy/edu/ie3/util/geo/GeoUtilsTest.groovy index 5401cb2c..cd08d850 100644 --- a/src/test/groovy/edu/ie3/util/geo/GeoUtilsTest.groovy +++ b/src/test/groovy/edu/ie3/util/geo/GeoUtilsTest.groovy @@ -74,10 +74,10 @@ class GeoUtilsTest extends Specification { DeprecatedGeoUtils.xyToPoint(52d, 9d) ] def coordinateDistances = [ - new CoordinateDistance(basePoint, points[0]), - new CoordinateDistance(basePoint, points[1]), - new CoordinateDistance(basePoint, points[2]), - new CoordinateDistance(basePoint, points[3]) + new DeprecatedCoordinateDistance(basePoint, points[0]), + new DeprecatedCoordinateDistance(basePoint, points[1]), + new DeprecatedCoordinateDistance(basePoint, points[2]), + new DeprecatedCoordinateDistance(basePoint, points[3]) ] expect: DeprecatedGeoUtils.getCoordinateDistances(basePoint, points) == new TreeSet(coordinateDistances) From 755c75d28705e438d1cef00a6397eff7e8da01b6 Mon Sep 17 00:00:00 2001 From: t-ober <63147366+t-ober@users.noreply.github.com> Date: Mon, 6 Dec 2021 12:12:13 +0100 Subject: [PATCH 05/17] safe line string test --- .../scala/edu/ie3/util/geo/GeoUtils.scala | 64 +++++++++---------- .../edu/ie3/util/geo/RichGeometries.scala | 11 ++-- .../edu/ie3/util/geo/GeoUtilsTest.groovy | 2 +- .../scala/edu/ie3/util/geo/GeoUtilsSpec.scala | 23 +++++++ 4 files changed, 60 insertions(+), 40 deletions(-) create mode 100644 src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala diff --git a/src/main/scala/edu/ie3/util/geo/GeoUtils.scala b/src/main/scala/edu/ie3/util/geo/GeoUtils.scala index d305a32f..8d15f32c 100644 --- a/src/main/scala/edu/ie3/util/geo/GeoUtils.scala +++ b/src/main/scala/edu/ie3/util/geo/GeoUtils.scala @@ -25,11 +25,9 @@ import tech.units.indriya.quantity.Quantities import tech.units.indriya.unit.Units.{METRE, RADIAN} import java.lang.Math.{atan2, cos, sin, sqrt, toRadians} -import java.util import javax.measure.Quantity import javax.measure.quantity.Length import scala.collection.immutable.{SortedSet, TreeSet} -import scala.jdk.CollectionConverters.IterableHasAsJava import scala.math.pow import scala.util.{Failure, Success, Try} @@ -40,6 +38,37 @@ object GeoUtils { val EARTH_RADIUS: ComparableQuantity[Length] = Quantities.getQuantity(6378137.0, METRE); + + /** Convert a given [[LineString]] with at least two points into a 'safe to be + * compared' [[LineString]] This is done by removing duplicates in the points + * in the provided linestring as well as a small change of the start + * coordinate if the linestring only consists of two coordinates. For details + * on the bug inside [[LineString]] that is addressed here, see + * https://github.com/locationtech/jts/issues/531 + * + * @param lineString + * the linestring that should be checked and maybe converted to a 'safe to + * be compared' linestring + * @return + * a 'safe to be compared' linestring + */ + def buildSafeLineString(lineString: LineString): LineString = + if (lineString.getCoordinates.length == 2) + buildSafeLineStringBetweenPoints( + lineString.getStartPoint, + lineString.getEndPoint + ) + // rebuild line with unique points + else { + val uniqueCoords: Array[Coordinate] = lineString.getCoordinates.distinct + if (uniqueCoords.length == 1) + buildSafeLineStringBetweenPoints( + lineString.getStartPoint, + lineString.getEndPoint + ) + else DEFAULT_GEOMETRY_FACTORY.createLineString(uniqueCoords) + } + /** Build an instance of [[LineString]] between two points that is safe to be * compared even if the provided two points consist of exactly the same * coordinates. This is done by increasing the coordinate of the provided @@ -116,35 +145,6 @@ object GeoUtils { new Coordinate(p1X, p1Y, p1Z) } - /** Convert a given [[LineString]] with at least two points into a 'safe to be - * compared' [[LineString]] This is done by removing duplicates in the points - * in the provided linestring as well as a small change of the start - * coordinate if the linestring only consists of two coordinates. For details - * on the bug inside [[LineString]] that is addressed here, see - * https://github.com/locationtech/jts/issues/531 - * - * @param lineString - * the linestring that should be checked and maybe converted to a 'safe to - * be compared' linestring - * @return - * a 'safe to be compared' linestring - */ - def buildSafeLineString(lineString: LineString): LineString = - if (lineString.getCoordinates.length == 2) - buildSafeLineStringBetweenPoints( - lineString.getStartPoint, - lineString.getEndPoint - ) - // rebuild line with unique points - else { - val uniqueCoords: Array[Coordinate] = lineString.getCoordinates.distinct - if (uniqueCoords.length == 1) - buildSafeLineStringBetweenPoints( - lineString.getStartPoint, - lineString.getEndPoint - ) - else DEFAULT_GEOMETRY_FACTORY.createLineString(uniqueCoords) - } /** Calculates and orders the coordinate distances from a base coordinate to a * list of coordinates @@ -295,7 +295,7 @@ object GeoUtils { ): Polygon = { val coordinates = (0 to 359) .map(deg => { - val angle = deg.toRadians + val angle = deg.toFloat.toRadians val x = center.x + radius.getValue.doubleValue() * cos(angle) val y = center.y + radius.getValue.doubleValue() * sin(angle) new Coordinate(x, y) diff --git a/src/main/scala/edu/ie3/util/geo/RichGeometries.scala b/src/main/scala/edu/ie3/util/geo/RichGeometries.scala index eae27655..a374662f 100644 --- a/src/main/scala/edu/ie3/util/geo/RichGeometries.scala +++ b/src/main/scala/edu/ie3/util/geo/RichGeometries.scala @@ -6,12 +6,9 @@ package edu.ie3.util.geo import edu.ie3.util.exceptions.GeoException -import edu.ie3.util.geo.GeoUtils.{ - DEFAULT_GEOMETRY_FACTORY, - buildPolygon, - calcHaversine -} +import edu.ie3.util.geo.GeoUtils.{DEFAULT_GEOMETRY_FACTORY, buildPolygon, calcHaversine} import edu.ie3.util.quantities.PowerSystemUnits.KILOMETRE +import edu.ie3.util.quantities.QuantityUtils.RichQuantityDouble import org.locationtech.jts.geom.{Coordinate, LineString, Point, Polygon} import tech.units.indriya.ComparableQuantity import tech.units.indriya.quantity.Quantities @@ -24,7 +21,7 @@ object RichGeometries { implicit class RichCoordinate(coordinate: Coordinate) { - /** Calculates the great circle disteance between two coordinates + /** Calculates the great circle distance between two coordinates * * @param coordinateB * coordinate b @@ -64,7 +61,7 @@ object RichGeometries { .add(coordinate.haversineDistance(b)) .getValue .doubleValue - if (abs(1 - (distancePassingMe / distance)) < epsilon) true + if (abs(1 - (distancePassingMe / distance.getValue.doubleValue())) < epsilon) true else false } diff --git a/src/test/groovy/edu/ie3/util/geo/GeoUtilsTest.groovy b/src/test/groovy/edu/ie3/util/geo/GeoUtilsTest.groovy index cd08d850..7bc20d01 100644 --- a/src/test/groovy/edu/ie3/util/geo/GeoUtilsTest.groovy +++ b/src/test/groovy/edu/ie3/util/geo/GeoUtilsTest.groovy @@ -96,7 +96,7 @@ class GeoUtilsTest extends Specification { for (int cnt = 0; cnt < coordinates.length; cnt++) { coordinates[cnt] == actualCoordinates[cnt] } - +from where: lineString | coordinates "{ \"type\": \"LineString\", \"coordinates\": [[7.411111, 51.49228],[7.411111, 51.49228]]}" | [ diff --git a/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala b/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala new file mode 100644 index 00000000..8a3bc477 --- /dev/null +++ b/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala @@ -0,0 +1,23 @@ +package edu.ie3.util.geo + +import edu.ie3.util.geo.GeoUtils.DEFAULT_GEOMETRY_FACTORY +import edu.ie3.util.geo.RichGeometries.RichCoordinate +import org.locationtech.jts.geom.Coordinate +import org.scalatest.matchers.should.Matchers +import org.scalatest.wordspec.AnyWordSpecLike + +class GeoUtilsSpec extends Matchers with AnyWordSpecLike{ + + "The GeoUtils" should { + "build a line string with two equal geo coordinates, avoiding the known bug in jts geometry" in { + val coordinate = new Coordinate(7.411111, 51.49228) + val coordinateWithCorrection = new Coordinate(7.4111110000001, 51.4922800000001) + val fromCoords = GeoUtils.buildSafeLineStringBetweenCoords(coordinate, coordinate).getCoordinates + val fromPoints = GeoUtils.buildSafeLineStringBetweenPoints(coordinate.toPoint, coordinate.toPoint).getCoordinates + val fromLineString = GeoUtils.buildSafeLineString(DEFAULT_GEOMETRY_FACTORY.createLineString(Array[Coordinate](coordinate, coordinate))) + fromCoords(0) == coordinateWithCorrection + fromCoords(1) == coordinate + } + } + +} From 9fda4483428239f9f6fdc4264554278e37e5004d Mon Sep 17 00:00:00 2001 From: t-ober <63147366+t-ober@users.noreply.github.com> Date: Wed, 8 Dec 2021 13:43:14 +0100 Subject: [PATCH 06/17] moar tests --- .../edu/ie3/util/geo/CoordinateDistance.scala | 46 ++--- .../scala/edu/ie3/util/geo/GeoUtils.scala | 108 +++++++---- .../edu/ie3/util/geo/RichGeometries.scala | 12 +- .../edu/ie3/util/geo/GeoUtilsTest.groovy | 4 +- .../scala/edu/ie3/util/geo/GeoUtilsSpec.scala | 181 +++++++++++++++++- .../edu/ie3/util/geo/RichGeometriesSpec.scala | 36 ++++ 6 files changed, 319 insertions(+), 68 deletions(-) create mode 100644 src/test/scala/edu/ie3/util/geo/RichGeometriesSpec.scala diff --git a/src/main/scala/edu/ie3/util/geo/CoordinateDistance.scala b/src/main/scala/edu/ie3/util/geo/CoordinateDistance.scala index b1190fd1..cfdda498 100644 --- a/src/main/scala/edu/ie3/util/geo/CoordinateDistance.scala +++ b/src/main/scala/edu/ie3/util/geo/CoordinateDistance.scala @@ -13,48 +13,48 @@ import tech.units.indriya.ComparableQuantity /** Wraps two coordinates with the distance between the first one and the second * one, can be compared by distance to another CoordinateDistance - * @param coordinateA + * @param from * The first coordinate - * @param coordinateB + * @param to * The second coordinate * @param distance * The distance from A to B */ class CoordinateDistance private ( - val coordinateA: Point, - val coordinateB: Point, + val from: Point, + val to: Point, val distance: ComparableQuantity[Length] ) extends Comparable[CoordinateDistance] { /** Calculates the distance from the first to the second coordinate using * [[GeoUtils.calcHaversine(double, double, double, double)]] * - * @param coordinateA + * @param from * The first coordinate - * @param coordinateB + * @param to * The second coordinate */ - def this(coordinateA: Point, coordinateB: Point) = { + def this(from: Point, to: Point) = { this( - coordinateA, - coordinateB, + from, + to, GeoUtils.calcHaversine( - coordinateA.getY, - coordinateA.getX, - coordinateB.getY, - coordinateB.getX + from.getY, + from.getX, + to.getY, + to.getX ) ) } - /** @return The first coordinate */ - def getCoordinateA: Point = { - coordinateA + /** The first coordinate */ + def getFrom: Point = { + from } - /** @return The second coordinate */ - def getCoordinateB: Point = { - coordinateB + /** The second coordinate */ + def getTo: Point = { + to } /** @return @@ -80,17 +80,17 @@ class CoordinateDistance private ( override def equals(other: Any): Boolean = other match { case other: CoordinateDistance if this == other => true case other: CoordinateDistance => - coordinateA.equals(other.coordinateA) && coordinateB.equals( - other.coordinateB + from.equals(other.from) && to.equals( + other.to ) && distance == other.distance case _ => false } override def hashCode: Int = { - Objects.hash(coordinateA, coordinateB, distance) + Objects.hash(from, to, distance) } override def toString: String = { - "CoordinateDistance{" + "coordinateA=" + coordinateA + ", coordinateB=" + coordinateB + ", distance=" + distance + '}' + "CoordinateDistance{" + "coordinateA=" + from + ", coordinateB=" + to + ", distance=" + distance + '}' } } diff --git a/src/main/scala/edu/ie3/util/geo/GeoUtils.scala b/src/main/scala/edu/ie3/util/geo/GeoUtils.scala index 8d15f32c..8fef07dc 100644 --- a/src/main/scala/edu/ie3/util/geo/GeoUtils.scala +++ b/src/main/scala/edu/ie3/util/geo/GeoUtils.scala @@ -6,7 +6,6 @@ package edu.ie3.util.geo import edu.ie3.util.exceptions.GeoException -import edu.ie3.util.quantities.PowerSystemUnits.KILOMETRE import org.apache.commons.lang3.ArrayUtils import org.locationtech.jts.algorithm.ConvexHull import org.locationtech.jts.geom.impl.CoordinateArraySequence @@ -24,7 +23,7 @@ import tech.units.indriya.ComparableQuantity import tech.units.indriya.quantity.Quantities import tech.units.indriya.unit.Units.{METRE, RADIAN} -import java.lang.Math.{atan2, cos, sin, sqrt, toRadians} +import java.lang.Math._ import javax.measure.Quantity import javax.measure.quantity.Length import scala.collection.immutable.{SortedSet, TreeSet} @@ -38,20 +37,19 @@ object GeoUtils { val EARTH_RADIUS: ComparableQuantity[Length] = Quantities.getQuantity(6378137.0, METRE); - /** Convert a given [[LineString]] with at least two points into a 'safe to be - * compared' [[LineString]] This is done by removing duplicates in the points - * in the provided linestring as well as a small change of the start - * coordinate if the linestring only consists of two coordinates. For details - * on the bug inside [[LineString]] that is addressed here, see - * https://github.com/locationtech/jts/issues/531 - * - * @param lineString - * the linestring that should be checked and maybe converted to a 'safe to - * be compared' linestring - * @return - * a 'safe to be compared' linestring - */ + * compared' [[LineString]] This is done by removing duplicates in the points + * in the provided linestring as well as a small change of the start + * coordinate if the linestring only consists of two coordinates. For details + * on the bug inside [[LineString]] that is addressed here, see + * https://github.com/locationtech/jts/issues/531 + * + * @param lineString + * the linestring that should be checked and maybe converted to a 'safe to + * be compared' linestring + * @return + * a 'safe to be compared' linestring + */ def buildSafeLineString(lineString: LineString): LineString = if (lineString.getCoordinates.length == 2) buildSafeLineStringBetweenPoints( @@ -145,7 +143,6 @@ object GeoUtils { new Coordinate(p1X, p1Y, p1Z) } - /** Calculates and orders the coordinate distances from a base coordinate to a * list of coordinates * @@ -185,7 +182,7 @@ object GeoUtils { latB: Double, longB: Double ): ComparableQuantity[Length] = { - val r = EARTH_RADIUS.to(KILOMETRE) + val r = EARTH_RADIUS val dLat = Quantities.getQuantity(toRadians(latB - latA), RADIAN) val dLon = Quantities.getQuantity(toRadians(longB - longA), RADIAN) val a = pow(sin(dLat.getValue.doubleValue / 2), 2) + cos( @@ -195,20 +192,25 @@ object GeoUtils { r.multiply(c) } - // Fixme: Might be slightly off as it does not account for earth's curvature - /** Builds a convex hull from a set of coordinates. + /** Builds a convex hull from a set of latitude/longitude coordinates. + * Projects them on a 2D-surface, calculates the convex hull and reverses the + * projection of resulting coordinates to latitude and longitude values. * * @param coordinates * the coordinates to consider * @return * a Try of the resulting polygon */ - def buildConvexHull(coordinates: Set[Coordinate]): Try[Polygon] = { + def buildConvexHull(coordinates: List[Coordinate]): Try[Polygon] = { + val projectedCoordinates = coordinates.map(equalAreaProjection) new ConvexHull( - coordinates.toArray, + projectedCoordinates.toArray, DEFAULT_GEOMETRY_FACTORY ).getConvexHull match { - case polygon: Polygon => Success(polygon) + case projectedPolygon: Polygon => + val coordinates = + projectedPolygon.getCoordinates.map(reverseEqualAreaProjection) + Success(buildPolygon(coordinates.toList)) case _: LineString => Failure( new GeoException( @@ -231,6 +233,14 @@ object GeoUtils { } } + def buildPoint(lat: Double, long: Double): Point = { + buildPoint(buildCoordinate(lat, long)) + } + + def buildPoint(coordinate: Coordinate): Point = { + DEFAULT_GEOMETRY_FACTORY.createPoint(coordinate) + } + /** Build a coordinate from latitude and longitude values. * * @param lat @@ -244,7 +254,9 @@ object GeoUtils { new Coordinate(long, lat) } - /** Builds a polygon from a List of coordinates + /** Builds a polygon from a List of coordinates. To build a Polygon the + * coordinates have to form a closed ring which means that the first and last + * coordinate have to be the same coordinate. * * @param coordinates * the coordinates for building the polygon @@ -258,10 +270,11 @@ object GeoUtils { new Polygon(linearRing, Array[LinearRing](), DEFAULT_GEOMETRY_FACTORY) } - /** Credits to Joe Kington + /** Sinusoidal equal-area projection of latitude and longitude values to a 2d + * surface. This is an approximation gets worse when lines become very long. + * + * Credits to Joe Kington * (https://stackoverflow.com/questions/4681737/how-to-calculate-the-area-of-a-polygon-on-the-earths-surface-using-python#:~:text=Basically%2C%20you%20just%20multiply%20the,the%20cosine%20of%20the%20latitude.) - * Sinusoidal equal-area projection of latitude and longitude values to a 2d - * surface. This is an approximation gets worse when lines become very long * * @param lat * latitude to project @@ -270,18 +283,38 @@ object GeoUtils { * @return * a projected Coordinate with values in metre */ - def equalAreaProjection(lat: Double, long: Double): Coordinate = { + def equalAreaProjection(coordinate: Coordinate): Coordinate = { + val lat = coordinate.getY + val long = coordinate.getX val latDist = - Math.PI * (EARTH_RADIUS.to(METRE).getValue.doubleValue() / 180.0) + PI * (EARTH_RADIUS.to(METRE).getValue.doubleValue() / 180d) val y = lat * latDist val x = long * latDist * cos(toRadians(lat)) new Coordinate(x, y) } + /** Reverses the [[equalAreaProjection()]] and returns a coordinate on earths + * surface + * + * @param coordinate + * the projected coordinate + * @return + * the latitude longitude based coordinate + */ + def reverseEqualAreaProjection(coordinate: Coordinate): Coordinate = { + val latDist = PI * (EARTH_RADIUS.to(METRE).getValue.doubleValue() / 180d) + val lat = coordinate.y / latDist + val long = coordinate.x / (latDist * cos(toRadians(lat))) + buildCoordinate(lat, long) + } + /** Draws a circle with a radius of the provided distance around the provided * center coordinates and returns the result as a drawable polygon (one point * per degree) * + * Source: https://www.movable-type.co.uk/scripts/latlong.html "Destination + * point given distance and bearing from start point" + * * @param center * coordinate of the circle's center * @param radius @@ -293,12 +326,19 @@ object GeoUtils { center: Coordinate, radius: Quantity[Length] ): Polygon = { - val coordinates = (0 to 359) - .map(deg => { - val angle = deg.toFloat.toRadians - val x = center.x + radius.getValue.doubleValue() * cos(angle) - val y = center.y + radius.getValue.doubleValue() * sin(angle) - new Coordinate(x, y) + val centerLat = toRadians(center.y) + val centerLon = toRadians(center.x) + val d = radius.divide(EARTH_RADIUS).getValue.doubleValue + val coordinates = (0 to 360) + .map(angle => { + val bearing = toRadians(angle) + val latRad = + asin(sin(centerLat) * cos(d) + cos(centerLat) * sin(d) * cos(bearing)) + val lonRad: Double = centerLon + atan2( + sin(bearing) * sin(d) * cos(centerLat), + cos(d) - sin(centerLat) * sin(latRad) + ) + new Coordinate(lonRad.toDegrees, latRad.toDegrees) }) .toList buildPolygon(coordinates) diff --git a/src/main/scala/edu/ie3/util/geo/RichGeometries.scala b/src/main/scala/edu/ie3/util/geo/RichGeometries.scala index a374662f..0d74539f 100644 --- a/src/main/scala/edu/ie3/util/geo/RichGeometries.scala +++ b/src/main/scala/edu/ie3/util/geo/RichGeometries.scala @@ -6,7 +6,11 @@ package edu.ie3.util.geo import edu.ie3.util.exceptions.GeoException -import edu.ie3.util.geo.GeoUtils.{DEFAULT_GEOMETRY_FACTORY, buildPolygon, calcHaversine} +import edu.ie3.util.geo.GeoUtils.{ + DEFAULT_GEOMETRY_FACTORY, + buildPolygon, + calcHaversine +} import edu.ie3.util.quantities.PowerSystemUnits.KILOMETRE import edu.ie3.util.quantities.QuantityUtils.RichQuantityDouble import org.locationtech.jts.geom.{Coordinate, LineString, Point, Polygon} @@ -61,7 +65,9 @@ object RichGeometries { .add(coordinate.haversineDistance(b)) .getValue .doubleValue - if (abs(1 - (distancePassingMe / distance.getValue.doubleValue())) < epsilon) true + if ( + abs(1 - (distancePassingMe / distance.getValue.doubleValue())) < epsilon + ) true else false } @@ -126,7 +132,7 @@ object RichGeometries { */ def equalAreaProjection: Polygon = { val projectedCoordinates = polygon.getCoordinates.map(coordinate => - GeoUtils.equalAreaProjection(coordinate.y, coordinate.getX) + GeoUtils.equalAreaProjection(coordinate) ) buildPolygon(projectedCoordinates.toList) } diff --git a/src/test/groovy/edu/ie3/util/geo/GeoUtilsTest.groovy b/src/test/groovy/edu/ie3/util/geo/GeoUtilsTest.groovy index 7bc20d01..af8551dc 100644 --- a/src/test/groovy/edu/ie3/util/geo/GeoUtilsTest.groovy +++ b/src/test/groovy/edu/ie3/util/geo/GeoUtilsTest.groovy @@ -45,7 +45,7 @@ class GeoUtilsTest extends Specification { ComparableQuantity actual = DeprecatedGeoUtils.calcHaversine(start.lat, start.lon, end.lat, end.lon) then: - Math.abs(actual.subtract(expected).to(METRE).value.doubleValue()) < tolerance.value.doubleValue() + Math.abs(actual.subtract(expected).to(METRE).value.doubleValue()) < tolerance .value.doubleValue() } def "Total length of LineString is correctly calculated"() { @@ -96,7 +96,7 @@ class GeoUtilsTest extends Specification { for (int cnt = 0; cnt < coordinates.length; cnt++) { coordinates[cnt] == actualCoordinates[cnt] } -from + from where: lineString | coordinates "{ \"type\": \"LineString\", \"coordinates\": [[7.411111, 51.49228],[7.411111, 51.49228]]}" | [ diff --git a/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala b/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala index 8a3bc477..08faa5bf 100644 --- a/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala +++ b/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala @@ -1,23 +1,192 @@ +/* + * © 2021. TU Dortmund University, + * Institute of Energy Systems, Energy Efficiency and Energy Economics, + * Research group Distribution grid planning and operation +*/ package edu.ie3.util.geo -import edu.ie3.util.geo.GeoUtils.DEFAULT_GEOMETRY_FACTORY +import edu.ie3.util.geo.GeoUtils.{ + DEFAULT_GEOMETRY_FACTORY, + buildCirclePolygon, + buildConvexHull, + buildPoint, + buildPolygon, + calcHaversine, + calcOrderedCoordinateDistances, + equalAreaProjection, + reverseEqualAreaProjection +} import edu.ie3.util.geo.RichGeometries.RichCoordinate +import edu.ie3.util.quantities.QuantityMatchers.equalWithTolerance import org.locationtech.jts.geom.Coordinate import org.scalatest.matchers.should.Matchers import org.scalatest.wordspec.AnyWordSpecLike +import tech.units.indriya.quantity.Quantities +import tech.units.indriya.unit.Units.METRE + +import scala.math.abs +import scala.util.{Failure, Success} -class GeoUtilsSpec extends Matchers with AnyWordSpecLike{ +class GeoUtilsSpec extends Matchers with AnyWordSpecLike { "The GeoUtils" should { "build a line string with two equal geo coordinates, avoiding the known bug in jts geometry" in { val coordinate = new Coordinate(7.411111, 51.49228) - val coordinateWithCorrection = new Coordinate(7.4111110000001, 51.4922800000001) - val fromCoords = GeoUtils.buildSafeLineStringBetweenCoords(coordinate, coordinate).getCoordinates - val fromPoints = GeoUtils.buildSafeLineStringBetweenPoints(coordinate.toPoint, coordinate.toPoint).getCoordinates - val fromLineString = GeoUtils.buildSafeLineString(DEFAULT_GEOMETRY_FACTORY.createLineString(Array[Coordinate](coordinate, coordinate))) + val coordinateWithCorrection = + new Coordinate(7.4111110000001, 51.4922800000001) + val fromCoords = GeoUtils + .buildSafeLineStringBetweenCoords(coordinate, coordinate) + .getCoordinates + val fromPoints = GeoUtils + .buildSafeLineStringBetweenPoints( + coordinate.toPoint, + coordinate.toPoint + ) + .getCoordinates + val fromLineString = GeoUtils + .buildSafeLineString( + DEFAULT_GEOMETRY_FACTORY.createLineString( + Array[Coordinate](coordinate, coordinate) + ) + ) + .getCoordinates fromCoords(0) == coordinateWithCorrection fromCoords(1) == coordinate + fromPoints(0) == coordinateWithCorrection + fromPoints(1) == coordinate + fromLineString(0) == coordinateWithCorrection + fromLineString(1) == coordinate + } + + "strip a line string to its distinct elements when building a safe one" in { + val coordinateA = new Coordinate(7.411111, 51.49228) + val coordinateB = new Coordinate(7.411111, 52.49228) + val safeLineString = GeoUtils.buildSafeLineString( + DEFAULT_GEOMETRY_FACTORY.createLineString( + Array[Coordinate](coordinateA, coordinateA, coordinateB, coordinateB) + ) + ) + val strippedCoordinates = safeLineString.getCoordinates + strippedCoordinates.length shouldBe 2 + strippedCoordinates(0) == coordinateA + strippedCoordinates(1) == coordinateB + } + + "calculate ordered coordinate distances correctly" in { + val base = buildPoint(51.53651260695586, 9.934043825403268) + val nearest = buildPoint(51.553166628160746, 9.968032777276226) + val middle = buildPoint(51.44074066208236, 10.039100585737868) + val furthest = buildPoint(51.6273949370414, 10.341911247878777) + val actual = + calcOrderedCoordinateDistances(base, List(middle, furthest, nearest)) + val list = actual.toList + list.size shouldBe 3 + list.head.to shouldBe nearest + list(1).to shouldBe middle + list(2).to shouldBe furthest } + + "calculate the haversine distance correctly" in { + val latA = 37.87532764735112 + val longA = -122.25311279296875 + val latB = 37.87934174490509 + val longB = -122.2537350654602 + val expected = Quantities.getQuantity(450.18011568984845, METRE) + val actual = calcHaversine(latA, longA, latB, longB) + actual should equalWithTolerance(expected) + } + + "build a convex hull correctly" in { + val topLeft = new Coordinate(7, 50) + val betweenTlTr = new Coordinate(7.5, 50) + val topRight = new Coordinate(8, 50) + val bottomRight = new Coordinate(8, 48) + val betweenBlBr = new Coordinate(7.5, 48) + val bottomLeft = new Coordinate(7, 48) + val coordinates = List( + topLeft, + betweenTlTr, + topRight, + bottomRight, + betweenBlBr, + bottomLeft + ) + val tryConvexHull = buildConvexHull(coordinates) + tryConvexHull match { + case Failure(exception) => throw exception + case Success(polygon) => + val hullCoordinates = polygon.getCoordinates + // contains all corner coordinates and filters out all colinear ones + hullCoordinates.size shouldBe 5 + hullCoordinates(0).equals2D(topLeft, 1e-10) + hullCoordinates(1).equals(topRight, 1e-10) shouldBe true + hullCoordinates(2).equals(bottomRight, 1e-10) shouldBe true + hullCoordinates(3).equals(bottomLeft, 1e-10) shouldBe true + hullCoordinates(4).equals(topLeft, 1e-10) shouldBe true + } + } + + "build a coordinate correctly" in { + val lat = 51.49860455457335 + val long = 7.468448342940863 + val actual = GeoUtils.buildCoordinate(lat, long) + actual.getX shouldBe long + actual.getY shouldBe lat + } + + "build a polygon correctly" in { + val coordinateA = new Coordinate(7.468448342940863, 51.49860455457335) + val coordinateB = new Coordinate(7.521007845835815, 51.50450661471354) + val coordinateC = new Coordinate(7.5598606548385545, 51.456498140367934) + val actual = + buildPolygon(List(coordinateA, coordinateB, coordinateC, coordinateA)) + actual.getCoordinates shouldBe Array( + coordinateA, + coordinateB, + coordinateC, + coordinateA + ) + } + + "do an equal area projection of a coordinate correctly" in { + val actual = equalAreaProjection( + new Coordinate(8.748497631269068, 51.72341137638795) + ) + actual.getX shouldBe 603277.0126920443 +- 1e-12 + actual.getY shouldBe 5757823.816510521 +- 1e-12 + } + + "reverse the equal area projection correctly" in { + val original = new Coordinate(7.468448342940863, 51.49860455457335) + val reversed = reverseEqualAreaProjection(equalAreaProjection(original)) + original.getX shouldBe reversed.getX +- 1e-12 + original.getY shouldBe reversed.getY +- 1e-12 + } + + "build a circle polygon correctly" in { + val center = new Coordinate(51.49860455457335, 7.468448342940863) + val radius = Quantities.getQuantity(50, METRE) + val actualCoordinates = buildCirclePolygon(center, radius).getCoordinates + implicit val precision: Double = 1e-6 + actualCoordinates.length shouldBe 361 + abs(actualCoordinates(90).y - center.y) shouldBe 0.0 +- precision + actualCoordinates(90).haversineDistance(center) should equalWithTolerance( + radius + ) + abs(actualCoordinates(180).x - center.x) shouldBe 0.0 +- precision + actualCoordinates(180).haversineDistance( + center + ) should equalWithTolerance(radius) + abs(actualCoordinates(270).y - center.y) shouldBe 0.0 +- precision + actualCoordinates(270).haversineDistance( + center + ) should equalWithTolerance(radius) + abs(actualCoordinates(360).x - center.x) shouldBe 0.0 +- precision + actualCoordinates(360).haversineDistance( + center + ) should equalWithTolerance(radius) + } + } } diff --git a/src/test/scala/edu/ie3/util/geo/RichGeometriesSpec.scala b/src/test/scala/edu/ie3/util/geo/RichGeometriesSpec.scala new file mode 100644 index 00000000..0cd77bc6 --- /dev/null +++ b/src/test/scala/edu/ie3/util/geo/RichGeometriesSpec.scala @@ -0,0 +1,36 @@ +/* + * © 2021. TU Dortmund University, + * Institute of Energy Systems, Energy Efficiency and Energy Economics, + * Research group Distribution grid planning and operation +*/ +package edu.ie3.util.geo + +import edu.ie3.util.geo.GeoUtils.DEFAULT_GEOMETRY_FACTORY +import edu.ie3.util.geo.RichGeometries.RichPolygon +import edu.ie3.util.quantities.QuantityUtils.RichQuantityDouble +import org.locationtech.jts.geom.Coordinate +import org.scalatest.matchers.should.Matchers +import org.scalatest.wordspec.AnyWordSpecLike + +class RichGeometriesSpec extends Matchers with AnyWordSpecLike { + + "A rich polygon" should { + + "calculate area on earth correctly" in { + val coordinateA = new Coordinate(8.748497631269068, 51.72341137638795) + val coordinateB = new Coordinate(8.76167264195022, 51.723225286136866) + val coordinateC = new Coordinate(8.76240220280227, 51.715568337479546) + val coordinateD = new Coordinate(8.74832596989211, 51.71546198184121) + val polygon = DEFAULT_GEOMETRY_FACTORY.createPolygon( + Array(coordinateA, coordinateB, coordinateC, coordinateD, coordinateA) + ) + val actual = 813431.49.asSquareMetre + polygon.calcAreaOnEarth + .divide(actual) + .getValue + .doubleValue() shouldBe 1d +- 0.01 + } + + } + +} From a58dd73b00dbfa7e5dc1d035d80c2a0d996017ef Mon Sep 17 00:00:00 2001 From: t-ober <63147366+t-ober@users.noreply.github.com> Date: Wed, 8 Dec 2021 14:50:49 +0100 Subject: [PATCH 07/17] rename RichGeometries --- .../edu/ie3/util/geo/RichGeometries.scala | 28 ++++++++++++------- .../scala/edu/ie3/util/osm/OsmEntities.scala | 2 +- .../scala/edu/ie3/util/geo/GeoUtilsSpec.scala | 2 +- .../edu/ie3/util/geo/RichGeometriesSpec.scala | 2 +- 4 files changed, 21 insertions(+), 13 deletions(-) diff --git a/src/main/scala/edu/ie3/util/geo/RichGeometries.scala b/src/main/scala/edu/ie3/util/geo/RichGeometries.scala index 0d74539f..99574ad5 100644 --- a/src/main/scala/edu/ie3/util/geo/RichGeometries.scala +++ b/src/main/scala/edu/ie3/util/geo/RichGeometries.scala @@ -9,7 +9,8 @@ import edu.ie3.util.exceptions.GeoException import edu.ie3.util.geo.GeoUtils.{ DEFAULT_GEOMETRY_FACTORY, buildPolygon, - calcHaversine + calcHaversine, + reverseEqualAreaProjection } import edu.ie3.util.quantities.PowerSystemUnits.KILOMETRE import edu.ie3.util.quantities.QuantityUtils.RichQuantityDouble @@ -23,7 +24,7 @@ import scala.util.{Failure, Success, Try} object RichGeometries { - implicit class RichCoordinate(coordinate: Coordinate) { + implicit class GeoCoordinate(coordinate: Coordinate) { /** Calculates the great circle distance between two coordinates * @@ -76,7 +77,7 @@ object RichGeometries { } } - implicit class RichLineString(lineString: LineString) { + implicit class GeoLineString(lineString: LineString) { /** Compute length of a [[LineString]] on earth's surface. * @@ -97,12 +98,19 @@ object RichGeometries { } } - implicit class RichPolygon(polygon: Polygon) { + implicit class GeoPolygon(polygon: Polygon) { - // Fixme: Is not precisely right because earths curvature - use equal area projection ? + // Fixme: Is this correct or should the function calculate the specific coordinates shared by the polygons + /** Calculates intersection between polygons + * + * @param polygonB + * polygon with which to calculate the intersection + * @return + */ def intersect(polygonB: Polygon): Try[Set[Coordinate]] = { - - Try(polygon.intersection(polygonB)) match { + Try( + polygon.equalAreaProjection.intersection(polygonB.equalAreaProjection) + ) match { case Failure(exception) => Failure( new GeoException( @@ -110,7 +118,8 @@ object RichGeometries { exception ) ) - case Success(geometry) => Success(geometry.getCoordinates.toSet) + case Success(geometry) => + Success(geometry.getCoordinates.map(reverseEqualAreaProjection).toSet) } } @@ -137,8 +146,7 @@ object RichGeometries { buildPolygon(projectedCoordinates.toList) } - // Fixme: Is not precisely right because earths curvatur - + // Fixme: Should we do an equal area projection here as well? /** Checks whether the polygon contains the coordinate. Uses "covers()" * insted of "contains()" so borders are included. * diff --git a/src/main/scala/edu/ie3/util/osm/OsmEntities.scala b/src/main/scala/edu/ie3/util/osm/OsmEntities.scala index 80b6ccf4..41e3ce84 100644 --- a/src/main/scala/edu/ie3/util/osm/OsmEntities.scala +++ b/src/main/scala/edu/ie3/util/osm/OsmEntities.scala @@ -6,7 +6,7 @@ package edu.ie3.util.osm import edu.ie3.util.geo.GeoUtils.buildPolygon -import edu.ie3.util.geo.RichGeometries.RichPolygon +import edu.ie3.util.geo.RichGeometries.GeoPolygon import org.locationtech.jts.geom.{Coordinate, Polygon} import tech.units.indriya.ComparableQuantity diff --git a/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala b/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala index 08faa5bf..2c7e6669 100644 --- a/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala +++ b/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala @@ -16,7 +16,7 @@ import edu.ie3.util.geo.GeoUtils.{ equalAreaProjection, reverseEqualAreaProjection } -import edu.ie3.util.geo.RichGeometries.RichCoordinate +import edu.ie3.util.geo.RichGeometries.GeoCoordinate import edu.ie3.util.quantities.QuantityMatchers.equalWithTolerance import org.locationtech.jts.geom.Coordinate import org.scalatest.matchers.should.Matchers diff --git a/src/test/scala/edu/ie3/util/geo/RichGeometriesSpec.scala b/src/test/scala/edu/ie3/util/geo/RichGeometriesSpec.scala index 0cd77bc6..a3c9a103 100644 --- a/src/test/scala/edu/ie3/util/geo/RichGeometriesSpec.scala +++ b/src/test/scala/edu/ie3/util/geo/RichGeometriesSpec.scala @@ -6,7 +6,7 @@ package edu.ie3.util.geo import edu.ie3.util.geo.GeoUtils.DEFAULT_GEOMETRY_FACTORY -import edu.ie3.util.geo.RichGeometries.RichPolygon +import edu.ie3.util.geo.RichGeometries.GeoPolygon import edu.ie3.util.quantities.QuantityUtils.RichQuantityDouble import org.locationtech.jts.geom.Coordinate import org.scalatest.matchers.should.Matchers From 35dfae32abf0610b64865bdf740ed7166fd30f76 Mon Sep 17 00:00:00 2001 From: t-ober <63147366+t-ober@users.noreply.github.com> Date: Wed, 8 Dec 2021 14:57:04 +0100 Subject: [PATCH 08/17] remove old java code --- .../geo/DeprecatedCoordinateDistance.java | 105 -- .../edu/ie3/util/geo/DeprecatedGeoUtils.java | 1034 ----------------- .../DeprecatedCoordinateDistanceTest.groovy | 42 - .../edu/ie3/util/geo/GeoUtilsTest.groovy | 191 --- 4 files changed, 1372 deletions(-) delete mode 100644 src/main/java/edu/ie3/util/geo/DeprecatedCoordinateDistance.java delete mode 100644 src/main/java/edu/ie3/util/geo/DeprecatedGeoUtils.java delete mode 100644 src/test/groovy/edu/ie3/util/geo/DeprecatedCoordinateDistanceTest.groovy delete mode 100644 src/test/groovy/edu/ie3/util/geo/GeoUtilsTest.groovy diff --git a/src/main/java/edu/ie3/util/geo/DeprecatedCoordinateDistance.java b/src/main/java/edu/ie3/util/geo/DeprecatedCoordinateDistance.java deleted file mode 100644 index ff671b4c..00000000 --- a/src/main/java/edu/ie3/util/geo/DeprecatedCoordinateDistance.java +++ /dev/null @@ -1,105 +0,0 @@ -/* - * © 2021. TU Dortmund University, - * Institute of Energy Systems, Energy Efficiency and Energy Economics, - * Research group Distribution grid planning and operation -*/ -package edu.ie3.util.geo; - -import java.util.Objects; -import javax.measure.quantity.Length; -import tech.units.indriya.ComparableQuantity; - -/** - * Wraps two coordinates with the distance between the first one and the second one, can be compared - * by distance to another CoordinateDistance - */ -public class DeprecatedCoordinateDistance implements Comparable { - private final org.locationtech.jts.geom.Point coordinateA; - private final org.locationtech.jts.geom.Point coordinateB; - private final ComparableQuantity distance; - - /** - * Calculates the distance from the first to the second coordinate using {@link - * DeprecatedGeoUtils#calcHaversine(double, double, double, double)} - * - * @param coordinateA The first coordinate - * @param coordinateB The second coordinate - */ - public DeprecatedCoordinateDistance( - org.locationtech.jts.geom.Point coordinateA, org.locationtech.jts.geom.Point coordinateB) { - this( - coordinateA, - coordinateB, - DeprecatedGeoUtils.calcHaversine( - coordinateA.getY(), coordinateA.getX(), coordinateB.getY(), coordinateB.getX())); - } - - /** - * @param coordinateA The first coordinate - * @param coordinateB The second coordinate - * @param distance The distance from A to B - */ - private DeprecatedCoordinateDistance( - org.locationtech.jts.geom.Point coordinateA, - org.locationtech.jts.geom.Point coordinateB, - ComparableQuantity distance) { - this.coordinateA = coordinateA; - this.coordinateB = coordinateB; - this.distance = distance; - } - - /** @return The first coordinate */ - public org.locationtech.jts.geom.Point getCoordinateA() { - return coordinateA; - } - - /** @return The second coordinate */ - public org.locationtech.jts.geom.Point getCoordinateB() { - return coordinateB; - } - - /** @return The distance from the first coordinate to the second coordinate in km */ - public ComparableQuantity getDistance() { - return distance; - } - - /** - * Compares two coordinate distances on the length of the distance alone, thus having a natural - * ordering that is inconsistent with equals - * - * @param that the distance to compare - * @return a number lower than 0 if this has a lower distance than that, 0 if they are the same, a - * number higher than 0 if that has a lower distance - */ - @Override - public int compareTo(DeprecatedCoordinateDistance that) { - return this.distance.compareTo(that.distance); - } - - @Override - public boolean equals(Object o) { - if (this == o) return true; - if (o == null || getClass() != o.getClass()) return false; - DeprecatedCoordinateDistance that = (DeprecatedCoordinateDistance) o; - return coordinateA.equals(that.coordinateA) - && coordinateB.equals(that.coordinateB) - && distance.equals(that.distance); - } - - @Override - public int hashCode() { - return Objects.hash(coordinateA, coordinateB, distance); - } - - @Override - public String toString() { - return "CoordinateDistance{" - + "coordinateA=" - + coordinateA - + ", coordinateB=" - + coordinateB - + ", distance=" - + distance - + '}'; - } -} diff --git a/src/main/java/edu/ie3/util/geo/DeprecatedGeoUtils.java b/src/main/java/edu/ie3/util/geo/DeprecatedGeoUtils.java deleted file mode 100644 index a136aae2..00000000 --- a/src/main/java/edu/ie3/util/geo/DeprecatedGeoUtils.java +++ /dev/null @@ -1,1034 +0,0 @@ -/* - * © 2020. TU Dortmund University, - * Institute of Energy Systems, Energy Efficiency and Energy Economics, - * Research group Distribution grid planning and operation -*/ -package edu.ie3.util.geo; - -import static edu.ie3.util.quantities.PowerSystemUnits.*; - -import com.google.common.collect.Lists; -import edu.ie3.util.copy.DeepCopy; -import edu.ie3.util.exceptions.GeoPreparationException; -import java.awt.*; -import java.util.*; -import java.util.List; -import java.util.stream.Collectors; -import javax.measure.Quantity; -import javax.measure.quantity.Angle; -import javax.measure.quantity.Area; -import javax.measure.quantity.Length; -import net.morbz.osmonaut.geometry.Bounds; -import net.morbz.osmonaut.geometry.Polygon; -import net.morbz.osmonaut.osm.*; -import org.apache.commons.lang3.ArrayUtils; -import org.locationtech.jts.geom.Coordinate; -import org.locationtech.jts.geom.GeometryFactory; -import org.locationtech.jts.geom.LineString; -import org.locationtech.jts.geom.PrecisionModel; -import org.locationtech.jts.geom.impl.CoordinateArraySequence; -import org.slf4j.Logger; -import org.slf4j.LoggerFactory; -import tech.units.indriya.ComparableQuantity; -import tech.units.indriya.quantity.Quantities; - -/** Functionality to deal with geographical and geometric information */ -public class DeprecatedGeoUtils { - private static final Logger logger = LoggerFactory.getLogger(DeprecatedGeoUtils.class); - - public static final ComparableQuantity EARTH_RADIUS = - Quantities.getQuantity(6378137.0, METRE); - - /** Offer a default geometry factory for the WGS84 coordinate system */ - public static final GeometryFactory DEFAULT_GEOMETRY_FACTORY = - new GeometryFactory(new PrecisionModel(), 4326); - - protected DeprecatedGeoUtils() { - throw new IllegalStateException("Utility classes cannot be instantiated"); - } - - /** - * Calculates the distance in km between two lat/long points using the haversine formula - * - * @param lat1 Latitude value of the first coordinate - * @param lng1 Longitude value of the first coordinate - * @param lat2 Latitude value of the second coordinate - * @param lng2 Longitude value of the second coordinate - * @return The distance between both coordinates in {@link - * edu.ie3.util.quantities.PowerSystemUnits#KILOMETRE} - */ - public static ComparableQuantity calcHaversine( - double lat1, double lng1, double lat2, double lng2) { - - // average radius of the earth in km - ComparableQuantity r = EARTH_RADIUS.to(KILOMETRE); - ComparableQuantity dLat = Quantities.getQuantity(Math.toRadians(lat2 - lat1), RADIAN); - ComparableQuantity dLon = Quantities.getQuantity(Math.toRadians(lng2 - lng1), RADIAN); - double a = - Math.sin(dLat.getValue().doubleValue() / 2) * Math.sin(dLat.getValue().doubleValue() / 2) - + Math.cos(Math.toRadians(lat1)) - * Math.cos(Math.toRadians(lat2)) - * Math.sin(dLon.getValue().doubleValue() / 2) - * Math.sin(dLon.getValue().doubleValue() / 2); - double c = 2 * Math.atan2(Math.sqrt(a), Math.sqrt(1 - a)); - return r.multiply(c); - } - - /** - * Calculates the total length of a LineString through building the sum of the distances between - * all points of LineString using {@link #calcHaversine(double, double, double, double)} - * - * @param lineString the line string which length shall be calculated - * @return The total length of the line string - */ - public static ComparableQuantity totalLengthOfLineString(LineString lineString) { - ComparableQuantity y = Quantities.getQuantity(0, KILOMETRE); - for (int i = 0; i < lineString.getNumPoints() - 1; i++) { - y = - y.add( - calcHaversine( - lineString.getCoordinateN(i).getX(), - lineString.getCoordinateN(i).getY(), - lineString.getCoordinateN(i + 1).getX(), - lineString.getCoordinateN(i + 1).getY())); - } - return y; - } - - /** - * Calculates and sorts the distances between a base coordinate and other given coordinates using - * {@link #calcHaversine(double, double, double, double)} - * - * @param baseCoordinate the base point - * @param coordinates the points to calculate the distance from the base point for - * @return a sorted set of distances between the base and other coordinates - */ - public static SortedSet getCoordinateDistances( - org.locationtech.jts.geom.Point baseCoordinate, - Collection coordinates) { - return coordinates.stream() - .map(coordinate -> new DeprecatedCoordinateDistance(baseCoordinate, coordinate)) - .collect(Collectors.toCollection(TreeSet::new)); - } - - /** - * Build an instance of {@link LineString} between two points that is safe to be compared even if - * the provided two points consist of exactly the same coordinates. This is done by increasing the - * coordinate of the provided Point {@code p1} by a small amount to make it different from Point - * {@code p2}. For details on the bug inside {@link LineString} that is addressed here, see - * https://github.com/locationtech/jts/issues/531 - * - * @param p1 start point of the linestring - * @param p2 end point of the linestring - * @return a {@link LineString} between the provided points - */ - public static LineString buildSafeLineStringBetweenPoints( - final org.locationtech.jts.geom.Point p1, final org.locationtech.jts.geom.Point p2) { - final org.locationtech.jts.geom.Point safePoint1 = p1.equals(p2) ? buildSafePoint(p1) : p1; - return DEFAULT_GEOMETRY_FACTORY.createLineString( - ArrayUtils.addAll(safePoint1.getCoordinates(), p2.getCoordinates())); - } - - /** - * Adapt the provided point as described in {@link #buildSafeCoord(Coordinate)} and return a new, - * adapted instance of {@link Point} - * - * @param p1 the point that should be adapted - * @return the adapted point with a slightly changed coordinate - */ - private static org.locationtech.jts.geom.Point buildSafePoint( - org.locationtech.jts.geom.Point p1) { - - Coordinate[] safeCoord = new Coordinate[] {buildSafeCoord(p1.getCoordinate())}; - CoordinateArraySequence safeCoordSeq = new CoordinateArraySequence(safeCoord); - - return new org.locationtech.jts.geom.Point(safeCoordSeq, p1.getFactory()); - } - - /** - * Adapted {@link Coordinate#x}, {@link Coordinate#y} and {@link Coordinate#z} of the provided - * {@link Coordinate} by 1e-13 and return a new, adapted instance of {@link Coordinate} - * - * @param coord the coordinate that should be adapted - * @return the adapted coordinate with slightly changed x,y,z values - */ - private static Coordinate buildSafeCoord(Coordinate coord) { - - double modVal = 1e-13; - double p1X = coord.getX() + modVal; - double p1Y = coord.getY() + modVal; - double p1Z = coord.getZ() + modVal; - - return new Coordinate(p1X, p1Y, p1Z); - } - - /** - * Build an instance of {@link LineString} between two coordinates that is safe to be compared - * even if the provided two coordinates are exactly the same coordinates. This is done by - * increasing the coordinate of the provided Point {@code c1} by a small amount to make it - * different from Point {@code c2}. For details on the bug inside {@link LineString} that is - * addressed here, see https://github.com/locationtech/jts/issues/531 - * - * @param c1 start coordinate of the linestring - * @param c2 end coordinate of the linestring - * @return A safely build line string - */ - public static LineString buildSafeLineStringBetweenCoords( - final Coordinate c1, final Coordinate c2) { - final Coordinate safeCoord1 = c1.equals(c2) ? buildSafeCoord(c1) : c1; - return DEFAULT_GEOMETRY_FACTORY.createLineString( - ArrayUtils.addAll(new Coordinate[] {safeCoord1}, c2)); - } - - /** - * Convert a given {@link LineString} with at least two points into a 'safe to be compared' {@link - * LineString} This is done by removing duplicates in the points in the provided linestring as - * well as a small change of the start coordinate if the linestring only consists of two - * coordinates. For details on the bug inside {@link LineString} that is addressed here, see - * https://github.com/locationtech/jts/issues/531 - * - * @param lineString the linestring that should be checked and maybe converted to a 'safe to be - * compared' linestring - * @return a 'safe to be compared' linestring - */ - public static LineString buildSafeLineString(LineString lineString) { - if (lineString.getCoordinates().length == 2) { - return buildSafeLineStringBetweenPoints(lineString.getStartPoint(), lineString.getEndPoint()); - } else { - // rebuild line with unique points - /* Please note, that using a simple HashSet here to obtain uniqueness is harmful, as it not necessarily maintains - * the order of coordinates (but most likely will). Additionally, the behaviour of a HashSet might change with the - * JVM (version) you use to execute the code. */ - Coordinate[] uniqueCoords = - Arrays.stream(lineString.getCoordinates()).distinct().toArray(Coordinate[]::new); - return uniqueCoords.length == 1 - ? buildSafeLineStringBetweenPoints(lineString.getStartPoint(), lineString.getEndPoint()) - : DEFAULT_GEOMETRY_FACTORY.createLineString(uniqueCoords); - } - } - - /** - * Calculates the geo position as a Point from a given Latlon (net.morbz.osmonaut.osm). Uses the - * WGS84 reference system. - * - * @param latLon Latlon from which the geo position shall be calculated - * @return calculated Point from the given Latlon - */ - public static org.locationtech.jts.geom.Point latlonToPoint(LatLon latLon) { - return xyToPoint(latLon.getLon(), latLon.getLat()); - } - - /** - * Wraps XY values in a JTS geometry point - * - * @param x longitude value - * @param y latitude value - * @return JTS geometry Point - */ - public static org.locationtech.jts.geom.Point xyToPoint(double x, double y) { - Coordinate coordinate = new Coordinate(x, y, 0); - return DEFAULT_GEOMETRY_FACTORY.createPoint(coordinate); - } - - /** - * This method takes all {@link RelationMember}s and joins the given ways to have closed ways. - * This only works, if the ways are connected all in a line. - * - * @param relation {@link Relation} relation to treat - * @return Deep copy of the {@code relation} with closed ways - * @deprecated This method is currently not under test and has to be revised thoroughly - */ - @Deprecated - public static Relation buildClosedWays(Relation relation) throws GeoPreparationException { - /* Copy relation and empty the Members */ - Relation closedRelation = DeepCopy.copy(relation); - closedRelation - .getMembers() - .removeAll( - closedRelation.getMembers().stream() - .filter(e -> e.getEntity() instanceof Way) - .collect(Collectors.toList())); - - List ways = - relation.getMembers().stream() - .filter(e -> e.getEntity() instanceof Way) - .map(e -> (Way) e.getEntity()) - .collect(Collectors.toList()); - - /* Get an idea, of which ways do have intersections and where */ - HashMap> intersections = new HashMap<>(); - List comparativeWays = new LinkedList<>(ways); - for (Way way1 : ways) { - List nodes1 = way1.getNodes(); - - for (Way way2 : comparativeWays) { - if (way1.equals(way2)) continue; - List nodes2 = way2.getNodes(); - Set sharedNodes = - nodes1.stream().filter(nodes2::contains).collect(Collectors.toSet()); - for (Node node : sharedNodes) { - intersections.putIfAbsent(node, new HashSet<>()); - intersections.get(node).add(way1); - intersections.get(node).add(way2); - } - } - - // TODO: When there are no shared nodes, iterate through both ways once again and find - // matching nodes based on a radius search - comparativeWays.remove(way1); - } - - if (intersections.values().stream().anyMatch(s -> s.size() > 2)) - throw new GeoPreparationException( - "There is at least one intersection apparent with more then two ways, which would result in a concave hull curve. Exiting here."); - - /* As long as there are untreated ways iterate through them an build closed ways */ - while (!intersections.isEmpty()) { - /* Take the first intersection node and it's adjacent ways. Remove those ways, to mark them as travelled */ - Node firstIntersection = intersections.keySet().iterator().next(); - Way closedWay = null; - - do { - if (closedWay != null && intersections.get(firstIntersection).isEmpty()) { - logger.warn( - "There are no more intersections left, but the way is not closed right now. Adding the first node once again."); - closedWay.getNodes().add(closedWay.getNodes().get(0)); - continue; - } - Way currentWay = intersections.get(firstIntersection).iterator().next(); - intersections.get(firstIntersection).remove(currentWay); - - /* Find the next way in order to get the two intersections on the current way */ - List candidateNodes = - intersections.entrySet().stream() - .filter(e -> e.getValue().contains(currentWay)) - .map(Map.Entry::getKey) - .collect(Collectors.toList()); - Node secondIntersection; - if (candidateNodes.size() > 1) - throw new GeoPreparationException( - "There is an intersection, where more then two ways meet. This is not supported, yet."); - else if (candidateNodes.isEmpty()) { - /* There is no more way. Close the way by adding the start node ones again */ - if (closedWay == null) { - logger.warn("There is only one way in this relation."); - closedWay = currentWay; - } else if (!currentWay.getNodes().contains(closedWay.getNodes().get(0))) { - throw new GeoPreparationException("Ran into an dead end."); - } - secondIntersection = closedWay.getNodes().get(0); - } else { - secondIntersection = candidateNodes.get(0); - intersections.get(secondIntersection).remove(currentWay); - } - - int[] nodePos = - new int[] { - currentWay.getNodes().indexOf(firstIntersection), - currentWay.getNodes().indexOf(secondIntersection) - }; - LinkedList nodesToAdd; - if (nodePos[0] <= nodePos[1]) { - /* Next way can be added in ascending order */ - nodesToAdd = new LinkedList<>(currentWay.getNodes().subList(nodePos[0], nodePos[1] + 1)); - } else { - /* The intersection is at the end of the next way. Next way can be added in descending order */ - nodesToAdd = new LinkedList<>(currentWay.getNodes().subList(nodePos[1], nodePos[0] + 1)); - Collections.reverse(nodesToAdd); - } - - /* If there is no buffer way, yet. Create a new one based on the part way identified first */ - if (closedWay == null) - closedWay = new Way(currentWay.getId(), currentWay.getTags(), nodesToAdd); - else closedWay.getNodes().addAll(nodesToAdd); - - /* Go one step ahead */ - firstIntersection = secondIntersection; - } while (!closedWay.isClosed()); - - /* Remove all travelled intersections */ - intersections.entrySet().removeIf(e -> e.getValue().isEmpty()); - - /* Add the closed way to the relation */ - closedRelation.getMembers().add(new RelationMember(closedWay, "outer")); - } - - return closedRelation; - } - - /** - * This method is manly inspired by the following internet sources:
- * https://en.wikipedia.org/wiki/Convex_hull_algorithms
- * https://en.wikipedia.org/wiki/Chan%27s_algorithm - * - * @param points The {@link Set} of {@link Point}s that shall be enclosed by the convex hull - * @param precision Prescision to use for "dirty casting" - * @param algorithm {@link ConvexHullAlgorithm} to use. CAUTION {@link ConvexHullAlgorithm#CHAN} - * currently not working! - * @return The {@link Polygon} of the convex hull - * @throws GeoPreparationException If points got missing, no boundary can be found, too few points - * are provided or anything else goes wrong - * @throws NullPointerException If there is null - * @deprecated ATTENTION! DON'T REMOVE: This method is currently not under test and has to be - * revised thoroughly - */ - @Deprecated - public static Polygon buildConvexHull( - Set points, int precision, ConvexHullAlgorithm algorithm) - throws GeoPreparationException { - /* "Cast" the points to java.awt.Point */ - Set candidatePoints = - points.stream() - .map( - p -> - new Point( - (int) (p.getX() * Math.pow(10, precision)), - (int) (p.getY() * Math.pow(10, precision)))) - .collect(Collectors.toSet()); - - if (candidatePoints.size() != points.size()) - throw new GeoPreparationException( - "Some points got lost while converting to java.awt classes."); - - /* - * Reduce computational complexity by removing all nodes that certainly do not belong the the convex hull - * (Akl–Toussaint heuristic). Find the bounding coordinates and those four points having - * 1) the minimum lat value, - * 2) the minimum lon value, - * 3) the maximum lat value, - * 4) the maximum lon value - * denoting a rectangle with the nodes certainly not in the convex hull. */ - OptionalDouble possiblyXMin = candidatePoints.stream().mapToDouble(Point::getX).min(); - OptionalDouble possiblyXMax = candidatePoints.stream().mapToDouble(Point::getX).max(); - OptionalDouble possiblyYMin = candidatePoints.stream().mapToDouble(Point::getY).min(); - OptionalDouble possiblyYMax = candidatePoints.stream().mapToDouble(Point::getY).max(); - if (!possiblyXMax.isPresent()) - throw new GeoPreparationException("Unable to get boundary rectangle of the node set."); - - double[] xBounds = new double[] {possiblyXMin.getAsDouble(), possiblyXMax.getAsDouble()}; - double[] yBounds = new double[] {possiblyYMin.getAsDouble(), possiblyYMax.getAsDouble()}; - - HashSet edgePoints = new HashSet<>(); - candidatePoints.stream() - .filter(p -> p.getX() == xBounds[0]) - .findFirst() - .ifPresent(edgePoints::add); - candidatePoints.stream() - .filter(p -> p.getX() == xBounds[1]) - .findFirst() - .ifPresent(edgePoints::add); - candidatePoints.stream() - .filter(p -> p.getY() == yBounds[0]) - .findFirst() - .ifPresent(edgePoints::add); - candidatePoints.stream() - .filter(p -> p.getY() == yBounds[1]) - .findFirst() - .ifPresent(edgePoints::add); - int[] x = edgePoints.stream().mapToInt(p -> p.x).toArray(); - int[] y = edgePoints.stream().mapToInt(p -> p.y).toArray(); - java.awt.Polygon filterQuadrilateral = new java.awt.Polygon(x, y, x.length); - - candidatePoints.removeIf(p -> filterQuadrilateral.contains(p) && !edgePoints.contains(p)); - - int pointCount = candidatePoints.size(); - - LinkedList hullPoints = new LinkedList<>(); - - double maxIterations = Math.log(Math.log(candidatePoints.size())) - 1; - if (maxIterations < 1 || algorithm.equals(ConvexHullAlgorithm.GRAHAM)) { - if (candidatePoints.size() < 3) - throw new GeoPreparationException("Cannot build a convex hull for less than three nodes."); - hullPoints.addAll(GrahamScan.getConvexHull(new ArrayList<>(candidatePoints))); - } else { - /* TODO: Currently not working! */ - /* Chan's algorithm (https://en.wikipedia.org/wiki/Chan%27s_algorithm) */ - Point point0 = new Point(-Integer.MAX_VALUE, 0); - Point point1 = - edgePoints.stream() - .filter(p -> p.getY() == yBounds[0]) - .findFirst() - .orElseThrow(() -> new RuntimeException("Cannot find a second node.")); - hullPoints.addLast(point1); - - /* Limit the amount of iterations to at most log(log(nodeCount)) */ - outloop: - for (int iterCount = 0; iterCount < maxIterations; iterCount++) { - double m = Math.pow(2, Math.pow(2, iterCount + 1d)); - - /* Split the candidate nodes into floor(pointCount/m) with about m points each subsets and calculate those convex hulls */ - int subsetCount = (int) Math.ceil((double) pointCount / m); - int pointsPerSubset = (int) Math.ceil((double) pointCount / subsetCount); - - LinkedList> subSets = new LinkedList<>(); - subSets.addLast(new ArrayList<>()); - for (Point point : candidatePoints) { - if (subSets.getLast().size() >= pointsPerSubset) subSets.addLast(new ArrayList<>()); - subSets.getLast().add(point); - } - - if (subSets.stream().anyMatch(s -> s.size() < 3)) { - throw new GeoPreparationException( - "The algorithm lead to a segmentation with less than three nodes, wherefore no convex hull may be built."); - } - - /* Build convex hull of each subset using Graham scan */ - LinkedList> subHulls = new LinkedList<>(); - for (List subPoints : subSets) { - subHulls.addLast(GrahamScan.getConvexHull(subPoints)); - } - - /* Perform jarvis binary search to build the convex hull of all points based on the sub hulls */ - for (int cntSearch = 0; cntSearch < m; cntSearch++) { - /* Collect all points beloging to the subhulls */ - List subHullPoints = new ArrayList<>(); - - /* Perform Jarvis binary search: - * Find the point in each convex hull maximising the angle between the last point of the overall hull - * and that point */ - Point pA; - Point pB; - if (hullPoints.size() == 1) { - pB = hullPoints.getLast(); - pA = point0; - } else { - pB = hullPoints.getLast(); - pA = hullPoints.get(hullPoints.size() - 2); - } - - double previousAngle = Math.atan2(pB.getY() - pA.getY(), pB.getX() - pA.getX()); - - for (List subHull : subHulls) { - /* Calculate the difference between the two straights */ - Point bestPoint = null; - double maxAngleDifference = Double.MIN_VALUE; - for (Point point : subHull) { - double angleDifference = - previousAngle - Math.atan2(pB.getY() - point.getY(), pB.getX() - point.getX()); - if (angleDifference > maxAngleDifference) { - maxAngleDifference = angleDifference; - bestPoint = point; - } - } - subHullPoints.add(bestPoint); - } - - /* Calculate the difference between the two straights */ - Point nextPoint = null; - double maxAngleDifference = Double.MIN_VALUE; - for (Point point : subHullPoints) { - double angleDifference = - previousAngle - Math.atan2(pB.getY() - point.getY(), pB.getX() - point.getX()); - if (angleDifference > maxAngleDifference) { - maxAngleDifference = angleDifference; - nextPoint = point; - } - } - if (nextPoint == null) - throw new GeoPreparationException("Cannot find a next node to add to the convex hull."); - - hullPoints.addLast(nextPoint); - - /* The hull is closed, therefore exit the search and give back the hull */ - if (hullPoints.getFirst().equals(hullPoints.getLast())) break outloop; - } - } - - throw new GeoPreparationException( - "Maximum amount of iterations reached, but no convex hull found..."); - } - - /* The hull is completed. Cast Point back to Node and build a polygon */ - int castCount = -1; - List hullNodes = new LinkedList<>(); - for (Point point : hullPoints) { - hullNodes.add( - new Node( - castCount++, - null, - new LatLon( - point.getY() / Math.pow(10, precision), point.getX() / Math.pow(10, precision)))); - } - Way hullWay = new Way(0, null, hullNodes); - if (!hullWay.isClosed()) { - logger.warn("Closed way manually."); - hullWay.getNodes().add(hullWay.getNodes().get(0)); - } - return new Polygon(hullWay); - } - - /** - * Calculates the intersection of two {@link Polygon}s {@code a} and {@code b}. - * - * @param a First {@link Polygon} - * @param b Second {@link Polygon} - * @return A {@link Polygon} if the intersection exists and null if not - * @deprecated ATTENTION! DON'T REMOVE: This method is currently not under test and has to be - * revised thoroughly - */ - @Deprecated - public static Polygon getIntersection(Polygon a, Polygon b) { - List sharedCoords = - a.getCoords().stream().filter(b::contains).collect(Collectors.toList()); - List additionalCoords = - b.getCoords().stream().filter(a::contains).collect(Collectors.toList()); - - if (sharedCoords.isEmpty() && additionalCoords.isEmpty()) return null; - else { - /* - * If the first node of the additional nodes has a higher distance to the last one of the shared nodes as to - * the first one, then inverse the order of the additional nodes - */ - if (sharedCoords.isEmpty()) sharedCoords.addAll(additionalCoords); - else if (!additionalCoords.isEmpty() - && calcHaversine( - additionalCoords.get(0).getLat(), - additionalCoords.get(0).getLon(), - sharedCoords.get(sharedCoords.size() - 1).getLat(), - sharedCoords.get(sharedCoords.size() - 1).getLon()) - .getValue() - .doubleValue() - > calcHaversine( - additionalCoords.get(0).getLat(), - additionalCoords.get(0).getLon(), - sharedCoords.get(0).getLat(), - sharedCoords.get(0).getLon()) - .getValue() - .doubleValue()) { - additionalCoords = Lists.reverse(additionalCoords); - } - - additionalCoords.removeAll(sharedCoords); // Remove duplicates - sharedCoords.addAll(additionalCoords); - - return new Polygon(sharedCoords); - } - } - - /** - * Calculates the area, which is surrounded by a closed way by the help of {@link - * DeprecatedGeoUtils#calcArea(Polygon)} - * - * @param w Closed way, that surrounds the area - * @return The covered area in {@link edu.ie3.util.quantities.PowerSystemUnits#SQUARE_METRE} - * @throws GeoPreparationException If some serious shit happens - * @deprecated ATTENTION! DON'T REMOVE: This method is currently not under test and has to be - * revised thoroughly - */ - @Deprecated - public static Quantity calcArea(Way w) throws GeoPreparationException { - if (!w.isClosed()) - throw new GeoPreparationException( - "Cannot determine the area covered by a way, that is not closed"); - return calcArea(new Polygon(w)); - } - - /** - * Calculate the area of an {@link Polygon} by adding subareas spanned between the longitude axis - * and the line segments on the polygon - * - * @param p {@link Polygon} whos area may be calculated - * @return The spanned area in {@link edu.ie3.util.quantities.PowerSystemUnits#SQUARE_METRE} - * @throws GeoPreparationException If some serious shit happens - * @deprecated ATTENTION! DON'T REMOVE: This method is currently not under test and has to be - * revised thoroughly - */ - @Deprecated - public static Quantity calcArea(Polygon p) throws GeoPreparationException { - /* Get the boundary of the Polygon */ - Bounds bounds = p.getBounds(); - double latMaxGlobal = bounds.getMaxLat(); - double lonMinGlobal = bounds.getMinLon(); - - /* Order the points in clockwise direction starting with the one with the highest y coordinate.. - * The "way" has to be closed (containing the start coordinate twice). But it may not be the "original" start - * coordinate, but the one with the highest latitutde. Otherwise the partial areas would be added incorrect. - * Therefore remove the duplicate coordinates (original start coordinate) and add the one with highest latitude - * once again.*/ - List coords = new LinkedList<>(p.getCoords()); - Map duplicateCoords = - coords.stream() - .filter(coord -> Collections.frequency(coords, coord) > 1) - .distinct() - .collect( - Collectors.toMap(coord -> coord, coord -> Collections.frequency(coords, coord))); - // Only remove one of the duplicate coords - Iterator it = coords.iterator(); - while (it.hasNext()) { - LatLon coord = it.next(); - if (duplicateCoords.containsKey(coord) && duplicateCoords.get(coord) > 1) { - duplicateCoords.put(coord, duplicateCoords.get(coord) - 1); - it.remove(); - } - } - Optional optStartCoord = - coords.stream().filter(c -> c.getLat() == latMaxGlobal).findFirst(); - if (!optStartCoord.isPresent()) - throw new GeoPreparationException( - "Did not find a suitable coordinate although a defined maximum latitude has been found..."); - - LinkedList orderedCoords = new LinkedList<>(); - LatLon startCoord = optStartCoord.get(); - int idxStart = coords.indexOf(startCoord); - int idxNext = idxStart + 1 == coords.size() ? 0 : idxStart + 1; - if (coords.get(idxNext).getLon() > coords.get(idxStart).getLon()) { - /* Order is already clockwise */ - orderedCoords.addAll(coords.subList(idxStart, coords.size())); - orderedCoords.addAll(coords.subList(0, idxStart)); - } else { - /* Order is counterclockwise */ - orderedCoords.addAll(Lists.reverse(coords.subList(0, idxNext))); - orderedCoords.addAll(Lists.reverse(coords.subList(idxNext, coords.size()))); - } - orderedCoords.addLast(orderedCoords.get(0)); - - /* Calculate the area by summing up the partial areas. - * Take two points. The distance between the mean of their longitudes and the polygons global longitude - * as well as the distance of both latitudes does do form the partial areas. - * Those of the right hand side (distance of lats > 0) do form positive areas, whereas the left hand side - * is subtracted. */ - Quantity area = Quantities.getQuantity(0d, SQUARE_METRE); - /* Go through the coordinates and calculate the partial areas */ - int idxPrevious = orderedCoords.size() - 1; - for (int idx = 0; idx < orderedCoords.size(); idx++) { - LatLon coord = orderedCoords.get(idx); - LatLon coordPrev = orderedCoords.get(idxPrevious); - - double maxLat = Math.max(coord.getLat(), coordPrev.getLat()); - double minLat = Math.min(coord.getLat(), coordPrev.getLat()); - double maxLon = Math.max(coord.getLon(), coordPrev.getLon()); - double minLon = Math.min(coord.getLon(), coordPrev.getLon()); - - /* This partial area obviously would be zero. */ - if (maxLat == minLat || maxLon == minLon) { - idxPrevious = idx; - continue; - } - - double meanLon = (maxLon + minLon) / 2; - - Quantity dX = calcHaversine(maxLat, meanLon, maxLat, lonMinGlobal); - Quantity dY = calcHaversine(maxLat, meanLon, minLat, meanLon); - Quantity partialArea = dX.multiply(dY).asType(Area.class).to(SQUARE_METRE); - if (coord.getLat() < coordPrev.getLat()) { - /* Right hand side*/ - area = area.add(partialArea); - } else { - /* Left hand side */ - area = area.subtract(partialArea); - } - - idxPrevious = idx; - } - - return area; - } - - /** - * Checks if Node c is between Node a and b - * - * @param a Node A - * @param b Node B - * @param c Node C - * @return True if node c is between node a and node b - * @deprecated ATTENTION! DON'T REMOVE: This method is currently not under test and has to be - * revised thoroughly - */ - @Deprecated - public static boolean isBetween(Node a, Node b, Node c) { - double crossProduct; - double dotProduct; - double squaredLengthBA; - double epsilon = 0.000000000001; // epsilon to check if a,b and c are aligned - // TODO: as the size of epsilon is crucial for the functionality, it should be available as an - // option - // lon = x - // lat = y - - crossProduct = - ((c.getLatlon().getLat() - a.getLatlon().getLat()) - * (b.getLatlon().getLon() - a.getLatlon().getLon()) - - (c.getLatlon().getLon() - a.getLatlon().getLon()) - * (b.getLatlon().getLat() - a.getLatlon().getLat())); - if (Math.abs(crossProduct) > epsilon) return false; - - dotProduct = - (c.getLatlon().getLon() - a.getLatlon().getLon()) - * (b.getLatlon().getLon() - a.getLatlon().getLon()) - + (c.getLatlon().getLat() - a.getLatlon().getLat()) - * (b.getLatlon().getLat() - a.getLatlon().getLat()); - - if (dotProduct < 0) return false; - - squaredLengthBA = - (b.getLatlon().getLon() - a.getLatlon().getLon()) - * (b.getLatlon().getLon() - a.getLatlon().getLon()) - + (b.getLatlon().getLat() - a.getLatlon().getLat()) - * (b.getLatlon().getLat() - a.getLatlon().getLat()); - - return dotProduct > squaredLengthBA; - } - - /** - * Calculates the area of a polygon in geo coordinates to an area in square kilometers NOTE: It - * may be possible, that (compared to the real building area), the size of the building area may - * be overestimated. To take this into account an optional correction factor might be used. - * - * @param geoArea: the area of the building based on geo coordinates - * @param cor: the optional correction factor - * @deprecated ATTENTION! DON'T REMOVE: This method is currently not under test and has to be - * revised thoroughly - */ - @Deprecated - public static Quantity calcGeo2qmNew(double geoArea, Quantity cor) { - double width = 51.5; - double length = 7.401; - - double square = Math.sqrt(geoArea); - double width2 = (width + square) / 180 * Math.PI; - double length2 = length / 180 * Math.PI; - double width3 = width / 180 * Math.PI; - double length3 = (length + square) / 180 * Math.PI; - width = width / 180 * Math.PI; - length = length / 180 * Math.PI; - - double e1 = - Math.acos( - Math.sin(width) * Math.sin(width2) - + Math.cos(width) * Math.cos(width2) * Math.cos(length2 - length)) - * EARTH_RADIUS.getValue().doubleValue(); - double e2 = - Math.acos( - Math.sin(width) * Math.sin(width3) - + Math.cos(width) * Math.cos(width3) * Math.cos(length3 - length)) - * EARTH_RADIUS.getValue().doubleValue(); - - /* (e1 * e2) - cor */ - return Quantities.getQuantity(e1, METRE) - .multiply(Quantities.getQuantity(e2, METRE)) - .asType(Area.class) - .subtract(cor); - } - - /** - * @deprecated ATTENTION! DON'T REMOVE: This method is currently not under test and has to be - * revised thoroughly - */ - @Deprecated - public static boolean isInsideLanduse(LatLon node, List landUses) { - for (Way landUse : landUses) { - if (rayCasting(new Polygon(landUse), node)) return true; - } - return false; - } - - /** - * @deprecated ATTENTION! DON'T REMOVE: This method is currently not under test and has to be - * revised thoroughly - */ - @Deprecated - public static boolean rayCasting(Polygon shape, LatLon node) { - boolean inside = false; - - // get lat lon from shape nodes - List shapeNodes = shape.getCoords(); - - for (int i = 1; i < shapeNodes.size(); i++) { - if (intersects(shapeNodes.get(i - 1), shapeNodes.get((i)), node)) inside = !inside; - } - return inside; - } - - /** - * @deprecated ATTENTION! DON'T REMOVE: This method is currently not under test and has to be - * revised thoroughly - */ - @Deprecated - private static boolean intersects(LatLon aIn, LatLon bIn, LatLon nIn) { - - // convert LatLons to arrays - double[] a = {aIn.getLon(), aIn.getLat()}; - double[] b = {bIn.getLon(), bIn.getLat()}; - double[] n = {nIn.getLon(), nIn.getLat()}; - - if (a[1] > b[1]) return intersects(bIn, aIn, nIn); - - if (n[1] == a[1] || n[1] == b[1]) n[1] += 0.0001; - - if (n[1] > b[1] || n[1] < a[1] || n[0] > Math.max(a[0], b[0])) return false; - - if (n[0] < Math.min(a[0], b[0])) return true; - - double red = (n[1] - a[1]) / (n[0] - a[0]); - double blue = (b[1] - a[1]) / (b[0] - a[0]); - return red >= blue; - } - - /** - * The algorithm assumes the usual mathematical convention that positive y points upwards. In - * computer systems where positive y is downwards (most of them) the easiest thing to do is list - * the vertices counter-clockwise using the "positive y down" coordinates. The two effects then - * cancel out to produce a positive area. - * - * @param building - * @return polygon area A - * @deprecated ATTENTION! DON'T REMOVE: This method is currently not under test and has to be - * revised thoroughly - */ - @Deprecated - public static double calculateBuildingArea(Way building) { - double area = 0.0; - - // extract number of polygon points - int j = - building.getNodes().size() - - 1; // -1 because the last vertex is the 'previous' one to the first - // x = long - // y = lat - for (int i = 0; i < building.getNodes().size(); i++) { - area = - area - + ((building.getNodes().get(j).getLatlon().getLon() - + building.getNodes().get(i).getLatlon().getLon()) - * (building.getNodes().get(j).getLatlon().getLat() - - building.getNodes().get(i).getLatlon().getLat())); - j = i; - } - - area = Math.abs(area); - - return area / 2; - } - - /** - * Draws a circle with a radius of the provided distance around the provided center coordinates - * and returns the result as a drawable polygon (one point per degree) - * - * @param center coordinate of the circle's center - * @param radius radius of the circle - * @return a polygon without the center, but with all points of the circle - */ - public static Polygon radiusWithCircleAsPolygon(LatLon center, Quantity radius) { - - List circlePoints = radiusWithCircle(center, radius); - - List circlePointsLatLon = - circlePoints.stream() - .map(point -> new LatLon(point.getLat(), point.getLon())) - .collect(Collectors.toList()); - - return new Polygon(circlePointsLatLon); - } - - /** - * Draws a circle with a radius of the provided distance around the provided center coordinates - * and returns the result as list of all circle points (one per degree) - * - * @param center - * @param radius - * @return a list with all coordinates of the circle points - * @deprecated ATTENTION! DON'T REMOVE: This method is currently not under test and has to be - * revised thoroughly - */ - @Deprecated - public static List radiusWithCircle(LatLon center, Quantity radius) { - - double lat1 = Math.toRadians(center.getLat()); - double lon1 = Math.toRadians(center.getLon()); - double d = (radius.divide(EARTH_RADIUS)).getValue().doubleValue(); - - List locs = new ArrayList<>(); - - for (int x = 0; x <= 360; x++) { - - double brng = Math.toRadians(x); - - double latRad = - Math.asin(Math.sin(lat1) * Math.cos(d) + Math.cos(lat1) * Math.sin(d) * Math.cos(brng)); - double lonRad = - lon1 - + Math.atan2( - Math.sin(brng) * Math.sin(d) * Math.cos(lat1), - Math.cos(d) - Math.sin(lat1) * Math.sin(latRad)); - - locs.add(new LatLon(Math.toDegrees(latRad), Math.toDegrees(lonRad))); - } - - return locs; - } - - /** - * Chains together several ways which are close together (e.g. as part of a relation) to one, - * correctly ordered way. This is a computational expensive method which might be worth - * refactoring! - * - * @param waysToChain the ways from a relation - * @param radius the radius that should be considered for search for the next connection point - * between two ways - * @param wayId the id of the newly created way - * @return A concatenated {@link Way} built from a collection of {@link Way}s - * @deprecated ATTENTION! DON'T REMOVE: This method is currently not under test and has to be - * revised thoroughly - */ - @Deprecated - public static Way wayFromWays(List waysToChain, Quantity radius, int wayId) { - - LinkedList waysCopy = new LinkedList<>(waysToChain); - - HashMap nodeToWayMap = new HashMap<>(); - for (Way way : waysToChain) { - way.getNodes().forEach(node -> nodeToWayMap.put(node, way)); - } - - LinkedList nodesList = new LinkedList<>(waysToChain.get(0).getNodes()); - waysCopy.remove(waysToChain.get(0)); - - while (!waysCopy.isEmpty()) { - Node lastNode = nodesList.getLast(); - - Polygon circle = radiusWithCircleAsPolygon(lastNode.getLatlon(), radius); - - double distance = Double.POSITIVE_INFINITY; - Way nextWay = null; - Node nextLastNode = null; - for (Map.Entry entry : nodeToWayMap.entrySet()) { - if (rayCasting(circle, entry.getKey().getLatlon())) { - double tempDistance; - tempDistance = - calcHaversine( - lastNode.getLatlon().getLat(), - lastNode.getLatlon().getLon(), - entry.getKey().getLatlon().getLat(), - entry.getKey().getLatlon().getLon()) - .to(KILOMETRE) - .getValue() - .doubleValue(); - if (tempDistance < distance) { - distance = tempDistance; - nextWay = entry.getValue(); - nextLastNode = entry.getKey(); - } - } - } - - if (nextWay != null) { - if (nextWay.getNodes().indexOf(nextLastNode) == nextWay.getNodes().size() - 1) { - Collections.reverse(nextWay.getNodes()); - } - nodesList.addAll(nextWay.getNodes()); - - logger.debug("Removing way with id {}", nextWay.getId()); - } - - nodeToWayMap.values().removeAll(Collections.singleton(nextWay)); - waysCopy.remove(nextWay); - } - - return new Way(wayId, new Tags(), nodesList); - } - - public enum ConvexHullAlgorithm { - CHAN, - GRAHAM - } -} diff --git a/src/test/groovy/edu/ie3/util/geo/DeprecatedCoordinateDistanceTest.groovy b/src/test/groovy/edu/ie3/util/geo/DeprecatedCoordinateDistanceTest.groovy deleted file mode 100644 index 691e988b..00000000 --- a/src/test/groovy/edu/ie3/util/geo/DeprecatedCoordinateDistanceTest.groovy +++ /dev/null @@ -1,42 +0,0 @@ -/* - * © 2021. TU Dortmund University, - * Institute of Energy Systems, Energy Efficiency and Energy Economics, - * Research group Distribution grid planning and operation - */ -package edu.ie3.util.geo - -import spock.lang.Specification - -class DeprecatedCoordinateDistanceTest extends Specification { - def "The constructor without a distance parameter calculates the distance as expected"() { - given: - def pointA = DeprecatedGeoUtils.xyToPoint(49d, 7d) - def pointB = DeprecatedGeoUtils.xyToPoint(50d, 7d) - def expectedDistance = DeprecatedGeoUtils.calcHaversine(pointA.y, pointA.x, pointB.y, pointB.x) - - when: - def coordinateDistance = new DeprecatedCoordinateDistance(pointA, pointB) - def expectedCoordinateDistance = new DeprecatedCoordinateDistance(coordinateDistance.coordinateA, coordinateDistance.coordinateB, expectedDistance) - - then: - coordinateDistance.distance == expectedDistance - // this equals can not be replaced with the == operator as this would cause equals not to be called - coordinateDistance.equals(expectedCoordinateDistance) - } - - def "CoordinateDistances are sortable using their distance field"() { - given: - def basePoint = DeprecatedGeoUtils.xyToPoint(49d, 7d) - def distA = new DeprecatedCoordinateDistance(basePoint, DeprecatedGeoUtils.xyToPoint(50d, 7d)) - def distB = new DeprecatedCoordinateDistance(basePoint, DeprecatedGeoUtils.xyToPoint(50d, 7.1d)) - def distC = new DeprecatedCoordinateDistance(basePoint, DeprecatedGeoUtils.xyToPoint(49d, 7.1d)) - def distD = new DeprecatedCoordinateDistance(basePoint, DeprecatedGeoUtils.xyToPoint(52d, 9d)) - def coordinateDistances = [distA, distB, distC, distD] - - when: - def sortedDistances = coordinateDistances.toSorted() - - then: - sortedDistances == [distC, distA, distB, distD].toList() - } -} diff --git a/src/test/groovy/edu/ie3/util/geo/GeoUtilsTest.groovy b/src/test/groovy/edu/ie3/util/geo/GeoUtilsTest.groovy deleted file mode 100644 index af8551dc..00000000 --- a/src/test/groovy/edu/ie3/util/geo/GeoUtilsTest.groovy +++ /dev/null @@ -1,191 +0,0 @@ -/* - * © 2020. TU Dortmund University, - * Institute of Energy Systems, Energy Efficiency and Energy Economics, - * Research group Distribution grid planning and operation - */ -package edu.ie3.util.geo - -import edu.ie3.util.quantities.PowerSystemUnits -import edu.ie3.util.quantities.QuantityUtil -import org.locationtech.jts.geom.Coordinate -import org.locationtech.jts.geom.LineString -import org.locationtech.jts.io.geojson.GeoJsonReader -import spock.lang.Shared - -import static edu.ie3.util.quantities.PowerSystemUnits.METRE - -import net.morbz.osmonaut.geometry.Polygon -import net.morbz.osmonaut.osm.LatLon -import spock.lang.Specification - -import javax.measure.Quantity -import javax.measure.quantity.Length -import tech.units.indriya.ComparableQuantity -import tech.units.indriya.quantity.Quantities - - - -class GeoUtilsTest extends Specification { - - @Shared - GeoJsonReader geoJsonReader - - def setupSpec() { - geoJsonReader = new GeoJsonReader() - } - - def "Test haversine (distance between two points given lat/lon)"() { - given: - LatLon start = new LatLon(37.87532764735112, -122.25311279296875) - LatLon end = new LatLon(37.87934174490509, -122.2537350654602) - ComparableQuantity tolerance = Quantities.getQuantity(1d, METRE) - ComparableQuantity expected = Quantities.getQuantity(450.18011568984845, METRE) - - when: - ComparableQuantity actual = DeprecatedGeoUtils.calcHaversine(start.lat, start.lon, end.lat, end.lon) - - then: - Math.abs(actual.subtract(expected).to(METRE).value.doubleValue()) < tolerance .value.doubleValue() - } - - def "Total length of LineString is correctly calculated"() { - given: - def lineString = DeprecatedGeoUtils.DEFAULT_GEOMETRY_FACTORY.createLineString([ - new Coordinate(22.69962d, 11.13038d, 0), - new Coordinate(20.84247d, 28.14743d, 0), - new Coordinate(24.21942d, 12.04265d, 0) - ] as Coordinate[]) - - when: - ComparableQuantity y = DeprecatedGeoUtils.totalLengthOfLineString(lineString) - - then: - QuantityUtil.isEquivalentAbs(y, Quantities.getQuantity(3463.37, PowerSystemUnits.KILOMETRE), 10) - // Value from Google Maps, error range of +-10 km - } - - def "The GridAndGeoUtils should get the CoordinateDistances between a base point and a collection of other points correctly"() { - given: - def basePoint = DeprecatedGeoUtils.xyToPoint(49d, 7d) - def points = [ - DeprecatedGeoUtils.xyToPoint(50d, 7d), - DeprecatedGeoUtils.xyToPoint(50d, 7.1d), - DeprecatedGeoUtils.xyToPoint(49d, 7.1d), - DeprecatedGeoUtils.xyToPoint(52d, 9d) - ] - def coordinateDistances = [ - new DeprecatedCoordinateDistance(basePoint, points[0]), - new DeprecatedCoordinateDistance(basePoint, points[1]), - new DeprecatedCoordinateDistance(basePoint, points[2]), - new DeprecatedCoordinateDistance(basePoint, points[3]) - ] - expect: - DeprecatedGeoUtils.getCoordinateDistances(basePoint, points) == new TreeSet(coordinateDistances) - } - - def "The GridAndGeoUtils should build a line string with two exact equal geo coordinates correctly avoiding the known bug in jts geometry"() { - given: - def line = geoJsonReader.read(lineString) as LineString - - when: - def safeLineString = DeprecatedGeoUtils.buildSafeLineString(line) - def actualCoordinates = safeLineString.coordinates - - then: - coordinates.length == actualCoordinates.length - for (int cnt = 0; cnt < coordinates.length; cnt++) { - coordinates[cnt] == actualCoordinates[cnt] - } - from - where: - lineString | coordinates - "{ \"type\": \"LineString\", \"coordinates\": [[7.411111, 51.49228],[7.411111, 51.49228]]}" | [ - new Coordinate(7.4111110000001, 51.4922800000001), - new Coordinate(7.411111, 51.49228) - ] as Coordinate[] - "{ \"type\": \"LineString\", \"coordinates\": [[7.411111, 51.49228],[7.411111, 51.49228],[7.411111, 51.49228],[7.411111, 51.49228]]}" | [ - new Coordinate(7.4111110000001, 51.4922800000001), - new Coordinate(7.411111, 51.49228) - ] as Coordinate[] - "{ \"type\": \"LineString\", \"coordinates\": [[7.411111, 51.49228],[7.411111, 51.49228],[7.311111, 51.49228],[7.511111, 51.49228]]}" | [ - new Coordinate(7.411111, 51.49228), - new Coordinate(7.311111, 51.49228), - new Coordinate(7.511111, 51.49228) - ] as Coordinate[] - } - - def "The GridAngGeoUtils maintain the correct order of coordinates, when overhauling a given LineString"() { - /* Remark: This test might even NOT fail, if the method is implemented incorrectly (utilizing a HashSet to - * maintain uniqueness). For detailed explanation cf. comment in method's implementation. */ - given: - def coordinates = [ - new Coordinate(51.49292, 7.41197), - new Coordinate(51.49333, 7.41183), - new Coordinate(51.49341, 7.41189), - new Coordinate(51.49391, 7.41172), - new Coordinate(51.49404, 7.41279) - ] as Coordinate[] - def lineString = DeprecatedGeoUtils.DEFAULT_GEOMETRY_FACTORY.createLineString(coordinates) - - when: - def actual = DeprecatedGeoUtils.buildSafeLineString(lineString).coordinates - - then: - coordinates.length == actual.length - for (int cnt = 0; cnt < coordinates.length; cnt++) { - coordinates[cnt] == actual[cnt] - } - } - - def "The GridAndGeoUtils should only modify a provided Coordinate as least as possible"() { - given: - def coord = new Coordinate(1, 1, 0) - - expect: - DeprecatedGeoUtils.buildSafeCoord(coord) == new Coordinate(1.0000000000001, 1.0000000000001, 1.0E-13) - } - - def "The GridAndGeoUtils should build a safe instance of a LineString between two provided coordinates correctly"() { - - expect: - DeprecatedGeoUtils.buildSafeLineStringBetweenCoords(coordA, coordB) == resLineString - // do not change or remove the following line, it is NOT equal to the line above in this case! - DeprecatedGeoUtils.buildSafeLineStringBetweenCoords(coordA, coordB).equals(resLineString) - - where: - coordA | coordB || resLineString - new Coordinate(1, 1, 0) | new Coordinate(1, 1, 0) || DeprecatedGeoUtils.DEFAULT_GEOMETRY_FACTORY.createLineString([ - new Coordinate(1.0000000000001, 1.0000000000001, 1.0E-13), - new Coordinate(1, 1, 0) - ] as Coordinate[]) - new Coordinate(1, 1, 0) | new Coordinate(2, 2, 0) || DeprecatedGeoUtils.DEFAULT_GEOMETRY_FACTORY.createLineString([ - new Coordinate(1, 1, 0), - new Coordinate(2, 2, 0) - ] as Coordinate[]) - - } - - def "Test radius with circle as polygon"() { - given: - LatLon center = new LatLon(52.02083574, 7.40110716) - Quantity radius = Quantities.getQuantity(50, METRE) - - when: - Polygon poly = DeprecatedGeoUtils.radiusWithCircleAsPolygon(center, radius) - List circlePoints = poly.getCoords() - - then: - // polygon should contain a center that is the provided center - Math.round(poly.center.lat * 100000000) / 100000000 == Math.round(center.lat * 100000000) / 100000000 - Math.round(poly.center.lon * 100000000) / 100000000 == Math.round(center.lon * 100000000) / 100000000 - - // number of expected circle points - circlePoints.size() == 361 - // rounded distance should be 50 meters - circlePoints.forEach({ point -> - Double distance = DeprecatedGeoUtils.calcHaversine(center.lat, center.lon, point.lat, point.lon).to(METRE).value.doubleValue() - Math.round(distance) == 50 - }) - - } -} From f8e5c177af15e5e670fa7831517b225b8f506f3e Mon Sep 17 00:00:00 2001 From: t-ober <63147366+t-ober@users.noreply.github.com> Date: Thu, 9 Dec 2021 11:41:46 +0100 Subject: [PATCH 09/17] adjust test --- src/main/scala/edu/ie3/util/geo/GeoUtils.scala | 8 +++----- .../scala/edu/ie3/util/geo/GeoUtilsSpec.scala | 15 +++++++++------ 2 files changed, 12 insertions(+), 11 deletions(-) diff --git a/src/main/scala/edu/ie3/util/geo/GeoUtils.scala b/src/main/scala/edu/ie3/util/geo/GeoUtils.scala index 8fef07dc..a2be9686 100644 --- a/src/main/scala/edu/ie3/util/geo/GeoUtils.scala +++ b/src/main/scala/edu/ie3/util/geo/GeoUtils.scala @@ -201,7 +201,7 @@ object GeoUtils { * @return * a Try of the resulting polygon */ - def buildConvexHull(coordinates: List[Coordinate]): Try[Polygon] = { + def buildConvexHull(coordinates: Set[Coordinate]): Try[Polygon] = { val projectedCoordinates = coordinates.map(equalAreaProjection) new ConvexHull( projectedCoordinates.toArray, @@ -276,10 +276,8 @@ object GeoUtils { * Credits to Joe Kington * (https://stackoverflow.com/questions/4681737/how-to-calculate-the-area-of-a-polygon-on-the-earths-surface-using-python#:~:text=Basically%2C%20you%20just%20multiply%20the,the%20cosine%20of%20the%20latitude.) * - * @param lat - * latitude to project - * @param long - * longitude to project + * @param coordinate + * the coordinate to project * @return * a projected Coordinate with values in metre */ diff --git a/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala b/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala index 2c7e6669..131ab71e 100644 --- a/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala +++ b/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala @@ -103,7 +103,7 @@ class GeoUtilsSpec extends Matchers with AnyWordSpecLike { val bottomRight = new Coordinate(8, 48) val betweenBlBr = new Coordinate(7.5, 48) val bottomLeft = new Coordinate(7, 48) - val coordinates = List( + val coordinates = Set( topLeft, betweenTlTr, topRight, @@ -116,13 +116,16 @@ class GeoUtilsSpec extends Matchers with AnyWordSpecLike { case Failure(exception) => throw exception case Success(polygon) => val hullCoordinates = polygon.getCoordinates + val cornerCoordinates = + List(topLeft, topRight, bottomRight, bottomLeft) // contains all corner coordinates and filters out all colinear ones + hullCoordinates.foreach(coord => println(coord.toString)) hullCoordinates.size shouldBe 5 - hullCoordinates(0).equals2D(topLeft, 1e-10) - hullCoordinates(1).equals(topRight, 1e-10) shouldBe true - hullCoordinates(2).equals(bottomRight, 1e-10) shouldBe true - hullCoordinates(3).equals(bottomLeft, 1e-10) shouldBe true - hullCoordinates(4).equals(topLeft, 1e-10) shouldBe true + hullCoordinates.foreach(hullCoordinate => + cornerCoordinates.exists(cornerCoordinate => + hullCoordinate.equals2D(cornerCoordinate, 1e-10) + ) + ) } } From 3de10c53d2e117383e265d1f8c9593dc30d5b560 Mon Sep 17 00:00:00 2001 From: t-ober <63147366+t-ober@users.noreply.github.com> Date: Thu, 9 Dec 2021 14:27:51 +0100 Subject: [PATCH 10/17] exten rich geometry tests --- .../scala/edu/ie3/util/geo/GeoUtils.scala | 8 +- .../edu/ie3/util/geo/RichGeometries.scala | 22 +-- .../scala/edu/ie3/util/osm/OsmEntities.scala | 2 +- .../scala/edu/ie3/util/geo/GeoUtilsSpec.scala | 24 +++- .../edu/ie3/util/geo/RichGeometriesSpec.scala | 131 +++++++++++++++++- 5 files changed, 160 insertions(+), 27 deletions(-) diff --git a/src/main/scala/edu/ie3/util/geo/GeoUtils.scala b/src/main/scala/edu/ie3/util/geo/GeoUtils.scala index a2be9686..0cc1b10f 100644 --- a/src/main/scala/edu/ie3/util/geo/GeoUtils.scala +++ b/src/main/scala/edu/ie3/util/geo/GeoUtils.scala @@ -155,7 +155,7 @@ object GeoUtils { */ def calcOrderedCoordinateDistances( baseCoordinate: Point, - coordinates: List[Point] + coordinates: Array[Point] ): SortedSet[CoordinateDistance] = { val coordinateDistances = coordinates.map(coordinate => new CoordinateDistance(baseCoordinate, coordinate) @@ -210,7 +210,7 @@ object GeoUtils { case projectedPolygon: Polygon => val coordinates = projectedPolygon.getCoordinates.map(reverseEqualAreaProjection) - Success(buildPolygon(coordinates.toList)) + Success(buildPolygon(coordinates)) case _: LineString => Failure( new GeoException( @@ -263,7 +263,7 @@ object GeoUtils { * @return * a [[Polygon]] */ - def buildPolygon(coordinates: List[Coordinate]): Polygon = { + def buildPolygon(coordinates: Array[Coordinate]): Polygon = { val arrayCoordinates = new CoordinateArraySequence(coordinates.toArray) val linearRing = new LinearRing(arrayCoordinates, DEFAULT_GEOMETRY_FACTORY) @@ -338,7 +338,7 @@ object GeoUtils { ) new Coordinate(lonRad.toDegrees, latRad.toDegrees) }) - .toList + .toArray buildPolygon(coordinates) } diff --git a/src/main/scala/edu/ie3/util/geo/RichGeometries.scala b/src/main/scala/edu/ie3/util/geo/RichGeometries.scala index 99574ad5..8b06b1c5 100644 --- a/src/main/scala/edu/ie3/util/geo/RichGeometries.scala +++ b/src/main/scala/edu/ie3/util/geo/RichGeometries.scala @@ -12,11 +12,11 @@ import edu.ie3.util.geo.GeoUtils.{ calcHaversine, reverseEqualAreaProjection } -import edu.ie3.util.quantities.PowerSystemUnits.KILOMETRE import edu.ie3.util.quantities.QuantityUtils.RichQuantityDouble import org.locationtech.jts.geom.{Coordinate, LineString, Point, Polygon} import tech.units.indriya.ComparableQuantity import tech.units.indriya.quantity.Quantities +import tech.units.indriya.unit.Units.METRE import javax.measure.quantity.{Area, Length} import scala.math.abs @@ -44,21 +44,23 @@ object RichGeometries { ) } - /** Checks if the coordinate lies between two coordinates a and b. + /** Checks if the coordinate lies between two coordinates a and b by + * comparing the distances between a and b with the sum of distances + * between the coordinate and a and the coordinate and b * * @param a * coordinate a * @param b * coordinate b * @param epsilon - * permitted deviation + * permitted relative deviation * @return * whether or not the coordinate lies between */ def isBetween( a: Coordinate, b: Coordinate, - epsilon: Double = 1e-11 + epsilon: Double = 1e-12 ): Boolean = { val distance = a.haversineDistance(b) val distancePassingMe = a @@ -88,7 +90,7 @@ object RichGeometries { val coordinates = lineString.getCoordinates.toVector val coordinateSize = coordinates.size coordinates.zipWithIndex - .foldLeft(Quantities.getQuantity(0, KILOMETRE))((acc, current) => { + .foldLeft(Quantities.getQuantity(0, METRE))((acc, current) => { if (current._2 < coordinateSize - 1) { val currentCoordinate = current._1 val nextCoordinate = coordinates(current._2 + 1) @@ -107,9 +109,9 @@ object RichGeometries { * polygon with which to calculate the intersection * @return */ - def intersect(polygonB: Polygon): Try[Set[Coordinate]] = { + def intersect(polygonB: Polygon): Try[Polygon] = { Try( - polygon.equalAreaProjection.intersection(polygonB.equalAreaProjection) + polygon.intersection(polygonB.equalAreaProjection) ) match { case Failure(exception) => Failure( @@ -118,8 +120,7 @@ object RichGeometries { exception ) ) - case Success(geometry) => - Success(geometry.getCoordinates.map(reverseEqualAreaProjection).toSet) + case Success(polygon: Polygon) => Success(polygon) } } @@ -143,10 +144,9 @@ object RichGeometries { val projectedCoordinates = polygon.getCoordinates.map(coordinate => GeoUtils.equalAreaProjection(coordinate) ) - buildPolygon(projectedCoordinates.toList) + buildPolygon(projectedCoordinates) } - // Fixme: Should we do an equal area projection here as well? /** Checks whether the polygon contains the coordinate. Uses "covers()" * insted of "contains()" so borders are included. * diff --git a/src/main/scala/edu/ie3/util/osm/OsmEntities.scala b/src/main/scala/edu/ie3/util/osm/OsmEntities.scala index 41e3ce84..d56251cc 100644 --- a/src/main/scala/edu/ie3/util/osm/OsmEntities.scala +++ b/src/main/scala/edu/ie3/util/osm/OsmEntities.scala @@ -133,7 +133,7 @@ object OsmEntities { } def toPolygon: Polygon = { - buildPolygon(getCoordinates) + buildPolygon(getCoordinates.toArray) } def calculateArea: ComparableQuantity[Area] = diff --git a/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala b/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala index 131ab71e..c2bb5cf1 100644 --- a/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala +++ b/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala @@ -9,6 +9,7 @@ import edu.ie3.util.geo.GeoUtils.{ DEFAULT_GEOMETRY_FACTORY, buildCirclePolygon, buildConvexHull, + buildCoordinate, buildPoint, buildPolygon, calcHaversine, @@ -16,7 +17,7 @@ import edu.ie3.util.geo.GeoUtils.{ equalAreaProjection, reverseEqualAreaProjection } -import edu.ie3.util.geo.RichGeometries.GeoCoordinate +import edu.ie3.util.geo.RichGeometries.{GeoCoordinate, GeoPolygon} import edu.ie3.util.quantities.QuantityMatchers.equalWithTolerance import org.locationtech.jts.geom.Coordinate import org.scalatest.matchers.should.Matchers @@ -78,7 +79,7 @@ class GeoUtilsSpec extends Matchers with AnyWordSpecLike { val middle = buildPoint(51.44074066208236, 10.039100585737868) val furthest = buildPoint(51.6273949370414, 10.341911247878777) val actual = - calcOrderedCoordinateDistances(base, List(middle, furthest, nearest)) + calcOrderedCoordinateDistances(base, Array(middle, furthest, nearest)) val list = actual.toList list.size shouldBe 3 list.head.to shouldBe nearest @@ -119,12 +120,11 @@ class GeoUtilsSpec extends Matchers with AnyWordSpecLike { val cornerCoordinates = List(topLeft, topRight, bottomRight, bottomLeft) // contains all corner coordinates and filters out all colinear ones - hullCoordinates.foreach(coord => println(coord.toString)) hullCoordinates.size shouldBe 5 hullCoordinates.foreach(hullCoordinate => cornerCoordinates.exists(cornerCoordinate => hullCoordinate.equals2D(cornerCoordinate, 1e-10) - ) + ) shouldBe true ) } } @@ -142,7 +142,7 @@ class GeoUtilsSpec extends Matchers with AnyWordSpecLike { val coordinateB = new Coordinate(7.521007845835815, 51.50450661471354) val coordinateC = new Coordinate(7.5598606548385545, 51.456498140367934) val actual = - buildPolygon(List(coordinateA, coordinateB, coordinateC, coordinateA)) + buildPolygon(Array(coordinateA, coordinateB, coordinateC, coordinateA)) actual.getCoordinates shouldBe Array( coordinateA, coordinateB, @@ -190,6 +190,18 @@ class GeoUtilsSpec extends Matchers with AnyWordSpecLike { ) should equalWithTolerance(radius) } - } + "checks if the polygon contains a coordinate" in { + val topLeft = new Coordinate(7, 50) + val topRight = new Coordinate(8, 50) + val bottomRight = new Coordinate(8, 48) + val bottomLeft = new Coordinate(7, 48) + val polygon = + buildPolygon(Array(topLeft, topRight, bottomRight, bottomLeft, topLeft)) + polygon.containsCoordinate(new Coordinate(7.5, 50)) shouldBe true + polygon.containsCoordinate(new Coordinate(7.5, 49)) shouldBe true + polygon.containsCoordinate(new Coordinate(7.5, 50.1)) shouldBe false + polygon.containsCoordinate(new Coordinate(8.1, 50)) shouldBe false + } + } } diff --git a/src/test/scala/edu/ie3/util/geo/RichGeometriesSpec.scala b/src/test/scala/edu/ie3/util/geo/RichGeometriesSpec.scala index a3c9a103..5d0e136f 100644 --- a/src/test/scala/edu/ie3/util/geo/RichGeometriesSpec.scala +++ b/src/test/scala/edu/ie3/util/geo/RichGeometriesSpec.scala @@ -5,16 +5,139 @@ */ package edu.ie3.util.geo -import edu.ie3.util.geo.GeoUtils.DEFAULT_GEOMETRY_FACTORY -import edu.ie3.util.geo.RichGeometries.GeoPolygon +import edu.ie3.util.geo.GeoUtils.{ + DEFAULT_GEOMETRY_FACTORY, + buildCoordinate, + buildPolygon, + calcHaversine +} +import edu.ie3.util.geo.RichGeometries.{ + GeoCoordinate, + GeoLineString, + GeoPolygon +} +import edu.ie3.util.quantities.QuantityMatchers.equalWithTolerance import edu.ie3.util.quantities.QuantityUtils.RichQuantityDouble import org.locationtech.jts.geom.Coordinate import org.scalatest.matchers.should.Matchers import org.scalatest.wordspec.AnyWordSpecLike +import tech.units.indriya.quantity.Quantities +import tech.units.indriya.unit.Units.METRE + +import scala.util.{Failure, Success} class RichGeometriesSpec extends Matchers with AnyWordSpecLike { - "A rich polygon" should { + "A rich GeoCordinate" should { + + "calculate the haversine distance correctly" in { + val latA = 37.87532764735112 + val longA = -122.25311279296875 + val latB = 37.87934174490509 + val longB = -122.2537350654602 + val coordinateA = buildCoordinate(latA, longA) + val coordinateB = buildCoordinate(latB, longB) + val expected = Quantities.getQuantity(450.18011568984845, METRE) + val actual = coordinateA.haversineDistance(coordinateB) + actual should equalWithTolerance(expected) + } + + "check if a coordinate is between two others" in { + val coordinateA = buildCoordinate(10, 5) + val coordinateB = buildCoordinate(20, 5) + val coordinateBetween = buildCoordinate(15, 5) + val coordinateNotExactlyBetween = buildCoordinate(16, 5.005) + val coordinateNotBetween = buildCoordinate(15, 6) + + coordinateBetween.isBetween(coordinateA, coordinateB) shouldBe true + coordinateBetween.isBetween(coordinateB, coordinateA) shouldBe true + coordinateNotExactlyBetween.isBetween( + coordinateA, + coordinateB + ) shouldBe false + coordinateNotExactlyBetween.isBetween( + coordinateA, + coordinateB, + 1e-5 + ) shouldBe true + coordinateNotBetween.isBetween(coordinateB, coordinateA) shouldBe false + } + + "transform to a point correctly" in { + val coordinate = buildCoordinate(10, 5) + val actual = coordinate.toPoint + actual.getCoordinate.x shouldBe coordinate.getX +- 1e-12 + actual.getCoordinate.y shouldBe coordinate.getY +- 1e-12 + } + } + + "A GeoLinestring" should { + + "calculate the aggregated length correctly" in { + val coordinateA = buildCoordinate(10, 5) + val coordinateB = buildCoordinate(11, 6) + val coordinateC = buildCoordinate(12, 7) + val lineString = DEFAULT_GEOMETRY_FACTORY.createLineString( + Array(coordinateA, coordinateB, coordinateC) + ) + lineString.haversineLength should equalWithTolerance( + coordinateA + .haversineDistance(coordinateB) + .add(coordinateB.haversineDistance(coordinateC)) + ) + } + } + + "A GeoPolygon" should { + + "calculate the intersection between polygons correctly" in { + val aBottomLeft = buildCoordinate(0, 1) + val aTopLeft = buildCoordinate(1, 1) + val aTopRight = buildCoordinate(1, 3) + val aBottomRight = buildCoordinate(0, 3) + + val bBottomLeft = buildCoordinate(0, 2) + val bTopLeft = buildCoordinate(0.5, 2) + val bTopRight = buildCoordinate(0.5, 4) + val bBottomRight = buildCoordinate(0, 4) + + val polygonA = buildPolygon( + Array( + aBottomLeft, + aTopLeft, + aTopRight, + aBottomRight, + aBottomLeft + ) + ) + + val polygonB = buildPolygon( + Array( + bBottomLeft, + bTopLeft, + bTopRight, + bBottomRight, + bBottomLeft + ) + ) + + val expectedIntersection = Array( + buildCoordinate(0.5, 3), + bBottomLeft, + bTopLeft, + aBottomRight + ) + + polygonA.intersect(polygonB) match { + case Failure(exception) => throw exception + case Success(polygon) => + polygon.getCoordinates.foreach(intCoordinate => + expectedIntersection.exists(expected => + expected.equals2D(intCoordinate) + ) shouldBe true + ) + } + } "calculate area on earth correctly" in { val coordinateA = new Coordinate(8.748497631269068, 51.72341137638795) @@ -30,7 +153,5 @@ class RichGeometriesSpec extends Matchers with AnyWordSpecLike { .getValue .doubleValue() shouldBe 1d +- 0.01 } - } - } From 44f8190e1df6adfe9e609f539e06139dfa03f685 Mon Sep 17 00:00:00 2001 From: t-ober <63147366+t-ober@users.noreply.github.com> Date: Thu, 9 Dec 2021 14:36:51 +0100 Subject: [PATCH 11/17] neglect coordinate projection --- src/main/scala/edu/ie3/util/geo/GeoUtils.scala | 8 ++------ src/main/scala/edu/ie3/util/geo/RichGeometries.scala | 3 +-- src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala | 3 +-- 3 files changed, 4 insertions(+), 10 deletions(-) diff --git a/src/main/scala/edu/ie3/util/geo/GeoUtils.scala b/src/main/scala/edu/ie3/util/geo/GeoUtils.scala index 0cc1b10f..a6b973d0 100644 --- a/src/main/scala/edu/ie3/util/geo/GeoUtils.scala +++ b/src/main/scala/edu/ie3/util/geo/GeoUtils.scala @@ -193,8 +193,6 @@ object GeoUtils { } /** Builds a convex hull from a set of latitude/longitude coordinates. - * Projects them on a 2D-surface, calculates the convex hull and reverses the - * projection of resulting coordinates to latitude and longitude values. * * @param coordinates * the coordinates to consider @@ -202,15 +200,13 @@ object GeoUtils { * a Try of the resulting polygon */ def buildConvexHull(coordinates: Set[Coordinate]): Try[Polygon] = { - val projectedCoordinates = coordinates.map(equalAreaProjection) + val projectedCoordinates = coordinates new ConvexHull( projectedCoordinates.toArray, DEFAULT_GEOMETRY_FACTORY ).getConvexHull match { case projectedPolygon: Polygon => - val coordinates = - projectedPolygon.getCoordinates.map(reverseEqualAreaProjection) - Success(buildPolygon(coordinates)) + Success(buildPolygon(projectedPolygon.getCoordinates)) case _: LineString => Failure( new GeoException( diff --git a/src/main/scala/edu/ie3/util/geo/RichGeometries.scala b/src/main/scala/edu/ie3/util/geo/RichGeometries.scala index 8b06b1c5..95cab330 100644 --- a/src/main/scala/edu/ie3/util/geo/RichGeometries.scala +++ b/src/main/scala/edu/ie3/util/geo/RichGeometries.scala @@ -9,8 +9,7 @@ import edu.ie3.util.exceptions.GeoException import edu.ie3.util.geo.GeoUtils.{ DEFAULT_GEOMETRY_FACTORY, buildPolygon, - calcHaversine, - reverseEqualAreaProjection + calcHaversine } import edu.ie3.util.quantities.QuantityUtils.RichQuantityDouble import org.locationtech.jts.geom.{Coordinate, LineString, Point, Polygon} diff --git a/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala b/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala index c2bb5cf1..08fc4ce3 100644 --- a/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala +++ b/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala @@ -123,7 +123,7 @@ class GeoUtilsSpec extends Matchers with AnyWordSpecLike { hullCoordinates.size shouldBe 5 hullCoordinates.foreach(hullCoordinate => cornerCoordinates.exists(cornerCoordinate => - hullCoordinate.equals2D(cornerCoordinate, 1e-10) + hullCoordinate.equals2D(cornerCoordinate) ) shouldBe true ) } @@ -201,7 +201,6 @@ class GeoUtilsSpec extends Matchers with AnyWordSpecLike { polygon.containsCoordinate(new Coordinate(7.5, 49)) shouldBe true polygon.containsCoordinate(new Coordinate(7.5, 50.1)) shouldBe false polygon.containsCoordinate(new Coordinate(8.1, 50)) shouldBe false - } } } From d330b50cb08f466548c799a687015e5fe60f57cf Mon Sep 17 00:00:00 2001 From: t-ober <63147366+t-ober@users.noreply.github.com> Date: Thu, 9 Dec 2021 14:43:21 +0100 Subject: [PATCH 12/17] adapt changelog --- CHANGELOG.md | 8 ++++++-- 1 file changed, 6 insertions(+), 2 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 7a9d197a..6569c2c3 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -6,8 +6,12 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ## [Unreleased/Snapshot] ### Added -- Added scala support -- Added `RichQuantityDouble` as double type enrichment to enable easy quantity conversions [#133](https://github.com/ie3-institute/PowerSystemUtils/issues/133) +- Added scala support +- Added `RichQuantityDouble` as double type enrichment to enable easy quantity conversions [#133](https://github.com/ie3-institute/PowerSystemUtils/issues/133) +- Added implicit classes for `loactiontec.jts` Geometries that represent geographical geometries with functionality before located in `GeoUtils` [#163] (https://github.com/ie3-institute/PowerSystemUtils/issues/163) + +### Changed +- Refactored `GeoUtils`, moved them to the scala package and tailored them toward the `loactiontec.jts` Geometries used in the `OsmModel` [#163] (https://github.com/ie3-institute/PowerSystemUtils/issues/163) ## [1.5.3] ### Fixed From 4d995ef0a257e0bfecff08baace830d002521446 Mon Sep 17 00:00:00 2001 From: t-ober <63147366+t-ober@users.noreply.github.com> Date: Fri, 7 Jan 2022 08:32:05 +0100 Subject: [PATCH 13/17] Update CHANGELOG.md Co-authored-by: Chris Kittl <44838605+ckittl@users.noreply.github.com> --- CHANGELOG.md | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 6569c2c3..d8212f9b 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -11,7 +11,8 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 - Added implicit classes for `loactiontec.jts` Geometries that represent geographical geometries with functionality before located in `GeoUtils` [#163] (https://github.com/ie3-institute/PowerSystemUtils/issues/163) ### Changed -- Refactored `GeoUtils`, moved them to the scala package and tailored them toward the `loactiontec.jts` Geometries used in the `OsmModel` [#163] (https://github.com/ie3-institute/PowerSystemUtils/issues/163) +- Refactored `GeoUtils`, moved them to the scala package and tailored them toward the `loactiontec.jts` Geometries used in the `OsmContainer` [#163] (https://github.com/ie3-institute/PowerSystemUtils/issues/163) + ## [1.5.3] ### Fixed From c4d69657c024faf3a5eb8dd6864a8a7753d84348 Mon Sep 17 00:00:00 2001 From: t-ober <63147366+t-ober@users.noreply.github.com> Date: Fri, 7 Jan 2022 09:34:43 +0100 Subject: [PATCH 14/17] Update src/main/scala/edu/ie3/util/osm/OsmUtils.scala Co-authored-by: Chris Kittl <44838605+ckittl@users.noreply.github.com> --- src/main/scala/edu/ie3/util/osm/OsmUtils.scala | 8 ++------ 1 file changed, 2 insertions(+), 6 deletions(-) diff --git a/src/main/scala/edu/ie3/util/osm/OsmUtils.scala b/src/main/scala/edu/ie3/util/osm/OsmUtils.scala index f2bb4f72..8728064e 100644 --- a/src/main/scala/edu/ie3/util/osm/OsmUtils.scala +++ b/src/main/scala/edu/ie3/util/osm/OsmUtils.scala @@ -10,11 +10,7 @@ import org.locationtech.jts.geom.Point object OsmUtils { - def isInsideLandUse(coordinate: Point, landuses: List[ClosedWay]): Boolean = { - landuses.foreach(landuse => - if (landuse.toPolygon.convexHull().contains(coordinate)) return true - ) - false - } +def isInsideLandUse(coordinate: Point, landuses: List[ClosedWay]): Boolean = + landuses.exists(_.toPolygon.convexHull().contains(coordinate)) } From 62622b37652adbed4123b27268ce6b3a7cdc276f Mon Sep 17 00:00:00 2001 From: t-ober <63147366+t-ober@users.noreply.github.com> Date: Fri, 7 Jan 2022 09:47:42 +0100 Subject: [PATCH 15/17] address review comments --- .../edu/ie3/util/geo/CoordinateDistance.scala | 38 +++++++-------- .../scala/edu/ie3/util/geo/GeoUtils.scala | 47 +++++++++---------- .../edu/ie3/util/geo/RichGeometries.scala | 31 +++++------- .../scala/edu/ie3/util/osm/OsmEntities.scala | 2 +- .../scala/edu/ie3/util/osm/OsmUtils.scala | 11 ++++- .../scala/edu/ie3/util/geo/GeoUtilsSpec.scala | 2 +- .../edu/ie3/util/geo/RichGeometriesSpec.scala | 6 +-- 7 files changed, 63 insertions(+), 74 deletions(-) diff --git a/src/main/scala/edu/ie3/util/geo/CoordinateDistance.scala b/src/main/scala/edu/ie3/util/geo/CoordinateDistance.scala index cfdda498..67c89d1c 100644 --- a/src/main/scala/edu/ie3/util/geo/CoordinateDistance.scala +++ b/src/main/scala/edu/ie3/util/geo/CoordinateDistance.scala @@ -26,27 +26,6 @@ class CoordinateDistance private ( val distance: ComparableQuantity[Length] ) extends Comparable[CoordinateDistance] { - /** Calculates the distance from the first to the second coordinate using - * [[GeoUtils.calcHaversine(double, double, double, double)]] - * - * @param from - * The first coordinate - * @param to - * The second coordinate - */ - def this(from: Point, to: Point) = { - this( - from, - to, - GeoUtils.calcHaversine( - from.getY, - from.getX, - to.getY, - to.getX - ) - ) - } - /** The first coordinate */ def getFrom: Point = { from @@ -94,3 +73,20 @@ class CoordinateDistance private ( "CoordinateDistance{" + "coordinateA=" + from + ", coordinateB=" + to + ", distance=" + distance + '}' } } + +object CoordinateDistance { + + def apply(from: Point, to: Point): CoordinateDistance = { + new CoordinateDistance( + from, + to, + GeoUtils.calcHaversine( + from.getY, + from.getX, + to.getY, + to.getX + ) + ) + } + +} diff --git a/src/main/scala/edu/ie3/util/geo/GeoUtils.scala b/src/main/scala/edu/ie3/util/geo/GeoUtils.scala index a6b973d0..974ad0ec 100644 --- a/src/main/scala/edu/ie3/util/geo/GeoUtils.scala +++ b/src/main/scala/edu/ie3/util/geo/GeoUtils.scala @@ -6,7 +6,7 @@ package edu.ie3.util.geo import edu.ie3.util.exceptions.GeoException -import org.apache.commons.lang3.ArrayUtils +import edu.ie3.util.geo.RichGeometries.RichCoordinate import org.locationtech.jts.algorithm.ConvexHull import org.locationtech.jts.geom.impl.CoordinateArraySequence import org.locationtech.jts.geom.{ @@ -82,10 +82,7 @@ object GeoUtils { * a [[LineString]] between the provided points */ def buildSafeLineStringBetweenPoints(p1: Point, p2: Point): LineString = { - val safePoint1 = if (p1 == p2) buildSafePoint(p1) else p1 - DEFAULT_GEOMETRY_FACTORY.createLineString( - safePoint1.getCoordinates ++ p2.getCoordinates - ) + buildSafeLineStringBetweenCoords(p1.getCoordinate, p2.getCoordinate) } /** Build an instance of [[LineString]] between two coordinates that is safe @@ -108,24 +105,10 @@ object GeoUtils { ): LineString = { val safeCoord1: Coordinate = if (c1 == c2) buildSafeCoord(c1) else c1 DEFAULT_GEOMETRY_FACTORY.createLineString( - ArrayUtils.addAll(Array[Coordinate](safeCoord1), c2) + Array(safeCoord1, c2) ) } - /** Adapt the provided point as described in buildSafeCoord ( Coordinate )} - * and return a new, adapted instance of [[Point]] - * - * @param p1 - * the point that should be adapted - * @return - * the adapted point with a slightly changed coordinate - */ - private def buildSafePoint(p1: Point): Point = { - val safeCoord = buildSafeCoord(p1.getCoordinate) - val safeCoordSeq = new CoordinateArraySequence(Array[Coordinate](safeCoord)) - new Point(safeCoordSeq, p1.getFactory) - } - /** Adapted [[Coordinate]] x, [[Coordinate]] y and [[Coordinate]] z of the * provided [[Coordinate]] by 1e-13 and return a new, adapted instance of * [[Coordinate]] @@ -143,6 +126,16 @@ object GeoUtils { new Coordinate(p1X, p1Y, p1Z) } + @deprecated( + "This method will be dropped and replaced by #calcOrderedCoordinateDistances" + ) + def getCoordinateDistances( + baseCoordinate: Point, + coordinates: Array[Point] + ): SortedSet[CoordinateDistance] = { + calcOrderedCoordinateDistances(baseCoordinate, coordinates) + } + /** Calculates and orders the coordinate distances from a base coordinate to a * list of coordinates * @@ -158,7 +151,7 @@ object GeoUtils { coordinates: Array[Point] ): SortedSet[CoordinateDistance] = { val coordinateDistances = coordinates.map(coordinate => - new CoordinateDistance(baseCoordinate, coordinate) + CoordinateDistance(baseCoordinate, coordinate) ) TreeSet(coordinateDistances: _*) } @@ -200,9 +193,8 @@ object GeoUtils { * a Try of the resulting polygon */ def buildConvexHull(coordinates: Set[Coordinate]): Try[Polygon] = { - val projectedCoordinates = coordinates new ConvexHull( - projectedCoordinates.toArray, + coordinates.toArray, DEFAULT_GEOMETRY_FACTORY ).getConvexHull match { case projectedPolygon: Polygon => @@ -230,11 +222,14 @@ object GeoUtils { } def buildPoint(lat: Double, long: Double): Point = { - buildPoint(buildCoordinate(lat, long)) + buildCoordinate(lat, long).toPoint } + @deprecated( + "Will be dropped. Please use implementation in rich class GeoCoordinate." + ) def buildPoint(coordinate: Coordinate): Point = { - DEFAULT_GEOMETRY_FACTORY.createPoint(coordinate) + coordinate.toPoint } /** Build a coordinate from latitude and longitude values. @@ -260,7 +255,7 @@ object GeoUtils { * a [[Polygon]] */ def buildPolygon(coordinates: Array[Coordinate]): Polygon = { - val arrayCoordinates = new CoordinateArraySequence(coordinates.toArray) + val arrayCoordinates = new CoordinateArraySequence(coordinates) val linearRing = new LinearRing(arrayCoordinates, DEFAULT_GEOMETRY_FACTORY) new Polygon(linearRing, Array[LinearRing](), DEFAULT_GEOMETRY_FACTORY) diff --git a/src/main/scala/edu/ie3/util/geo/RichGeometries.scala b/src/main/scala/edu/ie3/util/geo/RichGeometries.scala index 95cab330..a972711a 100644 --- a/src/main/scala/edu/ie3/util/geo/RichGeometries.scala +++ b/src/main/scala/edu/ie3/util/geo/RichGeometries.scala @@ -14,8 +14,6 @@ import edu.ie3.util.geo.GeoUtils.{ import edu.ie3.util.quantities.QuantityUtils.RichQuantityDouble import org.locationtech.jts.geom.{Coordinate, LineString, Point, Polygon} import tech.units.indriya.ComparableQuantity -import tech.units.indriya.quantity.Quantities -import tech.units.indriya.unit.Units.METRE import javax.measure.quantity.{Area, Length} import scala.math.abs @@ -23,7 +21,7 @@ import scala.util.{Failure, Success, Try} object RichGeometries { - implicit class GeoCoordinate(coordinate: Coordinate) { + implicit class RichCoordinate(coordinate: Coordinate) { /** Calculates the great circle distance between two coordinates * @@ -67,10 +65,7 @@ object RichGeometries { .add(coordinate.haversineDistance(b)) .getValue .doubleValue - if ( - abs(1 - (distancePassingMe / distance.getValue.doubleValue())) < epsilon - ) true - else false + abs(1 - (distancePassingMe / distance.getValue.doubleValue())) < epsilon } def toPoint: Point = { @@ -78,7 +73,7 @@ object RichGeometries { } } - implicit class GeoLineString(lineString: LineString) { + implicit class RichLineString(lineString: LineString) { /** Compute length of a [[LineString]] on earth's surface. * @@ -87,21 +82,17 @@ object RichGeometries { */ def haversineLength: ComparableQuantity[Length] = { val coordinates = lineString.getCoordinates.toVector - val coordinateSize = coordinates.size - coordinates.zipWithIndex - .foldLeft(Quantities.getQuantity(0, METRE))((acc, current) => { - if (current._2 < coordinateSize - 1) { - val currentCoordinate = current._1 - val nextCoordinate = coordinates(current._2 + 1) - acc.add(currentCoordinate.haversineDistance(nextCoordinate)) - } else acc - }) + val coordinatePairs = coordinates.init.zip(coordinates.tail) + coordinatePairs + .map { case (lhs, rhs) => + lhs.haversineDistance(rhs) + } + .reduce(_.add(_)) } } - implicit class GeoPolygon(polygon: Polygon) { + implicit class RichPolygon(polygon: Polygon) { - // Fixme: Is this correct or should the function calculate the specific coordinates shared by the polygons /** Calculates intersection between polygons * * @param polygonB @@ -110,7 +101,7 @@ object RichGeometries { */ def intersect(polygonB: Polygon): Try[Polygon] = { Try( - polygon.intersection(polygonB.equalAreaProjection) + polygon.intersection(polygonB) ) match { case Failure(exception) => Failure( diff --git a/src/main/scala/edu/ie3/util/osm/OsmEntities.scala b/src/main/scala/edu/ie3/util/osm/OsmEntities.scala index d56251cc..1bac833c 100644 --- a/src/main/scala/edu/ie3/util/osm/OsmEntities.scala +++ b/src/main/scala/edu/ie3/util/osm/OsmEntities.scala @@ -6,7 +6,7 @@ package edu.ie3.util.osm import edu.ie3.util.geo.GeoUtils.buildPolygon -import edu.ie3.util.geo.RichGeometries.GeoPolygon +import edu.ie3.util.geo.RichGeometries.RichPolygon import org.locationtech.jts.geom.{Coordinate, Polygon} import tech.units.indriya.ComparableQuantity diff --git a/src/main/scala/edu/ie3/util/osm/OsmUtils.scala b/src/main/scala/edu/ie3/util/osm/OsmUtils.scala index 8728064e..1d5f6584 100644 --- a/src/main/scala/edu/ie3/util/osm/OsmUtils.scala +++ b/src/main/scala/edu/ie3/util/osm/OsmUtils.scala @@ -5,12 +5,19 @@ */ package edu.ie3.util.osm +import edu.ie3.util.geo.RichGeometries.RichCoordinate import edu.ie3.util.osm.OsmEntities.ClosedWay -import org.locationtech.jts.geom.Point +import org.locationtech.jts.geom.{Coordinate, Point} object OsmUtils { -def isInsideLandUse(coordinate: Point, landuses: List[ClosedWay]): Boolean = + def isInsideLandUse(coordinate: Point, landuses: List[ClosedWay]): Boolean = landuses.exists(_.toPolygon.convexHull().contains(coordinate)) + def isInsideLandUse( + coordinate: Coordinate, + landuses: List[ClosedWay] + ): Boolean = + isInsideLandUse(coordinate.toPoint, landuses) + } diff --git a/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala b/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala index 08fc4ce3..214ad14f 100644 --- a/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala +++ b/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala @@ -17,7 +17,7 @@ import edu.ie3.util.geo.GeoUtils.{ equalAreaProjection, reverseEqualAreaProjection } -import edu.ie3.util.geo.RichGeometries.{GeoCoordinate, GeoPolygon} +import edu.ie3.util.geo.RichGeometries.{RichCoordinate, RichPolygon} import edu.ie3.util.quantities.QuantityMatchers.equalWithTolerance import org.locationtech.jts.geom.Coordinate import org.scalatest.matchers.should.Matchers diff --git a/src/test/scala/edu/ie3/util/geo/RichGeometriesSpec.scala b/src/test/scala/edu/ie3/util/geo/RichGeometriesSpec.scala index 5d0e136f..d98bf95b 100644 --- a/src/test/scala/edu/ie3/util/geo/RichGeometriesSpec.scala +++ b/src/test/scala/edu/ie3/util/geo/RichGeometriesSpec.scala @@ -12,9 +12,9 @@ import edu.ie3.util.geo.GeoUtils.{ calcHaversine } import edu.ie3.util.geo.RichGeometries.{ - GeoCoordinate, - GeoLineString, - GeoPolygon + RichCoordinate, + RichLineString, + RichPolygon } import edu.ie3.util.quantities.QuantityMatchers.equalWithTolerance import edu.ie3.util.quantities.QuantityUtils.RichQuantityDouble From ceaa9f4ff67b63cea1caad6659b055aaa2bc4f3a Mon Sep 17 00:00:00 2001 From: t-ober <63147366+t-ober@users.noreply.github.com> Date: Mon, 10 Jan 2022 08:39:42 +0100 Subject: [PATCH 16/17] fmt --- src/main/java/edu/ie3/util/geo/GeoUtils.java | 0 src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala | 1 - 2 files changed, 1 deletion(-) delete mode 100644 src/main/java/edu/ie3/util/geo/GeoUtils.java diff --git a/src/main/java/edu/ie3/util/geo/GeoUtils.java b/src/main/java/edu/ie3/util/geo/GeoUtils.java deleted file mode 100644 index e69de29b..00000000 diff --git a/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala b/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala index 214ad14f..235571a6 100644 --- a/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala +++ b/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala @@ -9,7 +9,6 @@ import edu.ie3.util.geo.GeoUtils.{ DEFAULT_GEOMETRY_FACTORY, buildCirclePolygon, buildConvexHull, - buildCoordinate, buildPoint, buildPolygon, calcHaversine, From df62d17722403992c201ea28effa3efe9ab7f79a Mon Sep 17 00:00:00 2001 From: t-ober <63147366+t-ober@users.noreply.github.com> Date: Mon, 10 Jan 2022 09:09:11 +0100 Subject: [PATCH 17/17] use wildcard imports --- src/main/scala/edu/ie3/util/geo/GeoUtils.scala | 11 +---------- .../edu/ie3/util/quantities/QuantityUtils.scala | 10 +--------- src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala | 12 +----------- .../scala/edu/ie3/util/geo/RichGeometriesSpec.scala | 3 +-- 4 files changed, 4 insertions(+), 32 deletions(-) diff --git a/src/main/scala/edu/ie3/util/geo/GeoUtils.scala b/src/main/scala/edu/ie3/util/geo/GeoUtils.scala index 974ad0ec..28c3d358 100644 --- a/src/main/scala/edu/ie3/util/geo/GeoUtils.scala +++ b/src/main/scala/edu/ie3/util/geo/GeoUtils.scala @@ -9,16 +9,7 @@ import edu.ie3.util.exceptions.GeoException import edu.ie3.util.geo.RichGeometries.RichCoordinate import org.locationtech.jts.algorithm.ConvexHull import org.locationtech.jts.geom.impl.CoordinateArraySequence -import org.locationtech.jts.geom.{ - Coordinate, - GeometryCollection, - GeometryFactory, - LineString, - LinearRing, - Point, - Polygon, - PrecisionModel -} +import org.locationtech.jts.geom._ import tech.units.indriya.ComparableQuantity import tech.units.indriya.quantity.Quantities import tech.units.indriya.unit.Units.{METRE, RADIAN} diff --git a/src/main/scala/edu/ie3/util/quantities/QuantityUtils.scala b/src/main/scala/edu/ie3/util/quantities/QuantityUtils.scala index 701a0a19..896ef0e5 100644 --- a/src/main/scala/edu/ie3/util/quantities/QuantityUtils.scala +++ b/src/main/scala/edu/ie3/util/quantities/QuantityUtils.scala @@ -24,15 +24,7 @@ import edu.ie3.util.quantities.interfaces.{ } import tech.units.indriya.ComparableQuantity import tech.units.indriya.quantity.Quantities -import tech.units.indriya.unit.Units.{ - AMPERE, - OHM, - PERCENT, - SIEMENS, - SQUARE_METRE, - VOLT, - WATT -} +import tech.units.indriya.unit.Units._ import javax.measure.MetricPrefix import javax.measure.quantity.{ diff --git a/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala b/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala index 235571a6..f19e6b4d 100644 --- a/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala +++ b/src/test/scala/edu/ie3/util/geo/GeoUtilsSpec.scala @@ -5,17 +5,7 @@ */ package edu.ie3.util.geo -import edu.ie3.util.geo.GeoUtils.{ - DEFAULT_GEOMETRY_FACTORY, - buildCirclePolygon, - buildConvexHull, - buildPoint, - buildPolygon, - calcHaversine, - calcOrderedCoordinateDistances, - equalAreaProjection, - reverseEqualAreaProjection -} +import edu.ie3.util.geo.GeoUtils._ import edu.ie3.util.geo.RichGeometries.{RichCoordinate, RichPolygon} import edu.ie3.util.quantities.QuantityMatchers.equalWithTolerance import org.locationtech.jts.geom.Coordinate diff --git a/src/test/scala/edu/ie3/util/geo/RichGeometriesSpec.scala b/src/test/scala/edu/ie3/util/geo/RichGeometriesSpec.scala index d98bf95b..e477220b 100644 --- a/src/test/scala/edu/ie3/util/geo/RichGeometriesSpec.scala +++ b/src/test/scala/edu/ie3/util/geo/RichGeometriesSpec.scala @@ -8,8 +8,7 @@ package edu.ie3.util.geo import edu.ie3.util.geo.GeoUtils.{ DEFAULT_GEOMETRY_FACTORY, buildCoordinate, - buildPolygon, - calcHaversine + buildPolygon } import edu.ie3.util.geo.RichGeometries.{ RichCoordinate,