isba_snow_agr.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 isba_snow_agr( HSNOW_ISBA, OMEB,                                &
                         pexns, pexna, pta, pqa, pzref, puref, pdircoszw, pvmod, &
                         pz0eff, pz0, pz0h, prr, psr,                            &
                         pemis, palb, ppsn, ppsng, ppsnv,                        &
                         prn, ph, ple, plei, pleg, plegi, plev, ples, pler,      &
                         pletr, pevap, psubl, pgflux, plvtt, plstt,              &
                         pustar,                                                 &
                         ples3l, plel3l, pevap3l,                                &
                         pswnet_v, pswnet_g, plwnet_v, plwnet_g, ph_v, ph_g,     &
                         plev_v_c, pletr_v_c, ples_v_c,                          &
                         pqs3l, palb3l,                                          &
                         prnsnow, phsnow, phpsnow,                               &
                         pswnetsnow, pswnetsnows, plwnetsnow,                    &
                         pgfluxsnow, pgsfcsnow, pustarsnow,                      &
                         pzgrndflux, pflsn_cor, pgrndflux, plesl,                &
                         pemisnow,                                               &
                         psnowtemp, pts_rad, pts, pri, pqs, phu,                 &
                         pcd, pcdn, pch, psnowhmass,                             &
                         prn_isba, ph_isba, pleg_isba, plegi_isba, plev_isba,    &
                         pletr_isba, pustar_isba, pler_isba, ple_isba,           &
                         plei_isba, pgflux_isba, pmeltadv, ptg,                  &
                         pemist, palbt, ple_flood, plei_flood,                   &
                         pffg, pffv, pff, ppalphan, ptc, omeb_litter, plelitter, &
                         plelitteri                                              )
!     ##########################################################################
!
!
!!****  *ISBA_SNOW_AGR* aggregates snow free and snow fluxes
!!
!!    PURPOSE
!!    -------
!     
!!**  METHOD
!!    ------
!
!!    EXTERNAL
!!    --------
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------ 
!!
!!      
!!    REFERENCE
!!    ---------
!!
!!    Noilhan and Planton (1989)
!!      
!!    AUTHOR
!!    ------
!!      V. Masson           * Meteo-France *
!!      (following A. Boone)
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    10/03/95 
!!      B. Decharme 01/2009  Floodplains 
!!      B. Decharme 01/2010  Effective surface temperature (for diag)
!!      B. Decharme 09/2012  Bug total sublimation flux: no PLESL
!!      B. Decharme 04/2013  Bug wrong radiative temperature
!!                           Sublimation diag flux
!!                           Qs for 3l or crocus (needed for coupling with atm)
!!      A. Boone    11/2014  MEB
!-------------------------------------------------------------------------------
!
!*       0.     DECLARATIONS
!               ------------
USE modd_surf_par,   ONLY : xundef
!
!
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb
!
IMPLICIT NONE
!
!*      0.1    declarations of arguments
!              -------------------------
!
!
!* general variables
!  -----------------
!
 CHARACTER(LEN=*),     INTENT(IN)  :: hsnow_isba ! 'DEF' = Default F-R snow scheme
!                                               !         (Douville et al. 1995)
!                                               ! '3-L' = 3-L snow scheme (option)
!                                               !         (Boone and Etchevers 2000)
!
LOGICAL,              INTENT(IN)  :: omeb       ! True = patch with multi-energy balance 
!                                               ! False = patch with classical ISBA 
!
!* surface and atmospheric parameters
!  ----------------------------------
!
REAL, DIMENSION(:), INTENT(IN)  :: pexns     ! Exner function at the surface
REAL, DIMENSION(:), INTENT(IN)  :: pexna     ! Exner function
REAL, DIMENSION(:), INTENT(IN)  :: pta       ! air temperature
REAL, DIMENSION(:), INTENT(IN)  :: pqa       ! air specific humidity
REAL, DIMENSION(:), INTENT(IN)  :: pzref     ! reference height of the first atmospheric level
REAL, DIMENSION(:), INTENT(IN)  :: puref     ! reference height of the wind
REAL, DIMENSION(:), INTENT(IN)  :: pdircoszw ! Cosinus of the angle between the normal to the surface and the vertical
REAL, DIMENSION(:), INTENT(IN)  :: pvmod     ! module of the horizontal wind
REAL, DIMENSION(:), INTENT(IN)  :: pz0eff    ! effective roughness length
REAL, DIMENSION(:), INTENT(IN)  :: pz0       ! roughness length for momentum
REAL, DIMENSION(:), INTENT(IN)  :: pz0h      ! roughness length for heat
REAL, DIMENSION(:), INTENT(IN)  :: prr       ! Rain rate (in kg/m2/s)
REAL, DIMENSION(:), INTENT(IN)  :: psr       ! Snow rate (in kg/m2/s)
!
!* surface parameters
!  ------------------
!
REAL, DIMENSION(:), INTENT(IN)  :: palb       ! albedo 
REAL, DIMENSION(:), INTENT(IN)  :: pemis      ! emissivity
!  'D95'     : represents aggregated (snow + flood + snow-flood-free) albedo and emissivity
!  '3-L'     : represents                    flood + snow-flood-free  albedo and emissivity
!  'MEB+3-L' : represents aggregated (snow + flood + snow-flood-free) albedo and emissivity
!
!
!* snow fractions
!  --------------
!
REAL, DIMENSION(:), INTENT(IN)  :: ppsn       ! fraction of the grid covered
!                                             ! by snow
REAL, DIMENSION(:), INTENT(IN)  :: ppsng      ! fraction of the the bare
!                                             ! ground covered by snow
REAL, DIMENSION(:), INTENT(IN)  :: ppsnv      ! fraction of the the veg.
!                                             ! covered by snow
REAL, DIMENSION(:), INTENT(IN)  :: ppalphan   ! fraction of the the explicit veg.
!                                             ! canopy buried by snow
!
!
!* ISBA-SNOW3L variables/parameters:
!  ---------------------------------
!
! Prognostic variables:
!
REAL, DIMENSION(:),   INTENT(IN) :: palb3l      ! Snow albedo
REAL, DIMENSION(:),   INTENT(IN) :: pqs3l       ! Surface humidity
! 
! Diagnostics:
!
REAL, DIMENSION(:), INTENT(IN)    :: pzgrndflux ! snow/soil-biomass interface flux (W/m2)
REAL, DIMENSION(:), INTENT(IN)    :: pflsn_cor  ! snow/soil-biomass correction flux (W/m2)
!
REAL, DIMENSION(:), INTENT(INOUT) :: pgrndflux  ! snow/soil-biomass interface flux (W/m2)
!
REAL, DIMENSION(:), INTENT(INOUT) :: phpsnow    ! heat release from rainfall (W/m2)
REAL, DIMENSION(:), INTENT(INOUT) :: psnowhmass ! snow heat content change from mass changes (J/m2)
REAL, DIMENSION(:), INTENT(INOUT) :: prnsnow    ! net radiative flux from snow (W/m2)
REAL, DIMENSION(:), INTENT(INOUT) :: pswnetsnow ! net shortwave snow radiative flux (W/m2)
REAL, DIMENSION(:), INTENT(INOUT) :: pswnetsnows! net shortwave snow radiative flux in sfc layer (W/m2)
REAL, DIMENSION(:), INTENT(INOUT) :: plwnetsnow ! net longwave snow radiative flux (W/m2)
REAL, DIMENSION(:), INTENT(INOUT) :: phsnow     ! sensible heat flux from snow (W/m2)
REAL, DIMENSION(:), INTENT(INOUT) :: pgfluxsnow ! net heat flux from snow (W/m2)
REAL, DIMENSION(:), INTENT(INOUT) :: pgsfcsnow  ! heat flux from snow sfc to sub sfc layers (W/m2)
REAL, DIMENSION(:), INTENT(INOUT) :: pswnet_v   ! net shortwave radiation of vegetation canopy 
REAL, DIMENSION(:), INTENT(INOUT) :: pswnet_g   ! net shortwave radiation of (below snow) surface
REAL, DIMENSION(:), INTENT(INOUT) :: plwnet_v   ! net longwave radiation of vegetation canopy 
REAL, DIMENSION(:), INTENT(INOUT) :: plwnet_g   ! net longwave radiation of (below snow) surface
REAL, DIMENSION(:), INTENT(IN)    :: pustarsnow ! friction velocity
REAL, DIMENSION(:), INTENT(OUT)   :: plesl      ! Evaporation (liquid) from wet snow (W/m2)
REAL, DIMENSION(:), INTENT(IN)    :: pemisnow   ! snow surface emissivity
REAL, DIMENSION(:), INTENT(OUT)   :: pts_rad    ! effective radiative temperature 
!                                                 of the natural surface (K)
REAL, DIMENSION(:), INTENT(OUT)   :: pts        ! effective surface temperature 
REAL, DIMENSION(:), INTENT(IN)    :: psnowtemp  ! snow layer sfc temperature (K)
REAL, DIMENSION(:), INTENT(IN)    :: ples3l     ! sublimation from ISBA-ES(3L)
REAL, DIMENSION(:), INTENT(IN)    :: plel3l     ! evaporation heat flux of water in the snow (W/m2)
REAL, DIMENSION(:), INTENT(INOUT) :: pevap3l    ! evaporation flux over snow from ISBA-ES (kg/m2/s)
REAL, DIMENSION(:), INTENT(IN)    :: plvtt, plstt    
!
!
! Prognostic variables:
!
REAL, DIMENSION(:),   INTENT(IN)  :: ptg              ! soil sfc layer average temperatures    (K)
REAL, DIMENSION(:),   INTENT(IN)  :: ptc              ! canopy air temperature                 (K)
!
!
!* diagnostic variables
!  --------------------
!
REAL, DIMENSION(:), INTENT(INOUT) :: pemist   ! total surface emissivity
REAL, DIMENSION(:), INTENT(INOUT) :: palbt    ! total surface albedo
!
!* surface fluxes
!  --------------
!
REAL, DIMENSION(:), INTENT(IN)    :: plev_v_c ! MEB: total evapotranspiration (no snow) from 
!                                             !  vegetation canopy overstory [W/m2]
REAL, DIMENSION(:), INTENT(IN)    :: ples_v_c ! MEB: total (intercepted) snow sublimation from 
!                                             !  vegetation canopy overstory [W/m2]
REAL, DIMENSION(:), INTENT(IN)    :: pletr_v_c! MEB: transpiration from overstory canopy 
!                                             !  vegetation [W/m2]
REAL, DIMENSION(:), INTENT(INOUT) :: prn      ! net radiation
REAL, DIMENSION(:), INTENT(INOUT) :: ph       ! sensible heat flux
REAL, DIMENSION(:), INTENT(INOUT) :: ph_v     ! sensible heat flux from explicit veg canopy
REAL, DIMENSION(:), INTENT(INOUT) :: ph_g     ! sensible heat flux from surface (below snow)
REAL, DIMENSION(:), INTENT(INOUT) :: ple      ! total latent heat flux
REAL, DIMENSION(:), INTENT(OUT)   :: plei     ! sublimation latent heat flux
REAL, DIMENSION(:), INTENT(INOUT) :: plegi    ! latent heat of sublimation over frozen soil
REAL, DIMENSION(:), INTENT(INOUT) :: pleg     ! latent heat of evaporation
REAL, DIMENSION(:), INTENT(IN)    :: plelitteri! sublimation of water in litter reservoir
REAL, DIMENSION(:), INTENT(IN)    :: plelitter !sublimation of water in litter reservoir
LOGICAL, INTENT(IN)               :: omeb_litter !True = litter option activated
!                                             ! over the ground
REAL, DIMENSION(:), INTENT(INOUT) :: plev     ! latent heat of evaporation
!                                             ! over the vegetation
REAL, DIMENSION(:), INTENT(INOUT) :: ples     ! latent heat of sublimation
!                                             ! over the snow
REAL, DIMENSION(:), INTENT(INOUT) :: pler     ! latent heat of the fraction
!                                             ! delta of water retained on the
!                                             ! foliage of the vegetation
REAL, DIMENSION(:), INTENT(INOUT) :: pletr    ! evapotranspiration of the rest
!                                             ! of the vegetation
REAL, DIMENSION(:), INTENT(INOUT) :: pevap    ! total evaporative flux (kg/m2/s)
REAL, DIMENSION(:), INTENT(INOUT) :: psubl    ! sublimation flux (kg/m2/s)
REAL, DIMENSION(:), INTENT(INOUT) :: pgflux   ! flux through the ground
REAL, DIMENSION(:), INTENT(INOUT) :: pustar   ! friction velocity
REAL, DIMENSION(:), INTENT(INOUT) :: pmeltadv ! advection heat flux from snowmelt (W/m2)
!
! The following surface fluxes are from snow-free portion of grid
! box when the ISBA-ES option is ON. Otherwise, they are equal
! to the same variables without the _ISBA extension.
!
REAL, DIMENSION(:), INTENT(OUT) :: prn_isba   ! net radiation
REAL, DIMENSION(:), INTENT(OUT) :: ph_isba    ! sensible heat flux
REAL, DIMENSION(:), INTENT(OUT) :: pleg_isba  ! latent heat of evaporation (ground)
REAL, DIMENSION(:), INTENT(OUT) :: plegi_isba ! latent heat of sublimation (ground)
REAL, DIMENSION(:), INTENT(OUT) :: plev_isba  ! latent heat of evaporation (vegetation)
REAL, DIMENSION(:), INTENT(OUT) :: pletr_isba ! latent heat of evaporation (transpiration)
REAL, DIMENSION(:), INTENT(OUT) :: pustar_isba! friction velocity
REAL, DIMENSION(:), INTENT(OUT) :: pler_isba  ! latent heat of evaporation (plant interception)
REAL, DIMENSION(:), INTENT(OUT) :: ple_isba   ! total latent heat flux 
REAL, DIMENSION(:), INTENT(OUT) :: plei_isba  ! sublimation latent heat flux 
REAL, DIMENSION(:), INTENT(OUT) :: pgflux_isba! flux through the ground
!
REAL, DIMENSION(:), INTENT(IN)    :: pffg,pffv,pff
REAL, DIMENSION(:), INTENT(INOUT) :: ple_flood, plei_flood ! Flood evaporation
!
REAL, DIMENSION(:), INTENT(INOUT) :: pri       ! grid-area Ridcharson number
REAL, DIMENSION(:), INTENT(INOUT) :: pqs       ! grid-area Surface humidity
REAL, DIMENSION(:), INTENT(INOUT) :: phu       ! grid-area near surface humidity
REAL, DIMENSION(:), INTENT(INOUT) :: pch       ! grid-area drag coefficient for heat
REAL, DIMENSION(:), INTENT(INOUT) :: pcd       ! grid-area drag coefficient for momentum
REAL, DIMENSION(:), INTENT(INOUT) :: pcdn      ! grid-area neutral drag coefficient for momentum
!
!
!*      0.2    declarations of local variables
!
REAL, DIMENSION(SIZE(PTA))       :: zwork
!
REAL(KIND=JPRB) :: zhook_handle
!
!-------------------------------------------------------------------------------
!
IF (lhook) CALL dr_hook('ISBA_SNOW_AGR',0,zhook_handle)
!
zwork(:) = 0.
!
IF(omeb)THEN
!
! Snow free (ground-based snow) diagnostics: canopy and ground blended (W m-2):
! NOTE that the effects of snow cover *fraction* are implicitly *included* in these fluxes 
! so do NOT multiply by snow fraction.

   prn_isba(:)    = pswnet_v(:) + pswnet_g(:) + plwnet_v(:) + plwnet_g(:)
   ph_isba(:)     = ph_v(:) + ph_g(:)
   IF (omeb_litter) THEN
    pleg_isba(:)   = plelitter(:)
    plegi_isba(:)  = plelitteri(:)
    pleg(:)        = plelitter(:)
    plegi(:)       = plelitteri(:)
   ELSE
    pleg_isba(:)   = pleg(:)
    plegi_isba(:)  = plegi(:)
   ENDIF
   plei_isba(:)   = plegi(:) + plei_flood(:) + ples(:) + ples_v_c(:)
   plev_isba(:)   = plev_v_c(:)
   pletr_isba(:)  = pletr_v_c(:) 
   pustar_isba(:) = pustar(:)                  ! NOTE for now, this is same as total Ustar (includes snow)
! LER does not include intercepted snow sublimation
   pler_isba(:)   = plev_v_c(:) - pletr_v_c(:) 
! LE includes intercepted snow sublimation
   ple_isba(:)    = pleg_isba(:) + plegi_isba(:) + plev_isba(:) + ple_flood(:) + ples_v_c(:) + plei_flood(:)
   pgflux_isba(:) = prn_isba(:) - ph_isba(:) - ple_isba(:)
!
   pemist(:)      = pemis(:)
!
! Effective surface temperature (for diag): for MEB:

   zwork(:)       =  ppalphan(:)*ppsn(:)
   pts(:)         = (1.0 - zwork(:))*ptc(:) + zwork(:)*psnowtemp(:)
!
! Total heat FLUX into snow/soil/vegetation surface:
!
   pgflux(:)      = prn(:) - ph(:) - ple(:) + phpsnow(:) 
!
ELSE
!
! * 2. Using an explicit snow scheme option with composite soil/veg ISBA:
!      ------------------------------------------------------------------
!
   IF(hsnow_isba == '3-L' .OR. hsnow_isba == 'CRO')THEN
!
!     Save fluxes from Force-Restore snow/explicit snow-free
!     portion of grid box (vegetation/soil):
!
      prn_isba(:)    = prn(:)
      ph_isba(:)     = ph(:)
      pleg_isba(:)   = pleg(:)
      plegi_isba(:)  = plegi(:)
      plev_isba(:)   = plev(:)
      pletr_isba(:)  = pletr(:)
      pustar_isba(:) = pustar(:)
      pler_isba(:)   = pler(:) 
      ple_isba(:)    = ple(:)
      pgflux_isba(:) = pgflux(:)
!  
      plei_isba(:)   = plegi(:)+plei_flood(:)+ples(:)
!
!     Effective surface temperature (for diag):
!
      pts(:)       = (1.-ppsn(:))*ptg(:)+ppsn(:)*psnowtemp(:)
!
!     Effective surface radiating temperature:
!
      palbt(:)    = palb(:)*(1.-ppsn(:)) + ppsn(:)*palb3l(:)
      pemist(:)    = pemis(:)*(1.-ppsn(:)) + ppsn(:)*pemisnow(:)
!  
      pts_rad(:)   = ( ((1.-ppsn(:))*pemis(:)*ptg(:)**4 +   &
                            ppsn(:) *pemisnow(:)*psnowtemp(:)**4     &
                            )/pemist(:) )**(0.25)  
!
!     Calculate actual fluxes from snow-free natural
!     portion of surface: NET flux from surface is the sum of
!     fluxes from snow free and snow covered portions
!     of natural portion of grid box when *ISBA-ES* in force.
!     when NOT in use, then these fluxes equal those above.
!
      prn(:)       = (1.-ppsn(:))  * prn(:)   + ppsn(:) * prnsnow(:)
      ph(:)        = (1.-ppsn(:))  * ph(:)    + ppsn(:) * phsnow(:)
!  
      pleg(:)      = (1.-ppsng(:)-pffg(:)) * pleg(:)  
      plegi(:)     = (1.-ppsng(:)-pffg(:)) * plegi(:)  
      plev(:)      = (1.-ppsnv(:)-pffv(:)) * plev(:)   
      pletr(:)     = (1.-ppsnv(:)-pffv(:)) * pletr(:)  
      pler(:)      = (1.-ppsnv(:)-pffv(:)) * pler(:)  
!
!     Total evapotranspiration flux (kg/m2/s):
!
      pevap(:)     = (plev(:) + pleg(:))/plvtt(:) + plegi(:)/plstt(:) + ple_flood(:)/plvtt(:) + &
                      plei_flood(:)/plstt(:) + ppsn(:) * pevap3l(:)
!
!     ISBA-ES/SNOW3L fluxes:
!
      ples(:)       =                           ppsn(:) * ples3l(:)
      plesl(:)      =                           ppsn(:) * plel3l(:)
      prnsnow(:)    =                           ppsn(:) * prnsnow(:)
      phsnow(:)     =                           ppsn(:) * phsnow(:)
      pgfluxsnow(:) =                           ppsn(:) * pgfluxsnow(:)
      psnowhmass(:) =                           ppsn(:) * psnowhmass(:)  ! (J m-2)
      phpsnow(:)    =                           ppsn(:) * phpsnow(:)
      pgsfcsnow(:)  =                           ppsn(:) * pgsfcsnow(:)
      pswnetsnow(:) =                           ppsn(:) * pswnetsnow(:)
      pswnetsnows(:)=                           ppsn(:) * pswnetsnows(:)
      pevap3l(:)    =                           ppsn(:) * pevap3l(:)
!
!     Total heat flux between snow and soil
!
      pgrndflux(:) =                            ppsn(:) * (pzgrndflux(:)+pflsn_cor(:))
      pmeltadv(:)  =                            ppsn(:) * pmeltadv(:)
!
!     Total evaporative flux (W/m2) :
!
      ple(:)       = pleg(:) + plev(:) + ples(:) + plesl(:) + plegi(:) + ple_flood(:) + plei_flood(:)
!
!     Total sublimation flux (W/m2) :
!
      plei(:)      = ples(:) + plegi(:) + plei_flood(:)
!
!     Total sublimation flux (kg/m2/s):
!
      psubl(:)     = plei(:)/plstt(:)
!
!     Total FLUX into snow/soil/vegetation surface:
!
      pgflux(:)    = prn(:) - ph(:) - ple(:) + phpsnow(:)  
!
!     surface humidity:
!
      pqs(:)       = (1.-ppsn(:))  * pqs(:)   + ppsn(:) * pqs3l(:)
!
!     near-surface humidity :
!  
      phu(:)       = (1.-ppsn(:))  * phu(:)   + ppsn(:)
!
!     Momentum fluxes:
!
      pustar(:)    = sqrt( (1.-ppsn(:))  * pustar(:)**2  + ppsn(:) * pustarsnow(:)**2 )
!
!     Richardson number and Drag coeff:
!
      CALL comput_ri_drag
!
   ELSE
!
      pts(:)  = ptg(:)
      pts_rad(:)  = ptg(:)
      palbt(:)  = palb(:)
      pemist(:)  = pemis(:)
!  
!     Total sublimation flux (W/m2) :
      plei(:)  = ples(:) + plegi(:) + plei_flood(:)
!
!     Total sublimation flux (kg/m2/s):
      psubl(:)  = plei(:)/plstt(:)
!
   ENDIF
!
ENDIF
!
IF (lhook) CALL dr_hook('ISBA_SNOW_AGR',1,zhook_handle)
!
!-------------------------------------------------------------------------------
 CONTAINS
!-------------------------------------------------------------------------------
!
SUBROUTINE comput_ri_drag
!
USE modd_surf_atm, ONLY : ldrag_coef_arp, lrrgust_arp,   &
                          xrrscale, xrrgamma, xutilgust  
!
USE modi_surface_ri
USE modi_surface_aero_cond
USE modi_surface_cd
USE modi_surface_cdch_1darp
USE modi_wind_threshold
!
!*      0.2    declarations of local variables
!
REAL, DIMENSION(SIZE(PTA)) :: zfp, zrrcor, zvmod, zac, zra
!
REAL(KIND=JPRB) :: zhook_handle
!
!-------------------------------------------------------------------------------
!
IF (lhook) CALL dr_hook('ISBA_SNOW_AGR:COMPUT_RI_DRAG',0,zhook_handle)
!
! * Richardson number
!
 CALL surface_ri(pts, pqs, pexns, pexna, pta, pqa,  &
                pzref, puref, pdircoszw, pvmod, pri)  
!
! * Wind check
!
zvmod = wind_threshold(pvmod,puref)
!
! * Drag coefficient for heat and momentum
!
IF (ldrag_coef_arp) THEN
   CALL surface_cdch_1darp(pzref, pz0eff, pz0h, zvmod, pta, ptg, &
                             pqa, pqs, pcd, pcdn, pch              )
ELSE
   CALL surface_aero_cond(pri, pzref, puref, zvmod, pz0, pz0h, zac, zra, pch)
   CALL surface_cd(pri, pzref, puref, pz0eff, pz0h, pcd, pcdn)
ENDIF
!
IF (lrrgust_arp) THEN
   zfp(:)=max(0.0,prr(:)+psr(:))
   zrrcor(:)=sqrt(1.0+((((zfp(:)/(zfp(:)+xrrscale))**xrrgamma)*xutilgust)**2) &
       /(pcd(:)*zvmod(:)**2))  
   pcd  = pcd  * zrrcor
   pch  = pch  * zrrcor
   pcdn = pcdn * zrrcor
ENDIF
!
IF (lhook) CALL dr_hook('ISBA_SNOW_AGR:COMPUT_RI_DRAG',1,zhook_handle)
!
END SUBROUTINE comput_ri_drag
!
!-------------------------------------------------------------------------------
!
END SUBROUTINE isba_snow_agr