prep_hor_teb_greenroof_field.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_teb_greenroof_field (DTCO, IG, I, UG, U, USS, TGR, TGRO, TGRPE, TGRP, &
                                         tg, top, &
                                         hprogram,hsurf,hatmfile,hatmfiletype,hpgdfile,hpgdfiletype,kpatch)
!     #################################################################################################
!
!!****  *PREP_HOR_TEB_GREENROOF_FIELD* - reads, interpolates and prepares an ISBA field for green roofs
!!
!!    PURPOSE
!!    -------
!
!!**  METHOD
!!    ------
!!     
!!    Based on "prep_hor_teb_garden_field"
!!
!!    REFERENCE
!!    ---------
!!      
!!
!!    AUTHOR
!!    ------
!!     V. Masson, A. Lemonsu & C. de Munck 
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    07/2011
!!------------------------------------------------------------------
!
!
USE modd_data_cover_n, ONLY : data_cover_t
USE modd_isba_grid_n, ONLY : isba_grid_t
USE modd_isba_n, ONLY : isba_t
USE modd_surf_atm_grid_n, ONLY : surf_atm_grid_t
USE modd_surf_atm_n, ONLY : surf_atm_t
USE modd_surf_atm_sso_n, ONLY : surf_atm_sso_t
USE modd_teb_greenroof_n, ONLY : teb_greenroof_t
USE modd_teb_greenroof_option_n, ONLY : teb_greenroof_options_t
USE modd_teb_greenroof_pgd_evol_n, ONLY : teb_greenroof_pgd_evol_t
USE modd_teb_greenroof_pgd_n, ONLY : teb_greenroof_pgd_t
USE modd_teb_grid_n, ONLY : teb_grid_t
USE modd_teb_option_n, ONLY : teb_options_t
!
USE modd_prep,               ONLY : cingrid_type, cinterp_type, xzs_ls,          &
                                    xlat_out, xlon_out, xx_out, xy_out,          &
                                    linterp, cmask
USE modd_prep_teb_greenroof, ONLY : xgrid_soil, ngrid_level,                     &
                                    xwsnow_gr, xrsnow_gr, xtsnow_gr,xlwcsnow_gr, &
                                    xagesnow_gr, xasnow_gr, lsnow_ideal_gr
USE modd_isba_par,           ONLY : xwgmin
USE modd_data_cover_par,     ONLY : nvegtype
USE modd_surf_par,           ONLY : xundef
!
USE modi_read_prep_teb_greenroof_conf
USE modi_read_prep_greenroof_snow
USE modi_prep_teb_greenroof_ascllv
USE modi_prep_teb_greenroof_grib
USE modi_prep_teb_greenroof_unif
USE modi_prep_teb_greenroof_buffer
USE modi_hor_interpol
USE modi_vegtype_grid_to_patch_grid
USE modi_prep_hor_snow_fields
USE modi_get_luout
USE modi_prep_teb_greenroof_extern
USE modi_put_on_all_vegtypes
!
USE modi_abor1_sfx
!
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(isba_t), INTENT(INOUT) :: i
TYPE(surf_atm_grid_t), INTENT(INOUT) :: ug
TYPE(surf_atm_t), INTENT(INOUT) :: u
TYPE(surf_atm_sso_t), INTENT(INOUT) :: uss
TYPE(teb_greenroof_t), INTENT(INOUT) :: tgr
TYPE(teb_greenroof_options_t), INTENT(INOUT) :: tgro
TYPE(teb_greenroof_pgd_evol_t), INTENT(INOUT) :: tgrpe
TYPE(teb_greenroof_pgd_t), INTENT(INOUT) :: tgrp
TYPE(teb_grid_t), INTENT(INOUT) :: tg
TYPE(teb_options_t), INTENT(INOUT) :: top
!
 CHARACTER(LEN=6),   INTENT(IN)  :: hprogram  ! program calling surf. schemes
 CHARACTER(LEN=7),   INTENT(IN)  :: hsurf     ! type of field
 CHARACTER(LEN=28),  INTENT(IN)  :: hatmfile    ! name of the Atmospheric file
 CHARACTER(LEN=6),   INTENT(IN)  :: hatmfiletype! type of the Atmospheric file
 CHARACTER(LEN=28),  INTENT(IN)  :: hpgdfile    ! name of the Atmospheric file
 CHARACTER(LEN=6),   INTENT(IN)  :: hpgdfiletype! type of the Atmospheric file
!
INTEGER,            INTENT(IN)  :: kpatch
!
!*      0.2    declarations of local variables
!
 CHARACTER(LEN=6)              :: yfiletype      ! type of input file
 CHARACTER(LEN=28)             :: yfile          ! name of file
 CHARACTER(LEN=6)              :: yfilepgdtype ! type of input file
 CHARACTER(LEN=28)             :: yfilepgd     ! name of file
 CHARACTER(LEN=6)              :: yfiletype_snow ! type of input file
 CHARACTER(LEN=28)             :: yfile_snow     ! name of file
 CHARACTER(LEN=6)              :: yfilepgdtype_snow ! type of input file
 CHARACTER(LEN=28)             :: yfilepgd_snow     ! name of file  
REAL, POINTER,     DIMENSION(:,:,:) :: zfieldin       ! field to interpolate horizontally
REAL, POINTER,     DIMENSION(:,:)   :: zfield         ! field to interpolate horizontally
REAL, ALLOCATABLE, DIMENSION(:,:,:) :: zfieldoutp ! field interpolated   horizontally
REAL, ALLOCATABLE, DIMENSION(:,:,:) :: zfieldoutv !
REAL, ALLOCATABLE, DIMENSION(:,:,:) :: zvegtype_patch ! vegtype for each patch
REAL, ALLOCATABLE, DIMENSION(:,:)   :: zw             ! work array (x, fine   soil grid, npatch)
REAL, ALLOCATABLE, DIMENSION(:,:)   :: zf             ! work array (x, output soil grid, npatch)
REAL, ALLOCATABLE, DIMENSION(:,:)   :: zdg            ! out T grid (x, output soil grid, npatch)
REAL, ALLOCATABLE, DIMENSION(:,:)   :: zpatch         ! work array for patches
REAL, ALLOCATABLE, DIMENSION(:)     :: zsg1snow, zsg2snow, zhistsnow
INTEGER                             :: iluout         ! output listing logical unit
!
LOGICAL                             :: gunif          ! flag for prescribed uniform field
LOGICAL                             :: gunif_snow     ! flag for prescribed uniform field
INTEGER                             :: jvegtype, jpatch    ! loop on vegtypes
INTEGER                             :: jlayer         ! loop on layers
INTEGER                             :: ini, inl, inp
INTEGER                             :: iwork          ! Work integer
REAL(KIND=JPRB) :: zhook_handle
!-------------------------------------------------------------------------------------
!
!
!*      1.     Reading of input file name and type
!
IF (lhook) CALL dr_hook('PREP_HOR_TEB_GREENROOF_FIELD',0,zhook_handle)
!
!
 CALL get_luout(hprogram,iluout)
!
 CALL read_prep_teb_greenroof_conf(hprogram,hsurf,yfile,yfiletype,yfilepgd,yfilepgdtype,&
                                  hatmfile,hatmfiletype,hpgdfile,hpgdfiletype,iluout,gunif)
!
 cmask = 'TOWN  '
!
ini=SIZE(tg%XLAT)
!
!-------------------------------------------------------------------------------------
!*      2.     Snow variables case
!
IF (hsurf=='SN_VEG ') THEN
  CALL read_prep_greenroof_snow(hprogram,tgr%CUR%TSNOW%SCHEME,tgr%CUR%TSNOW%NLAYER,yfile_snow,&
        yfiletype_snow,yfilepgd_snow,yfilepgdtype_snow,gunif_snow)
  IF(.NOT.gunif_snow.AND.len_trim(yfile_snow)==0.AND.len_trim(yfiletype_snow)==0)THEN
    !IF(LEN_TRIM(YFILE)/=0.AND.LEN_TRIM(YFILETYPE)/=0)THEN
    IF (yfiletype=='GRIB') THEN
      yfile_snow    =yfile
      yfiletype_snow=yfiletype
      yfilepgd_snow    =yfilepgd
      yfilepgdtype_snow=yfilepgdtype       
    ELSE            
      gunif_snow=.true.
      IF(all(xwsnow_gr==xundef))xwsnow_gr=0.0  
    ENDIF
  ENDIF          
  ALLOCATE(zsg1snow(SIZE(xwsnow_gr)))
  ALLOCATE(zsg2snow(SIZE(xwsnow_gr)))
  ALLOCATE(zhistsnow(SIZE(xwsnow_gr)))
  ALLOCATE(zpatch(SIZE(tgrp%XVEGTYPE,1),1))
  ALLOCATE(zvegtype_patch(SIZE(tgrp%XVEGTYPE,1),SIZE(tgrp%XVEGTYPE,2),1))
  !
  zpatch=1.
  zvegtype_patch(:,:,1) = tgrp%XVEGTYPE(:,:)
  CALL prep_hor_snow_fields(dtco, &
                            ig, u, &
                            hprogram,hsurf,                 &
                            yfile,yfiletype,                &
                            yfilepgd, yfilepgdtype,         &
                            iluout,gunif_snow, 1, kpatch,   &
                            SIZE(tg%XLAT),tgr%CUR%TSNOW, top%TTIME,        &
                            xwsnow_gr, xrsnow_gr, xtsnow_gr,&
                            xlwcsnow_gr, xasnow_gr,         &
                            lsnow_ideal_gr, zsg1snow,       &
                            zsg2snow, zhistsnow, xagesnow_gr,  &
                            tgrp%XVEGTYPE,zvegtype_patch, zpatch )
  DEALLOCATE(zsg1snow)
  DEALLOCATE(zsg2snow)
  DEALLOCATE(zhistsnow)
  DEALLOCATE(zpatch)
  DEALLOCATE(zvegtype_patch)
  IF (lhook) CALL dr_hook('PREP_HOR_TEB_GREENROOF_FIELD',1,zhook_handle)
  RETURN
END IF
!
!-------------------------------------------------------------------------------------
!
!*      3.     Reading of input  configuration (Grid and interpolation type)
!
IF (gunif) THEN
  CALL prep_teb_greenroof_unif(iluout,hsurf,zfieldin)
ELSE IF (yfiletype=='ASCLLV') THEN
  CALL prep_teb_greenroof_ascllv(dtco, ug, u, uss, &
                                 hprogram,hsurf,iluout,zfieldin)
ELSE IF (yfiletype=='GRIB  ') THEN
  CALL prep_teb_greenroof_grib(hprogram,hsurf,yfile,iluout,zfieldin)
ELSE IF (yfiletype=='MESONH' .OR. yfiletype=='ASCII ' .OR. yfiletype=='LFI   '.OR. yfiletype=='FA    ') THEN
   CALL prep_teb_greenroof_extern(dtco, i, u, &
                                  hprogram,hsurf,yfile,yfiletype,yfilepgd,yfilepgdtype,iluout,kpatch,zfieldin)
ELSE IF (yfiletype=='BUFFER') THEN
   CALL prep_teb_greenroof_buffer(hprogram,hsurf,iluout,zfieldin)
ELSE
   CALL abor1_sfx('PREP_HOR_TEB_GREENROOF_FIELD: data file type not supported : '//yfiletype)
END IF
!
!-------------------------------------------------------------------------------------
!
!*      5.     Horizontal interpolation for vegtype grid
!
inl = SIZE(zfieldin,2)
inp = SIZE(zfieldin,3)
!
ALLOCATE(zfieldoutp(ini,inl,inp))
ALLOCATE(zfield(SIZE(zfieldin,1),inl))
!
DO jpatch = 1, inp
  zfield=zfieldin(:,:,jpatch)
  IF (inp==nvegtype) linterp = (tgrp%XVEGTYPE(:,jpatch) > 0.)
  CALL hor_interpol(dtco, u, &
                    iluout,zfield,zfieldoutp(:,:,jpatch))
  linterp = .true.
END DO
!
DEALLOCATE(zfield)
!
ALLOCATE(zfieldoutv(ini,inl,nvegtype))
!
 CALL put_on_all_vegtypes(ini,inl,inp,nvegtype,zfieldoutp,zfieldoutv)
!
DEALLOCATE(zfieldoutp)
!
!-------------------------------------------------------------------------------------
!
!*      6.     Transformation from vegtype grid to patch grid
!
ALLOCATE(zw(ini,inl))
zw = 0.
DO jvegtype=1,nvegtype
  DO jlayer=1,SIZE(zw,2)
    zw(:,jlayer) = zw(:,jlayer) + tgrp%XVEGTYPE(:,jvegtype) * zfieldoutv(:,jlayer,jvegtype)
  END DO
END DO
!
!
!*      8.     Deallocations
!
DEALLOCATE(zfieldin )
DEALLOCATE(zfieldoutv) 
!
!-------------------------------------------------------------------------------------
!
!*      7.     Return to historical variable
!
!
SELECT CASE (hsurf)
  !
  !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  !
 CASE('WG     ') 
  ALLOCATE(zf(ini,tgro%NLAYER_GR))
  !
  !* interpolates on output levels
  CALL init_from_ref_grid(xgrid_soil,zw,tgrp%XDG,zf)
  !
  !* retrieves soil water content from soil relative humidity
  ALLOCATE(tgr%CUR%XWG(ini,tgro%NLAYER_GR))
  tgr%CUR%XWG(:,:) = tgrp%XWWILT + zf(:,:) * (tgrp%XWFC-tgrp%XWWILT)
  tgr%CUR%XWG(:,:) = max(min(tgr%CUR%XWG(:,:),tgrp%XWSAT),xwgmin)
  !
  WHERE(zf(:,:)==xundef)tgr%CUR%XWG(:,:)=xundef
  !
  DEALLOCATE(zf)
  !
  !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  !
 CASE('WGI    ')
  ALLOCATE(zf(ini,tgro%NLAYER_GR))
  !
  !* interpolates on output levels
  CALL init_from_ref_grid(xgrid_soil,zw,tgrp%XDG,zf)
  !
  !* retrieves soil ice content from soil relative humidity
  ALLOCATE(tgr%CUR%XWGI(ini,tgro%NLAYER_GR))
  tgr%CUR%XWGI(:,:) = zf(:,:) * tgrp%XWSAT
  tgr%CUR%XWGI(:,:) = max(min(tgr%CUR%XWGI(:,:),tgrp%XWSAT),0.)
  !
  WHERE(zf(:,:)==xundef)tgr%CUR%XWGI(:,:)=xundef
  !
  DEALLOCATE(zf)
  !
  !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  !
 CASE('TG     ') 
  iwork=tgro%NLAYER_GR
  ALLOCATE(tgr%CUR%XTG(ini,iwork))
  ALLOCATE(zdg(SIZE(tgrp%XDG,1),iwork))
  !* diffusion method, the soil grid is the same as for humidity
  zdg(:,:) = tgrp%XDG(:,:)
  CALL init_from_ref_grid(xgrid_soil,zw,zdg,tgr%CUR%XTG)
  DEALLOCATE(zdg)
  !
!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  !
 CASE('WR     ') 
  ALLOCATE(tgr%CUR%XWR(ini))
  tgr%CUR%XWR(:) = zw(:,1)
  !
  !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  !
 CASE('LAI    ') 
  !* LAI is updated only if present and pertinent (evolutive LAI) in input file
  WHERE (zw(:,1)/=xundef) tgrpe%CUR%XLAI(:) = zw(:,1)
  !
END SELECT
!
DEALLOCATE(zw)
!-------------------------------------------------------------------------------------
!
IF (lhook) CALL dr_hook('PREP_HOR_TEB_GREENROOF_FIELD',1,zhook_handle)
!
!-------------------------------------------------------------------------------------
!-------------------------------------------------------------------------------------
!
 CONTAINS
!
!-------------------------------------------------------------------------------------
!-------------------------------------------------------------------------------------
!
SUBROUTINE init_from_ref_grid(PGRID1,PT1,PD2,PT2)
!
USE modi_interp_grid_nat
!
REAL, DIMENSION(:,:), INTENT(IN)  :: pt1    ! variable profile
REAL, DIMENSION(:),   INTENT(IN)  :: pgrid1 ! normalized grid
REAL, DIMENSION(:,:), INTENT(IN)  :: pd2    ! output layer thickness
REAL, DIMENSION(:,:), INTENT(OUT) :: pt2    ! variable profile
!
INTEGER                                  :: ji, jl  ! loop counter
REAL, DIMENSION(SIZE(PT1,1),SIZE(PT1,2)) :: zd1 ! input grid
!
INTEGER :: ilayer1, ilayer2
REAL(KIND=JPRB) :: zhook_handle
!
!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
!
IF (lhook) CALL dr_hook('INIT_FROM_REF_GRID',0,zhook_handle)
!
IF (SIZE(pt1,2)==3) THEN
!
!* 1. case with only 3 input levels (typically coming from 'UNIF')
!     -----------------------------
!
       !surface layer (generally 0.01m imposed)
       pt2(:,1) = pt1(:,1) 
       !deep layers
       DO jl=2,tgro%NLAYER_GR
          pt2(:,jl) = pt1(:,3)
       END DO
       !if root layers
       DO ji=1,SIZE(pt1,1)
          DO jl=2,tgro%NLAYER_GR
             IF(tgrp%XROOTFRAC(ji,jl)<=1.0)THEN 
                pt2(ji,jl) = pt1(ji,2)
                EXIT
             ENDIF
          END DO
       END DO 
!
!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
ELSE
!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
!
!* 2. case with fine grid as input (general case)
!     ----------------------------
!
  DO jl=1,SIZE(pt1,2)
    zd1(:,jl) = pgrid1(jl)
  END DO
!
  CALL interp_grid_nat(zd1,pt1(:,:),pd2,pt2(:,:))
!
!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
END IF
!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
!  
IF (lhook) CALL dr_hook('INIT_FROM_REF_GRID',1,zhook_handle)
!
!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
END SUBROUTINE init_from_ref_grid
!-------------------------------------------------------------------------------------
!
END SUBROUTINE prep_hor_teb_greenroof_field