Kartierung
php Script Funktionen zur Koordinatentransformation
Für die hier gezeigten 8 Funktionen zur Koordinatentransformation wird die Existenz der folgenden Variablen vorausgesetzt
global $rw, $hw, $lp, $bp, $lw, $bw, $zone, $ew, $nw, $raster, $ew2, $nw2;
Diese haben folgende Bedeutung
| $rw, $hw | Rechts-, Hochwert im deutschen Gauss-Krüger-System |
| $lp, $bp | Geographische Länge und Breite auf dem Bessel-Ellipsoid (Potsdam Datum) in Grad mit Dezimalpunkt (nicht Komma) |
| $lw, $bw | Geographische Länge und Breite im auf dem Internationalen Ellipsoid (WGS84 Datum), westliche Länge und südliche Breite negativ |
| $zone, $ew, $nw | UTM Gitterzone (2 Ziffern für die Längenzone und ein Buchstabe für das Breitenband, Beispiel: 32U), 7 stelliger Ost- und Nordwert im UTM System |
| $raster, $ew2, $nw2 | UTMREF Raster (UTM Gitterzone und Kennung aus 2 Buchstaben für das 100 km × 100 km Feld, Beispiel: 32UMT), 5 stelliger Ost- und Nordwert |
Umrechnungen im WGS84 System
Die beiden folgenden Funktionen rechnen im WGS84 Kartenbezugssystem (world geodetic system 84) geographische Koordinaten in UTM Koordinaten um und umgekehrt.
Mit der Funktion geo2utm.php lassen sich geographische Längen und Breiten in UTM Zone, Ost- und Nordwert umrechnen:
function geo2utm()
{
/* Copyright (c) 2006, HELMUT H. HEIMEIER
Permission is hereby granted, free of charge, to any person obtaining a
copy of this software and associated documentation files (the "Software"),
to deal in the Software without restriction, including without limitation
the rights to use, copy, modify, merge, publish, distribute, sublicense,
and/or sell copies of the Software, and to permit persons to whom the
Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included
in all copies or substantial portions of the Software.*/
/* Die Funktion wandelt geographische Koordinaten in UTM Koordinaten
um. Geographische Länge lw und Breite bw müssen im WGS84 Datum
gegeben sein. Berechnet werden UTM Zone, Ostwert ew und Nordwert nw.*/
global $rw, $hw, $lp, $bp, $lw, $bw, $zone, $ew, $nw, $raster, $ew2, $nw2;
//Geographische Länge lw und Breite bw im WGS84 Datum
if ($lw =="" || $bw == "")
{
$zone = "";
$ew = "";
$nw = "";
return;
}
if ($lw <= -180 || $lw > 180 || $bw <= -80 || $bw >= 84)
{
echo("Werte nicht im Bereich des UTM Systems\n"+
"-180° <= LW < +180°, -80° < BW < 84° N");
$zone = "";
$ew = "";
$nw = "";
return;
}
$lw = doubleval($lw);
$bw = doubleval($bw);
//WGS84 Datum
//Große Halbachse a und Abplattung f
$a = 6378137.000;
$f = 3.35281068e-3;
$b_sel = 'CDEFGHJKLMNPQRSTUVWXX';
//Polkrümmungshalbmesser c
$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 Meridianbogenlänge
$e0 = $c*(pi()/180)*(1 - 3*$ex2/4 + 45*$ex4/64 - 175*$ex6/256 + 11025*$ex8/16384);
$e2 = $c*(- 3*$ex2/8 + 15*$ex4/32 - 525*$ex6/1024 + 2205*$ex8/4096);
$e4 = $c*(15*$ex4/256 - 105*$ex6/1024 + 2205*$ex8/16384);
$e6 = $c*(- 35*$ex6/3072 + 315*$ex8/12288);
//Längenzone lz und Breitenzone (Band) bz
$lzn = intval(($lw+180)/6) + 1;
$lz = "$lzn";
if ($lzn < 10) $lz = "0" + $lzn;
$bd = intval(1 + ($bw + 80)/8);
$bz = substr($b_sel,$bd-1,1);
//Geographische Breite in Radianten br
$br = $bw * pi()/180;
$tan1 = tan($br);
$tan2 = $tan1*$tan1;
$tan4 = $tan2*$tan2;
$cos1 = cos($br);
$cos2 = $cos1*$cos1;
$cos4 = $cos2*$cos2;
$cos3 = $cos2*$cos1;
$cos5 = $cos4*$cos1;
$etasq = $ex2*$cos2;
//Querkrümmungshalbmesser nd
$nd = $c/sqrt(1 + $etasq);
//Meridianbogenlänge g aus gegebener geographischer Breite bw
$g = ($e0*$bw) + ($e2*sin(2*$br)) + ($e4*sin(4*$br)) + ($e6*sin(6*$br));
//Längendifferenz dl zum Bezugsmeridian lh
$lh = ($lzn - 30)*6 - 3;
$dl = ($lw - $lh)*pi()/180;
$dl2 = $dl*$dl;
$dl4 = $dl2*$dl2;
$dl3 = $dl2*$dl;
$dl5 = $dl4*$dl;
//Maßstabsfaktor auf dem Bezugsmeridian bei UTM Koordinaten m = 0.9996
//Nordwert nw und Ostwert ew als Funktion von geographischer Breite und Länge
if ( $bw < 0 ) {
$nw = 10e6 + 0.9996*($g + $nd*$cos2*$tan1*$dl2/2 +
$nd*$cos4*$tan1*(5-$tan2+9*$etasq)*$dl4/24);
}
else {
$nw = 0.9996*($g + $nd*$cos2*$tan1*$dl2/2 +
$nd*$cos4*$tan1*(5-$tan2+9*$etasq)*$dl4/24);
}
$ew = 0.9996*( $nd*$cos1*$dl + $nd*$cos3*(1-$tan2+$etasq)*$dl3/6 +
$nd*$cos5*(5-18*$tan2+$tan4)*$dl5/120) + 500000;
$zone = $lz . $bz;
$nk = $nw - intval($nw);
if ($nk < 0.5) $nw = "" . intval($nw);
else $nw = "" . (intval($nw) + 1);
while (strlen($nw) < 7)
{
$nw = "0" . $nw;
}
$nk = $ew - intval($ew);
if ($nk < 0.5) $ew = "0" . intval($ew);
else $ew = "0" . intval($ew+1);
return;
}
Mit der Funktion utm2geo.php lassen sich UTM Zone, Ost- und Nordwert in geographische Längen und Breiten umrechnen:
function utm2geo()
{
/* Copyright (c) 2006, HELMUT H. HEIMEIER
Permission is hereby granted, free of charge, to any person obtaining a
copy of this software and associated documentation files (the "Software"),
to deal in the Software without restriction, including without limitation
the rights to use, copy, modify, merge, publish, distribute, sublicense,
and/or sell copies of the Software, and to permit persons to whom the
Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included
in all copies or substantial portions of the Software.*/
/* Die Funktion wandelt UTM Koordinaten in geographische Koordinaten
um. UTM Zone, Ostwert ew und Nordwert nw müssen gegeben sein.
Berechnet werden geographische Länge lw und Breite bw im WGS84 Datum.*/
global $rw, $hw, $lp, $bp, $lw, $bw, $zone, $ew, $nw, $raster, $ew2, $nw2;
// Längenzone zone, Ostwert ew und Nordwert nw im WGS84 Datum
if ($zone == "" || $ew == "" || $nw == "")
{
$zone = "";
$ew = "";
$nw = "";
return;
}
$band = substr($zone,2,1);
$z1 = intval($zone);
$ew1 = intval($ew);
$nw1 = intval($nw);
//WGS 84 Datum
//Große Halbachse a und Abplattung f
$a = 6378137.000;
$f = 3.35281068e-3;
//Polkrümmungshalbmesser c
$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);
//Entscheidung Nord-/Süd Halbkugel
if ($band >= "N"|| $band == "")
$m_nw = $nw1;
else
$m_nw = $nw1 - 10e6;
//Geographische Breite bf zur Meridianbogenlänge gf = m_nw
$sigma = ($m_nw/0.9996)/$e0;
$sigmr = $sigma*pi()/180;
$bf = $sigma + $f2*sin(2*$sigmr) + $f4*sin(4*$sigmr) + $f6*sin(6*$sigmr);
//Breite bf in Radianten
$br = $bf * pi()/180;
$tan1 = tan($br);
$tan2 = $tan1*$tan1;
$tan4 = $tan2*$tan2;
$cos1 = cos($br);
$cos2 = $cos1*$cos1;
$etasq = $ex2*$cos2;
//Querkrümmungshalbmesser nd
$nd = $c/sqrt(1 + $etasq);
$nd2 = $nd*$nd;
$nd4 = $nd2*$nd2;
$nd6 = $nd4*$nd2;
$nd3 = $nd2*$nd;
$nd5 = $nd4*$nd;
//Längendifferenz dl zum Bezugsmeridian lh
$lh = ($z1 - 30)*6 - 3;
$dy = ($ew1-500000)/0.9996;
$dy2 = $dy*$dy;
$dy4 = $dy2*$dy2;
$dy3 = $dy2*$dy;
$dy5 = $dy2*$dy3;
$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);
//Geographische Breite bw und Länge lw als Funktion von Ostwert ew
// und Nordwert nw
$bw = $bf + (180/pi()) * ($b2*$dy2 + $b4*$dy4 + $b6*$dy6);
$lw = $lh + (180/pi()) * ($l1*$dy + $l3*$dy3 + $l5*$dy5);
return;
}
Die beiden folgenden Funktionen rechnen im WGS84 Kartenbezugssystem (world geodetic system 84) zivile UTM Koordinaten in MGRS (military grid reference system), auch UTMREF Koordinaten genannt, um und umgekehrt.
Die Funktion utm2mgr.php ermittelt aus UTM Gitterzone und 7-stelligem Ost- und Nordwert das UTMREF Raster mit 5-stelligem Ost- und Nordwert:
function utm2mgr()
{
/* Copyright (c) 2006, HELMUT H. HEIMEIER
Permission is hereby granted, free of charge, to any person obtaining a
copy of this software and associated documentation files (the "Software"),
to deal in the Software without restriction, including without limitation
the rights to use, copy, modify, merge, publish, distribute, sublicense,
and/or sell copies of the Software, and to permit persons to whom the
Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included
in all copies or substantial portions of the Software.*/
/* Die Funktion wandelt zivile UTM Koordinaten in militärische Koordinaten
um. UTM Zone zone, Ostwert ew und Nordwert nw müssen gegeben sein.
Zurückgegeben wird das Rasterfeld raster sowie die aus den
letzten 5 Stellen von Ost- und Nordwert gebildete Koordinatenangabe
UTMREF.*/
global $rw, $hw, $lp, $bp, $lw, $bw, $zone, $ew, $nw, $raster, $ew2, $nw2;
// Längenzone zone, Ostwert ew und Nordwert nw im WGS84 Datum
if ($zone == "" || $ew == "" || $nw == "")
{
$zone = "";
$ew = "";
$nw = "";
return;
}
$z1 = substr($zone,0,2);
$z2 = substr($zone,2,1);
$ew1 = intval(substr($ew,0,2));
$nw1 = intval(substr($nw,0,2));
$ew2 = substr($ew,2,5);
$nw2 = substr($nw,2,5);
$m_east = 'ABCDEFGHJKLMNPQRSTUVWXYZ';
$m_north = 'ABCDEFGHJKLMNPQRSTUV';
if ($z1 < "01" || $z1 > "60" || $z2 < "C" ||$z2 > "X")
echo "$zone ist keine gültige UTM Zonenangabe";
$i = $z1 % 3;
if ($i == 1) $m_ce = $ew1 - 1;
if ($i == 2) $m_ce = $ew1 + 7;
if ($i == 0) $m_ce = $ew1 + 15;
$i = $z1 % 2;
if ($i == 1) $m_cn = 0;
else $m_cn = 5;
$i = $nw1;
while ($i-20 >= 0) $i = $i-20;
$m_cn = $m_cn + $i;
if ($m_cn > 19) $m_cn = $m_cn - 20;
$raster = $zone . substr($m_east,$m_ce,1) . substr($m_north,$m_cn,1);
return;
}
Die Funktion mgr2utm.php ermittelt aus UTMREF Raster und 5-stelligem Ost- und Nordwert die UTM Gitterzone und 7-stelligen Ost- und Nordwert:
function mgr2utm()
{
/* Copyright (c) 2006, HELMUT H. HEIMEIER
Permission is hereby granted, free of charge, to any person obtaining a
copy of this software and associated documentation files (the "Software"),
to deal in the Software without restriction, including without limitation
the rights to use, copy, modify, merge, publish, distribute, sublicense,
and/or sell copies of the Software, and to permit persons to whom the
Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included
in all copies or substantial portions of the Software.*/
/* Die Funktion wandelt militärische UTM Koordinaten (MGR oder
UTMREF) in zivile UTM Koordinaten um.
UTM Zone zone, raster und utmref müssen gegeben sein.
In zone muss die aus 2 Ziffern bestehende Längenzone enthaltens ein
gefolgt von der aus einem Buchstaben bestehenden Bandangabe.
In raster muss die aus 2 Buchstaben bestehende Kennung für das
100 km x 100 km Rasterfeld enthalten sein.
In UTMREF muss der 5 stellige Ostwert stehen gefolgt von einem blank
und dem 5 stelligen Nordwert.
Berechnet wird daraus der 7 stellige Ost- und Nordwert im zivilen
UTM System.*/
global $rw, $hw, $lp, $bp, $lw, $bw, $zone, $ew, $nw, $raster, $ew2, $nw2;
// Längenzone zone, Ostwert ew und Nordwert nw im WGS84 Datum
if ($raster == "" || $ew2 == "" || $nw2 == "")
{
$raster = "";
$ew2 = "";
$nw2 = "";
return;
}
$m_east_0 = "STUVWXYZ";
$m_east_1 = "ABCDEFGH";
$m_east_2 = "JKLMNPQR";
$m_north_0 = "FGHJKLMNPQRSTUVABCDE";
$m_north_1 = "ABCDEFGHJKLMNPQRSTUV";
$zone = substr($raster,0,3);
$r_east = substr($raster,3,1);
$r_north = substr($raster,4,1);
$i = intval(substr($raster,0,2)) % 3;
if ($i == 0) $m_ce = strpos($m_east_0,$r_east)+1;
if ($i == 1) $m_ce = strpos($m_east_1,$r_east)+1;
if ($i == 2) $m_ce = strpos($m_east_2,$r_east)+1;
$m_ce = "0" . $m_ce;
$ew = $m_ce . $ew2;
$i = intval(substr($raster,0,2)) % 2;
if ($i == 0) $m_cn = strpos($m_north_0,$r_north);
else $m_cn = strpos($m_north_1,$r_north);
$band = substr($raster,2,1);
if ($band >= "N"){
if ($band == "Q" && $m_cn < 10)
$m_cn = $m_cn + 20;
if ($band >= "R")
$m_cn = $m_cn + 20;
if ($band == "S" && $m_cn < 30)
$m_cn = $m_cn + 20;
if ($band >= "T")
$m_cn = $m_cn + 20;
if ($band == "U" && $m_cn < 50)
$m_cn = $m_cn + 20;
}
else {
if ($band == "C" && $m_cn < 10)
$m_cn = $m_cn + 20;
if ($band >= "D")
$m_cn = $m_cn + 20;
if ($band == "F" && $m_cn < 30)
$m_cn = $m_cn + 20;
if ($band >= "G")
$m_cn = $m_cn + 20;
if ($band == "H" && $m_cn < 50)
$m_cn = $m_cn + 20;
if ($band >= "J")
$m_cn = $m_cn + 20;
if ($band == "K" && $m_cn < 70)
$m_cn = $m_cn + 20;
if ($band >= "L")
$m_cn = $m_cn + 20;
}
if (strlen($m_cn) == 1)
$nw = "0" . $m_cn . $nw2;
else
$nw = "" .$m_cn . $nw2;
return;
}
Umrechnungen für das Bessel-Ellipsoid
Die beiden folgenden Funktionen rechnen für das Potsdam Kartendatum (Bessel Ellipsoid) geographische Koordinaten in Gauss-Krüger Koordinaten um und umgekehrt.
Mit der Funktion geo2gk.php lassen sich geographische Längen und Breiten für das Potsdam Kartendatum (Bessel-Ellipsoid) in deutsche Gauss-Krüger Koordinaten umrechnen:
function geo2gk()
{
/* Copyright (c) 2006, HELMUT H. HEIMEIER
Permission is hereby granted, free of charge, to any person obtaining a
copy of this software and associated documentation files (the "Software"),
to deal in the Software without restriction, including without limitation
the rights to use, copy, modify, merge, publish, distribute, sublicense,
and/or sell copies of the Software, and to permit persons to whom the
Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included
in all copies or substantial portions of the Software.*/
/* Die Funktion wandelt geographische Koordinaten in GK Koordinaten
um. Geographische Länge lp und Breite bp müssen im Potsdam Datum
gegeben sein. Berechnet werden Rechtswert rw und Hochwert hw.*/
global $rw, $hw, $lp, $bp, $lw, $bw, $zone, $ew, $nw, $raster, $ew2, $nw2;
//Geographische Länge lp und Breite bp im Potsdam Datum
if ($lp == "" || $bp == "")
{
$rw = "";
$hw = "";
return;
}
$lp = doubleval($lp);
$bp = doubleval($bp);
// Grenzen des Gauss-Krüger-Systems für Deutschland 46° N < bp < 55° N,
// 5° E < lp < 16° E
if ($bp < 46 || $bp > 56 || $lp < 5 || $lp > 16)
{
// echo("Werte außerhalb des für Deutschland definierten Gauss-Krüger-Systems\n"+
// " 5° E < LP < 16° E, 46° N < BP < 55° N");
$rw = " ";
$hw = " ";
return;
}
// Potsdam Datum
// Große Halbachse a und Abplattung f
$a = 6377397.155;
$f = 3.34277321e-3;
// Polkrümmungshalbmesser c
$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 Meridianbogenlänge
$e0 = $c*(pi()/180)*(1 - 3*$ex2/4 + 45*$ex4/64 - 175*$ex6/256 + 11025*$ex8/16384);
$e2 = $c*( - 3*$ex2/8 + 15*$ex4/32 - 525*$ex6/1024 + 2205*$ex8/4096);
$e4 = $c*(15*$ex4/256 - 105*$ex6/1024 + 2205*$ex8/16384);
$e6 = $c*( - 35*$ex6/3072 + 315*$ex8/12288);
// Breite in Radianten
$br = $bp * pi()/180;
$tan1 = tan($br);
$tan2 = $tan1*$tan1;
$tan4 = $tan2*$tan2;
$cos1 = cos($br);
$cos2 = $cos1*$cos1;
$cos4 = $cos2*$cos2;
$cos3 = $cos2*$cos1;
$cos5 = $cos4*$cos1;
$etasq = $ex2*$cos2;
// Querkrümmungshalbmesser nd
$nd = $c/sqrt(1 + $etasq);
// Meridianbogenlänge g aus gegebener geographischer Breite bp
$g = $e0*$bp + $e2*sin(2*$br) + $e4*sin(4*$br) + $e6*sin(6*$br);
// Längendifferenz dl zum Bezugsmeridian lh
$kz = intval(($lp+1.5)/3);
$lh = $kz*3;
$dl = ($lp - $lh)*pi()/180;
$dl2 = $dl*$dl;
$dl4 = $dl2*$dl2;
$dl3 = $dl2*$dl;
$dl5 = $dl4*$dl;
// Hochwert hw und Rechtswert rw als Funktion von geographischer Breite und Länge
$hw = ($g + $nd*$cos2*$tan1*$dl2/2 + $nd*$cos4*$tan1*(5-$tan2+9*$etasq)
*$dl4/24);
$rw = ($nd*$cos1*$dl + $nd*$cos3*(1-$tan2+$etasq)*$dl3/6 +
$nd*$cos5*(5-18*$tan2+$tan4)*$dl5/120 + $kz*1e6 + 500000);
$nk = $hw - intval($hw);
if ($nk < 0.5) $hw = intval($hw);
else $hw = intval($hw) + 1;
$nk = $rw - intval($rw);
if ($nk < 0.5) $rw = intval($rw);
else $rw = intval($rw+1);
return;
}
Mit der Funktion gk2geo.php lassen sich Gauss-Krüger Koordinaten in geographische Längen und Breiten des Potsdam Datums umrechnen:
function gk2geo()
{
/* Copyright (c) 2007, HELMUT H. HEIMEIER
Permission is hereby granted, free of charge, to any person obtaining a
copy of this software and associated documentation files (the "Software"),
to deal in the Software without restriction, including without limitation
the rights to use, copy, modify, merge, publish, distribute, sublicense,
and/or sell copies of the Software, and to permit persons to whom the
Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included
in all copies or substantial portions of the Software.*/
/* Die Funktion wandelt GK Koordinaten in geographische Koordinaten
um. Rechtswert rw und Hochwert hw müssen gegeben sein.
Berechnet werden geographische Länge lp und Breite bp
im Potsdam Datum.*/
global $rw, $hw, $lp, $bp, $lw, $bw, $zone, $ew, $nw, $raster, $ew2, $nw2;
// Rechtswert rw und Hochwert hw im Potsdam Datum
if ($rw == "" || $hw == "")
{
$lp = "";
$bp = "";
return;
}
$rw = intval($rw);
$hw = intval($hw);
// Potsdam Datum
// Große Halbachse a und Abplattung f
$a = 6377397.155;
$f = 3.34277321e-3;
// Polkrümmungshalbmesser c
$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 = $hw/$e0;
$sigmr = $sigma*pi()/180;
$bf = $sigma + $f2*sin(2*$sigmr) + $f4*sin(4*$sigmr) + $f6*sin(6*$sigmr);
// Breite bf in Radianten
$br = $bf * pi()/180;
$tan1 = tan($br);
$tan2 = $tan1*$tan1;
$tan4 = $tan2*$tan2;
$cos1 = cos($br);
$cos2 = $cos1*$cos1;
$etasq = $ex2*$cos2;
// Querkrümmungshalbmesser nd
$nd = $c/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 = intval($rw/1e6);
$lh = $kz*3;
$dy = $rw-($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
$bp = $bf + (180/pi()) * ($b2*$dy2 + $b4*$dy4 + $b6*$dy6);
$lp = $lh + (180/pi()) * ($l1*$dy + $l3*$dy3 + $l5*$dy5);
if ($lp < 5 || $lp > 16 || $bp < 46 || $bp > 56)
{
echo("RW und/oder HW ungültig für das deutsche Gauss-Krüger-System
\n");
$lp = "";
$bp = "";
}
return;
}
Funktionen zur Änderung des Kartenbezugssystems WGS84 - Potsdam
Mit den beiden folgenden Funktionen lässt sich das Kartenbezugsdatum vom WGS84 System auf das Potsdam Datum (Zentralpunkt Rauenberg) und umgekehrt verschieben. Bei der Transformation werden die Ellipsoidachsen parallel verschoben um dx = 587 m, dy = 16 m und dz = 393 m. Hierbei liegt der Koordinatenursprung im Erdmittelpunkt, die positive x-Achse schneidet Äquator und Nullmeridian, die positive y-Achse schneidet Äquator und 90 grad E Meridian und die positive z-Achse schneidet den Nordpol. Die shift-Parameter dx, dy und dz sind so gewählt, dass sich bei der Umformung übereinstimmung mit den von Garmin GPS-Geräten ermittelten Werten ergibt.
Man beachte, dass in der deutschen Landesvermessung im allgemeinen mit anderen Parametern gerechnet wird, nämlich mit dx = 631 m, dy = 23 m und dz = 451 m.
Die Funktion wgs2pot.php formt geographischen Längen und Breiten des WGS84 Datums in solche des Potsdam Datums um:
function wgs2pot()
{
/* Copyright (c) 2006, HELMUT H. HEIMEIER
Permission is hereby granted, free of charge, to any person obtaining a
copy of this software and associated documentation files (the "Software"),
to deal in the Software without restriction, including without limitation
the rights to use, copy, modify, merge, publish, distribute, sublicense,
and/or sell copies of the Software, and to permit persons to whom the
Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included
in all copies or substantial portions of the Software.*/
/* Die Funktion verschiebt das Kartenbezugssystem (map datum) vom
WGS84 Datum (World Geodetic System 84) zum in Deutschland
gebräuchlichen Potsdam-Datum. Geographische Länge lw und Breite
bw gemessen in grad auf dem WGS84 Ellipsoid müssen
gegeben sein. Ausgegeben werden geographische Länge lp
und Breite bp (in grad) auf dem Bessel-Ellipsoid.
Bei der Transformation werden die Ellipsoidachsen parallel
verschoben um dx = -587 m, dy = -16 m und dz = -393 m.*/
global $rw, $hw, $lp, $bp, $lw, $bw, $zone, $ew, $nw, $raster, $ew2, $nw2;
// Geographische Länge lw und Breite bw im WGS84 Datum
if ($lw == "" || $bw == "")
{
$lp = "";
$bp = "";
return;
}
$lw = doubleval($lw);
$bw = doubleval($bw);
// Quellsystem WGS 84 Datum
// Große Halbachse a und Abplattung fq
$a = 6378137.000;
$fq = 3.35281066e-3;
// Zielsystem Potsdam Datum
// Abplattung f
$f = $fq - 1.003748e-5;
// Parameter für datum shift
$dx = -587;
$dy = -16;
$dz = -393;
// Quadrat der ersten numerischen Exzentrizität in Quell- und Zielsystem
$e2q = (2*$fq-$fq*$fq);
$e2 = (2*$f-$f*$f);
// Breite und Länge in Radianten
$b1 = $bw * (pi()/180);
$l1 = $lw * (pi()/180);
// Querkrümmungshalbmesser nd
$nd = $a/sqrt(1 - $e2q*sin($b1)*sin($b1));
// Kartesische Koordinaten des Quellsystems WGS84
$xw = $nd*cos($b1)*cos($l1);
$yw = $nd*cos($b1)*sin($l1);
$zw = (1 - $e2q)*$nd*sin($b1);
// Kartesische Koordinaten des Zielsystems (datum shift) Potsdam
$x = $xw + $dx;
$y = $yw + $dy;
$z = $zw + $dz;
// Berechnung von Breite und Länge im Zielsystem
$rb = sqrt($x*$x + $y*$y);
$b2 = (180/pi()) * atan(($z/$rb)/(1-$e2));
if ($x > 0)
$l2 = (180/pi()) * atan($y/$x);
if ($x < 0 && $y > 0)
$l2 = (180/pi()) * atan($y/$x) + 180;
if ($x < 0 && $y < 0)
$l2 = (180/pi()) * atan($y/$x) - 180;
$lp = $l2;
$bp = $b2;
/* if ($lp < 5 || $lp > 16 || $bp < 46 || $bp > 56)
{
$lp = "";
$bp = "";
}*/
return;
}
Die Funktion pot2wgs.php formt geographische Längen und Breiten des Potsdam Datums in solche des WGS84 Datums um:
function pot2wgs()
{
/* Copyright (c) 2006, HELMUT H. HEIMEIER
Permission is hereby granted, free of charge, to any person obtaining a
copy of this software and associated documentation files (the "Software"),
to deal in the Software without restriction, including without limitation
the rights to use, copy, modify, merge, publish, distribute, sublicense,
and/or sell copies of the Software, and to permit persons to whom the
Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included
in all copies or substantial portions of the Software.*/
/* Die Funktion verschiebt das Kartenbezugssystem (map datum) vom in
Deutschland gebräuchlichen Potsdam-Datum zum WGS84 (World Geodetic
System 84) Datum. Geographische Länge lp und Breite bp gemessen in
grad auf dem Bessel-Ellipsoid müssen gegeben sein.
Ausgegeben werden geographische Länge lw und
Breite bw (in grad) auf dem WGS84-Ellipsoid.
Bei der Transformation werden die Ellipsoidachsen parallel
verschoben um dx = 587 m, dy = 16 m und dz = 393 m.*/
global $rw, $hw, $lp, $bp, $lw, $bw, $zone, $ew, $nw, $raster, $ew2, $nw2;
// Geographische Länge lp und Breite bp im Potsdam Datum
if ($lp == "" || $bp == "")
{
$lw = "";
$bw = "";
return;
}
$lp = doubleval($lp);
$bp = doubleval($bp);
// Quellsystem Potsdam Datum
// Große Halbachse a und Abplattung fq
$a = 6378137.000 - 739.845;
$fq = 3.35281066e-3 - 1.003748e-05;
// Zielsystem WGS84 Datum
// Abplattung f
$f = 3.35281066e-3;
// Parameter für datum shift
$dx = 587;
$dy = 16;
$dz = 393;
// Quadrat der ersten numerischen Exzentrizität in Quell- und Zielsystem
$e2q = (2*$fq-$fq*$fq);
$e2 = (2*$f-$f*$f);
// Breite und Länge in Radianten
$b1 = $bp * (pi()/180);
$l1 = $lp * (pi()/180);
// Querkrümmungshalbmesser nd
$nd = $a/sqrt(1 - $e2q*sin($b1)*sin($b1));
// Kartesische Koordinaten des Quellsystems Potsdam
$xp = $nd*cos($b1)*cos($l1);
$yp = $nd*cos($b1)*sin($l1);
$zp = (1 - $e2q)*$nd*sin($b1);
// Kartesische Koordinaten des Zielsystems (datum shift) WGS84
$x = $xp + $dx;
$y = $yp + $dy;
$z = $zp + $dz;
// Berechnung von Breite und Länge im Zielsystem
$rb = sqrt($x*$x + $y*$y);
$b2 = (180/pi()) * atan(($z/$rb)/(1-$e2));
if ($x > 0)
$l2 = (180/pi()) * atan($y/$x);
if ($x < 0 && $y > 0)
$l2 = (180/pi()) * atan($y/$x) + 180;
if ($x < 0 && $y < 0)
$l2 = (180/pi()) * atan($y/$x) - 180;
$lw = $l2;
$bw = $b2;
return;
}