cotwores.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 cotwores(PTSTEP, HPHOTO, OTR_ML, OSHADE,                   &
            pvegtype, ostressdef, pah, pbh, pf2i, pdmax,              &
            ppoi, pcsp, ptg, pf2, psw_rad, pra, pqa, pqsat, ple,      &
            ppsnv, pdelta, plai, prhoa, pzenith, pfzero, pepso,       &
            pgamm, pqdgamm, pgmes,  pgc, pqdgmes, pt1gmes, pt2gmes,   &
            pamax, pqdamax, pt1amax, pt2amax, pffv,                   &
            piacan_sunlit, piacan_shade, pfrac_sun, piacan,           &
            pabc, pan, panday, prs, panfm, pgpp, panf, presp_leaf     ) 
!   #########################################################################
!
!!****  *COTWORES*  
!!
!!    PURPOSE
!!    -------
!!
!!    Calculates net assimilation of CO2 and leaf conductance.
!!              
!!**  METHOD
!!    ------
!!    Calvet et al. 1998 Forr. Agri. Met. [from model of Jacobs(1994)]
!!
!!    EXTERNAL
!!    --------
!!    none
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!      
!!    USE MODD_CST
!!    USE MODD_CO2V_PAR
!!    USE MODI_COTWO
!!    USE MODI_CCETR
!!    USE MODE_THERMOS
!!
!!    REFERENCE
!!    ---------
!!
!!    Calvet et al. 1998 Forr. Agri. Met. 
!!      
!!    AUTHOR
!!    ------
!!
!!      A. Boone           * Meteo-France *
!!      (following Belair)
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    27/10/97 
!!      V. Masson and V. Rivailland 12/2003 modificatino of ISBA routines order
!!      L. Jarlan   27/10/04 : add of T2 to calculate soil respiration and use
!!                              of CRESPSL key to manage the calculation of soil
!!                              respiration
!!                             PAN et PPST in kgCO2 m-2 s-1 to be fully
!!                              compatible with vegetation growth module (lailoss.f90)
!!      P Le Moigne 09/2005 AGS modifs of L. Jarlan
!!      S. Lafont      03/09 : change units of EPSO GAMM ANDAY
!!      A.L. Gibelin   06/09 : suppress evolution of [CO2] in canopy
!!      A.L. Gibelin   06/09 : move calculations of some CO2 fluxes
!!      A.L. Gibelin   06/09 : add RESP_LEAF
!!      A.L. Gibelin   07/09 : ensure coherence between cotwores and cotworestress
!!      A.L. Gibelin   07/09 : Suppress PPST and PPSTF as outputs, and diagnose GPP
!!        S. Lafont    03/11 : Correct a bug fopr grassland below wilting point
!!      D. Carrer      04/11 : new radiative transfert 
!!      A. Boone       11/11 : add rsmax to MODD_ISBA_PAR
!!      B. Decharme    05/12 : Bug : flood fraction in COTWORES
!!                                   Optimization
!!      R. Alkama      04/12 : add 6 new tree vegtype (9 instead 3)
!!      C. Delire      01/14 : vertical profile of dark respiration for tropical forest 
!!                             (GTROP)   with Carrer radiative transfer (OTR_ML = T)               
!!Seferian & Delire  06/2015 : generalization of (i) linear water-stress reponse
!                              and (ii) exponential decrease of autothrophic respiration to all woody PFTs
!!      B. Decharme    07/15 : Suppress some numerical adjustement for F2 
!!
!-------------------------------------------------------------------------------
!
USE modd_csts,           ONLY : xmd, xtt, xlvtt
USE modd_isba_par,       ONLY : xrs_max, xdenom_min
USE modd_co2v_par,       ONLY : xparcf, xmco2, xdmax_ags,       &
                                xdmaxx, xdmaxn, xaw, xbw, xasw                              
USE modd_data_cover_par, ONLY : nvt_tebd, nvt_trbe, nvt_bone,   &
                                nvt_trbd, nvt_tebe, nvt_tene,   &
                                nvt_bobd, nvt_bond, nvt_shrb
USE modd_surf_par,       ONLY : xundef
!
USE modi_ccetr
USE modi_cotwo
!
!*       0.     DECLARATIONS
!               ------------
!
!
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb
!
IMPLICIT NONE
!
!*      0.1    declarations of arguments
!
!
REAL,                INTENT(IN)  :: ptstep      ! time step
 CHARACTER(LEN=*),    INTENT(IN)  :: hphoto      ! Kind of photosynthesis
!                                               ! 'NON'
!                                               ! 'AGS'
!                                               ! 'LAI'
!                                               ! 'AST'
!                                               ! 'LST'
LOGICAL,             INTENT(IN)  :: otr_ml      ! new TR
LOGICAL, DIMENSION(:),INTENT(IN) :: oshade
!
REAL, DIMENSION(:,:),INTENT(IN)  :: pvegtype
!                                PVEGTYPE  = type de vegetation (1 a 9)
!
LOGICAL,DIMENSION(:),INTENT(IN)  :: ostressdef
REAL,DIMENSION(:),   INTENT(IN)  :: pah, pbh, pf2i, pdmax
!                                    PAH, PBH  = coefficients for universal herbaceous
!                                                stress relation 
!                                    OSTRESSDEF   = water stress vegetation comportement 
!                                                (true:defensif false:offensif)
!                                    PF2I      = critical normalized soil water stress    
!                                    PDMAX     = maximum saturation deficit of 
!                                                  atmosphere tolerate by vegetation
!
REAL, DIMENSION(:),  INTENT(IN)  :: ppoi     ! Gaussian weights (as above)
!
REAL,DIMENSION(:),   INTENT(IN)  :: pcsp, ptg, pf2, psw_rad, pra 
!                                    PCSP  = atmospheric concentration of CO2
!                                    PTG   = updated leaf temperature
!                                    PF2   = normalized soil water stress factor
!                                    PSW_RAD = incident solar radiation
!                                    PRA   = aerodynamic resistance
!
REAL,DIMENSION(:),   INTENT(IN)  :: pqa, pqsat, ple, ppsnv, pdelta, plai, prhoa
!                                    PQA   = atmospheric mixing ratio
!                                    PQSAT = surface saturation mixing ratio
!                                    PLE   = evapotranspiration (kgH2O kgAir-1 m s-1)
!                                    PPSNV = snow cover fraction
!                                    PDELTA= fraction of the foliage covered
!                                        by intercepted water
!                                    PLAI  = leaf area index
!                                    PRHOA = air density
!
REAL,DIMENSION(:),    INTENT(IN)  :: pzenith
!                                    PZENITH = solar zenith angle needed 
!                                    for computation of diffusuion of solar
!                                    radiation: for CO2 model.
!
REAL,DIMENSION(:),    INTENT(IN)  :: pfzero, pepso, pgamm, pqdgamm, pgmes, pgc,     &
                                     pqdgmes, pt1gmes, pt2gmes, pamax, pqdamax,     &
                                     pt1amax, pt2amax      
!                                    PFZERO    = ideal value of F, no photorespiration or 
!                                                saturation deficit
!                                    PEPSO     = maximum initial quantum use efficiency 
!                                                (kgCO2 J-1 PAR)
!                                    PGAMM     = CO2 conpensation concentration (ppmv)
!                                    PQDGAMM   = Log of Q10 function for CO2 conpensation 
!                                                concentration
!                                    PGMES     = mesophyll conductance (m s-1)
!                                    PGC       = cuticular conductance (m s-1)
!                                    PQDGMES   = Log of Q10 function for mesophyll conductance 
!                                    PT1GMES   = reference temperature for computing 
!                                                compensation concentration function for 
!                                                mesophyll conductance: minimum temperature 
!                                    PT2GMES   = reference temperature for computing 
!                                                compensation concentration function for 
!                                                mesophyll conductance: maximum temperature
!                                    PAMAX     = leaf photosynthetic capacity (kg m-2 s-1)
!                                    PQDAMAX   = Log of Q10 function for leaf photosynthetic capacity
!                                    PT1AMAX   = reference temperature for computing 
!                                                compensation concentration function for leaf 
!                                                photosynthetic capacity: minimum temperature
!                                    PT2AMAX   = reference temperature for computing 
!                                                compensation concentration function for leaf 
!                                                photosynthetic capacity: maximum temperature
!
REAL, DIMENSION(:,:), INTENT(IN)    :: piacan_sunlit, piacan_shade, pfrac_sun
!
REAL, DIMENSION(:), INTENT(IN)      :: pffv ! Floodplain fraction over vegetation
!
REAL, DIMENSION(:,:), INTENT(INOUT) :: piacan ! PAR in the canopy at different gauss level
!
REAL,DIMENSION(:),  INTENT(INOUT) :: pabc, pan, panday, prs, panfm, pgpp
!                                    PABC  = Carrer radiative transfer: normalized heigh of considered layer (bottom=0, top=1)
!                                            Calvet radiative transfer: abcissa of the 3-points Gaussian quadrature 
!                                                (Goudriaan, Agric&For.Meteor, 38,1986)                                            
!                                    PAN   = Net assimilation of CO2 (kg_CO2/kg_air * m/s)
!                                    PANDAY= cumulated daily net assimilation of CO2 (kgCO2/m2/day)
!                                    PRS   = stomatal resistance
!                                    PANFM = maximum leaf assimilation
!                                    PGPP  = Gross Primary Production (kg_CO2/kg_air * m/s)
!
REAL,DIMENSION(:),    INTENT(OUT) :: panf
!                                    PANF  = total assimilation over canopy
!
REAL,DIMENSION(:),    INTENT(OUT) :: presp_leaf
!                                    PRESP_LEAF = dark respiration over canopy
!
!*      0.2    declarations of local variables
!
REAL, PARAMETER                :: zrs_min     = 1.e-4  ! minimum canopy resistance (s m-1)
!
INTEGER                     :: jint, jj ! index for loops
!
REAL, DIMENSION(SIZE(PLAI)) :: zconve1, ztspc, zia
!                                 ZTSPC = temperature conversion (K to C) 
!                                 ZIA   = absorbed PAR
REAL, DIMENSION(SIZE(PLAI)) :: zlai, zgmest, zfzero, zdmax
!                                 ZLAI = LAI 
!                                 ZFZERO  = ideal value of F, no photorespiration or 
!                                            saturation deficit
!                                 ZDMAX   = maximum saturation deficit of atmosphere
!                                           tolerate by vegetation
!
REAL, DIMENSION(SIZE(PLAI)) :: zgammt, zdsp, zanmax
!                                 ZGAMMT  = compensation point 
!                                 ZDSP    = saturation deficit of atmosphere 
!                                           verses the leaf surface (with correction)
!
REAL, DIMENSION(SIZE(PLAI)) :: zxmus, ztan, ztgs, zxia, zan0, zgs0, zxtgs, zrdk,zlaitop,ztrdk,zzlai  
!                                           ZXMUS = cosine of solar zenith angle
!                                           ZTAN  = canopy integrated net assimilation 
!                                           ZTGS  = canopy integrated  leaf conductance
!                                           ZXIA  = incident radiation after diffusion
!                                           ZAN0  = net asimilation at each interval
!                                                   in the canopy
!                                           ZGS0  = leaf conductance at each interval
!                                                   in the canopy        
!                                           ZXTGS = total canopy conductance
!                                           ZRDK  = dark respiration
!                                           ZLAITOP = LAI (thickness of canopy) above considered layer 
!                                           ZTRDK = canopy integrated dark respiration
!                                           ZZLAI = LAI, used for dark respiration profile
!
REAL, DIMENSION(SIZE(PLAI)) :: zan0_,zgs0_,zrdk_ ! parameters for shaded leaves
!
REAL, DIMENSION(SIZE(PLAI)) :: zepso
!                                           ZEPSO conversion of PEPSO in kgCO2/kgair m/s
!
REAL, DIMENSION(SIZE(PLAI)) :: zdmaxstar, zfzerostar, zfzeron, zgmestn  
!                                 ZDMAXSTAR  = maximum saturation deficit of atmosphere
!                                              tolerate by vegetation without soil water stress
!                                 ZFZEROSTAR = initial optimal ratio Ci/Cs for woody vegetation
!                                 ZFZERON    = minimum value for "fzero" in defensive woody strategy
!                                 ZGMESTN    = gmest value at zf2=zf2i in offensive woody strategy
!
!
REAL :: zabc, zweight
!                                           ZABC    = abscissa needed for integration
!                                                     of net assimilation and stomatal 
!                                                     conductance over canopy depth 
!                                                     (working scalar)
!
REAL, DIMENSION(SIZE(PLAI))    :: zwork !Work array
!
LOGICAL, DIMENSION(SIZE(PLAI)) :: gherb, gwood, gf2_inf_f2i, gtrop
!
INTEGER, DIMENSION(1)          :: idmax
!
REAL(KIND=JPRB) :: zhook_handle
!-------------------------------------------------------------------------------
!
! STOMATAL RESISTANCE: ENTRY VARIABLES TO CO2 ROUTINE:
!   CS        = CO2 concentration (kgCO2 kgair-1) cs
!   DSP       = specific humidity deficit (kgH2O kgair-1) ds
!   TSM       = surface temperature (C) ts
!   RG        = global radiation (W m-2) rg
!
! initialisation: convert from ppm to mg/m-3
!
IF (lhook) CALL dr_hook('COTWORES',0,zhook_handle)
!
zconve1(:) = xmco2*prhoa(:)/xmd
!
! initialisation: convert from K to C
!
ztspc(:)  = ptg(:) - xtt               
!
zlai(:)   = plai(:)
zgmest(:) = pgmes(:)
zfzero(:) = pfzero(:)
!
!GTROP = linear stress in case of tropical evergreen forest 
!        (with fixed f0=0.74 for Carrer rad. transf., f0=0.7 with Calvet rad. transf.)
gtrop(:) = (pvegtype(:,nvt_trbe) > 0.8) 
!
gherb(:) = (pvegtype(:,nvt_tebd) + pvegtype(:,nvt_trbe) + pvegtype(:,nvt_bone)   &
           +pvegtype(:,nvt_trbd) + pvegtype(:,nvt_tebe) + pvegtype(:,nvt_tene)   & 
           +pvegtype(:,nvt_bobd) + pvegtype(:,nvt_bond) + pvegtype(:,nvt_shrb)<0.5)
gwood(:) = (.NOT.gherb(:))
!
IF (hphoto=='AGS' .OR. hphoto=='LAI') THEN
  !
  !  Compute conductance and assimilation of CO2: 
  !
  zdmax = xdmax_ags
  !
  ! Add soil moisture stress effect to leaf conductance:
  !
  zgmest(:) = zgmest(:) * pf2(:)
  !
ELSEIF (hphoto=='AST' .OR. hphoto=='LST' .OR. hphoto=='NIT' .OR. hphoto=='NCB') THEN
  !
  WHERE (plai(:)==xundef) zlai(:)=0.0
  !
  !    See (Varlet-Granchet C., M. Chartier, G. Gosse,  and R. Bonhomme, 1981: 
  !    Rayonnement utilise pour la photosynthese des vegetaux en
  !    conditions naturelles: caracterisation et variations. 
  !    Oecol. Plant. 2(16), 189-202.)
  !
  !-------------------------------------
  ! Add soil moisture stress effect to leaf conductance:
  ! OFFENSIVE and DEFENSIVE water stress response 
  !  
  zdmax(:)  = pdmax(:)
  !
  gf2_inf_f2i(:) = (pf2(:)<pf2i(:))
  !
  ! -HERBACEOUS-
  !
  WHERE(gherb(:).AND.ostressdef(:))
    zdmax(:) = xdmaxn
  ENDWHERE
  WHERE(gherb(:).AND..NOT.ostressdef(:))
    zdmax(:) = xdmaxx
  ENDWHERE
  !
  ! PAH and PBH are original coefficients of Calvet 2000
  WHERE(gherb(:).AND.(.NOT.gf2_inf_f2i(:)))
    zdmaxstar(:) = exp((log(zgmest(:)*1000.)-pah(:))/pbh(:))/1000.
    zdmax(:) = zdmaxstar(:) - (zdmaxstar(:)-zdmax(:))*(1.-pf2(:))/(1.-pf2i(:))
  ENDWHERE
  !
  WHERE(gherb(:))
        zgmest(:) = exp(pah(:)+pbh(:)*log(zdmax(:)*1000.))/1000.
  ENDWHERE
  !
  WHERE (gherb(:).AND.gf2_inf_f2i(:).AND.ostressdef(:))
      zgmest(:) = zgmest(:) * pf2(:)/pf2i(:)
  ENDWHERE
  WHERE(gherb(:).AND.gf2_inf_f2i(:).AND.(.NOT.ostressdef(:)))
      zdmax(:) = zdmax(:) * pf2(:)/pf2i(:)
  ENDWHERE
  !
  ! to limit photosynthesis under wilting point
  WHERE (gherb(:).AND.(.NOT.ostressdef(:)).AND.zdmax(:)<=xdmaxn)
    zdmax(:)  = xdmaxn
    zgmest(:) = (exp(pah(:)+pbh(:)*log(xdmaxn*1000.))/1000.)*pf2(:)/pf2i(:)
  ENDWHERE
  !
  ! -WOODY but not tropical forest-
  !
  WHERE(gwood(:))
    zfzerostar(:) = ( xaw  - log(zgmest(:)*1000.) )/xbw
  ENDWHERE
  !
  WHERE (gwood(:).AND.ostressdef(:))
    zgmestn(:) = zgmest(:)
  ENDWHERE
  WHERE(gwood(:).AND.(.NOT.ostressdef(:)))
    zgmestn(:) = exp(xasw - xbw*zfzerostar(:))/1000.
  ENDWHERE
  !
  WHERE (gwood(:).AND.gf2_inf_f2i(:)) 
    zgmestn(:) = zgmestn(:)*pf2(:)/pf2i(:)
  ENDWHERE
  !
  WHERE(gwood(:))
    zwork(:) = max( xdenom_min, zgmestn(:) )
    zfzeron(:) = (xasw - log(zwork(:)*1000.))/xbw
  ENDWHERE
  !
  WHERE(gwood(:).AND.(.NOT.gf2_inf_f2i(:)).AND.ostressdef(:))
    zfzero(:) = zfzerostar(:)
    zfzero(:) = zfzero(:) - (zfzero(:)-zfzeron(:))*(1.-pf2(:))/(1.-pf2i(:))  
  ENDWHERE    
  WHERE(gwood(:).AND.(.NOT.gf2_inf_f2i(:)).AND.(.NOT.ostressdef(:)))
    zfzero(:) = zfzerostar(:)
    zgmest(:) = zgmest(:) - (zgmest(:)-zgmestn(:))*(1.-pf2(:))/(1.-pf2i(:))  
  ENDWHERE    
  !
  WHERE(gwood(:).AND.gf2_inf_f2i(:))
    zfzero(:) = min(.95, zfzeron(:))
    zgmest(:) = zgmestn(:)
  ENDWHERE  
  !
  ! -Tropical Forest-
  !
  WHERE(gtrop(:))
   zfzero(:) = pfzero(:)
   zgmest(:) = pgmes(:)*pf2(:)
  ENDWHERE
  !
ENDIF
!
!-------------------------
!
! compensation point (ppm): temperature response
!
!before optimization (with non log PQDGAMM) : 
!ZGAMMT(:) = PGAMM(:)*PQDGAMM(:)**(0.1*(ZTSPC(:)-25.0))
zwork(:) = (0.1*(ztspc(:)-25.0)) * pqdgamm(:)
zgammt(:) = pgamm(:) * exp(zwork(:))
!
! specific humidity deficit (kg kg-1)
!
zdsp(:)   = max( 0.0, pqsat(:) - pqa(:) - ple(:)*pra(:)/(prhoa*xlvtt) )
!
! cosine of solar zenith angle 
!
zxmus(:) = max(cos(pzenith(:)),0.01)
!
!
! Compute temperature response functions:
!
! kg/m2/s
!before optimization (with non log PQDAMAX) : 
!ZANMAX(:) = ( PAMAX(:)*PQDAMAX(:)**(0.1*(ZTSPC(:)-25.0)) ) / ...
zwork(:) = (0.1*(ztspc(:)-25.0)) * pqdamax(:)
zanmax(:) = ( pamax(:) * exp(zwork(:))  ) &
          / ( (1.0+exp(0.3*(pt1amax(:)-ztspc(:))))* (1.0+exp(0.3*(ztspc(:)-pt2amax(:)))) )
!
! m/s
!before optimization (with non log PQDGMES) : 
!ZGMEST(:) = ( ZGMEST(:)*PQDGMES(:)**(0.1*(ZTSPC(:)-25.0)) ) / ...
zwork(:) = (0.1*(ztspc(:)-25.0)) * pqdgmes(:)
zgmest(:) = ( zgmest(:) * exp(zwork(:)) ) &
          / ( (1.0+exp(0.3*(pt1gmes(:)-ztspc(:))))*  (1.0+exp(0.3*(ztspc(:)-pt2gmes(:)))) )  
!
!
! Integration over the canopy: SIZE(PABC) increments
! are used to approximate the integral.
!
ztan(:) = 0.0
ztgs(:) = 0.0
ztrdk(:)= 0.0
!
! Unit conversion
! ZANMAX and ZEPSO from kgCO2/m2/s to kgCO2/kgair m/s by dividing by RHOA (kgair/m3)
! ZGAMMT from ppm to kgCO2/kgair
zgammt(:)  = zgammt(:) * xmco2 / xmd * 1e-6
zanmax(:) = zanmax(:) / prhoa(:)
zepso(:)  = pepso(:)  / prhoa(:)
!
zia(:)     = psw_rad(:)*xparcf
!
DO jint = 1, SIZE(pabc)
  !
  !  Diffusion of incident radiation:
  !
  IF (otr_ml) THEN
    !
    zabc = 1.
    IF (jint.LT.SIZE(pabc)) zabc = pabc(jint+1)
    zweight = zabc - pabc(jint)
    zxia(:) = piacan_sunlit(:,jint)
    !
  ELSE
    !
    zabc = pabc(jint)
    zweight = ppoi(jint)
    !
    CALL ccetr(zxia,zia,zxmus,zabc,zlai)
    !
    ! PAR at different Gauss  level in micmolphot/m2/s
    !
    piacan(:,jint)= zxia(:)
    !
  ENDIF
  !
  ! Compute conductance and assimilation of CO2: 
  !
  !Extinction of respiration depends on LAI above only for tropical evergreen forest
  zlaitop(:) = 0.
  zzlai(:) = 1.
  IF (otr_ml) THEN         
    WHERE(gwood(:))  
      zlaitop(:) = (1.-(pabc(jint)+zabc)/2.)*zlai(:)
      zzlai(:) = zlai(:)
    ENDWHERE
  ENDIF
  CALL cotwo(pcsp, pf2, zxia, zdsp, zgammt,             &
             zfzero, zepso, zanmax, zgmest, pgc, zdmax, &  
             zan0, zgs0, zrdk, zlaitop, zzlai           )
  !
  IF (otr_ml) THEN
    !
    zxia(:) = piacan_shade(:,jint)
    CALL cotwo(pcsp, pf2, zxia, zdsp, zgammt,             &
               zfzero, zepso, zanmax, zgmest, pgc, zdmax, &  
               zan0_, zgs0_, zrdk_, zlaitop, zzlai        )
    !
    WHERE (oshade(:))
      !ponderate sum.
      zan0(:)=pfrac_sun(:,jint)*zan0(:)+(1.-pfrac_sun(:,jint))*zan0_(:)
      zrdk(:)=pfrac_sun(:,jint)*zrdk(:)+(1.-pfrac_sun(:,jint))*zrdk_(:)
      zgs0(:)=pfrac_sun(:,jint)*zgs0(:)+(1.-pfrac_sun(:,jint))*zgs0_(:)
    ENDWHERE
    !
  ENDIF
  !
  ! kgCO2/kgair m/s
  ztan(:) = ztan(:) + zan0(:)*zweight
  ztgs(:) = ztgs(:) + zgs0(:)*zweight
  ztrdk(:) = ztrdk(:) + zrdk(:)*zweight
  !
END DO
!
!
! Total assimilation
!
panf(:)= ztan(:)
!
! Net assimilation over canopy
!
pan(:) = (1.0-pdelta(:))*(1.0-ppsnv(:)-pffv(:))*panf(:)*zlai(:)
!
! Dark respiration over canopy (does not depend on radiation, 
! no need to integrate over vertical dimension)
!
presp_leaf(:) = (1.0-pdelta(:))*(1.0-ppsnv(:)-pffv(:))*ztrdk(:)*zlai(:)
!
! Gross primary production over canopy
!
pgpp(:) = pan(:) + presp_leaf(:)
!
! Cumulated daily net assimilation over canopy (kgCO2/m2/day)
!
panday(:) = panday(:) + pan(:) * ptstep * prhoa(:)
!
! Adjust maximum leaf assimilation:
!
panfm(:) = max( panf(:), panfm(:) )
!
! Total conductance over canopy 
!
zxtgs(:) = ztgs(:)*zlai(:)
!
! Canopy resistance from Ags:
!
prs(:) = min( 1.0/(zxtgs(:)+xdenom_min), xrs_max)
!
prs(:) = max( prs(:), zrs_min)
!
IF (lhook) CALL dr_hook('COTWORES',1,zhook_handle)
!
END SUBROUTINE cotwores