Rev 639 | Rev 664 | Go to most recent revision | Only display areas with differences | Ignore whitespace | Details | Blame | Last modification | View Log | RSS feed
Rev 639 | Rev 662 | ||
---|---|---|---|
1 | Redirecting to URL 'https://iland.boku.ac.at/svn/iland/tags/release_1.0/src/core/resourceunit.cpp': |
1 | Redirecting to URL 'https://iland.boku.ac.at/svn/iland/tags/release_1.0/src/core/resourceunit.cpp': |
2 | /** @class ResourceUnit
|
2 | /** @class ResourceUnit
|
3 | ResourceUnit is the spatial unit that encapsulates a forest stand and links to several environmental components
|
3 | ResourceUnit is the spatial unit that encapsulates a forest stand and links to several environmental components
|
4 | (Climate, Soil, Water, ...).
|
4 | (Climate, Soil, Water, ...).
|
5 | 5 | ||
6 | */
|
6 | */
|
7 | #include <QtCore>
|
7 | #include <QtCore>
|
8 | #include "global.h"
|
8 | #include "global.h"
|
9 | 9 | ||
10 | #include "resourceunit.h"
|
10 | #include "resourceunit.h"
|
11 | #include "resourceunitspecies.h"
|
11 | #include "resourceunitspecies.h"
|
12 | #include "speciesset.h"
|
12 | #include "speciesset.h"
|
13 | #include "species.h"
|
13 | #include "species.h"
|
14 | #include "production3pg.h"
|
14 | #include "production3pg.h"
|
15 | #include "model.h"
|
15 | #include "model.h"
|
16 | #include "climate.h"
|
16 | #include "climate.h"
|
17 | #include "watercycle.h"
|
17 | #include "watercycle.h"
|
18 | #include "snag.h"
|
18 | #include "snag.h"
|
19 | #include "soil.h"
|
19 | #include "soil.h"
|
20 | #include "helper.h"
|
20 | #include "helper.h"
|
21 | 21 | ||
22 | ResourceUnit::~ResourceUnit() |
22 | ResourceUnit::~ResourceUnit() |
23 | {
|
23 | {
|
24 | if (mWater) |
24 | if (mWater) |
25 | delete mWater; |
25 | delete mWater; |
26 | mWater = 0; |
26 | mWater = 0; |
27 | if (mSnag) |
27 | if (mSnag) |
28 | delete mSnag; |
28 | delete mSnag; |
29 | if (mSoil) |
29 | if (mSoil) |
30 | delete mSoil; |
30 | delete mSoil; |
31 | 31 | ||
32 | mSnag = 0; |
32 | mSnag = 0; |
33 | mSoil = 0; |
33 | mSoil = 0; |
34 | }
|
34 | }
|
35 | 35 | ||
36 | ResourceUnit::ResourceUnit(const int index) |
36 | ResourceUnit::ResourceUnit(const int index) |
37 | {
|
37 | {
|
38 | qDeleteAll(mRUSpecies); |
38 | qDeleteAll(mRUSpecies); |
39 | mSpeciesSet = 0; |
39 | mSpeciesSet = 0; |
40 | mClimate = 0; |
40 | mClimate = 0; |
41 | mPixelCount=0; |
41 | mPixelCount=0; |
42 | mStockedArea = 0; |
42 | mStockedArea = 0; |
43 | mStockedPixelCount = 0; |
43 | mStockedPixelCount = 0; |
44 | mIndex = index; |
44 | mIndex = index; |
45 | mSaplingHeightMap = 0; |
45 | mSaplingHeightMap = 0; |
46 | mEffectiveArea_perWLA = 0.; |
46 | mEffectiveArea_perWLA = 0.; |
47 | mWater = new WaterCycle(); |
47 | mWater = new WaterCycle(); |
48 | mSnag = 0; |
48 | mSnag = 0; |
49 | mSoil = 0; |
49 | mSoil = 0; |
50 | mID = 0; |
50 | mID = 0; |
51 | }
|
51 | }
|
52 | 52 | ||
53 | void ResourceUnit::setup() |
53 | void ResourceUnit::setup() |
54 | {
|
54 | {
|
55 | mWater->setup(this); |
55 | mWater->setup(this); |
56 | 56 | ||
57 | if (mSnag) |
57 | if (mSnag) |
58 | delete mSnag; |
58 | delete mSnag; |
59 | mSnag=0; |
59 | mSnag=0; |
60 | if (mSoil) |
60 | if (mSoil) |
61 | delete mSoil; |
61 | delete mSoil; |
62 | mSoil=0; |
62 | mSoil=0; |
63 | if (Model::settings().carbonCycleEnabled) { |
63 | if (Model::settings().carbonCycleEnabled) { |
64 | mSoil = new Soil(this); |
64 | mSoil = new Soil(this); |
65 | mSnag = new Snag; |
65 | mSnag = new Snag; |
66 | mSnag->setup(this); |
66 | mSnag->setup(this); |
67 | const XmlHelper &xml=GlobalSettings::instance()->settings(); |
67 | const XmlHelper &xml=GlobalSettings::instance()->settings(); |
68 | 68 | ||
69 | // setup contents of the soil of the RU; use values for C and N (kg/ha)
|
69 | // setup contents of the soil of the RU; use values for C and N (kg/ha)
|
70 | mSoil->setInitialState(CNPool(xml.valueDouble("model.site.youngLabileC", -1), |
70 | mSoil->setInitialState(CNPool(xml.valueDouble("model.site.youngLabileC", -1), |
71 | xml.valueDouble("model.site.youngLabileN", -1), |
71 | xml.valueDouble("model.site.youngLabileN", -1), |
72 | xml.valueDouble("model.site.youngLabileDecompRate", -1)), |
72 | xml.valueDouble("model.site.youngLabileDecompRate", -1)), |
73 | CNPool(xml.valueDouble("model.site.youngRefractoryC", -1), |
73 | CNPool(xml.valueDouble("model.site.youngRefractoryC", -1), |
74 | xml.valueDouble("model.site.youngRefractoryN", -1), |
74 | xml.valueDouble("model.site.youngRefractoryN", -1), |
75 | xml.valueDouble("model.site.youngRefractoryDecompRate", -1)), |
75 | xml.valueDouble("model.site.youngRefractoryDecompRate", -1)), |
76 | CNPair(xml.valueDouble("model.site.somC", -1), xml.valueDouble("model.site.somN", -1))); |
76 | CNPair(xml.valueDouble("model.site.somC", -1), xml.valueDouble("model.site.somN", -1))); |
77 | }
|
77 | }
|
78 | 78 | ||
79 | // setup variables
|
79 | // setup variables
|
80 | mUnitVariables.nitrogenAvailable = GlobalSettings::instance()->settings().valueDouble("model.site.availableNitrogen", 40); |
80 | mUnitVariables.nitrogenAvailable = GlobalSettings::instance()->settings().valueDouble("model.site.availableNitrogen", 40); |
81 | 81 | ||
82 | // if dynamic coupling of soil nitrogen is enabled, the calculate a starting value for available n.
|
82 | // if dynamic coupling of soil nitrogen is enabled, the calculate a starting value for available n.
|
83 | if (mSoil && Model::settings().useDynamicAvailableNitrogen && Model::settings().carbonCycleEnabled) { |
83 | if (mSoil && Model::settings().useDynamicAvailableNitrogen && Model::settings().carbonCycleEnabled) { |
84 | mSoil->setClimateFactor(1.); |
84 | mSoil->setClimateFactor(1.); |
85 | mSoil->calculateYear(); |
85 | mSoil->calculateYear(); |
86 | mUnitVariables.nitrogenAvailable = mSoil->availableNitrogen(); |
86 | mUnitVariables.nitrogenAvailable = mSoil->availableNitrogen(); |
87 | }
|
87 | }
|
88 | mAverageAging = 0.; |
88 | mAverageAging = 0.; |
89 | 89 | ||
90 | }
|
90 | }
|
91 | void ResourceUnit::setBoundingBox(const QRectF &bb) |
91 | void ResourceUnit::setBoundingBox(const QRectF &bb) |
92 | {
|
92 | {
|
93 | mBoundingBox = bb; |
93 | mBoundingBox = bb; |
94 | mCornerCoord = GlobalSettings::instance()->model()->grid()->indexAt(bb.topLeft()); |
94 | mCornerCoord = GlobalSettings::instance()->model()->grid()->indexAt(bb.topLeft()); |
95 | }
|
95 | }
|
96 | 96 | ||
97 | /// set species and setup the species-per-RU-data
|
97 | /// set species and setup the species-per-RU-data
|
98 | void ResourceUnit::setSpeciesSet(SpeciesSet *set) |
98 | void ResourceUnit::setSpeciesSet(SpeciesSet *set) |
99 | {
|
99 | {
|
100 | mSpeciesSet = set; |
100 | mSpeciesSet = set; |
101 | qDeleteAll(mRUSpecies); |
101 | qDeleteAll(mRUSpecies); |
102 | 102 | ||
103 | //mRUSpecies.resize(set->count()); // ensure that the vector space is not relocated
|
103 | //mRUSpecies.resize(set->count()); // ensure that the vector space is not relocated
|
104 | for (int i=0;i<set->count();i++) { |
104 | for (int i=0;i<set->count();i++) { |
105 | Species *s = const_cast<Species*>(mSpeciesSet->species(i)); |
105 | Species *s = const_cast<Species*>(mSpeciesSet->species(i)); |
106 | if (!s) |
106 | if (!s) |
107 | throw IException("ResourceUnit::setSpeciesSet: invalid index!"); |
107 | throw IException("ResourceUnit::setSpeciesSet: invalid index!"); |
108 | 108 | ||
109 | ResourceUnitSpecies *rus = new ResourceUnitSpecies(); |
109 | ResourceUnitSpecies *rus = new ResourceUnitSpecies(); |
110 | mRUSpecies.push_back(rus); |
110 | mRUSpecies.push_back(rus); |
111 | rus->setup(s, this); |
111 | rus->setup(s, this); |
112 | /* be careful: setup() is called with a pointer somewhere to the content of the mRUSpecies container.
|
112 | /* be careful: setup() is called with a pointer somewhere to the content of the mRUSpecies container.
|
113 | If the container memory is relocated (QVector), the pointer gets invalid!!!
|
113 | If the container memory is relocated (QVector), the pointer gets invalid!!!
|
114 | Therefore, a resize() is called before the loop (no resize()-operations during the loop)! */
|
114 | Therefore, a resize() is called before the loop (no resize()-operations during the loop)! */
|
115 | //mRUSpecies[i].setup(s,this); // setup this element
|
115 | //mRUSpecies[i].setup(s,this); // setup this element
|
116 | 116 | ||
117 | }
|
117 | }
|
118 | }
|
118 | }
|
119 | 119 | ||
120 | ResourceUnitSpecies &ResourceUnit::resourceUnitSpecies(const Species *species) |
120 | ResourceUnitSpecies &ResourceUnit::resourceUnitSpecies(const Species *species) |
121 | {
|
121 | {
|
122 | return *mRUSpecies[species->index()]; |
122 | return *mRUSpecies[species->index()]; |
123 | }
|
123 | }
|
124 | 124 | ||
125 | Tree &ResourceUnit::newTree() |
125 | Tree &ResourceUnit::newTree() |
126 | {
|
126 | {
|
127 | // start simple: just append to the vector...
|
127 | // start simple: just append to the vector...
|
128 | if (mTrees.isEmpty()) |
128 | if (mTrees.isEmpty()) |
129 | mTrees.reserve(100); // reserve a junk of memory for trees |
129 | mTrees.reserve(100); // reserve a junk of memory for trees |
130 | 130 | ||
131 | mTrees.append(Tree()); |
131 | mTrees.append(Tree()); |
132 | return mTrees.back(); |
132 | return mTrees.back(); |
133 | }
|
133 | }
|
134 | int ResourceUnit::newTreeIndex() |
134 | int ResourceUnit::newTreeIndex() |
135 | {
|
135 | {
|
136 | // start simple: just append to the vector...
|
136 | // start simple: just append to the vector...
|
137 | mTrees.append(Tree()); |
137 | mTrees.append(Tree()); |
138 | return mTrees.count()-1; |
138 | return mTrees.count()-1; |
139 | }
|
139 | }
|
140 | 140 | ||
141 | /// remove dead trees from tree list
|
141 | /// remove dead trees from tree list
|
142 | /// reduce size of vector if lots of space is free
|
142 | /// reduce size of vector if lots of space is free
|
143 | /// tests showed that this way of cleanup is very fast,
|
143 | /// tests showed that this way of cleanup is very fast,
|
144 | /// because no memory allocations are performed (simple memmove())
|
144 | /// because no memory allocations are performed (simple memmove())
|
145 | /// when trees are moved.
|
145 | /// when trees are moved.
|
146 | void ResourceUnit::cleanTreeList() |
146 | void ResourceUnit::cleanTreeList() |
147 | {
|
147 | {
|
148 | QVector<Tree>::iterator last=mTrees.end()-1; |
148 | QVector<Tree>::iterator last=mTrees.end()-1; |
149 | QVector<Tree>::iterator current = mTrees.begin(); |
149 | QVector<Tree>::iterator current = mTrees.begin(); |
150 | while (last>=current && (*last).isDead()) |
150 | while (last>=current && (*last).isDead()) |
151 | --last; |
151 | --last; |
152 | 152 | ||
153 | while (current<last) { |
153 | while (current<last) { |
154 | if ((*current).isDead()) { |
154 | if ((*current).isDead()) { |
155 | *current = *last; // copy data! |
155 | *current = *last; // copy data! |
156 | --last; // |
156 | --last; // |
157 | while (last>=current && (*last).isDead()) |
157 | while (last>=current && (*last).isDead()) |
158 | --last; |
158 | --last; |
159 | }
|
159 | }
|
160 | ++current; |
160 | ++current; |
161 | }
|
161 | }
|
162 | ++last; // last points now to the first dead tree |
162 | ++last; // last points now to the first dead tree |
163 | 163 | ||
164 | // free ressources
|
164 | // free ressources
|
165 | if (last!=mTrees.end()) { |
165 | if (last!=mTrees.end()) { |
166 | mTrees.erase(last, mTrees.end()); |
166 | mTrees.erase(last, mTrees.end()); |
167 | if (mTrees.capacity()>100) { |
167 | if (mTrees.capacity()>100) { |
168 | if (mTrees.count() / double(mTrees.capacity()) < 0.2) { |
168 | if (mTrees.count() / double(mTrees.capacity()) < 0.2) { |
169 | //int target_size = mTrees.count()*2;
|
169 | //int target_size = mTrees.count()*2;
|
170 | //qDebug() << "reduce size from "<<mTrees.capacity() << "to" << target_size;
|
170 | //qDebug() << "reduce size from "<<mTrees.capacity() << "to" << target_size;
|
171 | //mTrees.reserve(qMax(target_size, 100));
|
171 | //mTrees.reserve(qMax(target_size, 100));
|
172 | qDebug() << "reduce tree storage of RU" << index() << " from " << mTrees.capacity() << "to" << mTrees.count(); |
172 | qDebug() << "reduce tree storage of RU" << index() << " from " << mTrees.capacity() << "to" << mTrees.count(); |
173 | mTrees.squeeze(); |
173 | mTrees.squeeze(); |
174 | }
|
174 | }
|
175 | }
|
175 | }
|
176 | }
|
176 | }
|
177 | }
|
177 | }
|
178 | 178 | ||
179 | void ResourceUnit::newYear() |
179 | void ResourceUnit::newYear() |
180 | {
|
180 | {
|
181 | mAggregatedWLA = 0.; |
181 | mAggregatedWLA = 0.; |
182 | mAggregatedLA = 0.; |
182 | mAggregatedLA = 0.; |
183 | mAggregatedLR = 0.; |
183 | mAggregatedLR = 0.; |
184 | mEffectiveArea = 0.; |
184 | mEffectiveArea = 0.; |
185 | mPixelCount = mStockedPixelCount = 0; |
185 | mPixelCount = mStockedPixelCount = 0; |
186 | snagNewYear(); |
186 | snagNewYear(); |
187 | if (mSoil) |
187 | if (mSoil) |
188 | mSoil->newYear(); |
188 | mSoil->newYear(); |
189 | // clear statistics global and per species...
|
189 | // clear statistics global and per species...
|
190 | QList<ResourceUnitSpecies*>::const_iterator i; |
190 | QList<ResourceUnitSpecies*>::const_iterator i; |
191 | QList<ResourceUnitSpecies*>::const_iterator iend = mRUSpecies.constEnd(); |
191 | QList<ResourceUnitSpecies*>::const_iterator iend = mRUSpecies.constEnd(); |
192 | mStatistics.clear(); |
192 | mStatistics.clear(); |
193 | for (i=mRUSpecies.constBegin(); i!=iend; ++i) { |
193 | for (i=mRUSpecies.constBegin(); i!=iend; ++i) { |
194 | (*i)->statisticsDead().clear(); |
194 | (*i)->statisticsDead().clear(); |
195 | (*i)->statisticsMgmt().clear(); |
195 | (*i)->statisticsMgmt().clear(); |
- | 196 | (*i)->changeSapling().newYear(); |
|
196 | }
|
197 | }
|
197 | 198 | ||
198 | }
|
199 | }
|
199 | 200 | ||
200 | /** production() is the "stand-level" part of the biomass production (3PG).
|
201 | /** production() is the "stand-level" part of the biomass production (3PG).
|
201 | - The amount of radiation intercepted by the stand is calculated
|
202 | - The amount of radiation intercepted by the stand is calculated
|
202 | - the water cycle is calculated
|
203 | - the water cycle is calculated
|
203 | - statistics for each species are cleared
|
204 | - statistics for each species are cleared
|
204 | - The 3PG production for each species and ressource unit is called (calculates species-responses and NPP production)
|
205 | - The 3PG production for each species and ressource unit is called (calculates species-responses and NPP production)
|
205 | see also: http://iland.boku.ac.at/individual+tree+light+availability */
|
206 | see also: http://iland.boku.ac.at/individual+tree+light+availability */
|
206 | void ResourceUnit::production() |
207 | void ResourceUnit::production() |
207 | {
|
208 | {
|
208 | 209 | ||
209 | if (mAggregatedWLA==0 || mPixelCount==0) { |
210 | if (mAggregatedWLA==0 || mPixelCount==0) { |
210 | // nothing to do...
|
211 | // nothing to do...
|
211 | return; |
212 | return; |
212 | }
|
213 | }
|
213 | 214 | ||
214 | // the pixel counters are filled during the height-grid-calculations
|
215 | // the pixel counters are filled during the height-grid-calculations
|
215 | mStockedArea = 100. * mStockedPixelCount; // m2 (1 height grid pixel = 10x10m) |
216 | mStockedArea = 100. * mStockedPixelCount; // m2 (1 height grid pixel = 10x10m) |
216 | 217 | ||
217 | // calculate the leaf area index (LAI)
|
218 | // calculate the leaf area index (LAI)
|
218 | double LAI = mAggregatedLA / mStockedArea; |
219 | double LAI = mAggregatedLA / mStockedArea; |
219 | // calculate the intercepted radiation fraction using the law of Beer Lambert
|
220 | // calculate the intercepted radiation fraction using the law of Beer Lambert
|
220 | const double k = Model::settings().lightExtinctionCoefficient; |
221 | const double k = Model::settings().lightExtinctionCoefficient; |
221 | double interception_fraction = 1. - exp(-k * LAI); |
222 | double interception_fraction = 1. - exp(-k * LAI); |
222 | mEffectiveArea = mStockedArea * interception_fraction; // m2 |
223 | mEffectiveArea = mStockedArea * interception_fraction; // m2 |
223 | 224 | ||
224 | // calculate the total weighted leaf area on this RU:
|
225 | // calculate the total weighted leaf area on this RU:
|
225 | mLRI_modification = interception_fraction * mStockedArea / mAggregatedWLA; // p_WLA |
226 | mLRI_modification = interception_fraction * mStockedArea / mAggregatedWLA; // p_WLA |
226 | if (mLRI_modification == 0.) |
227 | if (mLRI_modification == 0.) |
227 | qDebug() << "lri modifaction==0!"; |
228 | qDebug() << "lri modifaction==0!"; |
228 | 229 | ||
229 | if (logLevelDebug()) { |
230 | if (logLevelDebug()) { |
230 | DBGMODE(qDebug() << QString("production: LAI: %1 (intercepted fraction: %2, stocked area: %4). LRI-Multiplier: %3") |
231 | DBGMODE(qDebug() << QString("production: LAI: %1 (intercepted fraction: %2, stocked area: %4). LRI-Multiplier: %3") |
231 | .arg(LAI) |
232 | .arg(LAI) |
232 | .arg(interception_fraction) |
233 | .arg(interception_fraction) |
233 | .arg(mLRI_modification) |
234 | .arg(mLRI_modification) |
234 | .arg(mStockedArea); |
235 | .arg(mStockedArea); |
235 | ); |
236 | ); |
236 | }
|
237 | }
|
237 | 238 | ||
238 | // calculate LAI fractions
|
239 | // calculate LAI fractions
|
239 | QList<ResourceUnitSpecies*>::const_iterator i; |
240 | QList<ResourceUnitSpecies*>::const_iterator i; |
240 | QList<ResourceUnitSpecies*>::const_iterator iend = mRUSpecies.constEnd(); |
241 | QList<ResourceUnitSpecies*>::const_iterator iend = mRUSpecies.constEnd(); |
241 | double ru_lai = leafAreaIndex(); |
242 | double ru_lai = leafAreaIndex(); |
242 | if (ru_lai < 1.) |
243 | if (ru_lai < 1.) |
243 | ru_lai = 1.; |
244 | ru_lai = 1.; |
244 | // note: LAIFactors are only 1 if sum of LAI is > 1. (see WaterCycle)
|
245 | // note: LAIFactors are only 1 if sum of LAI is > 1. (see WaterCycle)
|
245 | for (i=mRUSpecies.constBegin(); i!=iend; ++i) { |
246 | for (i=mRUSpecies.constBegin(); i!=iend; ++i) { |
246 | (*i)->setLAIfactor((*i)->statistics().leafAreaIndex() / ru_lai); |
247 | (*i)->setLAIfactor((*i)->statistics().leafAreaIndex() / ru_lai); |
247 | }
|
248 | }
|
248 | 249 | ||
249 | // soil water model - this determines soil water contents needed for response calculations
|
250 | // soil water model - this determines soil water contents needed for response calculations
|
250 | {
|
251 | {
|
251 | mWater->run(); |
252 | mWater->run(); |
252 | }
|
253 | }
|
253 | 254 | ||
254 | // invoke species specific calculation (3PG)
|
255 | // invoke species specific calculation (3PG)
|
255 | for (i=mRUSpecies.constBegin(); i!=iend; ++i) { |
256 | for (i=mRUSpecies.constBegin(); i!=iend; ++i) { |
256 | (*i)->calculate(); // CALCULATE 3PG |
257 | (*i)->calculate(); // CALCULATE 3PG |
257 | if (logLevelInfo() && (*i)->LAIfactor()>0) |
258 | if (logLevelInfo() && (*i)->LAIfactor()>0) |
258 | qDebug() << "ru" << mIndex << "species" << (*i)->species()->id() << "LAIfraction" << (*i)->LAIfactor() << "raw_gpp_m2" |
259 | qDebug() << "ru" << mIndex << "species" << (*i)->species()->id() << "LAIfraction" << (*i)->LAIfactor() << "raw_gpp_m2" |
259 | << (*i)->prod3PG().GPPperArea() << "area:" << productiveArea() << "gpp:" |
260 | << (*i)->prod3PG().GPPperArea() << "area:" << productiveArea() << "gpp:" |
260 | << productiveArea()*(*i)->prod3PG().GPPperArea() |
261 | << productiveArea()*(*i)->prod3PG().GPPperArea() |
261 | << "aging(lastyear):" << averageAging() << "f_env,yr:" << (*i)->prod3PG().fEnvYear(); |
262 | << "aging(lastyear):" << averageAging() << "f_env,yr:" << (*i)->prod3PG().fEnvYear(); |
262 | }
|
263 | }
|
263 | }
|
264 | }
|
264 | 265 | ||
265 | void ResourceUnit::calculateInterceptedArea() |
266 | void ResourceUnit::calculateInterceptedArea() |
266 | {
|
267 | {
|
267 | if (mAggregatedLR==0) { |
268 | if (mAggregatedLR==0) { |
268 | mEffectiveArea_perWLA = 0.; |
269 | mEffectiveArea_perWLA = 0.; |
269 | return; |
270 | return; |
270 | }
|
271 | }
|
271 | Q_ASSERT(mAggregatedLR>0.); |
272 | Q_ASSERT(mAggregatedLR>0.); |
272 | mEffectiveArea_perWLA = mEffectiveArea / mAggregatedLR; |
273 | mEffectiveArea_perWLA = mEffectiveArea / mAggregatedLR; |
273 | if (logLevelDebug()) qDebug() << "RU: aggregated lightresponse:" << mAggregatedLR << "eff.area./wla:" << mEffectiveArea_perWLA; |
274 | if (logLevelDebug()) qDebug() << "RU: aggregated lightresponse:" << mAggregatedLR << "eff.area./wla:" << mEffectiveArea_perWLA; |
274 | }
|
275 | }
|
275 | 276 | ||
276 | // function is called immediately before the growth of individuals
|
277 | // function is called immediately before the growth of individuals
|
277 | void ResourceUnit::beforeGrow() |
278 | void ResourceUnit::beforeGrow() |
278 | {
|
279 | {
|
279 | mAverageAging = 0.; |
280 | mAverageAging = 0.; |
280 | }
|
281 | }
|
281 | 282 | ||
282 | // function is called after finishing the indivdual growth / mortality.
|
283 | // function is called after finishing the indivdual growth / mortality.
|
283 | void ResourceUnit::afterGrow() |
284 | void ResourceUnit::afterGrow() |
284 | {
|
285 | {
|
285 | mAverageAging = leafArea()>0.?mAverageAging/leafArea():0; // calculate aging value (calls to addAverageAging() by individual trees) |
286 | mAverageAging = leafArea()>0.?mAverageAging/leafArea():0; // calculate aging value (calls to addAverageAging() by individual trees) |
286 | if (mAverageAging>0. && mAverageAging<0.00001) |
287 | if (mAverageAging>0. && mAverageAging<0.00001) |
287 | qDebug() << "ru" << mIndex << "aging <0.00001"; |
288 | qDebug() << "ru" << mIndex << "aging <0.00001"; |
288 | if (mAverageAging<0. || mAverageAging>1.) |
289 | if (mAverageAging<0. || mAverageAging>1.) |
289 | qDebug() << "Average aging invalid: (RU, LAI):" << index() << mStatistics.leafAreaIndex(); |
290 | qDebug() << "Average aging invalid: (RU, LAI):" << index() << mStatistics.leafAreaIndex(); |
290 | }
|
291 | }
|
291 | 292 | ||
292 | void ResourceUnit::yearEnd() |
293 | void ResourceUnit::yearEnd() |
293 | {
|
294 | {
|
294 | // calculate statistics for all tree species of the ressource unit
|
295 | // calculate statistics for all tree species of the ressource unit
|
295 | int c = mRUSpecies.count(); |
296 | int c = mRUSpecies.count(); |
296 | for (int i=0;i<c; i++) { |
297 | for (int i=0;i<c; i++) { |
297 | mRUSpecies[i]->statisticsDead().calculate(); // calculate the dead trees |
298 | mRUSpecies[i]->statisticsDead().calculate(); // calculate the dead trees |
298 | mRUSpecies[i]->statisticsMgmt().calculate(); // stats of removed trees |
299 | mRUSpecies[i]->statisticsMgmt().calculate(); // stats of removed trees |
299 | mRUSpecies[i]->updateGWL(); // get sum of dead trees (died + removed) |
300 | mRUSpecies[i]->updateGWL(); // get sum of dead trees (died + removed) |
300 | mRUSpecies[i]->statistics().calculate(); // calculate the living (and add removed volume to gwl) |
301 | mRUSpecies[i]->statistics().calculate(); // calculate the living (and add removed volume to gwl) |
301 | mStatistics.add(mRUSpecies[i]->statistics()); |
302 | mStatistics.add(mRUSpecies[i]->statistics()); |
302 | }
|
303 | }
|
303 | mStatistics.calculate(); // aggreagte on stand level |
304 | mStatistics.calculate(); // aggreagte on stand level |
304 | 305 | ||
305 | }
|
306 | }
|
306 | 307 | ||
307 | void ResourceUnit::addTreeAgingForAllTrees() |
308 | void ResourceUnit::addTreeAgingForAllTrees() |
308 | {
|
309 | {
|
309 | mAverageAging = 0.; |
310 | mAverageAging = 0.; |
310 | foreach(const Tree &t, mTrees) { |
311 | foreach(const Tree &t, mTrees) { |
311 | addTreeAging(t.leafArea(), t.species()->aging(t.height(), t.age())); |
312 | addTreeAging(t.leafArea(), t.species()->aging(t.height(), t.age())); |
312 | }
|
313 | }
|
313 | 314 | ||
314 | }
|
315 | }
|
315 | 316 | ||
316 | /// refresh of tree based statistics.
|
317 | /// refresh of tree based statistics.
|
317 | /// WARNING: this function is only called once (during startup).
|
318 | /// WARNING: this function is only called once (during startup).
|
318 | /// see function "yearEnd()" above!!!
|
319 | /// see function "yearEnd()" above!!!
|
319 | void ResourceUnit::createStandStatistics() |
320 | void ResourceUnit::createStandStatistics() |
320 | {
|
321 | {
|
321 | // clear statistics (ru-level and ru-species level)
|
322 | // clear statistics (ru-level and ru-species level)
|
322 | mStatistics.clear(); |
323 | mStatistics.clear(); |
323 | for (int i=0;i<mRUSpecies.count();i++) { |
324 | for (int i=0;i<mRUSpecies.count();i++) { |
324 | mRUSpecies[i]->statistics().clear(); |
325 | mRUSpecies[i]->statistics().clear(); |
325 | mRUSpecies[i]->statisticsDead().clear(); |
326 | mRUSpecies[i]->statisticsDead().clear(); |
326 | mRUSpecies[i]->statisticsMgmt().clear(); |
327 | mRUSpecies[i]->statisticsMgmt().clear(); |
327 | }
|
328 | }
|
328 | 329 | ||
329 | // add all trees to the statistics objects of the species
|
330 | // add all trees to the statistics objects of the species
|
330 | foreach(const Tree &t, mTrees) { |
331 | foreach(const Tree &t, mTrees) { |
331 | if (!t.isDead()) |
332 | if (!t.isDead()) |
332 | resourceUnitSpecies(t.species()).statistics().add(&t, 0); |
333 | resourceUnitSpecies(t.species()).statistics().add(&t, 0); |
333 | }
|
334 | }
|
334 | // summarize statistics for the whole resource unit
|
335 | // summarize statistics for the whole resource unit
|
335 | for (int i=0;i<mRUSpecies.count();i++) { |
336 | for (int i=0;i<mRUSpecies.count();i++) { |
336 | mRUSpecies[i]->statistics().calculate(); |
337 | mRUSpecies[i]->statistics().calculate(); |
337 | mStatistics.add(mRUSpecies[i]->statistics()); |
338 | mStatistics.add(mRUSpecies[i]->statistics()); |
338 | }
|
339 | }
|
339 | mStatistics.calculate(); |
340 | mStatistics.calculate(); |
340 | mAverageAging = mStatistics.leafAreaIndex()>0.?mAverageAging / (mStatistics.leafAreaIndex()*stockableArea()):0.; |
341 | mAverageAging = mStatistics.leafAreaIndex()>0.?mAverageAging / (mStatistics.leafAreaIndex()*stockableArea()):0.; |
341 | if (mAverageAging<0. || mAverageAging>1.) |
342 | if (mAverageAging<0. || mAverageAging>1.) |
342 | qDebug() << "Average aging invalid: (RU, LAI):" << index() << mStatistics.leafAreaIndex(); |
343 | qDebug() << "Average aging invalid: (RU, LAI):" << index() << mStatistics.leafAreaIndex(); |
343 | }
|
344 | }
|
344 | 345 | ||
345 | void ResourceUnit::setMaxSaplingHeightAt(const QPoint &position, const float height) |
346 | void ResourceUnit::setMaxSaplingHeightAt(const QPoint &position, const float height) |
346 | {
|
347 | {
|
347 | Q_ASSERT(mSaplingHeightMap); |
348 | Q_ASSERT(mSaplingHeightMap); |
348 | int pixel_index = cPxPerRU*(position.x()-mCornerCoord.x())+(position.y()-mCornerCoord.y()); |
349 | int pixel_index = cPxPerRU*(position.x()-mCornerCoord.x())+(position.y()-mCornerCoord.y()); |
349 | if (pixel_index<0 || pixel_index>=cPxPerRU*cPxPerRU) { |
350 | if (pixel_index<0 || pixel_index>=cPxPerRU*cPxPerRU) { |
350 | qDebug() << "setSaplingHeightAt-Error for position" << position << "for RU at" << boundingBox() << "with corner" << mCornerCoord; |
351 | qDebug() << "setSaplingHeightAt-Error for position" << position << "for RU at" << boundingBox() << "with corner" << mCornerCoord; |
351 | } else { |
352 | } else { |
352 | if (mSaplingHeightMap[pixel_index]<height) |
353 | if (mSaplingHeightMap[pixel_index]<height) |
353 | mSaplingHeightMap[pixel_index]=height; |
354 | mSaplingHeightMap[pixel_index]=height; |
354 | }
|
355 | }
|
355 | }
|
356 | }
|
356 | 357 | ||
357 | /// clear all saplings of all species on a given position (after recruitment)
|
358 | /// clear all saplings of all species on a given position (after recruitment)
|
358 | void ResourceUnit::clearSaplings(const QPoint &position) |
359 | void ResourceUnit::clearSaplings(const QPoint &position) |
359 | {
|
360 | {
|
360 | foreach(ResourceUnitSpecies* rus, mRUSpecies) |
361 | foreach(ResourceUnitSpecies* rus, mRUSpecies) |
361 | rus->clearSaplings(position); |
362 | rus->clearSaplings(position); |
- | 363 | ||
- | 364 | }
|
|
- | 365 | ||
- | 366 | /// kill all saplings within a given rect
|
|
- | 367 | void ResourceUnit::clearSaplings(const QRectF pixel_rect, const bool remove_from_soil) |
|
- | 368 | {
|
|
- | 369 | foreach(ResourceUnitSpecies* rus, mRUSpecies) { |
|
- | 370 | rus->changeSapling().clearSaplings(pixel_rect, remove_from_soil); |
|
- | 371 | }
|
|
362 | 372 | ||
363 | }
|
373 | }
|
364 | 374 | ||
365 | float ResourceUnit::saplingHeightForInit(const QPoint &position) const |
375 | float ResourceUnit::saplingHeightForInit(const QPoint &position) const |
366 | {
|
376 | {
|
367 | double maxh = 0.; |
377 | double maxh = 0.; |
368 | foreach(ResourceUnitSpecies* rus, mRUSpecies) |
378 | foreach(ResourceUnitSpecies* rus, mRUSpecies) |
369 | maxh = qMax(maxh, rus->sapling().heightAt(position)); |
379 | maxh = qMax(maxh, rus->sapling().heightAt(position)); |
370 | return maxh; |
380 | return maxh; |
371 | }
|
381 | }
|
372 | 382 | ||
373 | void ResourceUnit::calculateCarbonCycle() |
383 | void ResourceUnit::calculateCarbonCycle() |
374 | {
|
384 | {
|
375 | if (!snag()) |
385 | if (!snag()) |
376 | return; |
386 | return; |
377 | 387 | ||
378 | // (1) calculate the snag dynamics
|
388 | // (1) calculate the snag dynamics
|
379 | // because all carbon/nitrogen-flows from trees to the soil are routed through the snag-layer,
|
389 | // because all carbon/nitrogen-flows from trees to the soil are routed through the snag-layer,
|
380 | // all soil inputs (litter + deadwood) are collected in the Snag-object.
|
390 | // all soil inputs (litter + deadwood) are collected in the Snag-object.
|
381 | snag()->calculateYear(); |
391 | snag()->calculateYear(); |
382 | soil()->setClimateFactor( snag()->climateFactor() ); // the climate factor is only calculated once |
392 | soil()->setClimateFactor( snag()->climateFactor() ); // the climate factor is only calculated once |
383 | soil()->setSoilInput( snag()->labileFlux(), snag()->refractoryFlux()); |
393 | soil()->setSoilInput( snag()->labileFlux(), snag()->refractoryFlux()); |
384 | soil()->calculateYear(); // update the ICBM/2N model |
394 | soil()->calculateYear(); // update the ICBM/2N model |
385 | // use available nitrogen?
|
395 | // use available nitrogen?
|
386 | if (Model::settings().useDynamicAvailableNitrogen) |
396 | if (Model::settings().useDynamicAvailableNitrogen) |
387 | mUnitVariables.nitrogenAvailable = soil()->availableNitrogen(); |
397 | mUnitVariables.nitrogenAvailable = soil()->availableNitrogen(); |
388 | 398 | ||
389 | // debug output
|
399 | // debug output
|
390 | if (GlobalSettings::instance()->isDebugEnabled(GlobalSettings::dCarbonCycle) && !snag()->isEmpty()) { |
400 | if (GlobalSettings::instance()->isDebugEnabled(GlobalSettings::dCarbonCycle) && !snag()->isEmpty()) { |
391 | DebugList &out = GlobalSettings::instance()->debugList(index(), GlobalSettings::dCarbonCycle); |
401 | DebugList &out = GlobalSettings::instance()->debugList(index(), GlobalSettings::dCarbonCycle); |
392 | out << index() << id(); // resource unit index and id |
402 | out << index() << id(); // resource unit index and id |
393 | out << snag()->debugList(); // snag debug outs |
403 | out << snag()->debugList(); // snag debug outs |
394 | out << soil()->debugList(); // ICBM/2N debug outs |
404 | out << soil()->debugList(); // ICBM/2N debug outs |
395 | }
|
405 | }
|
396 | 406 | ||
397 | }
|
407 | }
|
398 | 408 | ||
399 | 409 | ||
400 | 410 |