From 42d8b865ffbe997f478be637f86c1eeb218d39ed Mon Sep 17 00:00:00 2001 From: Romain Gallet Date: Sat, 24 Feb 2018 23:58:10 +0000 Subject: [PATCH] Added midpoint --- README.md | 21 ++++++---- src/main/java/com/grum/geocalc/EarthCalc.java | 40 ++++++++++++++----- .../java/com/grum/geocalc/DistanceTest.java | 33 +++++++++++---- 3 files changed, 68 insertions(+), 26 deletions(-) diff --git a/README.md b/README.md index a0e5f03..0ec49b6 100644 --- a/README.md +++ b/README.md @@ -2,13 +2,15 @@ Java Geocalc ![alt text](https://api.travis-ci.org/grumlimited/geocalc.svg?branc ======= Geocalc is a simple java library aimed at doing arithmetics with Earth coordinates. -It is designed to be simple to embed in your existing applications and easy to use. Geocalc can: +It is designed to be simple to embed in your existing applications and easy to use. -1. Calculate the distance between two coordinates +Geocalc can: + +1. Calculate the distance between two coordinates (law of cosines, harvesine and vincenty) 2. Find a point at X distance from a standpoint, given a bearing -3. Calculate a rectangular area around a point +3. Calculate coordinates of a rectangular area around a point 4. Determine whether a Point is contained within that area -5. Calculate the azimuth and final bearings between two points +5. Calculate the azimuth, initial and final bearings between two points (vincenty) This library is being used on [www.rentbarometer.com](http://www.rentbarometer.com). @@ -17,14 +19,17 @@ This library implements in Java lots of ideas from [Movable-Type](http://www.mov ## Changelog ### 0.5.4 +* added `EarthCalc.midPoint(p1, p2)` +* renamed `EarthCalc.boundingArea()` to `EarthCalc.pointAt(...)` +* renamed `BoundingArea.boundingArea()` to `BoundingArea.around(...)` * Removed `get...()` out of `Point` and `BoundingArea` -* Added maven javadoc plugin +* added maven javadoc plugin ### 0.5.3 -* Changed constructors to default and private visibility -* Removed `get...()` keyword out of `EarthCalc` methods +* changed constructors to default and private visibility +* removed `get...()` keyword out of `EarthCalc` methods * `EarthCalc.getDistance()` is now `EarthCalc.gcdDistance()` -* Renamed `BoundingArea.isContainedWithin(...)` to `BoundingArea.contains(...)` +* renamed `BoundingArea.isContainedWithin(...)` to `BoundingArea.contains(...)` ### Embed diff --git a/src/main/java/com/grum/geocalc/EarthCalc.java b/src/main/java/com/grum/geocalc/EarthCalc.java index 620409c..61cf942 100644 --- a/src/main/java/com/grum/geocalc/EarthCalc.java +++ b/src/main/java/com/grum/geocalc/EarthCalc.java @@ -43,6 +43,30 @@ public class EarthCalc { public static final double EARTH_DIAMETER = 6371.01 * 1000; //meters + /** + * This is the half-way point along a great circle path between the two points. + * + * @param standPoint standPoint + * @param forePoint standPoint + * @return mid point + * @see + */ + public static Point midPoint(Point standPoint, Point forePoint) { + double λ1 = toRadians(standPoint.longitude); + double λ2 = toRadians(forePoint.longitude); + + double φ1 = toRadians(standPoint.latitude); + double φ2 = toRadians(forePoint.latitude); + + double Bx = cos(φ2) * cos(λ2 - λ1); + double By = cos(φ2) * sin(λ2 - λ1); + + double φ3 = atan2(sin(φ1) + sin(φ2), sqrt((cos(φ1) + Bx) * (cos(φ1) + Bx) + By * By)); + double λ3 = λ1 + atan2(By, cos(φ1) + Bx); + + return Point.at(Coordinate.fromRadians(φ3), Coordinate.fromRadians(λ3)); + } + /** * Returns the distance between two points at spherical law of cosines. * @@ -154,7 +178,6 @@ public static double vincentyDistance(Point standPoint, Point forePoint) { * @param standPoint The stand point * @param forePoint The fore point * @return (azimuth) bearing in degrees to the North - * * @see */ public static double vincentyBearing(Point standPoint, Point forePoint) { @@ -167,7 +190,6 @@ public static double vincentyBearing(Point standPoint, Point forePoint) { * @param standPoint The stand point * @param forePoint The fore point * @return (azimuth) bearing in degrees to the North - * * @see */ public static double getVincentyFinalBearing(Point standPoint, Point forePoint) { @@ -184,10 +206,9 @@ public static double getVincentyFinalBearing(Point standPoint, Point forePoint) * @param bearing Direction in degrees * @param distance distance in meters * @return forePoint coordinates - * * @see */ - public static Point pointRadialDistance(Point standPoint, double bearing, double distance) { + public static Point pointAt(Point standPoint, double bearing, double distance) { /** var φ2 = asin( sin(φ1)*cos(d/R) + cos(φ1)*sin(d/R)*cos(brng) ); var λ2 = λ1 + atan2(sin(brng)*sin(d/R)*cos(φ1), cos(d/R)-sin(φ1)*sin(φ2)); @@ -201,7 +222,7 @@ public static Point pointRadialDistance(Point standPoint, double bearing, double double lat2 = asin(sin(rlat1) * cos(rdistance) + cos(rlat1) * sin(rdistance) * cos(rbearing)); double lon2 = rlon1 + atan2(sin(rbearing) * sin(rdistance) * cos(rlat1), cos(rdistance) - sin(rlat1) * sin(lat2)); - return Point.at(new RadianCoordinate(lat2), new RadianCoordinate(lon2)); + return Point.at(Coordinate.fromRadians(lat2), Coordinate.fromRadians(lon2)); } /** @@ -231,16 +252,15 @@ public static double bearing(Point standPoint, Point forePoint) { * @param standPoint The centre of the area * @param distance Distance around standPoint, im meters * @return The area - * * @see */ - public static BoundingArea boundingArea(Point standPoint, double distance) { + public static BoundingArea around(Point standPoint, double distance) { //45 degrees going north-west - Point northWest = pointRadialDistance(standPoint, 45, distance); + Point northWest = pointAt(standPoint, 45, distance); //225 degrees going south-east - Point southEast = pointRadialDistance(standPoint, 225, distance); + Point southEast = pointAt(standPoint, 225, distance); return new BoundingArea(northWest, southEast); } @@ -251,7 +271,7 @@ private static class Vincenty { * initialBearing is the initial bearing, or forward azimuth (in reference to North point), in degrees * finalBearing is the final bearing (in direction p1→p2), in degrees */ - double distance, initialBearing, finalBearing; + final double distance, initialBearing, finalBearing; public Vincenty(double distance, double initialBearing, double finalBearing) { this.distance = distance; diff --git a/src/test/java/com/grum/geocalc/DistanceTest.java b/src/test/java/com/grum/geocalc/DistanceTest.java index 4fc5eb9..ec51aa9 100644 --- a/src/test/java/com/grum/geocalc/DistanceTest.java +++ b/src/test/java/com/grum/geocalc/DistanceTest.java @@ -135,7 +135,7 @@ public void testBoundingAreaDistance() { Coordinate lng = Coordinate.fromDegrees(-0.2912044); Point kew = Point.at(lat, lng); - BoundingArea area = EarthCalc.boundingArea(kew, 3000); + BoundingArea area = EarthCalc.around(kew, 3000); double northEastDistance = EarthCalc.gcdDistance(kew, area.northEast); logger.info("North East => " + northEastDistance); @@ -187,7 +187,7 @@ public void testBoundingAreaNorthPole() { Coordinate lng = Coordinate.fromDegrees(0); Point northPole = Point.at(lat, lng); - BoundingArea area = EarthCalc.boundingArea(northPole, 10000); + BoundingArea area = EarthCalc.around(northPole, 10000); logger.info("North East => " + area.northEast); logger.info("South West => " + area.southWest); @@ -205,7 +205,7 @@ public void testBoundingAreaNextToLondon() { Coordinate lng = Coordinate.fromDegrees(-0.1997387000000117); Point northPole = Point.at(lat, lng); - BoundingArea area = EarthCalc.boundingArea(northPole, 5); + BoundingArea area = EarthCalc.around(northPole, 5); logger.info("North East => " + area.northEast); logger.info("South West => " + area.southWest); @@ -223,15 +223,15 @@ public void testPointRadialDistanceZero() { Coordinate lng = Coordinate.fromDegrees(-0.2912044); Point kew = Point.at(lat, lng); - Point sameKew = EarthCalc.pointRadialDistance(kew, 45, 0); + Point sameKew = EarthCalc.pointAt(kew, 45, 0); assertEquals(lat.degrees(), sameKew.latitude, 1E-10); assertEquals(lng.degrees(), sameKew.longitude, 1E-10); - sameKew = EarthCalc.pointRadialDistance(kew, 90, 0); + sameKew = EarthCalc.pointAt(kew, 90, 0); assertEquals(lat.degrees(), sameKew.latitude, 1E-10); assertEquals(lng.degrees(), sameKew.longitude, 1E-10); - sameKew = EarthCalc.pointRadialDistance(kew, 180, 0); + sameKew = EarthCalc.pointAt(kew, 180, 0); assertEquals(lat.degrees(), sameKew.latitude, 1E-10); assertEquals(lng.degrees(), sameKew.longitude, 1E-10); } @@ -251,7 +251,7 @@ public void testPointRadialDistance() { double distance = EarthCalc.gcdDistance(kew, richmond); double bearing = EarthCalc.bearing(kew, richmond); - Point allegedRichmond = EarthCalc.pointRadialDistance(kew, bearing, distance); + Point allegedRichmond = EarthCalc.pointAt(kew, bearing, distance); assertEquals(richmond.latitude, allegedRichmond.latitude, 10E-5); assertEquals(richmond.longitude, allegedRichmond.longitude, 10E-5); @@ -315,8 +315,25 @@ public void testVincentyBearing() { lng = Coordinate.fromDegrees(-0.3035466); Point richmond = Point.at(lat, lng); - //comparing to results from ttp://www.movable-type.co.uk/scripts/latlong.html + //comparing to results from http://www.movable-type.co.uk/scripts/latlong.html assertEquals(EarthCalc.vincentyBearing(kew, richmond), new DMSCoordinate(198, 30, 19.58).degrees(), 10E-5); assertEquals(EarthCalc.getVincentyFinalBearing(kew, richmond), new DMSCoordinate(198, 29, 44.82).degrees(), 10E-5); } + + @Test + public void testMidPoint() { + //Kew + Coordinate lat = Coordinate.fromDegrees(51.4843774); + Coordinate lng = Coordinate.fromDegrees(-0.2912044); + Point kew = Point.at(lat, lng); + + //Richmond, London + lat = Coordinate.fromDegrees(51.4613418); + lng = Coordinate.fromDegrees(-0.3035466); + Point richmond = Point.at(lat, lng); + + //comparing to results from http://www.movable-type.co.uk/scripts/latlong.html + assertEquals(EarthCalc.midPoint(richmond, kew), Point.at(Coordinate.fromDegrees(51.47285976194266), Coordinate.fromDegrees(-0.2973770580524634))); + } + }