NCEPLIBS-ip 5.3.0
•All Data Structures Namespaces Files Functions Variables Pages
spffte.f
Go to the documentation of this file.
1C> @file
2C> @brief Perform multiple fast Fourier transforms.
3C>
4C> ### Program History Log
5C> Date | Programmer | Comments
6C> -----|------------|---------
7C> 1998-12-18 | Iredell | Initial.
8C> 2012-11-12 | Mirvis | fixing hard-wired types problem on Intel/Linux.
9C>
10C> @author Iredell @date 96-02-20
11
12C> This subprogram performs multiple fast Fourier transforms
13C> between complex amplitudes in Fourier space and real values
14C> in cyclic physical space.
15C>
16C> This subprogram must be invoked first with IDIR=0
17C> to initialize trigonemetric data. Use subprogram spfft1()
18C> to perform an FFT without previous initialization.
19C>
20C> This version invokes the IBM ESSL FFT.
21C>
22C> @note The restrictions on IMAX are that it must be a multiple
23C> of 1 to 25 factors of two, up to 2 factors of three,
24C> and up to 1 factor of five, seven and eleven.
25C>
26C> If IDIR=0, then W and G need not contain any valid data.
27C> The other parameters must be supplied and cannot change
28C> in succeeding calls until the next time it is called with IDIR=0.
29C>
30C> This subprogram is thread-safe.
31C>
32C> @param IMAX number of values in the cyclic physical space
33C> (see limitations on imax in remarks below.)
34C> @param INCW first dimension of the complex amplitude array
35C> (INCW >= IMAX/2+1)
36C> @param INCG first dimension of the real value array
37C> (INCG >= IMAX)
38C> @param KMAX number of transforms to perform
39C> @param[out] W complex amplitudes if IDIR>0
40C> @param[out] G real values if IDIR<0
41C> @param IDIR direction flag
42C> - IDIR=0 to initialize trigonometric data
43C> - IDIR>0 to transform from Fourier to physical space
44C> - IDIR<0 to transform from physical to Fourier space
45C> @param[out] AFFT auxiliary array if IDIR<>0
46C>
47C> @author Iredell @date 96-02-20
48 SUBROUTINE spffte(IMAX,INCW,INCG,KMAX,W,G,IDIR,AFFT)
49 IMPLICIT NONE
50 INTEGER,INTENT(IN):: IMAX,INCW,INCG,KMAX,IDIR
51 REAL,INTENT(INOUT):: W(2*INCW,KMAX)
52 REAL,INTENT(INOUT):: G(INCG,KMAX)
53 REAL(8),INTENT(INOUT):: AFFT(50000+4*IMAX)
54 REAL:: AFFTR(50000+4*IMAX)
55 INTEGER:: INIT,INC2X,INC2Y,N,M,ISIGN,NAUX1,NAUX2,NAUX3
56C ==EM== ^(4)
57 REAL:: SCALE
58 REAL :: AUX2(20000+2*IMAX),AUX3
59 INTEGER:: IACR,IARC
60
61 naux1=25000+2*imax
62 naux2=20000+2*imax
63 naux3=1
64 iacr=1
65 iarc=1+naux1
66 afftr=real(afft)
67
68C INITIALIZATION.
69C FILL AUXILIARY ARRAYS WITH TRIGONOMETRIC DATA
70 SELECT CASE(idir)
71 CASE(0)
72 init=1
73 inc2x=incw
74 inc2y=incg
75 n=imax
76 m=kmax
77 isign=-1
78 scale=1.
79 IF(digits(1.).LT.digits(1._8)) THEN
80 CALL scrft(init,w,inc2x,g,inc2y,n,m,isign,scale,
81 & afftr(iacr),naux1,aux2,naux2,aux3,naux3)
82 ELSE
83 CALL dcrft(init,w,inc2x,g,inc2y,n,m,isign,scale,
84 & afftr(iacr),naux1,aux2,naux2)
85 ENDIF
86 init=1
87 inc2x=incg
88 inc2y=incw
89 n=imax
90 m=kmax
91 isign=+1
92 scale=1./imax
93 IF(digits(1.).LT.digits(1._8)) THEN
94 CALL srcft(init,g,inc2x,w,inc2y,n,m,isign,scale,
95 & afftr(iarc),naux1,aux2,naux2,aux3,naux3)
96 ELSE
97 CALL drcft(init,g,inc2x,w,inc2y,n,m,isign,scale,
98 & afftr(iarc),naux1,aux2,naux2)
99 ENDIF
100
101C FOURIER TO PHYSICAL TRANSFORM.
102 CASE(1:)
103 init=0
104 inc2x=incw
105 inc2y=incg
106 n=imax
107 m=kmax
108 isign=-1
109 scale=1.
110 IF(digits(1.).LT.digits(1._8)) THEN
111 CALL scrft(init,w,inc2x,g,inc2y,n,m,isign,scale,
112 & afftr(iacr),naux1,aux2,naux2,aux3,naux3)
113 ELSE
114 CALL dcrft(init,w,inc2x,g,inc2y,n,m,isign,scale,
115 & afftr(iacr),naux1,aux2,naux2)
116 ENDIF
117
118C PHYSICAL TO FOURIER TRANSFORM.
119 CASE(:-1)
120 init=0
121 inc2x=incg
122 inc2y=incw
123 n=imax
124 m=kmax
125 isign=+1
126 scale=1./imax
127 IF(digits(1.).LT.digits(1._8)) THEN
128 CALL srcft(init,g,inc2x,w,inc2y,n,m,isign,scale,
129 & afftr(iarc),naux1,aux2,naux2,aux3,naux3)
130 ELSE
131 CALL drcft(init,g,inc2x,w,inc2y,n,m,isign,scale,
132 & afftr(iarc),naux1,aux2,naux2)
133 ENDIF
134 END SELECT
135 afft=real(afftr,kind=8)
136 END SUBROUTINE
drcft
subroutine drcft(init, x, ldx, y, ldy, n, m, isign, scale, table, n1, wrk, n2, z, nz)
Compute a set of m complex discrete n-point Fourier transforms of real data.
Definition fftpack.F:189
dcrft
subroutine dcrft(init, x, ldx, y, ldy, n, m, isign, scale, table, n1, wrk, n2, z, nz)
Computes a set of m real discrete n-point Fourier transforms of complex conjugate even data.
Definition fftpack.F:44
scrft
subroutine scrft(init, x, ldx, y, ldy, n, m, isign, scale, table, n1, wrk, n2, z, nz)
Compute a set of m real discrete n-point Fourier transforms of complex conjugate even data.
Definition fftpack.F:94
srcft
subroutine srcft(init, x, ldx, y, ldy, n, m, isign, scale, table, n1, wrk, n2, z, nz)
Compute a set of m complex discrete n-point Fourier transforms of real data.
Definition fftpack.F:248
spffte
subroutine spffte(imax, incw, incg, kmax, w, g, idir, afft)
This subprogram performs multiple fast Fourier transforms between complex amplitudes in Fourier space...
Definition spffte.f:49