Computes the Fourier coefficients of a complex periodic sequence.

Required Arguments

N — Length of the sequence to be transformed.   (Input)

SEQ — Complex array of length N containing the periodic sequence.   (Input)

COEF — Complex array of length N containing the Fourier coefficients.   (Output)

FORTRAN 90 Interface

Generic:          CALL FFTCF (N, SEQ, COEF)

Specific:         The specific interface names are S_FFTCF and D_FFTCF.

FORTRAN 77 Interface

Single:            CALL FFTCF (N, SEQ, COEF)

Double:          The double precision name is DFFTCF.

Description

The routine FFTCF computes the discrete complex Fourier transform of a complex vector of size N. The method used is a variant of the Cooley-Tukey algorithm, which is most efficient when N is a product of small prime factors. If N satisfies this condition, then the computational effort is proportional to N log N. This considerable savings has historically led people to refer to this algorithm as the “fast Fourier transform” or FFT.

Specifically, given an N-vector x, FFTCF returns in c = COEF

Furthermore, a vector of Euclidean norm S is mapped into a vector of norm

Finally, note that we can invert the Fourier transform as follows:

This formula reveals the fact that, after properly normalizing the Fourier coefficients, one has the coefficients for a trigonometric interpolating polynomial to the data. An unnormalized inverse is implemented in FFTCB. FFTCF is based on the complex FFT in FFTPACK. The package FFTPACK was developed by Paul Swarztrauber at the National Center for Atmospheric Research.

Comments

1.         Workspace may be explicitly provided, if desired, by use of F2TCF/DF2TCF. The reference is:

CALL F2TCF (N, SEQ, COEF, WFFTC, CPY)

The additional arguments are as follows:

WFFTC — Real array of length 4 * N + 15 initialized by FFTCI. The initialization depends on N.   (Input)

CPY — Real array of length 2 * N. (Workspace)

2.         The routine FFTCF is most efficient when N is the product of small primes.

3.         The arrays COEF and SEQ may be the same.

4.         If FFTCF/FFTCB is used repeatedly with the same value of N, then call FFTCI followed by repeated calls to F2TCF/F2TCB. This is more efficient than repeated calls to FFTCF/FFTCB.

Example

In this example, we input a pure exponential data vector and recover its Fourier series, which is a vector with all components zero except at the appropriate frequency where it has an N. Notice that the norm of the input vector is

but the norm of the output vector is N.

 

      USE FFTCF_INT

      USE CONST_INT

      USE UMACH_INT

 

      IMPLICIT   NONE

      INTEGER    N

      PARAMETER  (N=7)

!

      INTEGER    I, NOUT

      REAL       TWOPI

      COMPLEX    C, CEXP, COEF(N), H, SEQ(N)

      INTRINSIC  CEXP

!

      C     = (0.,1.)

      TWOPI = CONST('PI')

      TWOPI = 2.0 * TWOPI

!                                 Here we compute (2*pi*i/N)*3.

      H = (TWOPI*C/N)*3.

!                                 This loop fills out the data vector

!                                 with a pure exponential signal of

!                                 frequency 3.

      DO 10  I=1, N

         SEQ(I) = CEXP((I-1)*H)

   10 CONTINUE

!                                 Compute the Fourier transform of SEQ

      CALL FFTCF (N, SEQ, COEF)

!                                 Get output unit number and print

!                                 results

      CALL UMACH (2, NOUT)

      WRITE (NOUT,99998)

99998 FORMAT (9X, 'INDEX', 8X, 'SEQ', 15X, 'COEF')

      WRITE (NOUT,99999) (I, SEQ(I), COEF(I), I=1,N)

99999 FORMAT (1X, I11, 5X,'(',F5.2,',',F5.2,')', &

                      5X,'(',F5.2,',',F5.2,')')

      END

Output

 

INDEX        SEQ               COEF
  1     ( 1.00, 0.00)     ( 0.00, 0.00)
  2     (-0.90, 0.43)     ( 0.00, 0.00)
  3     ( 0.62,-0.78)     ( 0.00, 0.00)
  4     (-0.22, 0.97)     ( 7.00, 0.00)
  5     (-0.22,-0.97)     ( 0.00, 0.00)
  6     ( 0.62, 0.78)     ( 0.00, 0.00)
  7     (-0.90,-0.43)     ( 0.00, 0.00)

Visual Numerics, Inc. Visual Numerics - Developers of IMSL and PV-WAVE http://www.vni.com/ PHONE: 713.784.3131 FAX:713.781.9260