The NAG routine expects real and imaginary parts in separate arrays, whereas FFTPACK expects them in the same array, with corresponding real and imaginary values in adjacent elements. If the application can be changed to supply the input data in this format, so well and good. Otherwise you will have to have an extra work array in which to hold the input (and output) data in FFTPACK format. You would convert the supplied input data using code such as:
DOUBLE PRECISION X( N ), Y( N ), C( 2*N )
DO J = 1, N
I = 2*J
C( I - 1 ) = X( J )
C( I ) = Y( J )
END DO
or
DOUBLE PRECISION X( N ), Y( N ), C( 2, N )
DO J = 1, N
C( 1, J ) = X( J )
C( 2, J ) = Y( J )
END DO
where the X and Y arrays hold the supplied data, C is a work array, and N is the number of data points.
The work array passed to the FFT routine needs to be increased in size from N elements (for NAG) to 4*N+15 (for FFTPACK).
Replace the call
DOUBLE PRECISION X(N), Y(N), WORK(N)
CALL C06FCF( X, Y, N, WORK, IFAIL )
with
DOUBLE PRECISION C(2*N), WORK(4*N+15)
CALL PDA_DCFFTF( N, C, WORK )
The work array supplied to PDA_DCFFTF needs initialising before calling PDA_DCFFTF. This is done by calling PDA_DCFFTI:
DOUBLE PRECISION WORK( 4*N+15 )
CALL PDA_DCFFTI( N, WORK )
There is no need to re-initialise WORK if the value of N has not changed since the previous call to PDA_DCFFTI (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.
The Fourier coefficients created by FFTPACK are stored in a single array and are not normalised, whereas NAG stores them in two arrays and normalises them. You can either modify the way your application to use the FFTPACK format instead of the NAG format, or call the PDA_DC2NAG routine to convert the FFTPACK results into NAG format.
DOUBLE PRECISION X(N), Y(N), C(2*N)
CALL PDA_DC2NAG( N, C, X, Y )
where C is the output from PDA_DCFFTF, and X and Y hold the corresponding real and imaginary coefficients as returned by C06FCF.
PDA [1ex