;
; CANTERS LOW ERROR POLYCONIC C40
; ===============================
;
; Name:           Canters Low-error polyconic projection obtained through non-constrained
;                 optimization with 40 coefficients
; Quelle:         Canters, Small-scale Map Projection Design, Taylor & Francis, London 2002,
;                 p. 190, 191.
; 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.

; (C) Rolf Böhm 2006

; Benutzte Variablen
; ==================
;
; Die Variablennamen entsprechen weitgehend denen von Frank Canters

;
; Laufende Koordinaten
;
	_name	Canters~Low~Error~Polyconic~C40
;
	_var	phi
	_var	lambda
;
	_var	c10
	_var	c01
	_var	c20
	_var 	c11
	_var	c02
	_var	c30
	_var	c21
	_var	c12	
	_var	c03
	_var	c40
	_var	c31
	_var	c22
	_var	c13
	_var	c04
	_var	c50
	_var	c41
	_var	c32
	_var	c23
	_var	c14
	_var	c05
;
	_var	c'10
	_var	c'01
	_var	c'20
	_var 	c'11
	_var	c'02
	_var	c'30
	_var	c'21
	_var	c'12	
	_var	c'03
	_var	c'40
	_var	c'31
	_var	c'22
	_var	c'13
	_var	c'04
	_var	c'50
	_var	c'41
	_var	c'32
	_var	c'23
	_var	c'14
	_var	c'05
;
	_var	s0		; The y-Register
	_var	s1
	_var	s2
	_var	s3
	_var	s4
	_var	s5
;
; Konstanten der Transformation
;
	_var	lambda0		; Geogr. Länge des Bildmittelpunktes
	_var	scale		; Kartenmaßstabszahl (also 1000000, nicht 1/1000000)
;
; 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
; Koeffizienten initialisieren
	mov	c10	+1.0542
	mov	c01	-0.0173
	mov	c20	+0.0088
	mov	c11	-0.4166
	mov	c02	+0.0036
	mov	c30	-0.0449
	mov	c21	+0.0010
	mov	c12	-0.0432
	mov	c03	+0.0055
	mov	c40	+0.0001
	mov	c31	+0.0144
	mov	c22	-0.0007
	mov	c13	-0.0499
	mov	c04	+0.0053
	mov	c50	+0.0007
	mov	c41	-0.0006
	mov	c32	-0.0042
	mov	c23	-0.0020
	mov	c14	+0.0388
	mov	c05	+0.0017
;
	mov	c'10	-0.0009
	mov	c'01	+0.9744
	mov	c'20	+0.2607
	mov	c'11	-0.0040
	mov	c'02	+0.0712
	mov	c'30	+0.0005
	mov	c'21	-0.0708
	mov	c'12	-0.0002
	mov	c'03	-0.0694
	mov	c'40	-0.0081
	mov	c'31	+0.0001
	mov	c'22	-0.0146
	mov	c'13	+0.0004
	mov	c'04	-0.0207
	mov	c'50	-0.0001
	mov	c'41	+0.0047
	mov	c'32	+0.0003
	mov	c'23	-0.0298
	mov	c'14	+0.0014
	mov	c'05	+0.0491

; 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	
	clip	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
; -------------------
;
	mov	r0	1
	mov	r1	lambda
	mov	r2	lambda
	mov	r3	lambda
	mov	r4	lambda
	mov	r5	lambda

	mov	s0	1
	mov	s1	phi
	mov	s2	phi
	mov	s3	phi
	mov	s4	phi
	mov	s5	phi

	power	r2	2
	power	r3	3
	power	r4	4
	power	r5	5
	
	power	s2	2
	power	s3	3
	power	s4	4
	power	s5	5
;
; The X-Sum
;
	clr	x
	mov	a	c10	
	mul	a	r1
	add	x	a

	mov	a	c01
	mul	a	s1
	add	x	a

	mov	a	c20
	mul	a	r2
	add	x	a

	mov	a	c11
	mul	a	r1
	mul	a	s1
	add	x	a

	mov	a	c02
	mul	a	s2
	add	x	a

	mov	a	c30
	mul	a	r3
	add	x	a

	mov	a	c21
	mul	a	r2
	mul	a	s1
	add	x	a		

	mov	a	c12
	mul	a	r1
	mul	a	s2
	add	x	a

	mov	a	c03
	mul	a	s3
	add	x	a

	mov	a	c40
	mul	a	r4
	add	x	a	
	
	mov	a	c31
	mul	a	r3
	mul	a	s1
	add	x	a
	
	mov	a	c22
	mul	a	r2
	mul	a	s2
	add	x	a

	mov	a	c13
	mul	a	r1
	mul	a	s3
	add	x	a
	
	mov	a	c04
	mul	a	s4
	add	x	a

	mov	a	c50
	mul	a	r5	
	add	x	a

	mov	a	c41
	mul	a	r4
	mul	a	s1
	add	x	a

	mov	a	c32
	mul	a	r3
	mul	a	s2
	add	x	a

	mov	a	c23
	mul	a	r2
	mul	a	s3
	add	x	a

	mov	a	c14
	mul	a	r1
	mul	a	s4
	add	x	a

	mov	a	c05
	mul	a	s5
	add	x	a	; that's x
;
; The Y-Sum
;
	clr	y

	mov	a	c'10	
	mul	a	r1
	add	y	a

	mov	a	c'01
	mul	a	s1
	add	y	a

	mov	a	c'20
	mul	a	r2
	add	y	a

	mov	a	c'11
	mul	a	r1
	mul	a	s1
	add	y	a

	mov	a	c'02
	mul	a	s2
	add	y	a

	mov	a	c'30
	mul	a	r3
	add	y	a

	mov	a	c'21
	mul	a	r2
	mul	a	s1
	add	y	a		

	mov	a	c'12
	mul	a	r1
	mul	a	s2
	add	y	a

	mov	a	c'03
	mul	a	s3
	add	y	a

	mov	a	c'40
	mul	a	r4
	add	y	a	
	
	mov	a	c'31
	mul	a	r3
	mul	a	s1
	add	y	a
	
	mov	a	c'22
	mul	a	r2
	mul	a	s2
	add	y	a

	mov	a	c'13
	mul	a	r1
	mul	a	s3
	add	y	a
	
	mov	a	c'04
	mul	a	s4
	add	y	a

	mov	a	c'50
	mul	a	r5	
	add	y	a

	mov	a	c'41
	mul	a	r4
	mul	a	s1
	add	y	a

	mov	a	c'32
	mul	a	r3
	mul	a	s2
	add	y	a

	mov	a	c'23
	mul	a	r2
	mul	a	s3
	add	y	a

	mov	a	c'14
	mul	a	r1
	mul	a	s4
	add	y	a

	mov	a	c'05
	mul	a	s5
	add	y	a	; that's y
;
; 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
	_end