QSINF

Computes the coefficients of the sine Fourier transform with only odd wave numbers.

Required Arguments

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

SEQ — Array of length N containing the sequence. (Input)

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

FORTRAN 90 Interface

Generic: CALL QSINF (N, SEQ, COEF)

Specific: The specific interface names are S_QSINF and D_QSINF.

FORTRAN 77 Interface

Single: CALL QSINF (N, SEQ, COEF)

Double: The double precision name is DQSINF.

Description

The routine QSINF computes the discrete Fourier quarter sine transform of a real 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.

Specifically, given an N-vector s = SEQ, QSINF returns in c = COEF

Finally, note that the Fourier quarter sine transform has an (unnormalized) inverse, which is implemented in the IMSL routine QSINB. The routine QSINF is based on the quarter sine 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 Q2INF/DQ2INF. The reference is:

CALL Q2INF (N, SEQ, COEF, WQSIN)

The additional argument is:

WQSIN — Array of length 3 * N + 15 initialized by QSINI. The initialization depends on N. (Input)

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

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

4. If QSINF/QSINB is used repeatedly with the same value of N, then call QSINI followed by repeated calls to Q2INF/Q2INB. This is more efficient than repeated calls to QSINF/QSINB.

Example

In this example, we input a pure quarter sine wave as a data vector and recover its Fourier quarter sine series.

USE QSINF_INT

USE CONST_INT

USE UMACH_INT

IMPLICIT NONE

INTEGER N

PARAMETER (N=7)

!

INTEGER I, NOUT

REAL COEF(N), FLOAT, PI, SIN, SEQ(N)

INTRINSIC FLOAT, SIN

! Get output unit number

CALL UMACH (2, NOUT)

! Fill the data vector SEQ

! with a pure sine wave

PI = CONST('PI')

DO 10 I=1, N

SEQ(I) = SIN(FLOAT(I)*(PI/2.0)/FLOAT(N))

10 CONTINUE

! Compute the transform of SEQ

CALL QSINF (N, SEQ, COEF)

! Print results

WRITE (NOUT,99998)

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

99998 FORMAT (9X, 'INDEX', 6X, 'SEQ', 7X, 'COEF')

99999 FORMAT (1X, I11, 5X, F6.2, 5X, F6.2)

END

Output

INDEX SEQ COEF

1 0.22 7.00

2 0.43 0.00

3 0.62 0.00

4 0.78 0.00

5 0.90 0.00

6 0.97 0.00

7 1.00 0.00