
*	Maxon C++ Library:
*	Modul "math_doubtrans"
*	Jens Gelhar 19.03.93, 15.04.94

	xdef	_INIT_1_MatH_DoubTRans_o,_EXIT_1_MatH_DoubTRans_o
	xdef	_MathIeeeDoubTransBase,GetMathTrans
	xdef	_sin,_cos,_tan,_sinh,_cosh,_tanh,_asin,_acos,_atan
	xdef	_exp,_log,_log10,_sqrt,_pow,_atan2
	xdef	_modf,_frexp,_ldexp,_fmod

	xref	__OpenLib,_MathIeeeDoubBasBase,_ceil,_floor

Execbase	= 4

CloseLibrary	= -414
OpenLibrary	= -552

_INIT_1_MatH_DoubTRans_o:
	lea mathdoubtransname(pc),a1
	jsr __OpenLib
	move.l d0,_MathIeeeDoubTransBase
	rts

_EXIT_1_MatH_DoubTRans_o:
	move.l _MathIeeeDoubTransBase,d0
	beq.b ret
	move.l Execbase.w,a6
	move.l d0,a1
	jmp CloseLibrary(a6)

GetMathTrans:
	move.l _MathIeeeDoubTransBase,a6
	bne.b ret
	movem.l d0-d3/a0-a3,-(a7)
	bsr.b _INIT_1_MatH_DoubTRans_o
	movem.l (a7)+,d0-d3/a0-a3
	bra.b GetMathTrans
ret:	rts

_sin:	move.l a6,-(a7)
	movem.l 4+4(a7),d0/d1
	move.l _MathIeeeDoubTransBase,a6
	jsr -36(a6)
	move.l (a7)+,a6
	rts

_cos:	move.l a6,-(a7)
	movem.l 4+4(a7),d0/d1
	move.l _MathIeeeDoubTransBase,a6
	jsr -42(a6)
	move.l (a7)+,a6
	rts

_tan:	move.l a6,-(a7)
	movem.l 4+4(a7),d0/d1
	move.l _MathIeeeDoubTransBase,a6
	jsr -48(a6)
	move.l (a7)+,a6
	rts

_asin:	move.l a6,-(a7)
	movem.l 4+4(a7),d0/d1
	move.l _MathIeeeDoubTransBase,a6
	jsr -114(a6)
	move.l (a7)+,a6
	rts

_acos:	move.l a6,-(a7)
	movem.l 4+4(a7),d0/d1
	move.l _MathIeeeDoubTransBase,a6
	jsr -120(a6)
	move.l (a7)+,a6
	rts

_atan:	move.l a6,-(a7)
	movem.l 4+4(a7),d0/d1
	move.l _MathIeeeDoubTransBase,a6
	jsr -30(a6)
	move.l (a7)+,a6
	rts

_sinh:	move.l a6,-(a7)
	movem.l 4+4(a7),d0/d1
	move.l _MathIeeeDoubTransBase,a6
	jsr -60(a6)
	move.l (a7)+,a6
	rts

_cosh:	move.l a6,-(a7)
	movem.l 4+4(a7),d0/d1
	move.l _MathIeeeDoubTransBase,a6
	jsr -66(a6)
	move.l (a7)+,a6
	rts

_tanh:	move.l a6,-(a7)
	movem.l 4+4(a7),d0/d1
	move.l _MathIeeeDoubTransBase,a6
	jsr -72(a6)
	move.l (a7)+,a6
	rts

_exp:	move.l a6,-(a7)
	movem.l 4+4(a7),d0/d1
	move.l _MathIeeeDoubTransBase,a6
	jsr -78(a6)
	move.l (a7)+,a6
	rts

_log:	move.l a6,-(a7)
	movem.l 4+4(a7),d0/d1
	move.l _MathIeeeDoubTransBase,a6
	jsr -84(a6)
	move.l (a7)+,a6
	rts

_log10:	move.l a6,-(a7)
	movem.l 4+4(a7),d0/d1
	move.l _MathIeeeDoubTransBase,a6
	jsr -126(a6)
	move.l (a7)+,a6
	rts

_sqrt:	move.l a6,-(a7)
	movem.l 4+4(a7),d0/d1
	move.l _MathIeeeDoubTransBase,a6
	jsr -96(a6)
	move.l (a7)+,a6
	rts

_pow:	movem.l d2/d3/a6,-(a7)
	movem.l 3*4+4(a7),d0/d1;         = x
	move.l _MathIeeeDoubTransBase,a6
	jsr -84(a6)                     ; ln x
	movem.l 3*4+12(a7),d2/d3;        = y
	move.l _MathIeeeDoubBasBase,a6
	jsr -78(a6)                     ; y * ln x
	move.l _MathIeeeDoubTransBase,a6
	jsr -78(a6)                     ; exp(...)
	movem.l (a7)+,d2-d3/a6
	rts


pi_1	= $400921FB
pi_2	= $54442D0F

_atan2:	movem.l d2-d3/a0-a2/a6,-(a7)
	movem.l 6*4+4(a7),d0-d3
	tst.l d2
	bne.b at0
	; Division durch Null umgehen:
	move.l #$3FF921FB,d0            ; = pi/2
	move.l #$54442CF9,d1
	bra.b at1
at0:	move.l _MathIeeeDoubBasBase,a6
	jsr -84(a6)                     ; DIV
	move.l _MathIeeeDoubTransBase,a6
	jsr -30(a6)
at1:	tst.l 6*4+4(a7)         ; x < 0 ?
	bpl.b at3               ; dann Korrektur-Rechnung
	move.l #pi_1,d2
	move.l #pi_2,d3
	tst.l d0
	bmi.b at2
	; x = x'-pi
	move.l _MathIeeeDoubBasBase,a6
	jsr -72(a6)
	bra.b at3
at2:	; x = x'+pi
	move.l _MathIeeeDoubBasBase,a6
	jsr -84(a6)
at3:	movem.l (a7)+,d2-d3/a0-a2/a6
	rts


_modf:	; Stack-Parameter: double, double*
	movem.l d2-d4/a0,-(a7)
	movem.l 4*4+4(a7),d0/d1/a0
	move.l d0,d4
	beq.b modf_small
	bclr #31,d4
	swap d4
	lsr.w #4,d4             ; Exponent
	sub.w #$3FF,d4
	blo.b modf_small
	cmp.w #52,d4
	bhs.b modf_large
	move.l #$FFF00000,d2
	moveq  #$0,d3
	tst.w d4
	beq.b modf_end
modf_shift:	asr.l #1,d2
	roxr.l #1,d3
	subq.w #1,d4
	bne.b modf_shift
modf_end:	and.l d2,d0             ; d0|d1 ist Ganzzahl-Anteil
	and.l d3,d1
	cmp.l #0,a0
	beq.b modf_1
	movem.l d0/d1,(a0)
modf_1:	move.l d0,d2
	move.l d1,d3
	movem.l 4*4+4(a7),d0/d1
	move.l _MathIeeeDoubBasBase,a6
	jsr -72(a6)           ; sub
modf_ret:	movem.l (a7)+,d2-d4/a0
	rts
modf_small:	; Betrag des Parameters < 1:
	cmp.l #0,a0
	beq.b modf_ret
	clr.l (a0)+
	clr.l (a0)+
	bra.b modf_ret
modf_large:	; Betrag größer als 2^52:
	cmp.l #0,a0
	beq.b modf_2
	movem.l d0/d1,(a0)
modf_2:	moveq #0,d0
	moveq #0,d1
	bra.b modf_ret


_frexp:	; Parameter: double, int*
	movem.l d2/a0,-(a7)
	movem.l 2*4+4(a7),d0/d1/a0
	move.l d0,d2
	beq.b frexp_zero
	bclr #31,d2
	swap d2
	lsr.w #4,d2
	sub.w #$3FE,d2
	ext.l d2
	move.l d2,(a0)
	and.l #$000FFFFF,d0
	or.l  #$3FE00000,d0
frexp_ret:	movem.l (a7)+,d2/a0
	rts
frexp_zero:	clr.l (a0)
	bra.b frexp_ret


_ldexp:	; Parameter: double, int
	move.w 4(a7),d0
	beq.b ldexp_zero
	move.l 12(a7),d1
	bclr #15,d0
	lsr.w #4,d0
	ext.l d0
	add.l d1,d0
	cmp.l #$7FF,d0
	bhi.b ldexp_overflow
	swap d0
	asl.l #4,d0
	and.l #$800FFFFF,4(a7)
	or.l 4(a7),d0
	move.l 8(a7),d1
	rts
ldexp_overflow:	tst.l d1
	bmi.b ldexp_zero
	move.l #$7FFFFFFF,d0
	moveq #-1,d1
	rts
ldexp_zero:	moveq #0,d0
	moveq #0,d1
	rts

_fmod:	; x%y = y*frac(x/y)
	movem.l d2-d3/a0/a6,-(a7)
	movem.l 4*4+4(a7),d0-d3
	move.l _MathIeeeDoubBasBase,a6
	jsr -84(a6)                     ; x/y
	subq.l #8,a7
	move.l a7,-(a7)
	movem.l d0/d1,-(a7)
	bsr _modf                       ; Nachkommateil berechnen
	lea 12+8(a7),a7
	movem.l 4*4+12(a7),d2/d3
	move.l _MathIeeeDoubBasBase,a6
	jsr -78(a6)                     ; mit y multiplizieren
	movem.l (a7)+,d2-d3/a0/a6
	rts

mathdoubtransname:	dc.b 'mathieeedoubtrans.library',0

	SECTION D,DATA

_MathIeeeDoubTransBase:	blk.l 1

	end

