SUBROUTINE thermcell_env(ngrid,nlay,po,pt,pu,pv,pplay,  &
     &           pplev,zo,zh,zl,ztv,zthl,zu,zv,zpspsk,pqsat,lev_out)

!--------------------------------------------------------------
!thermcell_env: calcule les caracteristiques de l environnement
!necessaires au calcul des proprietes dans le thermique
!--------------------------------------------------------------


   USE thermcell_ini_mod, ONLY : prt_level,RLvCp,RKAPPA,RETV
   IMPLICIT NONE

! arguments

   integer,intent(in) :: ngrid,nlay,lev_out
   real,intent(in), dimension(ngrid,nlay) :: po,pt,pu,pv,pplay
   real,intent(in), dimension(ngrid,nlay+1) :: pplev
   real,intent(out), dimension(ngrid,nlay) :: zo,zl,zh,ztv,zthl
   real,intent(out), dimension(ngrid,nlay) :: zpspsk,zu,zv,pqsat
   
! Local

   integer ig,ll
   real dqsat_dT
   logical mask(ngrid,nlay)


!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
! Initialisations :
!------------------

   mask(:,:)=.true.

!
! calcul des caracteristiques de l environnement
   DO  ll=1,nlay
     DO ig=1,ngrid
        zo(ig,ll)=po(ig,ll)
        zl(ig,ll)=0.
        zh(ig,ll)=pt(ig,ll)
     enddo
   enddo

! Condensation :
!---------------
! Calcul de l'humidite a saturation et de la condensation

   call thermcell_qsat(ngrid*nlay,mask,pplev,pt,po,pqsat)
   do ll=1,nlay
      do ig=1,ngrid
         zl(ig,ll) = max(0.,po(ig,ll)-pqsat(ig,ll))
         zh(ig,ll) = pt(ig,ll)+RLvCp*zl(ig,ll)         !   T = Tl + Lv/Cp ql
         zo(ig,ll) = po(ig,ll)-zl(ig,ll)
      enddo
   enddo

!-----------------------------------------------------------------------
   if (prt_level.ge.1) print*,'0 OK convect8'

   do ll=1,nlay
      do ig=1,ngrid
          zpspsk(ig,ll)=(pplay(ig,ll)/100000.)**RKAPPA
          zu(ig,ll)=pu(ig,ll)
          zv(ig,ll)=pv(ig,ll)
!attention zh est maintenant le profil de T et plus le profil de theta !
! Quelle horreur ! A eviter.
!   T-> Theta
            ztv(ig,ll)=zh(ig,ll)/zpspsk(ig,ll)
!Theta_v
            ztv(ig,ll)=ztv(ig,ll)*(1.+RETV*(zo(ig,ll))-zl(ig,ll))
!Thetal
            zthl(ig,ll)=pt(ig,ll)/zpspsk(ig,ll)
!            
      enddo
   enddo
 
 RETURN
   END