;
; BARTHOLOMEW STAR
; ================
;
; Name:           Bartholomew Star
; Quelle:         Zitiert nach Carlos A. Furuti:
;                 http://www.progonos.com/furuti/MapProj/Normal/ProjInt/ProjStar/projStar.html:  
;                 "Reconstruction of a North pol tetrahedral projection, interrupting 
;                 started at 30° W" 
; Richtung:       Direkt-Transformation
;
; Take to compute Furutis Barthelemew Tetrahedral Reconstruction:
;
;   Sternpitzenzahl:      3
;   Grenzbreite:          23.5
;   Mittelounktslänge:   -30
;
; Das Programm erlaubt Südpolanzahlen von 3 bis 3600.
;
; *) Die beiden mit *) gekennzeichneten Befehle sorgen dafür, dass der Stern nicht
;    auf der Spitze steht, sondern auf "2 Zacken". 
;
; 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.

; Eigenentwicklung

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

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

;
; Laufende Koordinaten
;
	_name	Bartholomew~Star~(»Tetrahedral«)
	_var	phi		; Geographische Breite
	_var	lambda		; Geographische Länge
	_var	alpha		; Polarkoordinatenazimut
	_var	m		; Polarkoordinatenradius
;
; Transformationsparameter
;
	_var	lambda0		; Geogr. Länge des Bildmittelpunktes
	_var	scale		; Kartenmaßstabszahl (also 1000000, nicht 1/1000000)
	_var	Sternspitzenzahl
	_var	Grenzbreite
;
	_var	bon_da		; d Alpha in the Bonne (Werner)
	_var	azi_da		; d Alpha equaldistant azimuthal
;
; Konstanten
;
	_var	2pi
;
; Segmentberechnung
;
	_var	s-phase		; Phasenlage des Mittelmeridians des aktuellen Segments
	_var	s-off		; Restwinkel zwischen s-phase und aktueller Länge
	_var	s-width		; Segmentbreite 
	_var	s-halfw		; Halbe Segmentbreits
;
; 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	Sternspitzenzahl Anzahl~der~Sternspitzen~(3~...~3600,~Bartholomew:~3)
	input	Grenzbreite      Grenzbreite~(Pseudo-Äquator)~in~Grad~(-90~...~90,~Bartholomew:~23.5)
	input	lambda0	Mittelpunktslänge~in~Grad~(Bartholomew:~-30)
; Eingegebene Werte auf Min/Max bringen
	clip	scale			1	1E12
	clip	Sternspitzenzahl	3	3600
	clip	Grenzbreite		-90	90	
	clip	lambda0			-180	180
; In Bogenmaß
	mul	Grenzbreite	°(
; Reduktion nach auf Grenzparallel
	mov	r4	pi/2
	sub	r4	Grenzbreite	; r3: delta(Grenzbreite)
	mov	r3	r4
	sin	r3			; d-alpha Bonne
	mov	r2	r4
	div	r2	pi/2		; d-alpha Azimutal abstandstreu
	mov	r1	r2
	div	r1	r3		; d-alpha-B / d-alpha-A --> r1
; 2 Pi berechnen
	mov	2pi	pi
	mul	2pi	2
; Segmentierung vorbereiten	
	mov	s-width	pi
	mul	s-width	2
	div	s-width	Sternspitzenzahl
	mov	s-halfw	s-width
	div	s-halfw	2
; Programm ist initialisiert
	mov	initial	1
077$:

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

;
; Koordinaten übernehmen
; ----------------------
;
	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	°(
;
; Segmentierung
; -------------
;
	sub	lambda	s-halfw		; Verdrehung um 1/2 Zacken *)

	mov	s-phase	lambda		; Phase (des Zackens) ist Ganzzahlteil
	add	s-phase	2pi		; das vermeidet negative Winkel
	add	s-phase	s-halfw		; So rundet der fix: -0.5 ... 0.5 -> 0 ... 1
	div	s-phase	s-width		; Normalisieren
	fix	s-phase			; Ganzteil
	mul	s-phase	s-width		; Denormalisieren

	mov	s-off	lambda		; Restwinkel (im Zacken) ist gebrochener Teil
	add	s-off	2pi		; das vermeidet negative Winkel
	add	s-off	s-halfw		; So rundet der frac: -0.5 ... 0.5 -> 0 ... 1
	div	s-off	s-width		; Normalisieren
	frac	s-off			; Bruchteil	
	mul	s-off	s-width		; Denormalisieren

	mov	lambda	s-off		; Restwinkel ist neues Lambda
	sub	lambda	s-halfw		; (0 ... 1 -> -0.5 ... 0.5)
;
; Netzentwurf rechnen
; -------------------
;
; phi ... ganz einfach ...
; - - - - - - - - - - - - -
;
	mov	m	pi/2
	sub	m	phi
;
; lambda ... N und S unterscheiden ...
; - - - - - - - - - - - - - - - - - - -
;	
	cmplt	phi	Grenzbreite	south
north:
	mov	alpha	lambda
	jump	go
south:
	mov	delta	pi/2
	sub	delta	phi
	mov	bon_da	delta
	sin	bon_da			
	mov	azi_da	delta
	div	azi_da	pi/2
	mov	r5	bon_da
	div	r5	azi_da	

	mov	alpha	lambda
	mul	alpha	r5		; Reduktion im differentiellen Element
	mul	alpha	r1		; Reduktion auf Grenzbreite
go:
; Abschlussdrehungen
	add	alpha	s-phase		; Phasenverdrehung
	add	alpha	s-halfw		; Zurückdrehung um 1/2 Zacken *)
	sub	alpha	pi/2		; Restfehler 90 Grad korrigieren
;
; Polarkoordinaten in kartesische Koordinaten umrechnen
; -----------------------------------------------------
;
	mov	x	alpha	
	cos	x
	mul	x	m
	mov	y	alpha
	sin	y
	mul	y	m
;
; 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
; ---------------
;
	mov	x'	x
	mov	y'	y
	exit
	_end