;
; GYÖRFFY E1
; ==========
;
; Name:            Orthogonale pseudopolykonische Abbildung nach Györffy 2004
; Kurzbezeichnung: Györffy E1
; Quelle:          Györffy, Janos, Klinghammer, István, Budapest: Die kleinmaßstäbige Darstellung
;                  von Europa in Projektionen geringster Verzerrung - die Basis für einen Europa-Atlas
;                  Kartographische Nachrichten, Bonn-Bad Godesberg: Kirschbaum 6/2004.
; Richtung:        Direkt-Transformation
;
; 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.
;
; Ein Maßstab und der Mittelmeridian wird abgefragt. Der Mittelmeridian ist bei Györffy 25° östlicher Länge.

; Die Györffy-Entwürfe umfassen 3 Netze, die alle für den europäischen Raum 
; im Bereich -10° ... 60° geographischer Länge und 35° ... 70° geographischer Breite 
; optimier wurden:

; "E1":  5 Koeffizienten. Orthogonal pseudopolykonisch. Rechtwinklig. Verzerrungsmaß Ek=0.0082
; "E2":  8 Koeffizienten. Allgemein pseudopolykonisch.  Symmetrisch.  Verzerrungsmaß Ek=0.0076
; "E3": 17 Koeffizienten. Idealabbildung. Unsymmetrisch.              Verzerrungsmaß Ek=0.0070

; Zum Vergleich: Bonne Ek=0.0192, de I'sle Ek=0.0100.

; Entwurf (C) János Gyöffry 2004, Implementation: (C) Rolf Böhm 2005

; Benutzte Variablen
; ==================
;
; Die Variablennamen entsprechen weitgehend denen von Györffy.

;
; Laufende Koordinaten
;
	_name	Györffys~pseudopolykonische~Abbildung~E1~(orthogonal)
	_var	phi		; Geographische Breite
	_var	lambda		; Geographische Länge
	_var	lambda0		; Mittelmeridian
	_var	beta		; Poldistanz/geographisch, auch Schiefdistanz (sonst meist delta)
	_var	gamma		; Polarkoordinaten/Winkel (sonst meist alpha)
	_var	r		; Polarkoordinaten/Radius
	_var	c		; y-Verschiebung des (polykonischen) Poles
	_var	lambda^3	; lambda^3
	_var	beta^2		; beta^2
	_var	expo		; gamma-Exponent (-c1/r1)
	_var	t1		; temporäre Variablen
	_var	t2
	_var	t3
	_var 	t4
;
; Konstanten der Transformation
;
	_var	scale		; Kartenmaßstabszahl (also 1000000, nicht 1/1000000)
	_var	r1		; r1
	_var	r2		; r2
	_var	c1		; c1
	_var	f1		; f1
	_var	f2		; f2
	_var	........	; End of Symbol table
;
; 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	Mittelmeridian~in~Grad~(Gyöffry:~25) 
; Eingegebene Werte auf Min/Max bringen
	clip	scale	1	1E12
	clip	lambda0	-180	180
; Koeffinzienten
	mov	r1	1.164017
	mov	r2	0.001802
	mov	c1	-0.993478
	mov	f1	0.492496
	mov	f2	0.026293
	mov	expo	c1
	div	expo	r1
	neg	expo

; Programm ist initialisiert
	mov	initial	1
077$:

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

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

	sub	lambda	lambda0		; Mittelmeridian 25° Ost
	cmpgt	lambda	-180	10$
	add	lambda	360
10$:
	cmpgt	lambda	-180	30$
	add	lambda	360
30$:
	cmplt	lambda	180	40$
	sub	lambda	360
40$:
	cmplt	lambda	180	50$
	sub	lambda	360
50$:

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

; Netzentwurf rechnen
; -------------------

; Poldistanz beta berechnen

	mov	beta	pi/2
	sub	beta	phi
	
; Polarkoordinatenwinkel gamma berechnen - Formel S. 267 Mitte

	mov	lambda^3	lambda		; vorab ...
	power	lambda^3	3
	mov	beta^2		beta
	power 	beta^2		2	

	mov	t1	beta			; linker Ausdruck
	mov	t2	beta
	mul	t2	r2
	add	t2	r1
	div	t1	t2
	power	t1	expo

	mov	t2	lambda			; rechter Ausdruck, Zähler	
	mul	t2	f1
	mov	t3	lambda^3
	mul	t3	f2
	add	t2	t3

	mov	t3	beta			; rechter Ausdruck, Nenner
	mul	t3	r1	
	mov	t4	beta^2
	mul	t4	r2	
	add	t3	t4
	
	div	t2	t3			; rechter Ausdruck, Hauptdivision
	mul	t1	t2			; linker x rechter Ausdruck
	atan	t1
	mul	t1	2
	mov	gamma	t1	

; r berechnen

	mov	t1	beta
	mul	t1	r1
	mov	r	beta^2
	mul	r	r2
	add	r	t1

; c berechnen

	mov	c	beta
	mul	c 	c1
	add	c	r

; x berechnen

	mov	t1	gamma
	sin	t1
	mul	t1	r
	mov	x	t1

; y berechnen

	mov	t1	gamma
	cos	t1
	mul	t1	r
	mov	y	c
	sub	y	t1 		
;
; 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		; auf "Gesucht"
	mov	y'	y               ; auf "Gesucht"
	exit
	_end