| blob | 3829d650e6de657a1e5159f18533cfce0d7c8d3b |
1 !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
2 !
3 ! FFTPACK 5.0
4 ! Copyright (C) 1995-2004, Scientific Computing Division,
5 ! University Corporation for Atmospheric Research
6 ! Licensed under the GNU General Public License (GPL)
7 !
8 ! Authors: Paul N. Swarztrauber and Richard A. Valent
9 !
10 ! $Id: mradbg.f,v 1.2 2004/06/15 21:29:20 rodney Exp $
11 !
12 !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
14 SUBROUTINE MRADBG (M,IDO,IP,L1,IDL1,CC,C1,C2,IM1,IN1, &
15 & CH,CH2,IM2,IN2,WA)
16 REAL CH(IN2,IDO,L1,IP) ,CC(IN1,IDO,IP,L1) , &
17 & C1(IN1,IDO,L1,IP) ,C2(IN1,IDL1,IP), &
18 & CH2(IN2,IDL1,IP) ,WA(IDO)
19 !
20 M1D = (M-1)*IM1+1
21 M2S = 1-IM2
22 TPI=2.*4.*ATAN(1.0)
23 ARG = TPI/FLOAT(IP)
24 DCP = COS(ARG)
25 DSP = SIN(ARG)
26 IDP2 = IDO+2
27 NBD = (IDO-1)/2
28 IPP2 = IP+2
29 IPPH = (IP+1)/2
30 IF (IDO .LT. L1) GO TO 103
31 DO 102 K=1,L1
32 DO 101 I=1,IDO
33 M2 = M2S
34 DO 1001 M1=1,M1D,IM1
35 M2 = M2+IM2
36 CH(M2,I,K,1) = CC(M1,I,1,K)
37 1001 CONTINUE
38 101 CONTINUE
39 102 END DO
40 GO TO 106
41 103 DO 105 I=1,IDO
42 DO 104 K=1,L1
43 M2 = M2S
44 DO 1004 M1=1,M1D,IM1
45 M2 = M2+IM2
46 CH(M2,I,K,1) = CC(M1,I,1,K)
47 1004 CONTINUE
48 104 CONTINUE
49 105 END DO
50 106 DO 108 J=2,IPPH
51 JC = IPP2-J
52 J2 = J+J
53 DO 107 K=1,L1
54 M2 = M2S
55 DO 1007 M1=1,M1D,IM1
56 M2 = M2+IM2
57 CH(M2,1,K,J) = CC(M1,IDO,J2-2,K)+CC(M1,IDO,J2-2,K)
58 CH(M2,1,K,JC) = CC(M1,1,J2-1,K)+CC(M1,1,J2-1,K)
59 1007 CONTINUE
60 107 CONTINUE
61 108 END DO
62 IF (IDO .EQ. 1) GO TO 116
63 IF (NBD .LT. L1) GO TO 112
64 DO 111 J=2,IPPH
65 JC = IPP2-J
66 DO 110 K=1,L1
67 DO 109 I=3,IDO,2
68 IC = IDP2-I
69 M2 = M2S
70 DO 1009 M1=1,M1D,IM1
71 M2 = M2+IM2
72 CH(M2,I-1,K,J) = CC(M1,I-1,2*J-1,K)+CC(M1,IC-1,2*J-2,K)
73 CH(M2,I-1,K,JC) = CC(M1,I-1,2*J-1,K)-CC(M1,IC-1,2*J-2,K)
74 CH(M2,I,K,J) = CC(M1,I,2*J-1,K)-CC(M1,IC,2*J-2,K)
75 CH(M2,I,K,JC) = CC(M1,I,2*J-1,K)+CC(M1,IC,2*J-2,K)
76 1009 CONTINUE
77 109 CONTINUE
78 110 CONTINUE
79 111 END DO
80 GO TO 116
81 112 DO 115 J=2,IPPH
82 JC = IPP2-J
83 DO 114 I=3,IDO,2
84 IC = IDP2-I
85 DO 113 K=1,L1
86 M2 = M2S
87 DO 1013 M1=1,M1D,IM1
88 M2 = M2+IM2
89 CH(M2,I-1,K,J) = CC(M1,I-1,2*J-1,K)+CC(M1,IC-1,2*J-2,K)
90 CH(M2,I-1,K,JC) = CC(M1,I-1,2*J-1,K)-CC(M1,IC-1,2*J-2,K)
91 CH(M2,I,K,J) = CC(M1,I,2*J-1,K)-CC(M1,IC,2*J-2,K)
92 CH(M2,I,K,JC) = CC(M1,I,2*J-1,K)+CC(M1,IC,2*J-2,K)
93 1013 CONTINUE
94 113 CONTINUE
95 114 CONTINUE
96 115 END DO
97 116 AR1 = 1.
98 AI1 = 0.
99 DO 120 L=2,IPPH
100 LC = IPP2-L
101 AR1H = DCP*AR1-DSP*AI1
102 AI1 = DCP*AI1+DSP*AR1
103 AR1 = AR1H
104 DO 117 IK=1,IDL1
105 M2 = M2S
106 DO 1017 M1=1,M1D,IM1
107 M2 = M2+IM2
108 C2(M1,IK,L) = CH2(M2,IK,1)+AR1*CH2(M2,IK,2)
109 C2(M1,IK,LC) = AI1*CH2(M2,IK,IP)
110 1017 CONTINUE
111 117 CONTINUE
112 DC2 = AR1
113 DS2 = AI1
114 AR2 = AR1
115 AI2 = AI1
116 DO 119 J=3,IPPH
117 JC = IPP2-J
118 AR2H = DC2*AR2-DS2*AI2
119 AI2 = DC2*AI2+DS2*AR2
120 AR2 = AR2H
121 DO 118 IK=1,IDL1
122 M2 = M2S
123 DO 1018 M1=1,M1D,IM1
124 M2 = M2+IM2
125 C2(M1,IK,L) = C2(M1,IK,L)+AR2*CH2(M2,IK,J)
126 C2(M1,IK,LC) = C2(M1,IK,LC)+AI2*CH2(M2,IK,JC)
127 1018 CONTINUE
128 118 CONTINUE
129 119 CONTINUE
130 120 END DO
131 DO 122 J=2,IPPH
132 DO 121 IK=1,IDL1
133 M2 = M2S
134 DO 1021 M1=1,M1D,IM1
135 M2 = M2+IM2
136 CH2(M2,IK,1) = CH2(M2,IK,1)+CH2(M2,IK,J)
137 1021 CONTINUE
138 121 CONTINUE
139 122 END DO
140 DO 124 J=2,IPPH
141 JC = IPP2-J
142 DO 123 K=1,L1
143 M2 = M2S
144 DO 1023 M1=1,M1D,IM1
145 M2 = M2+IM2
146 CH(M2,1,K,J) = C1(M1,1,K,J)-C1(M1,1,K,JC)
147 CH(M2,1,K,JC) = C1(M1,1,K,J)+C1(M1,1,K,JC)
148 1023 CONTINUE
149 123 CONTINUE
150 124 END DO
151 IF (IDO .EQ. 1) GO TO 132
152 IF (NBD .LT. L1) GO TO 128
153 DO 127 J=2,IPPH
154 JC = IPP2-J
155 DO 126 K=1,L1
156 DO 125 I=3,IDO,2
157 M2 = M2S
158 DO 1025 M1=1,M1D,IM1
159 M2 = M2+IM2
160 CH(M2,I-1,K,J) = C1(M1,I-1,K,J)-C1(M1,I,K,JC)
161 CH(M2,I-1,K,JC) = C1(M1,I-1,K,J)+C1(M1,I,K,JC)
162 CH(M2,I,K,J) = C1(M1,I,K,J)+C1(M1,I-1,K,JC)
163 CH(M2,I,K,JC) = C1(M1,I,K,J)-C1(M1,I-1,K,JC)
164 1025 CONTINUE
165 125 CONTINUE
166 126 CONTINUE
167 127 END DO
168 GO TO 132
169 128 DO 131 J=2,IPPH
170 JC = IPP2-J
171 DO 130 I=3,IDO,2
172 DO 129 K=1,L1
173 M2 = M2S
174 DO 1029 M1=1,M1D,IM1
175 M2 = M2+IM2
176 CH(M2,I-1,K,J) = C1(M1,I-1,K,J)-C1(M1,I,K,JC)
177 CH(M2,I-1,K,JC) = C1(M1,I-1,K,J)+C1(M1,I,K,JC)
178 CH(M2,I,K,J) = C1(M1,I,K,J)+C1(M1,I-1,K,JC)
179 CH(M2,I,K,JC) = C1(M1,I,K,J)-C1(M1,I-1,K,JC)
180 1029 CONTINUE
181 129 CONTINUE
182 130 CONTINUE
183 131 END DO
184 132 CONTINUE
185 IF (IDO .EQ. 1) RETURN
186 DO 133 IK=1,IDL1
187 M2 = M2S
188 DO 1033 M1=1,M1D,IM1
189 M2 = M2+IM2
190 C2(M1,IK,1) = CH2(M2,IK,1)
191 1033 CONTINUE
192 133 END DO
193 DO 135 J=2,IP
194 DO 134 K=1,L1
195 M2 = M2S
196 DO 1034 M1=1,M1D,IM1
197 M2 = M2+IM2
198 C1(M1,1,K,J) = CH(M2,1,K,J)
199 1034 CONTINUE
200 134 CONTINUE
201 135 END DO
202 IF (NBD .GT. L1) GO TO 139
203 IS = -IDO
204 DO 138 J=2,IP
205 IS = IS+IDO
206 IDIJ = IS
207 DO 137 I=3,IDO,2
208 IDIJ = IDIJ+2
209 DO 136 K=1,L1
210 M2 = M2S
211 DO 1036 M1=1,M1D,IM1
212 M2 = M2+IM2
213 C1(M1,I-1,K,J) = WA(IDIJ-1)*CH(M2,I-1,K,J)-WA(IDIJ)* &
214 & CH(M2,I,K,J)
215 C1(M1,I,K,J) = WA(IDIJ-1)*CH(M2,I,K,J)+WA(IDIJ)* &
216 & CH(M2,I-1,K,J)
217 1036 CONTINUE
218 136 CONTINUE
219 137 CONTINUE
220 138 END DO
221 GO TO 143
222 139 IS = -IDO
223 DO 142 J=2,IP
224 IS = IS+IDO
225 DO 141 K=1,L1
226 IDIJ = IS
227 DO 140 I=3,IDO,2
228 IDIJ = IDIJ+2
229 M2 = M2S
230 DO 1040 M1=1,M1D,IM1
231 M2 = M2+IM2
232 C1(M1,I-1,K,J) = WA(IDIJ-1)*CH(M2,I-1,K,J)-WA(IDIJ)* &
233 & CH(M2,I,K,J)
234 C1(M1,I,K,J) = WA(IDIJ-1)*CH(M2,I,K,J)+WA(IDIJ)* &
235 & CH(M2,I-1,K,J)
236 1040 CONTINUE
237 140 CONTINUE
238 141 CONTINUE
239 142 END DO
240 143 RETURN
241 END