C....*...1.........2.........3.........4.........5.........6.........7.*.......8
C     FASTF   1/28/71
C
C     PURPOSE
C     COMPUTES ALL FOURIER COEFFICIENTS AND PERIODOGRAM ORDINATES FOR A
C     REAL TIME SERIES OF LENGTH N.  ON RETURN THE ARRANGEMENT IS:
C     C=A(0),A(1),B(1),...,A((N-2)/2),B((N-2)/2),A(N/2)   N COEFFICIENTS
C     P=P(W0),P(W1),P(W2),...,P(WN/2)                      NP ORDINATES
C
C     DPER is probably a better choice.
C
C     USAGE
C     CALL FASTF(X,N,C,P,NP)
C
C     SUBROUTINES CALLED
C     FFT, REALTR
C
C     ARGUMENTS
C     X - INPUT REAL TIME SERIES.
C     N - LENGTH OF X.  SEE RESTRICTIONS ON N BELOW.
C     C - VECTOR OF COMPUTED FOURIER COEFFICIENTS, DIMENSION C .GT. N+2
C     P - VECTOR OF COMPUTED PERIODOGRAM ORDINATES, DIMENSION P .GT. N+2
C     NP- COMPUTED NUMBER OF PERIODOGRAM ORDINATES
C
C     RESTRICTIONS ON N
C     FACTOR N INTO
C     N = N(1)*N(2)*N(3)*...*N(K)
C     WHERE EACH N(I) IS A PRIME NUMBER.
C     1. N MUST BE EVEN AND LARGER THAN 2.
C     2. THE LARGEST N(I) MUST BE LESS THAN OR EQUAL TO 23.
C     3. K MUST BE LESS THAN OR EQUAL TO 11.
C     4. THE PRODUCT OF THOSE PRIME NUMBERS WHICH OCCUR AN ODD NUMBER
C        OF TIMES IN THE FACTORIZATION OF N MUST BE LESS THAN OR EQUAL
C        TO 210.
C
C     EXAMPLE
C     10200 = 2*2*2*3*5*5*17
C     1. 10200 IS EVEN AND LARGER THAN 2.
C     2. 17 IS LESS THAN 23.
C     3. K=7 IS LESS THAN 11.
C     4. 2*3*17=102 IS LESS THAN 210.
C
C     REMARK
C     IF THE FOURIER COEFFICIENTS ARE NOT NEEDED THE USAGE
C     CALL FASTF(X,N,P,P,NP)
C     IS PERMISSIBLE.  IF X CAN BE DESTROYED THE USAGE
C     CALL FASTF(X,N,X,X,NP)
C     IS PERMISSIBLE.  X WILL CONTAIN THE PERIODOGRAM ORDINATES ON
C     RETURN.  DIMENSION X .GT. N+2 WITH THIS USAGE.
C
C     REFERENCE
C     R. C. SINGLETON. AN ALGORITHM FOR COMPUTING THE MIXED RADIX FAST
C     FOURIER TRANSFORM. IEEE TRANSACTIONS ON AUDIO AND ELECTROACOUSTICS
C     AU-17. JUNE 1969.
C
      SUBROUTINE FASTF(X,N,C,P,NP)
      implicit real*8 (a-h,o-z)
      save
      REAL*8 X(1),C(2),P(1)
      DO 10 I=1,N
10    C(I)=X(I)
      M=N/2
      NP=M+1
      ISN=2
      CALL FFT(C,C(2),M,M,M,ISN)
      CALL REALTR(C,C(2),M,ISN)
      A=1.d0/N
      B=N*.5d0
      C(1)=A*C(1)*.5d0
      P(1)=B*C(1)*C(1)*4.d0
      DO 20 I=2,M
      K=2*I-2
      C(K)=A*C(K+1)
      C(K+1)=A*C(K+2)
      P(I)=B*C(K)**2+B*C(K+1)**2
20    CONTINUE
      C(N)=A*C(N+1)*.5d0
      P(NP)=B*C(N)*C(N)*4.d0
      RETURN
      END
