get_teb_depths.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 get_teb_depths (DTCO,  &
                                 HFILE, HFILETYPE, HFILEPGD, HFILEPGDTYPE, &
                                 PD_ROOF, PD_ROAD, PD_WALL, PD_FLOOR,HDIR)
!     ##############################################################
!
!!**** *CONVERT_COVER* 
!!
!!    PURPOSE
!!    -------
!!
!!    METHOD
!!    ------
!!   
!!    EXTERNAL
!!    --------
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!
!!    REFERENCE
!!    ---------
!!
!!    AUTHOR
!!    ------
!!
!!    V. Masson        Meteo-France
!!
!!    MODIFICATION
!!    ------------
!!
!!    Original    01/2004
!     
!----------------------------------------------------------------------------
!
!*    0.     DECLARATION
!            -----------
!
USE modd_data_cover_n, ONLY : data_cover_t, data_cover_init
!
USE modd_surfex_mpi, ONLY : nrank, npio
USE modd_data_cover,     ONLY : xdata_d_roof, xdata_d_road, xdata_d_wall, xdata_d_floor
USE modd_data_cover_par, ONLY : ncover, ntype, ndata_roof_layer, ndata_road_layer, &
                                ndata_wall_layer, ndata_floor_layer
!
USE modi_open_aux_io_surf
USE modi_close_aux_io_surf
USE mode_read_surf_cov, ONLY : read_surf_cov
!
USE modi_read_surf
USE modi_av_pgd
USE modi_old_name
USE modi_thermal_layers_conf
USE modi_open_aux_io_surf
USE modi_close_aux_io_surf
USE modi_read_lecoclimap
USE modi_default_data_cover
!
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb
!
USE modi_abor1_sfx
!
IMPLICIT NONE
!
!*    0.1    Declaration of arguments
!            ------------------------
!
!
!
!
TYPE(data_cover_t), INTENT(INOUT) :: DTCO
!
 CHARACTER(LEN=28),  INTENT(IN)  :: HFILE ! type of input file
 CHARACTER(LEN=6),   INTENT(IN)  :: HFILETYPE ! type of input file
 CHARACTER(LEN=28),  INTENT(IN)  :: HFILEPGD  ! type of input file
 CHARACTER(LEN=6),   INTENT(IN)  :: HFILEPGDTYPE ! type of input file
!
REAL, DIMENSION(:,:), INTENT(OUT), OPTIONAL   :: PD_ROOF
REAL, DIMENSION(:,:), INTENT(OUT), OPTIONAL   :: PD_ROAD
REAL, DIMENSION(:,:), INTENT(OUT), OPTIONAL   :: PD_WALL
REAL, DIMENSION(:,:), INTENT(OUT), OPTIONAL   :: PD_FLOOR
!
 CHARACTER(LEN=1), INTENT(IN), OPTIONAL :: HDIR
!
!*    0.2    Declaration of local variables
!            ------------------------------
!
TYPE(data_cover_t) :: YDTCO
!
LOGICAL, DIMENSION(:), ALLOCATABLE   :: GCOVER ! flag to read the covers
REAL,    DIMENSION(:,:), ALLOCATABLE :: ZCOVER ! cover fractions
REAL,    DIMENSION(:,:), ALLOCATABLE :: ZD     ! depth of surface layers
REAL,    DIMENSION(:,:), ALLOCATABLE :: ZPAR_D ! depth of data_surface layers
!
REAL, DIMENSION(:,:), ALLOCATABLE :: ZDATA
!
INTEGER           :: IVERSION_PGD, IVERSION_PREP       ! surface version
INTEGER           :: IBUGFIX_PGD, IBUGFIX_PREP        ! surface bugfix version
INTEGER           :: IVERSION       ! surface version
INTEGER           :: IBUGFIX        ! surface bugfix version
 CHARACTER(LEN=1) :: YDIR
 CHARACTER(LEN=3)  :: YAREA          ! Area where field is to be averaged
 CHARACTER(LEN=5)  :: YSURF          ! Type of surface
 CHARACTER(LEN=12) :: YRECFM         ! Name of the article to be read
 CHARACTER(LEN=12) :: YRECFM0        ! Name of the article to be read
 CHARACTER(LEN=12) :: YRECFM1        ! Name of the article to be read
 CHARACTER(LEN=12) :: YRECFM2        ! Name of the article to be read
 CHARACTER(LEN=12) :: YRECFM3        ! Name of the article to be read 
INTEGER :: IRESP          ! reading return code
INTEGER :: ILAYER                   ! number of surface layers
INTEGER :: JL                   ! loop counter on surface layers
INTEGER :: IPAR_LAYER               ! number of data surface layers
INTEGER :: IDATA_LAYER              ! number of data surface layers from ecoclimap
INTEGER :: ILU                      ! number of points
LOGICAL           :: GECOCLIMAP, GECOSG
LOGICAL           :: GDATA, GDIM          ! T if depth is to be read in the file
LOGICAL           :: GREAD_EXT
REAL(KIND=JPRB) :: ZHOOK_HANDLE
!-------------------------------------------------------------------------------
!
!*    2.      SECONDARY VARIABLES
!             -------------------
!
!*    2.2     fields on artificial surfaces only
!             ----------------------------------
!
IF (lhook) CALL dr_hook('GET_TEB_DEPTHS',0,zhook_handle)
!
ydir = 'H'
IF (PRESENT(hdir)) ydir = hdir
!
 CALL open_aux_io_surf(hfilepgd,hfilepgdtype,'FULL  ')
yrecfm='VERSION'
 CALL read_surf(hfilepgdtype,yrecfm,iversion_pgd,iresp,hdir='-')
yrecfm='BUG'
 CALL read_surf(hfilepgdtype,yrecfm,ibugfix_pgd,iresp,hdir='-')
 CALL read_lecoclimap(hfilepgdtype,gecoclimap,gecosg,hdir='-')
 CALL close_aux_io_surf(hfilepgd,hfilepgdtype)
gdim = (iversion_pgd>8 .OR. iversion_pgd==8 .AND. ibugfix_pgd>0)
!
 CALL open_aux_io_surf(hfile,hfiletype,'FULL  ')
yrecfm='VERSION'
 CALL read_surf(hfiletype,yrecfm,iversion_prep,iresp,hdir='-')
yrecfm='BUG'
 CALL read_surf(hfiletype,yrecfm,ibugfix_prep,iresp,hdir='-')
 CALL close_aux_io_surf(hfile,hfiletype)
!
IF (gecosg) THEN
  ALLOCATE(zdata(sum(ntype),SIZE(xdata_d_roof,2)))
ELSE
  ALLOCATE(zdata(ncover,SIZE(xdata_d_roof,2)))
ENDIF
zdata(:,:) = 0.
!
 CALL data_cover_init(ydtco)
ALLOCATE(ydtco%XDATA_GARDEN(SIZE(zdata,1)))
ALLOCATE(ydtco%XDATA_TOWN  (SIZE(zdata,1)))
ALLOCATE(ydtco%XDATA_BLD   (SIZE(zdata,1)))
!
IF (SIZE(zdata,1)/=SIZE(xdata_d_roof,1)) THEN
  CALL default_data_cover(pdata_town=ydtco%XDATA_TOWN,pdata_garden=ydtco%XDATA_GARDEN,&
                          pdata_bld=ydtco%XDATA_BLD)
ELSE
  ydtco%XDATA_TOWN  (:) = dtco%XDATA_TOWN(:)
  ydtco%XDATA_GARDEN(:) = dtco%XDATA_GARDEN(:)
  ydtco%XDATA_BLD   (:) = dtco%XDATA_BLD(:)
ENDIF
!
IF (PRESENT(pd_roof)) THEN
  ysurf='ROOF '
  IF (gecosg) THEN
    zdata = xdata_d_roof
  ELSE
    CALL default_data_cover(pdata_town=ydtco%XDATA_TOWN,pdata_d_roof=zdata)
  ENDIF
  yrecfm0 = 'PAR_RF_LAYER'
  yrecfm1 = 'L_D_ROOF'
  yrecfm2 = 'D_D_ROOF'
  yrecfm3 = 'D_ROOF'
  idata_layer = ndata_roof_layer
  ilu     = SIZE(pd_roof,1)
  ilayer  = SIZE(pd_roof,2)
  yarea   = 'BLD'
END IF
IF (PRESENT(pd_wall)) THEN
  ysurf='WALL '
  IF (gecosg) THEN
    zdata = xdata_d_wall
  ELSE
    CALL default_data_cover(pdata_town=ydtco%XDATA_TOWN,pdata_d_wall=zdata)
  ENDIF
  yrecfm0 = 'PAR_WL_LAYER'
  yrecfm1 = 'L_D_WALL'
  yrecfm2 = 'D_D_WALL'
  yrecfm3 = 'D_WALL'
  idata_layer = ndata_wall_layer
  ilu     = SIZE(pd_wall,1)
  ilayer  = SIZE(pd_wall,2)
  yarea   = 'BLD'
END IF
IF (PRESENT(pd_road)) THEN
  ysurf='ROAD '
  IF (gecosg) THEN
    zdata = xdata_d_road
  ELSE
    CALL default_data_cover(pdata_town=ydtco%XDATA_TOWN,pdata_d_road=zdata)
  ENDIF
  yrecfm0 = 'PAR_RD_LAYER'
  yrecfm1 = 'L_D_ROAD'
  yrecfm2 = 'D_D_ROAD'
  yrecfm3 = 'D_ROAD'
  idata_layer = ndata_road_layer
  ilu     = SIZE(pd_road,1)
  ilayer  = SIZE(pd_road,2)
  yarea   = 'STR'
END IF
IF (PRESENT(pd_floor)) THEN
  ysurf='FLOOR'
  IF (gecosg) THEN
    zdata = xdata_d_floor
  ELSE
    WHERE (ydtco%XDATA_TOWN>0.) 
      zdata(:,1) = 0.01
      zdata(:,2) = 0.04
      zdata(:,3) = 0.10
    END WHERE
  ENDIF
  yrecfm0 = 'PAR_FL_LAYER'
  yrecfm1 = 'L_D_FLOOR'
  yrecfm2 = 'D_D_FLOOR'
  yrecfm3 = 'D_FLOOR'
  idata_layer = ndata_floor_layer
  ilu     = SIZE(pd_floor,1)
  ilayer  = SIZE(pd_floor,2)
  yarea   = 'BLD'
END IF
!
!* read if the depths description are written in the file
IF (iversion_pgd<7 .OR. (iversion_pgd==7 .AND. ibugfix_pgd<=2)) THEN
  gdata = .false.
ELSE
  !
  CALL open_aux_io_surf(hfilepgd,hfilepgdtype,'TOWN  ')
  !
  CALL read_surf(hfilepgdtype,yrecfm1,gdata,iresp,hdir='-')
  !* depths are read in the file
  IF (gdata) THEN
    !* gets number of data layers
    ipar_layer = 0
    CALL read_surf(hfilepgdtype,yrecfm0,ipar_layer,iresp,hdir='-')
    !* gets the data layers depths
    ALLOCATE(zd(ilu,ipar_layer)) 
    DO jl=1,ipar_layer
      WRITE(yrecfm,fmt='(A,I1)') trim(yrecfm2),jl
      CALL read_surf(hfilepgdtype,yrecfm,zd(:,jl),iresp,hdir=ydir)
    END DO
  ENDIF
  !
  CALL close_aux_io_surf(hfilepgd,hfilepgdtype)
  !  
END IF
!
!* depths are read in the file
IF (.NOT.gdata) THEN
  !
!* depths are deduced from the cover types
  ALLOCATE(zd(ilu,ilayer))
  ! 
  IF (iversion_prep>8 .OR. (iversion_prep==8 .AND. ibugfix_prep>=1)) THEN
    CALL open_aux_io_surf(hfile,hfiletype,'TOWN  ')
    CALL read_surf(hfiletype,'WRITE_EXT  ',gread_ext,iresp,hdir='-')
    IF (gread_ext) THEN
      DO jl=1,ilayer
        WRITE(yrecfm,fmt='(A,I1)') trim(yrecfm3),jl
        CALL read_surf(hfiletype,yrecfm,zd(:,jl),iresp,hdir=ydir)      
      END DO 
    ENDIF
    CALL close_aux_io_surf(hfile,hfiletype)
  ELSE
    gread_ext = .false.
  ENDIF
  !
  IF (.NOT.gread_ext) THEN
    !
    IF (gdim.AND.gecosg) THEN
      ALLOCATE(gcover(sum(ntype)))
    ELSE
      ALLOCATE(gcover(ncover))
    ENDIF
    !
    CALL open_aux_io_surf(hfilepgd,hfilepgdtype,'FULL  ')   
    !* reading of the cover to obtain the thickness of layers
    CALL old_name(hfilepgdtype,'COVER_LIST      ',yrecfm,'-')
    CALL read_surf(hfilepgdtype,yrecfm,gcover(:),iresp,hdir='-')
    !* reading of the cover fractions
    IF (ydir/='E'.OR.nrank==npio) THEN
      ALLOCATE(zcover(ilu,count(gcover)))
    ELSE
      ALLOCATE(zcover(0,0))
    ENDIF
    yrecfm='COVER'
    CALL read_surf_cov(hfilepgdtype,yrecfm,zcover(:,:),gcover,iresp,&
             hdir=ydir)
    CALL close_aux_io_surf(hfilepgd,hfilepgdtype)
    !
    ALLOCATE(zpar_d(ilu,idata_layer))
    IF (nrank==npio) THEN
      !* deduces the depths of each layer
      DO jl=1,idata_layer
        CALL av_pgd (ydtco,zpar_d(:,jl), zcover, zdata(:,jl),yarea,'ARI',gcover)
      END DO
    ENDIF
    CALL data_cover_init(ydtco)
    DEALLOCATE(zdata)
    DEALLOCATE(gcover,zcover)
    !
    IF (iversion_prep<7 .OR. (iversion_prep==7 .AND. ibugfix_prep<=2)) THEN
      !* ind version of TEB, the computational grid was equal to the data grid
      zd(:,:) = zpar_d(:,:)
    ELSEIF (nrank==npio) THEN
      !* recomputes the grid from the available data
      CALL thermal_layers_conf(ysurf,zpar_d,zd)
    END IF
    DEALLOCATE(zpar_d)
    !
  ENDIF  
  !  
ENDIF
!
IF (PRESENT(pd_roof )) pd_roof  = zd
IF (PRESENT(pd_wall )) pd_wall  = zd
IF (PRESENT(pd_road )) pd_road  = zd
IF (PRESENT(pd_floor)) pd_floor = zd
!
DEALLOCATE(zd)
!
IF (lhook) CALL dr_hook('GET_TEB_DEPTHS',1,zhook_handle)
!-------------------------------------------------------------------------------
!
END SUBROUTINE get_teb_depths