mirror/kissfft
1
0
Fork 0
mirror of https://github.com/mborgerding/kissfft.git synced 2025-05-27 21:20:27 -04:00
Code Issues Packages Projects Releases Wiki Activity

Added mucho comments.

Browse Source
This commit is contained in:
Mark Borgerding 2003-12-30 15:18:41 +00:00
parent ec3b64a62e
commit b468fd96d7
2 changed files with 68 additions and 4 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

4
tools/kfc.c
View File

@ -32,12 +32,12 @@ static const void * find_cached_fft(int nfft,int inverse)
cached_fft * prev=NULL;
while ( cur ) {
if ( cur->nfft == nfft && inverse == cur->inverse )
break;//found the right node
break;/*found the right node*/
prev = cur;
cur = (cached_fft*)prev->next;
}
if (cur== NULL) {
// no cached node found, need to create a new one
/* no cached node found, need to create a new one*/
kiss_fft_alloc(nfft,inverse,0,&len);
cur = (cached_fft*)malloc(sizeof(cached_fft) + len );
if (cur == NULL)

68
tools/kiss_fftnd.c
View File

@ -1,3 +1,5 @@
/*
Copyright (c) 2003, Mark Borgerding
@ -76,6 +78,68 @@ void * kiss_fftnd_alloc(int *dims,int ndims,int inverse_fft,void*mem,size_t*lenm
return st;
}
/*
This works by tackling one dimension at a time.
In effect,
Each stage starts out by reshaping the matrix into a DixSi 2d matrix.
A Di-sized fft is taken of each column, transposing the matrix as it goes.
Here's a 3-d example:
Take a 2x3x4 matrix, laid out in memory as a contiguous buffer
[ [ [ a b c d ] [ e f g h ] [ i j k l ] ]
[ [ m n o p ] [ q r s t ] [ u v w x ] ] ]
Stage 0 ( D=2): treat the buffer as a 2x12 matrix
[ [a b ... k l]
[m n ... w x] ]
FFT each column with size 2.
Transpose the matrix at the same time using kiss_fft_stride.
[ [ a+m a-m ]
[ b+n b-n]
...
[ k+w k-w ]
[ l+x l-x ] ]
Note fft([x y]) == [x+y x-y]
Stage 1 ( D=3) treats the buffer (the output of stage D=2) as an 3x8 matrix,
[ [ a+m a-m b+n b-n c+o c-o d+p d-p ]
[ e+q e-q f+r f-r g+s g-s h+t h-t ]
[ i+u i-u j+v j-v k+w k-w l+x l-x ] ]
And perform FFTs (size=3) on each of the columns as above, transposing
the matrix as it goes. The output of stage 1 is
(Legend: ap = [ a+m e+q i+u ]
am = [ a-m e-q i-u ] )
[ [ sum(ap) fft(ap)[0] fft(ap)[1] ]
[ sum(am) fft(am)[0] fft(am)[1] ]
[ sum(bp) fft(bp)[0] fft(bp)[1] ]
[ sum(bm) fft(bm)[0] fft(bm)[1] ]
[ sum(cp) fft(cp)[0] fft(cp)[1] ]
[ sum(cm) fft(cm)[0] fft(cm)[1] ]
[ sum(dp) fft(dp)[0] fft(dp)[1] ]
[ sum(dm) fft(dm)[0] fft(dm)[1] ] ]
Stage 2 ( D=4) treats this buffer as a 4*6 matrix,
[ [ sum(ap) fft(ap)[0] fft(ap)[1] sum(am) fft(am)[0] fft(am)[1] ]
[ sum(bp) fft(bp)[0] fft(bp)[1] sum(bm) fft(bm)[0] fft(bm)[1] ]
[ sum(cp) fft(cp)[0] fft(cp)[1] sum(cm) fft(cm)[0] fft(cm)[1] ]
[ sum(dp) fft(dp)[0] fft(dp)[1] sum(dm) fft(dm)[0] fft(dm)[1] ] ]
Then FFTs each column, transposing as it goes.
The resulting matrix is the 3d FFT of the 2x3x4 input matrix.
Note as a sanity check that the first element of the final
stage's output (DC term) is
sum( [ sum(ap) sum(bp) sum(cp) sum(dp) ] )
, i.e. the summation of all 24 input elements.
*/
void kiss_fftnd(const void * cfg,const kiss_fft_cpx *fin,kiss_fft_cpx *fout)
{
int i,k;
@ -92,8 +156,8 @@ void kiss_fftnd(const void * cfg,const kiss_fft_cpx *fin,kiss_fft_cpx *fout)
for ( k=0; k < st->ndims; ++k) {
int curdim = st->dims[k];
int stride = st->dimprod / curdim;
for (i=0;i<stride;++i)
kiss_fft_stride( st->states[k], bufin+i , bufout+i*curdim,stride );
for ( i=0 ; i<stride ; ++i )
kiss_fft_stride( st->states[k], bufin+i , bufout+i*curdim, stride );
/*toggle back and forth between the two buffers*/
if (bufout == st->tmpbuf){