share/fftpack5/fortran/mrftb1.f

   1 CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
   2 C
   3 C   FFTPACK 5.0
   4 C
   5 C   Authors:  Paul N. Swarztrauber and Richard A. Valent
   6 C
   7 C   $Id$
   8 C
   9 CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  10
  11       SUBROUTINE MRFTB1 (M,IM,N,IN,C,CH,WA,FAC)
  12       REAL       CH(M,*), C(IN,*), WA(N) ,FAC(15)
  13 C
  14       NF = FAC(2)
  15       NA = 0
  16       DO 10 K1=1,NF
  17       IP = FAC(K1+2)
  18       NA = 1-NA
  19       IF(IP .LE. 5) GO TO 10
  20       IF(K1 .EQ. NF) GO TO 10
  21       NA = 1-NA
  22    10 CONTINUE
  23       HALF = .5
  24       HALFM = -.5
  25       MODN = MOD(N,2)
  26       NL = N-2
  27       IF(MODN .NE. 0) NL = N-1
  28       IF (NA .EQ. 0) GO TO 120
  29       M2 = 1-IM
  30       DO 117 I=1,M
  31       M2 = M2+IM
  32       CH(I,1) = C(M2,1)
  33       CH(I,N) = C(M2,N)
  34   117 CONTINUE
  35       DO 118 J=2,NL,2
  36       M2 = 1-IM
  37       DO 118 I=1,M
  38          M2 = M2+IM
  39          CH(I,J) = HALF*C(M2,J)
  40          CH(I,J+1) = HALFM*C(M2,J+1)
  41   118 CONTINUE
  42       GO TO 124
  43   120 continue
  44       DO 122 J=2,NL,2
  45       M2 = 1-IM
  46       DO 122 I=1,M
  47          M2 = M2+IM
  48          C(M2,J) = HALF*C(M2,J)
  49          C(M2,J+1) = HALFM*C(M2,J+1)
  50   122 CONTINUE
  51   124 L1 = 1
  52       IW = 1
  53       DO 116 K1=1,NF
  54          IP = FAC(K1+2)
  55          L2 = IP*L1
  56          IDO = N/L2
  57          IDL1 = IDO*L1
  58          IF (IP .NE. 4) GO TO 103
  59          IX2 = IW+IDO
  60          IX3 = IX2+IDO
  61          IF (NA .NE. 0) GO TO 101
  62          CALL MRADB4 (M,IDO,L1,C,IM,IN,CH,1,M,WA(IW),WA(IX2),WA(IX3))
  63          GO TO 102
  64   101    CALL MRADB4 (M,IDO,L1,CH,1,M,C,IM,IN,WA(IW),WA(IX2),WA(IX3))
  65   102    NA = 1-NA
  66          GO TO 115
  67   103    IF (IP .NE. 2) GO TO 106
  68          IF (NA .NE. 0) GO TO 104
  69          CALL MRADB2 (M,IDO,L1,C,IM,IN,CH,1,M,WA(IW))
  70          GO TO 105
  71   104    CALL MRADB2 (M,IDO,L1,CH,1,M,C,IM,IN,WA(IW))
  72   105    NA = 1-NA
  73          GO TO 115
  74   106    IF (IP .NE. 3) GO TO 109
  75          IX2 = IW+IDO
  76          IF (NA .NE. 0) GO TO 107
  77          CALL MRADB3 (M,IDO,L1,C,IM,IN,CH,1,M,WA(IW),WA(IX2))
  78          GO TO 108
  79   107    CALL MRADB3 (M,IDO,L1,CH,1,M,C,IM,IN,WA(IW),WA(IX2))
  80   108    NA = 1-NA
  81          GO TO 115
  82   109    IF (IP .NE. 5) GO TO 112
  83          IX2 = IW+IDO
  84          IX3 = IX2+IDO
  85          IX4 = IX3+IDO
  86          IF (NA .NE. 0) GO TO 110
  87          CALL MRADB5 (M,IDO,L1,C,IM,IN,CH,1,M,WA(IW),WA(IX2),
  88      1                  WA(IX3),WA(IX4))
  89          GO TO 111
  90   110    CALL MRADB5 (M,IDO,L1,CH,1,M,C,IM,IN,WA(IW),WA(IX2),
  91      1                  WA(IX3),WA(IX4))
  92   111    NA = 1-NA
  93          GO TO 115
  94   112    IF (NA .NE. 0) GO TO 113
  95          CALL MRADBG (M,IDO,IP,L1,IDL1,C,C,C,IM,IN,CH,CH,1,
  96      1                            M,WA(IW))
  97          GO TO 114
  98   113    CALL MRADBG (M,IDO,IP,L1,IDL1,CH,CH,CH,1,M,C,C,IM,
  99      1                           IN,WA(IW))
 100   114    IF (IDO .EQ. 1) NA = 1-NA
 101   115    L1 = L2
 102          IW = IW+(IP-1)*IDO
 103   116 CONTINUE
 104       RETURN
 105       END