fapair.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 fapair(PABC, PFD_SKY, PIA, PLAI, PXMUS, PSSA_SUP, PSSA_INF, &
           PB_SUP, PB_INF, PALB_VEG, PALB_SOIL, OSHADE,            &
           PFAPR, PFAPR_BS, PLAI_EFF, PIACAN,                      &
           PIACAN_SHADE, PIACAN_SUNLIT, PFRAC_SUN                  )
!   #########################################################################
!
!!****  *FAPAIR*  
!!
!!    PURPOSE
!!    -------
!!    Calculates fraction of absorbed photosynthetic active radiation (FAPAR) and 
!!     fraction of absorbed near infrared (FAPIR) of vegetation and bare soil.
!!              
!!**  METHOD
!!    ------
!!    Carrer et al. 
!!
!!    EXTERNAL
!!    --------
!!    none
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!      
!!    USE MODD_SURF_PAR
!!    USE MODD_CSTS
!!    USE MODD_CO2V_PAR
!!    USE MODI_CCETR_PAIR
!!
!!    REFERENCE
!!    ---------
!!     Carrer et al, 2013 (doi:10.1002/jgrg20070)
!!      
!!    AUTHOR
!!    ------
!!      D. Carrer          * Meteo-France *
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    01/04/2011 
!!      Commented by C. Delire 07/13
!!      C. Delire   08/13 : moved calculation of diffuse fraction from here to radiative_transfert.F90
!!      A. Boone    02/17 : corrected computation of PFAPR_BS
!!
!!-------------------------------------------------------------------------------
USE modd_surf_par,       ONLY : xundef
USE modd_co2v_par,   ONLY : xk_sup, xk_inf, xxsi_sup, xxsi_inf ! clumping index parameters (Carrer et al 2.1.3)
!
USE modd_surfex_mpi, ONLY : nrank
USE modi_ccetr_pair  
!
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb
!
!*       0.     DECLARATIONS
!               ------------
!
IMPLICIT NONE
!
!*      0.1    declarations of arguments
!
REAL, DIMENSION(:), INTENT(IN) :: PABC    ! abscissa needed for integration
!                                         ! of net assimilation and stomatal
!                                         ! conductance over canopy depth
REAL, DIMENSION(:), INTENT(IN) :: PFD_SKY ! fraction of diffused radiation in sky
REAL, DIMENSION(:), INTENT(IN) :: PIA     ! incident PAR or NIR
REAL, DIMENSION(:), INTENT(IN) :: PLAI    ! leaf area index
REAL, DIMENSION(:), INTENT(IN) :: PXMUS   ! cosine of solar zenith angle
REAL,               INTENT(IN) :: PSSA_SUP, PSSA_INF  !single scatering albedo (PAR or NIR)
REAL, DIMENSION(:), INTENT(IN) :: PB_SUP, PB_INF      !sigma parameter in clumping (Table 1, eq4)
REAL, DIMENSION(:), INTENT(IN) :: PALB_VEG, PALB_SOIL
LOGICAL, DIMENSION(:), INTENT(IN) :: OSHADE   ! OSHADE = if 1 shading activated
!
REAL, DIMENSION(:), INTENT(OUT) :: PFAPR
REAL, DIMENSION(:), INTENT(OUT) :: PFAPR_BS
REAL, DIMENSION(:), OPTIONAL,   INTENT(OUT) :: PLAI_EFF
!
REAL, DIMENSION(:,:), OPTIONAL, INTENT(OUT) :: PIACAN ! PAR in the canopy at different gauss level
REAL, DIMENSION(:,:), OPTIONAL, INTENT(OUT) :: PIACAN_SHADE  ! PAR in the canopy at different gauss level
REAL, DIMENSION(:,:), OPTIONAL, INTENT(OUT) :: PIACAN_SUNLIT ! PAR in the canopy at different gauss level
REAL, DIMENSION(:,:), OPTIONAL, INTENT(OUT) :: PFRAC_SUN     ! fraction of sunlit leaves
!
!*      0.2    declarations of local variables
!
!
REAL, DIMENSION(SIZE(PLAI))   :: ZXIA, ZXIA_SUP, ZKMUSP_SUP, ZKMUSP_INF
REAL, DIMENSION(SIZE(PLAI))   :: ZB_DR_SUP, ZB_DR_INF, ZOMEGA_DR_SUP, ZOMEGA_DR_INF, &
                                 ZOMEGA_DF_SUP, ZOMEGA_DF_INF
!                                ZXIA  = abs. radiation of vegetation
REAL, DIMENSION(SIZE(PLAI))   :: ZTR, ZFD_VEG, ZFD_SUP, ZLAI_EFF0, ZLAI_EFF
!                                ZTR = transmittance
!                                ZFD_VEG, ZFD_SUP = fraction of radiation diffused by the considered medium (vegetation)     
!REAL, DIMENSION(SIZE(PLAI))  :: ZXIA_SUNLIT, ZXIA_SHADE, ZLAI_SUNLIT, ZLAI_SHADE
!                                ZXIA_SUNLIT = absorbed PAR of sunlit leaves
!                                ZXIA_SHADE = absorbed PAR of shaded leaves
!                                ZLAI_SUNLIT = LAI of sunlit leaves
!                                ZLAI_SHADE = LAI of shaded leaves
!REAL, DIMENSION(SIZE(PLAI))    :: ZRN_SUNLIT, ZRN_SHADE
REAL, DIMENSION(SIZE(PLAI),SIZE(PABC)) :: ZIACAN, ZIACAN_SUNLIT, ZIACAN_SHADE, ZFRAC_SUN
REAL                           :: ZABC, ZWEIGHT, ZCOEF, ZSUP, ZINF, &
                                  ZSSA_SUP, ZSSA_INF, ZB_DF_SUP, ZB_DF_INF
!                                 ZABC    = abscissa needed for integration
!                                            of net assimilation and stomatal 
!                                            conductance over canopy depth (working scalar)
!                                 ZSUP, ZINF = d_sup and d_inf from Table 1. Carrer et al.
!
INTEGER                        :: JINT, I ! index for loops
REAL(KIND=JPRB) :: ZHOOK_HANDLE
!-------------------------------------------------------------------------------
IF (lhook) CALL dr_hook('FAPAIR',0,zhook_handle)
!
! initialisation
!
ztr(:)      = 1.0
!
zxia_sup(:) = 0.
!
zfd_veg(:)  = 0.
zfd_sup(:)  = 0.
!
!ZXIA_SUNLIT(:) = 0.
!ZXIA_SHADE(:)  = 0.
!ZLAI_SUNLIT(:) = 0.
!ZLAI_SHADE(:)  = 0.
!
zlai_eff(:) = 0.
!
ziacan(:,:)        = 0.
ziacan_sunlit(:,:) = 0.
ziacan_shade(:,:)  = 0.
zfrac_sun(:,:)     = 0.
!
pfapr(:)    = 0.
pfapr_bs(:) = 0.
!PRN_SHADE(:) = 0.
!PRN_SUNLIT(:) = 0.
!
!
IF (pabc(SIZE(pabc)).GT.0.8) zfd_veg(:) = min(pfd_sky(:),1.)   
! set param sup / inf
!
zssa_sup = sqrt(1.-pssa_sup)
zssa_inf = sqrt(1.-pssa_inf)
!
zsup = - 0.461 * xxsi_sup + 3.8     ! d_sup Carrer et al Table 1
zinf = - 0.461 * xxsi_inf + 3.8     ! d_inf Carrer et al Table 1
!
DO i=1,SIZE(pia)
  IF (pia(i).NE.0.) THEN
    zkmusp_sup(i) = exp(-xk_sup*(acos(pxmus(i)))**zsup)
    zkmusp_inf(i) = exp(-xk_inf*(acos(pxmus(i)))**zinf)
    ! direct case
    ! Directional albedo of upper/lower layer
    zb_dr_sup(i) = 1.-(1.-zssa_sup)/(1.+2.*pxmus(i)*zssa_sup) ! Carrer et al. b_dr Table 1
    zb_dr_inf(i) = 1.-(1.-zssa_sup)/(1.+2.*pxmus(i)*zssa_inf)
    ! CLUMPING INDEX 
    zomega_dr_sup(i) = 1. / (1.+ pb_sup(i)*zkmusp_sup(i))     ! Carrer et al. eq. 4a
    zomega_dr_inf(i) = 1. / (1.+ pb_inf(i)*zkmusp_inf(i))
    ! diffuse case
    ! CLUMPING INDEX
    zomega_df_sup(i) = (1.+pb_sup(i)/2.)/(1.+pb_sup(i))       ! Carrer et al. eq. 4b
    zomega_df_inf(i) = (1.+pb_inf(i)/2.)/(1.+pb_inf(i))
  ENDIF
ENDDO
!
! Non-directional albedo of diffuse 
zb_df_sup = 1.-(1.-zssa_sup)/(1.+ zssa_sup)
zb_df_inf = 1.-(1.-zssa_inf)/(1.+ zssa_inf)
!
! Integration over the canopy: SIZE(PABC) increments
! are used to approximate the integral. And to calculate 
! absorded fluxes within the canopy and in the bare soil  
DO jint = SIZE(pabc),1,-1
!
  zabc = 1.                                ! normalized height unit of the layer above
  IF (jint.LT.SIZE(pabc)) zabc = pabc(jint+1)
  zweight = zabc - pabc(jint)
  !
  IF (pabc(jint).GT.0.8) THEN
    !  Compute transmittance of each level  
    CALL ccetr_pair (jint, pabc(jint), zabc, pia, pxmus, zb_dr_sup, &
                     zomega_dr_sup, zomega_df_sup, zb_df_sup, plai, &
                     palb_veg, palb_soil, pfd_sky, zfd_veg, ztr,    &
                     zxia, zlai_eff0              )
  ELSE
    CALL ccetr_pair (jint, pabc(jint), zabc, pia, pxmus, zb_dr_inf, &
                     zomega_dr_inf, zomega_df_inf, zb_df_inf, plai, &
                     palb_veg, palb_soil, pfd_sky, zfd_veg, ztr,    &
                     zxia, zlai_eff0              )
  ENDIF          
  ! 
  DO i=1,SIZE(pia)
    !
    zxia(i)        = max(0.,zxia(i))
    ziacan(i,jint) = max(0.,zxia(i)-zxia_sup(i))
    zxia_sup(i)    = zxia(i)
    !
    zlai_eff0(i) = max(0.,zlai_eff0(i))
    zlai_eff(i)  = zlai_eff(i) + zlai_eff0(i)
    !
    !calculate a FAPAR/FAPIR of the entire canopy
    pfapr(i)  = pfapr(i) + ziacan(i,jint)
    !
    !------------------------------------------------------
    ! If LSHADE=0 no shading, only sunlit leaves
    ! If LSHADE=1 shading
    ! PIACAN is used to calculate An of each level within the canopy in cotwores
    ! ZIACAN_SUNLIT used for net assimilation of a sunlit leave in COTWO
    ! ZIACAN_SHADE used  in A-gs for net assimilation of a shaded leave in COTWO
    IF (oshade(i)) THEN
      !
      !sunlit leaves
      !absorbed PAR of an equivalent canopy representative of the layer of leaves  eq. (8)
      zcoef = (1.0-zfd_sup(i))/ztr(i)+ zfd_sup(i)    
      ziacan_sunlit(i,jint) =             zcoef/(zweight*max(0.0001,plai(i)))*ziacan(i,jint)    
      !not sunlit leaves
      ziacan_shade(i,jint)  = max(0.,zfd_sup(i)/(zweight*max(0.0001,plai(i)))*ziacan(i,jint))
      !
      !ZXIA_SUNLIT(I) = ZXIA_SUNLIT(I) + ZWEIGHT*ZTR(I)      *ZIACAN_SUNLIT(I,JINT)
      !ZLAI_SUNLIT(I) = ZLAI_SUNLIT(I) + ZWEIGHT*ZTR(I)*ZCOEF*PLAI(I)
      !
      !ZXIA_SHADE(I)  = ZXIA_SHADE(I)  + ZWEIGHT*(1-ZTR(I))           *ZIACAN_SHADE(I,JINT)
      !ZLAI_SHADE(I)  = ZLAI_SHADE(I)  + ZWEIGHT*(1-ZTR(I))*ZFD_SUP(I)*PLAI(I)
      !
      zfrac_sun(i,jint) = ztr(i)  !fraction of sunlit leaves
      !      
    ELSE
      !
      ziacan_sunlit(i,jint) = max(0.,ziacan(i,jint)/(zweight*max(0.0001,plai(i))))
      !ZLAI_SUNLIT(I) = ZLAI_SUNLIT(I) + ZWEIGHT*PLAI(I)
      !
    ENDIF
    !
    zfd_sup(i) = zfd_veg(i)
    !
    ENDDO
  !
END DO
!
!
WHERE (pia(:).NE.0.)
  pfapr(:) = pfapr(:) / pia(:)
  pfapr_bs(:)= ztr(:)*(1.-palb_soil(:)*(1. - palb_veg(:)*(1.-ztr(:))))
END WHERE
!
!WHERE (ZLAI_SHADE(:) .NE.0.) ZRN_SHADE(:)  = ZXIA_SHADE(:) / ZLAI_SHADE(:)
!WHERE (ZLAI_SUNLIT(:).NE.0.) ZRN_SUNLIT(:) = ZXIA_SUNLIT(:)/ ZLAI_SUNLIT(:)
!
IF (PRESENT(plai_eff))      plai_eff      = zlai_eff
IF (PRESENT(piacan))        piacan        = ziacan
IF (PRESENT(piacan_sunlit)) piacan_sunlit = ziacan_sunlit
IF (PRESENT(piacan_shade))  piacan_shade  = ziacan_shade
IF (PRESENT(pfrac_sun))     pfrac_sun     = zfrac_sun
!
IF (lhook) CALL dr_hook('FAPAIR',1,zhook_handle)
!
END SUBROUTINE fapair