NAMD
Macros | Functions
pub3dfft.C File Reference
#include <math.h>
#include "pub3dfft.h"

Go to the source code of this file.

Macros

#define M_PI   3.14159265358979323846
 

Functions

int passb (int *nac, int *ido, int *ip, int *l1, int *idl1, double *cc, double *c1, double *c2, double *ch, double *ch2, double *wa)
 
int passb2 (int *ido, int *l1, double *cc, double *ch, double *wa1)
 
int passb3 (int *ido, int *l1, double *cc, double *ch, double *wa1, double *wa2)
 
int passb4 (int *ido, int *l1, double *cc, double *ch, double *wa1, double *wa2, double *wa3)
 
int passb5 (int *ido, int *l1, double *cc, double *ch, double *wa1, double *wa2, double *wa3, double *wa4)
 
int passf (int *nac, int *ido, int *ip, int *l1, int *idl1, double *cc, double *c1, double *c2, double *ch, double *ch2, double *wa)
 
int passf2 (int *ido, int *l1, double *cc, double *ch, double *wa1)
 
int passf3 (int *ido, int *l1, double *cc, double *ch, double *wa1, double *wa2)
 
int passf4 (int *ido, int *l1, double *cc, double *ch, double *wa1, double *wa2, double *wa3)
 
int passf5 (int *ido, int *l1, double *cc, double *ch, double *wa1, double *wa2, double *wa3, double *wa4)
 
int cfftb1 (int *n, double *c, double *ch, double *wa, int *ifac)
 
int cfftb (int *n, double *c, double *wsave)
 
int cfftf1 (int *n, double *c, double *ch, double *wa, int *ifac)
 
int cfftf (int *n, double *c, double *wsave)
 
int cffti1 (int *n, double *wa, int *ifac)
 
int cffti (int *n, double *wsave)
 
int pubz3di (int *n1, int *n2, int *n3, double *table, int *ntable)
 
int pubz3d (int *isign, int *n1, int *n2, int *n3, doublecomplex *w, int *ld1, int *ld2, double *table, int *ntable, doublecomplex *work)
 
int pubd3di (int n1, int n2, int n3, double *table, int ntable)
 
int pubdz3d (int isign, int n1, int n2, int n3, double *w, int ld1, int ld2, double *table, int ntable, double *work)
 
int pubzd3d (int isign, int n1, int n2, int n3, double *w, int ld1, int ld2, double *table, int ntable, double *work)
 

Macro Definition Documentation

◆ M_PI

#define M_PI   3.14159265358979323846

Copyright (c) 1995, 1996, 1997, 1998, 1999, 2000 by The Board of Trustees of the University of Illinois. All rights reserved.

Definition at line 15 of file pub3dfft.C.

Referenced by pubdz3d(), and pubzd3d().

Function Documentation

◆ cfftb()

int cfftb ( int *  n,
double *  c,
double *  wsave 
)

Definition at line 1401 of file pub3dfft.C.

References cfftb1().

Referenced by pubz3d().

1402 {
1403  static int iw1, iw2;
1404 
1405  /* Parameter adjustments */
1406  --c;
1407  --wsave;
1408 
1409  /* Function Body */
1410  if (*n == 1) {
1411  return 0;
1412  }
1413  iw1 = *n + *n + 1;
1414  iw2 = iw1 + *n + *n;
1415  cfftb1(n, &c[1], &wsave[1], &wsave[iw1], (int*) &wsave[iw2]);
1416  return 0;
1417 } /* cfftb_ */
cfftb1
int cfftb1(int *n, double *c, double *ch, double *wa, int *ifac)
Definition: pub3dfft.C:1277

◆ cfftb1()

int cfftb1 ( int *  n,
double *  c,
double *  ch,
double *  wa,
int *  ifac 
)

Definition at line 1277 of file pub3dfft.C.

References passb(), passb2(), passb3(), passb4(), and passb5().

Referenced by cfftb().

1279 {
1280  /* System generated locals */
1281  int i_1;
1282 
1283  /* Local variables */
1284  static int idot, i;
1285  static int k1, l1, l2, n2;
1286  static int na, nf, ip, iw, ix2, ix3, ix4, nac, ido, idl1;
1287 
1288  /* Parameter adjustments */
1289  --c;
1290  --ch;
1291  --wa;
1292  --ifac;
1293 
1294  /* Function Body */
1295  nf = ifac[2];
1296  na = 0;
1297  l1 = 1;
1298  iw = 1;
1299  i_1 = nf;
1300  for (k1 = 1; k1 <= i_1; ++k1) {
1301  ip = ifac[k1 + 2];
1302  l2 = ip * l1;
1303  ido = *n / l2;
1304  idot = ido + ido;
1305  idl1 = idot * l1;
1306  if (ip != 4) {
1307  goto L103;
1308  }
1309  ix2 = iw + idot;
1310  ix3 = ix2 + idot;
1311  if (na != 0) {
1312  goto L101;
1313  }
1314  passb4(&idot, &l1, &c[1], &ch[1], &wa[iw], &wa[ix2], &wa[ix3]);
1315  goto L102;
1316 L101:
1317  passb4(&idot, &l1, &ch[1], &c[1], &wa[iw], &wa[ix2], &wa[ix3]);
1318 L102:
1319  na = 1 - na;
1320  goto L115;
1321 L103:
1322  if (ip != 2) {
1323  goto L106;
1324  }
1325  if (na != 0) {
1326  goto L104;
1327  }
1328  passb2(&idot, &l1, &c[1], &ch[1], &wa[iw]);
1329  goto L105;
1330 L104:
1331  passb2(&idot, &l1, &ch[1], &c[1], &wa[iw]);
1332 L105:
1333  na = 1 - na;
1334  goto L115;
1335 L106:
1336  if (ip != 3) {
1337  goto L109;
1338  }
1339  ix2 = iw + idot;
1340  if (na != 0) {
1341  goto L107;
1342  }
1343  passb3(&idot, &l1, &c[1], &ch[1], &wa[iw], &wa[ix2]);
1344  goto L108;
1345 L107:
1346  passb3(&idot, &l1, &ch[1], &c[1], &wa[iw], &wa[ix2]);
1347 L108:
1348  na = 1 - na;
1349  goto L115;
1350 L109:
1351  if (ip != 5) {
1352  goto L112;
1353  }
1354  ix2 = iw + idot;
1355  ix3 = ix2 + idot;
1356  ix4 = ix3 + idot;
1357  if (na != 0) {
1358  goto L110;
1359  }
1360  passb5(&idot, &l1, &c[1], &ch[1], &wa[iw], &wa[ix2], &wa[ix3], &wa[
1361  ix4]);
1362  goto L111;
1363 L110:
1364  passb5(&idot, &l1, &ch[1], &c[1], &wa[iw], &wa[ix2], &wa[ix3], &wa[
1365  ix4]);
1366 L111:
1367  na = 1 - na;
1368  goto L115;
1369 L112:
1370  if (na != 0) {
1371  goto L113;
1372  }
1373  passb(&nac, &idot, &ip, &l1, &idl1, &c[1], &c[1], &c[1], &ch[1], &ch[
1374  1], &wa[iw]);
1375  goto L114;
1376 L113:
1377  passb(&nac, &idot, &ip, &l1, &idl1, &ch[1], &ch[1], &ch[1], &c[1], &
1378  c[1], &wa[iw]);
1379 L114:
1380  if (nac != 0) {
1381  na = 1 - na;
1382  }
1383 L115:
1384  l1 = l2;
1385  iw += (ip - 1) * idot;
1386 /* L116: */
1387  }
1388  if (na == 0) {
1389  return 0;
1390  }
1391  n2 = *n + *n;
1392  i_1 = n2;
1393  for (i = 1; i <= i_1; ++i) {
1394  c[i] = ch[i];
1395 /* L117: */
1396  }
1397  return 0;
1398 } /* cfftb1_ */
passb4
int passb4(int *ido, int *l1, double *cc, double *ch, double *wa1, double *wa2, double *wa3)
Definition: pub3dfft.C:408
passb2
int passb2(int *ido, int *l1, double *cc, double *ch, double *wa1)
Definition: pub3dfft.C:257
passb3
int passb3(int *ido, int *l1, double *cc, double *ch, double *wa1, double *wa2)
Definition: pub3dfft.C:319
passb5
int passb5(int *ido, int *l1, double *cc, double *ch, double *wa1, double *wa2, double *wa3, double *wa4)
Definition: pub3dfft.C:512
passb
int passb(int *nac, int *ido, int *ip, int *l1, int *idl1, double *cc, double *c1, double *c2, double *ch, double *ch2, double *wa)
Definition: pub3dfft.C:18

◆ cfftf()

int cfftf ( int *  n,
double *  c,
double *  wsave 
)

Definition at line 1544 of file pub3dfft.C.

References cfftf1().

Referenced by pubz3d().

1545 {
1546  static int iw1, iw2;
1547 
1548  /* Parameter adjustments */
1549  --c;
1550  --wsave;
1551 
1552  /* Function Body */
1553  if (*n == 1) {
1554  return 0;
1555  }
1556  iw1 = *n + *n + 1;
1557  iw2 = iw1 + *n + *n;
1558  cfftf1(n, &c[1], &wsave[1], &wsave[iw1], (int*) &wsave[iw2]);
1559  return 0;
1560 } /* cfftf_ */
cfftf1
int cfftf1(int *n, double *c, double *ch, double *wa, int *ifac)
Definition: pub3dfft.C:1420

◆ cfftf1()

int cfftf1 ( int *  n,
double *  c,
double *  ch,
double *  wa,
int *  ifac 
)

Definition at line 1420 of file pub3dfft.C.

References passf(), passf2(), passf3(), passf4(), and passf5().

Referenced by cfftf().

1422 {
1423  /* System generated locals */
1424  int i_1;
1425 
1426  /* Local variables */
1427  static int idot, i;
1428  static int k1, l1, l2, n2;
1429  static int na, nf, ip, iw, ix2, ix3, ix4, nac, ido, idl1;
1430 
1431  /* Parameter adjustments */
1432  --c;
1433  --ch;
1434  --wa;
1435  --ifac;
1436 
1437  /* Function Body */
1438  nf = ifac[2];
1439  na = 0;
1440  l1 = 1;
1441  iw = 1;
1442  i_1 = nf;
1443  for (k1 = 1; k1 <= i_1; ++k1) {
1444  ip = ifac[k1 + 2];
1445  l2 = ip * l1;
1446  ido = *n / l2;
1447  idot = ido + ido;
1448  idl1 = idot * l1;
1449  if (ip != 4) {
1450  goto L103;
1451  }
1452  ix2 = iw + idot;
1453  ix3 = ix2 + idot;
1454  if (na != 0) {
1455  goto L101;
1456  }
1457  passf4(&idot, &l1, &c[1], &ch[1], &wa[iw], &wa[ix2], &wa[ix3]);
1458  goto L102;
1459 L101:
1460  passf4(&idot, &l1, &ch[1], &c[1], &wa[iw], &wa[ix2], &wa[ix3]);
1461 L102:
1462  na = 1 - na;
1463  goto L115;
1464 L103:
1465  if (ip != 2) {
1466  goto L106;
1467  }
1468  if (na != 0) {
1469  goto L104;
1470  }
1471  passf2(&idot, &l1, &c[1], &ch[1], &wa[iw]);
1472  goto L105;
1473 L104:
1474  passf2(&idot, &l1, &ch[1], &c[1], &wa[iw]);
1475 L105:
1476  na = 1 - na;
1477  goto L115;
1478 L106:
1479  if (ip != 3) {
1480  goto L109;
1481  }
1482  ix2 = iw + idot;
1483  if (na != 0) {
1484  goto L107;
1485  }
1486  passf3(&idot, &l1, &c[1], &ch[1], &wa[iw], &wa[ix2]);
1487  goto L108;
1488 L107:
1489  passf3(&idot, &l1, &ch[1], &c[1], &wa[iw], &wa[ix2]);
1490 L108:
1491  na = 1 - na;
1492  goto L115;
1493 L109:
1494  if (ip != 5) {
1495  goto L112;
1496  }
1497  ix2 = iw + idot;
1498  ix3 = ix2 + idot;
1499  ix4 = ix3 + idot;
1500  if (na != 0) {
1501  goto L110;
1502  }
1503  passf5(&idot, &l1, &c[1], &ch[1], &wa[iw], &wa[ix2], &wa[ix3], &wa[
1504  ix4]);
1505  goto L111;
1506 L110:
1507  passf5(&idot, &l1, &ch[1], &c[1], &wa[iw], &wa[ix2], &wa[ix3], &wa[
1508  ix4]);
1509 L111:
1510  na = 1 - na;
1511  goto L115;
1512 L112:
1513  if (na != 0) {
1514  goto L113;
1515  }
1516  passf(&nac, &idot, &ip, &l1, &idl1, &c[1], &c[1], &c[1], &ch[1], &ch[
1517  1], &wa[iw]);
1518  goto L114;
1519 L113:
1520  passf(&nac, &idot, &ip, &l1, &idl1, &ch[1], &ch[1], &ch[1], &c[1], &
1521  c[1], &wa[iw]);
1522 L114:
1523  if (nac != 0) {
1524  na = 1 - na;
1525  }
1526 L115:
1527  l1 = l2;
1528  iw += (ip - 1) * idot;
1529 /* L116: */
1530  }
1531  if (na == 0) {
1532  return 0;
1533  }
1534  n2 = *n + *n;
1535  i_1 = n2;
1536  for (i = 1; i <= i_1; ++i) {
1537  c[i] = ch[i];
1538 /* L117: */
1539  }
1540  return 0;
1541 } /* cfftf1_ */
passf5
int passf5(int *ido, int *l1, double *cc, double *ch, double *wa1, double *wa2, double *wa3, double *wa4)
Definition: pub3dfft.C:1139
passf3
int passf3(int *ido, int *l1, double *cc, double *ch, double *wa1, double *wa2)
Definition: pub3dfft.C:946
passf
int passf(int *nac, int *ido, int *ip, int *l1, int *idl1, double *cc, double *c1, double *c2, double *ch, double *ch2, double *wa)
Definition: pub3dfft.C:647
passf4
int passf4(int *ido, int *l1, double *cc, double *ch, double *wa1, double *wa2, double *wa3)
Definition: pub3dfft.C:1035
passf2
int passf2(int *ido, int *l1, double *cc, double *ch, double *wa1)
Definition: pub3dfft.C:885

◆ cffti()

int cffti ( int *  n,
double *  wsave 
)

Definition at line 1676 of file pub3dfft.C.

References cffti1().

Referenced by pubz3di().

1677 {
1678  static int iw1, iw2;
1679 
1680  /* Parameter adjustments */
1681  --wsave;
1682 
1683  /* Function Body */
1684  if (*n == 1) {
1685  return 0;
1686  }
1687  iw1 = *n + *n + 1;
1688  iw2 = iw1 + *n + *n;
1689  cffti1(n, &wsave[iw1], (int*) &wsave[iw2]);
1690  return 0;
1691 } /* cffti_ */
cffti1
int cffti1(int *n, double *wa, int *ifac)
Definition: pub3dfft.C:1563

◆ cffti1()

int cffti1 ( int *  n,
double *  wa,
int *  ifac 
)

Definition at line 1563 of file pub3dfft.C.

Referenced by cffti().

1564 {
1565  /* Initialized data */
1566 
1567  static int ntryh[4] = { 3,4,2,5 };
1568 
1569  /* System generated locals */
1570  int i_1, i_2, i_3;
1571 
1572  /* Local variables */
1573  static double argh;
1574  static int idot, ntry, i, j;
1575  static double argld;
1576  static int i1, k1, l1, l2, ib;
1577  static double fi;
1578  static int ld, ii, nf, ip, nl, nq, nr;
1579  static double arg;
1580  static int ido, ipm;
1581  static double tpi;
1582 
1583  /* Parameter adjustments */
1584  --wa;
1585  --ifac;
1586 
1587  /* Function Body */
1588  nl = *n;
1589  nf = 0;
1590  j = 0;
1591 L101:
1592  ++j;
1593  if (j - 4 <= 0) {
1594  goto L102;
1595  } else {
1596  goto L103;
1597  }
1598 L102:
1599  ntry = ntryh[j - 1];
1600  goto L104;
1601 L103:
1602  ntry += 2;
1603 L104:
1604  nq = nl / ntry;
1605  nr = nl - ntry * nq;
1606  if (nr != 0) {
1607  goto L101;
1608  } else {
1609  goto L105;
1610  }
1611 L105:
1612  ++nf;
1613  ifac[nf + 2] = ntry;
1614  nl = nq;
1615  if (ntry != 2) {
1616  goto L107;
1617  }
1618  if (nf == 1) {
1619  goto L107;
1620  }
1621  i_1 = nf;
1622  for (i = 2; i <= i_1; ++i) {
1623  ib = nf - i + 2;
1624  ifac[ib + 2] = ifac[ib + 1];
1625 /* L106: */
1626  }
1627  ifac[3] = 2;
1628 L107:
1629  if (nl != 1) {
1630  goto L104;
1631  }
1632  ifac[1] = *n;
1633  ifac[2] = nf;
1634  tpi = 6.28318530717959;
1635  argh = tpi / (double) (*n);
1636  i = 2;
1637  l1 = 1;
1638  i_1 = nf;
1639  for (k1 = 1; k1 <= i_1; ++k1) {
1640  ip = ifac[k1 + 2];
1641  ld = 0;
1642  l2 = l1 * ip;
1643  ido = *n / l2;
1644  idot = ido + ido + 2;
1645  ipm = ip - 1;
1646  i_2 = ipm;
1647  for (j = 1; j <= i_2; ++j) {
1648  i1 = i;
1649  wa[i - 1] = 1.;
1650  wa[i] = 0.;
1651  ld += l1;
1652  fi = 0.;
1653  argld = (double) ld * argh;
1654  i_3 = idot;
1655  for (ii = 4; ii <= i_3; ii += 2) {
1656  i += 2;
1657  fi += 1.;
1658  arg = fi * argld;
1659  wa[i - 1] = cos(arg);
1660  wa[i] = sin(arg);
1661 /* L108: */
1662  }
1663  if (ip <= 5) {
1664  goto L109;
1665  }
1666  wa[i1 - 1] = wa[i - 1];
1667  wa[i1] = wa[i];
1668 L109:
1669  ;}
1670  l1 = l2;
1671 /* L110: */
1672  }
1673  return 0;
1674 } /* cffti1_ */

◆ passb()

int passb ( int *  nac,
int *  ido,
int *  ip,
int *  l1,
int *  idl1,
double *  cc,
double *  c1,
double *  c2,
double *  ch,
double *  ch2,
double *  wa 
)

Definition at line 18 of file pub3dfft.C.

Referenced by cfftb1().

21 {
22  /* System generated locals */
23  int ch_dim1, ch_dim2, ch_offset, cc_dim1, cc_dim2, cc_offset, c1_dim1,
24  c1_dim2, c1_offset, c2_dim1, c2_offset, ch2_dim1, ch2_offset,
25  i_1, i_2, i_3;
26 
27  /* Local variables */
28  static int idij, idlj, idot, ipph, i, j, k, l, jc, lc, ik, nt, idj,
29  idl, inc, idp;
30  static double wai, war;
31  static int ipp2;
32 
33  /* Parameter adjustments */
34  cc_dim1 = *ido;
35  cc_dim2 = *ip;
36  cc_offset = cc_dim1 * (cc_dim2 + 1) + 1;
37  cc -= cc_offset;
38  c1_dim1 = *ido;
39  c1_dim2 = *l1;
40  c1_offset = c1_dim1 * (c1_dim2 + 1) + 1;
41  c1 -= c1_offset;
42  c2_dim1 = *idl1;
43  c2_offset = c2_dim1 + 1;
44  c2 -= c2_offset;
45  ch_dim1 = *ido;
46  ch_dim2 = *l1;
47  ch_offset = ch_dim1 * (ch_dim2 + 1) + 1;
48  ch -= ch_offset;
49  ch2_dim1 = *idl1;
50  ch2_offset = ch2_dim1 + 1;
51  ch2 -= ch2_offset;
52  --wa;
53 
54  /* Function Body */
55  idot = *ido / 2;
56  nt = *ip * *idl1;
57  ipp2 = *ip + 2;
58  ipph = (*ip + 1) / 2;
59  idp = *ip * *ido;
60 
61  if (*ido < *l1) {
62  goto L106;
63  }
64  i_1 = ipph;
65  for (j = 2; j <= i_1; ++j) {
66  jc = ipp2 - j;
67  i_2 = *l1;
68  for (k = 1; k <= i_2; ++k) {
69  i_3 = *ido;
70  for (i = 1; i <= i_3; ++i) {
71  ch[i + (k + j * ch_dim2) * ch_dim1] = cc[i + (j + k * cc_dim2)
72  * cc_dim1] + cc[i + (jc + k * cc_dim2) * cc_dim1];
73  ch[i + (k + jc * ch_dim2) * ch_dim1] = cc[i + (j + k *
74  cc_dim2) * cc_dim1] - cc[i + (jc + k * cc_dim2) *
75  cc_dim1];
76 /* L101: */
77  }
78 /* L102: */
79  }
80 /* L103: */
81  }
82  i_1 = *l1;
83  for (k = 1; k <= i_1; ++k) {
84  i_2 = *ido;
85  for (i = 1; i <= i_2; ++i) {
86  ch[i + (k + ch_dim2) * ch_dim1] = cc[i + (k * cc_dim2 + 1) *
87  cc_dim1];
88 /* L104: */
89  }
90 /* L105: */
91  }
92  goto L112;
93 L106:
94  i_1 = ipph;
95  for (j = 2; j <= i_1; ++j) {
96  jc = ipp2 - j;
97  i_2 = *ido;
98  for (i = 1; i <= i_2; ++i) {
99  i_3 = *l1;
100  for (k = 1; k <= i_3; ++k) {
101  ch[i + (k + j * ch_dim2) * ch_dim1] = cc[i + (j + k * cc_dim2)
102  * cc_dim1] + cc[i + (jc + k * cc_dim2) * cc_dim1];
103  ch[i + (k + jc * ch_dim2) * ch_dim1] = cc[i + (j + k *
104  cc_dim2) * cc_dim1] - cc[i + (jc + k * cc_dim2) *
105  cc_dim1];
106 /* L107: */
107  }
108 /* L108: */
109  }
110 /* L109: */
111  }
112  i_1 = *ido;
113  for (i = 1; i <= i_1; ++i) {
114  i_2 = *l1;
115  for (k = 1; k <= i_2; ++k) {
116  ch[i + (k + ch_dim2) * ch_dim1] = cc[i + (k * cc_dim2 + 1) *
117  cc_dim1];
118 /* L110: */
119  }
120 /* L111: */
121  }
122 L112:
123  idl = 2 - *ido;
124  inc = 0;
125  i_1 = ipph;
126  for (l = 2; l <= i_1; ++l) {
127  lc = ipp2 - l;
128  idl += *ido;
129  i_2 = *idl1;
130  for (ik = 1; ik <= i_2; ++ik) {
131  c2[ik + l * c2_dim1] = ch2[ik + ch2_dim1] + wa[idl - 1] * ch2[ik
132  + (ch2_dim1 << 1)];
133  c2[ik + lc * c2_dim1] = wa[idl] * ch2[ik + *ip * ch2_dim1];
134 /* L113: */
135  }
136  idlj = idl;
137  inc += *ido;
138  i_2 = ipph;
139  for (j = 3; j <= i_2; ++j) {
140  jc = ipp2 - j;
141  idlj += inc;
142  if (idlj > idp) {
143  idlj -= idp;
144  }
145  war = wa[idlj - 1];
146  wai = wa[idlj];
147  i_3 = *idl1;
148  for (ik = 1; ik <= i_3; ++ik) {
149  c2[ik + l * c2_dim1] += war * ch2[ik + j * ch2_dim1];
150  c2[ik + lc * c2_dim1] += wai * ch2[ik + jc * ch2_dim1];
151 /* L114: */
152  }
153 /* L115: */
154  }
155 /* L116: */
156  }
157  i_1 = ipph;
158  for (j = 2; j <= i_1; ++j) {
159  i_2 = *idl1;
160  for (ik = 1; ik <= i_2; ++ik) {
161  ch2[ik + ch2_dim1] += ch2[ik + j * ch2_dim1];
162 /* L117: */
163  }
164 /* L118: */
165  }
166  i_1 = ipph;
167  for (j = 2; j <= i_1; ++j) {
168  jc = ipp2 - j;
169  i_2 = *idl1;
170  for (ik = 2; ik <= i_2; ik += 2) {
171  ch2[ik - 1 + j * ch2_dim1] = c2[ik - 1 + j * c2_dim1] - c2[ik +
172  jc * c2_dim1];
173  ch2[ik - 1 + jc * ch2_dim1] = c2[ik - 1 + j * c2_dim1] + c2[ik +
174  jc * c2_dim1];
175  ch2[ik + j * ch2_dim1] = c2[ik + j * c2_dim1] + c2[ik - 1 + jc *
176  c2_dim1];
177  ch2[ik + jc * ch2_dim1] = c2[ik + j * c2_dim1] - c2[ik - 1 + jc *
178  c2_dim1];
179 /* L119: */
180  }
181 /* L120: */
182  }
183  *nac = 1;
184  if (*ido == 2) {
185  return 0;
186  }
187  *nac = 0;
188  i_1 = *idl1;
189  for (ik = 1; ik <= i_1; ++ik) {
190  c2[ik + c2_dim1] = ch2[ik + ch2_dim1];
191 /* L121: */
192  }
193  i_1 = *ip;
194  for (j = 2; j <= i_1; ++j) {
195  i_2 = *l1;
196  for (k = 1; k <= i_2; ++k) {
197  c1[(k + j * c1_dim2) * c1_dim1 + 1] = ch[(k + j * ch_dim2) *
198  ch_dim1 + 1];
199  c1[(k + j * c1_dim2) * c1_dim1 + 2] = ch[(k + j * ch_dim2) *
200  ch_dim1 + 2];
201 /* L122: */
202  }
203 /* L123: */
204  }
205  if (idot > *l1) {
206  goto L127;
207  }
208  idij = 0;
209  i_1 = *ip;
210  for (j = 2; j <= i_1; ++j) {
211  idij += 2;
212  i_2 = *ido;
213  for (i = 4; i <= i_2; i += 2) {
214  idij += 2;
215  i_3 = *l1;
216  for (k = 1; k <= i_3; ++k) {
217  c1[i - 1 + (k + j * c1_dim2) * c1_dim1] = wa[idij - 1] * ch[i
218  - 1 + (k + j * ch_dim2) * ch_dim1] - wa[idij] * ch[i
219  + (k + j * ch_dim2) * ch_dim1];
220  c1[i + (k + j * c1_dim2) * c1_dim1] = wa[idij - 1] * ch[i + (
221  k + j * ch_dim2) * ch_dim1] + wa[idij] * ch[i - 1 + (
222  k + j * ch_dim2) * ch_dim1];
223 /* L124: */
224  }
225 /* L125: */
226  }
227 /* L126: */
228  }
229  return 0;
230 L127:
231  idj = 2 - *ido;
232  i_1 = *ip;
233  for (j = 2; j <= i_1; ++j) {
234  idj += *ido;
235  i_2 = *l1;
236  for (k = 1; k <= i_2; ++k) {
237  idij = idj;
238  i_3 = *ido;
239  for (i = 4; i <= i_3; i += 2) {
240  idij += 2;
241  c1[i - 1 + (k + j * c1_dim2) * c1_dim1] = wa[idij - 1] * ch[i
242  - 1 + (k + j * ch_dim2) * ch_dim1] - wa[idij] * ch[i
243  + (k + j * ch_dim2) * ch_dim1];
244  c1[i + (k + j * c1_dim2) * c1_dim1] = wa[idij - 1] * ch[i + (
245  k + j * ch_dim2) * ch_dim1] + wa[idij] * ch[i - 1 + (
246  k + j * ch_dim2) * ch_dim1];
247 /* L128: */
248  }
249 /* L129: */
250  }
251 /* L130: */
252  }
253  return 0;
254 } /* passb_ */

◆ passb2()

int passb2 ( int *  ido,
int *  l1,
double *  cc,
double *  ch,
double *  wa1 
)

Definition at line 257 of file pub3dfft.C.

Referenced by cfftb1().

259 {
260  /* System generated locals */
261  int cc_dim1, cc_offset, ch_dim1, ch_dim2, ch_offset, i_1, i_2;
262 
263  /* Local variables */
264  static int i, k;
265  static double ti2, tr2;
266 
267  /* Parameter adjustments */
268  cc_dim1 = *ido;
269  cc_offset = cc_dim1 * 3 + 1;
270  cc -= cc_offset;
271  ch_dim1 = *ido;
272  ch_dim2 = *l1;
273  ch_offset = ch_dim1 * (ch_dim2 + 1) + 1;
274  ch -= ch_offset;
275  --wa1;
276 
277  /* Function Body */
278  if (*ido > 2) {
279  goto L102;
280  }
281  i_1 = *l1;
282  for (k = 1; k <= i_1; ++k) {
283  ch[(k + ch_dim2) * ch_dim1 + 1] = cc[((k << 1) + 1) * cc_dim1 + 1] +
284  cc[((k << 1) + 2) * cc_dim1 + 1];
285  ch[(k + (ch_dim2 << 1)) * ch_dim1 + 1] = cc[((k << 1) + 1) * cc_dim1
286  + 1] - cc[((k << 1) + 2) * cc_dim1 + 1];
287  ch[(k + ch_dim2) * ch_dim1 + 2] = cc[((k << 1) + 1) * cc_dim1 + 2] +
288  cc[((k << 1) + 2) * cc_dim1 + 2];
289  ch[(k + (ch_dim2 << 1)) * ch_dim1 + 2] = cc[((k << 1) + 1) * cc_dim1
290  + 2] - cc[((k << 1) + 2) * cc_dim1 + 2];
291 /* L101: */
292  }
293  return 0;
294 L102:
295  i_1 = *l1;
296  for (k = 1; k <= i_1; ++k) {
297  i_2 = *ido;
298  for (i = 2; i <= i_2; i += 2) {
299  ch[i - 1 + (k + ch_dim2) * ch_dim1] = cc[i - 1 + ((k << 1) + 1) *
300  cc_dim1] + cc[i - 1 + ((k << 1) + 2) * cc_dim1];
301  tr2 = cc[i - 1 + ((k << 1) + 1) * cc_dim1] - cc[i - 1 + ((k << 1)
302  + 2) * cc_dim1];
303  ch[i + (k + ch_dim2) * ch_dim1] = cc[i + ((k << 1) + 1) * cc_dim1]
304  + cc[i + ((k << 1) + 2) * cc_dim1];
305  ti2 = cc[i + ((k << 1) + 1) * cc_dim1] - cc[i + ((k << 1) + 2) *
306  cc_dim1];
307  ch[i + (k + (ch_dim2 << 1)) * ch_dim1] = wa1[i - 1] * ti2 + wa1[i]
308  * tr2;
309  ch[i - 1 + (k + (ch_dim2 << 1)) * ch_dim1] = wa1[i - 1] * tr2 -
310  wa1[i] * ti2;
311 /* L103: */
312  }
313 /* L104: */
314  }
315  return 0;
316 } /* passb2_ */

◆ passb3()

int passb3 ( int *  ido,
int *  l1,
double *  cc,
double *  ch,
double *  wa1,
double *  wa2 
)

Definition at line 319 of file pub3dfft.C.

Referenced by cfftb1().

321 {
322  /* Initialized data */
323 
324  static double taur = -.5;
325  static double taui = .866025403784439;
326 
327  /* System generated locals */
328  int cc_dim1, cc_offset, ch_dim1, ch_dim2, ch_offset, i_1, i_2;
329 
330  /* Local variables */
331  static int i, k;
332  static double ci2, ci3, di2, di3, cr2, cr3, dr2, dr3, ti2, tr2;
333 
334  /* Parameter adjustments */
335  cc_dim1 = *ido;
336  cc_offset = (cc_dim1 << 2) + 1;
337  cc -= cc_offset;
338  ch_dim1 = *ido;
339  ch_dim2 = *l1;
340  ch_offset = ch_dim1 * (ch_dim2 + 1) + 1;
341  ch -= ch_offset;
342  --wa1;
343  --wa2;
344 
345  /* Function Body */
346  if (*ido != 2) {
347  goto L102;
348  }
349  i_1 = *l1;
350  for (k = 1; k <= i_1; ++k) {
351  tr2 = cc[(k * 3 + 2) * cc_dim1 + 1] + cc[(k * 3 + 3) * cc_dim1 + 1];
352  cr2 = cc[(k * 3 + 1) * cc_dim1 + 1] + taur * tr2;
353  ch[(k + ch_dim2) * ch_dim1 + 1] = cc[(k * 3 + 1) * cc_dim1 + 1] + tr2;
354 
355  ti2 = cc[(k * 3 + 2) * cc_dim1 + 2] + cc[(k * 3 + 3) * cc_dim1 + 2];
356  ci2 = cc[(k * 3 + 1) * cc_dim1 + 2] + taur * ti2;
357  ch[(k + ch_dim2) * ch_dim1 + 2] = cc[(k * 3 + 1) * cc_dim1 + 2] + ti2;
358 
359  cr3 = taui * (cc[(k * 3 + 2) * cc_dim1 + 1] - cc[(k * 3 + 3) *
360  cc_dim1 + 1]);
361  ci3 = taui * (cc[(k * 3 + 2) * cc_dim1 + 2] - cc[(k * 3 + 3) *
362  cc_dim1 + 2]);
363  ch[(k + (ch_dim2 << 1)) * ch_dim1 + 1] = cr2 - ci3;
364  ch[(k + ch_dim2 * 3) * ch_dim1 + 1] = cr2 + ci3;
365  ch[(k + (ch_dim2 << 1)) * ch_dim1 + 2] = ci2 + cr3;
366  ch[(k + ch_dim2 * 3) * ch_dim1 + 2] = ci2 - cr3;
367 /* L101: */
368  }
369  return 0;
370 L102:
371  i_1 = *l1;
372  for (k = 1; k <= i_1; ++k) {
373  i_2 = *ido;
374  for (i = 2; i <= i_2; i += 2) {
375  tr2 = cc[i - 1 + (k * 3 + 2) * cc_dim1] + cc[i - 1 + (k * 3 + 3) *
376  cc_dim1];
377  cr2 = cc[i - 1 + (k * 3 + 1) * cc_dim1] + taur * tr2;
378  ch[i - 1 + (k + ch_dim2) * ch_dim1] = cc[i - 1 + (k * 3 + 1) *
379  cc_dim1] + tr2;
380  ti2 = cc[i + (k * 3 + 2) * cc_dim1] + cc[i + (k * 3 + 3) *
381  cc_dim1];
382  ci2 = cc[i + (k * 3 + 1) * cc_dim1] + taur * ti2;
383  ch[i + (k + ch_dim2) * ch_dim1] = cc[i + (k * 3 + 1) * cc_dim1] +
384  ti2;
385  cr3 = taui * (cc[i - 1 + (k * 3 + 2) * cc_dim1] - cc[i - 1 + (k *
386  3 + 3) * cc_dim1]);
387  ci3 = taui * (cc[i + (k * 3 + 2) * cc_dim1] - cc[i + (k * 3 + 3) *
388  cc_dim1]);
389  dr2 = cr2 - ci3;
390  dr3 = cr2 + ci3;
391  di2 = ci2 + cr3;
392  di3 = ci2 - cr3;
393  ch[i + (k + (ch_dim2 << 1)) * ch_dim1] = wa1[i - 1] * di2 + wa1[i]
394  * dr2;
395  ch[i - 1 + (k + (ch_dim2 << 1)) * ch_dim1] = wa1[i - 1] * dr2 -
396  wa1[i] * di2;
397  ch[i + (k + ch_dim2 * 3) * ch_dim1] = wa2[i - 1] * di3 + wa2[i] *
398  dr3;
399  ch[i - 1 + (k + ch_dim2 * 3) * ch_dim1] = wa2[i - 1] * dr3 - wa2[
400  i] * di3;
401 /* L103: */
402  }
403 /* L104: */
404  }
405  return 0;
406 } /* passb3_ */

◆ passb4()

int passb4 ( int *  ido,
int *  l1,
double *  cc,
double *  ch,
double *  wa1,
double *  wa2,
double *  wa3 
)

Definition at line 408 of file pub3dfft.C.

Referenced by cfftb1().

410 {
411  /* System generated locals */
412  int cc_dim1, cc_offset, ch_dim1, ch_dim2, ch_offset, i_1, i_2;
413 
414  /* Local variables */
415  static int i, k;
416  static double ci2, ci3, ci4, cr2, cr3, cr4, ti1, ti2, ti3, ti4, tr1,
417  tr2, tr3, tr4;
418 
419  /* Parameter adjustments */
420  cc_dim1 = *ido;
421  cc_offset = cc_dim1 * 5 + 1;
422  cc -= cc_offset;
423  ch_dim1 = *ido;
424  ch_dim2 = *l1;
425  ch_offset = ch_dim1 * (ch_dim2 + 1) + 1;
426  ch -= ch_offset;
427  --wa1;
428  --wa2;
429  --wa3;
430 
431  /* Function Body */
432  if (*ido != 2) {
433  goto L102;
434  }
435  i_1 = *l1;
436  for (k = 1; k <= i_1; ++k) {
437  ti1 = cc[((k << 2) + 1) * cc_dim1 + 2] - cc[((k << 2) + 3) * cc_dim1
438  + 2];
439  ti2 = cc[((k << 2) + 1) * cc_dim1 + 2] + cc[((k << 2) + 3) * cc_dim1
440  + 2];
441  tr4 = cc[((k << 2) + 4) * cc_dim1 + 2] - cc[((k << 2) + 2) * cc_dim1
442  + 2];
443  ti3 = cc[((k << 2) + 2) * cc_dim1 + 2] + cc[((k << 2) + 4) * cc_dim1
444  + 2];
445  tr1 = cc[((k << 2) + 1) * cc_dim1 + 1] - cc[((k << 2) + 3) * cc_dim1
446  + 1];
447  tr2 = cc[((k << 2) + 1) * cc_dim1 + 1] + cc[((k << 2) + 3) * cc_dim1
448  + 1];
449  ti4 = cc[((k << 2) + 2) * cc_dim1 + 1] - cc[((k << 2) + 4) * cc_dim1
450  + 1];
451  tr3 = cc[((k << 2) + 2) * cc_dim1 + 1] + cc[((k << 2) + 4) * cc_dim1
452  + 1];
453  ch[(k + ch_dim2) * ch_dim1 + 1] = tr2 + tr3;
454  ch[(k + ch_dim2 * 3) * ch_dim1 + 1] = tr2 - tr3;
455  ch[(k + ch_dim2) * ch_dim1 + 2] = ti2 + ti3;
456  ch[(k + ch_dim2 * 3) * ch_dim1 + 2] = ti2 - ti3;
457  ch[(k + (ch_dim2 << 1)) * ch_dim1 + 1] = tr1 + tr4;
458  ch[(k + (ch_dim2 << 2)) * ch_dim1 + 1] = tr1 - tr4;
459  ch[(k + (ch_dim2 << 1)) * ch_dim1 + 2] = ti1 + ti4;
460  ch[(k + (ch_dim2 << 2)) * ch_dim1 + 2] = ti1 - ti4;
461 /* L101: */
462  }
463  return 0;
464 L102:
465  i_1 = *l1;
466  for (k = 1; k <= i_1; ++k) {
467  i_2 = *ido;
468  for (i = 2; i <= i_2; i += 2) {
469  ti1 = cc[i + ((k << 2) + 1) * cc_dim1] - cc[i + ((k << 2) + 3) *
470  cc_dim1];
471  ti2 = cc[i + ((k << 2) + 1) * cc_dim1] + cc[i + ((k << 2) + 3) *
472  cc_dim1];
473  ti3 = cc[i + ((k << 2) + 2) * cc_dim1] + cc[i + ((k << 2) + 4) *
474  cc_dim1];
475  tr4 = cc[i + ((k << 2) + 4) * cc_dim1] - cc[i + ((k << 2) + 2) *
476  cc_dim1];
477  tr1 = cc[i - 1 + ((k << 2) + 1) * cc_dim1] - cc[i - 1 + ((k << 2)
478  + 3) * cc_dim1];
479  tr2 = cc[i - 1 + ((k << 2) + 1) * cc_dim1] + cc[i - 1 + ((k << 2)
480  + 3) * cc_dim1];
481  ti4 = cc[i - 1 + ((k << 2) + 2) * cc_dim1] - cc[i - 1 + ((k << 2)
482  + 4) * cc_dim1];
483  tr3 = cc[i - 1 + ((k << 2) + 2) * cc_dim1] + cc[i - 1 + ((k << 2)
484  + 4) * cc_dim1];
485  ch[i - 1 + (k + ch_dim2) * ch_dim1] = tr2 + tr3;
486  cr3 = tr2 - tr3;
487  ch[i + (k + ch_dim2) * ch_dim1] = ti2 + ti3;
488  ci3 = ti2 - ti3;
489  cr2 = tr1 + tr4;
490  cr4 = tr1 - tr4;
491  ci2 = ti1 + ti4;
492  ci4 = ti1 - ti4;
493  ch[i - 1 + (k + (ch_dim2 << 1)) * ch_dim1] = wa1[i - 1] * cr2 -
494  wa1[i] * ci2;
495  ch[i + (k + (ch_dim2 << 1)) * ch_dim1] = wa1[i - 1] * ci2 + wa1[i]
496  * cr2;
497  ch[i - 1 + (k + ch_dim2 * 3) * ch_dim1] = wa2[i - 1] * cr3 - wa2[
498  i] * ci3;
499  ch[i + (k + ch_dim2 * 3) * ch_dim1] = wa2[i - 1] * ci3 + wa2[i] *
500  cr3;
501  ch[i - 1 + (k + (ch_dim2 << 2)) * ch_dim1] = wa3[i - 1] * cr4 -
502  wa3[i] * ci4;
503  ch[i + (k + (ch_dim2 << 2)) * ch_dim1] = wa3[i - 1] * ci4 + wa3[i]
504  * cr4;
505 /* L103: */
506  }
507 /* L104: */
508  }
509  return 0;
510 } /* passb4_ */

◆ passb5()

int passb5 ( int *  ido,
int *  l1,
double *  cc,
double *  ch,
double *  wa1,
double *  wa2,
double *  wa3,
double *  wa4 
)

Definition at line 512 of file pub3dfft.C.

Referenced by cfftb1().

515 {
516  /* Initialized data */
517 
518  static double tr11 = .309016994374947;
519  static double ti11 = .951056516295154;
520  static double tr12 = -.809016994374947;
521  static double ti12 = .587785252292473;
522 
523  /* System generated locals */
524  int cc_dim1, cc_offset, ch_dim1, ch_dim2, ch_offset, i_1, i_2;
525 
526  /* Local variables */
527  static int i, k;
528  static double ci2, ci3, ci4, ci5, di3, di4, di5, di2, cr2, cr3, cr5,
529  cr4, ti2, ti3, ti4, ti5, dr3, dr4, dr5, dr2, tr2, tr3, tr4, tr5;
530 
531  /* Parameter adjustments */
532  cc_dim1 = *ido;
533  cc_offset = cc_dim1 * 6 + 1;
534  cc -= cc_offset;
535  ch_dim1 = *ido;
536  ch_dim2 = *l1;
537  ch_offset = ch_dim1 * (ch_dim2 + 1) + 1;
538  ch -= ch_offset;
539  --wa1;
540  --wa2;
541  --wa3;
542  --wa4;
543 
544  /* Function Body */
545  if (*ido != 2) {
546  goto L102;
547  }
548  i_1 = *l1;
549  for (k = 1; k <= i_1; ++k) {
550  ti5 = cc[(k * 5 + 2) * cc_dim1 + 2] - cc[(k * 5 + 5) * cc_dim1 + 2];
551  ti2 = cc[(k * 5 + 2) * cc_dim1 + 2] + cc[(k * 5 + 5) * cc_dim1 + 2];
552  ti4 = cc[(k * 5 + 3) * cc_dim1 + 2] - cc[(k * 5 + 4) * cc_dim1 + 2];
553  ti3 = cc[(k * 5 + 3) * cc_dim1 + 2] + cc[(k * 5 + 4) * cc_dim1 + 2];
554  tr5 = cc[(k * 5 + 2) * cc_dim1 + 1] - cc[(k * 5 + 5) * cc_dim1 + 1];
555  tr2 = cc[(k * 5 + 2) * cc_dim1 + 1] + cc[(k * 5 + 5) * cc_dim1 + 1];
556  tr4 = cc[(k * 5 + 3) * cc_dim1 + 1] - cc[(k * 5 + 4) * cc_dim1 + 1];
557  tr3 = cc[(k * 5 + 3) * cc_dim1 + 1] + cc[(k * 5 + 4) * cc_dim1 + 1];
558  ch[(k + ch_dim2) * ch_dim1 + 1] = cc[(k * 5 + 1) * cc_dim1 + 1] + tr2
559  + tr3;
560  ch[(k + ch_dim2) * ch_dim1 + 2] = cc[(k * 5 + 1) * cc_dim1 + 2] + ti2
561  + ti3;
562  cr2 = cc[(k * 5 + 1) * cc_dim1 + 1] + tr11 * tr2 + tr12 * tr3;
563  ci2 = cc[(k * 5 + 1) * cc_dim1 + 2] + tr11 * ti2 + tr12 * ti3;
564  cr3 = cc[(k * 5 + 1) * cc_dim1 + 1] + tr12 * tr2 + tr11 * tr3;
565  ci3 = cc[(k * 5 + 1) * cc_dim1 + 2] + tr12 * ti2 + tr11 * ti3;
566  cr5 = ti11 * tr5 + ti12 * tr4;
567  ci5 = ti11 * ti5 + ti12 * ti4;
568  cr4 = ti12 * tr5 - ti11 * tr4;
569  ci4 = ti12 * ti5 - ti11 * ti4;
570  ch[(k + (ch_dim2 << 1)) * ch_dim1 + 1] = cr2 - ci5;
571  ch[(k + ch_dim2 * 5) * ch_dim1 + 1] = cr2 + ci5;
572  ch[(k + (ch_dim2 << 1)) * ch_dim1 + 2] = ci2 + cr5;
573  ch[(k + ch_dim2 * 3) * ch_dim1 + 2] = ci3 + cr4;
574  ch[(k + ch_dim2 * 3) * ch_dim1 + 1] = cr3 - ci4;
575  ch[(k + (ch_dim2 << 2)) * ch_dim1 + 1] = cr3 + ci4;
576  ch[(k + (ch_dim2 << 2)) * ch_dim1 + 2] = ci3 - cr4;
577  ch[(k + ch_dim2 * 5) * ch_dim1 + 2] = ci2 - cr5;
578 /* L101: */
579  }
580  return 0;
581 L102:
582  i_1 = *l1;
583  for (k = 1; k <= i_1; ++k) {
584  i_2 = *ido;
585  for (i = 2; i <= i_2; i += 2) {
586  ti5 = cc[i + (k * 5 + 2) * cc_dim1] - cc[i + (k * 5 + 5) *
587  cc_dim1];
588  ti2 = cc[i + (k * 5 + 2) * cc_dim1] + cc[i + (k * 5 + 5) *
589  cc_dim1];
590  ti4 = cc[i + (k * 5 + 3) * cc_dim1] - cc[i + (k * 5 + 4) *
591  cc_dim1];
592  ti3 = cc[i + (k * 5 + 3) * cc_dim1] + cc[i + (k * 5 + 4) *
593  cc_dim1];
594  tr5 = cc[i - 1 + (k * 5 + 2) * cc_dim1] - cc[i - 1 + (k * 5 + 5) *
595  cc_dim1];
596  tr2 = cc[i - 1 + (k * 5 + 2) * cc_dim1] + cc[i - 1 + (k * 5 + 5) *
597  cc_dim1];
598  tr4 = cc[i - 1 + (k * 5 + 3) * cc_dim1] - cc[i - 1 + (k * 5 + 4) *
599  cc_dim1];
600  tr3 = cc[i - 1 + (k * 5 + 3) * cc_dim1] + cc[i - 1 + (k * 5 + 4) *
601  cc_dim1];
602  ch[i - 1 + (k + ch_dim2) * ch_dim1] = cc[i - 1 + (k * 5 + 1) *
603  cc_dim1] + tr2 + tr3;
604  ch[i + (k + ch_dim2) * ch_dim1] = cc[i + (k * 5 + 1) * cc_dim1] +
605  ti2 + ti3;
606  cr2 = cc[i - 1 + (k * 5 + 1) * cc_dim1] + tr11 * tr2 + tr12 * tr3;
607 
608  ci2 = cc[i + (k * 5 + 1) * cc_dim1] + tr11 * ti2 + tr12 * ti3;
609  cr3 = cc[i - 1 + (k * 5 + 1) * cc_dim1] + tr12 * tr2 + tr11 * tr3;
610 
611  ci3 = cc[i + (k * 5 + 1) * cc_dim1] + tr12 * ti2 + tr11 * ti3;
612  cr5 = ti11 * tr5 + ti12 * tr4;
613  ci5 = ti11 * ti5 + ti12 * ti4;
614  cr4 = ti12 * tr5 - ti11 * tr4;
615  ci4 = ti12 * ti5 - ti11 * ti4;
616  dr3 = cr3 - ci4;
617  dr4 = cr3 + ci4;
618  di3 = ci3 + cr4;
619  di4 = ci3 - cr4;
620  dr5 = cr2 + ci5;
621  dr2 = cr2 - ci5;
622  di5 = ci2 - cr5;
623  di2 = ci2 + cr5;
624  ch[i - 1 + (k + (ch_dim2 << 1)) * ch_dim1] = wa1[i - 1] * dr2 -
625  wa1[i] * di2;
626  ch[i + (k + (ch_dim2 << 1)) * ch_dim1] = wa1[i - 1] * di2 + wa1[i]
627  * dr2;
628  ch[i - 1 + (k + ch_dim2 * 3) * ch_dim1] = wa2[i - 1] * dr3 - wa2[
629  i] * di3;
630  ch[i + (k + ch_dim2 * 3) * ch_dim1] = wa2[i - 1] * di3 + wa2[i] *
631  dr3;
632  ch[i - 1 + (k + (ch_dim2 << 2)) * ch_dim1] = wa3[i - 1] * dr4 -
633  wa3[i] * di4;
634  ch[i + (k + (ch_dim2 << 2)) * ch_dim1] = wa3[i - 1] * di4 + wa3[i]
635  * dr4;
636  ch[i - 1 + (k + ch_dim2 * 5) * ch_dim1] = wa4[i - 1] * dr5 - wa4[
637  i] * di5;
638  ch[i + (k + ch_dim2 * 5) * ch_dim1] = wa4[i - 1] * di5 + wa4[i] *
639  dr5;
640 /* L103: */
641  }
642 /* L104: */
643  }
644  return 0;
645 } /* passb5_ */

◆ passf()

int passf ( int *  nac,
int *  ido,
int *  ip,
int *  l1,
int *  idl1,
double *  cc,
double *  c1,
double *  c2,
double *  ch,
double *  ch2,
double *  wa 
)

Definition at line 647 of file pub3dfft.C.

Referenced by cfftf1().

650 {
651  /* System generated locals */
652  int ch_dim1, ch_dim2, ch_offset, cc_dim1, cc_dim2, cc_offset, c1_dim1,
653  c1_dim2, c1_offset, c2_dim1, c2_offset, ch2_dim1, ch2_offset,
654  i_1, i_2, i_3;
655 
656  /* Local variables */
657  static int idij, idlj, idot, ipph, i, j, k, l, jc, lc, ik, nt, idj,
658  idl, inc, idp;
659  static double wai, war;
660  static int ipp2;
661 
662  /* Parameter adjustments */
663  cc_dim1 = *ido;
664  cc_dim2 = *ip;
665  cc_offset = cc_dim1 * (cc_dim2 + 1) + 1;
666  cc -= cc_offset;
667  c1_dim1 = *ido;
668  c1_dim2 = *l1;
669  c1_offset = c1_dim1 * (c1_dim2 + 1) + 1;
670  c1 -= c1_offset;
671  c2_dim1 = *idl1;
672  c2_offset = c2_dim1 + 1;
673  c2 -= c2_offset;
674  ch_dim1 = *ido;
675  ch_dim2 = *l1;
676  ch_offset = ch_dim1 * (ch_dim2 + 1) + 1;
677  ch -= ch_offset;
678  ch2_dim1 = *idl1;
679  ch2_offset = ch2_dim1 + 1;
680  ch2 -= ch2_offset;
681  --wa;
682 
683  /* Function Body */
684  idot = *ido / 2;
685  nt = *ip * *idl1;
686  ipp2 = *ip + 2;
687  ipph = (*ip + 1) / 2;
688  idp = *ip * *ido;
689 
690  if (*ido < *l1) {
691  goto L106;
692  }
693  i_1 = ipph;
694  for (j = 2; j <= i_1; ++j) {
695  jc = ipp2 - j;
696  i_2 = *l1;
697  for (k = 1; k <= i_2; ++k) {
698  i_3 = *ido;
699  for (i = 1; i <= i_3; ++i) {
700  ch[i + (k + j * ch_dim2) * ch_dim1] = cc[i + (j + k * cc_dim2)
701  * cc_dim1] + cc[i + (jc + k * cc_dim2) * cc_dim1];
702  ch[i + (k + jc * ch_dim2) * ch_dim1] = cc[i + (j + k *
703  cc_dim2) * cc_dim1] - cc[i + (jc + k * cc_dim2) *
704  cc_dim1];
705 /* L101: */
706  }
707 /* L102: */
708  }
709 /* L103: */
710  }
711  i_1 = *l1;
712  for (k = 1; k <= i_1; ++k) {
713  i_2 = *ido;
714  for (i = 1; i <= i_2; ++i) {
715  ch[i + (k + ch_dim2) * ch_dim1] = cc[i + (k * cc_dim2 + 1) *
716  cc_dim1];
717 /* L104: */
718  }
719 /* L105: */
720  }
721  goto L112;
722 L106:
723  i_1 = ipph;
724  for (j = 2; j <= i_1; ++j) {
725  jc = ipp2 - j;
726  i_2 = *ido;
727  for (i = 1; i <= i_2; ++i) {
728  i_3 = *l1;
729  for (k = 1; k <= i_3; ++k) {
730  ch[i + (k + j * ch_dim2) * ch_dim1] = cc[i + (j + k * cc_dim2)
731  * cc_dim1] + cc[i + (jc + k * cc_dim2) * cc_dim1];
732  ch[i + (k + jc * ch_dim2) * ch_dim1] = cc[i + (j + k *
733  cc_dim2) * cc_dim1] - cc[i + (jc + k * cc_dim2) *
734  cc_dim1];
735 /* L107: */
736  }
737 /* L108: */
738  }
739 /* L109: */
740  }
741  i_1 = *ido;
742  for (i = 1; i <= i_1; ++i) {
743  i_2 = *l1;
744  for (k = 1; k <= i_2; ++k) {
745  ch[i + (k + ch_dim2) * ch_dim1] = cc[i + (k * cc_dim2 + 1) *
746  cc_dim1];
747 /* L110: */
748  }
749 /* L111: */
750  }
751 L112:
752  idl = 2 - *ido;
753  inc = 0;
754  i_1 = ipph;
755  for (l = 2; l <= i_1; ++l) {
756  lc = ipp2 - l;
757  idl += *ido;
758  i_2 = *idl1;
759  for (ik = 1; ik <= i_2; ++ik) {
760  c2[ik + l * c2_dim1] = ch2[ik + ch2_dim1] + wa[idl - 1] * ch2[ik
761  + (ch2_dim1 << 1)];
762  c2[ik + lc * c2_dim1] = -wa[idl] * ch2[ik + *ip * ch2_dim1];
763 /* L113: */
764  }
765  idlj = idl;
766  inc += *ido;
767  i_2 = ipph;
768  for (j = 3; j <= i_2; ++j) {
769  jc = ipp2 - j;
770  idlj += inc;
771  if (idlj > idp) {
772  idlj -= idp;
773  }
774  war = wa[idlj - 1];
775  wai = wa[idlj];
776  i_3 = *idl1;
777  for (ik = 1; ik <= i_3; ++ik) {
778  c2[ik + l * c2_dim1] += war * ch2[ik + j * ch2_dim1];
779  c2[ik + lc * c2_dim1] -= wai * ch2[ik + jc * ch2_dim1];
780 /* L114: */
781  }
782 /* L115: */
783  }
784 /* L116: */
785  }
786  i_1 = ipph;
787  for (j = 2; j <= i_1; ++j) {
788  i_2 = *idl1;
789  for (ik = 1; ik <= i_2; ++ik) {
790  ch2[ik + ch2_dim1] += ch2[ik + j * ch2_dim1];
791 /* L117: */
792  }
793 /* L118: */
794  }
795  i_1 = ipph;
796  for (j = 2; j <= i_1; ++j) {
797  jc = ipp2 - j;
798  i_2 = *idl1;
799  for (ik = 2; ik <= i_2; ik += 2) {
800  ch2[ik - 1 + j * ch2_dim1] = c2[ik - 1 + j * c2_dim1] - c2[ik +
801  jc * c2_dim1];
802  ch2[ik - 1 + jc * ch2_dim1] = c2[ik - 1 + j * c2_dim1] + c2[ik +
803  jc * c2_dim1];
804  ch2[ik + j * ch2_dim1] = c2[ik + j * c2_dim1] + c2[ik - 1 + jc *
805  c2_dim1];
806  ch2[ik + jc * ch2_dim1] = c2[ik + j * c2_dim1] - c2[ik - 1 + jc *
807  c2_dim1];
808 /* L119: */
809  }
810 /* L120: */
811  }
812  *nac = 1;
813  if (*ido == 2) {
814  return 0;
815  }
816  *nac = 0;
817  i_1 = *idl1;
818  for (ik = 1; ik <= i_1; ++ik) {
819  c2[ik + c2_dim1] = ch2[ik + ch2_dim1];
820 /* L121: */
821  }
822  i_1 = *ip;
823  for (j = 2; j <= i_1; ++j) {
824  i_2 = *l1;
825  for (k = 1; k <= i_2; ++k) {
826  c1[(k + j * c1_dim2) * c1_dim1 + 1] = ch[(k + j * ch_dim2) *
827  ch_dim1 + 1];
828  c1[(k + j * c1_dim2) * c1_dim1 + 2] = ch[(k + j * ch_dim2) *
829  ch_dim1 + 2];
830 /* L122: */
831  }
832 /* L123: */
833  }
834  if (idot > *l1) {
835  goto L127;
836  }
837  idij = 0;
838  i_1 = *ip;
839  for (j = 2; j <= i_1; ++j) {
840  idij += 2;
841  i_2 = *ido;
842  for (i = 4; i <= i_2; i += 2) {
843  idij += 2;
844  i_3 = *l1;
845  for (k = 1; k <= i_3; ++k) {
846  c1[i - 1 + (k + j * c1_dim2) * c1_dim1] = wa[idij - 1] * ch[i
847  - 1 + (k + j * ch_dim2) * ch_dim1] + wa[idij] * ch[i
848  + (k + j * ch_dim2) * ch_dim1];
849  c1[i + (k + j * c1_dim2) * c1_dim1] = wa[idij - 1] * ch[i + (
850  k + j * ch_dim2) * ch_dim1] - wa[idij] * ch[i - 1 + (
851  k + j * ch_dim2) * ch_dim1];
852 /* L124: */
853  }
854 /* L125: */
855  }
856 /* L126: */
857  }
858  return 0;
859 L127:
860  idj = 2 - *ido;
861  i_1 = *ip;
862  for (j = 2; j <= i_1; ++j) {
863  idj += *ido;
864  i_2 = *l1;
865  for (k = 1; k <= i_2; ++k) {
866  idij = idj;
867  i_3 = *ido;
868  for (i = 4; i <= i_3; i += 2) {
869  idij += 2;
870  c1[i - 1 + (k + j * c1_dim2) * c1_dim1] = wa[idij - 1] * ch[i
871  - 1 + (k + j * ch_dim2) * ch_dim1] + wa[idij] * ch[i
872  + (k + j * ch_dim2) * ch_dim1];
873  c1[i + (k + j * c1_dim2) * c1_dim1] = wa[idij - 1] * ch[i + (
874  k + j * ch_dim2) * ch_dim1] - wa[idij] * ch[i - 1 + (
875  k + j * ch_dim2) * ch_dim1];
876 /* L128: */
877  }
878 /* L129: */
879  }
880 /* L130: */
881  }
882  return 0;
883 } /* passf_ */

◆ passf2()

int passf2 ( int *  ido,
int *  l1,
double *  cc,
double *  ch,
double *  wa1 
)

Definition at line 885 of file pub3dfft.C.

Referenced by cfftf1().

887 {
888  /* System generated locals */
889  int cc_dim1, cc_offset, ch_dim1, ch_dim2, ch_offset, i_1, i_2;
890 
891  /* Local variables */
892  static int i, k;
893  static double ti2, tr2;
894 
895  /* Parameter adjustments */
896  cc_dim1 = *ido;
897  cc_offset = cc_dim1 * 3 + 1;
898  cc -= cc_offset;
899  ch_dim1 = *ido;
900  ch_dim2 = *l1;
901  ch_offset = ch_dim1 * (ch_dim2 + 1) + 1;
902  ch -= ch_offset;
903  --wa1;
904 
905  /* Function Body */
906  if (*ido > 2) {
907  goto L102;
908  }
909  i_1 = *l1;
910  for (k = 1; k <= i_1; ++k) {
911  ch[(k + ch_dim2) * ch_dim1 + 1] = cc[((k << 1) + 1) * cc_dim1 + 1] +
912  cc[((k << 1) + 2) * cc_dim1 + 1];
913  ch[(k + (ch_dim2 << 1)) * ch_dim1 + 1] = cc[((k << 1) + 1) * cc_dim1
914  + 1] - cc[((k << 1) + 2) * cc_dim1 + 1];
915  ch[(k + ch_dim2) * ch_dim1 + 2] = cc[((k << 1) + 1) * cc_dim1 + 2] +
916  cc[((k << 1) + 2) * cc_dim1 + 2];
917  ch[(k + (ch_dim2 << 1)) * ch_dim1 + 2] = cc[((k << 1) + 1) * cc_dim1
918  + 2] - cc[((k << 1) + 2) * cc_dim1 + 2];
919 /* L101: */
920  }
921  return 0;
922 L102:
923  i_1 = *l1;
924  for (k = 1; k <= i_1; ++k) {
925  i_2 = *ido;
926  for (i = 2; i <= i_2; i += 2) {
927  ch[i - 1 + (k + ch_dim2) * ch_dim1] = cc[i - 1 + ((k << 1) + 1) *
928  cc_dim1] + cc[i - 1 + ((k << 1) + 2) * cc_dim1];
929  tr2 = cc[i - 1 + ((k << 1) + 1) * cc_dim1] - cc[i - 1 + ((k << 1)
930  + 2) * cc_dim1];
931  ch[i + (k + ch_dim2) * ch_dim1] = cc[i + ((k << 1) + 1) * cc_dim1]
932  + cc[i + ((k << 1) + 2) * cc_dim1];
933  ti2 = cc[i + ((k << 1) + 1) * cc_dim1] - cc[i + ((k << 1) + 2) *
934  cc_dim1];
935  ch[i + (k + (ch_dim2 << 1)) * ch_dim1] = wa1[i - 1] * ti2 - wa1[i]
936  * tr2;
937  ch[i - 1 + (k + (ch_dim2 << 1)) * ch_dim1] = wa1[i - 1] * tr2 +
938  wa1[i] * ti2;
939 /* L103: */
940  }
941 /* L104: */
942  }
943  return 0;
944 } /* passf2_ */

◆ passf3()

int passf3 ( int *  ido,
int *  l1,
double *  cc,
double *  ch,
double *  wa1,
double *  wa2 
)

Definition at line 946 of file pub3dfft.C.

Referenced by cfftf1().

948 {
949  /* Initialized data */
950 
951  static double taur = -.5;
952  static double taui = -.866025403784439;
953 
954  /* System generated locals */
955  int cc_dim1, cc_offset, ch_dim1, ch_dim2, ch_offset, i_1, i_2;
956 
957  /* Local variables */
958  static int i, k;
959  static double ci2, ci3, di2, di3, cr2, cr3, dr2, dr3, ti2, tr2;
960 
961  /* Parameter adjustments */
962  cc_dim1 = *ido;
963  cc_offset = (cc_dim1 << 2) + 1;
964  cc -= cc_offset;
965  ch_dim1 = *ido;
966  ch_dim2 = *l1;
967  ch_offset = ch_dim1 * (ch_dim2 + 1) + 1;
968  ch -= ch_offset;
969  --wa1;
970  --wa2;
971 
972  /* Function Body */
973  if (*ido != 2) {
974  goto L102;
975  }
976  i_1 = *l1;
977  for (k = 1; k <= i_1; ++k) {
978  tr2 = cc[(k * 3 + 2) * cc_dim1 + 1] + cc[(k * 3 + 3) * cc_dim1 + 1];
979  cr2 = cc[(k * 3 + 1) * cc_dim1 + 1] + taur * tr2;
980  ch[(k + ch_dim2) * ch_dim1 + 1] = cc[(k * 3 + 1) * cc_dim1 + 1] + tr2;
981 
982  ti2 = cc[(k * 3 + 2) * cc_dim1 + 2] + cc[(k * 3 + 3) * cc_dim1 + 2];
983  ci2 = cc[(k * 3 + 1) * cc_dim1 + 2] + taur * ti2;
984  ch[(k + ch_dim2) * ch_dim1 + 2] = cc[(k * 3 + 1) * cc_dim1 + 2] + ti2;
985 
986  cr3 = taui * (cc[(k * 3 + 2) * cc_dim1 + 1] - cc[(k * 3 + 3) *
987  cc_dim1 + 1]);
988  ci3 = taui * (cc[(k * 3 + 2) * cc_dim1 + 2] - cc[(k * 3 + 3) *
989  cc_dim1 + 2]);
990  ch[(k + (ch_dim2 << 1)) * ch_dim1 + 1] = cr2 - ci3;
991  ch[(k + ch_dim2 * 3) * ch_dim1 + 1] = cr2 + ci3;
992  ch[(k + (ch_dim2 << 1)) * ch_dim1 + 2] = ci2 + cr3;
993  ch[(k + ch_dim2 * 3) * ch_dim1 + 2] = ci2 - cr3;
994 /* L101: */
995  }
996  return 0;
997 L102:
998  i_1 = *l1;
999  for (k = 1; k <= i_1; ++k) {
1000  i_2 = *ido;
1001  for (i = 2; i <= i_2; i += 2) {
1002  tr2 = cc[i - 1 + (k * 3 + 2) * cc_dim1] + cc[i - 1 + (k * 3 + 3) *
1003  cc_dim1];
1004  cr2 = cc[i - 1 + (k * 3 + 1) * cc_dim1] + taur * tr2;
1005  ch[i - 1 + (k + ch_dim2) * ch_dim1] = cc[i - 1 + (k * 3 + 1) *
1006  cc_dim1] + tr2;
1007  ti2 = cc[i + (k * 3 + 2) * cc_dim1] + cc[i + (k * 3 + 3) *
1008  cc_dim1];
1009  ci2 = cc[i + (k * 3 + 1) * cc_dim1] + taur * ti2;
1010  ch[i + (k + ch_dim2) * ch_dim1] = cc[i + (k * 3 + 1) * cc_dim1] +
1011  ti2;
1012  cr3 = taui * (cc[i - 1 + (k * 3 + 2) * cc_dim1] - cc[i - 1 + (k *
1013  3 + 3) * cc_dim1]);
1014  ci3 = taui * (cc[i + (k * 3 + 2) * cc_dim1] - cc[i + (k * 3 + 3) *
1015  cc_dim1]);
1016  dr2 = cr2 - ci3;
1017  dr3 = cr2 + ci3;
1018  di2 = ci2 + cr3;
1019  di3 = ci2 - cr3;
1020  ch[i + (k + (ch_dim2 << 1)) * ch_dim1] = wa1[i - 1] * di2 - wa1[i]
1021  * dr2;
1022  ch[i - 1 + (k + (ch_dim2 << 1)) * ch_dim1] = wa1[i - 1] * dr2 +
1023  wa1[i] * di2;
1024  ch[i + (k + ch_dim2 * 3) * ch_dim1] = wa2[i - 1] * di3 - wa2[i] *
1025  dr3;
1026  ch[i - 1 + (k + ch_dim2 * 3) * ch_dim1] = wa2[i - 1] * dr3 + wa2[
1027  i] * di3;
1028 /* L103: */
1029  }
1030 /* L104: */
1031  }
1032  return 0;
1033 } /* passf3_ */

◆ passf4()

int passf4 ( int *  ido,
int *  l1,
double *  cc,
double *  ch,
double *  wa1,
double *  wa2,
double *  wa3 
)

Definition at line 1035 of file pub3dfft.C.

Referenced by cfftf1().

1037 {
1038  /* System generated locals */
1039  int cc_dim1, cc_offset, ch_dim1, ch_dim2, ch_offset, i_1, i_2;
1040 
1041  /* Local variables */
1042  static int i, k;
1043  static double ci2, ci3, ci4, cr2, cr3, cr4, ti1, ti2, ti3, ti4, tr1,
1044  tr2, tr3, tr4;
1045 
1046  /* Parameter adjustments */
1047  cc_dim1 = *ido;
1048  cc_offset = cc_dim1 * 5 + 1;
1049  cc -= cc_offset;
1050  ch_dim1 = *ido;
1051  ch_dim2 = *l1;
1052  ch_offset = ch_dim1 * (ch_dim2 + 1) + 1;
1053  ch -= ch_offset;
1054  --wa1;
1055  --wa2;
1056  --wa3;
1057 
1058  /* Function Body */
1059  if (*ido != 2) {
1060  goto L102;
1061  }
1062  i_1 = *l1;
1063  for (k = 1; k <= i_1; ++k) {
1064  ti1 = cc[((k << 2) + 1) * cc_dim1 + 2] - cc[((k << 2) + 3) * cc_dim1
1065  + 2];
1066  ti2 = cc[((k << 2) + 1) * cc_dim1 + 2] + cc[((k << 2) + 3) * cc_dim1
1067  + 2];
1068  tr4 = cc[((k << 2) + 2) * cc_dim1 + 2] - cc[((k << 2) + 4) * cc_dim1
1069  + 2];
1070  ti3 = cc[((k << 2) + 2) * cc_dim1 + 2] + cc[((k << 2) + 4) * cc_dim1
1071  + 2];
1072  tr1 = cc[((k << 2) + 1) * cc_dim1 + 1] - cc[((k << 2) + 3) * cc_dim1
1073  + 1];
1074  tr2 = cc[((k << 2) + 1) * cc_dim1 + 1] + cc[((k << 2) + 3) * cc_dim1
1075  + 1];
1076  ti4 = cc[((k << 2) + 4) * cc_dim1 + 1] - cc[((k << 2) + 2) * cc_dim1
1077  + 1];
1078  tr3 = cc[((k << 2) + 2) * cc_dim1 + 1] + cc[((k << 2) + 4) * cc_dim1
1079  + 1];
1080  ch[(k + ch_dim2) * ch_dim1 + 1] = tr2 + tr3;
1081  ch[(k + ch_dim2 * 3) * ch_dim1 + 1] = tr2 - tr3;
1082  ch[(k + ch_dim2) * ch_dim1 + 2] = ti2 + ti3;
1083  ch[(k + ch_dim2 * 3) * ch_dim1 + 2] = ti2 - ti3;
1084  ch[(k + (ch_dim2 << 1)) * ch_dim1 + 1] = tr1 + tr4;
1085  ch[(k + (ch_dim2 << 2)) * ch_dim1 + 1] = tr1 - tr4;
1086  ch[(k + (ch_dim2 << 1)) * ch_dim1 + 2] = ti1 + ti4;
1087  ch[(k + (ch_dim2 << 2)) * ch_dim1 + 2] = ti1 - ti4;
1088 /* L101: */
1089  }
1090  return 0;
1091 L102:
1092  i_1 = *l1;
1093  for (k = 1; k <= i_1; ++k) {
1094  i_2 = *ido;
1095  for (i = 2; i <= i_2; i += 2) {
1096  ti1 = cc[i + ((k << 2) + 1) * cc_dim1] - cc[i + ((k << 2) + 3) *
1097  cc_dim1];
1098  ti2 = cc[i + ((k << 2) + 1) * cc_dim1] + cc[i + ((k << 2) + 3) *
1099  cc_dim1];
1100  ti3 = cc[i + ((k << 2) + 2) * cc_dim1] + cc[i + ((k << 2) + 4) *
1101  cc_dim1];
1102  tr4 = cc[i + ((k << 2) + 2) * cc_dim1] - cc[i + ((k << 2) + 4) *
1103  cc_dim1];
1104  tr1 = cc[i - 1 + ((k << 2) + 1) * cc_dim1] - cc[i - 1 + ((k << 2)
1105  + 3) * cc_dim1];
1106  tr2 = cc[i - 1 + ((k << 2) + 1) * cc_dim1] + cc[i - 1 + ((k << 2)
1107  + 3) * cc_dim1];
1108  ti4 = cc[i - 1 + ((k << 2) + 4) * cc_dim1] - cc[i - 1 + ((k << 2)
1109  + 2) * cc_dim1];
1110  tr3 = cc[i - 1 + ((k << 2) + 2) * cc_dim1] + cc[i - 1 + ((k << 2)
1111  + 4) * cc_dim1];
1112  ch[i - 1 + (k + ch_dim2) * ch_dim1] = tr2 + tr3;
1113  cr3 = tr2 - tr3;
1114  ch[i + (k + ch_dim2) * ch_dim1] = ti2 + ti3;
1115  ci3 = ti2 - ti3;
1116  cr2 = tr1 + tr4;
1117  cr4 = tr1 - tr4;
1118  ci2 = ti1 + ti4;
1119  ci4 = ti1 - ti4;
1120  ch[i - 1 + (k + (ch_dim2 << 1)) * ch_dim1] = wa1[i - 1] * cr2 +
1121  wa1[i] * ci2;
1122  ch[i + (k + (ch_dim2 << 1)) * ch_dim1] = wa1[i - 1] * ci2 - wa1[i]
1123  * cr2;
1124  ch[i - 1 + (k + ch_dim2 * 3) * ch_dim1] = wa2[i - 1] * cr3 + wa2[
1125  i] * ci3;
1126  ch[i + (k + ch_dim2 * 3) * ch_dim1] = wa2[i - 1] * ci3 - wa2[i] *
1127  cr3;
1128  ch[i - 1 + (k + (ch_dim2 << 2)) * ch_dim1] = wa3[i - 1] * cr4 +
1129  wa3[i] * ci4;
1130  ch[i + (k + (ch_dim2 << 2)) * ch_dim1] = wa3[i - 1] * ci4 - wa3[i]
1131  * cr4;
1132 /* L103: */
1133  }
1134 /* L104: */
1135  }
1136  return 0;
1137 } /* passf4_ */

◆ passf5()

int passf5 ( int *  ido,
int *  l1,
double *  cc,
double *  ch,
double *  wa1,
double *  wa2,
double *  wa3,
double *  wa4 
)

Definition at line 1139 of file pub3dfft.C.

Referenced by cfftf1().

1142 {
1143  /* Initialized data */
1144 
1145  static double tr11 = .309016994374947;
1146  static double ti11 = -.951056516295154;
1147  static double tr12 = -.809016994374947;
1148  static double ti12 = -.587785252292473;
1149 
1150  /* System generated locals */
1151  int cc_dim1, cc_offset, ch_dim1, ch_dim2, ch_offset, i_1, i_2;
1152 
1153  /* Local variables */
1154  static int i, k;
1155  static double ci2, ci3, ci4, ci5, di3, di4, di5, di2, cr2, cr3, cr5,
1156  cr4, ti2, ti3, ti4, ti5, dr3, dr4, dr5, dr2, tr2, tr3, tr4, tr5;
1157 
1158  /* Parameter adjustments */
1159  cc_dim1 = *ido;
1160  cc_offset = cc_dim1 * 6 + 1;
1161  cc -= cc_offset;
1162  ch_dim1 = *ido;
1163  ch_dim2 = *l1;
1164  ch_offset = ch_dim1 * (ch_dim2 + 1) + 1;
1165  ch -= ch_offset;
1166  --wa1;
1167  --wa2;
1168  --wa3;
1169  --wa4;
1170 
1171  /* Function Body */
1172  if (*ido != 2) {
1173  goto L102;
1174  }
1175  i_1 = *l1;
1176  for (k = 1; k <= i_1; ++k) {
1177  ti5 = cc[(k * 5 + 2) * cc_dim1 + 2] - cc[(k * 5 + 5) * cc_dim1 + 2];
1178  ti2 = cc[(k * 5 + 2) * cc_dim1 + 2] + cc[(k * 5 + 5) * cc_dim1 + 2];
1179  ti4 = cc[(k * 5 + 3) * cc_dim1 + 2] - cc[(k * 5 + 4) * cc_dim1 + 2];
1180  ti3 = cc[(k * 5 + 3) * cc_dim1 + 2] + cc[(k * 5 + 4) * cc_dim1 + 2];
1181  tr5 = cc[(k * 5 + 2) * cc_dim1 + 1] - cc[(k * 5 + 5) * cc_dim1 + 1];
1182  tr2 = cc[(k * 5 + 2) * cc_dim1 + 1] + cc[(k * 5 + 5) * cc_dim1 + 1];
1183  tr4 = cc[(k * 5 + 3) * cc_dim1 + 1] - cc[(k * 5 + 4) * cc_dim1 + 1];
1184  tr3 = cc[(k * 5 + 3) * cc_dim1 + 1] + cc[(k * 5 + 4) * cc_dim1 + 1];
1185  ch[(k + ch_dim2) * ch_dim1 + 1] = cc[(k * 5 + 1) * cc_dim1 + 1] + tr2
1186  + tr3;
1187  ch[(k + ch_dim2) * ch_dim1 + 2] = cc[(k * 5 + 1) * cc_dim1 + 2] + ti2
1188  + ti3;
1189  cr2 = cc[(k * 5 + 1) * cc_dim1 + 1] + tr11 * tr2 + tr12 * tr3;
1190  ci2 = cc[(k * 5 + 1) * cc_dim1 + 2] + tr11 * ti2 + tr12 * ti3;
1191  cr3 = cc[(k * 5 + 1) * cc_dim1 + 1] + tr12 * tr2 + tr11 * tr3;
1192  ci3 = cc[(k * 5 + 1) * cc_dim1 + 2] + tr12 * ti2 + tr11 * ti3;
1193  cr5 = ti11 * tr5 + ti12 * tr4;
1194  ci5 = ti11 * ti5 + ti12 * ti4;
1195  cr4 = ti12 * tr5 - ti11 * tr4;
1196  ci4 = ti12 * ti5 - ti11 * ti4;
1197  ch[(k + (ch_dim2 << 1)) * ch_dim1 + 1] = cr2 - ci5;
1198  ch[(k + ch_dim2 * 5) * ch_dim1 + 1] = cr2 + ci5;
1199  ch[(k + (ch_dim2 << 1)) * ch_dim1 + 2] = ci2 + cr5;
1200  ch[(k + ch_dim2 * 3) * ch_dim1 + 2] = ci3 + cr4;
1201  ch[(k + ch_dim2 * 3) * ch_dim1 + 1] = cr3 - ci4;
1202  ch[(k + (ch_dim2 << 2)) * ch_dim1 + 1] = cr3 + ci4;
1203  ch[(k + (ch_dim2 << 2)) * ch_dim1 + 2] = ci3 - cr4;
1204  ch[(k + ch_dim2 * 5) * ch_dim1 + 2] = ci2 - cr5;
1205 /* L101: */
1206  }
1207  return 0;
1208 L102:
1209  i_1 = *l1;
1210  for (k = 1; k <= i_1; ++k) {
1211  i_2 = *ido;
1212  for (i = 2; i <= i_2; i += 2) {
1213  ti5 = cc[i + (k * 5 + 2) * cc_dim1] - cc[i + (k * 5 + 5) *
1214  cc_dim1];
1215  ti2 = cc[i + (k * 5 + 2) * cc_dim1] + cc[i + (k * 5 + 5) *
1216  cc_dim1];
1217  ti4 = cc[i + (k * 5 + 3) * cc_dim1] - cc[i + (k * 5 + 4) *
1218  cc_dim1];
1219  ti3 = cc[i + (k * 5 + 3) * cc_dim1] + cc[i + (k * 5 + 4) *
1220  cc_dim1];
1221  tr5 = cc[i - 1 + (k * 5 + 2) * cc_dim1] - cc[i - 1 + (k * 5 + 5) *
1222  cc_dim1];
1223  tr2 = cc[i - 1 + (k * 5 + 2) * cc_dim1] + cc[i - 1 + (k * 5 + 5) *
1224  cc_dim1];
1225  tr4 = cc[i - 1 + (k * 5 + 3) * cc_dim1] - cc[i - 1 + (k * 5 + 4) *
1226  cc_dim1];
1227  tr3 = cc[i - 1 + (k * 5 + 3) * cc_dim1] + cc[i - 1 + (k * 5 + 4) *
1228  cc_dim1];
1229  ch[i - 1 + (k + ch_dim2) * ch_dim1] = cc[i - 1 + (k * 5 + 1) *
1230  cc_dim1] + tr2 + tr3;
1231  ch[i + (k + ch_dim2) * ch_dim1] = cc[i + (k * 5 + 1) * cc_dim1] +
1232  ti2 + ti3;
1233  cr2 = cc[i - 1 + (k * 5 + 1) * cc_dim1] + tr11 * tr2 + tr12 * tr3;
1234 
1235  ci2 = cc[i + (k * 5 + 1) * cc_dim1] + tr11 * ti2 + tr12 * ti3;
1236  cr3 = cc[i - 1 + (k * 5 + 1) * cc_dim1] + tr12 * tr2 + tr11 * tr3;
1237 
1238  ci3 = cc[i + (k * 5 + 1) * cc_dim1] + tr12 * ti2 + tr11 * ti3;
1239  cr5 = ti11 * tr5 + ti12 * tr4;
1240  ci5 = ti11 * ti5 + ti12 * ti4;
1241  cr4 = ti12 * tr5 - ti11 * tr4;
1242  ci4 = ti12 * ti5 - ti11 * ti4;
1243  dr3 = cr3 - ci4;
1244  dr4 = cr3 + ci4;
1245  di3 = ci3 + cr4;
1246  di4 = ci3 - cr4;
1247  dr5 = cr2 + ci5;
1248  dr2 = cr2 - ci5;
1249  di5 = ci2 - cr5;
1250  di2 = ci2 + cr5;
1251  ch[i - 1 + (k + (ch_dim2 << 1)) * ch_dim1] = wa1[i - 1] * dr2 +
1252  wa1[i] * di2;
1253  ch[i + (k + (ch_dim2 << 1)) * ch_dim1] = wa1[i - 1] * di2 - wa1[i]
1254  * dr2;
1255  ch[i - 1 + (k + ch_dim2 * 3) * ch_dim1] = wa2[i - 1] * dr3 + wa2[
1256  i] * di3;
1257  ch[i + (k + ch_dim2 * 3) * ch_dim1] = wa2[i - 1] * di3 - wa2[i] *
1258  dr3;
1259  ch[i - 1 + (k + (ch_dim2 << 2)) * ch_dim1] = wa3[i - 1] * dr4 +
1260  wa3[i] * di4;
1261  ch[i + (k + (ch_dim2 << 2)) * ch_dim1] = wa3[i - 1] * di4 - wa3[i]
1262  * dr4;
1263  ch[i - 1 + (k + ch_dim2 * 5) * ch_dim1] = wa4[i - 1] * dr5 + wa4[
1264  i] * di5;
1265  ch[i + (k + ch_dim2 * 5) * ch_dim1] = wa4[i - 1] * di5 - wa4[i] *
1266  dr5;
1267 /* L103: */
1268  }
1269 /* L104: */
1270  }
1271  return 0;
1272 } /* passf5_ */

◆ pubd3di()

int pubd3di ( int  n1,
int  n2,
int  n3,
double *  table,
int  ntable 
)

Definition at line 1851 of file pub3dfft.C.

References pubz3di().

1851  {
1852 
1853  int n1over2;
1854 
1855  n1over2 = n1 / 2;
1856  return pubz3di(&n1over2,&n2,&n3,table,&ntable);
1857 
1858 } /* pubd3di */
pubz3di
int pubz3di(int *n1, int *n2, int *n3, double *table, int *ntable)
Definition: pub3dfft.C:1706

◆ pubdz3d()

int pubdz3d ( int  isign,
int  n1,
int  n2,
int  n3,
double *  w,
int  ld1,
int  ld2,
double *  table,
int  ntable,
double *  work 
)

Definition at line 1863 of file pub3dfft.C.

References M_PI, and pubz3d().

1865  {
1866 
1867  int n1over2, ld1over2, rval;
1868  int i, j, j2, k, k2, i1, i2, imax;
1869  double *data, *data2;
1870  double TwoPiOverN, tmp1r, tmp1i, tmp2r, tmp2i;
1871 
1872  /* complex transform */
1873  n1over2 = n1 / 2;
1874  ld1over2 = ld1 / 2;
1875  rval = pubz3d(&isign, &n1over2, &n2, &n3, (doublecomplex*) w,
1876  &ld1over2, &ld2, table, &ntable, (doublecomplex*) work);
1877 
1878  /* rearrange data */
1879  TwoPiOverN = isign * 2.0 * M_PI / n1;
1880  imax = n1/4+1;
1881  for ( i=0; i<imax; ++i ) {
1882  work[2*i] = cos(i * TwoPiOverN);
1883  work[2*i+1] = sin(i * TwoPiOverN);
1884  }
1885  for ( k=0; k<n3; ++k ) {
1886  for ( j=0; j<n2; ++j ) {
1887  data = w + ld1*(ld2*k + j);
1888  data[n1] = data[0];
1889  data[n1+1] = data[1];
1890  }
1891  }
1892  for ( k=0; k<n3; ++k ) {
1893  k2 = k?(n3-k):0;
1894  for ( j=0; j<n2; ++j ) {
1895  j2 = j?(n2-j):0;
1896  data = w + ld1*(ld2*k + j);
1897  data2 = w + ld1*(ld2*k2 + j2);
1898  imax = n1/4;
1899  if ( (n1/2) & 1 ) imax += 1;
1900  else {
1901  if ( (2*j<n2) || (2*j==n2 && 2*k<=n3) ) imax +=1;
1902  if ( j==0 && 2*k>n3 ) imax -=1;
1903  }
1904  for ( i=0; i<imax; ++i ) {
1905  i1 = 2*i; i2 = n1-i1;
1906  tmp1r = data[i1] - data2[i2];
1907  tmp1i = data[i1+1] + data2[i2+1];
1908  tmp2r = tmp1r * work[i1+1] + tmp1i * work[i1];
1909  tmp2i = tmp1i * work[i1+1] - tmp1r * work[i1];
1910  tmp1r = data[i1] + data2[i2];
1911  tmp1i = data[i1+1] - data2[i2+1];
1912  data[i1] = 0.5 * ( tmp1r + tmp2r );
1913  data[i1+1] = 0.5 * ( tmp1i + tmp2i );
1914  data2[i2] = 0.5 * ( tmp1r - tmp2r );
1915  data2[i2+1] = 0.5 * ( tmp2i - tmp1i );
1916  }
1917  }
1918  }
1919 
1920  return rval;
1921 
1922 } /* pubdz3d */
M_PI
#define M_PI
Definition: pub3dfft.C:15
pubz3d
int pubz3d(int *isign, int *n1, int *n2, int *n3, doublecomplex *w, int *ld1, int *ld2, double *table, int *ntable, doublecomplex *work)
Definition: pub3dfft.C:1729

◆ pubz3d()

int pubz3d ( int *  isign,
int *  n1,
int *  n2,
int *  n3,
doublecomplex *  w,
int *  ld1,
int *  ld2,
double *  table,
int *  ntable,
doublecomplex *  work 
)

Definition at line 1729 of file pub3dfft.C.

References cfftb(), cfftf(), doublecomplex::i, and doublecomplex::r.

Referenced by pubdz3d(), and pubzd3d().

1732 {
1733  /* System generated locals */
1734  int w_dim1, w_dim2, w_offset, table_dim1, table_offset, i_1, i_2, i_3,
1735  i_4, i_5;
1736 
1737  /* Local variables */
1738  static int i, j, k;
1739 
1740  /* Parameter adjustments */
1741  w_dim1 = *ld1;
1742  w_dim2 = *ld2;
1743  w_offset = w_dim1 * (w_dim2 + 1) + 1;
1744  w -= w_offset;
1745  table_dim1 = *ntable;
1746  table_offset = table_dim1 + 1;
1747  table -= table_offset;
1748  --work;
1749 
1750  /* Function Body */
1751 /* ntable should be 4*max(n1,n2,n3) +15 */
1752 /* nwork should be max(n1,n2,n3) */
1753 
1754 /* transform along X first ... */
1755 
1756  i_1 = *n3;
1757  for (k = 1; k <= i_1; ++k) {
1758  i_2 = *n2;
1759  for (j = 1; j <= i_2; ++j) {
1760  i_3 = *n1;
1761  for (i = 1; i <= i_3; ++i) {
1762  i_4 = i;
1763  i_5 = i + (j + k * w_dim2) * w_dim1;
1764  work[i_4].r = w[i_5].r, work[i_4].i = w[i_5].i;
1765 /* L70: */
1766  }
1767  if (*isign == -1) {
1768  cfftf(n1, (double*) &work[1], &table[table_dim1 + 1]);
1769  }
1770  if (*isign == 1) {
1771  cfftb(n1, (double*) &work[1], &table[table_dim1 + 1]);
1772  }
1773  i_3 = *n1;
1774  for (i = 1; i <= i_3; ++i) {
1775  i_4 = i + (j + k * w_dim2) * w_dim1;
1776  i_5 = i;
1777  w[i_4].r = work[i_5].r, w[i_4].i = work[i_5].i;
1778 /* L80: */
1779  }
1780 /* L90: */
1781  }
1782 /* L100: */
1783  }
1784 
1785 /* transform along Y then ... */
1786 
1787  i_1 = *n3;
1788  for (k = 1; k <= i_1; ++k) {
1789  i_2 = *n1;
1790  for (i = 1; i <= i_2; ++i) {
1791  i_3 = *n2;
1792  for (j = 1; j <= i_3; ++j) {
1793  i_4 = j;
1794  i_5 = i + (j + k * w_dim2) * w_dim1;
1795  work[i_4].r = w[i_5].r, work[i_4].i = w[i_5].i;
1796 /* L170: */
1797  }
1798  if (*isign == -1) {
1799  cfftf(n2, (double*) &work[1], &table[(table_dim1 << 1) + 1]);
1800  }
1801  if (*isign == 1) {
1802  cfftb(n2, (double*) &work[1], &table[(table_dim1 << 1) + 1]);
1803  }
1804  i_3 = *n2;
1805  for (j = 1; j <= i_3; ++j) {
1806  i_4 = i + (j + k * w_dim2) * w_dim1;
1807  i_5 = j;
1808  w[i_4].r = work[i_5].r, w[i_4].i = work[i_5].i;
1809 /* L180: */
1810  }
1811 /* L190: */
1812  }
1813 /* L200: */
1814  }
1815 
1816 /* transform along Z finally ... */
1817 
1818  i_1 = *n1;
1819  for (i = 1; i <= i_1; ++i) {
1820  i_2 = *n2;
1821  for (j = 1; j <= i_2; ++j) {
1822  i_3 = *n3;
1823  for (k = 1; k <= i_3; ++k) {
1824  i_4 = k;
1825  i_5 = i + (j + k * w_dim2) * w_dim1;
1826  work[i_4].r = w[i_5].r, work[i_4].i = w[i_5].i;
1827 /* L270: */
1828  }
1829  if (*isign == -1) {
1830  cfftf(n3, (double*) &work[1], &table[table_dim1 * 3 + 1]);
1831  }
1832  if (*isign == 1) {
1833  cfftb(n3, (double*) &work[1], &table[table_dim1 * 3 + 1]);
1834  }
1835  i_3 = *n3;
1836  for (k = 1; k <= i_3; ++k) {
1837  i_4 = i + (j + k * w_dim2) * w_dim1;
1838  i_5 = k;
1839  w[i_4].r = work[i_5].r, w[i_4].i = work[i_5].i;
1840 /* L280: */
1841  }
1842 /* L290: */
1843  }
1844 /* L300: */
1845  }
1846  return 0;
1847 } /* pubz3d_ */
cfftb
int cfftb(int *n, double *c, double *wsave)
Definition: pub3dfft.C:1401
doublecomplex::r
double r
Definition: pub3dfft.h:10
cfftf
int cfftf(int *n, double *c, double *wsave)
Definition: pub3dfft.C:1544
doublecomplex::i
double i
Definition: pub3dfft.h:10

◆ pubz3di()

int pubz3di ( int *  n1,
int *  n2,
int *  n3,
double *  table,
int *  ntable 
)

Definition at line 1706 of file pub3dfft.C.

References cffti().

Referenced by pubd3di().

1708 {
1709  /* System generated locals */
1710  int table_dim1, table_offset;
1711 
1712  /* Local variables */
1713 
1714  /* Parameter adjustments */
1715  table_dim1 = *ntable;
1716  table_offset = table_dim1 + 1;
1717  table -= table_offset;
1718 
1719  /* Function Body */
1720 /* ntable should be 4*max(n1,n2,n3) +15 */
1721  cffti(n1, &table[table_dim1 + 1]);
1722  cffti(n2, &table[(table_dim1 << 1) + 1]);
1723  cffti(n3, &table[table_dim1 * 3 + 1]);
1724  return 0;
1725 } /* pubz3di_ */
cffti
int cffti(int *n, double *wsave)
Definition: pub3dfft.C:1676

◆ pubzd3d()

int pubzd3d ( int  isign,
int  n1,
int  n2,
int  n3,
double *  w,
int  ld1,
int  ld2,
double *  table,
int  ntable,
double *  work 
)

Definition at line 1927 of file pub3dfft.C.

References M_PI, and pubz3d().

1929  {
1930 
1931  int n1over2, ld1over2;
1932  int i, j, j2, k, k2, i1, i2, imax;
1933  double *data, *data2;
1934  double TwoPiOverN, tmp1r, tmp1i, tmp2r, tmp2i;
1935 
1936  /* rearrange data */
1937  TwoPiOverN = isign * 2.0 * M_PI / n1;
1938  imax = n1/4+1;
1939  for ( i=0; i<imax; ++i ) {
1940  work[2*i] = -cos(i * TwoPiOverN);
1941  work[2*i+1] = -sin(i * TwoPiOverN);
1942  }
1943  for ( k=0; k<n3; ++k ) {
1944  k2 = k?(n3-k):0;
1945  for ( j=0; j<n2; ++j ) {
1946  j2 = j?(n2-j):0;
1947  data = w + ld1*(ld2*k + j);
1948  data2 = w + ld1*(ld2*k2 + j2);
1949  imax = n1/4;
1950  if ( (n1/2) & 1 ) imax += 1;
1951  else {
1952  if ( (2*j<n2) || (2*j==n2 && 2*k<=n3) ) imax +=1;
1953  if ( j==0 && 2*k>n3 ) imax -=1;
1954  }
1955  for ( i=0; i<imax; ++i ) {
1956  i1 = 2*i; i2 = n1-i1;
1957  tmp1r = data[i1] - data2[i2];
1958  tmp1i = data[i1+1] + data2[i2+1];
1959  tmp2r = tmp1r * work[i1+1] + tmp1i * work[i1];
1960  tmp2i = tmp1i * work[i1+1] - tmp1r * work[i1];
1961  tmp1r = data[i1] + data2[i2];
1962  tmp1i = data[i1+1] - data2[i2+1];
1963  data[i1] = tmp1r + tmp2r;
1964  data[i1+1] = tmp1i + tmp2i;
1965  data2[i2] = tmp1r - tmp2r;
1966  data2[i2+1] = tmp2i - tmp1i;
1967  }
1968  }
1969  }
1970 
1971  /* complex transform */
1972  n1over2 = n1 / 2;
1973  ld1over2 = ld1 / 2;
1974  return pubz3d(&isign, &n1over2, &n2, &n3, (doublecomplex*) w,
1975  &ld1over2, &ld2, table, &ntable, (doublecomplex*) work);
1976 
1977 } /* pubzd3d */
M_PI
#define M_PI
Definition: pub3dfft.C:15
pubz3d
int pubz3d(int *isign, int *n1, int *n2, int *n3, doublecomplex *w, int *ld1, int *ld2, double *table, int *ntable, doublecomplex *work)
Definition: pub3dfft.C:1729
Generated by   doxygen 1.8.14