website/archive/binaries/mac/src/disembl/Tisean_3.0.1/source_f/slatec/radb3.f

   1 *DECK RADB3
   2       SUBROUTINE RADB3 (IDO, L1, CC, CH, WA1, WA2)
   3 C***BEGIN PROLOGUE  RADB3
   4 C***SUBSIDIARY
   5 C***PURPOSE  Calculate the fast Fourier transform of subvectors of
   6 C            length three.
   7 C***LIBRARY   SLATEC (FFTPACK)
   8 C***TYPE      SINGLE PRECISION (RADB3-S)
   9 C***AUTHOR  Swarztrauber, P. N., (NCAR)
  10 C***ROUTINES CALLED  (NONE)
  11 C***REVISION HISTORY  (YYMMDD)
  12 C   790601  DATE WRITTEN
  13 C   830401  Modified to use SLATEC library source file format.
  14 C   860115  Modified by Ron Boisvert to adhere to Fortran 77 by
  15 C           (a) changing dummy array size declarations (1) to (*),
  16 C           (b) changing definition of variable TAUI by using
  17 C               FORTRAN intrinsic function SQRT instead of a DATA
  18 C               statement.
  19 C   881128  Modified by Dick Valent to meet prologue standards.
  20 C   890831  Modified array declarations.  (WRB)
  21 C   891214  Prologue converted to Version 4.0 format.  (BAB)
  22 C   900402  Added TYPE section.  (WRB)
  23 C***END PROLOGUE  RADB3
  24       DIMENSION CC(IDO,3,*), CH(IDO,L1,3), WA1(*), WA2(*)
  25 C***FIRST EXECUTABLE STATEMENT  RADB3
  26       TAUR = -.5
  27       TAUI = .5*SQRT(3.)
  28       DO 101 K=1,L1
  29          TR2 = CC(IDO,2,K)+CC(IDO,2,K)
  30          CR2 = CC(1,1,K)+TAUR*TR2
  31          CH(1,K,1) = CC(1,1,K)+TR2
  32          CI3 = TAUI*(CC(1,3,K)+CC(1,3,K))
  33          CH(1,K,2) = CR2-CI3
  34          CH(1,K,3) = CR2+CI3
  35   101 CONTINUE
  36       IF (IDO .EQ. 1) RETURN
  37       IDP2 = IDO+2
  38       IF((IDO-1)/2.LT.L1) GO TO 104
  39       DO 103 K=1,L1
  40 CDIR$ IVDEP
  41          DO 102 I=3,IDO,2
  42             IC = IDP2-I
  43             TR2 = CC(I-1,3,K)+CC(IC-1,2,K)
  44             CR2 = CC(I-1,1,K)+TAUR*TR2
  45             CH(I-1,K,1) = CC(I-1,1,K)+TR2
  46             TI2 = CC(I,3,K)-CC(IC,2,K)
  47             CI2 = CC(I,1,K)+TAUR*TI2
  48             CH(I,K,1) = CC(I,1,K)+TI2
  49             CR3 = TAUI*(CC(I-1,3,K)-CC(IC-1,2,K))
  50             CI3 = TAUI*(CC(I,3,K)+CC(IC,2,K))
  51             DR2 = CR2-CI3
  52             DR3 = CR2+CI3
  53             DI2 = CI2+CR3
  54             DI3 = CI2-CR3
  55             CH(I-1,K,2) = WA1(I-2)*DR2-WA1(I-1)*DI2
  56             CH(I,K,2) = WA1(I-2)*DI2+WA1(I-1)*DR2
  57             CH(I-1,K,3) = WA2(I-2)*DR3-WA2(I-1)*DI3
  58             CH(I,K,3) = WA2(I-2)*DI3+WA2(I-1)*DR3
  59   102    CONTINUE
  60   103 CONTINUE
  61       RETURN
  62   104 DO 106 I=3,IDO,2
  63          IC = IDP2-I
  64 CDIR$ IVDEP
  65          DO 105 K=1,L1
  66             TR2 = CC(I-1,3,K)+CC(IC-1,2,K)
  67             CR2 = CC(I-1,1,K)+TAUR*TR2
  68             CH(I-1,K,1) = CC(I-1,1,K)+TR2
  69             TI2 = CC(I,3,K)-CC(IC,2,K)
  70             CI2 = CC(I,1,K)+TAUR*TI2
  71             CH(I,K,1) = CC(I,1,K)+TI2
  72             CR3 = TAUI*(CC(I-1,3,K)-CC(IC-1,2,K))
  73             CI3 = TAUI*(CC(I,3,K)+CC(IC,2,K))
  74             DR2 = CR2-CI3
  75             DR3 = CR2+CI3
  76             DI2 = CI2+CR3
  77             DI3 = CI2-CR3
  78             CH(I-1,K,2) = WA1(I-2)*DR2-WA1(I-1)*DI2
  79             CH(I,K,2) = WA1(I-2)*DI2+WA1(I-1)*DR2
  80             CH(I-1,K,3) = WA2(I-2)*DR3-WA2(I-1)*DI3
  81             CH(I,K,3) = WA2(I-2)*DI3+WA2(I-1)*DR3
  82   105    CONTINUE
  83   106 CONTINUE
  84       RETURN
  85       END