!*****************************************************************!
!* *!
!* SEA for 4C (FORESEE) Simulation Model *!
!* *!
!* *!
!* Subroutines: *!
!* *!
!* sea: control subroutine for sea *!
!* *!
!* sort_mansort: first sorting of mansort *!
!* sort_standsort: first sorting of standsort *!
!* sort_industrial *!
!* calculate_harvest_costs *!
!* calculate_assets *!
!* calculate_costs *!
!* calculate_npv *!
!* read_sea_prices *!
!* *!
!* 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 sea
use data_simul
use data_wpm
implicit none
character(150) pricesFile
! begin program
call setFlags
call allocate_in_output
call ini_input_sea
! read prices
pricesFile = trim(dirin)//'sea_prices.wpm'
! call fullPath(pricesFile, dirin)
call read_sea_prices(pricesFile)
! simulation
if ( associated(first_mansort) ) then
! first sorting
call sort_mansort
end if
if ( associated(first_standsort) ) then
call sort_standsort
end if
! sort 0.4 to the industrial wood
call sort_industrial
! harvest costs calculation
call calculate_harvest_costs
! timber selling assets calculation
call calculate_assets
! calculate rest costs
call calculate_costs
! calculate npv
call calculate_npv
end subroutine sea
!***************************************************************
! calculate timber grades from the mansort and standsort input
! input: data_wpm
! output: data_wpm
!***************************************************************
subroutine sort_mansort
use data_wpm
use data_simul
implicit none
integer i, index
real volume, pi, diam
character(4) act_typus
integer act_spec, act_year, set_year
pi = 3.1415926536 ! PI
i = nr_years
! set the first year, set an ima
act_year = first_mansort%mansort%year
set_year = first_manrec%manrec%year
! set the run pointer to the begin of the list
act_mansort => first_mansort
act_manrec => first_manrec
! check if actuelles management is not brushing or tending
if (trim(act_manrec%manrec%management) .eq. 'brushing') then
if( associated(act_manrec%next) ) then
act_manrec => act_manrec%next
endif
endif
if (trim(act_manrec%manrec%management) .eq. 'tending') then
if( associated(act_manrec%next) ) then
act_manrec => act_manrec%next
endif
endif
! check if last management year was some years ealier
do while (set_year < act_mansort%mansort%year)
act_mansort => act_mansort%next
end do
act_year = act_mansort%mansort%year
set_year = act_manrec%manrec%year
if( associated(act_manrec%next) ) then
act_manrec => act_manrec%next
end if
do while (associated(act_mansort))
if ( act_year <= act_manrec%manrec%year .and. &
act_manrec%manrec%year /= set_year ) then
! check the management in actual year
if (trim(act_manrec%manrec%management) .eq. 'brushing') then
if( associated(act_manrec%next) ) then
act_manrec => act_manrec%next
endif
endif
if (trim(act_manrec%manrec%management) .eq. 'tending') then
if( associated(act_manrec%next) ) then
act_manrec => act_manrec%next
endif
endif
if (trim(act_manrec%manrec%management) .ne. 'tending' .and. &
trim(act_manrec%manrec%management) .ne. 'brushing' ) then
! set next value for actual manrec year
set_year = act_manrec%manrec%year
endif
if( associated(act_manrec%next) ) then
act_manrec => act_manrec%next
end if
endif
! set species index
act_spec = act_mansort%mansort%spec
act_typus = act_mansort%mansort%typus
! calculate carbon for the actual line in mansort
volume = act_mansort%mansort%volume * act_mansort%mansort%number
select case (trim(act_typus))
case ('ste1', 'ste2')
diam = act_mansort%mansort%diam_wob
if (diam >=25 .and. diam < 30) index = 8
if (diam >=30 .and. diam < 35) index = 9
if (diam >=35) index = 10
case ('sg1', 'sg2')
diam = act_mansort%mansort%diam_wob
if (diam < 15) index = 3
if (diam >=15 .and. diam < 20) index = 4
if (diam >=20 .and. diam < 25) index = 5
if (diam >=25 .and. diam < 30) index = 6
if (diam >=30) index = 7
case ('in1', 'in2')
index = 2
case ('fue')
index = 1
end select
mansort_tg(act_spec, index, set_year) = mansort_tg(act_spec, index, set_year) + volume
! after using the mansort list item, go to the next
act_mansort => act_mansort%next
if (associated(act_mansort)) then
act_year = act_mansort%mansort%year
end if
end do
end subroutine sort_mansort
!***************************************************************
subroutine sort_industrial
use data_wpm
use data_simul
implicit none
integer i, j
real ind
! sort value*ind into the industrial wood
ind = 0.4
do i = 1, nr_spec
do j = 1, nr_timb_grades
select case (j>2)
case (.TRUE.)
mansort_tg(i, 2, :) = mansort_tg(i, j, :) * ind + mansort_tg(i, 2, :)
mansort_tg(i, j, :) = mansort_tg(i, j, :) * (1 - ind)
standsort_tg(i, 2, :) = standsort_tg(i, j, :) * ind + standsort_tg(i, 2, :)
standsort_tg(i, j, :) = standsort_tg(i, j, :) * (1 - ind)
end select
end do
end do
end subroutine sort_industrial
!****************************************************************
subroutine sort_standsort
use data_wpm
implicit none
integer j, index
real volume, pi, diam
character(4) act_typus
integer act_spec, act_year
pi = 3.1415926536 ! PI
j = 1
! set the first year
act_year = first_standsort%mansort%year
! set the run pointer to the begin of the list
act_standsort => first_standsort
! check if last management year was some years ealier
do while (act_year < act_standsort%mansort%year)
act_standsort => act_standsort%next
end do
do while (associated(act_standsort))
! check the management year
if ( act_year <= act_standsort%mansort%year )then
! set next value for actual standsort year
j = j+1
endif
! set species index
act_spec = act_standsort%mansort%spec
act_typus = act_standsort%mansort%typus
! calculate carbon for the actual line in standsort
volume = act_standsort%mansort%volume * act_standsort%mansort%number
select case (trim(act_typus))
case ('ste1', 'ste2')
diam = act_standsort%mansort%diam_wob
if (diam >=25 .and. diam < 30) index = 8
if (diam >=30 .and. diam < 35) index = 9
if (diam >=35) index = 10
case ('sg1', 'sg2')
diam = act_standsort%mansort%diam_wob
if (diam < 15) index = 3
if (diam >=15 .and. diam < 20) index = 4
if (diam >=20 .and. diam < 25) index = 5
if (diam >=25 .and. diam < 30) index = 6
if (diam >=30) index = 7
case ('in1', 'in2')
index = 2
case ('fue')
index = 1
end select
standsort_tg(act_spec, index, act_standsort%mansort%year) = standsort_tg(act_spec, index, act_standsort%mansort%year) + volume
! after using the standsort list item, go to the next
act_standsort => act_standsort%next
if (associated(act_standsort)) then
act_year = act_standsort%mansort%year
end if
end do
end subroutine sort_standsort
!*****************************************************************************
subroutine calculate_harvest_costs()
use data_wpm
use data_simul
implicit none
integer i, j
! calcultation of costs only implemented for monoculture stands
! differentiation between coniferous and deciduous trees
do i = 1, nr_spec
if (nr_spec.eq.2 .or. nr_spec.eq.3) then
do j = 1, nr_timb_grades
ms_costs(i,:) = hsystem(2) * mansort_tg(i,j,:) * chainsaw_prices(i,j) + &
hsystem(1) * mansort_tg(i,j,:) * harvester_prices(i,j) + &
ms_costs(i,:)
st_costs(i,:) = hsystem(2) * standsort_tg(i,j,:) * chainsaw_prices(i,j) + &
hsystem(1) * standsort_tg(i,j,:) * harvester_prices(i,j) + &
st_costs(i,:)
end do
else
do j = 1, nr_timb_grades
ms_costs(i,:) = mansort_tg(i,j,:) * chainsaw_prices(i,j) + ms_costs(i,:)
st_costs(i,:) = standsort_tg(i,j,:) * chainsaw_prices(i,j) + st_costs(i,:)
end do
endif
end do
end subroutine calculate_harvest_costs
!*****************************************************************************
subroutine calculate_assets()
use data_wpm
use data_simul
implicit none
integer i, j
do i = 1, nr_spec
do j = 1, nr_timb_grades
ms_assets(i,:) = mansort_tg(i,j,:) * net_prices(i,j) + &
ms_assets(i,:)
st_assets(i,:) = standsort_tg(i,j,:) * net_prices(i,j) + &
st_assets(i,:)
end do
end do
end subroutine calculate_assets
!*****************************************************************************
subroutine calculate_costs()
use data_wpm
use data_simul
use data_plant
implicit none
character(30) manag
integer i, act_year, spec
! sum of standsort
do i = 1, nr_spec
sum_costs(2,:) = st_assets(i,:) - st_costs(i,:) + sum_costs(2,:)
end do
! sum of mansort
do i = 1, nr_spec
sum_costs(3,:) = ms_assets(i,:) - ms_costs(i,:) + sum_costs(3,:)
end do
! silvicultural costs like tending etc.
act_manrec => first_manrec
do while (associated(act_manrec))
act_year = act_manrec%manrec%year
manag = trim(act_manrec%manrec%management)
! sum up silvicultural costs and subsidies
select case (trim(manag))
case ('tending')
! tending costs and subsidies
sum_costs(4, act_year) = - tending_prices(1) + sum_costs(4, act_year)
subsidy(2, act_year) = tending_prices(2) + subsidy(2, act_year)
case ('brushing')
! brushing costs and subsidies
sum_costs(4, act_year) = - brushing(1) + sum_costs(4, act_year)
subsidy(2, act_year) = brushing(2) + subsidy(2, act_year)
case ('felling')
! forest maintenance costs and subsidies
sum_costs(4,act_year) = - ext_for(2,1) + sum_costs(4,act_year)
sum_costs(5,act_year) = ext_for(2,2) + sum_costs(5,act_year)
case ('shelterwood system1')
! forest maintenance costs and subsidies
sum_costs(4,act_year) = - ext_for(1,1) + sum_costs(4,act_year)
sum_costs(4,act_year) = ext_for(1,2) + sum_costs(4,act_year)
sum_costs(4,act_year) = - ext_for(2,1) + sum_costs(4,act_year)
sum_costs(4,act_year) = ext_for(2,2) + sum_costs(4,act_year)
case ('shelterwood system2')
! brushing, forest maintenance, timber selling costs and subsidies
sum_costs(4,act_year) = - ext_for(1,1) + sum_costs(4,act_year)
sum_costs(4,act_year) = ext_for(1,2) + sum_costs(4,act_year)
sum_costs(4,act_year) = - ext_for(2,1) + sum_costs(4,act_year)
sum_costs(4,act_year) = ext_for(2,2) + sum_costs(4,act_year)
sum_costs(4, act_year) = - brushing(1) + sum_costs(4, act_year)
subsidy(2, act_year) = brushing(2) + subsidy(2, act_year)
case ('felling after shelterwood s.')
sum_costs(4,act_year) = - ext_for(2,1) + sum_costs(4,act_year)
sum_costs(4,act_year) = ext_for(2,2) + sum_costs(4,act_year)
case ('thinning')
! forest maintenance, timber selling
sum_costs(4,act_year) = - ext_for(1,1) + sum_costs(4,act_year)
sum_costs(4,act_year) = ext_for(1,2) + sum_costs(4,act_year)
sum_costs(4,act_year) = - ext_for(2,1) + sum_costs(4,act_year)
sum_costs(4,act_year) = ext_for(2,2) + sum_costs(4,act_year)
end select
act_manrec => act_manrec%next
end do
! planting
if (plant_year /= 0) then
select case (flag_plant)
case(8,7,6,5,4,33)
sum_costs(4,plant_year) = - sum(planting_prices * (npl_mix / 1000)) + sum_costs(4,plant_year)
subsidy(2,plant_year) = sum(planting_sub(1,plant_year)*npl_mix/npl_mix) + subsidy(2,plant_year)
! pine, beech, oak, spruce, birch
case (10)
spec = 3
case (11)
spec = 1
case (12)
spec = 4
case (13)
spec = 2
case (14)
spec = 5
end select
sum_costs(4,plant_year) = - planting_prices(spec) * numplant(spec)/1000 + sum_costs(4,plant_year)
subsidy(2,plant_year) = planting_sub(1,spec) + subsidy(2,plant_year)
! fence
sum_costs(4,plant_year) = - fence(1,spec) + sum_costs(4,plant_year)
sum_costs(5,plant_year) = fence(2,spec) + sum_costs(5,plant_year)
end if
! sum up subsidies
sum_costs(5,:) = subsidy(1,:) + subsidy(2,:) + fix(2) + sum_costs(5,:)
! sum up all except standsort
sum_costs(1,:) = sum_costs(3,:) + sum_costs(4,:) - fix(1) + sum_costs(5,:)
end subroutine calculate_costs
!*****************************************************************************
subroutine calculate_npv()
use data_wpm
use data_simul
use data_plant
implicit none
real, dimension(4, nr_years) :: rate
integer i, j
rate(:,:) = 0.
do i = 1, nr_years
do j = 1, 4
rate(j, i) = (1+int_rate(j))**i
npv(j,i) = (sum_costs(2,i) + sum_costs(3,i)) / rate(j, i)
npv(j+4,i) = sum(sum_costs(1,1:i)/rate(j,1:i))
npv(j+8,i) = npv(j+4,i) - npv(1,1) + npv(j,i)
end do
end do
end subroutine calculate_npv
!*****************************************************************************
subroutine read_sea_prices(pricesFile)
use data_wpm
use data_simul
implicit none
character(70) pricesFile
integer i, unit, ios
unit = getunit()
open(unit, FILE=pricesFile, STATUS='OLD', action='read')
! Headers
do i=1,5
read (unit, *)
enddo
read(unit, '(F6.2)',iostat=ios) fix(1)
read (unit, *)
read(unit, *)
read(unit, '(5(F6.2))',iostat=ios) &
planting_prices(1), &
planting_prices(2), &
planting_prices(3), &
planting_prices(4), &
planting_prices(5)
read(unit, *)
read(unit, '(5(F6.2))',iostat=ios) &
fence(1,1), &
fence(1,2), &
fence(1,3), &
fence(1,4), &
fence(1,5)
read(unit, *)
read(unit, '(F6.2)',iostat=ios) brushing(1)
read(unit, *)
read(unit, '(F6.2)',iostat=ios) tending_prices(1)
read(unit, *)
read(unit, *)
read(unit, '(2(F6.2))',iostat=ios) &
hsystem(1), &
hsystem(2)
read(unit, *)
read(unit, '(F6.2)',iostat=ios) dec_per
do i=1,4
read (unit, *)
enddo
do i=1, nr_timb_grades
read(unit, '(5(F6.2) )',iostat=ios) &
chainsaw_prices(1, i), &
chainsaw_prices(2, i), &
chainsaw_prices(3, i), &
chainsaw_prices(4, i), &
chainsaw_prices(5, i)
end do
read(unit, *)
read(unit, *)
do i=1, nr_timb_grades
read(unit, *) &
harvester_prices(1, i), &
harvester_prices(2, i), &
harvester_prices(3, i), &
harvester_prices(4, i), &
harvester_prices(5, i)
end do
read(unit, *)
read(unit, *)
do i=1, nr_timb_grades
read(unit, '(5(F6.2) )',iostat=ios) &
net_prices(1, i), &
net_prices(2, i), &
net_prices(3, i), &
net_prices(4, i), &
net_prices(5, i)
end do
read(unit, *)
read(unit, *)
read(unit, '(F6.2)',iostat=ios) ext_for(1,1)
read(unit, *)
read(unit, '(F6.2)',iostat=ios) ext_for(1,2)
read(unit, *)
read(unit, *)
read(unit, *)
read(unit, '(F6.2)',iostat=ios) fix(2)
read (unit, *)
read(unit, *)
read(unit, '(5(F6.2))',iostat=ios) &
planting_sub(1,1), &
planting_sub(1,2), &
planting_sub(1,3), &
planting_sub(1,4), &
planting_sub(1,5)
read(unit, *)
read(unit, '(5(F6.2))',iostat=ios) &
planting_sub(2,1), &
planting_sub(2,2), &
planting_sub(2,3), &
planting_sub(2,4), &
planting_sub(2,5)
read(unit, *)
read(unit, '(5(F6.2))',iostat=ios) &
fence(1,1), &
fence(1,2), &
fence(1,3), &
fence(1,4), &
fence(1,5)
read(unit, *)
read(unit, '(F6.2)',iostat=ios) brushing(2)
read(unit, *)
read(unit, '(F6.2)',iostat=ios) tending_prices(2)
read(unit, *)
read(unit, '(F6.2)',iostat=ios) ext_for(2,2)
read(unit, *)
read(unit, '(F6.2)',iostat=ios) ext_for(2,2)
read(unit, *)
read(unit, *) &
int_rate(2), &
int_rate(3), &
int_rate(4)
close(unit, status="keep")
end subroutine read_sea_prices