 
      SUBROUTINE LOLIP (NXG,XG,NYG,YG,NPI,XI,YI,NP,MP,MPM,NMAX,D)
C
C     THIS SUBROUTINE DETERMINES THE LOCAL INTERPOLATION POINTS FOR THE
C     GRID VERSION OF FRANKE'S METHOD OF SURFACE INTERPOLATION.
C     MINPTS POINTS ARE REQUIRED FOR EACH REGION.
C     IF FEWER THAN MINPTS POINTS ARE FOUND IN THE REGION, THE NEXT
C     CLOSEST POINTS (IN THE SUP NORM AFTER THE CURRENT RECTANGLE IS
C     TRANSFORMED ONTO (0,1)) ARE USED.  MINPTS IS SET TO 3, WHICH IS
C     THE RECOMMENDED VALUE, ALTHOUGH IT MAY BE ALTERED.
C
C     THE ARGUMENTS ARE AS FOLLOWS.
C
C        NXG  - INPUT.  NUMBER OF VERTICAL GRID LINES.
C        XG   - INPUT.  THE COORDINATES OF THE VERTICAL GRID LINES, IN
C                       INCREASING ORDER
C        NYG  - INPUT.  NUMBER OF HORIZONTAL GRID LINES.
C        YG   - INPUT.  THE COORDINATES OF THE HORIZONTAL GRID LINES,
C                       IN INCREASING ORDER.
C        NPI  - INPUT.  THE NUMBER OF DATA POINTS.
C        XI   - \
C        YI   - INPUT.  THE DATA POINTS (XI,YI), I=1,...,NPI.
C        FI   - /
C        NP   - OUTPUT.  AN ARRAY WHICH GIVES THE INITIAL SUBSCRIPT IN
C                       THE ARRAY MP AT WHICH THE SUBSCRIPTS FOR THE
C                       LOCAL INTERPOLATION POINTS ARE STORED.
C        MP   - OUTPUT.  AN ARRAY WHICH GIVES THE SUBSCRIPTS FOR THE
C                       LOCAL INTERPOLATION POINTS.
C        MPM  - INPUT.  DIMENSION OF THE ARRAY MP IN THE CALLING PROGRAM
C        NMAX - OUTPUT.  THE MAXIMUM NUMBER OF INTERPOLATION POINTS
C                        OVER ALL THE REGIONS.
C        D    - A WORK ARRAY OF DIMENSION AT LEAST NPI.
C
