diff options
Diffstat (limited to 'src/signalProcessing/fft/fft842.c')
| -rw-r--r-- | src/signalProcessing/fft/fft842.c | 173 |
1 files changed, 0 insertions, 173 deletions
diff --git a/src/signalProcessing/fft/fft842.c b/src/signalProcessing/fft/fft842.c deleted file mode 100644 index 162490fe..00000000 --- a/src/signalProcessing/fft/fft842.c +++ /dev/null @@ -1,173 +0,0 @@ -/* - * Scilab ( http://www.scilab.org/ ) - This file is part of Scilab - * Copyright (C) 2008 - INRIA - Allan SIMON - * - * This file must be used under the terms of the CeCILL. - * This source file is licensed as described in the file COPYING, which - * you should have received as part of this distribution. The terms - * are also available at - * http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt - * - */ - -#include "fft_internal.h" - - - -/* get binary log of integer argument; exact if n a power of 2 */ -static int fastlog2( int n) -{ - int log = -1; - while(n) { - log++; - n >>= 1; - } - return(log); -} - - - -void fft842 (doubleComplex* b, int size , int in) -{ - double fn; - doubleComplex temp ; - - int L[16],L1,L2,L3,L4,L5,L6,L7,L8,L9,L10,L11,L12,L13,L14,L15; - int j1,j2,j3,j4,j5,j6,j7,j8,j9,j10,j11,j12,j13,j14; - int i = 0, j, ij, ji, ij1, ji1; - int n2pow, n8pow, nthpo, ipass, nxtlt, lengt; - - n2pow = fastlog2( size ); - nthpo = size ; - fn = 1.0 / (double)nthpo; - - - if(in==FORWARD) - /* take conjugate */ - for(i=0;i< size ;i++) - { - b[i] = DoubleComplex ( zreals( b[i]) , - zimags (b[i])); - } - - - - n8pow = n2pow/3; - - if(n8pow) - { - /* if the size if a factor of a power of 8 we call r8tx */ - /* radix 8 iterations */ - for(ipass=1;ipass<=n8pow;ipass++) - { - nxtlt = 0x1 << (n2pow - 3*ipass); - lengt = 8*nxtlt; - - - - r8tx(nxtlt,nthpo,lengt, - b,b+nxtlt,b+2*nxtlt, - b+3*nxtlt,b+4*nxtlt,b+5*nxtlt, - b+6*nxtlt,b+7*nxtlt); - - - } - } - -/* if the size can be written this way 2^(3*n + 1) , then we call the radix 2 function - if can be written this way 2^(3*n + 1) the we call the radix 4 function */ - - if(n2pow%3 == 1) - { - /* radix 2 iteration needed */ - r2tx(nthpo,b,b+1); - } - - if(n2pow%3 == 2) - { - /* radix 4 iteration needed */ - - r4tx(nthpo,b,b+1,b+2,b+3); - } - - - - for(j=1;j<=15;j++) - { - L[j] = 1; - if(j-n2pow <= 0) L[j] = 0x1 << (n2pow + 1 - j); - } - /* this part can maybe be improved */ - - L15=L[1];L14=L[2];L13=L[3];L12=L[4];L11=L[5];L10=L[6];L9=L[7]; - L8=L[8];L7=L[9];L6=L[10];L5=L[11];L4=L[12];L3=L[13];L2=L[14];L1=L[15]; - - ij = 1; - -/* all the following instruction is to unscramble the output */ - for(j1=1;j1<=L1;j1++) - for(j2=j1;j2<=L2;j2+=L1) - for(j3=j2;j3<=L3;j3+=L2) - for(j4=j3;j4<=L4;j4+=L3) - for(j5=j4;j5<=L5;j5+=L4) - for(j6=j5;j6<=L6;j6+=L5) - for(j7=j6;j7<=L7;j7+=L6) - for(j8=j7;j8<=L8;j8+=L7) - for(j9=j8;j9<=L9;j9+=L8) - for(j10=j9;j10<=L10;j10+=L9) - for(j11=j10;j11<=L11;j11+=L10) - for(j12=j11;j12<=L12;j12+=L11) - for(j13=j12;j13<=L13;j13+=L12) - for(j14=j13;j14<=L14;j14+=L13) - for(ji=j14;ji<=L15;ji+=L14) - - { - ij1 = ij-1; - ji1 = ji-1; - - if(ij-ji<0) - { - temp = b[ij1]; - b[ij1] = b[ji1]; - b[ji1] = temp; - - /* - r = b[ij1].re; - b[ij1].re = b[ji1].re; - b[ji1].re = r; - fi = b[ij1].im; - b[ij1].im = b[ji1].im; - b[ji1].im = fi; - */ - } - ij++; - } - - if(in==FORWARD) /* take conjugates & unscramble outputs */ - for(i=0,j=size/2;j<size;i++,j++) - { - temp = DoubleComplex ( zreals ( b[j] ) ,- zimags( b[j] )); - b[j] = DoubleComplex ( zreals ( b[i] ) ,- zimags( b[i] )); - b[i] = DoubleComplex ( zreals ( temp ) , zimags( temp )); - - - /* r = b[j].re; fi = b[j].im; - b[j].re = b[i].re; b[j].im = -b[i].im; - b[i].re = r; b[i].im = -fi; - */ - } - - - -/* this code is two fix a problem of result order which was different from what scilab give */ - for ( i = 0 ; i < size /2 ; i++) - { - temp = b[i] ; - b[i] = b[i+(size/2)]; - b[i+(size/2)]= temp ; - - - - } - - -} |