\ Cmplx - Forth Complex Arithmetic Support Words for KFFT
\ Jerry Kallaus  02/14/89
\

INCLUDE? asm_fft? sp:fftcontrols

asm_fft? .IF  INCLUDE?  Z*  sp:fft.asm  .THEN

anew task-cmplx

\ ------------------ Basic support words  -------------------------

float_fft? .IF
: s+ COMPILE F+ ;  IMMEDIATE
: s- COMPILE F- ;  IMMEDIATE
: s* COMPILE F* ;  IMMEDIATE
: s/ COMPILE F/ ;  IMMEDIATE
: SNEGATE  COMPILE FNEGATE ; IMMEDIATE
: ZSCALE.DOWN DROP both ;

.ELSE
: S+  +  both ;
: S-  -  both ;
: S*  *  [ scale_fft negate ] literal ashift     ;
: S/  swap [ scale_fft ] literal ashift swap /   ;
: SNEGATE  NEGATE  both ;
scale_fft 2**  CONSTANT rone_fft
.THEN

float_fft? jforth2? AND .IF
pi   CONSTANT pi_fft
1.0  CONSTANT rone_fft
.THEN

float_fft? jforth2? NOT AND .IF
3.14159265+0 CONSTANT pi_fft
1+0          CONSTANT rone_fft        
: fsin   compile sin    ; immediate
: fcos   compile cos    ; immediate
: f.r    compile fp.rd  ; immediate
.THEN

float_fft? asm_fft? OR NOT .IF
: ZSCALE.DOWN  ?DUP IF negate dup>r ashift swap r> ashift swap THEN both ;
.THEN

fixasm_fft?  .IF
: ZSCALE.DOWN  compile  Z/2**N ; immediate
.THEN

asm_fft? NOT .IF
: 2CELL+  CELL+ CELL+  both ;
: 2CELL-  CELL- CELL-  both ;
: 2CELLS  CELLS 2*  both ;
: 4DUP    4 XDUP  both ;
: Z@      dup>r @ r> cell+ @ both ;
: Z!      dup>r cell+ ! r> ! both ;
.THEN

: 2CELL   8       both ;
: ZCELL   2CELL   both ;
: ZCELL+  2cell+  both ;
: ZCELL-  2cell-  both ;
: ZCELLS  2cells  both ;
: ZDROP   2drop   both ;
: ZDUP    2dup    both ;
: ZOVER   2over   both ;
: Z2DUP   4DUP    both ;

\ -----------------  Complex Arithmetic Stack Words  ---------------

fixasm_fft? NOT .IF
: Z+  ( a b c d --  a+c b+d )  rot s+ >r s+ r> both ;
: Z-  ( a b c d  --  a-c b-d)  rot swap s- >r s- r> both ;
: Z*  { a b c d  ---  ac-bd ad+bc }
            a c s* b d s* s-  a d s* b c s* s+  ;

\ : Z*  ( a b c d  --  ac-bd ad+bc ) 
\          2over 2over -rot s*  rot  s* s+ >r
\                      rot s* -rot  s* swap s-  r>  ;

: ZNEGATE ( a b -- -a -b ) 
    snegate swap snegate swap  both ;

.THEN

: CONJG  ( z -- conjugate z )  snegate  both ;
: ZI*    ( a b -- -b a , cmplex multipy by i )  snegate  swap  both ;

float_fft? .IF
: ZEXP ( z -- cosz sinz ) dup>r fcos r> fsin ;
.THEN

\ -----------------  Forth Complex Data Type Words  -------------------

\ Complex number definitions.
\ Convention for complex number on stack  is the imaginary part is on top.
\ Convention for complex number in memory is the real part is at lower addr.
\
: ZCONSTANT   create swap , ,  does>  Z@ ;
: ZVARIABLE   create 0 0 , ,   does>     ;
: ZARRAY  ( #elements -- )             ( ex: 1024 ZARRAY myarray         )
    create  2cells allot
    does>   swap 2cells + ;            ( i myarray         gets ith addr )

: ZPTR  ( addr-of-pointer -- )         ( ex: VARIABLE mypointer           )
    create ,                           ( mypointer ZPTR c-ptr             )
    does> @ @ swap 2cells + ;          ( 0 myarray mypointer !            )
                                       ( i c-ptr            gets ith addr )

rone_fft 0  ZCONSTANT  Z1