[Proj] Using Proj.4 to compress a 10-parameter Molodensky-Badekas datum shift into 7 parameters
ovv at hetnet.nl
Sun Oct 17 03:15:05 EST 2010
IN REPLY TO:
[Proj] Using Proj.4 to compress a 10-parameter Molodensky-Badekasdatum shift
into 7 parameters
Noel Zinn (cc) ndzinn at comcast.net
Sat Oct 16 09:42:42 EST 2010
> I haven't tested the geocentric -> geodetic code from PROJ4's geocent.c,
> it certainly looks good.
I contend that the correct output (delivered by my Matlab scripts and by
Blue Marble) is
=> 2464278.090 -5783490.396 974642.386 (to the same precision)
If my assertion is correct, then Proj.4 has a slight problem. It's not
blame, just a request for an independent verification of this particular
step (i) and (if verified) another possible improvement for Proj.4.
I do get exactly the same output as stated above, so up to that point there
is nothing wrong I guess.
Let's look at PROJ4's geocent. I took the code apart and tested both
Indeed, PROJ4's geocent.c has two functions. Which one is actually used,
depends on a compiler directive.
The first code is from Ralph Toms. This code is optimized for speed. It uses
a small fudge factor (AD_C). This introduces a small error, but it enhances
the error stability especially for extremely high altitudes.
The second code is an iterative procedure, which gives nearly exact results
(the iteration cut-off error is 1e-12).
Let's see. In the problem at hand my Molodenksy-Badekas routine arrives at
an X,Y,Z of:
X = 2464278.0897929626; Y = -5783490.396202258; Z = 974642.3859339222;
Converted to lat, lon, height on WGS84 with the exact Borkowski procedure
9.580277979760344, -66.08186260475033, 180.5140539117856
Converted with the iterative procedure from geocent.c gives in my test:
9.580277979760346, -66.08186260475033, 180.514053911902
Converted with the Toms procedure from geocent.c gives in my test:
9.580277981233836, -66.08186260475033, 180.51408141944557
The latter test differs by 5.3 micro-arcseconds in the latitude with the
Borkowski method and 27.5 micrometers in the height.
I believe PROJ4's geocent uses by default the somewhat inaccurate Toms
method for converting geocentric to geodetic.
But at least in my test it seems that the error is rather insignificantly
So the core of geocent.c seems OK.
Are the global structs filled correctly and are all parameters passed
Oscar van Vlijmen
More information about the Proj