(Comp.sources.hp48)
Item: 140 by kcooke@milton.u.washington.edu [Ken Cooke]
Subj: Machine Code FFT source code
Date: 14 Jul 1992
I'm glad to hear that someone liked the machine code FFT programs.
In response to email requests, here is the SASM source code. It can
be assembled with the HP tools. I didn't include FFT, because it is a
subset of IFFT (I commented the lines to remove).
I apologize for the shortage of higher-level comments (most of this is
still on paper). Also, I traded some readability to avoid subroutine calls
(hardware stack levels are very precious) and keep variables in registers.
A few of the ROM calls are "unsupported entry points", which may need
fixing in the future. With rumors of new 48's in the works, ones with
machine code Eqn Writers, re-written plotter user interfaces, and maybe
even faster Saturn processors...well, lets just say that it may be wise to
stick to entries listed in ENTRIES.A.
I hope this is of use to someone interested in ML math software. While I
like playing games on my 48 as much as the next guy, it would be nice to
see more math tools using the capabilities of sys-RPL and ML. (I do,
occasionally, use my 48 as a calculating device :) Unfortunately,
there is no documentation on using the low-level floating point routines
(MULTF, RADDF, etc) so a little hacking with Voyager or MLDL is required.
I would love to hear from other math-hackers out there...
Ken Cooke
N7VFE
kcooke@u.washington.edu
----------
* ANGLE.S WORKS FOR ANY SIZE VECTOR OR ARRAY
* Copyright 1992 Ken Cooke
* Takes a VECTOR/ARRAY in STK1. Returns a real VECTOR/ARRAY of
* same size, that contains the angle of each complex element,
* according to the current trig mode.
TITLE ANGLE
ASSEMBLE
NIBASC /HPHP48-A/
RPL
::
CK1NoBlame
CK&DISPATCH1
FOUR
::
* IF REAL ARRAY, EXIT UNCHANGED
DUP TYPECARRY? NOT?SEMI ( ORIG )
* MALLOC ANSWER ARRAY (SAME SIZE AS ORIG BUT REAL )
DUP MDIMS ITE TWO{}N ONE{}N %0 MAKEARRY ( ORIG ANS )
CODE
****************** UNLISTED ENTRY POINTS *****************
ATTN?Lp = #0CA88
ARG15 = #2B6D7
SetTrigMode = #2AEF6
************************ MCODE ANGLE ***********************
* ENTRY:
* STK2: ORIGINAL ARRY1 (COMPLEX)
* STK1: ANSWER ARRY (REAL)
* EXIT: STK1: ANSWER ARRAY (REAL)
NIBHEX 823 ;CLEAR SB AND XM IN HST
D1=D1+ 5 ;DROP STK1
D=D+1 A ;RETURN MEM
GOSBVL =SAVPTR ;SAVE REGS WITH ONLY orig ON STACK
A=DAT1 A ;A->ORIG
D1=D1- 5 ;POINT BACK TO ANS
C=DAT1 A ;C->ANS
R3=C.F A ;SAVE ->ANS IN R3
D0=C
D0=D0+ 15
D0=D0+ 10 ;D0->ANS DATA (OR DIM2 IF ARRAY)
D1=A ;D1->ORIG
D1=D1+ 15 ;D1->#DIMS
A=DAT1 A ;A=#DIMS
D1=D1+ 5 ;D1->DIM1
C=DAT1 A ;C=DIM1
D1=D1+ 5 ;D1->ORIG DATA (OR DIM2 IF ARRAY)
A=A-1 A ;CHECK FOR VECTOR
?A=0 A
GOYES ISVECT
***IF ARRAY, #ELEMENTS=DIM1*DIM2. ALSO, SKIP OVER DIM2 FIELD
A=DAT1 A ;A=DIM2
D1=D1+ 5 ;D1->ORIG DATA
D0=D0+ 5 ;D0->ANS DATA
SETHEX
GOSBVL =MUL# ;B=A*C
C=B A ;C=DIM1*DIM2
ISVECT C=C-1 A ;SUB 1 SO CAN USE BORROW FOR LOOPING
R4=C.F A ;SAVE #ELEMENTS-1 IN R4
ABSLOOP GOSBVL =ATTN?Lp ;ABORT IF ATTN (ON) WAS PRESSED
A=DAT1 W ;A=RE
D1=D1+ 16 ;D1->IM
CD0EX
RSTK=C ;PUSH D0
C=DAT1 W ;C=IM
D1=D1+ 16 ;D1->NEXT RE
SETDEC
GOSBVL =SPLTAC ;CONVERT BOTH TO LONGS
GOSBVL =SetTrigMode ;SET ST ACCORDING TO TRIG MODE FLAGS
GOSBVL =ARG15 ;AB=ANGLE (TRASHES A,B,C,D,D0,R0,R1)
GOSBVL =PACK ;A=ABSVAL
C=RSTK
D0=C ;POP D0
DAT0=A W ;STORE ABSVAL IN ANSWER ARRAY
D0=D0+ 16 ;POINT AT NEXT ANSWER
* CHECK FOR DONE
C=R4.F A ;GET COUNTER
SETHEX
C=C-1 A
R4=C.F A ;DEC COUNTER
GONC ABSLOOP ;LOOP UNTIL BORROW
***NOW ANSWER CONTAINS ABS VALUES OF ORIG
A=R3.F A
GOSBVL =GPOverWrALp ;GETPTR, REPLACE TOS WITH ANS, LOOP
ENDCODE
;
;
END_SRC
BEGIN_SRC ifft.s
******************************************************************
* MACHINE CODE IFFT by Ken Cooke
* Copyright 1992
******************************************************************
* RADIX-2 DIT IFFT
* Input: real or complex vector on STK1, with size=2^n
* Ouput: complex vector containing IFFT
* Note: uses the definition that FFT uses exp(-xxx) and IFFT uses exp(+xxx)
TITLE IFFT
ASSEMBLE
NIBASC /HPHP48-A/
RPL
::
CK1NoBlame
CK&DISPATCH1
FOUR
::
DUP MDIMS ( ARRY #N flag )
* CHECK FOR MATRIX (NcaseDIMERR = 37DF6 )
NOT NcaseDIMERR ( ARRY #N )
* CHECK FOR SIZE=2^N
DUP ( ARRY #N #N )
CODE
****************** CHECK THAT SIZE=POWER OF TWO **************
*CHKSIZE2 REPLACES BINT (ON TOS) WITH TRUE IF BINT IS A
* POWER OF 2, FALSE OTHERWISE.
*NOTE: WORKS UP TO 2^19, AND RETURNS TRUE FOR BINT=0
WrTLoop EQU #03B1A
WrFLoop EQU #03B06
C=DAT1 A ;C->BINT
CD1EX ;SAVE D1 IN C.A
D1=D1+ 5 ;D1->DATA
A=DAT1 A ;A.A=BINT
D1=C ;RESTORE D1
SETHEX
C=A A
A=-A-1 A ;A=ONES COMP
C=-C A ;C=TWOS COMP
A=A!C A ;NO ZEROS IF BINT=2^N
A=A+1 A ;WILL SET CARRY IF BINT=2^N
GOC DoTrue
GOVLNG WrFLoop ;OVERWRITE BINT WITH FALSE AND CONT
DoTrue GOVLNG WrTLoop ;OVERWRITE BINT WITH TRUE AND CONT
ENDCODE
* ( ARRY #N FLAG )
NcaseDIMERR ( IF FALSE THEN "INVALID DIMENSION" ERROR )
* IF SIZE=1, RETURN WITH ARRY ON STACK
DUP#1= ITE DROP ( ARRY #N )
::
* CREATE TEMP STORAGE (MUST BE 172 NIBS BIGGER THAN UNPACKED A FOR VARS)
FORTYTWO OVER #* 172 #+ MAKEHEX$ SWAP ( ARRY HEX$ #N )
ONE{}N C%0 MAKEARRY ( ARRY HEX$ ANSARRY )
* DO THE FFT
CODE
****************** UNLISTED ENTRY POINTS *****************
PUTAB1 = #2C066
GETCD1 = #2C017
aSIN15 = #2B6E0
ATTN?Lp = #0CA88
longpi = #2A45D
varnibs = 172
************************ MCODE FFT ***********************
* ENTRY:
* STK3: (orig) real/complex data array
* STK2: (temp) hex$ (to hold vars and unpacked reals during FFT)
* STK1: (ans) complex answer array
* EXIT: STK1: complex FFT array (ans)
NIBHEX 823 ;CLEAR SB AND XM IN HST
D1=D1+ 10 ;DROP temp AND ans
D=D+1 A ;RETURN MEM
D=D+1 A
GOSBVL =SAVPTR ;SAVE REGS WITH ONLY orig ON STACK
* SAVE NEEDED POINTERS IN SCRATCH REGS
D1=D1- 5 ;POINT BACK TO temp
A=DAT1 A ;A->temp
A=A+CON A,10 ;SKIP HEADER TO vars
R4=A.F A ;R4->vars
P= 0 ;*****REMEMBER:ALWAYS HAVE P=0 BEFORE LC(5)*****
LC(5) varnibs ;
SETHEX
A=A+C A ;SKIP OVER VARS
R1=A.F A ;R1->data
D1=D1+ 5 ;POINT TO orig
A=DAT1 A ;A->orig
D0=A
D0=D0+ 10 ;D0->OBJ PROLOG
* SET FLAG IF REAL, FOR USE WHEN UNPACKING
* USE D.15 INSTEAD OF ST TO AVOID TROUBLE
D=0 S ;USE D.15 AS FLAG (0=COMPLEX,1=REAL)
A=DAT0 A ;A=OB PROLOG
LC(5) #02977 ;C=COMPLEX PROLOG
?A=C A
GOYES CMPLEX
D=D+1 S ;IF REAL, D.15=1
CMPLEX D0=D0+ 10 ;D0->SIZE
A=DAT0 A ;A.A=SIZE
R2=A.F A ;STORE SIZE IN R2
D0=D0+ 5 ;D0->DATA
CD0EX
R0=C.F A ;R0->origdata
*FIND LOG2 OF SIZE
C=0 S ;C.S WILL HOLD SIZE2=LOGbase2(SIZE)
LOG2LP ASRB.F A
C=C+1 S
?ABIT=0 0
GOYES LOG2LP
R2=C.F S ;STORE SIZE2 IN R2.15
*************************** BITREV AND UNPACK ****************************
*NOW, CONVERT REALS IN orig INTO EXT.REALS IN data
* AT THE SAME TIME, REARRANGE INTO BIT-REVERSED ORDER
A=R0.F A
D0=A ;D0->origdata
D=0 A ;D HOLDS INDEX INTO orig (START AT 0)
BITREV CPEX 15 ;PUT BIT-COUNTER IN P
P=P-1 ;SUB ONE SO CAN USE BORROW FOR BITLOOP
C=D A ;C GETS INDEX INTO orig
B=0 A ;B WILL BE BIT-REV INDEX
*BIT REVERSAL
BITLOOP B=B+B A ;SHIFT B LEFT ONE BIT
?CBIT=0 0 ;TEST LSB
GOYES NOCARRY ;SKIP IF ZERO
B=B+1 A ;ELSE FEED INTO B
NOCARRY CSRB.F A ;SHIFT C RIGHT ONE BIT
P=P-1 ;DEC BIT-COUNTER (BORROW WHEN DONE)
GONC BITLOOP ;LOOP IF MORE
*MULT B (BITREV INDEX) BY 21 AND ADD TO ->data
A=B A ;COULD USE =MUL# BUT THIS IS FASTER/SHORTER
B=B+B A
B=B+B A ;B=4B
A=A+B A ;A=5B
B=B+B A
B=B+B A ;B=16B
B=A+B A ;B=21B
B=B+B A ;B=42B
A=R1.F A ;A->data
A=A+B A ;
D1=A ;D1 POINTS INTO data AT BITREV POSITION
A=DAT0 W ;GET RE(Ai)
D0=D0+ 16 ;POINT AT NEXT REAL
GOSBVL =SPLITA ;NO NEED TO SETDEC (DONE IN SPLITA)
GOSBVL =PUTAB1 ;PUT LONG AT DAT1, D1+=21
* IF INPUT ARRAY WAS REAL, USE IMAG=0
?D=0 S ;IF COMPLEX
GOYES GETIM ; GET IMAG PART
A=0 W ;ELSE USE IMAG=0
B=0 W
GONC NOGETIM ;GO ALWAYS
GETIM A=DAT0 W ;GET IM(Ai)
D0=D0+ 16
GOSBVL =SPLITA
NOGETIM GOSBVL =PUTAB1 ;STORE IMAG PART
SETHEX
D=D+1 A ;INC INDEX
C=R2 ;LOAD SIZE AND SIZE2 INTO C
?D=C A ;COMPARE INDEX TO SIZE
GOYES BITDONE
GOTO BITREV ;REPEAT IF NOT DONE
BITDONE
********** data NOW HOLDS LONG BITREV ARRAY **********
******************** BUTTERFLIES *******************
* ENTRY: temp CONTAINS BIT-REV LONG COMPLEX DATA
* A IS NO LONGER NEEDED
* Define template for var storage: (var names from Num.Recipes)
istep = 0
n = (istep)+5
mmax = (n)+5
data0 = (mmax)+5
wr = (data0)+5
wi = (wr)+21
tempr = (wi)+21
wtemp = (tempr)+21
theta = (wtemp)+21
wpi = (theta)+21
wpr = (wpi)+21
m = (wpr)+21
*FIRST, STORE SCRATCH REGS IN VARS
SETHEX
A=R4.F A ;A->istep
D1=A
D1=D1+ 5 ;D1->n
C=R2.F A ;C=size
C=C+C A ;C=2*size
DAT1=C A ;n=2*size
D1=D1+ 5 ;D1->mmax
P= 0
LC(5) #2 ;C=2
DAT1=C A ;mmax=2
D1=D1+ 5 ;D1->data0 ***STORE POINTER TO data[0] HERE
C=R1.F A ;***THIS CORRECTS FOR "FORTRAN INDEXING"
C=C-CON A,16 ;*** SINCE ALGORITHM ASSUMES INDEXES START AT 1
C=C-CON A,5 ;*** I.E. data[1] = data0+21 = FIRST SPOT
DAT1=C A
**REDUCED TO ONE SINE/MAINLOOP BY USING wtemp TO HOLD LAST SIN(theta)
** AND CALCULATING SIN(theta/2)
LC(5) wtemp ;INIT wtemp TO AVOID SIN(PI) SINCE NOT RETURN ZERO
C=C+A A ;A STILL =R4
D0=C ;D0->wtemp
A=0 W
B=0 W ;AB=0.0
GOSBVL =PUTAB0 ;STORE wtemp=0.0, D0->theta
**INIT theta TO -pi FOR FFT, +pi FOR IFFT
D1=(5) =longpi
GOSBVL =GETAB1 ;get pi from ROM
* SETDEC *********enable for FFT
* A=-A-1 S ;A=-PI *********enable for FFT
GOSBVL =PUTAB0 ;STORE IN theta
************************ BEGIN MAIN LOOP ***********************
MAINLOOP
C=R4.F A ;C->vars
D1=C
D1=D1+ 5 ;D1->n
C=DAT1 A ;C=n
D1=D1+ 5 ;POINT AT mmax
A=DAT1 A ;A=mmax
* MAINLOOP TEST
?A<C A ;IF (mmax>=n)
GOYES MAINCONT ; THEN DONE
GOTO MAINDONE
** A=mmax, D1->mmax
MAINCONT
SETHEX
A=A+A A ;A=2*mmax
D1=D1- 10 ;istep
DAT1=A A ;STORE istep=2*mmax
D1=D1+ 10
D1=D1+ 10 ;D1->wr
A=0 W
B=0 W
P= 14 ;CREATE AB=1.0
B=B+1 P
GOSBVL =PUTAB1 ;STORE wr=1.0, D1->wi
B=0 P
GOSBVL =PUTAB1 ;STORE wi=0.0, D1->tempr
D1=D1+ 16
D1=D1+ 5 ;D1->wtemp
CD1EX
RSTK=C ;PUSH ->wtemp
D0=C ;D0->wtemp
D1=C ;D1->wtemp
D1=D1+ 16
D1=D1+ 16
D1=D1+ 10 ;D1->wpi
SETDEC :***NOTE: SETDEC BEFORE ALL FP MATH (WILL STAY IN DEC)
GOSBVL =GETAB0 ;AB=wtemp, D0->theta
GOSBVL =PUTAB1 ;STORE wpi=wtemp, D1->wpr
GOSBVL =GETAB0 ;AB=theta, D0->wpi
GOSBVL =DIV2 ;AB=theta/2
D0=D0- 16
D0=D0- 5 ;D0->theta
GOSBVL =PUTAB0 ;STORE theta=theta/2
ST=1 9 ;SET UP FOR SIN (IN RADIANS)
ST=0 4
GOSBVL =aSIN15 ;**NOTE:TRASHES A,B,C,D,P,D0,R0,R1!!!
C=RSTK ;POP ->wtemp
D0=C ;D0->wtemp
GOSBVL =PUTAB0 ;STORE wtemp=SIN(theta)
C=B W
D=C W ;CD=AB
C=A W
GOSBVL =MULTF ;SQUARE wtemp
C=B W
D=C W ;CD=AB
C=A W
GOSBVL =RADDF ;DOUBLE AB (D0 TRASHED)
A=-A-1 S ;NEGATE AB
GOSBVL =PUTAB1 ;STORE wpr=-2*SQR(wtemp), D1->m
A=0 A
A=A+1 B ;A=1
DAT1=A A ;STORE m=1
************************ BEGIN OUTER LOOP ***********************
OUTLOOP SETHEX
A=R4.F A
D0=A
D0=D0+ 10 ;D0->mmax
P= 0
LC(5) m
C=C+A A
D1=C ;D1->m
A=DAT0 A ;A=mmax
C=DAT1 A ;C=m
*OUTLOOP TEST
?C<A A ;IF (m>=mmax)
GOYES NOTDONE ; THEN DONE
GOTO OUTDONE
NOTDONE RSTK=C ;PUSH i=m (i KEPT ON TOS)
A=R4.F A
D1=A ;D1->n-5
************************ BEGIN INNER LOOP ***********************
INLOOP
GOSBVL =ATTN?Lp ;IF ATTN(ON) HAS BEEN PRESSED, GETPTR AND LOOP
D1=D1+ 5 ;D1->n
A=DAT1 A ;A=n
C=RSTK ;C=i
*INLOOP TEST
?C<=A A ;IF (i>n)
GOYES INCONT ; THEN DONE
GOTO INDONE
INCONT RSTK=C ;PUSH i BACK
D1=D1+ 5 ;D1->mmax
A=DAT1 A ;A=mmax
SETHEX
A=A+C A ;A=j
B=A A
B=B+B A
B=B+B A
A=A+B A
B=B+B A
B=B+B A
A=A+B A ;A=21*j
B=C A
C=C+C A
C=C+C A
B=B+C A
C=C+C A
C=C+C A
B=C+B A ;B=21*i
D1=D1+ 5 ;D1->data0
C=DAT1 A ;C->data
A=A+C A ;A->data[j]
C=B+C A ;C->data[i]
RSTK=C ;PUSH ->data[i]
C=A A
RSTK=C ;PUSH ->data[j]
D1=D1+ 5 ;D1->wr
D0=C ;D0->data[j]
GOSBVL =GETAB1 ;AB=wr, D1->wi
GOSBVL =STAB2 ;SAVE wr IN R2/R3
GOSBVL =GETCD0 ;CD=data[j], D0->data[j+1]
SETDEC
GOSBVL =MULTF ;AB=wr*data[j]
GOSBVL =STAB0 ;SAVE IN R0/R1
GOSBVL =GETAB1 ;AB=wi, D1->tempr
GOSBVL =GETCD0 ;CD=data[j+1]
GOSBVL =MULTF ;AB=wi*data[j+1]
GOSBVL =RCCD0 ;CD=wr*data[j]
A=-A-1 S ;NEGATE AB
GOSBVL =RADDF ;AB=RESULT (D0 TRASHED)
GOSBVL =PUTAB1 ;tempr=RESULT, D1->tempi
C=RSTK
RSTK=C ;COPY ->data[j] FROM TOS
D0=C ;D0->data[j]
GOSBVL =GETCD0 ;CD=data[j], D0->data[j+1]
D1=D1- 16
D1=D1- 16
D1=D1- 10 ;D1->wi
GOSBVL =GETAB1 ;AB=wi
GOSBVL =MULTF
GOSBVL =STAB0 ;STORE wi*data[j] IN R0/R1
GOSBVL =GETCD0 ;CD=data[j+1]
GOSBVL =RCAB2 ;AB=wr
GOSBVL =MULTF ;AB=wr*data[j+1]
GOSBVL =RCCD0 ;CD=wi*data[j]
GOSBVL =RADDF ;AB=tempi (D0 TRASHED)
GOSBVL =STAB0 ;SAVE tempi IN R0/R1
C=RSTK
***STORE D1=>data[j]
D1=C ;D1->data[j]
C=RSTK
RSTK=C ;COPY ->data[i] FROM TOS
C=C+CON A,16
C=C+CON A,5 ;C->data[i+1]
RSTK=C ;PUSH ->data[i+1]
D0=C
GOSBVL =GETCD0 ;CD=data[i+1]
GOSBVL =STCD2 ;SAVE data[i+1] IN R2/R3
A=-A-1 S ;NEGATE tempi
GOSBVL =RADDF ;AB=data[i+1] - tempi
CD1EX
D1=C ;C->data[j]
D0=C
D0=D0+ 16
D0=D0+ 5 ;D0->data[j+1]
GOSBVL =PUTAB0 ;STORE data[j+1]
GOSBVL =RCAB0 ;AB=tempi
GOSBVL =RCCD2 ;CD=data[i+1]
GOSBVL =RADDF ;AB=tempi+(data[i+1])
C=RSTK ;POP ->data[i+1]
D0=C ;D0->data[i+1]
GOSBVL =PUTAB0 ;STORE data[i+1]
A=R4.F A
P= 0
LC(5) tempr
SETHEX
A=A+C A
D0=A ;D0->tempr
GOSBVL =GETAB0 ;AB=tempr
GOSBVL =STAB0 ;STORE tempr IN R0/R1
C=RSTK
RSTK=C ;COPY ->data[i] FROM TOS
D0=C ;D0->data[i]
GOSBVL =GETCD0 ;CD=data[i]
GOSBVL =STCD2 ;SAVE data[i] IN R2/R3
SETDEC
A=-A-1 S ;NEGATE tempr
GOSBVL =RADDF ;AB=data[i]-tempr
**USE D1=>data[j]
GOSBVL =PUTAB1 ;STORE data[j]
GOSBVL =RCAB0 ;AB=tempr
GOSBVL =RCCD2 ;CD=data[i]
GOSBVL =RADDF ;AB=data[i]+tempr
C=RSTK ;POP ->data[i]
D0=C ;D0->data[i]
GOSBVL =PUTAB0 ;STORE data[i]
A=R4.F A
D1=A ;D1->istep
A=DAT1 A ;A=istep
C=RSTK ;POP i
SETHEX
C=C+A A ;i=i+istep
RSTK=C ;PUSH NEW i
GOTO INLOOP ;NEXT ITERATION
************************ END INNER LOOP ***********************
INDONE
A=R4.F A
P= 0
LC(5) wr
SETHEX
C=C+A A
RSTK=C ;PUSH *wr
D0=C ;D0->wr
LC(5) wpr
C=C+A A
D1=C ;D1->wpr
SETDEC
GOSBVL =GETAB0 ;AB=wr
GOSBVL =STAB0 ;SAVE OLDwr IN R0/R1
GOSBVL =GETCD1 ;CD=wpr, D1->wpr+21
GOSBVL =MULTF ;AB=wr*wpr
GOSBVL =STAB2 ;SAVE IN R2/R3
GOSBVL =GETCD0 ;CD=wi
D1=D1- 16
D1=D1- 16
D1=D1- 10 ;D1->wpi
GOSBVL =GETAB1 ;AB=wpi, D1->wpr
GOSBVL =MULTF ;AB=wi*wpi
A=-A-1 S ;NEGATE AB
GOSBVL =RCCD2 ;CD=wr*wpr
GOSBVL =RADDF ;AB=wr*wpr - wi*wpi (D0 TRASHED)
GOSBVL =RCCD0 ;CD=OLDwr
GOSBVL =RADDF ;AB=RESULT (D0 TRASHED)
C=RSTK ;POP *wr
D0=C ;D0->wr
GOSBVL =PUTAB0 ;STORE NEW wr, D0->wi
D1=D1- 16
D1=D1- 5 ;D1->wpi
GOSBVL =GETAB1 ;AB=wpi, D1->wpr
GOSBVL =RCCD0 ;CD=OLDwr
GOSBVL =MULTF ;AB=OLDwr*wpi
GOSBVL =STAB0 ;SAVE IN R0/R1
GOSBVL =GETCD0 ;CD=wi
GOSBVL =STCD2 ;SAVE wi IN R2/R3
GOSBVL =GETAB1 ;AB=wpr, D1->m
GOSBVL =MULTF ;AB=wpr*wi
GOSBVL =RCCD0 ;CD=OLDwr*wpi
GOSBVL =RADDF ;AB=wpr*wi + OLDwr*wpi (D0 TRASHED)
GOSBVL =RCCD2 ;CD=wi
GOSBVL =RADDF ;AB=RESULT
C=R4.F A
D0=C
D0=D0+ 16
D0=D0+ 16
D0=D0+ 9 ;D0->wi
GOSBVL =PUTAB0 ;STORE NEW wi
C=DAT1 A ;C=m
SETHEX
C=C+1 A
C=C+1 A ;m=m+2
DAT1=C A ;STORE NEW m
GOTO OUTLOOP ;NEXT ITERATION
************************ END OUTER LOOP ***********************
OUTDONE SETHEX ;D0 STILL POINTS TO mmax, A=mmax
A=A+A A ;A=2*mmax
DAT0=A A ;STORE NEW mmax
GOTO MAINLOOP
************************ END MAIN LOOP ***********************
MAINDONE
* temp NOW HOLDS vars (FIRST 172 NIBS) AND LONG COMPLEX ARRAY OF ANSWERS
***************** PACK AND COPY INTO ans *****************
A=R4.F A
D1=A ;D1->vars
D1=D1+ 5 ;D1->n
A=DAT1 A ;A=n
R2=A.F A ;SAVE n IN R2 (COUNTER IN PACKLOOP)
***REMOVE FOR FFT
ASRB.F A ;A=size
************* CONVERT A.A INTO EXT.REAL ************
* USING HORNER'S RULE, 10100 = 2^4 + 2^2 = ((((1)2+0)2+1)2+0)2+0
* SINCE THE *2 AND ADDITION IS DONE IN DEC, RESULT IN DEC
* ENTER: BINT IN A.A (ANY)
* EXIT: EXT.REAL IN A/B (DEC)
* TRASHES: A,B,C,P
P= 0
LCHEX 19 ;REPEAT FOR EACH BIT IN A.A
B=0 W ;B WILL HOLD MANTISSA
BITLP B=B+B W ;B*=2
SETHEX
A=A+A A ;SHIFT MSB OF BINT INTO CARRY
SETDEC
GONC SKIP
B=B+1 W
SKIP C=C-1 B ;DEC BIT COUNTER
GONC BITLP
P= 5 ;A.A IS 0 FROM SHIFTING IN ZEROS
A=A-1 X
NORMLP A=A+1 X ;A ACCUM'S EXPONENT BY COUNTING DIGITS
BSRC ;ROTATE UNTIL MANT IS LEFT-JUSTIFIED
?B#0 WP
GOYES NORMLP
BSR W ;SHIFT MANT INTO PLACE (B.S=0)
A=0 S ;SIGN IS POSITIVE
* AB NOW HOLDS EXT.REAL
**************** END OF CONVERSION ****************
GOSBVL =STAB0 ;SAVE SIZE IN R0/R1
*******REMOVE FOR FFT
SETHEX
A=R4.F A ;**THIS IS THE LAST TIME WE USE R4(temp)
P= 0
LC(5) varnibs ;SKIP OVER VARS TO DATA
C=A+C A ;C->data
D0=C ;D0->data (SOURCE OF COPY)
D1=(5) =DSKTOP ;RECALL SAVED D1
C=DAT1 A
D1=C ;D1->STK3 (orig)
D1=D1- 10 ;D1->STK1 (ans)
C=DAT1 A
R3=C.F A ;SAVE ->ans IN R3
D1=C ;D1->ans
D1=D1+ 15 ;SKIP HEADER
D1=D1+ 10 ;D1->ansdata (DEST OF COPY)
* NOW PACK AND COPY
* NOTE: OVER/UNDERFLOW ERRORS ARE TRAPPED IN =PACK. orig IS SAVED AS TOS
PAKLOOP
GOSBVL =GETAB0 ;GET NEXT LONG
GOSBVL =RCCD0 ;RECALL SIZE *REMOVE FOR FFT
SETDEC
GOSBVL =DIVF ;DIVIDE BY SIZE *REMOVE FOR FFT
GOSBVL =PACK ;EXT.REAL->REAL
DAT1=A W ;PUT REAL IN PLACE
D1=D1+ 16 ;POINT AT NEXT REAL
C=R2.F A ;GET COUNTER
SETHEX
C=C-1 A ;DEC COUNTER
R2=C.F A
?C#0 A
GOYES PAKLOOP
* NOW RESULT IS PACKED INTO ans ARRAY
* AT LAST, REPLACE STK1 (ORIGINAL ARRAY) WITH RESULT ARRAY
A=R3.F A ;A->ans
GOVLNG =GPOverWrALp ;RESTORE REGS, PUT @A ON TOS, CONT RPL
ENDCODE
;
;
;
END_SRC
BEGIN_SRC absv.s
* ABSV.S WORKS FOR ANY SIZE VECTOR OR ARRAY
* Copyright 1992 Ken Cooke
* Takes a VECTOR/ARRAY in STK1. Returns a real VECTOR/ARRAY of
* same size, that contains the element-by-element absolute values.
TITLE ABSV
ASSEMBLE
NIBASC /HPHP48-A/
RPL
::
CK1NoBlame
CK&DISPATCH1
FOUR
::
* IF REAL ARRAY, EXIT UNCHANGED
DUP TYPECARRY? NOT?SEMI ( ORIG )
* MALLOC ANSWER ARRAY (SAME SIZE AS ORIG BUT REAL )
DUP MDIMS ITE TWO{}N ONE{}N %0 MAKEARRY ( ORIG ANS )
CODE
****************** UNLISTED ENTRY POINTS *****************
ATTN?Lp = #0CA88
SQR15 = #2B9F3
************************ MCODE ABSV ***********************
* ENTRY:
* STK2: ORIGINAL ARRY1 (COMPLEX)
* STK1: ANSWER ARRY (REAL)
* EXIT: STK1: ANSWER ARRAY (REAL)
NIBHEX 823 ;CLEAR SB AND XM IN HST
D1=D1+ 5 ;DROP STK1
D=D+1 A ;RETURN MEM
GOSBVL =SAVPTR ;SAVE REGS WITH ONLY orig ON STACK
* NOW IF USER HITS "ON" KEY, WILL LEAVE ORIGINAL STACK
A=DAT1 A ;A->ORIG
D1=D1- 5 ;POINT BACK TO ANS
C=DAT1 A ;C->ANS
R0=C.F A ;SAVE ->ANS IN R0
D0=C
D0=D0+ 15
D0=D0+ 10 ;D0->ANS DATA (OR DIM2 IF ARRAY)
D1=A ;D1->ORIG
D1=D1+ 15 ;D1->#DIMS
A=DAT1 A ;A=#DIMS
D1=D1+ 5 ;D1->DIM1
C=DAT1 A ;C=DIM1
D1=D1+ 5 ;D1->ORIG DATA (OR DIM2 IF ARRAY)
A=A-1 A ;CHECK FOR VECTOR
?A=0 A
GOYES ISVECT
***IF ARRAY, #ELEMENTS=DIM1*DIM2. ALSO, SKIP OVER DIM2 FIELD
A=DAT1 A ;A=DIM2
D1=D1+ 5 ;D1->ORIG DATA
D0=D0+ 5 ;D0->ANS DATA
SETHEX
GOSBVL =MUL# ;B=A*C
C=B A ;C=DIM1*DIM2
ISVECT C=C-1 A ;SUB 1 SO CAN USE BORROW FOR LOOPING
R4=C.F A ;SAVE #ELEMENTS-1 IN R4
ABSLOOP GOSBVL =ATTN?Lp ;ABORT IF ATTN (ON) WAS PRESSED
SETDEC
GOSUB SQDAT1 ;SQUARE RE PART
GOSBVL =STAB2 ;SAVE IN R2/R3
GOSUB SQDAT1 ;SQUARE IM PART
CD0EX
RSTK=C ;PUSH D0
GOSBVL =RCCD2 ;CD=SQR(RE)
GOSBVL =RADDF ;AB=SQR(RE)+SQR(IM) (D0 TRASHED)
GOSBVL =SQR15 ;AB=LONG ABSVAL
GOSBVL =PACK ;A=ABSVAL
C=RSTK
D0=C ;POP D0
DAT0=A W ;STORE ABSVAL IN ANSWER ARRAY
D0=D0+ 16 ;POINT AT NEXT ANSWER
* CHECK FOR DONE
C=R4.F A ;GET COUNTER
SETHEX
C=C-1 A
R4=C.F A ;DEC COUNTER
GONC ABSLOOP ;LOOP UNTIL BORROW
***NOW ANSWER CONTAINS ABS VALUES OF ORIG
GOVLNG =GPOverWrR0Lp ;GETPTR, REPLACE TOS WITH ANS, LOOP
************************ SQDAT1 ***********************************
* GET REAL FROM DAT1, UNPACK AND SQUARE IT, INC D1
SQDAT1 A=DAT1 W ;A=REAL
D1=D1+ 16 ;D1->NEXT
GOSBVL =SPLITA ;AB=REAL
C=B W
D=C W
C=A W ;COPY INTO CD
GOSBVL =MULTF ;AB=SQR(REAL)
RTN
ENDCODE
;
;
