prep_hor_snow_fields.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 prep_hor_snow_fields (DTCO, &
                                  ig, u, &
                                 hprogram,hsurf,              &
                                hfile,hfiletype,             &
                                hfilepgd,hfilepgdtype,       &
                                kluout,ounif,kpatch,         &
                                kteb_patch, &
                                kl,tpsnow, tptime,           &
                                punif_wsnow, punif_rsnow,    &
                                punif_tsnow, punif_lwcsnow,  &
                                punif_asnow, osnow_ideal,    &
                                punif_sg1snow, punif_sg2snow,&
                                punif_histsnow,punif_agesnow,&
                                pvegtype, pvegtype_patch,    &
                                ppatch, okey                 )  
!     #######################################################
!
!
!!****  *PREP_HOR_SNOW_FIELDS* - prepares all snow fields for one surface scheme.
!!
!!    PURPOSE
!!    -------
!
!!**  METHOD
!!    ------
!!
!!    REFERENCE
!!    ---------
!!      
!!
!!    AUTHOR
!!    ------
!!     V. Masson 
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    01/2004
!!      B. Decharme 10/2013, Phasage Arpege-Climat
!!      B. Decharme 04/2014, Init permsnow
!!------------------------------------------------------------------
!
!
!
!
!
!
USE modd_data_cover_n, ONLY : data_cover_t
!
USE modd_isba_grid_n, ONLY : isba_grid_t
USE modd_surf_atm_n, ONLY : surf_atm_t
!
USE modd_type_snow
USE modd_type_date_surf, ONLY : date_time
!
USE modd_surf_par,       ONLY : xundef
USE modd_snow_par,       ONLY : xaglamin, xaglamax
!
USE modd_data_cover_par, ONLY : nvegtype
!
USE modi_allocate_gr_snow
USE modi_prep_hor_snow_field
USE mode_snow3l
USE modi_open_aux_io_surf
USE modi_read_surf
USE modi_close_aux_io_surf
!
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb
!
IMPLICIT NONE
!
!*      0.1    declarations of arguments
!
!
!
TYPE(data_cover_t), INTENT(INOUT) :: dtco
!
TYPE(isba_grid_t), INTENT(INOUT) :: ig
TYPE(surf_atm_t), INTENT(INOUT) :: u
!
 CHARACTER(LEN=6),   INTENT(IN)  :: hprogram  ! program calling surf. schemes
 CHARACTER(LEN=7),   INTENT(IN)  :: hsurf     ! type of field
 CHARACTER(LEN=28),  INTENT(IN)  :: hfile     ! file name
 CHARACTER(LEN=6),   INTENT(IN)  :: hfiletype ! file type
 CHARACTER(LEN=28),  INTENT(IN)  :: hfilepgd     ! file name
 CHARACTER(LEN=6),   INTENT(IN)  :: hfilepgdtype ! file type
INTEGER,            INTENT(IN)  :: kluout    ! logical unit of output listing
LOGICAL,            INTENT(IN)  :: ounif     ! flag for prescribed uniform field
INTEGER,            INTENT(IN)  :: kpatch    ! patch number for output scheme
INTEGER,            INTENT(IN) :: kteb_patch
INTEGER,            INTENT(IN)  :: kl        ! number of points
TYPE(surf_snow)                 :: tpsnow    ! snow fields
TYPE(date_time),    INTENT(IN)  :: tptime    ! date and time
REAL, DIMENSION(:), INTENT(IN)  :: punif_wsnow ! prescribed snow content (kg/m2)
REAL, DIMENSION(:), INTENT(IN)  :: punif_rsnow ! prescribed density (kg/m3)
REAL, DIMENSION(:), INTENT(IN)  :: punif_tsnow ! prescribed temperature (K)
REAL, DIMENSION(:), INTENT(IN)  :: punif_lwcsnow ! prescribed snow liquid water content (kg/m3)
REAL,               INTENT(IN)  :: punif_asnow ! prescribed albedo (-)
LOGICAL,            INTENT(IN)  :: osnow_ideal
REAL, DIMENSION(:), INTENT(IN)  :: punif_sg1snow ! 
REAL, DIMENSION(:), INTENT(IN)  :: punif_sg2snow ! 
REAL, DIMENSION(:), INTENT(IN)  :: punif_histsnow ! 
REAL, DIMENSION(:), INTENT(IN)  :: punif_agesnow ! 

REAL,DIMENSION(:,:),    INTENT(IN ), OPTIONAL :: pvegtype       ! fraction of each vegtype
REAL,DIMENSION(:,:,:),  INTENT(IN ), OPTIONAL :: pvegtype_patch ! fraction of each vegtype per patch
REAL,DIMENSION(:,:),    INTENT(IN ), OPTIONAL :: ppatch         ! fraction of each patch
LOGICAL,                INTENT(OUT), OPTIONAL :: okey
!
!
!*      0.2    declarations of local variables
!
 CHARACTER(LEN=10)                   :: ysnsurf   ! type of field
REAL,ALLOCATABLE,DIMENSION(:,:,:)   :: zw        ! total snow content
REAL,ALLOCATABLE,DIMENSION(:,:)     :: zwrho     ! total snow content from rho profile alone
REAL,ALLOCATABLE,DIMENSION(:,:,:)   :: zd        ! total snow depth
REAL,ALLOCATABLE,DIMENSION(:,:,:)   :: zdepth    ! snow depth of each layer
REAL,ALLOCATABLE,DIMENSION(:,:)     :: zdtot     ! total snow depth
REAL,DIMENSION(KL,KPATCH)           :: zpatch    ! fraction of each patch
REAL,DIMENSION(:,:),   ALLOCATABLE  :: zvegtype          ! fraction of each patch
REAL,DIMENSION(:,:,:), ALLOCATABLE  :: zvegtype_patch    ! fraction of each vegtype per patch
!
INTEGER                             :: jpatch    ! loop counter on patches
INTEGER                             :: jlayer    ! loop counter on layers
INTEGER                             :: iversion  ! surface version
 CHARACTER(LEN=16)                   :: yrecfm    ! record name
INTEGER                             :: iresp     ! error return code
LOGICAL                             :: gglacier
!
REAL(KIND=JPRB) :: zhook_handle
!---------------------------------------------------------------------------
!
IF (lhook) CALL dr_hook('PREP_HOR_SNOW_FIELDS',0,zhook_handle)
!
IF (present(ppatch)) THEN
   zpatch = ppatch
ELSE
   zpatch = 1.
ENDIF
IF (present(pvegtype)) THEN
  ALLOCATE(zvegtype(kl,SIZE(pvegtype,2)))
  zvegtype = pvegtype
ELSE
  ALLOCATE(zvegtype(kl,nvegtype))
  zvegtype = 1.
ENDIF
IF (present(pvegtype_patch)) THEN
  ALLOCATE(zvegtype_patch(kl,SIZE(pvegtype_patch,2),kpatch))
  zvegtype_patch = pvegtype_patch
ELSE
  ALLOCATE(zvegtype_patch(kl,1,kpatch))
  zvegtype_patch = 1.
ENDIF
!
!*      1.     Allocation of output field
!
 CALL allocate_gr_snow(tpsnow,kl,kpatch)
!
!---------------------------------------------------------------------------
!
!*      2.     Find if PERMSNOW must be done
!
IF(present(okey))THEN
!  
  IF ( (hfiletype=='MESONH' .OR. hfiletype=='ASCII ' .OR. hfiletype=='LFI   '.OR. hfiletype=='FA    ') &
       .AND. (hsurf=='SN_VEG ')  ) THEN
!       
     CALL open_aux_io_surf(&
                           hfile,hfiletype,'FULL  ')
     yrecfm='VERSION'
     CALL read_surf(&
                    hfiletype,yrecfm,iversion,iresp)    
     CALL close_aux_io_surf(hfile,hfiletype)
!  
     IF(iversion>7)THEN       
       CALL open_aux_io_surf(&
                           hfile,hfiletype,'NATURE')
       yrecfm='GLACIER'
       CALL read_surf(&
                    hfiletype,yrecfm,gglacier,iresp)    
       CALL close_aux_io_surf(hfile,hfiletype)  
       IF(gglacier)okey=.false.
     ENDIF
!
  ENDIF
!  
  IF(osnow_ideal)okey=.false.
!
ENDIF
!
!---------------------------------------------------------------------------
!
!*      3.     Treatment of total snow content (kg/m2)
!
ALLOCATE(zw(kl,tpsnow%NLAYER,kpatch))
!
ysnsurf='WWW'//hsurf
 CALL prep_hor_snow_field(dtco, &
                         ig, u, &
                         hprogram, hfile, hfiletype, hfilepgd, hfilepgdtype,  &
                         kluout, ounif, ysnsurf, kpatch, kteb_patch, kl, tpsnow, tptime,  &
                         punif_wsnow, punif_rsnow, punif_tsnow, punif_lwcsnow,&
                         punif_asnow, osnow_ideal, punif_sg1snow,             &
                         punif_sg2snow, punif_histsnow,punif_agesnow,         &                      
                         pf=zw, pvegtype=zvegtype,                            &
                         pvegtype_patch=zvegtype_patch, ppatch=zpatch         )  
!
!----------------------------------------------------------------------------
!
!*      4.     Treatment of total snow depth
!
ALLOCATE(zd(kl,tpsnow%NLAYER,kpatch))
!
ysnsurf='DEP'//hsurf
 CALL prep_hor_snow_field(dtco, &
                         ig, u, &
                         hprogram, hfile, hfiletype, hfilepgd, hfilepgdtype,  &
                         kluout, ounif, ysnsurf, kpatch, kteb_patch, kl, tpsnow, tptime,  &
                         punif_wsnow, punif_rsnow, punif_tsnow, punif_lwcsnow,&
                         punif_asnow, osnow_ideal, punif_sg1snow,             &
                         punif_sg2snow, punif_histsnow,punif_agesnow,         &
                         pf=zd, pvegtype=zvegtype,                            &
                         pvegtype_patch=zvegtype_patch, ppatch=zpatch         )
!
!* snow layer thickness definition
!
ALLOCATE(zdepth(kl,tpsnow%NLAYER,kpatch))
!
IF (osnow_ideal) THEN
  zdepth(:,:,:) = zd(:,:,:)
ELSE
  IF (tpsnow%NLAYER==1) THEN
    DO jpatch=1,kpatch
      zdepth(:,1,jpatch) = zd(:,1,jpatch)
    END DO
  ELSEIF (tpsnow%SCHEME=='3-L') THEN
      zdepth(:,:,:)=zd(:,:,:)
  ELSEIF (tpsnow%SCHEME=='CRO') THEN
    ALLOCATE(zdtot(kl,kpatch))
    zdtot(:,:)=0.0
    DO jlayer=1,tpsnow%NLAYER
       zdtot(:,:)=zdtot(:,:)+zd(:,jlayer,:)
    END DO
    DO jpatch=1,kpatch
       CALL snow3lgrid(zdepth(:,:,jpatch),zdtot(:,jpatch))
    END DO
    DEALLOCATE(zdtot)
  ENDIF
ENDIF
!
!----------------------------------------------------------------------------
!
!*      4.     Snow density profile
!              --------------------
!
!* density profile
ysnsurf='RHO'//hsurf
 CALL prep_hor_snow_field(dtco, &
                         ig, u, &
                         hprogram,hfile,hfiletype,hfilepgd,hfilepgdtype,           &
                         kluout,ounif,ysnsurf, kpatch, kteb_patch, kl, tpsnow, tptime,         &
                         punif_wsnow, punif_rsnow, punif_tsnow, punif_lwcsnow,     &
                         punif_asnow, osnow_ideal, punif_sg1snow,                  &
                         punif_sg2snow, punif_histsnow,punif_agesnow,              &
                         pdepth=zdepth, pvegtype=zvegtype,                         &
                         pvegtype_patch=zvegtype_patch, ppatch=zpatch              )  
!
!----------------------------------------------------------------------------
!
!*      5.     Snow water content profile
!              --------------------------
!
IF (osnow_ideal) THEN
  !
  tpsnow%WSNOW(:,:,:) = zw(:,:,:)
  !
ELSE
  !
  ALLOCATE(zwrho(SIZE(tpsnow%WSNOW,1),kpatch))
  ALLOCATE(zdtot(SIZE(tpsnow%WSNOW,1),kpatch))
  zwrho(:,:) = 0.0
  zdtot(:,:) = 0.0
  !
  !* snow depth estimated from rho profile
  DO jpatch=1,kpatch
    DO jlayer=1,tpsnow%NLAYER
      WHERE (zpatch(:,jpatch)>0. .AND. tpsnow%RHO(:,jlayer,jpatch)/=xundef)
        zwrho(:,jpatch) = zwrho(:,jpatch) + tpsnow%RHO(:,jlayer,jpatch) * zdepth(:,jlayer,jpatch)
      ELSEWHERE
        zwrho(:,jpatch) = xundef
      END WHERE
    END DO
  END DO
  !
  !* modification of snow depth: coherence between rho profile, total snow and total depth
  DO jpatch=1,kpatch
    DO jlayer=1,tpsnow%NLAYER
      WHERE(zpatch(:,jpatch)>0. .AND. zwrho(:,jpatch)/=0. .AND. zwrho(:,jpatch)/=xundef .AND. zw(:,1,jpatch)>0.0)
        zdtot(:,jpatch) = zdtot(:,jpatch) + zdepth(:,jlayer,jpatch) * zw(:,1,jpatch) / zwrho(:,jpatch)
      ENDWHERE
    END DO
    CALL snow3lgrid(zdepth(:,:,jpatch),zdtot(:,jpatch))
  END DO
  !
  !* snow content profile for each grid level
  DO jpatch=1,kpatch
    DO jlayer=1,tpsnow%NLAYER
      WHERE(zpatch(:,jpatch)>0..AND.tpsnow%RHO(:,jlayer,jpatch)/=xundef.AND.zdtot(:,jpatch)>0.)
        tpsnow%WSNOW(:,jlayer,jpatch) = tpsnow%RHO(:,jlayer,jpatch) * zdepth(:,jlayer,jpatch)
      ELSEWHERE(zpatch(:,jpatch)>0..AND.(tpsnow%RHO(:,jlayer,jpatch)==xundef.OR.zdtot(:,jpatch)==0.0))
        tpsnow%WSNOW(:,jlayer,jpatch) = 0.0
      ELSEWHERE
        tpsnow%WSNOW(:,jlayer,jpatch) = xundef
      END WHERE
    END DO
  END DO
  !
  DEALLOCATE(zwrho)
  DEALLOCATE(zdtot)
  !
ENDIF
!
!----------------------------------------------------------------------------
!
!*      6.     Albedo, snow heat content, and age
!              ----------------------------------
!
!* albedo
ysnsurf='ALB'//hsurf
 CALL prep_hor_snow_field(dtco, &
                         ig, u, &
                         hprogram,hfile,hfiletype,hfilepgd,hfilepgdtype,         &
                         kluout,ounif,ysnsurf, kpatch, kteb_patch, kl, tpsnow, tptime,        &
                         punif_wsnow, punif_rsnow, punif_tsnow, punif_lwcsnow,    &
                         punif_asnow, osnow_ideal, punif_sg1snow,                 &
                         punif_sg2snow, punif_histsnow,punif_agesnow,             &
                         pdepth=zdepth, pvegtype=zvegtype,                        &
                         pvegtype_patch=zvegtype_patch, ppatch=zpatch             ) 
!
IF (tpsnow%SCHEME/='D95') THEN
  !
  !* heat in snowpack profile
  ysnsurf='HEA'//hsurf
  CALL prep_hor_snow_field(dtco, &
                         ig, u, &
                         hprogram,hfile,hfiletype,hfilepgd,hfilepgdtype,          &
                           kluout,ounif,ysnsurf, kpatch, kteb_patch, kl, tpsnow, tptime,        &
                           punif_wsnow, punif_rsnow, punif_tsnow, punif_lwcsnow,    &
                           punif_asnow, osnow_ideal, punif_sg1snow,                 &
                           punif_sg2snow, punif_histsnow,punif_agesnow,             &
                           pdepth=zdepth, pvegtype=zvegtype,                        &
                           pvegtype_patch=zvegtype_patch, ppatch=zpatch             )
  !
ENDIF
!
IF (tpsnow%SCHEME=='CRO'.OR. tpsnow%SCHEME=='3-L') THEN
  !
  !* age in snowpack profile
  ysnsurf='AGE'//hsurf
  CALL prep_hor_snow_field(dtco, &
                         ig, u, &
                         hprogram,hfile,hfiletype,hfilepgd,hfilepgdtype,        &
                         kluout,ounif,ysnsurf, kpatch, kteb_patch, kl, tpsnow, tptime,        &
                         punif_wsnow, punif_rsnow, punif_tsnow, punif_lwcsnow,    &
                         punif_asnow, osnow_ideal, punif_sg1snow,                 &
                         punif_sg2snow, punif_histsnow,punif_agesnow,             &
                         pdepth=zdepth, pvegtype=zvegtype,                        &
                         pvegtype_patch=zvegtype_patch, ppatch=zpatch             )   
  !
  WHERE(tpsnow%WSNOW(:,1,:)>0.0.AND.tpsnow%WSNOW(:,1,:)/=xundef.AND. &
        tpsnow%AGE(:,1,:)==0.0.AND.tpsnow%ALB(:,:)<xaglamin)
        tpsnow%ALB(:,:)=(xaglamin+xaglamax)/2.0
  ENDWHERE
  !
ENDIF
!
!----------------------------------------------------------------------------
!
!*      7.     Crocus specific parameters
!              --------------------------
!
IF (tpsnow%SCHEME=='CRO') THEN
  !
  ysnsurf='SG1'//hsurf
  CALL prep_hor_snow_field(dtco, &
                         ig, u, &
                         hprogram,hfile,hfiletype,hfilepgd,hfilepgdtype,        &
                         kluout,ounif,ysnsurf, kpatch, kteb_patch, kl, tpsnow, tptime,        &
                         punif_wsnow, punif_rsnow, punif_tsnow, punif_lwcsnow,    &
                         punif_asnow, osnow_ideal, punif_sg1snow,                 &
                         punif_sg2snow, punif_histsnow,punif_agesnow,             &
                         pdepth=zdepth, pvegtype=zvegtype,                        &
                         pvegtype_patch=zvegtype_patch, ppatch=zpatch             )   
  !
  ysnsurf='SG2'//hsurf
  CALL prep_hor_snow_field(dtco, &
                         ig, u, &
                         hprogram,hfile,hfiletype,hfilepgd,hfilepgdtype,        &
                         kluout,ounif,ysnsurf, kpatch, kteb_patch, kl, tpsnow, tptime,        &
                         punif_wsnow, punif_rsnow, punif_tsnow, punif_lwcsnow,    &
                         punif_asnow, osnow_ideal, punif_sg1snow,                 &
                         punif_sg2snow, punif_histsnow,punif_agesnow,             &
                         pdepth=zdepth, pvegtype=zvegtype,                        &
                         pvegtype_patch=zvegtype_patch, ppatch=zpatch             )   
  !
  ysnsurf='HIS'//hsurf
  CALL prep_hor_snow_field(dtco, &
                         ig, u, &
                         hprogram,hfile,hfiletype,hfilepgd,hfilepgdtype,        &
                         kluout,ounif,ysnsurf, kpatch, kteb_patch, kl, tpsnow, tptime,        &
                         punif_wsnow, punif_rsnow, punif_tsnow, punif_lwcsnow,    &
                         punif_asnow, osnow_ideal, punif_sg1snow,                 &
                         punif_sg2snow, punif_histsnow,punif_agesnow,             &
                         pdepth=zdepth, pvegtype=zvegtype,                        &
                         pvegtype_patch=zvegtype_patch, ppatch=zpatch             )   
  !
ENDIF
!
!*      8.     Deallocations
!
DEALLOCATE(zd      )
DEALLOCATE(zw      )
DEALLOCATE(zdepth  )
DEALLOCATE(zvegtype)
DEALLOCATE(zvegtype_patch)
!
IF (lhook) CALL dr_hook('PREP_HOR_SNOW_FIELDS',1,zhook_handle)
!
!----------------------------------------------------------------------------
!
END SUBROUTINE prep_hor_snow_fields