diff options
Diffstat (limited to '2.3-1/src/c/signalProcessing/fft/zfftma.c')
-rw-r--r-- | 2.3-1/src/c/signalProcessing/fft/zfftma.c | 156 |
1 files changed, 156 insertions, 0 deletions
diff --git a/2.3-1/src/c/signalProcessing/fft/zfftma.c b/2.3-1/src/c/signalProcessing/fft/zfftma.c new file mode 100644 index 00000000..48616eaf --- /dev/null +++ b/2.3-1/src/c/signalProcessing/fft/zfftma.c @@ -0,0 +1,156 @@ +/* + * 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 + * + */ +#define FFT842 1 +#define DFFT2 0 + +#include <stdlib.h> +#include <math.h> +#include "fft.h" +#include "lapack.h" +#include "fft_internal.h" + +void zfftma ( doubleComplex* in , int rows, int cols, doubleComplex* out) +{ + + int choosenAlgo = DFFT2 ; + + int size = rows*cols ; + int sizeTemp = 0; + + int rowsTemp = 0 ; + int colsTemp = 0 ; + + int ierr = 0 ; + int isn = -1; + int i = 0; + + int increment=1; + + + double* realIn = (double*) malloc ( sizeof (double) * (unsigned int) size ); + double* imagIn = (double*) malloc ( sizeof (double) * (unsigned int) size ); + doubleComplex* inCopy = (doubleComplex*) malloc ( sizeof (doubleComplex) * (unsigned int) size); + + doubleComplex* inTemp = (doubleComplex*) malloc ( sizeof (doubleComplex) * (unsigned int) size ); + + zimaga ( in , size , imagIn) ; + zreala ( in , size , realIn) ; + for(i=0;i<size;i++) inCopy[i]=in[i]; + + if ( rows == 1 || cols == 1 ) + { + sizeTemp = (int) pow(2.0, (int) (log(size + 0.5) / log(2.0))); + if ( size == sizeTemp ) + { + if ( size <= pow(2.0, 15.0)) + { + fft842 ( inCopy , size , 0 ); + choosenAlgo = FFT842 ; + } + else + { + dfft2 ( realIn , imagIn , 1 , size , 1 , isn , ierr ); + } + } + else + { + dfft2 ( realIn , imagIn , 1 , size , 1 , isn , ierr ); + } + } + else + { + rowsTemp = (int) pow(2.0, (int)(log(rows + 0.5) / log(2.0))) ; + colsTemp = (int) pow(2.0 ,(int)(log(cols + 0.5) / log(2.0))) ; + + if ( rows == rowsTemp) + { + if ( rows <= pow(2.0, 15.0)) + { + for ( i = 0 ; i < cols ; i++ ) + { + fft842 ( &inCopy[ rows*i] , rows , 0); + /* stock new inCopy in realIn and imagIn + if the second call don't call fft842 + ex : matrix 2*3 */ + zimaga ( inCopy , size , imagIn) ; + zreala ( inCopy , size , realIn) ; + } + } + else + { + dfft2 ( realIn, imagIn ,cols , rows , 1 , isn , ierr); + /* stock new realIn and imagIn in inCopy + if the second call call fft842 + ex : matrix 3*2 */ + inCopy=DoubleComplexMatrix(realIn,imagIn,size); + } + } + else + { + dfft2 ( realIn, imagIn ,cols , rows , 1 , isn , ierr); + /* stock new realIn and imagIn in inCopy + if the second call call fft842 + ex : matrix 3*2 */ + inCopy=DoubleComplexMatrix(realIn,imagIn,size); + } + /*second call*/ + if ( colsTemp == cols ) + { + if ( cols <= pow(2.0, 15.0) ) + { + /*compute the fft on each line of the matrix */ + for (i = 0 ; i < rows ; i++ ) + { + C2F(zcopy) ( &cols, inCopy + i, &rows, inTemp , &increment ); + + fft842( inTemp , cols , 0); + choosenAlgo = FFT842 ; + C2F(zcopy) ( &cols, inTemp , &increment, inCopy + i, &rows ); + + } + } + else + { + dfft2 ( realIn, imagIn, 1, cols, rows, isn, ierr); + } + } + else + { + dfft2 ( realIn, imagIn, 1, cols, rows, isn, ierr); + } + + } + + + + if ( choosenAlgo == FFT842 ) + { + for ( i = 0 ; i < size ; i++) + { + out[i] = DoubleComplex ( zreals(inCopy[i]) , zimags(inCopy[i]) ); + } + } + else + { + for ( i = 0 ; i < size ; i++) + { + out[i] = DoubleComplex ( realIn[i] , imagIn[i] ); + } + + } + + free(realIn); + free(imagIn); + free(inCopy); + free(inTemp); + +} |