[Proj] Using Proj.4 to compress a 10-parameter Molodensky-Badekas datum shift into 7 parameters
ovv at hetnet.nl
Sun Oct 17 10:16:49 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
Sun Oct 17 07:44:35 EST 2010
To help me understand whether we agree of not, let me decompose Mikael
Rittri's problem simply.
A 7-parameter datum transformation (geodetic <=> geodetic) has 3 parts: (1)
geodetic => geocentric, (2) matrix manipulation on 3D geocentric (ECEF)
coordinates involving translation, rotation and scale change, and then (3)
geocentric => geodetic (which is the only inexact piece of the process, but
you can get any precision you desire).
> 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;
So, we agree then? The Proj.4 geocentric => geocentric routine (cs2cs -f
"%.3f" +proj=geocent +towgs84=0,0,0,5.226,1.238,-2.381,-5.109 +to
+datum=WGS84 +proj=geocent) has a decimetric difference with your code, my
code and with Blue Marble.
I have verified that this difference is not due to the difference between
the single rotation matrix normally used (small angle approximation) and the
3 concatenated rotation matrices used by ESRI (for example), which is only
millimetric in quantity. Is it something going on under the hood in Proj.4,
like an embedded call to geodetic <=> geocentric for some reason? I don't
know. We could debate whether 19 centimeters in Y and 12 centimeters in Z
is significant or not. For a straight-up matrix manipulation, I think it
We do agree. Something is somewhere wrong. But what and where?
If I were to get a bug report stating that in a procedure spanning several
operations the result departs from the expected value, I would throw this
Every single step should be tested against independent results from others,
not just the end result of a concatenation of a handful of steps.
Several steps are clear, for instance:
- The right choice of rotation model is made.
- BTW: rotation order has no effect if you work with approximations in the
rotation matrices, like PROJ does.
- The datum transformation functions of PROJ seem to be correct.
- At least the core of the geocentric->geodetic conversion function is OK.
But where does it go wrong?
- Are all parameters passed correctly? Is somewhere a lot of precision lost
because at a higher level not enough decimals are set?
- Are the global structs filled correctly by initialization functions?
- Is the software used as it should be? Take cs2cs: you can do incredible
things with it, but next to nothing is documented (see its man page).
Oscar van Vlijmen
More information about the Proj