      SUBROUTINE GGPNT3(U,V,UG,VG,ILG,ILG1,ILAT,DLON,DLAT,
     1                  SLON,SLAT,USP,VSP,UNP,VNP,IFLD)
C
C     * MAR 25/13 - S.KHARIN
C
C     * CUBIC INTERPOLATION FROM A SCALAR FIELD UG(ILG1,ILAT) OR
C     * A VECTOR FIELD (UG(ILG1,ILAT), VG(ILG1,ILAT)) TO A SINGLE
C     * POINT U OR (U,V) WITH LATITUDE DLAT AND LONGITUDE DLON.
C     * DLAT IS MEASURED IN DEGREES FROM THE SOUTH POLE.
C     * DLON IS MEASURED IN DEGREES EASTWARD FROM THE GREENWICH.
C     * SLAT AND SLON ARE INPUT GRID LATITUDES AND LONGITUDES.
C     * USP,VSP,UNP,VNP ARE SCALAR OR VECTOR ESTIMATES AT THE POLES.
C     * IFLD=0 USE SCALAR FIELD UG.
C     * IFLD=1 USE A VECTOR COMPONENT UG.
C     * IFLD=2 USE BOTH VECTOR COMPONENTS (UG,VG).
C     * RETURN SPVAL IF ANY OF NEIGHBORS HAS SPVAL.
C
      IMPLICIT NONE
      REAL U,V,DLAT,DLON,USP,VSP,UNP,VNP
      INTEGER ILG,ILG1,ILAT,IFLD,I,J
      REAL UG(ILG1,ILAT),VG(ILG1,ILAT),SLON(ILG1),SLAT(ILAT)
      REAL G(4,4),X,Y
      INTEGER IX(4),IY(4)
      REAL SPVAL,SPVALT
      DATA SPVAL/1.E+38/,SPVALT/1.E+32/
C-----------------------------------------------------------------------
      IF(IFLD.LT.0.OR.IFLD.GT.2) CALL              XIT('GGPNT3',-1)
      IF(ILG.LT.4.AND.ILG.NE.1.OR.ILAT.LT.4) CALL  XIT('GGPNT3',-2)
C
C     * FIND 4 CLOSEST INPUT LONGITUDES (TAKE INTO ACCOUNT PERIODICITY)
C
      IF(ILG.EQ.1)THEN
C
C       * ZONAL AVERAGE CASE
C
        IX(1)=1
        IX(2)=1
        IX(3)=1
        IX(4)=1
        X=0.E0
      ELSE IF(DLON.LT.SLON(1))THEN
        IX(1)=ILG-1
        IX(2)=ILG
        IX(3)=1
        IX(4)=2
        X=(DLON+360.E0-SLON(IX(2)))/(SLON(IX(3))+360.E0-SLON(IX(2)))
      ELSE IF(DLON.LT.SLON(2))THEN
        IX(1)=ILG
        IX(2)=1
        IX(3)=2
        IX(4)=3
        X=(DLON-SLON(IX(2)))/(SLON(IX(3))-SLON(IX(2)))
      ELSE IF(DLON.GE.SLON(ILG))THEN
        IX(1)=ILG-1
        IX(2)=ILG
        IX(3)=1
        IX(4)=2
        X=(DLON-SLON(IX(2)))/(SLON(IX(3))+360.E0-SLON(IX(2)))
      ELSE IF(DLON.GE.SLON(ILG-1))THEN
        IX(1)=ILG-2
        IX(2)=ILG-1
        IX(3)=ILG
        IX(4)=1
        X=(DLON-SLON(IX(2)))/(SLON(IX(3))-SLON(IX(2)))
      ELSE
        DO I=1,ILG1
          IF(DLON.LT.SLON(I)) GO TO 110
        ENDDO
 110    CONTINUE
        IX(1)=I-2
        IX(2)=I-1
        IX(3)=I
        IX(4)=I+1
        X=(DLON-SLON(IX(2)))/(SLON(IX(3))-SLON(IX(2)))
      ENDIF
C
C     * FIND 4 CLOSEST INPUT LATITUDES (RESTRICT TO 0...ILAT+1)
C
      IF(DLAT.LT.SLAT(2)) THEN
        IY(1)=0
        IY(2)=1
        IY(3)=2
        IY(4)=3
      ELSE IF(DLAT.GE.SLAT(ILAT-1)) THEN
        IY(1)=ILAT-2
        IY(2)=ILAT-1
        IY(3)=ILAT
        IY(4)=ILAT+1
      ELSE
        DO J=1,ILAT
          IF(DLAT.LT.SLAT(J)) GO TO 120
        ENDDO
 120    CONTINUE
        IY(1)=J-2
        IY(2)=J-1
        IY(3)=J
        IY(4)=J+1
      ENDIF
      Y=(DLAT-SLAT(IY(2)))/(SLAT(IY(3))-SLAT(IY(2)))
C
C     * SELECT 4X4 ARRAY FROM UG FOR CUBIC INTERPOLATION
C     * (USE SOUTH OR NORTH POLE APPROXIMATIONS FOR J=0, ILAT+1)
C
      DO J=1,4
        DO I=1,4
          G(I,J)=SPVAL
          IF(IY(J).EQ.0)THEN
C
C           * SOUTH POLE
C
            IF(IFLD.EQ.0)THEN
              G(I,J)=USP
            ELSE
              IF(ABS(USP-SPVAL).GT.SPVALT)
     1             G(I,J)=USP*COSD(DLON)+VSP*SIND(DLON)
            ENDIF
          ELSE IF(IY(J).EQ.ILAT+1)THEN
C
C           * NORTH POLE
C
            IF(IFLD.EQ.0)THEN
              G(I,J)=UNP
            ELSE
              IF(ABS(UNP-SPVAL).GT.SPVALT)
     1             G(I,J)=UNP*COSD(DLON)+VNP*SIND(DLON)
            ENDIF
          ELSE
C
C           * AWAY FROM THE POLES
C
            G(I,J)=UG(IX(I),IY(J))
          ENDIF
          IF(ABS(G(I,J)-SPVAL).LT.SPVALT)THEN
            U=SPVAL
            RETURN
          ENDIF
        ENDDO
      ENDDO
C
C     * CUBIC INTERPOLATION OF U
C
      CALL GINTC2(U,G,4,4,2.E0+X,2.E0+Y)
      IF(IFLD.LE.1)RETURN
C
C     * SELECT 4X4 ARRAY FROM VG FOR CUBIC INTERPOLATION
C     * (USE SOUTH OR NORTH POLE APPROXIMATIONS FOR J=0, ILAT+1)
C
      DO J=1,4
        DO I=1,4
          IF(IY(J).EQ.0)THEN
C
C           * SOUTH POLE
C
            IF(ABS(USP-SPVAL).GT.SPVALT)
     1           G(I,J)=USP*SIND(DLON)-VSP*COSD(DLON)
          ELSE IF(IY(J).EQ.ILAT+1)THEN
C
C           * NORTH POLE
C
            IF(ABS(UNP-SPVAL).GT.SPVALT)
     1           G(I,J)=VNP*COSD(DLON)-UNP*SIND(DLON)
          ELSE
C
C           * AWAY FROM THE POLES
C
            G(I,J)=VG(IX(I),IY(J))
          ENDIF
          IF(ABS(G(I,J)-SPVAL).LT.SPVALT)THEN
            V=SPVAL
            RETURN
          ENDIF
        ENDDO
      ENDDO
C
C     * CUBIC INTERPOLATION OF V
C
      CALL GINTC2(V,G,4,4,2.E0+X,2.E0+Y)
      RETURN
      END
