      SUBROUTINE POLINT(XA,YA,X,Y)

C ADAPTED FROM NUMERICAL RECIPES, PRESS et al., p. 82, 
C WITH N = 3 AND YA(2) AS THE INITIAL ESTIMATE FOR Y.
C XA ARE KNOWN

C      INCLUDE  'PARAMETR.INC'

      REAL*8    XA(3), YA(3), C(3), D(3)
      REAL*8    HO, HP, W, DEN, X, Y, DY
      INTEGER*4 NS, I, M

      NS           = 1

      DO I = 1, 3
          C(I)     = YA(I)
          D(I)     = YA(I)
      END DO    

      Y            = YA(2)

      DO M = 1, 2

          DO I = 1, 3-M
              HO   = XA(I)-X
              HP   = XA(I+M)-X
              W    = C(I+1)-D(I)
              DEN  = HO-HP
              DEN  = W/DEN
              D(I) = HP*DEN
              C(I) = HO*DEN
          END DO

          IF (2*NS.LT.3-M) THEN
              DY   = C(NS+1)
          ELSE
              DY   = D(NS)
              NS   = NS-1
          END IF

          Y        = Y+DY

      END DO

      RETURN
      END


