Vimage Special
RTA: Einen Kartennetzentwurf programmieren
Dieses Spezialkapitel erläutert, wie man einen Kartennetzentwurf in RTA programmiert.
Vimage rechnet Kartennetzentwürfe in einem zweistufigen Verfahren:
• Vimage enthält Projection engines, die eine Koordinatentransformation lediglich abstrakt mit einer allgemeinen Funktion (x’,y’)=ƒ(x, y) ausführen.
• Die konkrete Funktion wird nun mit einem RT-Assemblerprogramm beschrieben. Das RT-Assemblerprogramm ist eine externe Datei. Dieses Programm wird innerhalb der Projection engine auf einem virtuellen Prozessor ausgeführt. Dabei wird die tatsächliche Formel, z. B. x’ = sin(x), y’ = y gerechnet.
Für die meisten Kartenprojektionen gibt es bereits derartige RTA-Programme. Wenn aber ein nur wenig bekanntes oder neu entwickeltes Kartennetz gerechnet werden soll, so kann man sich selbst ein RTA-Programm schreiben.
Zuerst ein Beispiel
Hier zunächst ein komplettes RT-Assemblerprogramm:
;
; FLÄCHENTREUER ZYLINDERENTWURF
; =============================
;
; Deklarationen
; -------------
;
; Deklaration des Programmnamen (empfohlen, aber nicht unbedingt erforderlich)
;
_name Flächentreuer~Zylinderentwurf
;
; Variabendeklarationen (empfohlen, aber nicht unbedingt erforderlich)
;
_var initial ; Wenn nicht 0, ist das Programm initialisiert
_var phi ; Geographische Breite
_var lambda ; Geographische Länge
_var scale ; Kartenmaßstabszahl (also z. B. 10000000)
;
; x, y, x', y', Cx', Cy', Rx’, Ry’, °(, (° u. a. sind vordefinierte Variablen
;
; Initialbereich (wird nur einmal durchlaufen)
; --------------------------------------------
;
tstne initial go$ ; Wenn initial nicht 0, dann Sprung nach go$
input scale Maßstabszahl ; Maßstabszahl scale mit Dialogtext "Maßstabszahl"
; abfragen
clip scale 1 1E12 ; Maßstabszahl auf 1 ... 10^12 begrenzen
; (Die Maßstabszahl darf nicht Null sein.)
mov initial 1 ; Variable initial auf 1 setzen
go$: ; Definition der Marke go$
;
; Laufbereich (wird für jedes Pixel einmal durchlaufen)
; -----------------------------------------------------
;
; Maßstab, Kartenmittelpunkt etc. einrechnen ("Einheitskugelreduktion")
;
sub x Cx' ; Bildmittelpunkt von x subtrahieren
div x Rx' ; x durch Erdradius dividieren
mul x scale ; x mit Kartenmaßstabszahl multiplizieren
sub y Cy' ; Bildmittelpunkt von y subtrahieren
div y Ry' ; y durch Erdradius dividieren
mul y scale ; y mit Kartenmaßstabszahl multiplizieren
;
; Eigentlicher Entwurf, dieser invers
;
mov lambda x ; Geographische Länge von x übernehmen
mov phi y ; Geographische Breite von y übernehmen
asin phi ; Inversen Sinus (Arkussinus) von phi berechnen
errjump out$ ; Wenn arcsin nicht definiert ist: Außerhalb
;
; In Gradmaß umrechnen
;
mul phi (° ; Bogen- in Gradmaß
mul lambda (° ; Bogen- in Gradmaß
;
; Ausserhalbtest
;
cmplt phi -90 out$ ; wenn phi < -90°, Sprung zur Marke out$
cmpgt phi 90 out$ ; wenn phi > 90°, Sprung zur Marke out$
cmplt lambda -180 out$ ; wenn lambda < -180°, Sprung nach out$
cmpgt lambda 180 out$ ; wenn lambda > 180°, Sprung nach out$
;
; Schlussarbeiten
; ---------------
;
; Reguläres Ende
;
mov x' lambda ; lambda auf die Variable x' schaffen
mov y' phi ; phi auf y' schaffen
exit ; Programmende
;
; Ausserhalb-Ende
;
out$: ; Markendefinition von out$
mov x' -9999 ; -9999 auf x' schaffen
mov y' -9999 ; -9999 auf y' schaffen
exit ; Programmende
_end ; Ende des Programmtextes
Über die Programmiersprache RTA
RTA ist eine richtige, turingvollständige Assemblersprache, etwa wie die Assemblersprachen, mit denen früher die Großrechner programmiert wurden – nur nicht so kompliziert. Ein Programm besteht aus Zeilen, jede Zeile enthält einen Befehl, z. B.
mul a 10
Dies ist ein Multiplikationsbefehl, beginnend mit dem Befehlscode mul, gefolgt von zwei Operanden a und 10. Der Befehl multipliziert den Wert einer Variablen a mit 10.
Indem nun aufeinanderfolgende Zeilen hintereinander abgearbeitet werden, entstehen komplette Formeln. Beispiel:
sub x Cx' ; Bildmittelpunkt Cx’ von x subtrahieren
div x Rx' ; x durch Erdradius Rx’ dividieren
mul x scale ; x mit Kartenmaßstabszahl scale multiplizieren
Auch Verzweigungen und Sprünge sind möglich, hierfür gibt es Befehlsfolgen wie z. B.
cmplt phi 180 m1 ; Wenn phi < 180 Sprung zur Marke m1
sub phi 360 ; Ziehe 360 von phi ab (wenn phi > 180)
m1: ; Das ist die Marke m1
Befehle, die mit Unterstrich beginnen, wie _var oder _name sind Pseudobefehle. Diese werden nicht bei der Programmabarbeitung abgearbeitet, sondern vor dem Programmstart.
Der Pseudobefehl
_name Flächentreuer~Zylinderentwurf
gibt z. B. dem Programm den Namen „Flächentreuer Zylinderentwurf“. Anmerkung: Das Leerzeichen wird in RT als Tilde „~“ (Taste AltGr/+) notiert. Ein Leerzeichen würde „Flächentreuer“ und „Zylinderentwurf“ als 2 getrennte Operanden kennzeichnen.
Pseudobefehle, wie
_var phi
deklarieren eine Variable phi, d. h. sie führen sie in das Programm ein. Dies ist nicht unbedingt erforderlich, weil der RTA unbekannte Variablen selbstständig deklariert. Es ist aber guter Programmierstil, alle benutzten Variabeln am Programmanfang zu deklarieren. So erfährt der Leser des Programms die Bedeutung der Variablen.
Insgesamt gibt es 97 verschiedene Befehle. Das ist im Prinzip alles. RTA ist eine ganz einfache Sprache, viel einfacher als Basic, PHP oder Java oder C. Eine genaue Beschreibung passt auf 7 Seiten und ist im Kapitel „Beschreibung der Programmiersprache RTA“ abgedruckt.
Die Koordinatenübergabe
Man erkennt bestimmte Variablen, wie x, y, Rx’, Cy’, die nicht deklariert werden. Dies sind vordefinierte globale Variablen, über die jedes Programm automatisch verfügt.
Nun ist es die Aufgabe eines jeden RT-Assemblerprogrammes die folgende Koordinatentransformation auszuführen:
(x’, y’) = ƒ (x, y).
Die Werte x und y werden auf den vordefinierten Variablen x und y bereitgestellt. Nach Programmbeendigung müssen die errechneten (transformierten) Koordinaten vom Programm auf den Variablen mit den Namen x’ und y’ eingetragen worden sein.
Implementierung der Formel
Formeln, mit denen ein Kartennetzentwurf berechnet wird, finden sich in der Literatur. Zusätzlich sind einige Nebenaufgaben zu lösen.
• Reduktion auf Einheitskugel und Kartenmitte: Zunächst ist zu berücksichtigen, dass die „Lehrbuchformeln“ der Kartennetzentwürfe aus der Literatur auf der Einheitskugel der Mathematik rechnen. Weiterhin ist es praktisch, das Koordinatenzentrum in den Kartenmittelpunkt (und nicht in eine Kartenecke) zu setzen. Der Ansatz hierfür lautet (bei Inverstransformation, s. u.)
Koordinate = ((Koordinate – Kartenmittelpunkt) / Erdradius) · Maßstabszahl
Praktisch erledigen dies die Befehle
sub x Cx' ; Bildmittelpunkt von x subtrahieren
div x Rx' ; x durch Erdradius dividieren
mul x scale ; x mit Kartenmaßstabszahl multiplizieren
sub y Cy' ; Bildmittelpunkt von y subtrahieren
div y Ry' ; y durch Erdradius dividieren
mul y scale ; y mit Kartenmaßstabszahl multiplizieren
Hier sind Cx’ und Cy’ die Koordinaten des Kartenmittelpunktes und Rx’ und Ry’ Erdradien je bereits im Zielbild-Geokoordinatenmaß, die auf vordefinierten globalen Variablen automatisch richtig bereitgestellt werden. (Anmerkung: Die Trennung in Rx’, Cx’ und Ry’, Cy’ hat methodische Gründe, die hier nicht weiter erörtert werden sollen.)
• Inverse Formel: Die „Lehrbuchformel“ z. B. für den flächentreuen Zylinderentwurf lautet:
x = sin(λ) y = φ
Die Projection engine der Hauptreihe, der Transformator „scannt“ die Kartenebene pixelweise „durch“ und berechnet zu jedem Pixel die Erdkoordinate, die gelesen werden muss. Hier ist es wenig hilfreich, Kartenkoordinaten zu berechnen. Vielmehr sind x, y ja bereits gegeben und λ und φ sind unbekannt. Die Netzentwurfsformeln müssen also invertiert werden, d. h. sie sind so umzustellen, dass gegebene und gesuchte Variablen ihre Plätze tauschen.
Im Falle unseres flächentreuen Zylinderentwurfes ist die Inversion leicht möglich. Die inversen Formeln lauten:
λ = arcsin (x) φ = y.
Für RTA ist der Arkussinus kein Problem, denn es gibt nicht nur für die Winkelfunktionen Befehle, sondern auch für deren Umkehrfunktionen. Die Befehlsfolge
mov lambda x ; Geographische Länge von x übernehmen
mov phi y ; Geographische Breite von y übernehmen
asin phi ; Inversen Sinus (Arkussinus) von phi berechnen
rechnet den Entwurf.
• Umrechnung auf das Bogenmaß: Winkelfunktionen und auch Arkusfunktionen rechnen immer im Bogenmaß und nicht im Gradmaß. Für die Umrechnung gilt:
Grad- in Bogenmaß: Wert mit 0,0174532925 multiplizieren,
Bogen- in Gradmaß: Wert mit 57,2957795130 multiplizieren.
Hierfür stellt der RTA-Assembler die beiden vordefinierte Variablen mit den Namen „(°“ und „°(“ bereit. Derartige Namen sind in RTA durchaus zulässig. Variablennamen können aus ganz beliebigen Druckzeichen gebildet werden. Es müssen keinesfalls etwa nur Buchstaben sein. Der Code zur Umrechnung von Bogenmaß in Gradmaß lautet:
mul phi (° ; Bogen- in Gradmaß
mul lambda (° ; Bogen- in Gradmaß
Nun sind die Ergebnisse fertig berechnet und können auf die Zielvariablen x’ und y’ übertragen werden. Das erledigen die beiden mov-Befehle
mov x' lambda ; lambda auf die Variable x' schaffen
mov y' phi ; phi auf die y' schaffen,
Ein abschließendes
exit ; Programmende
beendet das Programm.
Weitere Hinweise zum Programmieren
Um ein neues Netz zu implementieren ist es günstig, ein einfaches vorhandenes Programm von der Vimage-Installations-CD zu nehmen, und nur die eigentlichen Formelbefehle zu ändern. So geht eine Implementation am schnellsten.
Das Programm wird für jeden Bildpunkt einmal abgearbeitet. Diese Arbeitsweise heißt Parallelprozess. Es spart im Parallelprozess Rechenzeit, wenn Konstanten nur einmalig beim ersten Prozessordurchlauf berechnet werden. Dies ist die Aufgabe der Variablen „initial“ im Beispiel. Diese wird im ersten Durchlauf von 0 auf 1 gesetzt. Sofern sie nicht 0 ist, wird der Initialcode übersprungen.
In diesem Code des ersten Prozessordurchlaufes lassen sich auch günstig Dialoge (z. B. mit dem input-Befehl) anordnen. So wird der Kartenmaßstab scale abgefragt.
Zur Programmentwicklung gibt es den RT-Debugger im Extra-Menü. Er enthält verschiedene Schaltflächen, mit denen man sich u. a. Symboltabelle, Codetabelle, ein Laufzeitfehlerprotokoll etc. ansehen kann.
Mit dem Laufzeitfehlerprotokoll kann man sehr gut Laufzeitfehler im Parallelprozess finden. Zu jedem Befehl wird hier ein Fehlerzähler geführt, sowie der letzte aufgetretene Fehlercode aufgezeichnet. So lassen sich auch Fehler finden, die nur bei einigen wenigen Koordinatenwerten im Parallelprozess auftreten.
Weiterhin kann man die Projection engine Generator bei der Programmentwicklung einsetzen um Zwischenergebnisse im Parallelprozess zu testen. Man schaffe diese einstweilen einfach nach x’ und y’. Nach dem Generatorlauf stehen x’ und y’ im 1. und 2. Band des erzeugten Generatorbildes.
Die Direkttransformation
Die Inversion einer Netzentwurfsformel ist eines der Grundprobleme der Theorie der Kartennetze (sogar bis hin zur Differentialgeometrie) und nicht immer möglich.
Für den Fall, dass ein Netz gerechnet werden muss, für welches es keine Inversformel
(λ, φ) = f(x, y)
gibt, sondern nur eine Vorwärtsformel
(x, y) = f -1(λ, φ )
verfügt Vimage über spezielle Projection engines, z. B. Direkttransformator, Direktgenerator, und Direkttranslator.
Auf diesen kann man dann geeignete „Direktprogramme“ abarbeiten, die mit Vorwärtsformeln rechnen.
Direktprogramme rechnen aus Sicht von Vimage ebenfalls die gewöhnliche allgemeine Koordinatentransformation
(x’, y’) = f -1(x, y).
Sie erhalten allerdings nun die geographische Länge λ auf der vordefinierten Variablen x bereitgestellt und die geographische Breite φ auf der Variablen y. Am Programmende muss sich auf der Variablen x’ die Kartenkoordinate x und auf der Variablen y’ die Kartenkoordinate y befinden.
Achtung! Der Buchstabe „x“ hat hier eine zweifache Bedeutung. Beim Programmstart befindet sich die geographische Länge λ auf der Variablen x. Beim Programmende muss sich die x-Kartenkoordinate auf der Variablen x’ befinden. Analoges gilt für y.
Direkttransformationsprogramme müssen den gesamten Rechenablauf umgekehrt ausführen.
Sie dürfen auch nicht mit den gewöhnlichen Projektionsprogrammen verwechselt werden. Es wird empfohlen, ihre Dateinamen mit „direct_…“ beginnen zu lassen.
Das Direkttransformieren besitzt einige Besonderheiten. Bildverarbeitungstheoretisch erfolgt keine Bild-Neuabtastung (Resampling), sondern eine „Altbild-Übertragung“ („Dissampling“). Für das Dissampling braucht man etwas Erfahrung. Es gibt eine „dissamplingtypische“ Schrittweite. Bei zu großer Schrittweite können Bildlöcher entstehen. Siehe hierzu auch das Kapitel „RTA: Einen Kartennetzentwurf rechnen …“ .
(Auszug aus: Vimage-Buch, 4. Aufl., Kapitel 4.3)
-SP-89-00