soilgrid.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 soilgrid(PSOILGRID, PSOILDEPTH, PDG, KWG_LAYER  )

!     ##########################################################################
!
!!****  *SOILGRID*  
!!
!!    PURPOSE
!!    -------
!
!     Calculates the soil grid configuration using a reference grid
!     Also compute the root fraction
!         
!     
!!**  METHOD
!!    ------
!
!     Direct calculation
!
!!    EXTERNAL
!!    --------
!
!     None
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!
!!    REFERENCE
!!    ---------
!!
!!    Noilhan and Planton (1989)
!!    Belair (1995)
!!    Boone (2000)
!!    Boone et al. (2000)
!!    Habets et al. (2003)
!!    Decharme et al. (2011)
!!      
!!    AUTHOR
!!    ------
!!      A. Boone           * Meteo-France *
!!      new version :
!!      B. Decharme        * Meteo-France *
!!
!!    MODIFICATIONS
!!    -------------
!!      Original     12/04/03
!!      new version :10/08/2011
!!      modif       :   09/2012 soildepth can reach 12m (permafrost)
!!                              bug coef algo
!-------------------------------------------------------------------------------
!
!*       0.     DECLARATIONS
!               ------------
!
USE modd_surf_par, ONLY : xundef, nundef
USE modd_isba_par, ONLY : noptimlayer, xoptimgrid
!
USE modi_abor1_sfx
!
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb
!
IMPLICIT NONE
!
!*      0.1    declarations of arguments
!
REAL,    DIMENSION(:),     INTENT(IN)    :: psoilgrid   ! reference soil grid          (m)
REAL,    DIMENSION(:),   INTENT(IN)    :: psoildepth  ! total soil depth             (m)   
REAL,    DIMENSION(:,:), INTENT(OUT)   :: pdg         ! depth of base of soil layers (m)
INTEGER, DIMENSION(:),   INTENT(OUT)   :: kwg_layer   ! last layers for soil moisture
!
!*      0.2    declarations of local variables
!
REAL,DIMENSION(SIZE(PDG,1)) :: zref
!
REAL               :: zwork
!
INTEGER            :: ini,inl
INTEGER            :: jj,jl
!
LOGICAL            :: loptimgrid
!
REAL(KIND=JPRB) :: zhook_handle
!-------------------------------------------------------------------------------
!        0.     Initialization
!               --------------
!
IF (lhook) CALL dr_hook('SOILGRID',0,zhook_handle)
!
ini    = SIZE(pdg,1)
inl    = SIZE(pdg,2)
!
kwg_layer(:) = 0
zref(:) = xundef
!
!-------------------------------------------------------------------------------
!
!*       1.     Grid configuration
!               ------------------
!
!*       1.1    Check consistency
!               -----------------
!
loptimgrid=.false.
IF(inl==noptimlayer) THEN
  IF( all(psoilgrid(1:inl)==xoptimgrid(1:noptimlayer)) ) loptimgrid=.true.
ENDIF
!
!*       1.2    Assign soil layer depths if ECOCLIMAP
!               -------------------------------------
!
IF(loptimgrid)THEN
  !
  !Optimized ECOCLIMAP soil grid
  CALL optimsoilgrid
  !
ELSE
  !
  WHERE(psoildepth(:)/=xundef)
    pdg(:,1)=min(0.01,psoilgrid(1))
  ELSEWHERE
    pdg(:,1)=xundef
  ENDWHERE
  !
    DO jj=1,ini 
      !
      IF( psoildepth(jj)==xundef )THEN
        !
        pdg(jj,:) = xundef              
        kwg_layer(jj) = nundef
        !
      ELSE
        !
        DO jl=2,inl
          !
          pdg(jj,jl) = psoilgrid(jl)
          !
          IF ( psoilgrid(jl)-psoilgrid(jl-1)<=0.3 ) THEN       
            zwork = abs(psoilgrid(jl)-psoildepth(jj))
            IF(zwork<=zref(jj))THEN
              kwg_layer(jj) = jl
              zref(jj) = zwork
            ENDIF
          ELSEIF(psoildepth(jj)>=(psoilgrid(jl)*0.3+psoilgrid(jl-1)*0.7))THEN
            kwg_layer(jj) = jl
          ENDIF
        ENDDO
        !
      ENDIF
    ENDDO
  !
ENDIF
!
IF(any(kwg_layer(:)==0))THEN
  CALL abor1_sfx('SOILGRID: WITH CISBA=DIF NWG_LAYER MUST BE DEFINED FOR EACH POINT')
ENDIF
!
IF (lhook) CALL dr_hook('SOILGRID',1,zhook_handle)
!
!-------------------------------------------------------------------------------
CONTAINS
!-------------------------------------------------------------------------------
!
SUBROUTINE optimsoilgrid
!
USE modd_reprod_oper, ONLY : cdgdif
!
IMPLICIT NONE
!
! declarations of local variables
!
INTEGER, PARAMETER                  :: ndlim = 13
!
REAL, DIMENSION(NDLIM), PARAMETER   :: zdlim = &
      (/1.25,1.50,1.75,2.00,2.25,2.50,2.75,3.00,3.50,4.00,4.50,5.00,5.50/)
!
REAL,DIMENSION(SIZE(PDG,1)) :: zdg_water
!
LOGICAL            :: lwork
REAL(KIND=JPRB)    :: zhook_handle
!
!-------------------------------------------------------------------------------
! init
!
IF (lhook) CALL dr_hook('SOILGRID:OPTIMSOILGRID',0,zhook_handle)
!
zdg_water(:) = xundef
!
!-------------------------------------------------------------------------------
!a. Arranged depth
!
IF(cdgdif=='ROOT')THEN
!
     DO jj=1,ini
        IF(psoildepth(jj)<=1.1)THEN
          zdg_water(jj)=min(1.0,psoildepth(jj))
        ELSEIF(psoildepth(jj)>1.1.AND.psoildepth(jj)<=1.25)THEN
          zdg_water(jj)=1.25          
        ELSEIF(psoildepth(jj)>5.50.AND.psoildepth(jj)<=8.00)THEN
              zdg_water(jj)=8.00
        ELSEIF(psoildepth(jj)>8.00.AND.psoildepth(jj)<xundef)THEN
              zdg_water(jj)=12.00 ! Permafrost case
        ELSE
          DO jl=1,ndlim-1         
             IF(psoildepth(jj)>zdlim(jl).AND.psoildepth(jj)<=zdlim(jl+1))THEN
               zdg_water(jj)=merge(zdlim(jl),zdlim(jl+1),psoildepth(jj)<(0.8*zdlim(jl)+0.2*zdlim(jl+1)))
             ENDIF
          ENDDO
        ENDIF
     ENDDO
!
ELSE
!
     DO jj=1,ini
        IF(psoildepth(jj)<1.25)THEN
          zdg_water(jj)=min(1.0,psoildepth(jj))
        ELSEIF(psoildepth(jj)>=5.50.AND.psoildepth(jj)<6.50)THEN
              zdg_water(jj)=5.50
        ELSEIF(psoildepth(jj)>=6.50.AND.psoildepth(jj)<10.50)THEN
              zdg_water(jj)=8.00
        ELSEIF(psoildepth(jj)>=10.50.AND.psoildepth(jj)<xundef)THEN
              zdg_water(jj)=12.00 ! Permafrost case
        ELSE
          DO jl=1,ndlim-1         
             IF(psoildepth(jj)>=zdlim(jl).AND.psoildepth(jj)<zdlim(jl+1))THEN
               zdg_water(jj)=merge(zdlim(jl),zdlim(jl+1),psoildepth(jj)<(0.4*zdlim(jl)+0.6*zdlim(jl+1)))
             ENDIF
          ENDDO
        ENDIF
     ENDDO
!
ENDIF
!
!-------------------------------------------------------------------------------
!b. General cases
!
  DO jj=1,ini
    !
    IF(psoildepth(jj)==xundef)THEN
      !
      pdg(jj,:) = xundef      
      kwg_layer(jj) = nundef    
      !
    ELSE
      !
      pdg(jj,:) = psoilgrid(:)
      !      
      lwork=(zdg_water(jj)<=1.0.OR.&
             zdg_water(jj)==1.5.OR.&
             zdg_water(jj)==2.0.OR.&
             zdg_water(jj)==3.0.OR.&
             zdg_water(jj)==5.0.OR.&
             zdg_water(jj)==8.0.OR.&
             zdg_water(jj)==12.0   )
      ! 
      IF (lwork) THEN    
        DO jl=2,inl      
          zwork = abs(psoilgrid(jl)-zdg_water(jj))                 
          IF(zwork<=zref(jj))THEN
            kwg_layer(jj)=jl
            zref(jj)=zwork
          ENDIF
        ENDDO           
      ENDIF
      !
    ENDIF
    !
  ENDDO
!
!-------------------------------------------------------------------------------
!c. Particular cases
!
WHERE (zdg_water(:)==1.25)
  kwg_layer(:) = 9
  pdg(:,9) = zdg_water(:)
ELSEWHERE (zdg_water(:)==1.75 .OR. zdg_water(:)==2.25)
  kwg_layer(:) = 10
  pdg(:,10) = zdg_water(:)
  WHERE (zdg_water(:)==1.75) pdg(:,9) = 1.25
ELSEWHERE (zdg_water(:)==2.50 .OR. zdg_water(:)==2.75 .OR. zdg_water(:)==3.50)
  kwg_layer(:) = 11
  pdg(:,11) = zdg_water(:)
ELSEWHERE (zdg_water(:)==4.00 .OR. zdg_water(:)==4.50 .OR. zdg_water(:)==5.50)
  kwg_layer(:) = 12
  pdg(:,12) = zdg_water(:)
ENDWHERE
!
IF (lhook) CALL dr_hook('SOILGRID:OPTIMSOILGRID',1,zhook_handle)
!
END SUBROUTINE optimsoilgrid
!
!-------------------------------------------------------------------------------
END SUBROUTINE soilgrid