[Proj] Using Proj.4 to compress a 10-parameter Molodensky-Badekas datum shift into 7 parameters

OvV_HN ovv at hetnet.nl
Sun Oct 17 03:15:05 EST 2010


[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,
> but
> 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 
geocentric->geodetic functions.
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 mailing list