      subroutine extrap2 (val,xtrp)
c
c=======================================================================
c                                                                    ===
c  This subroutine provides extrapolated values for those values so  ===
c  indicated.                                                        ===
c                                                                    ===
c  ------                                                            ===
c  Input:                                                            ===
c  ------                                                            ===
c                                                                    ===
c  VAL    Accepted & questionable values.        (real array)        ===
c  XTRP   flags to indicate extrapolation.       (logical array)     ===
c            [.true.]   an extrapolated value.                       ===
c            [.false.]  an accepted value.                           ===
c                                                                    ===
c  -------                                                           ===
c  Output:                                                           ===
c  -------                                                           ===
c                                                                    ===
c  VAL    Accepted & corrected values.           (real array)        ===
c                                                                    ===
c  For the second order pressure gradient, this routine chooses a    ===
c  horizontal interpolation scheme to  minimize pressure gradient    ===
c  errors.   This scheme varies depending upon the number of good    ===
c  (nonextrapolated) points.                                         ===
c                                                                    ===
c  * For 0 good points:  The first value is arbitrarily transmitted  ===
c                        to the other three points (No gradient).    ===
c                                                                    ===
c  * For 1 good point:  That value is transmitted to the other       ===
c                       three points (No gradient).                  ===
c                                                                    ===
c  * For 4 good points:  nothing is done.                            ===
c                                                                    ===
c  * For 2 or 3 good points:  The orientation needs to be            ===
c                             determined.                            ===
c                                                                    ===
c  The vectors of data (val) and extrapolation flags (xtrp) come in  ===
c  the following orientation with the following index properties:    ===
c                                                                    ===
c     3--4                     sum of all indices = 10               ===
c     |  |                                                           ===
c     1--2         sum of opposite corner indices = 5                ===
c                                                                    ===
c  The orientations can be decided using the sum of the indices of   ===
c  the good points INDSUM. Denoting good points with a "g" and bad   ===
c  points with a "b", the horizontal interpolations go as follows:   ===
c                                                                    ===
c     g--b     >>  2 good points, opposite corners:                  ===
c     |  |           Bad points get the mean of the good points      ===
c     b--g           (evenly distribute diagonal gradient in ordinal ===
c                    directions). This is only configuration with    ===
c                    INDSUM=5.                                       ===
c                                                                    ===
c     g----g   >>  2 good points, adjacent:                          ===
c     |    |         Transmit good values according to figure        ===
c    \|/  \|/        (preserve known gradient in given direction).   ===
c     b----b                                                         ===
c                If N is the index for one bad point then:           ===
c                                                                    ===
c                         index of other bad point = 10-INDSUM-N     ===
c                index of good point adjacent to N = 5-(10-INDSUM-N) ===
c                                                  = INDSUM+N-5      ===
c                                                                    ===
c     +1       >>  3 good points :                                   ===
c      g--b          Use second order extrapolation (weights at      ===
c      |  |          right). Best estimate.                          ===
c      g--g                                                          ===
c     -1  +1     If N is the index for the bad point then:           ===
c                                                                    ===
c                                       N = 10-INDSUM                ===
c                index of opposite corner = 5-N                      ===
c                                         = INDSUM-5                 ===
c                                                                    ===
c  Note:  The cases of 2 good points, opposite corners and 3 good    ===
c         points are written as a value plus a "small correction".   ===
c         In each case, the value and correction pair are chosen to  ===
c         give the smallest relative correction.                     ===
c                                                                    ===
c  Calls:  none                                                      ===
c                                                                    ===
c=======================================================================
c
c-----------------------------------------------------------------------
c  Define global data.
c-----------------------------------------------------------------------
c
#include <cdefs.h>
#include <pconst.h>
c
c-----------------------------------------------------------------------
c  Define local data.
c-----------------------------------------------------------------------
c
      integer nval
      parameter (nval=4)
c
      logical xtrp(nval)
      integer indsum,n,nbad,ngood1,no_acc
      FLOAT
     *      avg_val,diff1,diff2,val1,val2
      FLOAT
     *      val(nval)
c
c=======================================================================
c  Begin executable code.
c=======================================================================
c
c-----------------------------------------------------------------------
c  Interpolate values depending on the number of acceptible values.
c-----------------------------------------------------------------------
c
      indsum=0
      no_acc=0
      avg_val=c0
      do 10 n=1,nval
        if(.not.xtrp(n)) then
          indsum=indsum+n
          no_acc=no_acc+1
          avg_val=avg_val+val(n)
        endif
  10  continue
c
      if((indsum.eq.5).or.(no_acc.eq.3)) then
c
        if(indsum.eq.5) then
          n=1
          do 20 while (xtrp(n))
            n=n+1
  20      continue
          val1=val(n)
          diff1=p5*(val(5-n)-val(n))
          val2=val(5-n)
          diff2=p5*(val(n)-val(5-n))
        else
          nbad=10-indsum
          if (mod(nbad,2).eq.1) then
            ngood1=nbad+1
          else
            ngood1=nbad-1
          endif
          val1=val(ngood1)
          diff1=val(5-ngood1)-val(indsum-5)
          val2=val(5-ngood1)
          diff2=val(ngood1)-val(indsum-5)
        endif
c
        if(abs(val2*diff1).le.abs(val1*diff2)) then
          avg_val=val1+diff1
        else
          avg_val=val2+diff2
        endif
c
      endif
c
      if((no_acc.eq.1).or.(no_acc.eq.3).or.(indsum.eq.5)) then
        do 30 n=1,nval
          if(xtrp(n)) val(n)=avg_val
  30    continue
      elseif(no_acc.eq.2) then
        do 40 n=1,nval
          if(xtrp(n)) val(n)=val(indsum+n-5)
  40    continue
      elseif(no_acc.eq.0) then
        do 50 n=2,nval
          val(n)=val(1)
  50    continue
      endif
c
      return
      end
