movl	%edx,%esi		andl	$0x38000000,%esi		addl	%esi,%edx		addb	%bl,%bl			cmc				notl	%eax			adcl	$0,%eax			notl	%edx			adcl	$0,%edx			movl	$0x3fff,%ecx		adcl	$0,%ecx			orl	$0x80000000,%edx		xchgl	(%edi),%eax	#swap arguments	xchgl	4(%edi),%edx		xchgl	8(%edi),%ecx		movl	%ecx,%ebx	#registers = x-C	addw	$0xc006,%bx		jle	log52			xchgl	%ecx,%ebx		shll	%cl,%edx		shrl	%cl,%edx		xchgl	%ecx,%ebx		jmp	log54		#        subttl   arctangent of x#        page# take arctangent of x# input:        x in cx:edx:eax and ST# output:       atan(x) in cx:edx:eax	.globl	epatan	.balign	16atncon:				#coefficients for polynomial	.word	0,0x45d1,0x5d17,0xd174,0x1745,0xf45d	.word	0,0xe38e,0x8e38,0x38e3,0xe38e,0x0e38	.word	0,0xdb6e,0x6db6,0xb6db,0xdb6d,0xedb6	.word	0,0x999a,0x9999,0x9999,0x9999,0x1999	.word	0,0x5555,0x5555,0x5555,0x5555,0xd555	.word	0,0x0000,0x0000,0x0000,0x0000,0x8000	.word	-1	.balign	16attab:				#interpolation table for n/32	.word	0x0000,0x0000,0x0000,0x0000,0x0000,1 # atan(0.000000)	.word	0x2542,0x4bb1,0xaddd,0xffea,0x3ff9,0 # atan(0.031250)	.word	0x4e37,0x67ef,0xddb9,0xffaa,0x3ffa,0 # atan(0.062500)	.word	0x7460,0x1788,0xc130,0xbf70,0x3ffb,0 # atan(0.093750)	.word	0x6e33,0x617b,0xd4d5,0xfead,0x3ffb,0 # atan(0.125000)	.word	0x62c4,0x3313,0x7746,0x9eb7,0x3ffc,0 # atan(0.156250)	.word	0x1135,0x72d8,0xda5e,0xbdcb,0x3ffc,0 # atan(0.187500)	.word	0x1023,0x9305,0xba94,0xdc86,0x3ffc,0 # atan(0.218750)	.word	0xeb16,0x6406,0xafc9,0xfadb,0x3ffc,0 # atan(0.250000)	.word	0xc130,0x5f8b,0xad18,0x8c5f,0x3ffd,0 # atan(0.281250)	.word	0x5e6a,0x3f5e,0xb9b8,0x9b13,0x3ffd,0 # atan(0.312500)	.word	0x4eda,0x8e6a,0x6cca,0xa985,0x3ffd,0 # atan(0.343750)	.word	0x8473,0x26f7,0xca0f,0xb7b0,0x3ffd,0 # atan(0.375000)	.word	0x2b6d,0x50d9,0x69ca,0xc592,0x3ffd,0 # atan(0.406250)	.word	0xe5b6,0x611f,0x761e,0xd327,0x3ffd,0 # atan(0.437500)	.word	0x7c68,0x764f,0xa64a,0xe06d,0x3ffd,0 # atan(0.468750)	.word	0x7b46,0x0dda,0x382b,0xed63,0x3ffd,0 # atan(0.500000)	.word	0xbe5c,0xa0a0,0xe85a,0xfa06,0x3ffd,0 # atan(0.531250)	.word	0x7e2a,0xd986,0xf4a6,0x832b,0x3ffe,0 # atan(0.562500)	.word	0x47b5,0xde95,0xecdf,0x892a,0x3ffe,0 # atan(0.593750)	.word	0x9f9c,0xf7f5,0x5d5e,0x8f00,0x3ffe,0 # atan(0.625000)	.word	0x86f6,0x8471,0x72c9,0x94ac,0x3ffe,0 # atan(0.656250)	.word	0xda20,0x71bd,0x80e6,0x9a2f,0x3ffe,0 # atan(0.687500)	.word	0xa1ed,0xf4b7,0xfdc4,0x9f89,0x3ffe,0 # atan(0.718750)	.word	0x0924,0x34f7,0x7d19,0xa4bc,0x3ffe,0 # atan(0.750000)	.word	0xf5c8,0x4830,0xabdc,0xa9c7,0x3ffe,0 # atan(0.781250)	.word	0xc080,0xb4d8,0x4c38,0xaeac,0x3ffe,0 # atan(0.812500)	.word	0x3691,0x1f04,0x31c9,0xb36b,0x3ffe,0 # atan(0.843750)	.word	0x9e74,0xc231,0x3e2b,0xb805,0x3ffe,0 # atan(0.875000)	.word	0xf281,0xe98a,0x5dea,0xbc7b,0x3ffe,0 # atan(0.906250)	.word	0x6640,0xac52,0x85b8,0xc0ce,0x3ffe,0 # atan(0.937500)	.word	0xbd54,0xbf8f,0xaffa,0xc4ff,0x3ffe,0 # atan(0.968750)atnw8:	movl	8(%edi),%ecx	#return 3pi/4	andl	$0x80000000,%ecx		movw	$0x4000,%cx		movl	$0x96cbe3f9,%edx		movl	$0x990e91a7,%eax		orb	$precis,h_stat		call	exprop			ret			atnzer:	movb	11(%edi),%ch	#return signed zero	movb	$1,%cl			shll	$16,%ecx		xorl	%edx,%edx		xorl	%eax,%eax		ret			atnw11:	orl	%ecx,%ecx	# INF / INF	js	atnw8			orb	$precis,h_stat		call	exprop			jmp	atn90		atnw3:	cmpw	$0x7fff,%cx	# INF / y	jz	atnw11		atnw9:	orb	$precis,h_stat		call	exprop		atn82:	movl	$0x3fff,%ecx	#INF returns +- pi/2	jmp	atn91		atnw5:	orl	%ecx,%ecx	# y / INF or 0 / x	jns	atnzer			orb	$precis,h_stat		call	exprop		atn84:	movl	$0x4000,%ecx	#here return +-pi	jmp	atn91		atnw:					call	retnan		#proper return for illegal arguments	testb	$1,10(%edi)	#weed out x = 0	jnz	atnw5			btl	$16,%ecx	#weed out y = 0	jc	atnw9			cmpw	$0x7fff,8(%edi)	#weed out x = INF	jz	atnw3			cmpw	$0x7fff,%cx	#weed out y = INF	jz	atnw5			call	tnormal		#normalize arguments	xchgl	(%edi),%eax		xchgl	4(%edi),%edx		xchgl	8(%edi),%ecx		call	tnormal			xchgl	(%edi),%eax		xchgl	4(%edi),%edx		xchgl	8(%edi),%ecx		jmp	atn2		atn90:	orl	%ecx,%ecx		js	atnw8			movl	$0x3ffe,%ecx	#atan(1) return +-pi/4atn91:	movl	$0xc90fdaa2,%edx		movl	$0x2168c235,%eax		jmp	atn23		atn89:	movb	s_X+5(%ebp),%bl		testb	$1,%bl			jz	atn82			orb	%bl,%bl			js	atn84			call	epdiv1			jmp	atn23		epatan:					shldl	$16,%ecx,%ebx	#save sign bits	movb	11(%edi),%bl		movw	%bx,s_X+4(%ebp)		shldl	$16,%ecx,%ebx	#weed out special arguments	orb	10(%edi),%bl		jnz	atnw		atn2:	#PREEXS  precis          #set precision exception	orb	$precis,h_stat		testb	$precis,h_ctrl		jnz	LL8			orw	$0x8080,h_stat		movl	$exret,-4(%ebp)	LL8:					cmpw	8(%edi),%cx	#compare x and y	jl	atn5			jnz	atn7			cmpl	4(%edi),%edx		jc	atn5			jnz	atn7			cmpl	(%edi),%eax		jz	atn90			jc	atn5						atn7:	xchgl	(%edi),%eax	#switch arguments, set flag	xchgl	4(%edi),%edx		xchgl	8(%edi),%ecx		incb	s_X+5(%ebp)	atn5:	movl	%ecx,%ebx	#underflow here means either	addw	$tBIAS,%bx	#   atn(INF) = +-pi/2 or	subw	8(%edi),%bx	#   atn(0) = 0 or pi	cmpw	$0x3f40,%bx		jl	atn89		atn6:	call	epdiv1		#divide x / y	andl	$0x7fffffff,%ecx		#PUTEP   	movl	%eax,(%edi)		movl	%edx,4(%edi)		movl	%ecx,8(%edi)		xorl	%ebx,%ebx	#get interpolation interval, 0-31	subw	$tBIAS-6,%cx		cmpw	$9,%cx			jc	atn4			movw	$0,%cx		atn4:	shldl	%cl,%edx,%ebx		pushl	%ebx			shll	%cl,%edx	#calculate x - C = d	shrl	%cl,%edx		movw	8(%edi),%cx		js	atn55			jnz	atn51			orl	%eax,%eax		jnz	atn51			xorl	%ecx,%ecx		jmp	atn54		atn51:	decw	%cx			addl	%eax,%eax		adcl	%edx,%edx		jns	atn51		atn55:	pushl	%eax		#push d	pushl	%edx			pushl	%ecx							rorl	$5,%ebx		#calculate x * C	movl	4(%edi),%eax		mull	%ebx			movl	%eax,%ecx		xchgl	%ebx,%edx		movl	(%edi),%eax		mull	%edx			addl	%edx,%ecx		adcl	$0,%ebx			movl	%ebx,%edx	#get x * C into right registers	xchgl	%ecx,%eax		xorl	%ebx,%ebx		movl	$0x3fff,%ecx		subw	8(%edi),%cx	atn11:	shrl	$1,%edx		#normalize x * C to exponent 3fff	rcrl	$1,%eax			loop	atn11			orl	$0x80000000,%edx	#calculate 1 + x*C	movw	$0x3fff,%cx		#PUTEP   	movl	%eax,(%edi)		movl	%edx,4(%edi)		movl	%ecx,8(%edi)		popl	%ecx		#calculate d/(1 + x*C)	popl	%edx			popl	%eax			call	epdiv1			addb	%bl,%bl			adcl	$0,%eax			adcl	$0,%edx			#PUTEP   	movl	%eax,(%edi)		movl	%edx,4(%edi)		movl	%ecx,8(%edi)		call	square		#polynomial uses x^2	movl	$atncon,%ebx	#do the polynomial	call	sigmax			addw	$0x3ffe,%cx		call	epmul1		#multiply with xatn54:	popl	%ebx		#then add atan(C)       	shll	$2,%ebx			lea	attab(%ebx,%ebx,2),%ebx		movl	%cs:(%ebx),%esi		movl	%esi,(%edi)		movl	%cs:4(%ebx),%esi		movl	%esi,4(%edi)		movl	%cs:8(%ebx),%esi		movl	%esi,8(%edi)		call	epadd1			movb	s_X+5(%ebp),%bl	#adjust for actual quadrant	cmpb	$0x01,%bl		jz	atn23			cmpb	$0x80,%bl		jz	atn25			orl	$0x80000000,%ecx	atn25:	movw	$0x4000,8(%edi)		cmpb	$0x81,%bl		jz	atn26			movw	$0x3fff,8(%edi)	atn26:	movl	$0xc90fdaa2,4(%edi)		movl	$0x2168c235,0(%edi)		call	epadd1		atn23:	roll	$8,%ecx		#put in sign	movb	s_X+4(%ebp),%cl		rorl	$8,%ecx			ret			# routine to compute a Taylor polynomial# for speed and accuracy, we perform fixed-point scaled calculations# in    ebx = pointer to coefficient table#       ecx:edx:eax = EP argument# out   ecx:edx:eax = EP resultshf22:	movl	%edx,%eax	#big shift here	xorl	%edx,%edx		subw	$32,%cx			cmpw	$32,%cx			jc	shf9			xorl	%eax,%eax		jmp	shf19		sigmax:					pushl	%edi		#save registers	movl	%ecx,%esi		subw	$0x3ffe,%cx	#get the scaling factor	negl	%ecx			cmpw	$32,%cx			jnc	shf22		shf9:	shrdl	%cl,%edx,%eax	#shift x right to scale	shrl	%cl,%edx	shf19:	movl	%eax,s_A(%ebp)	#store x 	movl	%edx,s_A+4(%ebp)		movw	%cs:2(%ebx),%ax	#initial value of sum = first coefficient	movl	%cs:4(%ebx),%ecx		movl	%cs:8(%ebx),%edx		addl	$12,%ebx	sig5:	pushl	%ebx		#keep constant pointer on the stack	movl	%edx,%edi		mulw	s_A+6(%ebp)	#multiply with x01	movw	%dx,%si		#   put result in ebx:ecx:si	movl	%edi,%eax		mull	s_A+4(%ebp)		movl	%edx,%ebx		xchgl	%ecx,%eax		orl	%edi,%edi		jns	sig42			subl	s_A(%ebp),%ecx		sbbl	s_A+4(%ebp),%ebx	sig42:	pushl	%eax			mull	s_A+4(%ebp)		shrl	$16,%eax		addw	%ax,%si			adcl	%edx,%ecx		adcl	$0,%ebx			popl	%eax		#multiply with x23	shrl	$16,%eax		mulw	s_A+2(%ebp)		addw	%dx,%si			adcl	$0,%ecx			adcl	$0,%ebx			movl	%edi,%eax		mull	s_A(%ebp)		shrl	$16,%eax		addw	%si,%ax			adcl	%edx,%ecx		adcl	$0,%ebx			movl	%ebx,%edx	#sum now in edx:ecx:ax	popl	%ebx		#address of constant table in bx	orl	%esi,%esi	#we may have to complement the sum	jns	sig31			notl	%edx			notl	%ecx			negw	%ax			sbbl	$-1,%ecx		sbbl	$-1,%edx	sig31:	addw	%cs:2(%ebx),%ax	#add next coefficient to sum	adcl	%cs:4(%ebx),%ecx		adcl	%cs:8(%ebx),%edx		addl	$12,%ebx	#proceed to next element	cmpb	$-1,(%ebx)		jnz	sig5			movb	%ah,%bl		#restore registers	movl	%ecx,%eax		popl	%edi			movl	$1,%ecx		#normalize number	orl	%edx,%edx		js	sig18			decl	%ecx			addb	%bl,%bl			adcl	%eax,%eax		adcl	%edx,%edx	sig18:	ret			#routine to raise number to second power# in    x in ecx:edx:eax# out   x^2 in samesqu80:	ret			#return unchangedsquare:					orw	%cx,%cx		#may be called with zero	jz	squ80			pushl	%edi		#preserve registers	pushl	%ebx			pushl	%ebp			pushl	%ecx			movl	%edx,%esi	#move mantissa	movl	%eax,%edi		shrl	$16,%eax	#x23 * x23	mulw	%ax			movb	%dh,%bl			movl	%esi,%eax	#x01 * x01 -> bp:cx:bl	mull	%eax			movl	%eax,%ecx		movl	%edx,%ebp		movl	%esi,%eax	#double x01 * x23	mull	%edi			shrl	$16,%eax		addb	%al,%bl			adcl	%edx,%ecx		adcl	$0,%ebp			addb	%al,%bl			adcl	%edx,%ecx		adcl	$0,%ebp			movl	%ecx,%eax	#move mantissa to proper place	movl	%ebp,%edx		popl	%ecx		#double exponent	addl	%ecx,%ecx		subw	$tBIAS-1,%cx		addb	%bl,%bl		#round	adcl	$0,%eax			adcl	$0,%edx			popl	%ebp		#restore registers	popl	%ebx			popl	%edi			ret			# routine to reduce the trigonometric argument to between 0 and pi/4# returns octant number (0-7) in cl	.globl	reductred80:	orb	$4,h_stat+1	#argument too big	addl	$4,%esp			ret			red85:	call	getqn		#here INF	orb	$invop,h_stat		call	exproc		red99:	xorb	%bl,%bl			ret			red96:	andl	$0x8000ffff,%ecx		jmp	red99		redw:	call	chkarg		#weed out illegals	jc	red99			btl	$16,%ecx	#weed out INF, zero	jc	red99			cmpw	$1,%cx			jz	red96			orw	%cx,%cx			jnz	red85			btsl	$18,%ecx		cmpb	$7,rmcode(%ebp)		jz	red99			orb	$underf,h_stat		call	exproc			jmp	red99		reduct:					andb	$0x0fb,h_stat+1	#clear reduction bit	testl	$0x00030000,%ecx	#jump if special value	jnz	redw			cmpw	$tBIAS-1,%cx	#small x, do nothing	jc	red99			subw	$tBIAS-2,%cx	#get shift count -> cx	cmpw	$64,%cx		#   reduce only up to 63	jg	red80			pushl	%edi		#preserve register	movl	$0x2168c234,%edi		movl	$0xc90fdaa2,%esi		xorl	%ebx,%ebx	#divide by subtractionred61:	subb	$0x0c0,%bh	#   remainder -> edx:eax:bh	sbbl	%edi,%eax		sbbl	%esi,%edx		sbbb	$0,%ch			jnc	red63			xorb	%ch,%ch			addb	$0x0c0,%bh		adcl	%edi,%eax		adcl	%esi,%edx	red63:	cmc				adcb	%bl,%bl			decb	%cl			jle	red65			addb	%bh,%bh			adcl	%eax,%eax		adcl	%edx,%edx		adcb	$0,%ch			jmp	red61		red65:	andb	$7,%bl		#octant number -> bh	testb	$1,%bl		#if octant 1 or 3 return pi/4 - x	jz	red41			negb	%bh			addb	$0x0c0,%bh		notl	%eax			adcl	%edi,%eax		notl	%edx			adcl	%esi,%edx	red41:	movl	$0x3ffe,%ecx	#exponent of pi/4	popl	%edi		red37:	orl	%edx,%edx	#normalize so bit 31 = 1	js	red19			jnz	red25			xchgl	%eax,%edx		xchgb	%bh,%al			shll	$16,%eax	#fix 10/28/97#	shll	$16,%ebx		subl	$32,%ecx		jmp	red37		red25:	decl	%ecx			addb	%bh,%bh			adcl	%eax,%eax		adcl	%edx,%edx		jns	red25		red19:	ret			#__lib   endp#emuseg  ends#        end