home *** CD-ROM | disk | FTP | other *** search
-
- ; *******************************************************
- ; * *
- ; * Turbo Pascal Runtime Library Version 7.0 *
- ; * Real ArcTan *
- ; * *
- ; * Copyright (C) 1989-1993 Norbert Juffa *
- ; * *
- ; *******************************************************
-
- TITLE FPATN
-
-
- CODE SEGMENT BYTE PUBLIC
-
- ASSUME CS:CODE
-
- ; Externals
-
- EXTRN RealAdd:NEAR,RealPoly:NEAR,RealTrunc:NEAR,RealFloat:NEAR
- EXTRN RealDivRev:NEAR,RealMulNChk2:NEAR,CmpAbsValue:NEAR
- EXTRN ShortMulRev:NEAR
-
- ; Publics
-
- PUBLIC RArcTan
-
- ;-------------------------------------------------------------------------------
- ; RArcTan computes the arctangent of a TURBO-Pascal six byte floating-point
- ; number. No overflow is possible, since function values range from -pi/2 to
- ; pi/2. The argument is reduced such that |z| <= 1/16 by an interval scheme
- ; based on the identity Arctan (x) = Arctan (a) + Arctan ((x-a) / (1 + x*a)).
- ; The arctangent of the reduced argument is then computed using a fast poly-
- ; nomial approximation. The result is finally adjusted depending on which
- ; reduction interval the argument was in. This polynomial approximation to the
- ; arctangent is used:
- ;
- ; p(z):=((-1.420288085938e-1*x^2+1.999981270101e-1)*x^2
- ; -3.333333321307e-1)*x^2*x+x
- ;
- ; This approximation has a theoretical maximum relative error of 2.98e-13.
- ; Maximum observed error when evaluated in REAL arithmetic is 1.11e-12.
- ;
- ; INPUT: DX:BX:AX argument
- ;
- ; OUTPUT: DX:BX:AX arctan of argument
- ;
- ; DESTROYS: AX,BX,CX,DX,SI,DI,Flags
- ;-------------------------------------------------------------------------------
-
- RArcTanExt PROC FAR
- $atn_x: RET ; return zero
- RArcTanExt ENDP
-
- ALIGN 4
-
- RArcTan PROC FAR
- CMP AL, 06Dh ; x very small ?
- JB $atn_x ; yes, return x
- PUSH BP ; save TURBO-Pascal frame pointer
- MOV CH, 80h ; load mask for sign bit
- AND CH, DH ; isolate sign bit
- XOR DH, CH ; absolute value of x
- PUSH CX ; save sign bit
- MOV CX, 81h ; load
- XOR SI, SI ; floating
- MOV DI, SI ; constant 1
- CALL CmpAbsValue ; x <= 1 ?
- PUSHF ; save division flag
- JB $no_division ; x <= 1, no division
- CALL RealDivRev ; compute x = 1/x
- $no_division:MOV DI, DX ; save
- MOV SI, BX ; x
- PUSH AX ;
- ADD AL, 3 ; x * 8
- MOV CH, 0FFh ; signal rounding desired
- CALL RealTrunc ; n = Round (x * 8)
- MOV BP, AX ; save n
- POP CX ; DI:SI:CX = x
- OR AX, AX ; n = 0 ?
- MOV DX, DI ; reload x
- MOV BX, SI ; into
- XCHG AX, CX ; DX:BX:AX
- JZ $atn_appr ; n=0, evaluate approximating polynomial
- XCHG AX, CX ; AX=n, CX = LSW of x
- PUSH AX ; save n
- PUSH CX ; save LSW of x
- CWD ; convert n to longint (n <= 8)
- CALL RealFloat ; convert n to REAL
- POP CX ; get LSW of x
- SUB AL, 3 ; n/8
- MOV BP, AX ; save
- MOV AH, DH ; AX = packed n
- XCHG AX, BP ; BP = packed n
- XOR DH, 80h ; -n/8
- PUSH DI ; save
- PUSH SI ; x
- PUSH CX ; save x
- CALL RealAdd ; x - n/8
- POP CX ; get
- POP SI ; x
- POP DI ; from stack
- PUSH DX ; save x - n/8
- PUSH BX ; on
- PUSH AX ; stack
- MOV AX, BP ; get packed floating point n
- XOR DX, DX ; unpack
- MOV BX, DX ; floating point
- XCHG AH, DH ; n
- CALL ShortMulRev ; x * n
- MOV CX, 81h ; load
- XOR SI, SI ; real constant
- MOV DI, SI ; 1.0
- CALL RealAdd ; compute 1 + x*n/8
- POP CX ; get
- POP SI ; back
- POP DI ; x - n/8
- CALL RealDivRev ; compute (x - n/8) / (1 + x*n/8)
- POP BP ; get integer n
- MOV CX, BP ; compute
- ADD CX, CX ; offset into
- ADD BP, CX ; constant field
- ADD BP, BP ; as 6*n
- ADD BP, OFFSET ATN_CST; address of constant B for interval n
- $atn_appr: MOV CX, 3 ; approximation uses three coefficients
- XOR SI, SI ; polynomial of type P(x²)*x+x
- MOV DI,OFFSET AT_COEFF; pointer to first coefficient
- CALL RealPoly ; compute arctan (z) = z + z * P(z)
- CMP BP, OFFSET ATN_CST+6; first or second interval ?
- JB $first_intv ; first interval, no correction needed
- JNE $normal_add ; not second interval, no special add
- MOV CX, 03365h ; add eps for second interval
- MOV SI, 07B6Eh ; for a pseudo
- MOV DI, 05561h ; multiple precision
- CALL RealAdd ; addition
- $normal_add: MOV CX, CS:[BP-6] ; load value
- MOV SI, CS:[BP-4] ; of constant B appropriate
- MOV DI, CS:[BP-2] ; for reduction interval
- CALL RealAdd ; add B to arctan (z)
- $first_intv: POPF ; get flag for division
- JB $no_add ; if no division then done
- MOV CX, 02181h ; load
- MOV SI, 0DAA2h ; real constant
- MOV DI, 0490Fh ; 0.5*pi
- XOR DH, 80h ; - arctan (z)
- CALL RealAdd ; compute 0.5*pi - arctan (z)
- $no_add: POP CX ; get sign mask
- POP BP ; restore TURBO-Pascal frame pointer
- OR DH, CH ; make sign bit
- RET ; done
- ATN_CST DB 07Dh,000h,000h,0D4h,0ADh,07Eh ; arctan (1/8) - eps
- DB 07Eh,064h,0C9h,0AFh,0DBh,07Ah ; arctan (2/8)
- DB 07Fh,027h,00Fh,0CAh,0B0h,037h ; arctan (3/8)
- DB 07Fh,00Eh,02Bh,038h,063h,06Dh ; arctan (4/8)
- DB 080h,0F8h,05Eh,05Dh,000h,00Fh ; arctan (5/8)
- DB 080h,035h,019h,07Dh,0BCh,024h ; arctan (6/8)
- DB 080h,0C2h,02Bh,03Eh,005h,038h ; arctan (7/8)
- DB 080h,021h,0A2h,0DAh,00Fh,049h ; arctan (8/8)
- AT_COEFF DB 07Eh, 070h,091h ; -1.420288085938e-1
- DB 07Eh,014h,01Bh,04Fh,0CCh,04Ch ; 1.999981270101e-1
- DB 07Fh,056h,0A0h,0AAh,0AAh,0AAh ; -3.333333321307e-1
- RArcTan ENDP
-
- ALIGN 4
-
- CODE ENDS
-
- END