      subroutine getshp(xlat, iday, estozn, fgprf, aptoz, apprf, 
     2       lprint, delshp) 
c
c***********************************************************************
c
c     getshp 
c
c     march, 2001  by charlie wellemeyer SSAI
c
c     purpose
c        determine a priori ozone profile shape a given ozone, day and
c        latitude
c
c     modification
c        none
c
c     variables
c       name      type  i/o   description
c       ----      ----  ---   ------------
c
c      arguments
c       xlat       r*4    i    latitude
c       iday       i*4    i    day of year
c       estozn     r*4    i    build 1 total ozone estimate
c       fgprf(11)  r*4    i    build 1 standard profile
c       apprf(11)  r*4    o    3-D a priori profile   
c       delshp(11) r*4    o    a priori ozone profile shape (zero sum)
c
c      internal
c       l1         i*4         latitude index for lower interpolation node
c       l2         i*4         latitude index for upper interpolation node
c       fracl      r*4         interpolation factor for latitude
c       cntrs(18)  r*4         center latitudes of climatological bands
c       m1         i*4         month index for previous interpolation node
c       m2         i*4         month index for subsequent interpolation node
c       fracm      r*4         interpolation factor for month
c
c***********************************************************************
c
      implicit none
      real t234a, t234b
      integer*4 lyr, i
c
c     -- input parameters
c
      real xlat, estozn, fgprf(11)
      integer iday 
      logical lprint(30), dmpit
c
c     -- internal parameters
c
      real aptoz, apprf(11), stnapr(11), newprf(11), aploc,
     1     staploc, fgloc, dfloc, apfnloc
c
c     -- output parameters
c
      real delshp(11)
c
c  call stnprf to interpolate standard profile to a priori total ozone
c
      call stnprf(xlat, aptoz, lprint, stnapr)
c
c  add profile shape difference at aptoz to standard profile at estozn
c  then subtract build 1 profile from resultant profile to get shape diff.
c
      dmpit = .false.
      t234a=0.0
      t234b=0.0
      aploc=0.0
      staploc=0.0
      fgloc=0.0
      dfloc=0.0
      apfnloc=0.0
      do 100 lyr=1,11
        newprf(lyr) = apprf(lyr) - stnapr(lyr) + fgprf(lyr)
        delshp(lyr) = newprf(lyr) - fgprf(lyr)
        apprf(lyr)=newprf(lyr)
        aploc=aploc+apprf(lyr)
        staploc=staploc+stnapr(lyr)
        fgloc=fgloc+fgprf(lyr)
        dfloc=dfloc+newprf(lyr) - fgprf(lyr)
        apfnloc=apfnloc+newprf(lyr)
        if(lyr.ge.2.and.lyr.le.4) then
          t234a=t234a+stnapr(lyr)
          t234b=t234b+newprf(lyr)
        endif
  100 continue
c
c  redistribute ozone in layers 2-4 to avoid negative values
c
      if(abs(xlat).lt.30.0.or.estozn.lt.175.0) then
        do 200 lyr=2,4
          newprf(lyr)=t234b*stnapr(lyr)/t234a
          delshp(lyr) = newprf(lyr) - fgprf(lyr)
          apprf(lyr)=newprf(lyr)
  200   continue
      endif
c
      if (lprint(21)) write (6,1000) xlat,iday,estozn,aptoz,
     1 apprf,aploc,stnapr,staploc,(apprf(i)-stnapr(i),i=1,11),
     2 dfloc,fgprf,fgloc,(fgprf(i)+apprf(i)-stnapr(i),i=1,11),
     3 apfnloc,delshp
c
      return
c
c  format statements
c
1000  format (/'Subroutine getshp'/'Input: xlat =',f7.3,
     1 '  iday  =',i7,'  estozn = ',f7.1,'  aptoz = ',f7.1,/,
     2  ' apprf, stnapr, apprf - stnapr, fgprf, fgprf + apprf',
     3  1x,'- stnapr, and delshp',/,'  0    1    2    3    4',
     4  '    5    6    7    8    9   10    total',6(/11f6.1,f8.1))
c
      end
