This is C code for 2 dim FFTs for data length a power of two. thanks to John W. Herman for the pointer to the code ! this is also avialable at ftp://ftp.teleos.com/VISION-LIST-ARCHIVE/SHAREWARE/CODE/FFT/FFT.shar from the file: ** forward_fft(array, rows, cols) ** COMPLEX *array; ** int rows, cols; ** ** inverse_fft(array, rows, cols) ** COMPLEX *array; ** int rows, cols; ** ** These entry points compute the forward and inverse DFT's, respectively, ** of a single-precision COMPLEX array. ** ** The result is a COMPLEX array of the same size, returned in ** the same space as the input array. That is, the original array is ** overwritten and destroyed. ** ** Rows and columns must each be an integral power of 2. ** ** These routines return integer value ERROR if an error was detected, ** NO_ERROR otherwise ** ** This implementation of the DFT uses the transform pair defined as follows. ** ** Let there be two COMPLEX arrays each with n rows and m columns ** Index them as ** f(x,y): 0 <= x <= m - 1, 0 <= y <= n - 1 ** F(u,v): -m/2 <= u <= m/2 - 1, -n/2 <= v <= n/2 - 1 ** ** Then the forward and inverse transforms are related as ** ** Forward: ** ** F(u,v) = \sum_{x=0}^{m-1} \sum_{y=0}^{n-1} ** f(x,y) \exp{-2\pi i (ux/m + vy/n)} ** ** ** Inverse: ** ** f(x,y) = 1/(mn) \sum_{u=-m/2}^{m/2-1} \sum_{v=-n/2}^{n/2-1} ** F(u,v) \exp{2\pi i (ux/m + vy/n)} ** ** Therefore, the transforms have these properties: ** 1. \sum_x \sum_y f(x,y) = F(0,0) ** 2. m n \sum_x \sum_y |f(x,y)|^2 = \sum_u \sum_v |F(u,v)|^2