;
; ABSTANDSTREUER KEGELENTWURF IN NORMALER (POLARER) LAGE MIT 2 LÄNGENTREUEN PARALLELKREISEN
; =========================================================================================
;
; Nr. bei Wagner: 2
; Name:           Abstandstreuer Kegelentwurf mit 2 längentreuen Parallelkreisen in normaler Lage
; Variante:       Polversetzend
; Autor:
; Quelle:         Karlheinz Wagner, Kartographische Netzentwürfe, Leipzig 1949, S. 34ff.
; 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

; Polversetzend/Zentrierend: Bei Kegelentwürfen ist das Rechenzentrum (der Pol) meist 
; nicht der Kartenmittelpunkt. Dies ist ein polversetzendes Programm, d. h. man kann/muss eine 
; Distanz angeben, um die der Pol nach N versetzt wird. Das Gegenstück bilden die zentrierenden 
; Kegelentwurfsprogramme, die Pol so versetzen, dass der Mitte der Quellkarte zur Mitte der Zielkarte wird.

; Benutzte Variablen
; ==================
;
; Die Variablennamen entsprechen weitgehend denen von Karlheinz Wagner,

;
; Laufende Koordinaten
;
	_name	Abstandstreuer~Schnittkegelentwurf~(polversetzend)
	_var	phi		; Geographische Breite
	_var	lambda		; Geographische Länge
	_var	alpha		; Azimut ebene Polarkoord., auch Schiefazimut
	_var	delta		; Poldistanz/geographisch, auch Schiefdistanz
	_var	epsilon		; Delta0 minus Delta
	_var	m		; Radius ebene Polarkoordinaten
	_var	n		; Winkelreduktionsfaktor (Cos(Delta0))
;
; Konstanten der Transformation
;
	_var	phi1		; Geogr. Breite des Schnittparallels 1
	_var	phi2		; Geogr. Breite des Schnittparallels 2
	_var	delta0		; Poldistanz des Mittelparallels zwischen den Berührungsparallelen
	_var	delta1		; Poldistanz des 1. Schnittparallels
	_var	delta2		; Poldistanz des 1. Schnittparallels
	_var	delta0		; Poldistanz des Mittelparallels zwischen den Berührungsparallelen
	_var	epsilon0	; delta0 mius delta1 oder delta2
	_var	m0		; Halbmesser des delta0-Parallels
	_var	lambda0		; Geogr. Breite des Mittelmeridians
	_var	scale		; Kartenmaßstabszahl (also 1000000, nicht 1/1000000)
	_var	pol-y		; Polversatz nach Nord
	_var	2pi		; 360 Grad
	_var	sigma		; Hemisphärenvorzeichen
	_var	.a		; Temporäre Variable
	_var	.b		; Temporäre Variable
;
; 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	lambda0	Mittelpunktlänge~in~Grad
	input	pol-y	Pollage~nördlich~vom~Bildmittelpunkt~in~Geokoordinaten
; Eingegebene Werte auf Min/Max bringen
	clip	scale	1	1E12
	clip	phi1	-90	90
	clip	phi2	-90	90
	clip	lambda0	-180	180
	clip	pol-y	-1E6	1E6
; Sigma
	mov	sigma	phi1
	add	sigma	phi2
	sgn	sigma
	mul	phi1	sigma
	mul	phi2	sigma
; Konstanten berechnen
	mul	phi1	°(
	mul	phi2	°(
	mov	delta1	pi/2
	sub	delta1	phi1
	mov	delta2	pi/2
	sub	delta2	phi2
	clr	delta0
	add	delta0	delta1
	add	delta0	delta2
	div	delta0	2		; Delta0 fertig

	mov	epsilon0	delta0
	sub	epsilon0	delta1
	abs	epsilon0		; epsilon0 fertig

	mov	.a	delta0
	tan	.a
	mov	m0	epsilon0
	cot	m0
	mul	m0	.a
	mul	m0	epsilon0		; m0 fertig

	mov	.a	epsilon0
	sin	.a
	mov	n	delta0
	cos	n
	mul	n	.a
	div	n	epsilon0		; n fertig
; 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'
	sub	y	pol-y
	div	y	Ry'
	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	pi
	sub	alpha	pi
$153:
	power	x	2
	power	y	2
	clr	m
	add	m	x
	add	m	y
	root	m	2
;
; Eigentlicher abstandstreuer Entwurf, dieser invers
; --------------------------------------------------
;
	mov	lambda	alpha		; Geographische Länge
	div	lambda	n

	mov	epsilon	m0		; Geographische Breite
	sub	epsilon	m

	mov	delta	delta0
	sub	delta	epsilon

	mov	phi	pi/2
	sub	phi	delta
;
; Kontern, in Gradmaß umrechnen und Ausserhalbtest
; ------------------------------------------------
;
	mul	lambda	sigma		; Südhemisphäre:
	mul	phi	sigma		; ggf. kontern
	mul	lambda	(°
	mul	phi	(°
	cmplt	lambda	-180	out
	cmpgt	lambda	180	out
	cmplt	phi	-90	out
	cmpgt	phi	90	out
;
; 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