FFTCB

Computes the complex periodic sequence from its Fourier coefficients.

Required Arguments

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

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

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

FORTRAN 90 Interface

Generic:          CALL FFTCB (N, COEF, SEQ)

Specific:         The specific interface names are S_FFTCB and D_FFTCB.

FORTRAN 77 Interface

Single:            CALL FFTCB (N, COEF, SEQ)

Double:          The double precision name is DFFTCB.

Description

The routine FFTCB computes the inverse 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 c = COEF, FFTCB returns in s = SEQ

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

Finally, note that we can invert the inverse 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. FFTCB is based on the complex inverse 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 F2TCB/DF2TCB. The reference is:

CALL F2TCB (N, COEF, SEQ, 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 FFTCB 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 first compute the Fourier transform of the vector x, where xj = j for j = 1 to N. Note that the norm of x is (N[N + 1][2N + 1]/6) 1/2, and hence, the norm of the transformed vector

is N([N + 1][2N + 1]/6) 1/2. The vector

is used as input into FFTCB with the resulting output s = Nx, that is, sj = jN, for j = 1 to N.

 

      USE FFTCB_INT

      USE FFTCF_INT

      USE UMACH_INT

 

      IMPLICIT   NONE

      INTEGER    N

      PARAMETER  (N=7)

!

      INTEGER    I, NOUT

      COMPLEX    CMPLX, SEQ(N), COEF(N), X(N)

      INTRINSIC  CMPLX

!                                 This loop fills out the data vector

!                                 with X(I)=I, I=1,N

      DO 10  I=1, N

         X(I) = CMPLX(I,0)

   10 CONTINUE

!                                 Compute the forward transform of X

      CALL FFTCF (N, X, COEF)

!                                 Compute the backward transform of

!                                 COEF

      CALL FFTCB (N, COEF, SEQ)

!                                 Get output unit number

      CALL UMACH (2, NOUT)

!                                 Print results

      WRITE (NOUT,99998)

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

99998 FORMAT (5X, 'INDEX', 9X, 'INPUT', 9X, 'FORWARD TRANSFORM', 3X, &

            'BACKWARD TRANSFORM')

99999 FORMAT (1X, I7, 7X,'(',F5.2,',',F5.2,')', &

                     7X,'(',F5.2,',',F5.2,')', &

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

      END

Output

 

INDEX      INPUT         FORWARD TRANSFORM   BACKWARD TRANSFORM
 1       ( 1.00, 0.00)       (28.00, 0.00)       ( 7.00, 0.00)
 2       ( 2.00, 0.00)       (-3.50, 7.27)       (14.00, 0.00)
 3       ( 3.00, 0.00)       (-3.50, 2.79)       (21.00, 0.00)
 4       ( 4.00, 0.00)       (-3.50, 0.80)       (28.00, 0.00)
 5       ( 5.00, 0.00)       (-3.50,-0.80)       (35.00, 0.00)
 6       ( 6.00, 0.00)       (-3.50,-2.79)       (42.00, 0.00)
 7       ( 7.00, 0.00)       (-3.50,-7.27)       (49.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