NaN result for inverse transform on ellipsoid with very close points

[nnseva]

I've passed very close coordinates with difference at fifth digit after dp in seconds to invgeod and received a NaN result: echo 56d16'19.00001 43d58'35.00001 56d16'19.00002 43d58'35.00002 | invgeod +ellps=WGS84 - -2147483648d-2147483648'-nan" -2147483648d-2147483648'-nan" -nan

What I am expecting: Basing on passed values it should be 0,0,0

Really the invgeod filters input (see # define DTOL at the top of the geod_inv.c file), but it is not enough as we can see.

The problem has been initially found and reported to pyproj project v.1.2.9/trunk. It also has been reproduced on invgeod v.4.7.1. As I can see in the geod_inv.c source, it should still remain in trunk version.

[nnseva]

I suppose that the DTOL value increased to 1e-8 works the problem around. The zero result then is returned till about 0.15 meters between points in my case, I think this is good enough (results below are from pyproj compiled with fixed DTOL):

>>> geod.inv(43.9763413300,56.2720659176,43.9763419000,56.272067000) (0.0, 0.0, 0.0) >>> geod.inv(43.9763413300,56.2720659176,43.976342000,56.272067000) (19.00427040948362, -160.99572903328846, 0.13428452977133082)

BUT: The strong mathematical base should be used here to substantiate input filter algorithm instead of empiric like this.

[warmerdam]

Note that I'm not overly enthusiastic about working on the geodesic code in PROJ.4. There are better implementations elsewhere. See GeodesicCalculations

I'll leave this ticket open in case someone has a fix to offer, but I'm not planning to look into the issue.