Taking the Fourier transform

Actually taking the Fourier transform of a complex data structure is quite straightforward. The forward transform is performed by the Figaro function `fft' and the reverse transform is performed by `bfft'. The only parameters for `fft' and `bfft' are the names of the input and output structures.

For example, to calculate the Fourier transform of the spectrum held in the non-complex structure `myspect', it should be enough to do:

   ICL> r2cmplx myspect cmplxspect
   ICL> fft  cmplxspect cmplxspect

which results in `cmplxspect' containing the `fft' of `myspect'. If the power spectrum of `myspect' is required it can be obtained by `cmplx2m'. For example, it may be plotted using:

   ICL> cmplx2m cmplxspect modulus
   ICL> splot   modulus reset accept

Actually, the modulus is the square root of the power spectrum; the statement that the power spectrum is available represents a slightly cavalier attitude towards the proper usage of terms. However, if the power spectrum is only needed in order to get a feel for the frequency distribution of the original data, then the modulus will generally do just as well.

Often, a power spectrum needs to be plotted logarithmically to produce a sensible plot, so it is quite acceptable to try

   ICL> ilog  modulus modulus
   ICL> splot modulus reset accept

The FFT routines generate data with the low frequency terms in the lowest numbered elements in each dimension of the data. `fft' (the Figaro program) re-orders the transformed data so that the zero frequency term is in the centre of each dimension, and the resulting data now goes from -N to +N (where N is the Nyquist frequency). The plots of the modulus made as described above will show a plot with an axis labeled from -1 to +1, (the unit being the Nyquist frequency), that is symmetrical about the centre and which peaks (for most input data) in the centre. New axis structures are created by `fft' to reflect the new axis values. The reason for this re-ordering of the data is that it seems to be easier to visualise complex filters, particularly in 2-dimensions, when the low frequency terms are collected in the middle of the data rather than scattered into the corners.

The reverse FFT is performed by `bfft'. Note that the precise definition of `forward transform' and `reverse transform' differ between FFT implementations. The PDA routines, for example, do not introduce any scaling factors between data that is first forward and then reverse transformed and the original data. This means that the sequence

   ICL> r2cmplx myspect cmplxspect
   ICL> fft  cmplxspect cmplxspect
   ICL> bfft cmplxspect cmplxspect
   ICL> cmplx2r cmplxspect newspect

should generate a `newspect' that is exactly the same as the original `myspect', except that any axis information will have been lost.