      PROGRAM FFT

C     Do an FFT of a real time series
C     (Should be modfied to do a real FFT ...)

      PARAMETER ( MAXNT = 200000 )

      REAL    T( MAXNT )
      COMPLEX S( MAXNT )

C     *** Read in data ***

      DO IT = 1, MAXNT
         READ( *, *, END = 1000 ) T( IT ), STEMP
         S( IT ) = STEMP   ! real to complex
         NT = IT

         IF ( NT .GT. MAXNT ) STOP 'Too many points'
      END DO

 1000 TMAX = T( NT )

C     *** compute spectrum ***

C     *** Check NT is a power of 2 ***

      IF ( NT .LE. 0 )
     &   STOP 'FATAL ERROR: NT must be positive'

      N1 = 2**( INT( LOG10( REAL( NT ) ) / 0.30104 ) +1 )

      IF ( N1 .NE. NT ) THEN
         WRITE( *, * ) 'NT = ', NT
         STOP 'FATAL ERROR: NT must be a power of 2'
      END IF

      CALL CFFT( S, NT, 1 )
      CALL CFFTSC( S, NT )

C     *** write out the spectrum (w/out dc component) ***

      DO IT = 2, NT / 2 + 1
         DELTAF = 1.0 / TMAX
         F = ( IT - 1 ) * DELTAF
         WRITE( *, * ) F, S( IT )
      END DO

      STOP
      END
