home *** CD-ROM | disk | FTP | other *** search
/ Frostbyte's 1980s DOS Shareware Collection / floppyshareware.zip / floppyshareware / DOOG / PCSSP1.ZIP / ITRPAPSM.ZIP / FORIT.FOR < prev    next >
Text File  |  1985-11-29  |  3KB  |  95 lines
  1. C
  2. C     ..................................................................
  3. C
  4. C        SUBROUTINE FORIT
  5. C
  6. C        PURPOSE
  7. C           FOURIER ANALYSIS OF A PERIODICALLY TABULATED FUNCTION.
  8. C           COMPUTES THE COEFFICIENTS OF THE DESIRED NUMBER OF TERMS
  9. C           IN THE FOURIER SERIES F(X)=A(0)+SUM(A(K)COS KX+B(K)SIN KX)
  10. C           WHERE K=1,2,...,M TO APPROXIMATE A GIVEN SET OF
  11. C           PERIODICALLY TABULATED VALUES OF A FUNCTION.
  12. C
  13. C        USAGE
  14. C           CALL FORIT(FNT,N,M,A,B,IER)
  15. C
  16. C        DESCRIPTION OF PARAMETERS
  17. C           FNT-VECTOR OF TABULATED FUNCTION VALUES OF LENGTH 2N+1
  18. C           N  -DEFINES THE INTERVAL SUCH THAT 2N+1 POINTS ARE TAKEN
  19. C               OVER THE INTERVAL (0,2PI). THE SPACING IS THUS 2PI/2N+1
  20. C           M  -MAXIMUM ORDER OF HARMONICS TO BE FITTED
  21. C           A  -RESULTANT VECTOR OF FOURIER COSINE COEFFICIENTS OF
  22. C               LENGTH M+1
  23. C               A SUB 0, A SUB 1,..., A SUB M
  24. C           B  -RESULTANT VECTOR OF FOURIER SINE COEFFICIENTS OF
  25. C               LENGTH M+1
  26. C               B SUB 0, B SUB 1,..., B SUB M
  27. C           IER-RESULTANT ERROR CODE WHERE
  28. C               IER=0  NO ERROR
  29. C               IER=1  N NOT GREATER OR EQUAL TO M
  30. C               IER=2  M LESS THAN 0
  31. C
  32. C        REMARKS
  33. C           M MUST BE GREATER THAN OR EQUAL TO ZERO
  34. C           N MUST BE GREATER THAN OR EQUAL TO M
  35. C           THE FIRST ELEMENT OF VECTOR B IS ZERO IN ALL CASES
  36. C
  37. C        SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED
  38. C           NONE
  39. C
  40. C        METHOD
  41. C           USES RECURSIVE TECHNIQUE DESCRIBED IN A. RALSTON, H. WILF,
  42. C           'MATHEMATICAL METHODS FOR DIGITAL COMPUTERS', JOHN WILEY
  43. C           AND SONS, NEW YORK, 1960, CHAPTER 24. THE METHOD OF INDEXING
  44. C           THROUGH THE PROCEDURE HAS BEEN MODIFIED TO SIMPLIFY THE
  45. C           COMPUTATION.
  46. C
  47. C     ..................................................................
  48. C
  49.       SUBROUTINE FORIT(FNT,N,M,A,B,IER)
  50.       DIMENSION A(1),B(1),FNT(1)
  51. C
  52. C        CHECK FOR PARAMETER ERRORS
  53. C
  54.       IER=0
  55.    20 IF(M) 30,40,40
  56.    30 IER=2
  57.       RETURN
  58.    40 IF(M-N) 60,60,50
  59.    50 IER=1
  60.       RETURN
  61. C
  62. C        COMPUTE AND PRESET CONSTANTS
  63. C
  64.    60 AN=N
  65.       COEF=2.0/(2.0*AN+1.0)
  66.       CONST=3.141593*COEF
  67.       S1=SIN(CONST)
  68.       C1=COS(CONST)
  69.       C=1.0
  70.       S=0.0
  71.       J=1
  72.       FNTZ=FNT(1)
  73.    70 U2=0.0
  74.       U1=0.0
  75.       I=2*N+1
  76. C
  77. C        FORM FOURIER COEFFICIENTS RECURSIVELY
  78. C
  79.    75 U0=FNT(I)+2.0*C*U1-U2
  80.       U2=U1
  81.       U1=U0
  82.       I=I-1
  83.       IF(I-1) 80,80,75
  84.    80 A(J)=COEF*(FNTZ+C*U1-U2)
  85.       B(J)=COEF*S*U1
  86.       IF(J-(M+1)) 90,100,100
  87.    90 Q=C1*C-S1*S
  88.       S=C1*S+S1*C
  89.       C=Q
  90.       J=J+1
  91.       GO TO 70
  92.   100 A(1)=A(1)*0.5
  93.       RETURN
  94.       END
  95.