Rev 522 | Rev 524 | Go to most recent revision | Only display areas with differences | Regard whitespace | Details | Blame | Last modification | View Log | RSS feed
Rev 522 | Rev 523 | ||
---|---|---|---|
1 | Redirecting to URL 'https://iland.boku.ac.at/svn/iland/tags/release_1.0/src/core/snag.cpp': |
1 | Redirecting to URL 'https://iland.boku.ac.at/svn/iland/tags/release_1.0/src/core/snag.cpp': |
2 | #include "snag.h"
|
2 | #include "snag.h"
|
3 | #include "tree.h"
|
3 | #include "tree.h"
|
4 | #include "species.h"
|
4 | #include "species.h"
|
5 | #include "globalsettings.h"
|
5 | #include "globalsettings.h"
|
6 | #include "expression.h"
|
6 | #include "expression.h"
|
7 | // for calculation of climate decomposition
|
7 | // for calculation of climate decomposition
|
8 | #include "resourceunit.h"
|
8 | #include "resourceunit.h"
|
9 | #include "watercycle.h"
|
9 | #include "watercycle.h"
|
10 | #include "climate.h"
|
10 | #include "climate.h"
|
11 | 11 | ||
12 | /** @class Snag
|
12 | /** @class Snag
|
13 | Snag deals with carbon / nitrogen fluxes from the forest until the reach soil pools.
|
13 | Snag deals with carbon / nitrogen fluxes from the forest until the reach soil pools.
|
14 | Snag lives on the level of the ResourceUnit; carbon fluxes from trees enter Snag, and parts of the biomass of snags
|
14 | Snag lives on the level of the ResourceUnit; carbon fluxes from trees enter Snag, and parts of the biomass of snags
|
15 | is subsequently forwarded to the soil sub model.
|
15 | is subsequently forwarded to the soil sub model.
|
16 | Carbon is stored in three classes (depending on the size)
|
16 | Carbon is stored in three classes (depending on the size)
|
17 | 17 | ||
18 | */
|
18 | */
|
19 | // static variables
|
19 | // static variables
|
20 | double Snag::mDBHLower = 0.; |
20 | double Snag::mDBHLower = 0.; |
21 | double Snag::mDBHHigher = 0.; |
21 | double Snag::mDBHHigher = 0.; |
22 | double Snag::mCarbonThreshold[] = {0., 0., 0.}; |
22 | double Snag::mCarbonThreshold[] = {0., 0., 0.}; |
23 | 23 | ||
24 | double CNPool::biomassCFraction = biomassCFraction; // get global from globalsettings.h |
24 | double CNPool::biomassCFraction = biomassCFraction; // get global from globalsettings.h |
25 | 25 | ||
26 | void Snag::setupThresholds(const double lower, const double upper) |
26 | void Snag::setupThresholds(const double lower, const double upper) |
27 | {
|
27 | {
|
28 | if (mDBHLower == lower) |
28 | if (mDBHLower == lower) |
29 | return; |
29 | return; |
30 | mDBHLower = lower; |
30 | mDBHLower = lower; |
31 | mDBHHigher = upper; |
31 | mDBHHigher = upper; |
32 | mCarbonThreshold[0] = lower / 2.; |
32 | mCarbonThreshold[0] = lower / 2.; |
33 | mCarbonThreshold[1] = lower + (upper - lower)/2.; |
33 | mCarbonThreshold[1] = lower + (upper - lower)/2.; |
34 | mCarbonThreshold[2] = upper + (upper - lower)/2.; |
34 | mCarbonThreshold[2] = upper + (upper - lower)/2.; |
35 | //# threshold levels for emptying out the dbh-snag-classes
|
35 | //# threshold levels for emptying out the dbh-snag-classes
|
36 | //# derived from Psme woody allometry, converted to C, with a threshold level set to 10%
|
36 | //# derived from Psme woody allometry, converted to C, with a threshold level set to 10%
|
37 | //# values in kg!
|
37 | //# values in kg!
|
38 | for (int i=0;i<3;i++) |
38 | for (int i=0;i<3;i++) |
39 | mCarbonThreshold[i] = 0.10568*pow(mCarbonThreshold[i],2.4247)*0.5*0.1; |
39 | mCarbonThreshold[i] = 0.10568*pow(mCarbonThreshold[i],2.4247)*0.5*0.1; |
40 | }
|
40 | }
|
41 | 41 | ||
42 | 42 | ||
43 | // species specific soil paramters
|
43 | // species specific soil paramters
|
44 | struct SoilParameters
|
44 | struct SoilParameters
|
45 | {
|
45 | {
|
46 | SoilParameters(): kyl(0.15), kyr(0.0807), ksw(0.015), cnFoliage(75.), cnFineroot(40.), cnWood(300.), halfLife(15.) {} |
46 | SoilParameters(): kyl(0.15), kyr(0.0807), ksw(0.015), cnFoliage(75.), cnFineroot(40.), cnWood(300.), halfLife(15.) {} |
47 | double kyl; // litter decomposition rate |
47 | double kyl; // litter decomposition rate |
48 | double kyr; // downed woody debris (dwd) decomposition rate |
48 | double kyr; // downed woody debris (dwd) decomposition rate |
49 | double ksw; // standing woody debris (swd) decomposition rate |
49 | double ksw; // standing woody debris (swd) decomposition rate |
50 | double cnFoliage; // C/N foliage litter |
50 | double cnFoliage; // C/N foliage litter |
51 | double cnFineroot; // C/N ratio fine root |
51 | double cnFineroot; // C/N ratio fine root |
52 | double cnWood; // C/N Wood: used for brances, stem and coarse root |
52 | double cnWood; // C/N Wood: used for brances, stem and coarse root |
53 | double halfLife; // half-life period of the given species (for calculation of SWD->DWD transition) |
53 | double halfLife; // half-life period of the given species (for calculation of SWD->DWD transition) |
54 | } soilparams; |
54 | } soilparams; |
55 | 55 | ||
56 | 56 | ||
57 | Snag::Snag() |
57 | Snag::Snag() |
58 | {
|
58 | {
|
59 | mRU = 0; |
59 | mRU = 0; |
60 | CNPool::setCFraction(biomassCFraction); |
60 | CNPool::setCFraction(biomassCFraction); |
61 | Snag::setupThresholds( |
61 | Snag::setupThresholds(20., 60.); |
62 | }
|
62 | }
|
63 | 63 | ||
64 | void Snag::setup( const ResourceUnit *ru) |
64 | void Snag::setup( const ResourceUnit *ru) |
65 | {
|
65 | {
|
66 | mRU = ru; |
66 | mRU = ru; |
67 | mClimateFactor = 0.; |
67 | mClimateFactor = 0.; |
68 | // branches
|
68 | // branches
|
69 | mBranchCounter=0; |
69 | mBranchCounter=0; |
70 | for (int i=0;i<3;i++) { |
70 | for (int i=0;i<3;i++) { |
71 | mTimeSinceDeath[i] = 0.; |
71 | mTimeSinceDeath[i] = 0.; |
72 | mNumberOfSnags[i] = 0.; |
72 | mNumberOfSnags[i] = 0.; |
73 | mAvgDbh[i] = 0.; |
73 | mAvgDbh[i] = 0.; |
74 | mAvgHeight[i] = 0.; |
74 | mAvgHeight[i] = 0.; |
75 | mAvgVolume[i] = 0.; |
75 | mAvgVolume[i] = 0.; |
76 | mKSW[i] = 0.; |
76 | mKSW[i] = 0.; |
77 | mCurrentKSW[i] = 0.; |
77 | mCurrentKSW[i] = 0.; |
78 | }
|
78 | }
|
79 | mTotalSnagCarbon = 0.; |
79 | mTotalSnagCarbon = 0.; |
80 | }
|
80 | }
|
81 | 81 | ||
82 | // debug outputs
|
82 | // debug outputs
|
83 | QList<QVariant> Snag::debugList() |
83 | QList<QVariant> Snag::debugList() |
84 | {
|
84 | {
|
85 | // list columns
|
85 | // list columns
|
86 | // for three pools
|
86 | // for three pools
|
87 | QList<QVariant> list; |
87 | QList<QVariant> list; |
88 | 88 | ||
89 |
|
89 | // totals
|
- | 90 | list << mTotalSnagCarbon << mTotalIn.C << mTotalToAtm.C << mSWDtoSoil.C << mSWDtoSoil.N; |
|
90 | // fluxes to labile soil pool and to refractory soil pool
|
91 | // fluxes to labile soil pool and to refractory soil pool
|
91 | list << mLabileFlux.C << mLabileFlux.N << mRefractoryFlux.C << mRefractoryFlux.N << mSWDtoSoil.C << mSWDtoSoil.N; |
92 | list << mLabileFlux.C << mLabileFlux.N << mRefractoryFlux.C << mRefractoryFlux.N << mSWDtoSoil.C << mSWDtoSoil.N; |
92 | 93 | ||
93 | for (int i=0;i<3;i++) { |
94 | for (int i=0;i<3;i++) { |
94 | // pools "swdx_c", "swdx_n", "swdx_count", "swdx_tsd", "toswdx_c", "toswdx_n"
|
95 | // pools "swdx_c", "swdx_n", "swdx_count", "swdx_tsd", "toswdx_c", "toswdx_n"
|
95 | list << mSWD[i].C << mSWD[i].N << mNumberOfSnags[i] << mTimeSinceDeath[i] << mToSWD[i].C << mToSWD[i].N; |
96 | list << mSWD[i].C << mSWD[i].N << mNumberOfSnags[i] << mTimeSinceDeath[i] << mToSWD[i].C << mToSWD[i].N; |
96 | }
|
97 | }
|
97 | 98 | ||
98 | // branch pools (5 yrs)
|
99 | // branch pools (5 yrs)
|
99 | list << mBranches[mBranchCounter].C << mBranches[mBranchCounter].N |
100 | list << mBranches[mBranchCounter].C << mBranches[mBranchCounter].N |
100 | << mBranches[(mBranchCounter+1)%5].C << mBranches[(mBranchCounter+1)%5].N |
101 | << mBranches[(mBranchCounter+1)%5].C << mBranches[(mBranchCounter+1)%5].N |
101 | << mBranches[(mBranchCounter+2)%5].C << mBranches[(mBranchCounter+2)%5].N |
102 | << mBranches[(mBranchCounter+2)%5].C << mBranches[(mBranchCounter+2)%5].N |
102 | << mBranches[(mBranchCounter+3)%5].C << mBranches[(mBranchCounter+3)%5].N |
103 | << mBranches[(mBranchCounter+3)%5].C << mBranches[(mBranchCounter+3)%5].N |
103 | << mBranches[(mBranchCounter+4)%5].C << mBranches[(mBranchCounter+4)%5].N; |
104 | << mBranches[(mBranchCounter+4)%5].C << mBranches[(mBranchCounter+4)%5].N; |
104 | return list; |
105 | return list; |
105 | }
|
106 | }
|
106 | 107 | ||
107 | void Snag::newYear() |
108 | void Snag::newYear() |
108 | {
|
109 | {
|
109 | for (int i=0;i<3;i++) { |
110 | for (int i=0;i<3;i++) { |
110 | mToSWD[i].clear(); // clear transfer pools to standing-woody-debris |
111 | mToSWD[i].clear(); // clear transfer pools to standing-woody-debris |
111 | mCurrentKSW[i] = 0.; |
112 | mCurrentKSW[i] = 0.; |
112 | }
|
113 | }
|
113 | mLabileFlux.clear(); |
114 | mLabileFlux.clear(); |
114 | mRefractoryFlux.clear(); |
115 | mRefractoryFlux.clear(); |
115 | mTotalToAtm.clear(); |
116 | mTotalToAtm.clear(); |
116 | mTotalToExtern.clear(); |
117 | mTotalToExtern.clear(); |
117 | mTotalIn.clear(); |
118 | mTotalIn.clear(); |
118 | mSWDtoSoil.clear(); |
119 | mSWDtoSoil.clear(); |
119 | }
|
120 | }
|
120 | 121 | ||
121 | /// calculate the dynamic climate modifier for decomposition 're'
|
122 | /// calculate the dynamic climate modifier for decomposition 're'
|
122 | /// calculation is done on the level of ResourceUnit because the water content per day is needed.
|
123 | /// calculation is done on the level of ResourceUnit because the water content per day is needed.
|
123 | double Snag::calculateClimateFactors() |
124 | double Snag::calculateClimateFactors() |
124 | {
|
125 | {
|
125 | double rel_wc; |
126 | double rel_wc; |
126 | double ft, fw; |
127 | double ft, fw; |
127 | double f_sum = 0.; |
128 | double f_sum = 0.; |
128 | for (const ClimateDay *day=mRU->climate()->begin(); day!=mRU->climate()->end(); ++day) |
129 | for (const ClimateDay *day=mRU->climate()->begin(); day!=mRU->climate()->end(); ++day) |
129 | {
|
130 | {
|
130 | rel_wc = mRU->waterCycle()->relContent(day->day); // relative water content |
131 | rel_wc = mRU->waterCycle()->relContent(day->day)*100.; // relative water content in per cent of the day |
131 | ft = exp(308.56*(1./56.02-1./((273.+day->temperature)-227.13))); // empirical variable Q10 model of Lloyd and Taylor (1994), see also Adair et al. (2008) |
132 | ft = exp(308.56*(1./56.02-1./((273.+day->temperature)-227.13))); // empirical variable Q10 model of Lloyd and Taylor (1994), see also Adair et al. (2008) |
132 | fw = pow(1.-exp(-0.2*rel_wc),5.); // # see Standcarb for the 'stable soil' pool |
133 | fw = pow(1.-exp(-0.2*rel_wc),5.); // # see Standcarb for the 'stable soil' pool |
133 | f_sum += ft*fw; |
134 | f_sum += ft*fw; |
134 | }
|
135 | }
|
135 | // the climate factor is defined as the arithmentic annual mean value
|
136 | // the climate factor is defined as the arithmentic annual mean value
|
136 | mClimateFactor = f_sum / double(mRU->climate()->daysOfYear()); |
137 | mClimateFactor = f_sum / double(mRU->climate()->daysOfYear()); |
137 | return mClimateFactor; |
138 | return mClimateFactor; |
138 | }
|
139 | }
|
139 | 140 | ||
140 | /// do the yearly calculation
|
141 | /// do the yearly calculation
|
141 | /// see http://iland.boku.ac.at/snag+dynamics
|
142 | /// see http://iland.boku.ac.at/snag+dynamics
|
142 | void Snag::processYear() |
143 | void Snag::processYear() |
143 | {
|
144 | {
|
144 | mSWDtoSoil.clear(); |
145 | mSWDtoSoil.clear(); |
145 | if (isEmpty()) // nothing to do |
146 | if (isEmpty()) // nothing to do |
146 | return; |
147 | return; |
147 | 148 | ||
148 | // process branches: every year one of the five baskets is emptied and transfered to the refractory soil pool
|
149 | // process branches: every year one of the five baskets is emptied and transfered to the refractory soil pool
|
149 | mRefractoryFlux+=mBranches[mBranchCounter]; |
150 | mRefractoryFlux+=mBranches[mBranchCounter]; |
- | 151 | mSWDtoSoil += mBranches[mBranchCounter]; |
|
150 | mBranches[mBranchCounter].clear(); |
152 | mBranches[mBranchCounter].clear(); |
151 | mBranchCounter= (mBranchCounter+1) % 5; // increase index, roll over to 0. |
153 | mBranchCounter= (mBranchCounter+1) % 5; // increase index, roll over to 0. |
152 | 154 | ||
153 | // process standing snags.
|
155 | // process standing snags.
|
154 | // the input of the current year is in the mToSWD-Pools
|
156 | // the input of the current year is in the mToSWD-Pools
|
155 | 157 | ||
156 | const double climate_factor_re = calculateClimateFactors(); |
158 | const double climate_factor_re = calculateClimateFactors(); |
157 | for (int i=0;i<3;i++) { |
159 | for (int i=0;i<3;i++) { |
158 | 160 | ||
159 | // update the swd-pool with this years' input
|
161 | // update the swd-pool with this years' input
|
160 | if (!mToSWD[i].isEmpty()) { |
162 | if (!mToSWD[i].isEmpty()) { |
161 | // update decay rate (apply average yearly input to the state parameters)
|
163 | // update decay rate (apply average yearly input to the state parameters)
|
162 | mKSW[i] = mKSW[i]*(mSWD[i].C/(mSWD[i].C+mToSWD[i].C)) + mCurrentKSW[i]*(mToSWD[i].C/(mSWD[i].C+mToSWD[i].C)); |
164 | mKSW[i] = mKSW[i]*(mSWD[i].C/(mSWD[i].C+mToSWD[i].C)) + mCurrentKSW[i]*(mToSWD[i].C/(mSWD[i].C+mToSWD[i].C)); |
163 | //move content to the SWD pool
|
165 | //move content to the SWD pool
|
164 | mSWD[i] += mToSWD[i]; |
166 | mSWD[i] += mToSWD[i]; |
165 | }
|
167 | }
|
166 | 168 | ||
167 | if (mSWD[i].C > 0) { |
169 | if (mSWD[i].C > 0) { |
168 | // reduce the Carbon (note: the N stays, thus the CN ratio changes)
|
170 | // reduce the Carbon (note: the N stays, thus the CN ratio changes)
|
169 | // use the decay rate that is derived as a weighted average of all standing woody debris
|
171 | // use the decay rate that is derived as a weighted average of all standing woody debris
|
170 |
|
172 | double survive_rate = exp(-mKSW[i] *climate_factor_re * 1. ); // 1: timestep |
- | 173 | mTotalToAtm.C += mSWD[i].C * (1. - survive_rate); |
|
- | 174 | mSWD[i].C *= survive_rate; |
|
171 | 175 | ||
172 | // transition to downed woody debris
|
176 | // transition to downed woody debris
|
173 | // update: use negative exponential decay, species parameter: half-life
|
177 | // update: use negative exponential decay, species parameter: half-life
|
174 | // modified for the climatic effect on decomposition, i.e. if decomp is slower, snags stand longer and vice versa
|
178 | // modified for the climatic effect on decomposition, i.e. if decomp is slower, snags stand longer and vice versa
|
175 | // this is loosely oriented on Standcarb2 (http://andrewsforest.oregonstate.edu/pubs/webdocs/models/standcarb2.htm),
|
179 | // this is loosely oriented on Standcarb2 (http://andrewsforest.oregonstate.edu/pubs/webdocs/models/standcarb2.htm),
|
176 | // where lag times for cohort transitions are linearly modified with re although here individual good or bad years have
|
180 | // where lag times for cohort transitions are linearly modified with re although here individual good or bad years have
|
177 | // an immediate effect, the average climatic influence should come through (and it is inherently transient)
|
181 | // an immediate effect, the average climatic influence should come through (and it is inherently transient)
|
178 | // note that swd.hl is species-specific, and thus a weighted average over the species in the input (=mortality)
|
182 | // note that swd.hl is species-specific, and thus a weighted average over the species in the input (=mortality)
|
179 | // needs to be calculated, followed by a weighted update of the previous swd.hl.
|
183 | // needs to be calculated, followed by a weighted update of the previous swd.hl.
|
180 | // As weights here we use stem number, as the processes here pertain individual snags
|
184 | // As weights here we use stem number, as the processes here pertain individual snags
|
181 | // calculate the transition probability of SWD to downed dead wood
|
185 | // calculate the transition probability of SWD to downed dead wood
|
182 | 186 | ||
183 | double half_life = mHalfLife[i] / climate_factor_re; |
187 | double half_life = mHalfLife[i] / climate_factor_re; |
184 | double rate = -M_LN2 / half_life; // M_LN2: math. constant |
188 | double rate = -M_LN2 / half_life; // M_LN2: math. constant |
185 | 189 | ||
186 | // higher decay rate for the class with smallest diameters
|
190 | // higher decay rate for the class with smallest diameters
|
187 | if (i==0) |
191 | if (i==0) |
188 | rate*=2.; |
192 | rate*=2.; |
189 | 193 | ||
190 | double transfer = exp(rate); |
194 | double transfer = 1. - exp(rate); |
191 | 195 | ||
192 | // calculate flow to soil pool...
|
196 | // calculate flow to soil pool...
|
193 | mSWDtoSoil += mSWD[i] * transfer; |
197 | mSWDtoSoil += mSWD[i] * transfer; |
194 | mRefractoryFlux += mSWD[i] * transfer; |
198 | mRefractoryFlux += mSWD[i] * transfer; |
195 | mSWD[i] *= (1.-transfer); // reduce pool |
199 | mSWD[i] *= (1.-transfer); // reduce pool |
196 | // calculate the stem number of remaining snags
|
200 | // calculate the stem number of remaining snags
|
197 | mNumberOfSnags[i] = mNumberOfSnags[i] * (1. - transfer); |
201 | mNumberOfSnags[i] = mNumberOfSnags[i] * (1. - transfer); |
- | 202 | ||
- | 203 | mTimeSinceDeath[i] += 1.; |
|
198 | // if stems<0.5, empty the whole cohort into DWD, i.e. release the last bit of C and N and clear the stats
|
204 | // if stems<0.5, empty the whole cohort into DWD, i.e. release the last bit of C and N and clear the stats
|
199 | // also, if the Carbon of an average snag is less than 10% of the original average tree
|
205 | // also, if the Carbon of an average snag is less than 10% of the original average tree
|
200 | // (derived from allometries for the three diameter classes), the whole cohort is emptied out to DWD
|
206 | // (derived from allometries for the three diameter classes), the whole cohort is emptied out to DWD
|
201 | if (mNumberOfSnags[i] < 0.5 || mSWD[i].C / mNumberOfSnags[i] < mCarbonThreshold[i]) { |
207 | if (mNumberOfSnags[i] < 0.5 || mSWD[i].C / mNumberOfSnags[i] < mCarbonThreshold[i]) { |
202 | // clear the pool: add the rest to the soil, clear statistics of the pool
|
208 | // clear the pool: add the rest to the soil, clear statistics of the pool
|
203 | mRefractoryFlux += mSWD[i]; |
209 | mRefractoryFlux += mSWD[i]; |
204 | mSWDtoSoil += mSWD[i]; |
210 | mSWDtoSoil += mSWD[i]; |
205 | mSWD[i].clear(); |
211 | mSWD[i].clear(); |
206 | mAvgDbh[i] = 0.; |
212 | mAvgDbh[i] = 0.; |
207 | mAvgHeight[i] = 0.; |
213 | mAvgHeight[i] = 0.; |
208 | mAvgVolume[i] = 0.; |
214 | mAvgVolume[i] = 0.; |
209 | mKSW[i] = 0.; |
215 | mKSW[i] = 0.; |
210 | mCurrentKSW[i] = 0.; |
216 | mCurrentKSW[i] = 0.; |
211 | mHalfLife[i] = 0.; |
217 | mHalfLife[i] = 0.; |
212 | mTimeSinceDeath[i] = 0.; |
218 | mTimeSinceDeath[i] = 0.; |
213 | }
|
219 | }
|
214 | 220 | ||
215 | }
|
221 | }
|
216 | 222 | ||
217 | }
|
223 | }
|
218 | // total carbon in the snag-container on the RU *after* processing is the content of the
|
224 | // total carbon in the snag-container on the RU *after* processing is the content of the
|
219 | // standing woody debris pools + the branches
|
225 | // standing woody debris pools + the branches
|
220 | mTotalSnagCarbon = mSWD[0].C + mSWD[1].C + mSWD[2].C + |
226 | mTotalSnagCarbon = mSWD[0].C + mSWD[1].C + mSWD[2].C + |
221 | mBranches[0].C + mBranches[1].C + mBranches[2].C + mBranches[3].C + mBranches[4].C; |
227 | mBranches[0].C + mBranches[1].C + mBranches[2].C + mBranches[3].C + mBranches[4].C; |
222 | }
|
228 | }
|
223 | 229 | ||
224 | /// foliage and fineroot litter is transferred during tree growth.
|
230 | /// foliage and fineroot litter is transferred during tree growth.
|
225 | void Snag::addTurnoverLitter(const Tree *tree, const double litter_foliage, const double litter_fineroot) |
231 | void Snag::addTurnoverLitter(const Tree *tree, const double litter_foliage, const double litter_fineroot) |
226 | {
|
232 | {
|
227 | const SoilParameters &soil_params = soilparams; // to be updated |
233 | const SoilParameters &soil_params = soilparams; // to be updated |
228 | mLabileFlux.addBiomass(litter_foliage, soil_params.cnFoliage); |
234 | mLabileFlux.addBiomass(litter_foliage, soil_params.cnFoliage); |
229 | mLabileFlux.addBiomass(litter_fineroot, soil_params.cnFineroot); |
235 | mLabileFlux.addBiomass(litter_fineroot, soil_params.cnFineroot); |
230 | }
|
236 | }
|
231 | 237 | ||
232 | /// after the death of the tree the five biomass compartments are processed.
|
238 | /// after the death of the tree the five biomass compartments are processed.
|
233 | void Snag::addMortality(const Tree *tree) |
239 | void Snag::addMortality(const Tree *tree) |
234 | {
|
240 | {
|
235 | const SoilParameters &soil_params = soilparams; // to be updated |
241 | const SoilParameters &soil_params = soilparams; // to be updated |
236 | 242 | ||
237 | // immediate flows: 100% of foliage, 100% of fine roots: they go to the labile pool
|
243 | // immediate flows: 100% of foliage, 100% of fine roots: they go to the labile pool
|
238 | // 100% of coarse root biomass: that goes to the refractory pool
|
244 | // 100% of coarse root biomass: that goes to the refractory pool
|
239 | mLabileFlux.addBiomass(tree->biomassFoliage(), soil_params.cnFoliage); |
245 | mLabileFlux.addBiomass(tree->biomassFoliage(), soil_params.cnFoliage); |
240 | mLabileFlux.addBiomass(tree->biomassFineRoot(), soil_params.cnFineroot); |
246 | mLabileFlux.addBiomass(tree->biomassFineRoot(), soil_params.cnFineroot); |
241 | mRefractoryFlux.addBiomass(tree->biomassCoarseRoot(), soil_params.cnWood); |
247 | mRefractoryFlux.addBiomass(tree->biomassCoarseRoot(), soil_params.cnWood); |
242 | 248 | ||
243 | // branches are equally distributed over five years:
|
249 | // branches are equally distributed over five years:
|
244 | double biomass_branch = tree->biomassBranch() * 0.2; |
250 | double biomass_branch = tree->biomassBranch() * 0.2; |
245 | for (int i=0;i<5; i++) |
251 | for (int i=0;i<5; i++) |
246 | mBranches[i].addBiomass(biomass_branch, soil_params.cnWood); |
252 | mBranches[i].addBiomass(biomass_branch, soil_params.cnWood); |
247 | 253 | ||
- | 254 | // just for book-keeping: keep track of all inputs into branches / swd
|
|
- | 255 | mTotalIn.addBiomass(tree->biomassBranch() + tree->biomassStem(), soil_params.cnWood); |
|
248 | // stem biomass is transferred to the standing woody debris pool (SWD), increase stem number of pool
|
256 | // stem biomass is transferred to the standing woody debris pool (SWD), increase stem number of pool
|
249 | int pi = poolIndex(tree->dbh()); // get right transfer pool |
257 | int pi = poolIndex(tree->dbh()); // get right transfer pool |
250 | CNPool &swd = mToSWD[pi]; |
- | |
251 | swd.addBiomass(tree->biomassStem(), soil_params.cnWood); |
- | |
252 | 258 | ||
253 | // update statistics - stemnumber-weighted averages
|
259 | // update statistics - stemnumber-weighted averages
|
254 | // note: here the calculations are repeated for every died trees (i.e. consecutive weighting ... but delivers the same results)
|
260 | // note: here the calculations are repeated for every died trees (i.e. consecutive weighting ... but delivers the same results)
|
255 | double p_old = mNumberOfSnags[pi] / (mNumberOfSnags[pi] + 1); // weighting factor for state vars (based on stem numbers) |
261 | double p_old = mNumberOfSnags[pi] / (mNumberOfSnags[pi] + 1); // weighting factor for state vars (based on stem numbers) |
256 | double p_new = 1. / (mNumberOfSnags[pi] + 1); // weighting factor for added tree (p_old + p_new = 1). |
262 | double p_new = 1. / (mNumberOfSnags[pi] + 1); // weighting factor for added tree (p_old + p_new = 1). |
257 | mAvgDbh[pi] = mAvgDbh[pi]*p_old + tree->dbh()*p_new; |
263 | mAvgDbh[pi] = mAvgDbh[pi]*p_old + tree->dbh()*p_new; |
258 | mAvgHeight[pi] = mAvgHeight[pi]*p_old + tree->height()*p_new; |
264 | mAvgHeight[pi] = mAvgHeight[pi]*p_old + tree->height()*p_new; |
259 | mAvgVolume[pi] = mAvgVolume[pi]*p_old + tree->volume()*p_new; |
265 | mAvgVolume[pi] = mAvgVolume[pi]*p_old + tree->volume()*p_new; |
260 | mTimeSinceDeath[pi] = mTimeSinceDeath[pi]*p_old + 1.*p_new; |
266 | mTimeSinceDeath[pi] = mTimeSinceDeath[pi]*p_old + 1.*p_new; |
261 | mHalfLife[pi] = mHalfLife[pi]*p_old + soil_params.halfLife * p_new; |
267 | mHalfLife[pi] = mHalfLife[pi]*p_old + soil_params.halfLife * p_new; |
262 | 268 | ||
263 | // average the decay rate (ksw); this is done based on the carbon content
|
269 | // average the decay rate (ksw); this is done based on the carbon content
|
264 | // aggregate all trees that die in the current year (and save weighted decay rates to CurrentKSW)
|
270 | // aggregate all trees that die in the current year (and save weighted decay rates to CurrentKSW)
|
265 | if (tree->biomassStem()==0) |
271 | if (tree->biomassStem()==0) |
266 | throw IException("Snag::addMortality: tree without stem biomass!!"); |
272 | throw IException("Snag::addMortality: tree without stem biomass!!"); |
267 | p_old = mToSWD[pi].C / (mToSWD[pi].C + tree->biomassStem()* biomassCFraction); |
273 | p_old = mToSWD[pi].C / (mToSWD[pi].C + tree->biomassStem()* biomassCFraction); |
268 | p_new =tree->biomassStem()* biomassCFraction / (mToSWD[pi].C + tree->biomassStem()* biomassCFraction); |
274 | p_new =tree->biomassStem()* biomassCFraction / (mToSWD[pi].C + tree->biomassStem()* biomassCFraction); |
269 | mCurrentKSW[pi] = mCurrentKSW[pi]*p_old + soil_params.ksw * p_new; |
275 | mCurrentKSW[pi] = mCurrentKSW[pi]*p_old + soil_params.ksw * p_new; |
270 | mNumberOfSnags[pi]++; |
276 | mNumberOfSnags[pi]++; |
- | 277 | ||
- | 278 | // finally add the biomass
|
|
- | 279 | CNPool &swd = mToSWD[pi]; |
|
- | 280 | swd.addBiomass(tree->biomassStem(), soil_params.cnWood); |
|
271 | }
|
281 | }
|
272 | 282 | ||
273 | /// add residual biomass of 'tree' after harvesting.
|
283 | /// add residual biomass of 'tree' after harvesting.
|
274 | /// remove_{stem, branch, foliage}_fraction: percentage of biomass compartment that is *removed* by the harvest operation (i.e.: not to stay in the system)
|
284 | /// remove_{stem, branch, foliage}_fraction: percentage of biomass compartment that is *removed* by the harvest operation (i.e.: not to stay in the system)
|
275 | /// the harvested biomass is collected.
|
285 | /// the harvested biomass is collected.
|
276 | void Snag::addHarvest(const Tree* tree, const double remove_stem_fraction, const double remove_branch_fraction, const double remove_foliage_fraction ) |
286 | void Snag::addHarvest(const Tree* tree, const double remove_stem_fraction, const double remove_branch_fraction, const double remove_foliage_fraction ) |
277 | {
|
287 | {
|
278 | const SoilParameters &soil_params = soilparams; // to be updated |
288 | const SoilParameters &soil_params = soilparams; // to be updated |
279 | 289 | ||
280 | // immediate flows: 100% of residual foliage, 100% of fine roots: they go to the labile pool
|
290 | // immediate flows: 100% of residual foliage, 100% of fine roots: they go to the labile pool
|
281 | // 100% of coarse root biomass: that goes to the refractory pool
|
291 | // 100% of coarse root biomass: that goes to the refractory pool
|
282 | mLabileFlux.addBiomass(tree->biomassFoliage() * (1. - remove_foliage_fraction), soil_params.cnFoliage); |
292 | mLabileFlux.addBiomass(tree->biomassFoliage() * (1. - remove_foliage_fraction), soil_params.cnFoliage); |
283 | mLabileFlux.addBiomass(tree->biomassFineRoot(), soil_params.cnFineroot); |
293 | mLabileFlux.addBiomass(tree->biomassFineRoot(), soil_params.cnFineroot); |
284 | mRefractoryFlux.addBiomass(tree->biomassCoarseRoot(), soil_params.cnWood); |
294 | mRefractoryFlux.addBiomass(tree->biomassCoarseRoot(), soil_params.cnWood); |
285 | 295 | ||
286 | // residual branches are equally distributed over five years:
|
296 | // residual branches are equally distributed over five years:
|
287 | for (int i=0;i<5; i++) |
297 | for (int i=0;i<5; i++) |
288 | mBranches[i].addBiomass(tree->biomassBranch() * remove_branch_fraction * 0.2, soil_params.cnWood); |
298 | mBranches[i].addBiomass(tree->biomassBranch() * (1. - remove_branch_fraction) * 0.2, soil_params.cnWood); |
289 | 299 | ||
- | 300 | mTotalToExtern.addBiomass(tree->biomassBranch()*remove_branch_fraction + tree->biomassStem()*remove_stem_fraction, soil_params.cnWood); |
|
290 | // stem biomass is transferred to the standing woody debris pool (SWD), increase stem number of pool
|
301 | // stem biomass is transferred to the standing woody debris pool (SWD), increase stem number of pool
|
291 | // TODO: what to do with harvest and stems??? I think harvested stems (that are not removed) should
|
302 | // TODO: what to do with harvest and stems??? I think harvested stems (that are not removed) should
|
292 | // go directly to the DWD??
|
303 | // go directly to the DWD??
|
293 | CNPool &swd = mToSWD[poolIndex(tree->dbh())]; // get right transfer pool |
304 | CNPool &swd = mToSWD[poolIndex(tree->dbh())]; // get right transfer pool |
294 | swd.addBiomass(tree->biomassStem() * remove_stem_fraction, soil_params.cnWood); |
305 | swd.addBiomass(tree->biomassStem() * remove_stem_fraction, soil_params.cnWood); |
295 | if (remove_stem_fraction < 1.) |
306 | if (remove_stem_fraction < 1.) |
296 | mNumberOfSnags[poolIndex(tree->dbh())]++; |
307 | mNumberOfSnags[poolIndex(tree->dbh())]++; |
297 | }
|
308 | }
|
298 | 309 | ||
299 | 310 |