LfCalc.f90

Go to the documentation of this file.

!==============================================================================
!
!  Name   :  Lfcalc
!
!  Purpose:
!
!  Remarks:  contains the SUBROUTINE
!
!            Lfcalc (<LocalDynOneTimeStep<LocalDynOneTimeStep<DistrLocalDynThisTime<
!                      TimeLoop<TreeMig)
!
!
!==============================================================================
! design          : H. Lischke, N. Zimmermann
! author(s)       : H. Lischke, N. Zimmermann
! implementation  : H. Lischke, N. Zimmermann
! cleaner         : T.J. Loeffler
! copyright       : (c) 1999, 03 by H. Lischke
!==============================================================================
!===================================================================== 
!===============================================================
 
SUBROUTINE lfcalc(nspc, heightcl, LAIsum, LAIsqsum, numin0)
 
    ! <CALL modules for TYPE definitions and +- fixed variables>---------------------------------------------------------------
    use all_par, only: maxspc, maxhc, maxlc, pltsiz, facvar, &
                       ltfrq, ltcum, ltprf, la, &
                       speccl                                ! Use control, stand and species parameters
    IMPLICIT NONE
 
    ! <Passed variables>--------------------------------------------------------------------------------
    INTEGER, INTENT(in) :: heightcl, &                 ! height  class                            [=height]
                           nspc                      ! number of species                         [=spec]
    REAL, INTENT(in)    :: numin0(0:maxhc, 1:maxspc)  ! temp var originating from numNInH           [=..]
    REAL, INTENT(inout) :: LAIsqsum, &                ! sq. sum of LA                        [=LAIsqrSum]
                           LAIsum                     ! leave area (LA)                         [=LAISum]
 
    ! <Local variables>---------------------------------------------------------------------------------
    INTEGER :: lightcl, &                              ! lightclass index variable                    [=l]
               i                                     ! looping variable  for species
    REAL    :: ntinh, &                      ! # tree in height class                    [=nInH]
               newla, newla2, &                      ! new lai  and squared lai          [=rNewLeafArea]
               treelai, &                      ! lai  per tree                     [=rNewLeafArea]
               sigma, &                      ! std?                                     [=sigma]
               sum, &                      ! running sum                                [=sum]
               disbot(maxlc), &                      ! Distr. of light (bottom of class)?    [=distrBot]
               distop, &                      ! Distr. of light (top of class)?       [=distrTop]
               value                                 ! Out-variable from NormDist
 
    ! <Here the SUBROUTINE starts>**********************************************************************
    !----- By species, the number of trees is calculated and evaluated per height class
    byspecies: do i = 1, nspc
        !----- The current height class is completely included, i.e. there is also competition between trees in this height class
        ntinh = numin0(heightcl, i)
        !----- The leaf area and the squared leaf area are calculated
        if (ntinh /= 0.0) then
            treelai = speccl(i)%lain(heightcl)/pltsiz
            newla = ntinh*treelai
            newla2 = ntinh*treelai**2
            laisum = laisum + newla
            laisqsum = laisqsum + newla2
        end if
    end do byspecies
 
    !----- LAIsum is the mean, sigma the standard deviation of the leaf area Normal distribution
    !----- facvar describes how narrow the Normal distributions is.
    !----- 0: uniform spatial distribution, 1: random spat. distribution, >1 clumped spatial distribution.
    !----- It is parameterized in All_PAR
 
    if ((facvar <= 0.0) .or. (laisqsum  <  1.e-7)) then
        sigma = 0.001         !----- sigma has to be larger 0
    else
        sigma = sqrt(laisqsum)*facvar
    end if
 
    ltfrq = 0.0
    !----- if there are (nearly) no leaves, then there is full light...
    if (abs(laisum)  <  1.e-10) then
        ltfrq(maxlc) = 1.0
        sum = 1.0
    else
        !----- ... otherwise for each light class from full light  (lightcl=maxlc ) to complete darkness (lightcl=1)
        !----- the (cumulative) normal distribution function with parameters LAIsum and sigma is evaluated
        !----- at the leaf area (la(light)) required to produce the light value. The difference between
        !----- two values of the distribution function is the frequency of these leaf areas, resp. light values.
        sum = 0.0
        call normdist(value, laisum, sigma, la(1))
        ltfrq(1) = 1 - value
        sum = sum + 1 - value
        call normdist(value, laisum, sigma, la(maxlc - 1))
        ltfrq(maxlc) = value
        sum = sum + value
        distop = value
 
        byltcls1: do lightcl = maxlc - 1, 2, -1
            call normdist(value, laisum, sigma, la(lightcl - 1))
            disbot(lightcl) = value
            ltfrq(lightcl) = disbot(lightcl) - distop
            sum = sum + ltfrq(lightcl)
            distop = disbot(lightcl)
        end do byltcls1
    end if
 
    !----- If the leaf area is so large that even the la value set for the lowest light class is much lower,
    !----- and all lftrq's are 0, then the lowest light class frequency is set to 1
    if (abs(sum) < 1.e-20) then
        ltfrq(1) = 1.0
        sum = 1.0
    end if
 
    !----- The light distribution is normalized, because it is truncated a 0
    ltfrq(maxlc) = ltfrq(maxlc)/sum
    ltprf(heightcl, maxlc) = ltfrq(maxlc)
    ltcum(maxlc) = ltfrq(maxlc)
 
    byltcls2: do lightcl = maxlc - 1, 1, -1
        ltfrq(lightcl) = ltfrq(lightcl)/sum
        ltprf(heightcl, lightcl) = ltfrq(lightcl)
        ltcum(lightcl) = ltcum(lightcl + 1) + ltfrq(lightcl)
    end do byltcls2
 
end SUBROUTINE lfcalc
!done