C06FAF is the NAG routine for finding the FFT of a one-dimensional sequence of real data values. The routine performs a forward transform, storing the FFT as a ``Hermitian'' sequence in which only half of the real and imaginary terms are kept. The inverse transform is obtained by calling C06FBF, which finds the FFT of a one-dimensional Hermitian sequence.
The following steps are involved in replacing C06FAF calls with equivalent FFTPACK calls:
The work array passed to the FFT routine needs to be increased in size from N elements (for NAG) to 2*N+15 (for FFTPACK).
Replace the call
DOUBLE PRECISION X(N), WORK(N)
CALL C06FAF( X, N, WORK, IFAIL )
with
DOUBLE PRECISION X(N), WORK(2*N+15)
CALL PDA_DRFFTF( N, X, WORK )
The work array supplied to PDA_DRFFTF needs initialising before calling PDA_DRFFTF. This is done by calling PDA_DRFFTI:
DOUBLE PRECISION WORK(2*N+15)
CALL PDA_DRFFTI( N, WORK )
There is no need to re-initialise WORK if the value of N has not changed since the previous call to PDA_DRFFTI (and if the contents of the work array have not been altered). No harm will occur (except for significant slowing down of execution) if the WORK array is unnecessarily re-initialised, but it is a good idea to include some logic to prevent this.
Compared to the Fourier coefficients created by NAG, those created by FFTPACK are stored in a different order in the output array and are normalised differently. You can either modify your application to use the FFTPACK format throughout, or call the PDA_DR2NAG routine to convert the FFTPACK results into NAG format.
DOUBLE PRECISION X(N)
CALL PDA_DR2NAG( N, X )
where X is the output from PDA_DRFFTF. On return, X holds a NAG-style Hermitian sequence.
PDA [1ex