Einfache Geofunktionen in Delphi

V0.0 rev004 05.01.2005 (fs)

Letztlich geht es in diesem Artikel darum, GPS-Informationen in vorhandene Karten einzuzeichnen. Beispielsweise, um eine momentane Position in der Karte anzuzeigen oder eine aufgezeichnete GPS-Spur in eine Karte einzuzeichnen.

Dazu sind zwei Schritte nötig: Zunächst muss man wissen, welchen Ausschnitt der Erde die Karte zeigt. Dazu muss man beispielsweise die Weltkoordinaten eines Eckpunktes der Karte sowie ihre Ausdehnung in Nord-Süd- und Ost-West-Richtung kennen. Ich nenne das, die Karte kalibrieren.

(Die Angelegenheit ist kompliziert, wenn man sie exakt erfassen will. Die Erde ist annähernd eine Kugel, und die Gestalt ihrer Oberfläche verzerrt sich, sobald man sie auf eine Ebene projeziert. Glücklicherweise ist die Erde recht groß und unser Horizont recht klein. In einem Stadtplan oder einer Wanderkarte ist der Fehler der Verzerrung geringer als ein Pixel auf dem Bildschirm. Und das ist doch unser Maßstab, oder?)

Im zweiten Schritt kann man dann die Weltkoordinaten, die uns das GPS liefert, in die Karte einzeichnen.

Das Koordinatensystem der Erde

Längengrade

Wenn wir uns die Erdkugel als Männlein vorstellen, bei dem der Nordpol oben und der Südpol unten ist, dann ist ein Längengrad (Meridian) ein Umfang der Erde der Länge nach, der durch Nordpol und Südpol geht. Einmal rundrum sind 360 Grad; ein Grad ist unterteilt in 60 Minuten mit je 60 Sekunden. Der Nullmeridian geht durch Greenwich in England. Häufig, etwas beim GPS, zählt man vom Nullmeridian jeweils 180 Grad nach Ost und nach West. Berlin etwa, liegt auf cirka 13 Grad östlicher Länge.

Breitengrade

Ein Breitengrad ist die "Dicke" (also Breite) des Erdmännchens an einer bestimmten Stelle. Und wie bei einem Männchen, misst man die Dicke entlang einer Ebene senkrecht zur Achse. Die dickste Stelle ist der Äquator. Der Einfachheit halben nehmen wir die Erde als genaue Kugel, also ist der Äquator ebenso lang wie alle Längengrade. Der Äquator liegt auf Breite Null, der Nordpol bei 90 Grad Nord (oder +90 Grad) und der Südpol bei 90 Grad Süd (oder -90 Grad).

Mittlerer Erdradius: 6371,299 km

Daraus errechnet sich der mittlere Erdumfang (2* Pi * Radius): 40032,052 km. 40000 km, das passt, das haben wir schon mal irgendwo gehört.

Auf einer Karte ist Norden oben und Osten rechts. Ein Längengrad verläuft also in Richtung der Y-Achse, ein Breitengrad in Richtung der X-Achse.

Entfernungsberechnung

Zunächst wollen wir die Entfernung entlang der Grade, und dann zwischen zwei beliebeigen Punkten ausrechnen.

Entlang der Längengrade

Weil die Längengrade alle den gleichen Umfang haben, ist entlang derer die Entfernungsmessung besonders einfach. Zunächst wird die Differenz der Breitenwerte bestimmt, und dann mit einem 360stel des Erdumfangs multipliziert. In der einfachsten Form sieht das dann so aus:

function WorldDistanceY(Y1, Y2: real):real;
const WorldUmfang = 40032.052;
begin
 WorldDistanceY:=WorldUmfang/360*(Y2-Y1);
end;

Das funktioniert so ganz prima, nur akzeptiert die Funktion auch Eingabewerte ausserhalb des Wertebereichs klaglos. Delphi stellt uns den wunderbaren Mechanismus der Exceptions zur Verfügung, um diesen Fall abzufangen. Das sieht dann so aus:

function WorldDistanceY(Y1, Y2: real):real;
const WorldUmfang = 40032.052;
begin
 if (Y1 < -90) or (Y1 > 90) or (Y2 < -90) or (Y2 > 90)
 then raise Exception.create('Breitenwerte müssen zwischen -90 und +90 liegen');
 WorldDistanceY:=WorldUmfang/360*(Y2-Y1);
end;

Nun steht die Plausibilitätsprüfung an. Das Projekt GeoDemo1 eignet sich dazu, verschiedene Werte durchzurechnen. Dabei testen wir Werte, deren Ergebis wir kennen. Bei Null Unterschied soll auch die Entfernung 0 sein. Und für 90 Grad soll sie ein Viertel Erdumfang betragen, also ca. 10000 km. Dabei muss es egal sein, wo der Anfangspunkt liegt.

Die letzte Zeile lässt uns kurz die Augenbrauen heben. Aha! Das Vorzeichen gibt an, in welche Richtung die Bewegung erfolgt. Minus heißt nach Süden, Plus nach Norden.

Entlang der Breitengrade

Entfernungen entlang der Breitengrade zu messen, ist etwas schwieriger. Zunächst einmal braucht man eine zusätzliche Information, nämlich um welchen Breitengrad es sich handelt, da sie ja alle verschieden lang sind. Hat man daraus die Länge des gewünschten Breitengrades ermittelt, so kann man die Entfernung wie gehabt errechnen.

Zunächst gilt es also, den Radius des gesuchten Breitengrades zu ermitteln. Die folgende Zeichnung verdeutlicht die Verhältnisse.

Die rote Linie r_E geht vom Erdmittelpunkt bis zu irgendeinem Punkt auf dem Breitengrad. Sie ist der mittlere Erdradius. Der Winkel Beta ist die Breite des Breitengrades. Die grüne Linie r_B ist der Radius des Breitengrades. Aus dem Geometrieunterricht wissen wir, dass sich der Winkel Beta auch zwischen r_B und r_E findet. Wieder der Geometrieunterricht:

cos(Beta) = rB/rE ( = Ankathete / Hypotenuse) oder

rB = cos(Beta)*rE

Bevor wir aber den Algoritmus in Delphi verwirklichen können, müssen wir uns vor Augen führen, dass die Winkelfunktionen in Delphi Winkel nicht in Grad, sondern im Bogenmaß Radiant (rad) erwarten. Die Umrechnung ist einfach:

Grad = rad*180/Pi   oder:

rad = Grad*Pi/180

Ins Unreine geschrieben wird daraus:

LaengeInRad:=LaengeInGrad*Pi/180;
BreitengradRadius := cos(LaengeInRad*PI/180)*WorldRadius;
WorldDistanceX:=Pi*BreitengradRadius/180*(X2-X1);

Zwischenvariable, die nur einmal genutzt werden, sollte man vermeiden (In den dokumentierenden Kommentaren hingegen sind sie gut aufgehoben). Also lässt sich einiges kürzen, und die endgültige Implementierung sieht so aus:

function WorldDistanceX(Y1, X1, X2: real):real;
const WorldRadius = 6371.299;
      Pi = 3.1415926535;
begin
 if (X1 < -180) or (X1 > 180) or (X2 < -180) or (X2 > 180)
 then raise Exception.create('Längenwerte müssen zwischen -180 und +180 liegen');
 if (Y1 < -90) or (Y1 > 90)
 then raise Exception.create('Breitenwerte müssen zwischen -90 und +90 liegen');
 WorldDistanceX:=Pi*cos(Y1*PI/180)*WorldRadius/180*(X2-X1);
end;

Auch hier wieder die Plausibilitätsprüfung. Dafür habe ich das Projekt GeoDemo2 erstellt. Bei Null Längen-Unterschied soll auch die Entfernung 0 sein, und zwar unabhängig von der Breite. Und für 90 Grad Längenunterschied soll sie am Äquator (0 Grad Breite) ein Viertel Erdumfang betragen, am Pol (90 Grad Breite) jedoch 0.

Auch hier gibt wieder das Vorzeichen die Bewegungsrichtung an: Minus bedeutet: Go West.

Fazit

Mit diesen Funktionen sollte man nun eigentlich eine GPS-Spur visuell darstellen können. Das soll demnächst Inhalt eines neuen Artikels sein.