Subversion Repositories public iLand

Rev

Rev 1202 | Rev 1217 | 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
 
534 werner 21
/** @class ResourceUnit
22
  ResourceUnit is the spatial unit that encapsulates a forest stand and links to several environmental components
23
  (Climate, Soil, Water, ...).
697 werner 24
  @ingroup core
25
  A resource unit has a size of (currently) 100x100m. Many processes in iLand operate on the level of a ResourceUnit.
26
  Each resource unit has the same Climate and other properties (e.g. available nitrogen).
27
  Proceses on this level are, inter alia, NPP Production (see Production3PG), water calculations (WaterCycle), the modeling
28
  of dead trees (Snag) and soil processes (Soil).
534 werner 29
 
30
  */
31
#include <QtCore>
32
#include "global.h"
33
 
34
#include "resourceunit.h"
35
#include "resourceunitspecies.h"
36
#include "speciesset.h"
37
#include "species.h"
38
#include "production3pg.h"
39
#include "model.h"
40
#include "climate.h"
41
#include "watercycle.h"
42
#include "snag.h"
43
#include "soil.h"
44
#include "helper.h"
45
 
46
ResourceUnit::~ResourceUnit()
47
{
48
    if (mWater)
49
        delete mWater;
50
    mWater = 0;
51
    if (mSnag)
52
        delete mSnag;
53
    if (mSoil)
54
        delete mSoil;
55
 
738 werner 56
    qDeleteAll(mRUSpecies);
57
 
1159 werner 58
    if (mSaplings)
59
        delete[] mSaplings;
60
 
534 werner 61
    mSnag = 0;
62
    mSoil = 0;
1159 werner 63
    mSaplings = 0;
534 werner 64
}
65
 
66
ResourceUnit::ResourceUnit(const int index)
67
{
68
    qDeleteAll(mRUSpecies);
69
    mSpeciesSet = 0;
70
    mClimate = 0;
71
    mPixelCount=0;
72
    mStockedArea = 0;
73
    mStockedPixelCount = 0;
1157 werner 74
    mStockableArea = 0;
1024 werner 75
    mAggregatedWLA = 0.;
76
    mAggregatedLA = 0.;
77
    mAggregatedLR = 0.;
78
    mEffectiveArea = 0.;
79
    mLRI_modification = 0.;
534 werner 80
    mIndex = index;
81
    mSaplingHeightMap = 0;
82
    mEffectiveArea_perWLA = 0.;
83
    mWater = new WaterCycle();
84
    mSnag = 0;
85
    mSoil = 0;
1159 werner 86
    mSaplings = 0;
569 werner 87
    mID = 0;
534 werner 88
}
89
 
90
void ResourceUnit::setup()
91
{
92
    mWater->setup(this);
93
 
94
    if (mSnag)
95
        delete mSnag;
96
    mSnag=0;
97
    if (mSoil)
98
        delete mSoil;
99
    mSoil=0;
100
    if (Model::settings().carbonCycleEnabled) {
591 werner 101
        mSoil = new Soil(this);
534 werner 102
        mSnag = new Snag;
103
        mSnag->setup(this);
104
        const XmlHelper &xml=GlobalSettings::instance()->settings();
105
 
106
        // setup contents of the soil of the RU; use values for C and N (kg/ha)
107
        mSoil->setInitialState(CNPool(xml.valueDouble("model.site.youngLabileC", -1),
108
                                      xml.valueDouble("model.site.youngLabileN", -1),
109
                                      xml.valueDouble("model.site.youngLabileDecompRate", -1)),
110
                               CNPool(xml.valueDouble("model.site.youngRefractoryC", -1),
111
                                      xml.valueDouble("model.site.youngRefractoryN", -1),
112
                                      xml.valueDouble("model.site.youngRefractoryDecompRate", -1)),
113
                               CNPair(xml.valueDouble("model.site.somC", -1), xml.valueDouble("model.site.somN", -1)));
114
    }
115
 
1159 werner 116
    if (mSaplings)
117
        delete mSaplings;
118
    if (Model::settings().regenerationEnabled) {
119
        mSaplings = new SaplingCell[cPxPerHectare];
120
    }
121
 
534 werner 122
    // setup variables
123
    mUnitVariables.nitrogenAvailable = GlobalSettings::instance()->settings().valueDouble("model.site.availableNitrogen", 40);
124
 
895 werner 125
    // if dynamic coupling of soil nitrogen is enabled, a starting value for available N is calculated
534 werner 126
    if (mSoil && Model::settings().useDynamicAvailableNitrogen && Model::settings().carbonCycleEnabled) {
127
        mSoil->setClimateFactor(1.);
128
        mSoil->calculateYear();
895 werner 129
        mUnitVariables.nitrogenAvailable = soil()->availableNitrogen();
534 werner 130
    }
664 werner 131
    mHasDeadTrees = false;
534 werner 132
    mAverageAging = 0.;
133
 
134
}
135
void ResourceUnit::setBoundingBox(const QRectF &bb)
136
{
137
    mBoundingBox = bb;
1118 werner 138
    mCornerOffset = GlobalSettings::instance()->model()->grid()->indexAt(bb.topLeft());
534 werner 139
}
140
 
1203 werner 141
/// return the sapling cell at given LIF-coordinates
142
SaplingCell *ResourceUnit::saplingCell(const QPoint &lifCoords) const
143
{
144
    // LIF-Coordinates are global, we here need (RU-)local coordinates
145
    int ix = lifCoords.x() % cPxPerRU;
146
    int iy = lifCoords.y() % cPxPerRU;
147
    int i = iy*cPxPerRU+ix;
148
    Q_ASSERT(i>=0 && i<cPxPerHectare);
149
    return &mSaplings[i];
150
}
151
 
534 werner 152
/// set species and setup the species-per-RU-data
153
void ResourceUnit::setSpeciesSet(SpeciesSet *set)
154
{
155
    mSpeciesSet = set;
156
    qDeleteAll(mRUSpecies);
157
 
158
    //mRUSpecies.resize(set->count()); // ensure that the vector space is not relocated
159
    for (int i=0;i<set->count();i++) {
160
        Species *s = const_cast<Species*>(mSpeciesSet->species(i));
161
        if (!s)
162
            throw IException("ResourceUnit::setSpeciesSet: invalid index!");
163
 
164
        ResourceUnitSpecies *rus = new ResourceUnitSpecies();
165
        mRUSpecies.push_back(rus);
166
        rus->setup(s, this);
167
        /* be careful: setup() is called with a pointer somewhere to the content of the mRUSpecies container.
168
           If the container memory is relocated (QVector), the pointer gets invalid!!!
169
           Therefore, a resize() is called before the loop (no resize()-operations during the loop)! */
170
        //mRUSpecies[i].setup(s,this); // setup this element
171
 
172
    }
173
}
174
 
175
ResourceUnitSpecies &ResourceUnit::resourceUnitSpecies(const Species *species)
176
{
177
    return *mRUSpecies[species->index()];
178
}
179
 
1040 werner 180
const ResourceUnitSpecies *ResourceUnit::constResourceUnitSpecies(const Species *species) const
181
{
182
    return mRUSpecies[species->index()];
183
}
184
 
534 werner 185
Tree &ResourceUnit::newTree()
186
{
187
    // start simple: just append to the vector...
188
    if (mTrees.isEmpty())
189
        mTrees.reserve(100); // reserve a junk of memory for trees
190
 
191
    mTrees.append(Tree());
192
    return mTrees.back();
193
}
194
int ResourceUnit::newTreeIndex()
195
{
734 werner 196
    newTree();
197
    return mTrees.count()-1; // return index of the last tree
534 werner 198
}
199
 
200
/// remove dead trees from tree list
201
/// reduce size of vector if lots of space is free
202
/// tests showed that this way of cleanup is very fast,
203
/// because no memory allocations are performed (simple memmove())
204
/// when trees are moved.
205
void ResourceUnit::cleanTreeList()
206
{
664 werner 207
    if (!mHasDeadTrees)
208
        return;
209
 
534 werner 210
    QVector<Tree>::iterator last=mTrees.end()-1;
211
    QVector<Tree>::iterator current = mTrees.begin();
212
    while (last>=current && (*last).isDead())
213
        --last;
214
 
215
    while (current<last) {
216
        if ((*current).isDead()) {
217
            *current = *last; // copy data!
218
            --last; //
219
            while (last>=current && (*last).isDead())
220
                --last;
221
        }
222
        ++current;
223
    }
224
    ++last; // last points now to the first dead tree
225
 
226
    // free ressources
227
    if (last!=mTrees.end()) {
228
        mTrees.erase(last, mTrees.end());
229
        if (mTrees.capacity()>100) {
230
            if (mTrees.count() / double(mTrees.capacity()) < 0.2) {
231
                //int target_size = mTrees.count()*2;
232
                //qDebug() << "reduce size from "<<mTrees.capacity() << "to" << target_size;
233
                //mTrees.reserve(qMax(target_size, 100));
664 werner 234
                if (logLevelDebug())
235
                    qDebug() << "reduce tree storage of RU" << index() << " from " << mTrees.capacity() << "to" << mTrees.count();
534 werner 236
                mTrees.squeeze();
237
            }
238
        }
239
    }
664 werner 240
    mHasDeadTrees = false; // reset flag
534 werner 241
}
242
 
243
void ResourceUnit::newYear()
244
{
245
    mAggregatedWLA = 0.;
246
    mAggregatedLA = 0.;
247
    mAggregatedLR = 0.;
248
    mEffectiveArea = 0.;
249
    mPixelCount = mStockedPixelCount = 0;
250
    snagNewYear();
609 werner 251
    if (mSoil)
252
        mSoil->newYear();
534 werner 253
    // clear statistics global and per species...
254
    QList<ResourceUnitSpecies*>::const_iterator i;
255
    QList<ResourceUnitSpecies*>::const_iterator iend = mRUSpecies.constEnd();
256
    mStatistics.clear();
257
    for (i=mRUSpecies.constBegin(); i!=iend; ++i) {
258
        (*i)->statisticsDead().clear();
259
        (*i)->statisticsMgmt().clear();
260
    }
261
 
262
}
263
 
264
/** production() is the "stand-level" part of the biomass production (3PG).
265
    - The amount of radiation intercepted by the stand is calculated
266
    - the water cycle is calculated
267
    - statistics for each species are cleared
268
    - The 3PG production for each species and ressource unit is called (calculates species-responses and NPP production)
269
    see also: http://iland.boku.ac.at/individual+tree+light+availability */
270
void ResourceUnit::production()
271
{
272
 
1107 werner 273
    if (mAggregatedWLA==0. || mPixelCount==0) {
936 werner 274
        // clear statistics of resourceunitspecies
275
        for ( QList<ResourceUnitSpecies*>::const_iterator i=mRUSpecies.constBegin(); i!=mRUSpecies.constEnd(); ++i)
276
            (*i)->statistics().clear();
277
        mEffectiveArea = 0.;
278
        mStockedArea = 0.;
534 werner 279
        return;
280
    }
281
 
282
    // the pixel counters are filled during the height-grid-calculations
1184 werner 283
    mStockedArea = cHeightPerRU*cHeightPerRU * mStockedPixelCount; // m2 (1 height grid pixel = 10x10m)
1107 werner 284
    if (leafAreaIndex()<3.) {
285
        // estimate stocked area based on crown projections
286
        double crown_area = 0.;
287
        for (int i=0;i<mTrees.count();++i)
288
            crown_area += mTrees.at(i).isDead() ? 0. : mTrees.at(i).stamp()->reader()->crownArea();
534 werner 289
 
1157 werner 290
        if (logLevelDebug())
291
            qDebug() << "crown area: lai" << leafAreaIndex() << "stocked area (pixels)" << mStockedArea << " area (crown)" << crown_area;
292
        if (leafAreaIndex()<1.) {
293
            mStockedArea = std::min(crown_area, mStockedArea);
1107 werner 294
        } else {
1184 werner 295
            // for LAI between 1 and 3:
296
            // interpolate between sum of crown area of trees (at LAI=1) and the pixel-based value (at LAI=3 and above)
1157 werner 297
            double px_frac = (leafAreaIndex()-1.)/2.; // 0 at LAI=1, 1 at LAI=3
298
            mStockedArea = mStockedArea * px_frac + std::min(crown_area, mStockedArea) * (1. - px_frac);
1107 werner 299
        }
300
        if (mStockedArea==0.)
301
            return;
302
    }
303
 
534 werner 304
    // calculate the leaf area index (LAI)
305
    double LAI = mAggregatedLA / mStockedArea;
306
    // calculate the intercepted radiation fraction using the law of Beer Lambert
307
    const double k = Model::settings().lightExtinctionCoefficient;
308
    double interception_fraction = 1. - exp(-k * LAI);
309
    mEffectiveArea = mStockedArea * interception_fraction; // m2
310
 
311
    // calculate the total weighted leaf area on this RU:
312
    mLRI_modification = interception_fraction *  mStockedArea / mAggregatedWLA; // p_WLA
313
    if (mLRI_modification == 0.)
314
        qDebug() << "lri modifaction==0!";
315
 
611 werner 316
    if (logLevelDebug()) {
534 werner 317
    DBGMODE(qDebug() << QString("production: LAI: %1 (intercepted fraction: %2, stocked area: %4). LRI-Multiplier: %3")
318
            .arg(LAI)
319
            .arg(interception_fraction)
320
            .arg(mLRI_modification)
321
            .arg(mStockedArea);
322
    );
611 werner 323
    }
534 werner 324
 
325
    // calculate LAI fractions
326
    QList<ResourceUnitSpecies*>::const_iterator i;
327
    QList<ResourceUnitSpecies*>::const_iterator iend = mRUSpecies.constEnd();
328
    double ru_lai = leafAreaIndex();
329
    if (ru_lai < 1.)
330
        ru_lai = 1.;
331
    // note: LAIFactors are only 1 if sum of LAI is > 1. (see WaterCycle)
332
    for (i=mRUSpecies.constBegin(); i!=iend; ++i) {
720 werner 333
        double lai_factor = (*i)->statistics().leafAreaIndex() / ru_lai;
1157 werner 334
 
335
        //DBGMODE(
336
        if (lai_factor > 1.) {
337
                        const ResourceUnitSpecies* rus=*i;
338
                        qDebug() << "LAI factor > 1: species ru-index:" << rus->species()->name() << rus->ru()->index();
339
                    }
340
        //);
720 werner 341
        (*i)->setLAIfactor( lai_factor );
534 werner 342
    }
343
 
344
    // soil water model - this determines soil water contents needed for response calculations
345
    {
346
    mWater->run();
347
    }
348
 
349
    // invoke species specific calculation (3PG)
350
    for (i=mRUSpecies.constBegin(); i!=iend; ++i) {
1157 werner 351
        //DBGMODE(
352
        if ((*i)->LAIfactor() > 1.) {
353
                    const ResourceUnitSpecies* rus=*i;
354
                    qDebug() << "LAI factor > 1: species ru-index value:" << rus->species()->name() << rus->ru()->index() << rus->LAIfactor();
355
                    }
356
        //);
534 werner 357
        (*i)->calculate(); // CALCULATE 3PG
1196 werner 358
 
359
        // debug output related to production
360
        if (GlobalSettings::instance()->isDebugEnabled(GlobalSettings::dStandGPP) && (*i)->LAIfactor()>0.) {
361
            DebugList &out = GlobalSettings::instance()->debugList(index(), GlobalSettings::dStandGPP);
362
            out << (*i)->species()->id() << index() << id();
363
            out << (*i)->LAIfactor() << (*i)->prod3PG().GPPperArea() << productiveArea()*(*i)->LAIfactor()*(*i)->prod3PG().GPPperArea() << averageAging() << (*i)->prod3PG().fEnvYear() ;
364
 
365
        }
534 werner 366
    }
367
}
368
 
369
void ResourceUnit::calculateInterceptedArea()
370
{
371
    if (mAggregatedLR==0) {
372
        mEffectiveArea_perWLA = 0.;
373
        return;
374
    }
375
    Q_ASSERT(mAggregatedLR>0.);
376
    mEffectiveArea_perWLA = mEffectiveArea / mAggregatedLR;
377
    if (logLevelDebug()) qDebug() << "RU: aggregated lightresponse:" << mAggregatedLR  << "eff.area./wla:" << mEffectiveArea_perWLA;
378
}
379
 
380
// function is called immediately before the growth of individuals
381
void ResourceUnit::beforeGrow()
382
{
383
    mAverageAging = 0.;
384
}
385
 
386
// function is called after finishing the indivdual growth / mortality.
387
void ResourceUnit::afterGrow()
388
{
389
    mAverageAging = leafArea()>0.?mAverageAging/leafArea():0; // calculate aging value (calls to addAverageAging() by individual trees)
390
    if (mAverageAging>0. && mAverageAging<0.00001)
391
        qDebug() << "ru" << mIndex << "aging <0.00001";
392
    if (mAverageAging<0. || mAverageAging>1.)
393
        qDebug() << "Average aging invalid: (RU, LAI):" << index() << mStatistics.leafAreaIndex();
394
}
395
 
396
void ResourceUnit::yearEnd()
397
{
398
    // calculate statistics for all tree species of the ressource unit
399
    int c = mRUSpecies.count();
400
    for (int i=0;i<c; i++) {
401
        mRUSpecies[i]->statisticsDead().calculate(); // calculate the dead trees
402
        mRUSpecies[i]->statisticsMgmt().calculate(); // stats of removed trees
403
        mRUSpecies[i]->updateGWL(); // get sum of dead trees (died + removed)
404
        mRUSpecies[i]->statistics().calculate(); // calculate the living (and add removed volume to gwl)
405
        mStatistics.add(mRUSpecies[i]->statistics());
406
    }
407
    mStatistics.calculate(); // aggreagte on stand level
408
 
1157 werner 409
    // update carbon flows
410
    if (soil() && GlobalSettings::instance()->model()->settings().carbonCycleEnabled) {
411
        double area_factor = stockableArea() / cRUArea; //conversion factor
412
        mUnitVariables.carbonUptake = statistics().npp() * biomassCFraction;
413
        mUnitVariables.carbonUptake += statistics().nppSaplings() * biomassCFraction;
414
 
415
        double to_atm = snag()->fluxToAtmosphere().C / area_factor; // from snags, kgC/ha
416
        to_atm += soil()->fluxToAtmosphere().C *cRUArea/10.; // soil: t/ha -> t/m2 -> kg/ha
417
        mUnitVariables.carbonToAtm = to_atm;
418
 
419
        double to_dist = snag()->fluxToDisturbance().C / area_factor;
420
        to_dist += soil()->fluxToDisturbance().C * cRUArea/10.;
421
        double to_harvest = snag()->fluxToExtern().C / area_factor;
422
 
423
        mUnitVariables.NEP = mUnitVariables.carbonUptake - to_atm - to_dist - to_harvest; // kgC/ha
424
 
425
        // incremental values....
426
        mUnitVariables.cumCarbonUptake += mUnitVariables.carbonUptake;
427
        mUnitVariables.cumCarbonToAtm += mUnitVariables.carbonToAtm;
428
        mUnitVariables.cumNEP += mUnitVariables.NEP;
429
 
430
    }
431
 
534 werner 432
}
433
 
434
void ResourceUnit::addTreeAgingForAllTrees()
435
{
436
    mAverageAging = 0.;
437
    foreach(const Tree &t, mTrees) {
438
        addTreeAging(t.leafArea(), t.species()->aging(t.height(), t.age()));
439
    }
440
 
441
}
442
 
443
/// refresh of tree based statistics.
444
/// WARNING: this function is only called once (during startup).
445
/// see function "yearEnd()" above!!!
446
void ResourceUnit::createStandStatistics()
447
{
448
    // clear statistics (ru-level and ru-species level)
449
    mStatistics.clear();
450
    for (int i=0;i<mRUSpecies.count();i++) {
451
        mRUSpecies[i]->statistics().clear();
452
        mRUSpecies[i]->statisticsDead().clear();
453
        mRUSpecies[i]->statisticsMgmt().clear();
1178 werner 454
        mRUSpecies[i]->saplingStat().clearStatistics();
534 werner 455
    }
456
 
457
    // add all trees to the statistics objects of the species
458
    foreach(const Tree &t, mTrees) {
459
        if (!t.isDead())
460
            resourceUnitSpecies(t.species()).statistics().add(&t, 0);
461
    }
1178 werner 462
    // summarise sapling stats
463
    GlobalSettings::instance()->model()->saplings()->calculateInitialStatistics(this);
464
 
534 werner 465
    // summarize statistics for the whole resource unit
466
    for (int i=0;i<mRUSpecies.count();i++) {
1178 werner 467
        mRUSpecies[i]->saplingStat().calculate(mRUSpecies[i]->species(), this);
468
        mRUSpecies[i]->statistics().add(&mRUSpecies[i]->saplingStat());
534 werner 469
        mRUSpecies[i]->statistics().calculate();
470
        mStatistics.add(mRUSpecies[i]->statistics());
471
    }
472
    mStatistics.calculate();
575 werner 473
    mAverageAging = mStatistics.leafAreaIndex()>0.?mAverageAging / (mStatistics.leafAreaIndex()*stockableArea()):0.;
534 werner 474
    if (mAverageAging<0. || mAverageAging>1.)
475
        qDebug() << "Average aging invalid: (RU, LAI):" << index() << mStatistics.leafAreaIndex();
1178 werner 476
 
534 werner 477
}
478
 
720 werner 479
/** recreate statistics. This is necessary after events that changed the structure
480
    of the stand *after* the growth of trees (where stand statistics are updated).
481
    An example is after disturbances.  */
1157 werner 482
void ResourceUnit::recreateStandStatistics(bool recalculate_stats)
720 werner 483
{
1202 werner 484
    // when called after disturbances (recalculate_stats=false), we
485
    // clear only the tree-specific variables in the stats (i.e. we keep NPP, and regen carbon),
486
    // and then re-add all trees (since TreeGrowthData is NULL no NPP is available).
487
    // The statistics are not summarised here, because this happens for all resource units
488
    // in the yearEnd function of RU.
720 werner 489
    for (int i=0;i<mRUSpecies.count();i++) {
1202 werner 490
        if (recalculate_stats)
491
            mRUSpecies[i]->statistics().clear();
492
        else
493
            mRUSpecies[i]->statistics().clearOnlyTrees();
720 werner 494
    }
495
    foreach(const Tree &t, mTrees) {
496
        resourceUnitSpecies(t.species()).statistics().add(&t, 0);
497
    }
1157 werner 498
 
499
    if (recalculate_stats) {
500
        for (int i=0;i<mRUSpecies.count();i++) {
501
            mRUSpecies[i]->statistics().calculate();
502
        }
937 werner 503
    }
720 werner 504
}
505
 
824 werner 506
 
534 werner 507
 
508
 
509
void ResourceUnit::calculateCarbonCycle()
510
{
511
    if (!snag())
512
        return;
513
 
514
    // (1) calculate the snag dynamics
515
    // because all carbon/nitrogen-flows from trees to the soil are routed through the snag-layer,
516
    // all soil inputs (litter + deadwood) are collected in the Snag-object.
517
    snag()->calculateYear();
518
    soil()->setClimateFactor( snag()->climateFactor() ); // the climate factor is only calculated once
519
    soil()->setSoilInput( snag()->labileFlux(), snag()->refractoryFlux());
520
    soil()->calculateYear(); // update the ICBM/2N model
521
    // use available nitrogen?
522
    if (Model::settings().useDynamicAvailableNitrogen)
523
        mUnitVariables.nitrogenAvailable = soil()->availableNitrogen();
524
 
525
    // debug output
526
    if (GlobalSettings::instance()->isDebugEnabled(GlobalSettings::dCarbonCycle) && !snag()->isEmpty()) {
527
        DebugList &out = GlobalSettings::instance()->debugList(index(), GlobalSettings::dCarbonCycle);
605 werner 528
        out << index() << id(); // resource unit index and id
534 werner 529
        out << snag()->debugList(); // snag debug outs
530
        out << soil()->debugList(); // ICBM/2N debug outs
531
    }
532
 
533
}
600 werner 534
 
535