snow_cover_1layer.F90

Go to the documentation of this file.

!SFX_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
!SFX_LIC This is part of the SURFEX software governed by the CeCILL-C licence
!SFX_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt  
!SFX_LIC for details. version 1.
!     #########
    SUBROUTINE snow_cover_1layer(PTSTEP, PANSMIN, PANSMAX, PTODRY,         &
                                   prhosmin, prhosmax, prhofold, oall_melt,  &
                                   pdrain_time, pwcrn, pz0sn, pz0hsn,        &
                                   ptsnow, pasnow, prsnow, pwsnow, pts_snow, &
                                   pesnow,                                   &
                                   ptg, ptg_coefa, ptg_coefb,                &
                                   pabs_sw, plw1, plw2,                      &
                                   pta, pqa, pvmod, pps, prhoa, psr,         &
                                   pzref, puref,                             &
                                   prnsnow, phsnow, plesnow, pgsnow, pmelt,  &
                                   pdqs_snow, pabs_lw                        )  
!   ##########################################################################
!
!!****  *SNOW_COVER_1LAYER*  
!!
!!    PURPOSE
!!    -------
!
!     One layer snow mantel scheme
!         
!     
!!**  METHOD
!     ------
!
!
! The temperature equation is written as:
!
!              b T+ = y
!
!
!!    EXTERNAL
!!    --------
!!
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!
!!    MODD_CST
!!
!!      
!!    REFERENCE
!!    ---------
!!
!!      
!!    AUTHOR
!!    ------
!!
!!      V. Masson           * Meteo-France *
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    08/09/98 
!!      J. Escobar 24/10/2012 : BUF PGI10.X , rewrite some 1 line WHERE statement
!!      V. Masson  13/09/2013 : implicitation of coupling with roof below
!-------------------------------------------------------------------------------
!
!*       0.     DECLARATIONS
!               ------------
!
USE modd_csts,       ONLY : xtt, xci, xrholi, xrholw, xcpd, xlstt, xlmtt, xday, xcondi
USE modd_snow_par,   ONLY : xemissn
USE modd_surf_par,   ONLY : xundef
!
USE mode_thermos
!
USE modi_surface_ri
USE modi_surface_aero_cond
!
!
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb
!
IMPLICIT NONE
!
!*      0.1    declarations of arguments
!
!
REAL,                 INTENT(IN)    :: ptstep   ! time step
REAL,                 INTENT(IN)    :: pansmin  ! minimum snow albedo
REAL,                 INTENT(IN)    :: pansmax  ! maximum snow albedo
REAL,                 INTENT(IN)    :: ptodry   ! snow albedo decreasing constant
REAL,                 INTENT(IN)    :: prhosmin ! minimum snow density
REAL,                 INTENT(IN)    :: prhosmax ! maximum snow density
REAL,                 INTENT(IN)    :: prhofold ! snow density increasing constant
LOGICAL,              INTENT(IN)    :: oall_melt! T --> all snow runs off if
                                                ! lower surf. temperature is
                                                ! positive
REAL,                 INTENT(IN)    :: pdrain_time ! drainage folding time (days)
REAL,                 INTENT(IN)    :: pwcrn    ! critical snow amount necessary
                                                ! to cover the considered surface
REAL,                 INTENT(IN)    :: pz0sn    ! snow roughness length for momentum
REAL,                 INTENT(IN)    :: pz0hsn   ! snow roughness length for heat
REAL, DIMENSION(:), INTENT(INOUT) :: pwsnow   ! snow reservoir (kg/m2)
REAL, DIMENSION(:), INTENT(INOUT) :: ptsnow   ! snow temperature
REAL, DIMENSION(:), INTENT(INOUT) :: pasnow   ! snow albedo
REAL, DIMENSION(:), INTENT(INOUT) :: prsnow   ! snow density
REAL, DIMENSION(:), INTENT(INOUT) :: pts_snow ! snow surface temperature
REAL, DIMENSION(:), INTENT(INOUT) :: pesnow   ! snow emissivity
REAL, DIMENSION(:), INTENT(IN)    :: ptg      ! underlying ground temperature
REAL, DIMENSION(:), INTENT(IN)    :: ptg_coefa! underlying ground temperature
REAL, DIMENSION(:), INTENT(IN)    :: ptg_coefb! implicit terms
REAL, DIMENSION(:), INTENT(IN)    :: pabs_sw  ! absorbed SW energy (Wm-2)
REAL, DIMENSION(:), INTENT(IN)    :: plw1     ! LW coef independant of TSNOW
                                              ! (Wm-2)     usually equal to:
                                              !      emis_snow * LW_down
                                              !
REAL, DIMENSION(:), INTENT(IN)    :: plw2     ! LW coef dependant   of TSNOW
                                              ! (Wm-2 K-4) usually equal to:
                                              ! -1 * emis_snow * stefan_constant
                                              !
REAL, DIMENSION(:), INTENT(IN)    :: pta      ! temperature at the lowest level
REAL, DIMENSION(:), INTENT(IN)    :: pqa      ! specific humidity
                                                ! at the lowest level
REAL, DIMENSION(:), INTENT(IN)    :: pvmod    ! module of the horizontal wind
REAL, DIMENSION(:), INTENT(IN)    :: pps      ! pressure at the surface
REAL, DIMENSION(:), INTENT(IN)    :: prhoa    ! air density
                                                ! at the lowest level
REAL, DIMENSION(:), INTENT(IN)    :: psr      ! snow rate
REAL, DIMENSION(:), INTENT(IN)    :: pzref    ! reference height of the first
                                              ! atmospheric level (temperature)
REAL, DIMENSION(:), INTENT(IN)    :: puref    ! reference height of the first
                                              ! atmospheric level (wind)
REAL, DIMENSION(:), INTENT(OUT)   :: prnsnow  ! net radiation over snow
REAL, DIMENSION(:), INTENT(OUT)   :: phsnow   ! sensible heat flux over snow
REAL, DIMENSION(:), INTENT(OUT)   :: plesnow  ! latent heat flux over snow
REAL, DIMENSION(:), INTENT(OUT)   :: pgsnow   ! flux under the snow
REAL, DIMENSION(:), INTENT(OUT)   :: pmelt    ! snow melting rate (kg/m2/s)
REAL, DIMENSION(:), INTENT(OUT)   :: pdqs_snow! heat storage inside snow
REAL, DIMENSION(:), INTENT(OUT)   :: pabs_lw ! absorbed LW rad by snow (W/m2)
!
!
!*      0.2    declarations of local variables
!
REAL :: zexpl = 0.
REAL :: zimpl = 1.
!
REAL, DIMENSION(SIZE(PWSNOW)) :: zexns, zexna, zdircoszw
REAL, DIMENSION(SIZE(PWSNOW)) :: zz0      ! roughness length for momentum
REAL, DIMENSION(SIZE(PWSNOW)) :: zz0h     ! roughness length forheat
!
REAL, DIMENSION(SIZE(PWSNOW)) :: zri      ! Richardson number
REAL, DIMENSION(SIZE(PWSNOW)) :: zac      ! aerodynamical conductance
REAL, DIMENSION(SIZE(PWSNOW)) :: zra      ! aerodynamical resistance
REAL, DIMENSION(SIZE(PWSNOW)) :: zch      ! drag coefficient for heat
REAL, DIMENSION(SIZE(PWSNOW)) :: zb, zy   ! coefficients in Ts eq.
REAL, DIMENSION(SIZE(PWSNOW)) :: zwsnow   ! snow before evolution
REAL, DIMENSION(SIZE(PWSNOW)) :: zsnow_hc ! snow heat capacity
REAL, DIMENSION(SIZE(PWSNOW)) :: zsnow_tc ! snow thermal conductivity
REAL, DIMENSION(SIZE(PWSNOW)) :: zsnow_d  ! snow depth
REAL, DIMENSION(SIZE(PWSNOW)) :: zmelt    ! snow melting rate (kg/m3/s)
REAL, DIMENSION(SIZE(PWSNOW)) :: zts_snow ! snow surface temperature
                                          ! at previous time-step
REAL, DIMENSION(SIZE(PWSNOW)) :: zqsat    ! specific humidity
!                                         ! for ice
REAL, DIMENSION(SIZE(PWSNOW)) :: zdqsat   ! d(specific humidity)/dT
!                                         ! for ice
!
REAL, DIMENSION(SIZE(PWSNOW)) :: zsr1, zsr2   ! norm. snow precip.
!
LOGICAL, DIMENSION(SIZE(PWSNOW)) :: gsnowmask ! where snow is
!                                             ! at previuos time-step
LOGICAL, DIMENSION(SIZE(PWSNOW)) :: gfluxmask ! where fluxes can
!                                             ! be computed at
!                                             ! new time-step
!                                             ! i.e. snow occurence
!                                             ! at previous time-step
!                                             ! OR snow fall
INTEGER, DIMENSION(SIZE(PWSNOW)) :: jsnowmask1, jsnowmask2, jsnowmask3 ! where snow is or not
!                                                                      ! at previuos time-step
INTEGER, DIMENSION(SIZE(PWSNOW)) :: jfluxmask ! where fluxes can
!                                             ! be computed at
!                                             ! new time-step
!                                             ! i.e. snow occurence
!                                             ! at previous time-step
!                                             ! OR snow fall
!
REAL :: zwsnow_min = 0.1 ! minimum value of snow content (kg/m2) for prognostic
!                        ! computations
!
REAL, DIMENSION(SIZE(PWSNOW)) :: zei_snow  ! internal energy of snow
REAL, DIMENSION(SIZE(PWSNOW)) :: zpei_snow ! internal energy of snow at t+
REAL, DIMENSION(SIZE(PWSNOW)) :: zwork1
REAL, DIMENSION(SIZE(PWSNOW)) :: zdqsati, zqsati
!
INTEGER                         :: jj, ji, jcompt_snow1, jcompt_snow2, jcompt_snow3, jcompt_flux
REAL(KIND=JPRB) :: zhook_handle
!-------------------------------------------------------------------------------
!
!
IF (lhook) CALL dr_hook('SNOW_COVER_1LAYER',0,zhook_handle)
zb(:)=0.
zy(:)=0.
zmelt(:) = 0.
pmelt(:) = 0.
prnsnow(:) = 0.
phsnow(:) = 0.
plesnow(:) = 0.
pgsnow(:) = 0.
!RJ: workaround to prevent decomposition unstable xundef masks for Tx_LWA_SN_RD fields
!RJ: in TEB_DIAGNOSTICS.nc during TEB_GARDEN_GREENROOF_BEM_3L_IRRIG_* tests
!RJ: problem with decomposition handling somewhere else
#ifdef RJ_PFIX
pabs_lw(:) = 0.0
#else
pabs_lw(:) = xundef
#endif
!
!* snow reservoir before evolution
!
zwsnow(:) = pwsnow(:)
zts_snow(:) = min(xtt,ptg(:))
!
zsnow_d(:) = 0.
zsnow_tc(:) = 0.
zsnow_hc(:) = 0.
!
!-------------------------------------------------------------------------------
!
!*      1.1    most useful masks
!              -----------------
!
gsnowmask(:)=.false.
gfluxmask(:)=.false.
jsnowmask1(:)=0.
jsnowmask2(:)=0.
jsnowmask3(:)=0.
jfluxmask(:)=0.

  !*      1.2    drag
!              ----
!
!*      1.2.1  defaults
!              --------
!
!* variation of temperature with altitude is neglected
!
zexns(:) = 1.
zexna(:) = 1.
!
!* slope is neglected in drag computation
!
zdircoszw(:) = 1.
!
!* roughness length are imposed:
!
zz0(:) = pz0sn
zz0h(:) = pz0hsn

!
!
!*      1.1    most useful masks
!              -----------------
!* snow occurence at previous time-step
!
!* snow occurence during the time-step for fluxes computation
!
jcompt_snow1=0
jcompt_snow2=0
jcompt_snow3=0
jcompt_flux=0
DO jj=1,SIZE(zwsnow)
  IF (zwsnow(jj)>0.) THEN
    gsnowmask(jj)=.true.
    !* surface temperature
    zts_snow(jj)=pts_snow(jj)
    gfluxmask(jj)=.true.
    !gsnowmask=t
    jcompt_snow1=jcompt_snow1+1
    jsnowmask1(jcompt_snow1) = jj
    !gfluxmask=t
    jcompt_flux=jcompt_flux+1
    jfluxmask(jcompt_flux) = jj
    IF (zwsnow(jj)>=zwsnow_min) THEN
      !second snow mask
      jcompt_snow3=jcompt_snow3+1
      jsnowmask3(jcompt_snow3)=jj
    ELSE
      !lower limit of snow cover for prognostic computations
      !0.1 kg/m2 of snow water content
      ptsnow(jj)=min(ptg(jj),xtt)
    ENDIF
  ELSE
    ptsnow(jj)=min(ptg(jj),xtt)
    !gsnowmask=false
    jcompt_snow2=jcompt_snow2+1
    jsnowmask2(jcompt_snow2) = jj
    IF (psr(jj)>0.) THEN
      gfluxmask(jj)=.true.
      jcompt_flux=jcompt_flux+1
      jfluxmask(jcompt_flux) = jj
    ENDIF
  ENDIF
ENDDO
!
!-------------------------------------------------------------------------------
!
!*      1.2    drag
!              ----
!
!*      1.2.2   qsat (Tsnow)
!              ------------
!
zqsat(:) = qsati(zts_snow(:), pps(:) )
!
!*      1.2.3  Richardson number
!              -----------------
!
!* snow is present on all the considered surface.
!* computation occurs where snow is and/or falls.
!
 CALL surface_ri(zts_snow, zqsat, zexns, zexna, pta, pqa, &
                  pzref, puref, zdircoszw, pvmod, zri      )  
!
!*      1.2.4  Aerodynamical conductance
!              -------------------------
!
 CALL surface_aero_cond(zri, pzref, puref, pvmod, zz0, zz0h, zac, zra, zch)
!
!-------------------------------------------------------------------------------
!
!*      2.     snow thermal characteristics
!              ----------------------------
!cdir nodep
DO jj=1,jcompt_snow1
  !
  ji = jsnowmask1(jj)
  !
  !*      2.1    snow heat capacity
  zsnow_hc(ji) = prsnow(ji) * xci * xrholi / xrholw
!*      2.2    snow depth
  zsnow_d(ji) = zwsnow(ji) / prsnow(ji)
!*      2.3    snow thermal conductivity
  zsnow_tc(ji) = xcondi * (prsnow(ji)/xrholw)**1.885
!*      2.4    internal energy of snow
  zei_snow(ji) = zsnow_hc(ji)*zsnow_d(ji)*ptsnow(ji)
  !
ENDDO
!
!cdir nodep
DO jj=1,jcompt_snow2
  !
  ji = jsnowmask2(jj)
  !
  !*      2.1    snow heat capacity
  zsnow_hc(ji) = prhosmin * xci * xrholi / xrholw
!*      2.2    snow depth
  zsnow_d(ji) = ptstep * psr(ji) / prhosmin
!*      2.3    snow thermal conductivity
  zsnow_tc(ji) = xcondi * (prhosmin /xrholw)**1.885
!*      2.4    internal energy of snow
  zei_snow(ji) = 0.
!
ENDDO
!
!-------------------------------------------------------------------------------
!
!*      3.     Snow temperature evolution
!              --------------------------
!
!*      3.5    dqsat/ dT (Tsnow)
!              -----------------
!
zdqsati = dqsati(zts_snow(:),pps(:),zqsat(:))
WHERE (gsnowmask(:) .AND. zwsnow(:)>=zwsnow_min)
  zdqsat(:) = zdqsati(:)
END WHERE
!
!*      3.1    coefficients from Temperature tendency
!              --------------------------------------
!
!cdir nodep
DO jj=1,jcompt_snow3
!
  ji=jsnowmask3(jj)

  zwork1(ji) = zsnow_d(ji) * zsnow_hc(ji) / ptstep
!
  zb(ji) = zb(ji) + zwork1(ji)
!
!*      3.2    coefficients from solar radiation
!          ---------------------------------
!    
  zy(ji) = zy(ji) + zwork1(ji) * ptsnow(ji) + pabs_sw(ji)
!
!
!*      3.3    coefficients from infra-red radiation
!              -------------------------------------
!
  zwork1(ji) = plw2(ji) * ptsnow(ji)**3
!
  zb(ji) = zb(ji) - 4 * zimpl * zwork1(ji)
!
  zy(ji) = zy(ji) + plw1(ji) + zwork1(ji) * (zexpl-3.*zimpl) * ptsnow(ji)
!
!
!*      3.4    coefficients from sensible heat flux
!              ------------------------------------
!
 zwork1(ji) = xcpd * prhoa(ji) * zac(ji)
! 
  zb(ji) = zb(ji) + zwork1(ji) *   zimpl
!
  zy(ji) = zy(ji) - zwork1(ji) * ( zexpl * ptsnow(ji) - pta(ji) )
!
!
!*      3.6    coefficients from latent heat flux
!              ----------------------------------
!
  zwork1(ji) =  xlstt * prhoa(ji) * zac(ji)
!
  zb(ji) = zb(ji) + zwork1(ji) *  zimpl * zdqsat(ji)
!
  zy(ji) = zy(ji) - zwork1(ji) * (  zqsat(ji) - pqa(ji) - zimpl * zdqsat(ji)*ptsnow(ji) )
!
!*      3.7    coefficients from conduction flux at snow base
!              ----------------------------------------------
!
  zwork1(ji) = zsnow_tc(ji)/(0.5*zsnow_d(ji))
!
  zb(ji) = zb(ji) + zwork1(ji) *  zimpl / ( 1. + zwork1(ji)*ptg_coefa(ji) )
!
  zy(ji) = zy(ji) - zwork1(ji) * (zexpl * ptsnow(ji) - ptg_coefb(ji)) &
                   / ( 1. + zwork1(ji)*ptg_coefa(ji) )
!
!*      3.8    guess of snow temperature before accumulation and melting
!              ---------------------------------------------------------
!
  ptsnow(ji) = zy(ji) / zb(ji)
!
ENDDO
!
!-------------------------------------------------------------------------------
!
!*      4.     Snow melt
!              ---------
!
!*      4.1    melting
!              -------
!
!cdir nodep
DO jj=1,jcompt_snow1
!
  ji = jsnowmask1(jj)
!
  zmelt(ji)  = max( ptsnow(ji) - xtt , 0. ) * zsnow_hc(ji) /  xlmtt / ptstep
!
  zmelt(ji)  = min( zmelt(ji) , zwsnow(ji) / zsnow_d(ji) / ptstep )
!
  ptsnow(ji) = min( ptsnow(ji) , xtt )
!
ENDDO
!
!*      4.2    run-off of all snow if lower surface temperature is positive
!              ------------------------------------------------------------
!
!* this option is used when snow is located on sloping roofs for example.
!
IF (oall_melt) THEN
  WHERE ( gsnowmask(:) .AND. ptg(:)>xtt .AND. zwsnow(:)>=zwsnow_min )
    pmelt(:) = pmelt(:) + zwsnow(:) / ptstep
  END WHERE
END IF
!
!*      4.3    output melting in kg/m2/s
!              -------------------------
!
pmelt(:) = zmelt(:) * zsnow_d(:)
!
!-------------------------------------------------------------------------------
!
!*      5.     fluxes
!              ------
!
!*      5.3    qsat (Tsnow)
!              ------------
!
zqsati = qsati(ptsnow(:),pps(:))
WHERE (gfluxmask(:)) 
   zqsat(:) = zqsati(:)
END WHERE
!
!*      5.1    net radiation (with Ts lin. extrapolation)
!              -------------
!
!cdir nodep
DO jj = 1, jcompt_flux
!
  ji = jfluxmask(jj)
!
  pabs_lw(ji) =  plw1(ji) + plw2(ji) * ptsnow(ji)**4
!
  prnsnow(ji) = pabs_sw(ji) + pabs_lw(ji)
!
!
!*      5.2    sensible heat flux
!              ------------------
!
  phsnow(ji) = xcpd * prhoa(ji) * zac(ji) * ( ptsnow(ji) - pta(ji) )
!
!
!*      5.4    latent heat flux
!              ----------------
!
  plesnow(ji) = xlstt * prhoa(ji) * zac(ji) * ( zqsat(ji) - pqa(ji) )
  !
!
!*      5.5    Conduction heat flux
!              --------------------
!
  !PGSNOW(JI) = ZSNOW_TC(JI)/(0.5*ZSNOW_D(JI)) * ( PTSNOW(JI) - PTG(JI) )
  pgsnow(ji) = zsnow_tc(ji)/(0.5*zsnow_d(ji)) * ( ptsnow(ji) - ptg_coefb(ji) ) &
             / ( 1. + zsnow_tc(ji)/(0.5*zsnow_d(ji))*ptg_coefa(ji) )
!
!
!*      5.6    If ground T>0 C, Melting is estimated from conduction heat flux
!              ---------------------------------------------------------------
!
  IF (ptg(ji)>xtt)  pmelt(ji) = max(pmelt(ji), -pgsnow(ji)/xlmtt)
!
ENDDO
!
!-------------------------------------------------------------------------------
!
!*      6.     reservoir evolution
!              -------------------
!
!cdir nodep
DO jj = 1, SIZE(pwsnow)
!
!*      6.1    snow fall
!              ---------
!
  pwsnow(jj) = pwsnow(jj) + ptstep * psr(jj)
!
!
!*      6.2    sublimation
!              -----------
!
  plesnow(jj) = min( plesnow(jj), xlstt*pwsnow(jj)/ptstep )
!
  pwsnow(jj)  = max( pwsnow(jj) - ptstep * plesnow(jj)/xlstt , 0.)
!
  IF ( pwsnow(jj)<1.e-8 * ptstep ) pwsnow(jj) = 0.
!
!*      6.3    melting
!              -------
!
  pmelt(jj) = min( pmelt(jj), pwsnow(jj)/ptstep )
!
  pwsnow(jj)= max( pwsnow(jj) - ptstep * pmelt(jj) , 0.)
!
  IF ( pwsnow(jj)<1.e-8 * ptstep ) pwsnow(jj) = 0.
!
  IF (pwsnow(jj)==0.) pgsnow(jj) = max( pgsnow(jj), - pmelt(jj)*xlmtt )
!
ENDDO
!
!*      6.4    time dependent drainage
!              -----------------------
!
IF (pdrain_time>0.) THEN
  WHERE ( pwsnow(:)>0.)
    pwsnow(:) = pwsnow(:) * exp(-ptstep/pdrain_time/xday)
  END WHERE
END IF
!
!*      6.5    melting of last 1mm of snow depth
!              ---------------------------------
!
WHERE ( pwsnow(:)<zwsnow_min .AND. pmelt(:)>0. .AND. psr(:)==0. )
  pmelt(:) = pmelt(:) + pwsnow(:) / ptstep
  pwsnow(:)=0.
END WHERE
!
WHERE ( pwsnow(:)<1.e-8 * ptstep ) 
   pwsnow(:) = 0.
END WHERE
!
!-------------------------------------------------------------------------------
!
!*      7.     albedo evolution
!              ----------------
!
!*      7.1    If melting occurs or not
!              -----------------------
!
!
!cdir nodep
DO jj=1,jcompt_snow1
!
  ji = jsnowmask1(jj)
!
  IF (pmelt(ji) > 0. ) THEN
!
    pasnow(ji) = (pasnow(ji)-pansmin)*exp(-prhofold*ptstep/xday) + pansmin   &
                + psr(ji)*ptstep/pwcrn*pansmax  
!
  ELSEIF (pmelt(ji)==0.) THEN 
    pasnow(ji) = pasnow(ji) - ptodry*ptstep/xday                          &
                + psr(ji)*ptstep/pwcrn*pansmax  
!
  ENDIF
!
ENDDO
!
!-------------------------------------------------------------------------------
!
!*      8.     density evolution
!              -----------------
!
!*      8.1    old snow
!              --------
!
!cdir nodep
DO jj = 1, jcompt_snow1
! 
  ji = jsnowmask1(jj)
!
  IF (pwsnow(ji)>0. ) THEN
!
    zsr1(ji) = max( pwsnow(ji) , psr(ji) * ptstep )
!
    prsnow(ji) = (prsnow(ji)-prhosmax)*exp(-prhofold*ptstep/xday) + prhosmax
    prsnow(ji) = ( (zsr1(ji)-psr(ji)*ptstep) * prsnow(ji)    &
                  + (psr(ji)*ptstep) * prhosmin ) / zsr1(ji)  
  ENDIF
!
ENDDO
!
!*      8.2    fresh snow
!              ----------
!
!cdir nodep
DO jj=1,SIZE(pwsnow)
  IF (  pwsnow(jj)>0. ) THEN
    pasnow(jj) = max(pasnow(jj),pansmin)
    pasnow(jj) = min(pasnow(jj),pansmax)
    IF (zwsnow(jj)==0.) THEN
      pasnow(jj) = pansmax
      pesnow(jj) = xemissn
      prsnow(jj) = prhosmin
    ENDIF
  ENDIF
    ENDDO
!
!-------------------------------------------------------------------------------
!
!*      9.     fresh snow accumulation (if more than 1mm of snow depth)
!              -----------------------
!
!cdir nodep
DO jj=1,jcompt_snow3
!
  ji = jsnowmask3(jj)
!
  IF (psr(ji)>0. .AND. pwsnow(ji)>0.) THEN
!
    zsr2(ji) = min( pwsnow(ji) , psr(ji) * ptstep )
!
    ptsnow(ji) =( ( pwsnow(ji) - zsr2(ji) ) *      ptsnow(ji)        &
              +                zsr2(ji)   * min( pta(ji) ,xtt ))&
                /(   pwsnow(ji) )  
  ENDIF
!
ENDDO
!
!-------------------------------------------------------------------------------
!
!*     10.     Surface temperature
!              -------------------
!
!* note that if the relation between snow pack temperature and its
!  surface temperature is modified, think to modify it also in
!  subroutine init_snow_lw.f90
!
WHERE (gsnowmask(:) )
  pts_snow(:) = ptsnow(:)
END WHERE
!
!-------------------------------------------------------------------------------
!
!*     11.     bogus values
!              ------------
!
!*     11.1    snow characteristics where snow IS present at current time-step
!              ---------------------------------------------------------------
!
WHERE (pwsnow(:)==0.)
  ptsnow(:) = xundef
  prsnow(:) = xundef
  pasnow(:) = xundef
  pts_snow(:) = xundef
  pesnow(:) = xundef
END WHERE
!
!
!-------------------------------------------------------------------------------
!
!*     12.     Heat storage inside snow pack
!
WHERE (gsnowmask(:))
  zpei_snow(:) = zsnow_hc(:)*zsnow_d(:)*ptsnow(:)
ELSEWHERE
  zpei_snow(:) = 0.
END WHERE
pdqs_snow(:) = (zpei_snow(:)-zei_snow(:))/ptstep
!
IF (lhook) CALL dr_hook('SNOW_COVER_1LAYER',1,zhook_handle)

!------------------------------------------------------------------------------- !
END SUBROUTINE snow_cover_1layer