subroutine meanrho
c
c=======================================================================
c                                                                    ===
c  This routine computes the background density field RHOBAR at      ===
c  each grid point.                                                  ===
c                                                                    ===
c=======================================================================
c
#include <cdefs.h>
#include <param.h>
#include <pconst.h>
#include <onedim.h>
#include <options.h>
#include <hybrid.h>
#include <rhomean.h>
c
c-----------------------------------------------------------------------
c  Define local data.
c-----------------------------------------------------------------------
c
      integer i,j,k
      FLOAT
     *      der1,derkm,frsd,sal,tem,z
      FLOAT
     *      d2s(mprof),d2t(mprof),depth(imt,km),rho(imt,km),
     *      s(imt,km),t(imt,km),wk(mprof)
c
      parameter (der1=c1e30,derkm=c1e30)
c
c-----------------------------------------------------------------------
c  Begin executable code.
c-----------------------------------------------------------------------
c
c  Calculate background density field
c
      if (iflag(2).eq.1) then
        call spline(tinit(1,3),tinit(1,1),nprof,der1,derkm,d2t,wk)
        call spline(tinit(1,3),tinit(1,2),nprof,der1,derkm,d2s,wk)
      end if
c
      do 30 j=1,jmt
        call depthslab(j,tgrid,depth)
        do 10 i=1,imt
        do 10 k=1,km
          z=depth(i,k)
          if (iflag(2).eq.1) then
            call splint(tinit(1,3),tinit(1,1),d2t,nprof,z,tem,frsd)
            call splint(tinit(1,3),tinit(1,2),d2s,nprof,z,sal,frsd)
           else
            call lintrp(nprof,tinit(1,3),tinit(1,1),1,z,tem)
            call lintrp(nprof,tinit(1,3),tinit(1,2),1,z,sal)
          end if
          t(i,k)=tem
          s(i,k)=sal-smean
  10    continue
        call state(t,s,depth,rho)
        do 20 k=1,km
        do 20 i=1,imt
          rhobar(i,j,k)=rho(i,k)
  20    continue
  30  continue
c
      return
      end