PLaSK library
Loading...
Searching...
No Matches
rfft2b.c
Go to the documentation of this file.
1/* rfft2b.f -- translated by f2c (version 20100827).
2 You must link the resulting object file with libf2c:
3 on Microsoft Windows system, link with libf2c.lib;
4 on Linux or Unix systems, link with .../path/to/libf2c.a -lm
5 or, if you install libf2c.a in a standard place, with -lf2c -lm
6 -- in that order, at the end of the command line, as in
7 cc *.o -lf2c -lm
8 Source for libf2c is in /netlib/f2c/libf2c.zip, e.g.,
9
10 http://www.netlib.org/f2c/libf2c.zip
11*/
12
13#include "f2c.h"
14
15/* Table of constant values */
16
17static integer c__6 = 6;
18static integer c__8 = 8;
19static integer c_n6 = -6;
20static integer c__1 = 1;
21static integer c_n5 = -5;
22
23/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
24/* * * */
25/* * copyright (c) 2011 by UCAR * */
26/* * * */
27/* * University Corporation for Atmospheric Research * */
28/* * * */
29/* * all rights reserved * */
30/* * * */
31/* * FFTPACK version 5.1 * */
32/* * * */
33/* * A Fortran Package of Fast Fourier * */
34/* * * */
35/* * Subroutines and Example Programs * */
36/* * * */
37/* * by * */
38/* * * */
39/* * Paul Swarztrauber and Dick Valent * */
40/* * * */
41/* * of * */
42/* * * */
43/* * the National Center for Atmospheric Research * */
44/* * * */
45/* * Boulder, Colorado (80307) U.S.A. * */
46/* * * */
47/* * which is sponsored by * */
48/* * * */
49/* * the National Science Foundation * */
50/* * * */
51/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
52
53/* Subroutine */ int rfft2b_(integer *ldim, integer *l, integer *m,
54 doublereal *r__, doublereal *wsave, integer *lensav, doublereal *work,
55 integer *lenwrk, integer *ier)
56{
57 /* System generated locals */
58 integer r_dim1, r_offset, i__1, i__2, i__3;
59
60 /* Builtin functions */
61 double log(doublereal);
62
63 /* Local variables */
64 integer i__, j;
65 extern /* Subroutine */ int r2w_(integer *, integer *, integer *, integer
66 *, doublereal *, doublereal *), w2r_(integer *, integer *,
67 integer *, integer *, doublereal *, doublereal *);
68 integer ldh, ldw, ldx, ier1, modl, modm, mmsav, lwsav, mwsav;
69 extern /* Subroutine */ int cfftmb_(integer *, integer *, integer *,
70 integer *, doublereal *, integer *, doublereal *, integer *,
71 doublereal *, integer *, integer *), rfftmb_(integer *, integer *,
72 integer *, integer *, doublereal *, integer *, doublereal *,
73 integer *, doublereal *, integer *, integer *), xerfft_(char *,
74 integer *, ftnlen);
75
76
77
78/* INITIALIZE IER */
79
80 /* Parameter adjustments */
81 r_dim1 = *ldim;
82 r_offset = 1 + r_dim1;
83 r__ -= r_offset;
84 --wsave;
85 --work;
86
87 /* Function Body */
88 *ier = 0;
89
90/* VERIFY LENSAV */
91
92 lwsav = *l + (integer) (log((doublereal) (*l)) / log(2.)) + 4;
93 mwsav = (*m << 1) + (integer) (log((doublereal) (*m)) / log(2.)) + 4;
94 mmsav = *m + (integer) (log((doublereal) (*m)) / log(2.)) + 4;
95 modl = *l % 2;
96 modm = *m % 2;
97
98 if (*lensav < lwsav + mwsav + mmsav) {
99 *ier = 2;
100 xerfft_("RFFT2F", &c__6, (ftnlen)6);
101 goto L100;
102 }
103
104/* VERIFY LENWRK */
105
106 if (*lenwrk < (*l + 1) * *m) {
107 *ier = 3;
108 xerfft_("RFFT2F", &c__8, (ftnlen)6);
109 goto L100;
110 }
111
112/* VERIFY LDIM IS AS BIG AS L */
113
114 if (*ldim < *l) {
115 *ier = 5;
116 xerfft_("RFFT2F", &c_n6, (ftnlen)6);
117 goto L100;
118 }
119
120/* TRANSFORM SECOND DIMENSION OF ARRAY */
121
122 i__1 = ((*m + 1) / 2 << 1) - 1;
123 for (j = 2; j <= i__1; ++j) {
124 r__[j * r_dim1 + 1] += r__[j * r_dim1 + 1];
125 }
126 i__1 = *m;
127 for (j = 3; j <= i__1; j += 2) {
128 r__[j * r_dim1 + 1] = -r__[j * r_dim1 + 1];
129 }
130 i__1 = *m * *ldim;
131 rfftmb_(&c__1, &c__1, m, ldim, &r__[r_offset], &i__1, &wsave[lwsav +
132 mwsav + 1], &mmsav, &work[1], lenwrk, &ier1);
133 ldh = (*l + 1) / 2;
134 if (ldh > 1) {
135 ldw = ldh + ldh;
136
137/* R AND WORK ARE SWITCHED BECAUSE THE THE FIRST DIMENSION */
138/* OF THE INPUT TO COMPLEX CFFTMF MUST BE EVEN. */
139
140 r2w_(ldim, &ldw, l, m, &r__[r_offset], &work[1]);
141 i__1 = ldh - 1;
142 i__2 = ldh * *m;
143 i__3 = *l * *m;
144 cfftmb_(&i__1, &c__1, m, &ldh, &work[2], &i__2, &wsave[lwsav + 1], &
145 mwsav, &r__[r_offset], &i__3, &ier1);
146 if (ier1 != 0) {
147 *ier = 20;
148 xerfft_("RFFT2B", &c_n5, (ftnlen)6);
149 goto L100;
150 }
151 w2r_(ldim, &ldw, l, m, &r__[r_offset], &work[1]);
152 }
153
154 if (modl == 0) {
155 i__1 = ((*m + 1) / 2 << 1) - 1;
156 for (j = 2; j <= i__1; ++j) {
157 r__[*l + j * r_dim1] += r__[*l + j * r_dim1];
158 }
159 i__1 = *m;
160 for (j = 3; j <= i__1; j += 2) {
161 r__[*l + j * r_dim1] = -r__[*l + j * r_dim1];
162 }
163 i__1 = *m * *ldim;
164 rfftmb_(&c__1, &c__1, m, ldim, &r__[*l + r_dim1], &i__1, &wsave[lwsav
165 + mwsav + 1], &mmsav, &work[1], lenwrk, &ier1);
166 }
167
168/* PRINT*, 'BACKWARD TRANSFORM IN THE J DIRECTION' */
169/* DO I=1,L */
170/* PRINT*, (R(I,J),J=1,M) */
171/* END DO */
172
173/* TRANSFORM FIRST DIMENSION OF ARRAY */
174
175 ldx = ((*l + 1) / 2 << 1) - 1;
176 i__1 = ldx;
177 for (i__ = 2; i__ <= i__1; ++i__) {
178 i__2 = *m;
179 for (j = 1; j <= i__2; ++j) {
180 r__[i__ + j * r_dim1] += r__[i__ + j * r_dim1];
181 }
182 }
183 i__1 = *m;
184 for (j = 1; j <= i__1; ++j) {
185 i__2 = ldx;
186 for (i__ = 3; i__ <= i__2; i__ += 2) {
187 r__[i__ + j * r_dim1] = -r__[i__ + j * r_dim1];
188 }
189 }
190 i__1 = *m * *ldim;
191 i__2 = *l + (integer) (log((doublereal) (*l)) / log(2.)) + 4;
192 rfftmb_(m, ldim, l, &c__1, &r__[r_offset], &i__1, &wsave[1], &i__2, &work[
193 1], lenwrk, &ier1);
194
195
196/* PRINT*, 'BACKWARD TRANSFORM IN THE I DIRECTION' */
197/* DO I=1,L */
198/* PRINT*, (R(I,J),J=1,M) */
199/* END DO */
200
201 if (ier1 != 0) {
202 *ier = 20;
203 xerfft_("RFFT2F", &c_n5, (ftnlen)6);
204 goto L100;
205 }
206
207 if (ier1 != 0) {
208 *ier = 20;
209 xerfft_("RFFT2F", &c_n5, (ftnlen)6);
210 goto L100;
211 }
212
213L100:
214
215 return 0;
216} /* rfft2b_ */
217