[Proj] Using Proj.4 to compress a 10-parameter Molodensky-Badekas datum shift into 7 parameters
Clifford J Mugnier
cjmce at lsu.edu
Tue Oct 12 09:37:26 EST 2010
The "evaluation point" is the Datum Origin, a rather important item. In regard to the "unexplained" difference in ellipsoid height results, my guess is the approximations made by the inverse algorithm from Geocentric coordinates.
Personally, I prefer the elegance of incorporating the Datum Origin.
I have also been told by Dr. Muneendra Kumar that Badekas had nothing to do with Molodensky's original method. I believe they were in school together at the same time, and Prof. Molodensky preceeded them on another continent by some years if not decades.
Clifford J. Mugnier, C.P., C.M.S.
Chief of Geodesy,
Center for GeoInformatics
Department of Civil Engineering
Patrick F. Taylor Hall 3223A
LOUISIANA STATE UNIVERSITY
Baton Rouge, LA 70803
Voice and Facsimile: (225) 578-8536 [Academic]
Voice and Facsimile: (225) 578-4578 [Research]
Cell: (225) 238-8975 [Academic & Research]
Honorary Life Member of the
Louisiana Society of Professional Surveyors
Fellow Emeritus of the ASPRS
Member of the Americas Petroleum Survey Group
________________________________
From: proj-bounces at lists.maptools.org on behalf of Mikael Rittri
Sent: Tue 12-Oct-10 03:45
To: PROJ.4 and general Projections Discussions
Subject: [Proj] Using Proj.4 to compress a 10-parameter Molodensky-Badekas datum shift into 7 parameters
A Molodensky-Badekas datum shift is represented by 10 parameters.
But it can be compressed to a 7-parameter representation.
As Noel Zinn has explained (http://lists.maptools.org/pipermail/proj/2010-October/005424.html),
the 10-parameter representation has advantages over the 7-parameter
one, for people who derive, analyze and compare datum shifts.
But for people who just use a datum shift to transform coordinates,
the 7-parameter representation gives the same results. In this
limited sense, the compression is lossless, which is useful if
your software (like Proj.4) does not support 10-parameter M-B.
I learned how to do the compression from Melita Kennedy (thanks
again!). And I have now realized that one can persuade Proj.4 to
do most of the work of the compression process. Since I thought
this is rather nice, I'd like to explain it by an example.
Example: M-B datum shift EPSG:1096, "La Canoa to WGS 84 (2)",
which has the same parameters as EPSG:1771, "La Canoa to REGVEN (1)".
Source ellipsoid: International 1924.
Target ellipsoid: WGS 84.
dX = -270.933 m
dY = 115.599 m
dZ = -360.226 m
RX = -5.266 arc sec
RY = -1.238 arc sec
RZ = 2.381 arc sec
dS = -5.109 ppm
eX = 2464351.59 m
eY = -5783466.61 m
eZ = 974809.81 m
where eX,eY,eZ are the geocentric coordinates of the so-called
evaluation point.
In principle, the equivalent 7-parameter datum shift should have
the same rotations RX, RY, RZ and the same scale difference dS.
But there is a catch: the sign convention for the rotations.
M-B uses the same sign convention as Coordinate Frame Rotation,
while cs2cs uses the opposite convention of the Position Vector
Transform, so we must negate the three rotations. So, in cs2cs
syntax, the 7-parameter datum shift we seek has the form
+towgs84=?,?,?,5.266,1.238,-2.381,-5.109
To find the question marks, we have to do three steps.
i) Make a 7-parameter "something" by temporarily replacing
the question marks by zeros. It is not really a datum
shift, but we might call it a 7-parameter compressor.
Use it in a cs2cs command like this
>cs2cs -f "%.3f" +proj=geocent +towgs84=0,0,0,5.226,1.238,-2.381,-5.109 +to +datum=WGS84 +proj=geocent
As input, use the evaluation point eX,eY,eZ of the M-B datum shift:
2464351.59 -5783466.61 974809.81
We get the output
2464278.090 -5783490.207 974643.507
ii) Now subtract the output row from the input row, and we get
73.500 23.597 166.303
iii) Then add the original dX,dY,dZ from the M-B datum shift,
which were -270.933,115.599,-360.226, and we get
-197.433 139.196 -193.923
The last three values are the question marks in our 7-parameter
datum shift.
To test it, we use the test point given in section 2.4.4.1,
pages 114-115, EPSG Guidance Note 7.2, http://www.epsg.org/guides/docs/G7-2.pdf
>cs2cs +ellps=intl +proj=longlat +towgs84=-197.433,139.196,-193.923,5.226,1.238,-2.381,-5.109 +to +datum=WGS84 +proj=longlat
66d4'48.091"W 9d35'0.386"N 201.46
66d4'54.705"W 9d34'49.001"N 180.499
According to GN 7.2, the result should have been
66°04'54.705"W 9°34'49.001"N 180.51 m
so we have achieved a lossless compression - apart from the
unexplained difference of about 1.1 cm in the ellipoidal height.
Regards,
Mikael Rittri
Carmenta AB
Sweden
www.carmenta.com
_______________________________________________
Proj mailing list
Proj at lists.maptools.org
http://lists.maptools.org/mailman/listinfo/proj
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://lists.maptools.org/pipermail/proj/attachments/20101012/9d6198fe/attachment.htm
More information about the Proj
mailing list