Problem bei Umrechnung von Pulkovo 1942 in WGS84 (Geodäsie/Vermessung)
Hallo,
ich habe ein Problem bei einer Koordinatentransformation.
Ich habe Koordinaten im GaussKrüger Format (Rechtswert/Hochwert) im Quellkoordinatensystem Pulkovo 1942 (Krassovsky). Ich würde diese gerne ins WGS84 transformieren.
Zum Vorgehen:
1.) Umwandeln von GaussKrüger Koordinaten in Geo Koordinaten
2.) Umwandeln von Geographischen Koordinaten in Geozentrische Koordinaten
3.) Speziallfall Molodensky-Badekas Tranformation
4.) Umwandeln von Geozentrischen Koordinaten in Geographische Koordinaten
Ich nehme als Referenz folgende Online-Transformation:
https://epsg.io/transform#s_srs=2399-1675&t_srs=4326&x=5389306.9210650&y=5823166.2710200
Meine Daten:
Rechtswert: 5389306.921065
Hochwert: 5823166.27102
Was ich erhalten möchte:
Längengrad: 13.3672138°
Breitengrad: 52.5249447°
1.)
a = 6378245;
f = 1/293.3;
c = a / (1 - f);
// Quadrat der zweiten numerischen Exzentrizität
ex2 = (2 * f - f * f) / ((1 - f) * (1 - f));
ex4 = ex2 * ex2;
ex6 = ex4 * ex2;
ex8 = ex4 * ex4;
// Koeffizienten zur Berechnung der geographischen Breite aus gegebener
// Meridianbogenlänge
e0 = c * (pi / 180) * (1 - 3 * ex2 / 4 + 45 * ex4 / 64 - 175 * ex6 / 256 + 11025 * ex8 / 16384);
f2 = (180 / pi) * (3 * ex2 / 8 - 3 * ex4 / 16 + 213 * ex6 / 2048 - 255 * ex8 / 4096);
f4 = (180 / pi) * (21 * ex4 / 256 - 21 * ex6 / 256 + 533 * ex8 / 8192);
f6 = (180 / pi) * (151 * ex6 / 6144 - 453 * ex8 / 12288);
// Geographische Breite bf zur Meridianbogenlänge gf = hw
sigma = (northing) / e0;
sigmr = sigma * pi / 180;
bf = sigma + f2 * Math.sin(2 * sigmr) + f4 * Math.sin(4 * sigmr) + f6 * Math.sin(6 * sigmr);
// Breite bf in Radianten
br = bf * pi / 180;
tan1 = Math.tan(br);
tan2 = tan1 * tan1;
tan4 = tan2 * tan2;
cos1 = Math.cos(br);
cos2 = cos1 * cos1;
etasq = ex2 * cos2;
// Querkrümmungshalbmesser nd
nd = c / Math.sqrt(1 + etasq);
nd2 = nd * nd;
nd4 = nd2 * nd2;
nd6 = nd4 * nd2;
nd3 = nd2 * nd;
nd5 = nd4 * nd;
// Längendifferenz dl zum Bezugsmeridian lh
kz = (int)(easting / 1e6);
lh = kz * 3;
dy = (easting - (kz * 1e6 + 500000));
dy2 = dy * dy;
dy4 = dy2 * dy2;
dy3 = dy2 * dy;
dy5 = dy4 * dy;
dy6 = dy3 * dy3;
b2 = -tan1 * (1 + etasq) / (2 * nd2);
b4 = tan1 * (5 + 3 * tan2 + 6 * etasq * (1 - tan2)) / (24 * nd4);
b6 = -tan1 * (61 + 90 * tan2 + 45 * tan4) / (720 * nd6);
l1 = 1 / (nd * cos1);
l3 = -(1 + 2 * tan2 + etasq) / (6 * nd3 * cos1);
l5 = (5 + 28 * tan2 + 24 * tan4) / (120 * nd5 * cos1);
// Geographischer Breite bp und Länge lp als Funktion von Rechts- und
// Hochwert
latitude = bf + (180 / pi) * (b2 * dy2 + b4 * dy4 + b6 * dy6);
longitude = lh + (180 / pi) * (l1 * dy + l3 * dy3 + l5 * dy5)
Nach 1.) erhalte ich folgende Werte:
Längengrad: 13.369036382735441
Breitengrad: 52.525467788632824
2.) Umrechnung in geozentrische Koordinaten:
pi = Math.PI;
a = 6378245;
f = 1/293.3;
e2 = (2 * f - f * f);
longitudeRAD = longitude * (pi / 180);
latitudeRAD = latitude * (pi / 180);
// Querkrümmungshalbmesser
double nd = a /sqrt(1 - e2 * sin(latitudeRAD)²);
// Umrechnung in geozentrische Koordinaten
xp = (nd) * cos(latitudeRAD) * cos(longitudeRAD);
yp = (nd) * cos(latitudeRAD) * sin(longitudeRAD);
zp = ((1 - e2wgs84) * nd) * sin(latitudeRAD);
Nach 2.) erhalte ich folgende geozentrische Koordinaten:
xp: 3783403.403351258
yp: 899173.4481951335
zp: 5038677.955342111
3.) Transformation von Pulkovo 1942 zu WGS84
Transformationsparameter:
dx = 24.0;
dy = -123.0;
dz = -94.0;
rx = -0.02 * (pi / 180);
ry = 0.25 * (pi / 180);
rz = 0.13 * (pi / 180);
scale = 1.1;
M = 1 + scale*10e-6;
// 7-Parameter Transformation
x = xp + (M * (dx - (rz * dy) + (ry * dz)));
y = yp + (M * ((rz * dx) + dy - (rx * dz)));
z = zp + (M * ((-ry * dx) + (rx * dy) + dz));
4.) Umrechnung geozentrische Koordinaten in geographische Koordinaten
rb = Math.sqrt(Math.pow(x, 2) + Math.pow(y, 2));
latitutdeWGS84 = (180 / pi) * Math.atan((z / rb) / (1 - e2wgs84));
longitudeWGS84 = (180 / pi) * Math.atan(y / x);
Hier meine Frage:
Laut https://epsg.io/1032-method multipliziere ich die Quellkoordinaten mit den Rotationen und M und verschiebe dann um dx,dy,dz.
Mache ich es so, erhalte ich die folgenden Werte:
52.27755746434763°, 13.447067545543165°
Wenn ich wie oben in 3.) beschrieben dx,dy und dz rotiere und mit M multiplizieren und dann auf xp,yp,zp addiere erhalte ich folgende Werte:
52.52501509125434°, 13.367192233116334°
Ersteres macht mehr Sinn und ist ja auch so definiert. Die Werte sehen umgedreht allerdings besser aus.
Die unteren Werte haben eine Distanz von 8m zum erwarteten Punkt. Ich hätte es gerne noch genauer und vermute ich habe irgendwo einen Denkfehler.
Ich wäre über jegliche Anregung dankbar!
Problem bei Umrechnung von Pulkovo 1942 in WGS84
Hi,
Ich wäre über jegliche Anregung dankbar!
Manchmal gibt es Unterschiede in den Konventionen der Drehwinkel. Ich würde im übrigen meinen, dass die angegebenen Werte in Sekunden sind und nicht in Grad! Das gewünschte Ergebnis konnte ich mit GeoTra aber rekonstruieren, vgl. den angefügten Screenshot.
Im übrigen gibt es aktuellere Transformationsparameter von 2001, vgl. Tabelle 3 im Artikel Möglichkeiten zur Datumstransformation zwischen Koordinatenreferenzsystemen
/Micha
--
applied-geodesy.org - OpenSource Least-Squares Adjustment Software for Geodetic Sciences
Gelöst - Problem bei Umrechnung von Pulkovo 1942 in WGS84
Hi,
Manchmal gibt es Unterschiede in den Konventionen der Drehwinkel. Ich würde im übrigen meinen, dass die angegebenen Werte in Sekunden sind und nicht in Grad!
Oh mann, manchmal sieht man den Wald vor lauter Bäumen nicht, danke!
Rechne die Parameter nun von arcsec in rad um und erhalte die gewünschten Werten.
Im übrigen gibt es aktuellere Transformationsparameter von 2001, vgl. Tabelle 3 im Artikel Möglichkeiten zur Datumstransformation zwischen Koordinatenreferenzsystemen
Ja, danke. Ich wollte voerst die Referenz erfüllen und dadurch die dort hinterlegten Parameter.
Danke nochmal.
/Micha