cosqb1.f 1.9 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657
  1. *DECK COSQB1
  2. SUBROUTINE COSQB1 (N, X, W, XH)
  3. C***BEGIN PROLOGUE COSQB1
  4. C***SUBSIDIARY
  5. C***PURPOSE Compute the unnormalized inverse of COSQF1.
  6. C***LIBRARY SLATEC (FFTPACK)
  7. C***CATEGORY J1A3
  8. C***TYPE SINGLE PRECISION (COSQB1-S)
  9. C***KEYWORDS FFTPACK, FOURIER TRANSFORM
  10. C***AUTHOR Swarztrauber, P. N., (NCAR)
  11. C***DESCRIPTION
  12. C
  13. C Subroutine COSQB1 computes the fast Fourier transform of quarter
  14. C wave data. That is, COSQB1 computes a sequence from its
  15. C representation in terms of a cosine series with odd wave numbers.
  16. C The transform is defined below at output parameter X.
  17. C
  18. C***REFERENCES P. N. Swarztrauber, Vectorizing the FFTs, in Parallel
  19. C Computations (G. Rodrigue, ed.), Academic Press,
  20. C 1982, pp. 51-83.
  21. C***ROUTINES CALLED RFFTB
  22. C***REVISION HISTORY (YYMMDD)
  23. C 790601 DATE WRITTEN
  24. C 830401 Modified to use SLATEC library source file format.
  25. C 860115 Modified by Ron Boisvert to adhere to Fortran 77 by
  26. C changing dummy array size declarations (1) to (*).
  27. C 881128 Modified by Dick Valent to meet prologue standards.
  28. C 891214 Prologue converted to Version 4.0 format. (BAB)
  29. C 920501 Reformatted the REFERENCES section. (WRB)
  30. C***END PROLOGUE COSQB1
  31. DIMENSION X(*), W(*), XH(*)
  32. C***FIRST EXECUTABLE STATEMENT COSQB1
  33. NS2 = (N+1)/2
  34. NP2 = N+2
  35. DO 101 I=3,N,2
  36. XIM1 = X(I-1)+X(I)
  37. X(I) = X(I)-X(I-1)
  38. X(I-1) = XIM1
  39. 101 CONTINUE
  40. X(1) = X(1)+X(1)
  41. MODN = MOD(N,2)
  42. IF (MODN .EQ. 0) X(N) = X(N)+X(N)
  43. CALL RFFTB (N,X,XH)
  44. DO 102 K=2,NS2
  45. KC = NP2-K
  46. XH(K) = W(K-1)*X(KC)+W(KC-1)*X(K)
  47. XH(KC) = W(K-1)*X(K)-W(KC-1)*X(KC)
  48. 102 CONTINUE
  49. IF (MODN .EQ. 0) X(NS2+1) = W(NS2)*(X(NS2+1)+X(NS2+1))
  50. DO 103 K=2,NS2
  51. KC = NP2-K
  52. X(K) = XH(K)+XH(KC)
  53. X(KC) = XH(K)-XH(KC)
  54. 103 CONTINUE
  55. X(1) = X(1)+X(1)
  56. RETURN
  57. END