!*****************************************************************!
!* *!
!* 4C (FORESEE) Simulation Model *!
!* *!
!* *!
!* Subroutine *!
!* target thinning - *!
!* thinning routine with given values of stem number per *!
!* thinning year as target values *!
!* rtargetm i given inN/ha *!
! or target value 0 < rtargetm <1 for *!
!* basal area reduction *!
!* *!
!* Copyright (C) 1996-2018 *!
!* Potsdam Institute for Climate Impact Reserach (PIK) *!
!* Authors and contributors see AUTHOR file *!
!* This file is part of 4C and is licensed under BSD-2-Clause *!
!* See LICENSE file or under: *!
!* http://www.https://opensource.org/licenses/BSD-2-Clause *!
!* Contact: *!
!* https://gitlab.pik-potsdam.de/foresee/4C *!
!* *!
!*****************************************************************!
SUBROUTINE target_thinning_bas(i)
use data_stand
use data_manag
use data_simul
use data_species
use data_par
implicit none
real :: rtargetm=0. ! target value of stem biomass
integer :: target_help = 0.
real :: dbhmin=0, &
dbhmin_us = 0, &
wpa=0, & ! Weibull parameter
wpb=0, & ! -"-
wpc=0, & ! -"-
d63=0, &
help=0, &
pequal, &
tdbh=0, &
bas_area=0, &
bas_help=0., &
bas_area_us = 0., &
bas_area_test=0., &
target_help1=0, &
dbh_h =0, &
db_l = 0., &
db_u = 0., &
d_est=0., &
w_kb=0., &
stembiom, &
stembiom_us, &
stembiom_re, &
stembiom_all, &
diff, &
mdiam, &
mdiam_us
integer :: nrmin, &
nrmin_us, &
lowtree, &
undertree, &
flagth, &
taxnr, &
counth, &
min_id, &
max_id, &
ih1,ih2,ncoh, &
coun1
integer :: flag_bas
! auxilarity for thinning routine 4: selective thinning
integer :: count, i, &
idum,third, ipot, isel, ih
integer,dimension(0:anz_coh) :: cohl
integer, dimension(anz_coh) :: id_pot
real :: h1, h2 , tar_h
real,external :: gasdev
real:: ran0
! reacalculation of target to kg DW/patch
flag_bas=0
h1 = 0.
h2 = 0.
count = 0
cohl = -1
flagth = 0
coun1 = 0
help=0.
lowtree=0
undertree = 0
anz_tree_dbh = 0
bas_area = 0.
bas_area_us = 0.
! stem biomass of overstorey
stembiom = 0.
! stem biomass of understorey
stembiom_us = 0.
stembiom_all = 0.
tar_h = 300.
if (time.eq.73.and. ip .eq.87) then
stembiom = 0
end if
taxnr = thin_spec(i)
mdiam = 0.
mdiam_us = 0.
! calculation of mean diameter (correspondung to med_diam) and basal area of stand
zeig => pt%first
DO
IF (.NOT. ASSOCIATED(zeig)) EXIT
! Modification for V Kint: no test for diameter
IF((zeig%coh%ntreea>0).and.zeig%coh%species.eq.taxnr.and.zeig%coh%underst.eq.0) THEN
! overstorey
stembiom = stembiom + (zeig%coh%x_sap + zeig%coh%x_hrt)*zeig%coh%ntreea
help = help + zeig%coh%ntreea*(zeig%coh%diam**2)
bas_area = bas_area + zeig%coh%ntreea*(zeig%coh%diam**2)*pi/4
if( zeig%coh%diam>0) then
anz_tree_dbh = anz_tree_dbh + zeig%coh%ntreea
mdiam = mdiam + zeig%coh%ntreea * (zeig%coh%diam**2)
end if
! Trees with DBH=0 for populations and per species; Baeume mit DBH =0 fuer Bestand und pro Spezie
ELSE IF( (zeig%coh%ntreea>0).and.zeig%coh%species.eq.taxnr.and.zeig%coh%underst.eq.1) THEN
! seedings/regeneration
stembiom_re = stembiom_re + (zeig%coh%x_sap + zeig%coh%x_hrt)*zeig%coh%ntreea
lowtree = lowtree + zeig%coh%ntreea
ELSE if((zeig%coh%ntreea>0).and.zeig%coh%species.eq.taxnr.and.zeig%coh%underst.eq.2.and.zeig%coh%underst.eq.2) THEN
! understorey
stembiom_us = stembiom_us + (zeig%coh%x_sap + zeig%coh%x_hrt)*zeig%coh%ntreea
mdiam_us = mdiam_us + zeig%coh%ntreea * (zeig%coh%diam**2)
undertree = undertree + zeig%coh%ntreea
bas_area_us = bas_area_us + zeig%coh%ntreea*(zeig%coh%diam**2)*pi/4
ENDIF
zeig => zeig%next
ENDDO
! mean diameter for over and understorey
stembiom_all = stembiom + stembiom_us
if(anz_tree_dbh.ne.0) mdiam = sqrt(mdiam/real(anz_tree_dbh))
if(undertree.ne.0) mdiam_us = sqrt(mdiam_us/undertree)
third = nint(anz_tree_dbh*0.333333)
anz_tree_ha = nint(anz_tree_dbh*10000./kpatchsize)
anz_tree = anz_tree_dbh + undertree
IF(anz_tree>0)THEN
if(lowtree<anz_tree) help = sqrt(help/(anz_tree-lowtree))
ENDIF
! new basal area management
if(target_mass(i).le.1) then
flag_bas=1
if(thin_stor(i).eq.0) then
target_help = bas_area
else
target_help = bas_area_us
end if
else
if(thin_stor(i).eq.0) then
target_help = anz_tree_dbh
else
target_help = undertree
end if
end if
! tending
if(thin_tysp(i).eq.4) then
target_help = bas_area_us
end if
if(target_mass(i).gt.1) then
rtargetm = target_mass(i)*kpatchsize/10000.
else
rtargetm = target_mass(i)
end if
! cuttting
if (rtargetm.eq.0.)then
zeig => pt%first
DO
IF (.NOT. ASSOCIATED(zeig)) EXIT
if(zeig%coh%species.eq.taxnr.and. zeig%coh%underst.eq.thin_stor(i)) then
zeig%coh%ntreem = zeig%coh%ntreea
zeig%coh%ntreea = 0
zeig%coh%nta = 0
end if
zeig=> zeig%next
END DO
!tending of regeneration
else if(thin_tysp(i).eq.4) then
rtargetm = target_mass(i) * target_help
min_id = 1000
max_id = 0.
zeig=>pt%first
do
if(.not.associated(zeig)) exit
if(zeig%coh%species.eq.taxnr.and. zeig%coh%underst.eq.1 .or. zeig%coh%underst.eq.2) then
ih1 = zeig%coh%ident
if(ih1.lt.min_id) min_id = ih1
ih2 = zeig%coh%ident
if (ih2.gt.max_id) max_id = ih2
end if
zeig=> zeig%next
end do
target_help1 = 0.
do
call random_number(pequal)
ncoh = min_id +(max_id-min_id)*pequal
zeig=>pt%first
do
if(.not.associated(zeig)) exit
if(zeig%coh%species.eq.taxnr.and. zeig%coh%underst.eq.1 .or.zeig%coh%underst.eq.2 .and. zeig%coh%ident.eq.ncoh ) then
zeig%coh%ntreea = zeig%coh%ntreea - 1
zeig%coh%nta = zeig%coh%nta-1
zeig%coh%ntreem = zeig%coh%ntreem +1
if (flag_bas.eq.1) then
target_help=target_help - zeig%coh%diam*zeig%coh%diam*pi/4
else
target_help = target_help - zeig%coh%ntreea
end if
exit
end if
zeig=>zeig%next
end do
diff = target_help - rtargetm
if(diff.lt.0.01) exit
end do
! different thinnings from below and above
else IF ( rtargetm .ne. 0.) then
if(target_mass(i).lt.1.) then
rtargetm = target_mass(i) * target_help
! No management if rtargetm=1
else if (rtargetm.eq.1) then
return
endif
select case(thin_tysp(i))
case(1)
! medium lower thinning
d_est = 1.02
w_kb = 2.5
case(2)
! strong lower thinning
d_est = 1.03
w_kb = 1.5
case(3)
! High thinning
d_est = 1.04
w_kb = 1.2
end select
! calculation of Weibull-Parameter
call min_dbh_overs(nrmin,dbhmin,taxnr)
call min_dbh_unders(nrmin_us,dbhmin_us, taxnr)
! write (7777,*) time, nrmin,nrmin_us, dbhmin, dbhmin_us
bas_help = bas_area
if (thin_stor(i).eq.0) then
wpa = dbhmin
else
wpa = dbhmin_us
end if
if (thin_stor(i).eq.0) then
d63 = mdiam*d_est
else
d63 = mdiam_us * d_est
end if
wpb = (d63 - wpa)/ w_kb
wpc = 2
wpc = 0.8
if ((thin_tysp(i).ne.4) .and. rtargetm.lt.target_help) then
!selection of trees for thinning
DO ! 11 begin selecting
IF (thin_stor(i) .eq. 2 .and. nrmin_us .eq. -1) EXIT ! exit thinning loop if there are no trees with dbh in understory an understory shall be thinned
coun1=coun1+1
call random_number(pequal)
tdbh = wpa + wpb*(-log(1.-pequal))**(1./wpc)
flagth = 0
! list of potential thinned tree chorts
counth = 0
id_pot = 0
ipot = 1
zeig => pt%first
DO ! 2 list preparing
IF (.NOT. ASSOCIATED(zeig)) EXIT
if(zeig%coh%notviable.eqv. .TRUE.) then
if(flag_mort.eq.0) then
id_pot(ipot)=zeig%coh%ident
ipot=ipot + 1
endif
else if(zeig%coh%ntreea.gt.0.and.zeig%coh%species.eq.taxnr.and. zeig%coh%underst.eq.thin_stor(i)) then
dbh_h = zeig%coh%diam
db_l = dbh_h - 0.1*dbh_h
db_u = dbh_h + 0.1*dbh_h
counth = counth +1
if (tdbh.ge.db_l.and.tdbh.le.db_u.and. zeig%coh%ntreea.ne. 0) then
id_pot(ipot) = zeig%coh%ident
ipot = ipot + 1
end if
if(counth.gt. 100000) exit
end if
zeig=> zeig%next
END DO ! 2 list of potential thinned tree cohorts
! selecting one equal distributed tree from the list of cohorts
if ((ipot-1).ge.1) then
if((ipot-1).eq.1) then
isel = 1
else
call random_number(pequal)
pequal = ran0(idum)
isel = int(pequal*(ipot-1)) +1
end if
ih = id_pot(isel)
zeig => pt%first
DO ! 3 removing one tree
IF (.NOT. ASSOCIATED(zeig)) EXIT
if(zeig%coh%ident.eq. ih.and.zeig%coh%ntreea.ne. 0 ) then
if(zeig%coh%notviable.eqv..TRUE.) then
if(flag_mort.eq.0) then
zeig%coh%ntreem = zeig%coh%ntreea
zeig%coh%ntreea=0
zeig%coh%nta=0
if (flag_bas.eq.1) then
target_help = target_help - (zeig%coh%diam**2)*pi/4
else
target_help = target_help -1
endif
end if
else
zeig%coh%ntreea = zeig%coh%ntreea -1
zeig%coh%nta = zeig%coh%nta -1
zeig%coh%ntreem = zeig%coh%ntreem +1
if(flag_bas.eq.1) then
! write (8888,*) 'target_help', time, target_help, rtargetm
target_help = target_help - (zeig%coh%diam**2)*pi/4
else
target_help = target_help -1
end if
exit
! end if
end if
end if
zeig =>zeig%next
ENDDO ! 3 thinning of one tree
end if
diff = target_help - rtargetm
if(diff.le.0.1) exit
if(coun1.gt.100000) exit
ENDDO ! 11 total thinning
end if ! thintype 1,2,3
END IF ! all thinnings and tending
zeig=>pt%first
do
if(.not.associated(zeig)) exit
if(zeig%coh%ntreea>0.and.zeig%coh%species.eq.taxnr) then
bas_area_test = bas_area_test + zeig%coh%ntreea*(zeig%coh%diam**2)*pi/4
end if
zeig=>zeig%next
end do
! adding biomasses to litter pools depending on stage of stand
zeig=>pt%first
do
if(.not.associated(zeig)) exit
if(zeig%coh%ntreem>0.and.zeig%coh%species.eq.taxnr) then
! all parts without stems of trees are input for litter
h1 = zeig%coh%ntreem*zeig%coh%x_fol
h2 = h2 + h1
zeig%coh%litC_fol = zeig%coh%litC_fol + zeig%coh%ntreem*(1.-spar(taxnr)%psf)*zeig%coh%x_fol*cpart
zeig%coh%litN_fol = zeig%coh%litN_fol + zeig%coh%ntreem*((1.-spar(taxnr)%psf)*zeig%coh%x_fol*cpart)/spar(taxnr)%cnr_fol
zeig%coh%litC_frt = zeig%coh%litC_frt + zeig%coh%ntreem*zeig%coh%x_frt*cpart
zeig%coh%litN_frt = zeig%coh%litN_frt + zeig%coh%ntreem*zeig%coh%x_frt*cpart/spar(taxnr)%cnr_frt
zeig%coh%litC_tb = zeig%coh%litC_tb + zeig%coh%ntreem*zeig%coh%x_tb*cpart
zeig%coh%litN_tb = zeig%coh%litN_tb + zeig%coh%ntreem*zeig%coh%x_tb*cpart/spar(taxnr)%cnr_tbc
zeig%coh%litC_crt = zeig%coh%litC_crt + zeig%coh%ntreem*zeig%coh%x_crt*cpart
zeig%coh%litN_crt = zeig%coh%litN_crt + zeig%coh%ntreem*zeig%coh%x_crt*cpart/spar(taxnr)%cnr_crt
endif
zeig=>zeig%next
enddo
call class_man
END SUBROUTINE target_thinning_bas