cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c        this file contains the subroutines which correspond       c
c        to those in slgcm.f in the 3d model, except the routine   c
c        called slgcst.                                            c
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
      subroutine fipole(win,fi,vin,angin,spole)
      include "slgcm.slgprm"
      dimension  win(imp3,jmp3), fi(imp3,jm)
c
      logical spole
c
      gl1st = angin
      if (spole) then
         ph1st = -pio2 + vin
         if (ph1st .le. -pio2 + dph*0.05) then
             fip = win(1,3)
            go to 101
         endif
      else
         ph1st =  pio2 - vin
         if (ph1st .ge. pio2 - dph*0.05) then
            fip = win(1,jmp2)
            go to 101
         endif
      endif
c
      xgl1st = (gl1st+dlm)*rdlm
      xph1st = (ph1st+pio2)*rdph + 3
c
c     if (ipr .le. 1 .or. ipr .ge. imp2 .or. jpr .le. 1 .or.
c    *    jpr .ge. jmp2) print *,' in fipole ipr=',ipr,' jpr=', jpr
c     print *,' in fipole ipr=',ipr,' jpr=', jpr
c
      if (.not. tricub) then
         ipr    = xgl1st + 0.5
         jpr    = xph1st + 0.5
         btpr   = xph1st-jpr
         alpr   = xgl1st-ipr
         if (ipr .eq. 1) ipr = imp1
         iprp1  = ipr + 1
         jprp1  = jpr + 1
         iprm1  = ipr - 1
         jprm1  = jpr - 1
         tem1   = 1.0 - alpr
         tem2   = 1.0 + alpr
         c1a =  - alpr * tem1 * 0.5
         c2a =    tem1 * tem2
         c3a =    alpr * tem2 * 0.5
c
         tem1   = 1.0 - btpr
         tem2   = 1.0 + btpr
         c1b =  - btpr * tem1 * 0.5
         c2b =    tem1 * tem2
         c3b =    btpr * tem2 * 0.5
c
         fip =  c1a*(c1b*win(iprm1,jprm1 ) + c2b*win(iprm1,jpr )
     *             + c3b*win(iprm1,jprp1 ))
     *       +  c2a*(c1b*win(ipr  ,jprm1 ) + c2b*win(ipr  ,jpr )
     *             + c3b*win(ipr  ,jprp1 ))
     *       +  c3a*(c1b*win(iprp1,jprm1 ) + c2b*win(iprp1,jpr )
     *             + c3b*win(iprp1,jprp1 ))
      else
         ipr    = xgl1st
         jpr    = xph1st
         btpr   = xph1st-jpr
         alpr   = xgl1st-ipr
         if (ipr .eq. 1) ipr = imp1
         iprp1  = ipr + 1
         jprp1  = jpr + 1
         iprm1  = ipr - 1
         jprm1  = jpr - 1
         alprsq = alpr*alpr
         btprsq = btpr*btpr
         iprp2  = ipr + 2
         jprp2  = jpr + 2
         twomal = (2.-alpr) * 0.5
         twombt = (2.-btpr) * 0.5
         onmasq = 1. - alprsq
         onmbsq = 1. - btprsq
c
c
         c1a = onmasq*twomal
         c2a = (alprsq+alpr)*twomal
         c3a = third*(alprsq-alpr)*twomal
         c4a = sixth*onmasq*(-alpr)
c
         c1b = onmbsq*twombt
         c2b = (btprsq+btpr)*twombt
         c3b = third*(btprsq-btpr)*twombt
         c4b = sixth*onmbsq*(-btpr)
c
         fip  = c1a*(c1b*win(ipr  ,jpr   ) + c2b*win(ipr  ,jprp1 )
     *             + c3b*win(ipr  ,jprm1 ) + c4b*win(ipr  ,jprp2 ))
     *        + c2a*(c1b*win(iprp1,jpr   ) + c2b*win(iprp1,jprp1 )
     *             + c3b*win(iprp1,jprm1 ) + c4b*win(iprp1,jprp2 ))
     *        + c3a*(c1b*win(iprm1,jpr   ) + c2b*win(iprm1,jprp1 )
     *             + c3b*win(iprm1,jprm1 ) + c4b*win(iprm1,jprp2 ))
     *        + c4a*(c1b*win(iprp2,jpr   ) + c2b*win(iprp2,jprp1 )
     *             + c3b*win(iprp2,jprm1 ) + c4b*win(iprp2,jprp2 ))
c
c     else
c        fip = alpr*(btpr*win(iprp1,jprp1) + (1.0-btpr)*win(iprp1,jpr))
c    *   + (1-alpr)*(btpr*win(ipr  ,jprp1) + (1.0-btpr)*win(ipr,  jpr))
      endif
  101 continue
c
      if (spole) then
         do 310 i=1,imp3
         fi(i,1) = fip
  310    continue
      else
         do 320 i=1,imp3
         fi(i,jm) = fip
  320    continue
      endif
c
      return
      end
        subroutine vpole (imx, dlm, twopi, twooim, pifac
     *,                  vin, clam, slam, vout, angout)
c
        dimension  vin(imx),clam(imx),slam(imx)
        dimension  wka(1000),wkb(1000)
c
        sumvcs = 0.0
        sumvsn = 0.0
        do 10 i=1,imx
        wka(i) = vin(i) * clam(i)
        wkb(i) = vin(i) * slam(i)
   10   continue
        do 20 i=1,imx
        sumvcs    = sumvcs + wka(i)
        sumvsn    = sumvsn + wkb(i)
   20   continue
        sumvcs    = sumvcs * twooim
        sumvsn    = sumvsn * twooim
c
c
        if (sumvcs .ne. 0.0) then
           del   = atan2(sumvsn,sumvcs) + pifac
       elseif (sumvsn .gt. 0.0) then
           del    = twopi * 0.25 + pifac
       else
           del    = - twopi * 0.25 + pifac
       endif
       vout  = sqrt(sumvcs*sumvcs+sumvsn*sumvsn)
c
c  redefine lam* such that dlm < lam* < (twopi+dlm).
c
        ratio  = (del-dlm)/twopi
        iratio = ratio
        if (ratio .lt. 0) iratio = iratio-1
        result = (ratio-iratio)*twopi
        angout = result + dlm
  100   continue
c
        return
        end
      subroutine filpu(ufill)
c
c*********************************************************************
c  fill in the 'ghost' values at u - points j=1 and j=jm
c*********************************************************************
c
      include "slgcm.slgprm"
      include "slgcm.slgclm"
c
      dimension ufill(imp3,jm)

      annp = 0.0
      bnnp = 0.0
      ansp = 0.0
      bnsp = 0.0
      do 20 i=1,im
      annp = annp + ufill(i,jmm1) * clm(i,1)
      bnnp = bnnp + ufill(i,jmm1) * slm(i,1)
      ansp = ansp + ufill(i,2)     * clm(i,1)
      bnsp = bnsp + ufill(i,2)     * slm(i,1)
 20   continue
      annp = annp * twooim
      bnnp = bnnp * twooim
      ansp = ansp * twooim
      bnsp = bnsp * twooim
      do 30 i=1,im
      ufill(i,1 )        = ansp * clm(i,1) + bnsp * slm(i,1)
      ufill(i,jm)        = annp * clm(i,1) + bnnp * slm(i,1)
  30  continue
c
      ufill(imp1,1) = ufill(1,1)
      ufill(imp2,1) = ufill(2,1)
      ufill(imp3,1) = ufill(3,1)
      ufill(imp1,jm) = ufill(1,jm)
      ufill(imp2,jm) = ufill(2,jm)
      ufill(imp3,jm) = ufill(3,jm)
c

      return
      end
        subroutine filpv(vfill)
c
c*********************************************************************
c  fill in the 'ghost' values at v - points j=jm
c*********************************************************************
c
      include "slgcm.slgprm"
c
        dimension vfill(imp3,jm)
c
        do 10 i=1,imb2
        vfill(i,jm)      = - vfill(i+imb2,jmm1)
        vfill(i+imb2,jm) = - vfill(i,jmm1)
   10   continue
c
        vfill(imp1,jm) = vfill(1,jm)
        vfill(imp2,jm) = vfill(2,jm)
        vfill(imp3,jm) = vfill(3,jm)
c
        return
        end
        subroutine fillv(vfill)
c
c*********************************************************************
c  fill in the 'ghost' values at v - points j=1,2 and j=jmp2,jmp3
c*********************************************************************
c
      include "slgcm.slgprm"
c
        dimension vfill(imp3,jmp3)
c
c
        do 10 i=1,imb2
        vfill(i,1   ) = - vfill(i+imb2,4   )
        vfill(i,2   ) = - vfill(i+imb2,3   )
        vfill(i,jmp2) = - vfill(i+imb2,jmp1)
        vfill(i,jmp3) = - vfill(i+imb2,jm  )
c
        vfill(i+imb2,1   ) = - vfill(i,4   )
        vfill(i+imb2,2   ) = - vfill(i,3   )
        vfill(i+imb2,jmp2) = - vfill(i,jmp1)
        vfill(i+imb2,jmp3) = - vfill(i,jm  )
   10   continue
c
        vfill(imp1,1   ) = vfill(1,1   )
        vfill(imp2,1   ) = vfill(2,1   )
        vfill(imp3,1   ) = vfill(3,1   )
        vfill(imp1,2   ) = vfill(1,2   )
        vfill(imp2,2   ) = vfill(2,2   )
        vfill(imp3,2   ) = vfill(3,2   )
        vfill(imp1,jmp2) = vfill(1,jmp2)
        vfill(imp2,jmp2) = vfill(2,jmp2)
        vfill(imp3,jmp2) = vfill(3,jmp2)
        vfill(imp1,jmp3) = vfill(1,jmp3)
        vfill(imp2,jmp3) = vfill(2,jmp3)
        vfill(imp3,jmp3) = vfill(3,jmp3)
c
        return
        end
      subroutine fillu(ufill)
c
c*********************************************************************
c  fill in the 'ghost' values at u - points j=1,2,3 and j=jmp2,jmp3
c*********************************************************************
c
      include "slgcm.slgprm"
c
      dimension ufill(imp3,jmp3)
c
      do 10 i=1,imb2
      ufill(i,   1) = -ufill(i+imb2,5)
      ufill(i,   2) = -ufill(i+imb2,4)
      ufill(i,jmp3) = -ufill(i+imb2,jmp1)
c
      ufill(i+imb2,   1) = -ufill(i,5)
      ufill(i+imb2,   2) = -ufill(i,4)
      ufill(i+imb2,jmp3) = -ufill(i,jmp1)
   10 continue
c
      ufill(imp1,1   ) = ufill(1,1   )
      ufill(imp2,1   ) = ufill(2,1   )
      ufill(imp3,1   ) = ufill(3,1   )
      ufill(imp1,2   ) = ufill(1,2   )
      ufill(imp2,2   ) = ufill(2,2   )
      ufill(imp3,2   ) = ufill(3,2   )
c
      ufill(imp1,jmp3) = ufill(1,jmp3)
      ufill(imp2,jmp3) = ufill(2,jmp3)
      ufill(imp3,jmp3) = ufill(3,jmp3)
c
      return
      end
      subroutine fillp(fi)
c
c*********************************************************************
c  fill in the 'ghost' values at phi - points j=1,2 and j=jmp3
c*********************************************************************
c
      include "slgcm.slgprm"
c
      dimension fi(imp3,jmp3)
c
      do 10 i=1,imb2
      fi(i,1     ) = fi(i+imb2,5   )
      fi(i,2   ) = fi(i+imb2,4   )
      fi(i,jmp3) = fi(i+imb2,jmp1)
c
      fi(i+imb2,1)        = fi(i,5)
      fi(i+imb2,2)     = fi(i,4)
      fi(i+imb2,jmp3) = fi(i,jmp1)
   10 continue
c
      fi(imp1,1)     = fi(1,1   )
      fi(imp2,1)     = fi(2 ,1  )
      fi(imp3,1)     = fi(3  ,1 )
      fi(imp1,2   )  = fi(1,2  )
      fi(imp2,2   )  = fi(2,2  )
      fi(imp3,2   )  = fi(3,2  )
      fi(imp1,jmp3) = fi(1,jmp3)
      fi(imp2,jmp3) = fi(2,jmp3)
      fi(imp3,jmp3) = fi(3,jmp3)
c
      return
      end
