;
; ROBINSONPROJEKTION POLYNOMISCH
; ==============================
;
; Name:           Robinsonprojektion
; Quelle:         Esque, Alexandre; Bergada, Alex 2004: 
;                 Cataleg de projeccions cartografiques per la mapamundis. [Barcelona]
;                 http://mercator.upc.es/Cataleg/PROJECTE.HTM
; Richtung:       Direkt-Transformation
;
; Von der Robinson-Projektion gibt es zwei bekannte Ansätze, einmal tabellierte
; Koeffizienten, einmal Polynome. Beide Varianten sind sehr ähnlich, weichen aber
; in höhener Breiten etwas (grob ca. 1 Grad) voneinander ab.
;
; Dies ist Version mit Polynomen.
;
; Das Programm übernimmt die Koordinaten eines Punktes (x/y) und transformiert
; diese in einen Punkt (x'/y').
;
; x/y sind ebene Zielpunktkoordinaten, x'/y' geben die geogr. Breite und Länge
; der Position auf der Quell-Erdkugel, auf der der Zielpunkt gelesen werden
; kann.

; (C) Dr.-Ing. Rolf Böhm 2005

; Benutzte Variablen
; ==================

;
; Laufende Koordinaten
;
	_name	Robinsonprojektion~(polynomisch)
	_var	phi		; Geographische Breite
	_var	lambda		; Geographische Länge
	_var	oldphi
;
; Konstanten der Transformation
;
	_var	lambda0		; Geogr. Länge des Bildmittelpunktes
	_var	scale		; Kartenmaßstabszahl (also 1000000, nicht 1/1000000)
	_var	p0		; Polynomterme 0 bis 5
	_var	p1
	_var	p2
	_var	p3
	_var	p4
	_var	p5
;
; x, y, x', y', Cx', Cy', °(, (°, pi, pi/2 etc. sind vordefinierte globale Konstanten
;
; Initialisierung
; ===============
;
	tstne	initial	077$
; Dialog
	pause	Hinweis:~Dieses~Programm~rechnet~eine~Vorwärtstransformation.\\Es~muss~mit~einer~direkt~arbeitenden~Projection~engine~abgearbeitet~werden.
	input	scale	Maßstabszahl
	input	lambda0	Mittelpunktslänge~in~Grad
; Eingegebene Werte auf Min/Max bringen
	clip	scale	1	1E12
	clip	lambda0	-180	180
; Konstanten berechnen
	mov	cosphi0	37.5
	cos	cosphi0
; Programm ist initialisiert
	mov	initial	1
077$:

;
; SIMD-Laufbereich
; ================

;
; Eigentlicher Entwurf, dieser direkt!
; ------------------------------------
;
	mov	lambda	x		; Geographische Länge
	mov	phi	y		; Geographische Breite

	sub	lambda	lambda0	
	cmod	lambda	-180	180

; x, y haben eine Doppelbedeutung:
; - einmal die verdefinierten RTA-Input-Koordinaten, also eigentlich Phi und Lambda,
; - dann aber auch die ebenen Kartenkoordinaten, die errechnet werden

;
; Umrechnung in Bogenmaß
; ----------------------
;
	mul	phi	°(
	mul	lambda	°(

; Netzentwurf rechnen
; -------------------
;
	cmpeq	oldphi	phi	1101$
	; die Potenzen von phi berechnen	
	mov	p0	1
	mov	p1	phi
	mov	p2	phi
	mov	p3	phi
	mov	p4	phi
	mov	p5	phi
	power	p2	2
	power	p3	3
	power	p4	4
	power	p5	5
	; diese mit den Koeffizienten multiplizieren
	mul	p0	0.8507
	mul	p1	0.9642
	mul	p2	-0.1450
	mul	p3	-0.0013
	mul	p4	-0.0104
	mul	p5	-0.0129
	mov	oldphi	phi
1101$:
	; x berechnen
	clr	x	
	add	x	p0
	add	x	p2
	add	x	p4
	mul	x	lambda
	; y berechnen
	clr	y	
	add	y	p1
	add	y	p3
	add	y	p5	
	;
;
; Maßstab, Kartenmittelpunkt etc. einrechnen
; ------------------------------------------
;
	mul	x	Rx'		; Erdradius
	div	x	scale		; Kartenmaßstab
	add	x	Cx'

	mul	y	Ry'
	div	y	scale
	add	y	Cy'
;
; Schlussarbeiten
; ---------------
;
111$:
	mov	x'	x
	mov	y'	y
	exit
out:
	mov	x'	-9999
	mov	y'	-9999
	exit
	_end