5

I have a number of points which make getOrthodromicDistance method to fail with exception in geotools lib, while these points are valid lat lon points:

Point which throws the exception (lat,lon):

val p1= (5.318765,-75.786109)
val p2= (-6.32907,106.09254)

eg exception: No convergence for points 75°47,2'W 06°19,7'S and 106°05,6'E 05°19,1'N. java.lang.ArithmeticException: No convergence for points 75°47,2'W 06°19,7'S and 106°05,6'E 05°19,1'N. at org.geotools.referencing.GeodeticCalculator.computeDirection(GeodeticCalculator.java:1073)

Code used in Scala:

  def latlonDistance(p1:(Double,Double), p2:(Double,Double)):Double={
      val world= new GeodeticCalculator()
      world.setStartingGeographicPoint(p1._2, p2._1)
      world.setDestinationGeographicPoint(p2._2, p1._1)
      world.getOrthodromicDistance
   }

Note: My point format passed in latlonDistance is (lat,lon) as mentioned above while, setStartingGeographicPoint, setDestinationGeographicPoint need (lon,lat) order.

Version used:

        <dependency>
          <groupId>org.geotools</groupId>
          <artifactId>gt-referencing</artifactId>
          <version>13.2</version>
        </dependency>

In python works as expected:

>>> from geopy.distance import vincenty
>>> pt1= [5.318765,-75.786109]
>>> pt2= [-6.32907,106.09254]
>>> vincenty(pt1 , pt2)
Distance(19791.6883647)

It is the orthodromicDistance method in org.geotools.referencing.datum.DefaultEllipsoid which does not converge. Any workarounds?

skonto
  • 91
  • 1
  • 7

1 Answers1

3

The problem is that this is not a straightforward calculation as the Vincenty algorithm is an iterative process and some sets of points don't necessarily converge (within the limit set).

There are two possible solutions 1 - edit the GeodeticCalculator to increase the number of possible iterations from 12 to 15, this works in this case but I can't guarantee it in others. Or 2 use another algorithm, following links from this question's answers I found the GeographicLib library at Sourceforge and used that instead with your points. It is written by the author (@cffk) of the paper linked to in the other answer.

For your points it gives a very reasonable looking 20004 Km.

Ian Turton
  • 10,018
  • 1
  • 28
  • 47
  • I've raised a bug report (https://osgeo-org.atlassian.net/browse/GEOT-5191) if you want to watch it – Ian Turton Aug 12 '15 at 14:46
  • Yes i investigated it as well GeographicLib library is a candidate solution it works fine because it has a different algorithm implemented, one that converges for any pair of points while Vincenty by design does not. I will watch the bug thnx! – skonto Aug 13 '15 at 15:30
  • Please find a set of 130 pair of points which also could serve as a test case [bad points](https://www.dropbox.com/s/39jaqaitckzjcja/badpoints.txt?dl=0) – skonto Aug 13 '15 at 15:39
  • 1
    I've submitted this [pull request](https://github.com/geotools/geotools/pull/945) for geotools to address this issue. Please try it out. (Currently, the build checks are failing; but I think that's unrelated to the patch.) – cffk Aug 23 '15 at 13:33
  • 1
    This will no longer be a problem for GeoTools 14+ thanks to @cffk – Ian Turton Aug 30 '15 at 16:29