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),SIZE(PDG,3)) :: zref
!
REAL               :: zwork
!
INTEGER            :: ini,inl,ipatch
INTEGER            :: jj,jl,jpatch
!
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)
ipatch = SIZE(pdg,3)
!
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 jpatch=1,ipatch
    DO jj=1,ini 
      !
      IF( psoildepth(jj,jpatch)==xundef )THEN
        !
        pdg(jj,:,jpatch) = xundef              
        kwg_layer(jj,  jpatch) = nundef
        !
      ELSE
        !
        DO jl=2,inl
          !
          pdg(jj,jl,jpatch) = psoilgrid(jl)
          !
          IF ( psoilgrid(jl)-psoilgrid(jl-1)<=0.3 ) THEN       
            zwork = abs(psoilgrid(jl)-psoildepth(jj,jpatch))
            IF(zwork<=zref(jj,jpatch))THEN
              kwg_layer(jj,jpatch) = jl
              zref(jj,jpatch) = zwork
            ENDIF
          ELSEIF(psoildepth(jj,jpatch)>=(psoilgrid(jl)*0.3+psoilgrid(jl-1)*0.7))THEN
            kwg_layer(jj,jpatch) = jl
          ENDIF
        ENDDO
        !
      ENDIF
    ENDDO
  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),SIZE(PDG,3)) :: 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 jpatch=1,ipatch
     DO jj=1,ini
        IF(psoildepth(jj,jpatch)<=1.1)THEN
          zdg_water(jj,jpatch)=min(1.0,psoildepth(jj,jpatch))
        ELSEIF(psoildepth(jj,jpatch)>1.1.AND.psoildepth(jj,jpatch)<=1.25)THEN
          zdg_water(jj,jpatch)=1.25          
        ELSEIF(psoildepth(jj,jpatch)>5.50.AND.psoildepth(jj,jpatch)<=8.00)THEN
              zdg_water(jj,jpatch)=8.00
        ELSEIF(psoildepth(jj,jpatch)>8.00.AND.psoildepth(jj,jpatch)<xundef)THEN
              zdg_water(jj,jpatch)=12.00 ! Permafrost case
        ELSE
          DO jl=1,ndlim-1         
             IF(psoildepth(jj,jpatch)>zdlim(jl).AND.psoildepth(jj,jpatch)<=zdlim(jl+1))THEN
               zdg_water(jj,jpatch)=merge(zdlim(jl),zdlim(jl+1),psoildepth(jj,jpatch)<(0.8*zdlim(jl)+0.2*zdlim(jl+1)))
             ENDIF
          ENDDO
        ENDIF
     ENDDO
  ENDDO
!
ELSE
!
  DO jpatch=1,ipatch
     DO jj=1,ini
        IF(psoildepth(jj,jpatch)<1.25)THEN
          zdg_water(jj,jpatch)=min(1.0,psoildepth(jj,jpatch))
        ELSEIF(psoildepth(jj,jpatch)>=5.50.AND.psoildepth(jj,jpatch)<6.50)THEN
              zdg_water(jj,jpatch)=5.50
        ELSEIF(psoildepth(jj,jpatch)>=6.50.AND.psoildepth(jj,jpatch)<10.50)THEN
              zdg_water(jj,jpatch)=8.00
        ELSEIF(psoildepth(jj,jpatch)>=10.50.AND.psoildepth(jj,jpatch)<xundef)THEN
              zdg_water(jj,jpatch)=12.00 ! Permafrost case
        ELSE
          DO jl=1,ndlim-1         
             IF(psoildepth(jj,jpatch)>=zdlim(jl).AND.psoildepth(jj,jpatch)<zdlim(jl+1))THEN
               zdg_water(jj,jpatch)=merge(zdlim(jl),zdlim(jl+1),psoildepth(jj,jpatch)<(0.4*zdlim(jl)+0.6*zdlim(jl+1)))
             ENDIF
          ENDDO
        ENDIF
     ENDDO
  ENDDO
!
ENDIF
!
!-------------------------------------------------------------------------------
!b. General cases
!
DO jpatch=1,ipatch
  DO jj=1,ini
    !
    IF(psoildepth(jj,jpatch)==xundef)THEN
      !
      pdg(jj,:,jpatch) = xundef      
      kwg_layer(jj,  jpatch) = nundef    
      !
    ELSE
      !
      pdg(jj,:,jpatch) = psoilgrid(:)
      !      
      lwork=(zdg_water(jj,jpatch)<=1.0.OR.&
             zdg_water(jj,jpatch)==1.5.OR.&
             zdg_water(jj,jpatch)==2.0.OR.&
             zdg_water(jj,jpatch)==3.0.OR.&
             zdg_water(jj,jpatch)==5.0.OR.&
             zdg_water(jj,jpatch)==8.0.OR.&
             zdg_water(jj,jpatch)==12.0   )
      ! 
      IF (lwork) THEN    
        DO jl=2,inl      
          zwork = abs(psoilgrid(jl)-zdg_water(jj,jpatch))                 
          IF(zwork<=zref(jj,jpatch))THEN
            kwg_layer(jj,jpatch)=jl
            zref(jj,jpatch)=zwork
          ENDIF
        ENDDO           
      ENDIF
      !
    ENDIF
    !
  ENDDO
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