init_dst.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 init_dst (DSTK, U, &
                     HPROGRAM,  &    ! Program calling unit
                  KSIZE_P, & ! Number of nature points in a patch
                  KR_P, &    ! Mask from patch --> nature vectors
                  PVEGTYPE_PATCH  ) ! fraction (in a nature point) of a vegtype for a patch

!
USE modd_dst_n, ONLY : dst_t
USE modd_surf_atm_n, ONLY : surf_atm_t
!
USE modd_dst_surf
USE modd_data_cover_par, ONLY : nvt_no, nvt_rock
!
USE modi_get_luout
USE modi_get_vegtype_2_patch_mask
USE modi_abor1_sfx
!
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb
!
IMPLICIT NONE
!
!PASSED VARIABLES
!
TYPE(dst_t), INTENT(INOUT) :: DSTK
TYPE(surf_atm_t), INTENT(INOUT) :: U
!
 CHARACTER(LEN=6), INTENT(IN)      :: HPROGRAM              !Passing unit
!
INTEGER, INTENT(IN) :: KSIZE_P
INTEGER, DIMENSION(:), INTENT(IN) :: KR_P
REAL, DIMENSION(:,:), INTENT(IN) :: PVEGTYPE_PATCH
!
!LOCAL VARIABLES
 CHARACTER(LEN=4)    :: CRGUNIT               ! type of log-normal geometric mean radius
INTEGER             :: JVEG                  ! Counter for vegetation classes
INTEGER             :: JVEG_IN               ! Vegetation index
INTEGER             :: JPATCH                ! Counter for patches
INTEGER             :: JMODE                 ! Counter for dust modes
INTEGER             :: JMODE_IDX             ! Index for dust modes
INTEGER             :: ILUOUT
INTEGER             :: ISIZE_LARGEST_DST
REAL(KIND=JPRB) :: ZHOOK_HANDLE

!get output listing unit
IF (lhook) CALL dr_hook('INIT_DST',0,zhook_handle)
 CALL get_luout(hprogram,iluout)
!
!Allocate memory for the real values which will be used by the model
ALLOCATE(dstk%XEMISRADIUS_DST(ndstmde))
ALLOCATE(dstk%XEMISSIG_DST   (ndstmde))
ALLOCATE(dstk%XMSS_FRC_SRC   (ndstmde))
!
!Get initial size distributions. This is cut and pasted
!from dead routine dstpsd.F90
!Check for different source parameterizations
IF(cemisparam_dst.EQ."She84")THEN
  crgunit   = 'MASS'
  xemisradius_ini_dst(:) = 0.5d6 * (/ 0.0111e-6, 2.524e-6, 42.10e-6 /)  ! [um]  Mass median radius She84 p. 75 Table 1
  xemissig_ini_dst(:) = (/ 1.89 , 2.0  , 2.13 /)             ! [frc] Geometric standard deviation She84 p. 75 Table 1
  xmss_frc_src_ini(:) = (/2.6e-6, 0.781, 0.219/)             ! [frc] Mass fraction She84 p. 75 Table 1
ELSEIF(cemisparam_dst.EQ."PaG77")THEN
  crgunit   = 'MASS'
  xemisradius_ini_dst(:) = 0.5d6 * (/0.27e-6  ,  5.6e-6  ,  57.6e-6 /) ! [um] Mass median radius PaG77 p. 2080 Table 1 
  xemissig_ini_dst(:) = (/ 1.88, 2.2  , 1.62 /)              ! [frc] Geometric standard deviation PaG77 p. 2080 Table 1
  xmss_frc_src_ini(:) = (/0.036, 0.957, 0.007/)              ! [frc] Mass fraction BSM96 p. 73 Table 2 (ad hoc)  
ELSEIF(cemisparam_dst.EQ."Dal87") THEN 
  ! D'Almeida, 1987 as default
  crgunit   = 'MASS'
  xemisradius_ini_dst(:) = 0.5d6 * (/ 0.832e-6 ,  4.82e-6 , 19.38e-6 /) ! [um] Mass median radius BSM96 p. 73 Table 2
  xemissig_ini_dst(:) = (/ 2.10, 1.90 , 1.60 /)              ! [frc] Geometric standard deviation BSM96 p. 73 Table 2
  xmss_frc_src_ini(:) = (/0.036, 0.957, 0.007/)              ! [frc] Mass fraction BSM96 p. 73 Table 2  
ELSEIF (cemisparam_dst.EQ."alf98".OR.cemisparam_dst.eq."EXPLI") THEN  !  Alfaro et al 1998 as default
  IF (cemisparam_dst.EQ."alf98") xflx_mss_fdg_fct = 6e-4
  IF (cemisparam_dst.EQ."EXPLI") xflx_mss_fdg_fct = 3.5e-4
  crgunit   = 'MASS'
  xemisradius_ini_dst(:) = 0.5*(/ 1.5, 6.7, 14.2 /)      ! [um] Mass median radius BSM96 p. 73 Table 2
  xemissig_ini_dst(:) = (/1.70, 1.60, 1.50/)          ! [frc] Geometric standard deviation BSM96 p. 73 Table 2
  xmss_frc_src_ini(:) = (/0.01, 0.19, 0.8 /)          ! [frc] Mass fraction BSM96 p. 73 Table 2  
ELSEIF (cemisparam_dst.EQ."AMMA ") THEN ! Default: New distribution from AMMA
  xflx_mss_fdg_fct = 105.e-4
  crgunit   = 'NUMB' 
  xemisradius_ini_dst(:) = 0.5*(/ 0.078, 0.641, 5.00 /)  ! [um] Number median radius 
  xemissig_ini_dst(:) = (/ 1.75,  1.76, 1.70/)        ! [frc] Geometric standard deviation 
  xmss_frc_src_ini(:) = (/0.008, 0.092, 0.99/)        ! [frc] Mass fraction  
ELSEIF (cemisparam_dst.EQ."CRUM ") THEN ! Default: New distribution from AMMA
  xflx_mss_fdg_fct = 10.e-4
  crgunit   = 'NUMB' 
  xemisradius_ini_dst(:) = 0.5*(/ 0.078, 0.641, 5.00 /)  ! [um] Number median radius 
  xemissig_ini_dst(:) = (/  1.75,   1.76, 1.70   /)   ! [frc] Geometric standard deviation   
  xmss_frc_src_ini(:) = (/0.0005, 0.0029, 0.9966/)    ! [frc] Mass fraction  
ELSE
  WRITE(iluout,*) " FATAL ERROR "
  WRITE(iluout,*) " YOU MUST DECIDE THE EMISSIUON PARAMETERIZATION, YOU USES "
  WRITE(iluout,*) " CEMISPARAM_DST = ",cemisparam_dst," AND IT IS NOT DEFINED "
  WRITE(iluout,*) " see init_dstn.f90 to see what dust parameterization is available. "
  CALL abor1_sfx("INIT_DST: CEMISPARAM_DST not defined ")
ENDIF
!
DO jmode=1,ndstmde
  jmode_idx=jorder_dst(jmode)
  !
  dstk%XEMISSIG_DST   (jmode) = xemissig_ini_dst(jmode_idx)
  dstk%XEMISRADIUS_DST(jmode) = xemisradius_ini_dst(jmode_idx)
  dstk%XMSS_FRC_SRC   (jmode) = xmss_frc_src_ini(jmode_idx)
  !
  !Get emisradius, and at the same time convert to number median radius
  IF (crgunit=='MASS') &
    dstk%XEMISRADIUS_DST(jmode) = dstk%XEMISRADIUS_DST(jmode) * exp(-3.d0 * (log(dstk%XEMISSIG_DST(jmode)))**2)  
  !
ENDDO
!
!Normalize the sum of the emissions to 1 so that all dust is
!put in one mode or the other
IF(sum(dstk%XMSS_FRC_SRC(:)).LT.1.) dstk%XMSS_FRC_SRC(:) = dstk%XMSS_FRC_SRC(:) / sum(dstk%XMSS_FRC_SRC(:))
!
!Allocate memory
!ALLOCATE(NVEGNO_DST)
!Set the number of classes that can emit dust (fxm: set this elsewhere)
nvegno_dst = 2
!
!Allocate memory for the vegtype-translator
ALLOCATE(dstk%NVT_DST(nvegno_dst))
!
!Set the dust/vegtype translator vector
dstk%NVT_DST(1)  = nvt_no
dstk%NVT_DST(2)  = nvt_rock
!
!Allocate memory for roughness lengths of erodible surfaces
ALLOCATE(dstk%Z0_EROD_DST(nvegno_dst))
!
!Set the roughness lengths corresponding to erodible surfaces
!Smooth roughness length is given to 1.d-5 (dstmbl.f90)
dstk%Z0_EROD_DST(1) = 30.d-6    !m (30 um) 
dstk%Z0_EROD_DST(2) = 200.d-6   !m (200 um) 
!
!Allocate memory for dust emitter surface vectors in patch vectors
IF (.NOT.ASSOCIATED(dstk%NSIZE_PATCH_DST)) ALLOCATE(dstk%NSIZE_PATCH_DST(nvegno_dst))
!
DO jveg = 1,nvegno_dst
  !Count all the points in the patch where you have dust emitter vegetation
  dstk%NSIZE_PATCH_DST(jveg) = count(pvegtype_patch(:,dstk%NVT_DST(jveg)) > 0.) 
ENDDO
!
!Find the largest dust emitter vector in any patch
!ALLOCATE (NSIZE_LARGEST_DST)
isize_largest_dst = 0
DO jveg = 1,nvegno_dst
  isize_largest_dst = max(isize_largest_dst,dstk%NSIZE_PATCH_DST(jveg))
ENDDO
!
!Allocate memory for KR_PATCH_DST mask translate from patch vector to dust vector
ALLOCATE(dstk%NR_PATCH_DST(isize_largest_dst,nvegno_dst))
!
!Initialize the mask array
dstk%NR_PATCH_DST(:,:)=0
!
!Get values from the dust emitter vegetation mask
DO jveg=1,nvegno_dst
  jveg_in = dstk%NVT_DST(jveg)          ! Get the real vegtype index
#ifdef RJ_OFIX
  CALL get_vegtype_2_patch_mask(iluout,    &
           dstk%NSIZE_PATCH_DST(jveg),             &!I Size of dust emitter vector
           ksize_p,                   &!I Size of patch vector
!RJ: attempt to make this call generic
           kr_p,&!I Mask from patch to nature
           pvegtype_patch,                           &!I Fraction of vegtype of nature point within jpatch 
           dstk%NR_PATCH_DST(:dstk%NSIZE_PATCH_DST(jveg),jveg),  &!O Part of mask array to fill with values
           jveg_in                                  &!I Index of vegtype in question
             )  
#else
  CALL get_vegtype_2_patch_mask(iluout,    &
           dstk%NSIZE_PATCH_DST(jveg),             &!I Size of dust emitter vector
           ksize_p,                   &!I Size of patch vector
           kr_p,                              &!I Mask from patch to nature
           pvegtype_patch,                           &!I Fraction of vegtype of nature point within jpatch 
           dstk%NR_PATCH_DST(:dstk%NSIZE_PATCH_DST(jveg),jveg),  &!O Part of mask array to fill with values
           jveg_in                                  &!I Index of vegtype in question
             )  
#endif
ENDDO    !Loop on veg-types
!
IF (lhook) CALL dr_hook('INIT_DST',1,zhook_handle)
!
END SUBROUTINE init_dst