SUBROUTINE SMOOFT(Y,N,PTS)
      SAVE
      PARAMETER(MMAX=4096)
      DIMENSION Y(MMAX), X(MMAX)
      M=2
      NMIN=N+2.*PTS
1     IF (M.LT.NMIN) THEN
        M=2*M
          GO TO 1
      END IF
      IF(M.GT.MMAX) PAUSE 'MMAX too small'
      CONST=(PTS/M)**2
C
C  Remove a slope from the data ( to be restored after the smoothing)
C
      DO 5 I = 1, N
          X(I) = FLOAT(I)
    5 CONTINUE
      CALL FITLINE(N, X, Y, SLOPE, ZPOINT )
      Y1=Y(1)
      YN=Y(N)
      RN1=1./(N-1.)
      DO 11 J=1,N
COMMENT  Y(J)=Y(J)-RN1*(Y1*(N-J)+YN*(J-1))
         X(J) = Y(J) - ZPOINT - SLOPE * FLOAT(J)
   11 CONTINUE
      IF (N+1 .LE. M) THEN
        DO 12 J=N+1,M
          X(J)=0.
   12   CONTINUE
      END IF
      MO2=M/2
      CALL REALFT(X, MO2, 1 )
      X(1)=X(1)/MO2
      FAC=1.
      DO 13 J=1,MO2-1
        K=2*J+1
        IF(FAC.NE.0.)THEN
          FAC=AMAX1(0.,(1.-CONST*J**2)/MO2)
          X(K)=FAC*X(K)
          X(K+1)=FAC*X(K+1)
        ELSE
          X(K)=0.
          X(K+1)=0.
        END IF
13    CONTINUE
      FAC=AMAX1(0.,(1.-0.25*PTS**2)/MO2)
      X(2)=FAC*X(2)
      CALL REALFT(X,MO2,-1)
C
C   Restore the linear trend in the data.
C
      DO 14 J=1,N
COMMENT  Y(J)=RN1*(Y1*(N-J)+YN*(J-1))+Y(J)
         Y(J) = X(J) + ZPOINT + SLOPE * FLOAT(J)
14    CONTINUE
      RETURN
      END
      SUBROUTINE REALFT(DATA,N,ISIGN)
      SAVE
      REAL*8 WR,WI,WPR,WPI,WTEMP,THETA
      DIMENSION DATA(*)
      THETA=6.28318530717959D0/2.0D0/DBLE(N)
      C1=0.5
      IF (ISIGN.EQ.1) THEN
        C2=-0.5
        CALL FOUR1(DATA,N,+1)
      ELSE
        C2=0.5
        THETA=-THETA
      ENDIF
      WPR=-2.0D0*DSIN(0.5D0*THETA)**2
      WPI=DSIN(THETA)
      WR=1.0D0+WPR
      WI=WPI
      N2P3=2*N+3
      DO 11 I=2,N/2+1
        I1=2*I-1
        I2=I1+1
        I3=N2P3-I2
        I4=I3+1
        WRS=SNGL(WR)
        WIS=SNGL(WI)
        H1R=C1*(DATA(I1)+DATA(I3))
        H1I=C1*(DATA(I2)-DATA(I4))
        H2R=-C2*(DATA(I2)+DATA(I4))
        H2I=C2*(DATA(I1)-DATA(I3))
        DATA(I1)=H1R+WRS*H2R-WIS*H2I
        DATA(I2)=H1I+WRS*H2I+WIS*H2R
        DATA(I3)=H1R-WRS*H2R+WIS*H2I
        DATA(I4)=-H1I+WRS*H2I+WIS*H2R
        WTEMP=WR
        WR=WR*WPR-WI*WPI+WR
        WI=WI*WPR+WTEMP*WPI+WI
11    CONTINUE
      IF (ISIGN.EQ.1) THEN
        H1R=DATA(1)
        DATA(1)=H1R+DATA(2)
        DATA(2)=H1R-DATA(2)
      ELSE
        H1R=DATA(1)
        DATA(1)=C1*(H1R+DATA(2))
        DATA(2)=C1*(H1R-DATA(2))
        CALL FOUR1(DATA,N,-1)
      ENDIF
      RETURN
      END
      SUBROUTINE FOUR1(DATA,NN,ISIGN)
      SAVE
      REAL*8 WR,WI,WPR,WPI,WTEMP,THETA
      DIMENSION DATA(*)
      N=2*NN
      J=1
      DO 11 I=1,N,2
        IF(J.GT.I)THEN
          TEMPR=DATA(J)
          TEMPI=DATA(J+1)
          DATA(J)=DATA(I)
          DATA(J+1)=DATA(I+1)
          DATA(I)=TEMPR
          DATA(I+1)=TEMPI
        ENDIF
        M=N/2
1       IF ((M.GE.2).AND.(J.GT.M)) THEN
          J=J-M
          M=M/2
        GO TO 1
        ENDIF
        J=J+M
11    CONTINUE
      MMAX=2
2     IF (N.GT.MMAX) THEN
        ISTEP=2*MMAX
        THETA=6.28318530717959D0/(ISIGN*MMAX)
        WPR=-2.D0*DSIN(0.5D0*THETA)**2
        WPI=DSIN(THETA)
        WR=1.D0
        WI=0.D0
        DO 13 M=1,MMAX,2
          DO 12 I=M,N,ISTEP
            J=I+MMAX
            TEMPR=SNGL(WR)*DATA(J)-SNGL(WI)*DATA(J+1)
            TEMPI=SNGL(WR)*DATA(J+1)+SNGL(WI)*DATA(J)
            DATA(J)=DATA(I)-TEMPR
            DATA(J+1)=DATA(I+1)-TEMPI
            DATA(I)=DATA(I)+TEMPR
            DATA(I+1)=DATA(I+1)+TEMPI
12        CONTINUE
          WTEMP=WR
          WR=WR*WPR-WI*WPI+WR
          WI=WI*WPR+WTEMP*WPI+WI
13      CONTINUE
        MMAX=ISTEP
      GO TO 2
      ENDIF
      RETURN
      END