; ; FLÄCHENTREUER KEGELENTWURF IN SCHIEFER (FREIER) LAGE MIT 2 LÄNGENTREUEN PARALLELKREISEN ; ======================================================================================= ; ; Nr. bei Wagner: 4 ; Name: Flächentreuer Kegelentwurf auf dem Schnittkegel oder ; Flächentreuer Kegelentwurf mit 2 längentreuen Parallelkreisen in normaler Lage ; (Sog. Albersscher Kegelentwurf oder Albersscher flächentreuer Kegelentwurf) ; Variante: Zentrierend ; Autor: ; Quelle: Karlheinz Wagner, Kartographische Netzentwürfe, Leipzig 1949, S. 47ff. ; Richtung: Inverse Transformation ; ; Literatur: ; Wagner: Kartographische Netzentwürfe, Leipzig: Bibliographisches Institut 1949 ; Fiala: Mathematische Kartographie, Berlin: Verlag Technik 1957 ; (C) Rolf Böhm Bad Schandau 2004, 2008 ; Polversetzend/(Selbst-)Zentrierend: Bei Kegelentwürfen ist das Rechenzentrum (der Pol) meist ; nicht der Kartenmittelpunkt. Dies ist ein (selbst)zentrierendes Programm, welches den Pol ; so versetzt, dass der Mitte der Quellkarte zur Mitte der Zielkarte wird. Das Gegenstück bilden ; polversetzende Polgramme, bei denen der Polversatz "zu Fuss" abgefragt wird. ; ; Anmerkung 1: Die Berührungsparallel oder die Schnittparallele haben k e i n e n ; Einfluss auf die Zentrierung. Entscheidend ist der Quellbildmittelparallel ; ; Anmerkung 2: Auf die X-Lagerung hat die Selbstzentrierung ebenfalls keinen Einfluss, ; dies ist, wie bei allen anderen Programmen, mit lambda0 einzustellen. ; ; Anmerkung 3: Die Selbstzentrierung ist ganz schön aufwändig. Immerhin muss ; dafür zusätzlich zu der Inversformel eine Vorwärtsformel implementiert ; werden (wenngleich diese auch nur ein einziges Mal gerechnet wird) ; ; Anmerkung 4: Schiefachsige Variante. Der außerhalbtest wird hier weggelassen, um evtl. ; erwünschte Randüberschreitungen nicht abzuschneiden. ; ; Benutzte Variablen ; ================== ; ; Die Variablennamen entsprechen weitgehend denen von Karlheinz Wagner, ; ; Laufende Koordinaten ; _name Flächentreuer~Schnittkegelentwurf~(zentrierend)~in~freier~Lage _var phi ; Geographische Breite _var lambda ; Geographische Länge _var alpha ; Azimut ebene Polarkoord. _var m ; Radius ebene Polarkoord. _var n ; n _var delta ; Poldistanz/geographisch _var phi1 ; Geogr. Breite des 1. Berührungsparallels _var phi2 ; Geogr. Breite des 2. Berührungsparallels _var delta1 ; Poldistanz des 1. Berührungsparallels _var delta2 ; Poldistanz des 2. Berührungsparallels _var C ; Die additive Konstante in Formel (4) Wagner S. 53 unten) _var sigma ; Hemisphärenvorzeichen _var t1 ; Temporäre Variable 1 _var lambda0 ; Geogr. Breite des Mittelmeridians _var poly ; Polverschiebung _var 2pi _var scale ; Kartenmaßstabszahl (also 1000000, nicht 1/1000000) ; ; Für Transformation in freie Lage ; _var phi0 _var sinphi0 _var cosphi0 _var phi' ; Kartographische Breite _var lambda' ; Kartographische Länge in 4 Varianten _var lambda" _var lambda- _var lambda= _var d11 ; Differenzen _var d12 _var d21 _var d22 _var dx ; Differenzen »gekreuzt« _var d= ; Differenzen »even« ; ; x, y, x', y', Cx', Cy', Rx', Ry', °(, (°, pi, pi/2 etc. sind vordefinierte globale Konstanten ; ; Initialisierung ; =============== ; tstne initial 077$ ; Dialog input scale Maßstabszahl input phi1 1.~Schnittparallelbreite~in~Grad input phi2 2.~Schnittparallelbreite~in~Grad input phi0 Polbreite~in~Grad input lambda0 Mittelpunktlänge~in~Grad ; Eingegebene Werte auf Min/Max bringen clip scale 1 1E12 clip phi1 -90 90 clip phi2 -90 90 clip phi0 -90 90 clip lambda0 -180 180 ; SinCosPhi0 mul phi0 °( mov sinphi0 phi0 sin sinphi0 mov cosphi0 phi0 cos cosphi0 ; Sigma mov sigma phi1 add sigma phi2 sgn sigma mul phi1 sigma mul phi2 sigma ; Konstanten berechnen mul phi1 °( mov delta1 pi/2 sub delta1 phi1 ; Delta1 fertig mul phi2 °( mov delta2 pi/2 sub delta2 phi2 ; Delta2 fertig mov t1 delta1 cos t1 mov n delta2 cos n add n t1 div n 2 ; n fertig (Formel unter (4) Wagner S. 53 unten) mov C 4 mov t1 n power t1 2 div C t1 mov t1 delta1 div t1 2 sin t1 power t1 2 mul C t1 mov t1 delta2 div t1 2 sin t1 power t1 2 mul C t1 ; C fertig (Summand in der Formel (4)) mov 2pi pi add 2pi pi ; 2pi fertig ; ; Das "Fangen" des Bildes: Es wird ein poly berechnet, das den Pol so verschieben wird, dass ; die Mitte des "Erdbildes" in Kartenmitte zu liegen kommt. ; Hierzu wird der Quellbildmittelpunkt komplett vorwärts transformiert ... ; mov t1 Cy ; "Phi" von Cy abholen. "Lambda" alias Cx sei immer 0 ... mul t1 sigma ; ;;;;; mul t1 °( ; "Phi" in Bogenmaß umrechnen ; mov t2 pi/2 ; Eigentlicher Entwurf: Phi (in t1) im m (alias poly) umrechnen sub t2 t1 ; pi/2 - phi div t2 2 sin t2 power t2 2 mul t2 4 div t2 n add t2 C root t2 2 mov poly t2 ; ; m in y umrechnen. Aus x=0 folgt, da m=sqrt(x²+y²) m=y! div poly scale ; Maßstab einrechnen mul poly Ry' ; Erdradius einrechnen mul poly sigma ; hemisphärenabhängig proof poly ; Programm ist initialisiert mov initial 1 077$: ; ; SIMD-Laufbereich ; ================ ; ; Maßstab, Kartenmittelpunkt etc. einrechnen ; ------------------------------------------ ; sub x Cx' ; Bildmittelpunkt div x Rx' ; Erdradius mul x scale ; Kartenmaßstab sub y Cy' ; Bildmittelpunkt sub y poly ; Erdradius div y Ry' ; Kartenmaßstab mul y scale mul x sigma mul y sigma ; ; Kartesische Koordinaten x/y in Polarkoordinaten m/alpha umwandeln ; ----------------------------------------------------------------- ; mov alpha x ; kartesische Koord. x y in Polarkoordinaten m alpha div alpha y atan alpha tstgt y 151$ ; Der atan ist doppeldeutig!: Wenn y negativ dann add alpha pi ; 180° addieren 151$: neg alpha ; Jetzt ist alpha berechnet, liegt aber noch falsch add alpha pi cmplt alpha pi 153$ sub alpha 2pi 153$: power x 2 power y 2 clr m add m x add m y root m 2 ; ; Eigentlicher flächentreuer Entwurf, dieser invers ; ------------------------------------------------- ; ; div 0 0 ; Testhalt mov lambda alpha ; Geographische Länge div lambda n mov t1 m ; Geographische Breite power t1 2 sub t1 C mul t1 n div t1 4 root t1 2 errjump out asin t1 errjump out mul t1 2 mov delta t1 mov phi pi/2 sub phi delta ; Uff! Das ist die Formel (4) Wagner S. 53 invertiert ; ; SÜdhalbkugel erforderlichenfalls kontern ; ---------------------------------------- ; mul lambda sigma ; Südhemisphäre: mul phi sigma ; ggf. kontern ; ; Transformation in schiefe (freie) Lage: Umrechnung der polständigen phi/lambda in schiefe phi'/lambda'-Koordinatenpaare ; ----------------------------------------------------------------------------------------------------------------------- ; ; Das Problem ist, dass die Rechnung über arcsin läuft und dieser doppeldeutig ist. ; Die Lösung: Die Sache wird ein zweitesmal mit einer Formal angegangen, die einen ; arccos liefert. Dieser ist nun zwar auch doppeldeutig. Es ist aber nun diejenige Lösung ; richtig, die sowohl arcsin als auch arccos liefern. Wegen Rundungsunsicherheiten wird ; nicht auf Identität, sondern auf Differenzminimum getestet. ; Vorab: sin/cos berechnen mov sinl lambda sin sinl ; Sin(Lambda) mov cosl lambda cos cosl ; Cos(Lambda) mov sinp phi sin sinp ; Sin(Phi) mov cosp phi cos cosp ; Cos(Phi) ; phi' berechnen mov t1 sinp ; phi' berechnen (Formeln Fiala, S. 79f.) mul t1 sinphi0 mov t2 cosp mul t2 cosphi0 mul t2 cosl add t1 t2 mov phi' t1 asin phi' ; Lambda (1): lambda' berechnen: Mittels Arcsin mov t1 cosp neg t1 mul t1 sinl mov t2 phi' cos t2 div t1 t2 mov lambda' t1 asin lambda' ; Lambda (2): lambda- berechnen: Die zweite Lösung, die der Arcsin liefert tstlt lambda' 201$ mov lambda- pi sub lambda- lambda' jump 203$ 201$: clr lambda- sub lambda- pi sub lambda- lambda' 203$: ; Lambda (3): lambda" berechnen: Mittels Arccos mov t1 cosphi0 neg t1 mul t1 sinp mov t2 sinphi0 mul t2 cosp mul t2 cosl mov t3 phi' cos t3 add t1 t2 div t1 t3 mov lambda" t1 acos lambda" ; Lambda (4): lambda= berechnen: Die zweite Lösung, die der arccos liefert mov lambda= lambda" ; Gleichheitszeichen lies als Doppelstrich quer neg lambda= ; Lambda (5): Ein Lambda aus demjenigen Paar ist richtig, welches die kleinste Differenz aufweist ; Differenzen für 2. Entscheidungsstufe berechnen mov d11 lambda' sub d11 lambda" abs d11 mov d12 lambda' sub d12 lambda= abs d12 mov d21 lambda- sub d21 lambda" abs d21 mov d22 lambda- sub d22 lambda= abs d22 ; Differenzen für 1. Entscheidungsstufe berechnen mov dx d12 ; "gekreuzt" mul dx d21 mov d= d11 ; "even" mul d= d22 ; 1. Entscheidungsstufe: gekreuzt oder even cmpgt dx d= =richtig ; 2. Entscheidungsstufe: Einzellösung xrichtig: cmpgt d12 d21 21richtig jump 12richtig =richtig: cmpgt d11 d22 22richtig jump 11richtig ; Lösungen zusammensuchen 21richtig: mov lambda' lambda- jump ok 22richtig: mov lambda' lambda- 12richtig: 11richtig: ok: ; ; In Gradmaß umrechnen ; -------------------- ; mul lambda' (° mul phi' (° neg lambda' ; ** Empirisch bemerkt: Bild sonst seitenverkehrt ; ; Schlussarbeiten ; --------------- ; mov x' lambda' mov y' phi' ; Lambda kreisen lassen add x' lambda0 cmod x' -180 180 exit out: mov x' -9999 mov y' -9999 exit _end