Rev 1111 | Rev 1162 | Go to most recent revision | Details | Compare with Previous | Last modification | View Log | RSS feed
Rev | Author | Line No. | Line |
---|---|---|---|
1 | |||
671 | werner | 2 | /******************************************************************************************** |
3 | ** iLand - an individual based forest landscape and disturbance model |
||
4 | ** http://iland.boku.ac.at |
||
5 | ** Copyright (C) 2009- Werner Rammer, Rupert Seidl |
||
6 | ** |
||
7 | ** This program is free software: you can redistribute it and/or modify |
||
8 | ** it under the terms of the GNU General Public License as published by |
||
9 | ** the Free Software Foundation, either version 3 of the License, or |
||
10 | ** (at your option) any later version. |
||
11 | ** |
||
12 | ** This program is distributed in the hope that it will be useful, |
||
13 | ** but WITHOUT ANY WARRANTY; without even the implied warranty of |
||
14 | ** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the |
||
15 | ** GNU General Public License for more details. |
||
16 | ** |
||
17 | ** You should have received a copy of the GNU General Public License |
||
18 | ** along with this program. If not, see <http://www.gnu.org/licenses/>. |
||
19 | ********************************************************************************************/ |
||
20 | |||
595 | werner | 21 | #include "sapling.h" |
22 | #include "model.h" |
||
23 | #include "species.h" |
||
24 | #include "resourceunit.h" |
||
25 | #include "resourceunitspecies.h" |
||
26 | #include "tree.h" |
||
27 | |||
28 | /** @class Sapling |
||
697 | werner | 29 | @ingroup core |
595 | werner | 30 | Sapling stores saplings per species and resource unit and computes sapling growth (before recruitment). |
31 | http://iland.boku.ac.at/sapling+growth+and+competition |
||
32 | Saplings are established in a separate step (@sa Regeneration). If sapling reach a height of 4m, they are recruited and become "real" iLand-trees. |
||
33 | Within the regeneration layer, a cohort-approach is applied. |
||
34 | |||
35 | */ |
||
36 | |||
37 | double Sapling::mRecruitmentVariation = 0.1; // +/- 10% |
||
1063 | werner | 38 | double Sapling::mBrowsingPressure = 0.; |
595 | werner | 39 | |
40 | Sapling::Sapling() |
||
41 | { |
||
42 | mRUS = 0; |
||
43 | clearStatistics(); |
||
663 | werner | 44 | mAdded = 0; |
595 | werner | 45 | } |
46 | |||
663 | werner | 47 | // reset statistics, called at newYear |
48 | void Sapling::clearStatistics() |
||
49 | { |
||
50 | // mAdded: removed |
||
51 | mRecruited=mDied=mLiving=0; |
||
52 | mSumDbhDied=0.; |
||
53 | mAvgHeight=0.; |
||
54 | mAvgAge=0.; |
||
55 | mAvgDeltaHPot=mAvgHRealized=0.; |
||
56 | } |
||
57 | |||
1063 | werner | 58 | void Sapling::updateBrowsingPressure() |
59 | { |
||
60 | if (GlobalSettings::instance()->settings().valueBool("model.settings.browsing.enabled")) |
||
61 | Sapling::mBrowsingPressure = GlobalSettings::instance()->settings().valueDouble("model.settings.browsing.browsingPressure"); |
||
62 | else |
||
63 | Sapling::mBrowsingPressure = 0.; |
||
64 | } |
||
65 | |||
595 | werner | 66 | /// get the *represented* (Reineke's Law) number of trees (N/ha) |
67 | double Sapling::livingStemNumber(double &rAvgDbh, double &rAvgHeight, double &rAvgAge) const |
||
68 | { |
||
69 | double total = 0.; |
||
70 | double dbh_sum = 0.; |
||
71 | double h_sum = 0.; |
||
72 | double age_sum = 0.; |
||
73 | const SaplingGrowthParameters &p = mRUS->species()->saplingGrowthParameters(); |
||
1111 | werner | 74 | for (QVector<SaplingTreeOld>::const_iterator it = mSaplingTrees.constBegin(); it!=mSaplingTrees.constEnd(); ++it) { |
595 | werner | 75 | float dbh = it->height / p.hdSapling * 100.f; |
76 | if (dbh<1.) // minimum size: 1cm |
||
77 | continue; |
||
78 | double n = p.representedStemNumber(dbh); // one cohort on the pixel represents that number of trees |
||
79 | dbh_sum += n*dbh; |
||
80 | h_sum += n*it->height; |
||
81 | age_sum += n*it->age.age; |
||
82 | total += n; |
||
83 | } |
||
84 | if (total>0.) { |
||
85 | dbh_sum /= total; |
||
86 | h_sum /= total; |
||
87 | age_sum /= total; |
||
88 | } |
||
89 | rAvgDbh = dbh_sum; |
||
90 | rAvgHeight = h_sum; |
||
91 | rAvgAge = age_sum; |
||
92 | return total; |
||
93 | } |
||
94 | |||
967 | werner | 95 | double Sapling::representedStemNumber(float height) const |
96 | { |
||
97 | const SaplingGrowthParameters &p = mRUS->species()->saplingGrowthParameters(); |
||
98 | float dbh = height / p.hdSapling * 100.f; |
||
99 | double n = p.representedStemNumber(dbh); |
||
100 | return n; |
||
101 | } |
||
102 | |||
595 | werner | 103 | /// maintenance function to clear dead/recruited saplings from storage |
104 | void Sapling::cleanupStorage() |
||
105 | { |
||
1111 | werner | 106 | QVector<SaplingTreeOld>::iterator forw=mSaplingTrees.begin(); |
107 | QVector<SaplingTreeOld>::iterator back; |
||
595 | werner | 108 | |
109 | // seek last valid |
||
110 | for (back=mSaplingTrees.end()-1; back>=mSaplingTrees.begin(); --back) |
||
111 | if ((*back).isValid()) |
||
112 | break; |
||
113 | |||
114 | if (back<mSaplingTrees.begin()) { |
||
115 | mSaplingTrees.clear(); // no valid trees available |
||
116 | return; |
||
117 | } |
||
118 | |||
119 | while (forw < back) { |
||
120 | if (!(*forw).isValid()) { |
||
121 | *forw = *back; // copy (fill gap) |
||
122 | while (back>forw) // seek next valid |
||
123 | if ((*--back).isValid()) |
||
124 | break; |
||
125 | } |
||
126 | ++forw; |
||
127 | } |
||
128 | if (back != mSaplingTrees.end()-1) { |
||
129 | // free resources... |
||
130 | mSaplingTrees.erase(back+1, mSaplingTrees.end()); |
||
131 | } |
||
132 | } |
||
133 | |||
134 | // not a very good way of checking if sapling is present |
||
135 | // maybe better: use also a (local) maximum sapling height grid |
||
136 | // maybe better: use a bitset: |
||
137 | // position: index of pixel on LIF (absolute index) |
||
138 | bool Sapling::hasSapling(const QPoint &position) const |
||
139 | { |
||
140 | const QPoint &offset = mRUS->ru()->cornerPointOffset(); |
||
141 | int index = (position.x()- offset.x())*cPxPerRU + (position.y() - offset.y()); |
||
802 | werner | 142 | if (index<0) |
143 | qDebug() << "Sapling error"; |
||
595 | werner | 144 | return mSapBitset[index]; |
145 | /* |
||
146 | float *target = GlobalSettings::instance()->model()->grid()->ptr(position.x(), position.y()); |
||
147 | QVector<SaplingTree>::const_iterator it; |
||
148 | for (it = mSaplingTrees.constBegin(); it!=mSaplingTrees.constEnd(); ++it) { |
||
149 | if (it->pixel==target) |
||
150 | return true; |
||
151 | } |
||
152 | return false; |
||
153 | */ |
||
154 | } |
||
155 | |||
600 | werner | 156 | /// retrieve the height of the sapling at the location 'position' (given in LIF-coordinates) |
157 | /// this is quite expensive and only done for initialization |
||
158 | double Sapling::heightAt(const QPoint &position) const |
||
159 | { |
||
160 | if (!hasSapling(position)) |
||
161 | return 0.; |
||
162 | // ok, we'll have to search through all saplings |
||
1111 | werner | 163 | QVector<SaplingTreeOld>::const_iterator it; |
600 | werner | 164 | float *lif_ptr = GlobalSettings::instance()->model()->grid()->ptr(position.x(), position.y()); |
165 | for (it = mSaplingTrees.constBegin(); it!=mSaplingTrees.constEnd(); ++it) { |
||
166 | if (it->isValid() && it->pixel == lif_ptr) |
||
167 | return it->height; |
||
168 | } |
||
169 | return 0.; |
||
595 | werner | 170 | |
600 | werner | 171 | } |
172 | |||
173 | |||
941 | werner | 174 | void Sapling::setBit(const QPoint &pos_index, bool value) |
695 | werner | 175 | { |
176 | int index = (pos_index.x() - mRUS->ru()->cornerPointOffset().x()) * cPxPerRU +(pos_index.y() - mRUS->ru()->cornerPointOffset().y()); |
||
941 | werner | 177 | mSapBitset.set(index,value); // set bit: now there is a sapling there |
695 | werner | 178 | } |
179 | |||
901 | werner | 180 | /// add a sapling at the given position (index on the LIF grid, i.e. 2x2m) |
951 | werner | 181 | int Sapling::addSapling(const QPoint &pos_lif, const float height, const int age) |
595 | werner | 182 | { |
183 | // adds a sapling... |
||
1111 | werner | 184 | mSaplingTrees.push_back(SaplingTreeOld()); |
185 | SaplingTreeOld &t = mSaplingTrees.back(); |
||
911 | werner | 186 | t.height = height; // default is 5cm height |
951 | werner | 187 | t.age.age = age; |
595 | werner | 188 | Grid<float> &lif_map = *GlobalSettings::instance()->model()->grid(); |
189 | t.pixel = lif_map.ptr(pos_lif.x(), pos_lif.y()); |
||
941 | werner | 190 | setBit(pos_lif, true); |
595 | werner | 191 | mAdded++; |
911 | werner | 192 | return mSaplingTrees.count()-1; // index of the newly added tree. |
595 | werner | 193 | } |
194 | |||
195 | /// clear saplings on a given position (after recruitment) |
||
196 | void Sapling::clearSaplings(const QPoint &position) |
||
197 | { |
||
198 | float *target = GlobalSettings::instance()->model()->grid()->ptr(position.x(), position.y()); |
||
1111 | werner | 199 | QVector<SaplingTreeOld>::const_iterator it; |
595 | werner | 200 | for (it = mSaplingTrees.constBegin(); it!=mSaplingTrees.constEnd(); ++it) { |
201 | if (it->pixel==target) { |
||
202 | // trick: use a const iterator to avoid a deep copy of the vector; then do an ugly const_cast to actually write the data |
||
662 | werner | 203 | //const SaplingTree &t = *it; |
204 | //const_cast<SaplingTree&>(t).pixel=0; |
||
1111 | werner | 205 | clearSapling(const_cast<SaplingTreeOld&>(*it), false); // kill sapling and move carbon to soil |
595 | werner | 206 | } |
207 | } |
||
941 | werner | 208 | setBit(position, false); // clear bit: now there is no sapling on this position |
209 | //int index = (position.x() - mRUS->ru()->cornerPointOffset().x()) * cPxPerRU +(position.y() - mRUS->ru()->cornerPointOffset().y()); |
||
210 | //mSapBitset.set(index,false); // clear bit: now there is no sapling on this position |
||
595 | werner | 211 | |
212 | } |
||
213 | |||
662 | werner | 214 | /// clear saplings within a given rectangle |
215 | void Sapling::clearSaplings(const QRectF &rectangle, const bool remove_biomass) |
||
216 | { |
||
1111 | werner | 217 | QVector<SaplingTreeOld>::const_iterator it; |
662 | werner | 218 | FloatGrid *grid = GlobalSettings::instance()->model()->grid(); |
219 | for (it = mSaplingTrees.constBegin(); it!=mSaplingTrees.constEnd(); ++it) { |
||
902 | werner | 220 | if (rectangle.contains(grid->cellCenterPoint(it->coords()))) { |
1111 | werner | 221 | clearSapling(const_cast<SaplingTreeOld&>(*it), remove_biomass); |
662 | werner | 222 | } |
223 | } |
||
224 | } |
||
225 | |||
1111 | werner | 226 | void Sapling::clearSapling(SaplingTreeOld &tree, const bool remove) |
662 | werner | 227 | { |
941 | werner | 228 | QPoint p=tree.coords(); |
662 | werner | 229 | tree.pixel=0; |
941 | werner | 230 | setBit(p, false); // no tree left |
662 | werner | 231 | if (!remove) { |
232 | // killing of saplings: |
||
233 | // if remove=false, then remember dbh/number of trees (used later in calculateGrowth() to estimate carbon flow) |
||
234 | mDied++; |
||
235 | mSumDbhDied+=tree.height / mRUS->species()->saplingGrowthParameters().hdSapling * 100.; |
||
236 | } |
||
237 | |||
238 | } |
||
239 | |||
911 | werner | 240 | // |
241 | void Sapling::clearSapling(int index, const bool remove) |
||
242 | { |
||
243 | Q_ASSERT(index < mSaplingTrees.count()); |
||
244 | clearSapling(mSaplingTrees[index], remove); |
||
245 | } |
||
246 | |||
595 | werner | 247 | /// growth function for an indivudal sapling. |
248 | /// returns true, if sapling survives, false if sapling dies or is recruited to iLand. |
||
249 | /// see also http://iland.boku.ac.at/recruitment |
||
1111 | werner | 250 | bool Sapling::growSapling(SaplingTreeOld &tree, const double f_env_yr, Species* species) |
595 | werner | 251 | { |
252 | QPoint p=GlobalSettings::instance()->model()->grid()->indexOf(tree.pixel); |
||
1068 | werner | 253 | //GlobalSettings::instance()->model()->heightGrid()[Grid::index5(tree.pixel-GlobalSettings::instance()->model()->grid()->begin())]; |
595 | werner | 254 | |
255 | // (1) calculate height growth potential for the tree (uses linerization of expressions...) |
||
802 | werner | 256 | double h_pot = species->saplingGrowthParameters().heightGrowthPotential.calculate(tree.height); // TODO check if this can be source of crashes (race condition) |
595 | werner | 257 | double delta_h_pot = h_pot - tree.height; |
258 | |||
259 | // (2) reduce height growth potential with species growth response f_env_yr and with light state (i.e. LIF-value) of home-pixel. |
||
260 | double lif_value = *tree.pixel; |
||
261 | double h_height_grid = GlobalSettings::instance()->model()->heightGrid()->valueAtIndex(p.x()/cPxPerHeight, p.y()/cPxPerHeight).height; |
||
1157 | werner | 262 | if (h_height_grid==0.) |
595 | werner | 263 | throw IException(QString("growSapling: height grid at %1/%2 has value 0").arg(p.x()).arg(p.y())); |
264 | |||
816 | werner | 265 | double rel_height = tree.height / h_height_grid; |
266 | |||
595 | werner | 267 | double lif_corrected = mRUS->species()->speciesSet()->LRIcorrection(lif_value, rel_height); // correction based on height |
683 | werner | 268 | |
595 | werner | 269 | double lr = mRUS->species()->lightResponse(lif_corrected); // species specific light response (LUI, light utilization index) |
270 | |||
271 | double delta_h_factor = f_env_yr * lr; // relative growth |
||
272 | |||
273 | if (h_pot<0. || delta_h_pot<0. || lif_corrected<0. || lif_corrected>1. || delta_h_factor<0. || delta_h_factor>1. ) |
||
274 | qDebug() << "invalid values in Sapling::growSapling"; |
||
275 | |||
1063 | werner | 276 | // check browsing |
277 | if (mBrowsingPressure>0. && tree.height<=2.f) { |
||
278 | double p = mRUS->species()->saplingGrowthParameters().browsingProbability; |
||
279 | // calculate modifed annual browsing probability via odds-ratios |
||
280 | // odds = p/(1-p) -> odds_mod = odds * browsingPressure -> p_mod = odds_mod /( 1 + odds_mod) === p*pressure/(1-p+p*pressure) |
||
281 | double p_browse = p*mBrowsingPressure / (1. - p + p*mBrowsingPressure); |
||
282 | if (drandom() < p_browse) { |
||
283 | delta_h_factor = 0.; |
||
284 | } |
||
285 | } |
||
286 | |||
595 | werner | 287 | // check mortality of saplings |
288 | if (delta_h_factor < species->saplingGrowthParameters().stressThreshold) { |
||
289 | tree.age.stress_years++; |
||
290 | if (tree.age.stress_years > species->saplingGrowthParameters().maxStressYears) { |
||
291 | // sapling dies... |
||
662 | werner | 292 | clearSapling(tree, false); // false: put carbon to the soil |
595 | werner | 293 | return false; |
294 | } |
||
295 | } else { |
||
296 | tree.age.stress_years=0; // reset stress counter |
||
297 | } |
||
298 | DBG_IF(delta_h_pot*delta_h_factor < 0.f || delta_h_pot*delta_h_factor > 2., "Sapling::growSapling", "inplausible height growth."); |
||
299 | |||
300 | // grow |
||
301 | tree.height += delta_h_pot * delta_h_factor; |
||
302 | tree.age.age++; // increase age of sapling by 1 |
||
303 | |||
304 | // recruitment? |
||
305 | if (tree.height > 4.f) { |
||
306 | mRecruited++; |
||
307 | |||
308 | ResourceUnit *ru = const_cast<ResourceUnit*> (mRUS->ru()); |
||
309 | float dbh = tree.height / species->saplingGrowthParameters().hdSapling * 100.f; |
||
310 | // the number of trees to create (result is in trees per pixel) |
||
311 | double n_trees = species->saplingGrowthParameters().representedStemNumber(dbh); |
||
312 | int to_establish = (int) n_trees; |
||
863 | werner | 313 | |
595 | werner | 314 | // if n_trees is not an integer, choose randomly if we should add a tree. |
315 | // e.g.: n_trees = 2.3 -> add 2 trees with 70% probability, and add 3 trees with p=30%. |
||
707 | werner | 316 | if (drandom() < (n_trees-to_establish) || to_establish==0) |
595 | werner | 317 | to_establish++; |
318 | |||
319 | // add a new tree |
||
320 | for (int i=0;i<to_establish;i++) { |
||
321 | Tree &bigtree = ru->newTree(); |
||
322 | bigtree.setPosition(p); |
||
323 | // add variation: add +/-10% to dbh and *independently* to height. |
||
707 | werner | 324 | bigtree.setDbh(dbh * nrandom(1. - mRecruitmentVariation, 1. + mRecruitmentVariation)); |
325 | bigtree.setHeight(tree.height * nrandom(1. - mRecruitmentVariation, 1. + mRecruitmentVariation)); |
||
595 | werner | 326 | bigtree.setSpecies( species ); |
327 | bigtree.setAge(tree.age.age,tree.height); |
||
328 | bigtree.setRU(ru); |
||
329 | bigtree.setup(); |
||
855 | werner | 330 | const Tree *t = &bigtree; |
331 | mRUS->statistics().add(t, 0); // count the newly created trees already in the stats |
||
595 | werner | 332 | } |
333 | // clear all regeneration from this pixel (including this tree) |
||
662 | werner | 334 | clearSapling(tree, true); // remove this tree (but do not move biomass to soil) |
335 | ru->clearSaplings(p); // remove all other saplings on the same pixel |
||
855 | werner | 336 | |
595 | werner | 337 | return false; |
338 | } |
||
339 | // book keeping (only for survivors) |
||
340 | mLiving++; |
||
341 | mAvgHeight+=tree.height; |
||
342 | mAvgAge+=tree.age.age; |
||
343 | mAvgDeltaHPot+=delta_h_pot; |
||
344 | mAvgHRealized += delta_h_pot * delta_h_factor; |
||
345 | return true; |
||
346 | |||
347 | } |
||
348 | |||
608 | werner | 349 | |
662 | werner | 350 | /** main growth function for saplings. |
351 | Statistics are cleared at the beginning of the year. |
||
352 | */ |
||
595 | werner | 353 | void Sapling::calculateGrowth() |
354 | { |
||
355 | Q_ASSERT(mRUS); |
||
821 | werner | 356 | if (mSaplingTrees.count()==0) |
595 | werner | 357 | return; |
358 | |||
359 | ResourceUnit *ru = const_cast<ResourceUnit*> (mRUS->ru() ); |
||
360 | Species *species = const_cast<Species*>(mRUS->species()); |
||
361 | |||
362 | // calculate necessary growth modifier (this is done only once per year) |
||
363 | mRUS->calculate(true); // calculate the 3pg module (this is done only if that did not happen up to now); true: call comes from regeneration |
||
364 | double f_env_yr = mRUS->prod3PG().fEnvYear(); |
||
365 | |||
366 | mLiving=0; |
||
1111 | werner | 367 | QVector<SaplingTreeOld>::const_iterator it; |
595 | werner | 368 | for (it = mSaplingTrees.constBegin(); it!=mSaplingTrees.constEnd(); ++it) { |
1111 | werner | 369 | const SaplingTreeOld &tree = *it; |
595 | werner | 370 | if (tree.height<0) |
371 | qDebug() << "Sapling::calculateGrowth(): h<0"; |
||
821 | werner | 372 | // if sapling is still living check execute growth routine |
595 | werner | 373 | if (tree.isValid()) { |
821 | werner | 374 | // growing (increases mLiving if tree did not die, mDied otherwise) |
1111 | werner | 375 | if (growSapling(const_cast<SaplingTreeOld&>(tree), f_env_yr, species)) { |
595 | werner | 376 | // set the sapling height to the maximum value on the current pixel |
902 | werner | 377 | ru->setMaxSaplingHeightAt(tree.coords(),tree.height); |
595 | werner | 378 | } |
379 | } |
||
380 | } |
||
381 | if (mLiving) { |
||
382 | mAvgHeight /= double(mLiving); |
||
383 | mAvgAge /= double(mLiving); |
||
384 | mAvgDeltaHPot /= double(mLiving); |
||
385 | mAvgHRealized /= double(mLiving); |
||
386 | } |
||
387 | // calculate carbon balance |
||
608 | werner | 388 | CNPair old_state = mCarbonLiving; |
595 | werner | 389 | mCarbonLiving.clear(); |
608 | werner | 390 | |
595 | werner | 391 | CNPair dead_wood, dead_fine; // pools for mortality |
392 | // average dbh |
||
393 | if (mLiving) { |
||
1157 | werner | 394 | // calculate the avg dbh and number of stems |
595 | werner | 395 | double avg_dbh = mAvgHeight / species->saplingGrowthParameters().hdSapling * 100.; |
396 | double n = mLiving * species->saplingGrowthParameters().representedStemNumber( avg_dbh ); |
||
397 | // woody parts: stem, branchse and coarse roots |
||
398 | double woody_bm = species->biomassWoody(avg_dbh) + species->biomassBranch(avg_dbh) + species->biomassRoot(avg_dbh); |
||
399 | double foliage = species->biomassFoliage(avg_dbh); |
||
400 | double fineroot = foliage*species->finerootFoliageRatio(); |
||
401 | |||
402 | mCarbonLiving.addBiomass( woody_bm*n, species->cnWood() ); |
||
403 | mCarbonLiving.addBiomass( foliage*n, species->cnFoliage() ); |
||
404 | mCarbonLiving.addBiomass( fineroot*n, species->cnFineroot() ); |
||
608 | werner | 405 | |
595 | werner | 406 | // turnover |
802 | werner | 407 | if (mRUS->ru()->snag()) |
408 | mRUS->ru()->snag()->addTurnoverLitter(species, foliage*species->turnoverLeaf(), fineroot*species->turnoverRoot()); |
||
595 | werner | 409 | |
410 | // calculate the "mortality from competition", i.e. carbon that stems from reduction of stem numbers |
||
411 | // from Reinekes formula. |
||
412 | // |
||
413 | if (avg_dbh>1.) { |
||
414 | double avg_dbh_before = (mAvgHeight - mAvgHRealized) / species->saplingGrowthParameters().hdSapling * 100.; |
||
415 | double n_before = mLiving * species->saplingGrowthParameters().representedStemNumber( qMax(1.,avg_dbh_before) ); |
||
416 | if (n<n_before) { |
||
417 | dead_wood.addBiomass( woody_bm * (n_before-n), species->cnWood() ); |
||
418 | dead_fine.addBiomass( foliage * (n_before-n), species->cnFoliage() ); |
||
419 | dead_fine.addBiomass( fineroot * (n_before-n), species->cnFineroot() ); |
||
420 | } |
||
421 | } |
||
422 | |||
423 | } |
||
424 | if (mDied) { |
||
425 | double avg_dbh_dead = mSumDbhDied / double(mDied); |
||
426 | double n = mDied * species->saplingGrowthParameters().representedStemNumber( avg_dbh_dead ); |
||
427 | // woody parts: stem, branchse and coarse roots |
||
428 | |||
429 | dead_wood.addBiomass( ( species->biomassWoody(avg_dbh_dead) + species->biomassBranch(avg_dbh_dead) + species->biomassRoot(avg_dbh_dead)) * n, species->cnWood() ); |
||
430 | double foliage = species->biomassFoliage(avg_dbh_dead)*n; |
||
431 | |||
432 | dead_fine.addBiomass( foliage, species->cnFoliage() ); |
||
433 | dead_fine.addBiomass( foliage*species->finerootFoliageRatio(), species->cnFineroot() ); |
||
434 | } |
||
435 | if (!dead_wood.isEmpty() || !dead_fine.isEmpty()) |
||
802 | werner | 436 | if (mRUS->ru()->snag()) |
437 | mRUS->ru()->snag()->addToSoil(species, dead_wood, dead_fine); |
||
595 | werner | 438 | |
608 | werner | 439 | // calculate net growth: |
440 | // delta of stocks |
||
625 | werner | 441 | mCarbonGain = mCarbonLiving + dead_fine + dead_wood - old_state; |
442 | if (mCarbonGain.C < 0) |
||
443 | mCarbonGain.clear(); |
||
608 | werner | 444 | |
595 | werner | 445 | if (mSaplingTrees.count() > mLiving*1.3) |
446 | cleanupStorage(); |
||
447 | |||
448 | mRUS->statistics().add(this); |
||
615 | werner | 449 | GlobalSettings::instance()->systemStatistics()->saplingCount+=mLiving; |
450 | GlobalSettings::instance()->systemStatistics()->newSaplings+=mAdded; |
||
663 | werner | 451 | mAdded = 0; // reset |
452 | |||
595 | werner | 453 | //qDebug() << ru->index() << species->id()<< ": (living/avg.height):" << mLiving << mAvgHeight; |
454 | } |
||
455 | |||
456 | /// fill a grid with the maximum height of saplings per pixel (2x2m). |
||
457 | /// this function is used for visualization only |
||
458 | void Sapling::fillMaxHeightGrid(Grid<float> &grid) const |
||
459 | { |
||
1111 | werner | 460 | QVector<SaplingTreeOld>::const_iterator it; |
595 | werner | 461 | for (it = mSaplingTrees.begin(); it!=mSaplingTrees.end(); ++it) { |
462 | if (it->isValid()) { |
||
902 | werner | 463 | QPoint p=it->coords(); |
595 | werner | 464 | if (grid.valueAtIndex(p)<it->height) |
465 | grid.valueAtIndex(p) = it->height; |
||
466 | } |
||
467 | } |
||
468 | |||
469 | } |
||
600 | werner | 470 | |
608 | werner | 471 | |
472 | |||
695 | werner | 473 |