Rev 1221 | Details | Compare with Previous | Last modification | View Log | RSS feed
| Rev | Author | Line No. | Line |
|---|---|---|---|
| 671 | werner | 1 | /******************************************************************************************** |
| 2 | ** iLand - an individual based forest landscape and disturbance model |
||
| 3 | ** http://iland.boku.ac.at |
||
| 4 | ** Copyright (C) 2009- Werner Rammer, Rupert Seidl |
||
| 5 | ** |
||
| 6 | ** This program is free software: you can redistribute it and/or modify |
||
| 7 | ** it under the terms of the GNU General Public License as published by |
||
| 8 | ** the Free Software Foundation, either version 3 of the License, or |
||
| 9 | ** (at your option) any later version. |
||
| 10 | ** |
||
| 11 | ** This program is distributed in the hope that it will be useful, |
||
| 12 | ** but WITHOUT ANY WARRANTY; without even the implied warranty of |
||
| 13 | ** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the |
||
| 14 | ** GNU General Public License for more details. |
||
| 15 | ** |
||
| 16 | ** You should have received a copy of the GNU General Public License |
||
| 17 | ** along with this program. If not, see <http://www.gnu.org/licenses/>. |
||
| 18 | ********************************************************************************************/ |
||
| 19 | |||
| 534 | werner | 20 | #include "global.h" |
| 21 | #include "tests.h" |
||
| 22 | |||
| 23 | #include "helper.h" |
||
| 808 | werner | 24 | #include "debugtimer.h" |
| 534 | werner | 25 | #include "random.h" |
| 26 | #include "model.h" |
||
| 27 | #include "resourceunit.h" |
||
| 28 | #include "expressionwrapper.h" |
||
| 29 | #include "expression.h" |
||
| 30 | /// |
||
| 31 | #include "climate.h" |
||
| 32 | #include "species.h" |
||
| 33 | #include "speciesresponse.h" |
||
| 34 | #include "watercycle.h" |
||
| 35 | #include "csvfile.h" |
||
| 36 | #include "xmlhelper.h" |
||
| 37 | #include "environment.h" |
||
| 38 | #include "exception.h" |
||
| 39 | #include "seeddispersal.h" |
||
| 40 | #include "establishment.h" |
||
| 1111 | werner | 41 | #include "saplings.h" |
| 534 | werner | 42 | // |
| 43 | #include "standloader.h" |
||
| 44 | #include "soil.h" |
||
| 45 | |||
| 543 | werner | 46 | #include "mapgrid.h" |
| 544 | werner | 47 | #include "management.h" |
| 642 | werner | 48 | #include "dem.h" |
| 49 | #include "modelcontroller.h" |
||
| 543 | werner | 50 | |
| 652 | werner | 51 | #include "modules.h" |
| 52 | #include "../plugins/fire/fireplugin.h" |
||
| 53 | #include "../plugins/fire/firemodule.h" |
||
| 702 | werner | 54 | #include "../plugins/wind/windmodule.h" |
| 55 | #include "../plugins/wind/windplugin.h" |
||
| 642 | werner | 56 | #include <QInputDialog> |
| 780 | werner | 57 | #include <QtConcurrent/QtConcurrent> |
| 642 | werner | 58 | |
| 766 | werner | 59 | #include "spatialanalysis.h" |
| 60 | |||
| 808 | werner | 61 | #include "forestmanagementengine.h" |
| 876 | werner | 62 | #include "fmstp.h" |
| 808 | werner | 63 | |
| 534 | werner | 64 | Tests::Tests(QObject *wnd) |
| 65 | { |
||
| 66 | //QObject::connect(this, SIGNAL(needRepaint),wnd, SLOT(repaint)); |
||
| 67 | mParent=wnd; |
||
| 68 | |||
| 69 | } |
||
| 70 | |||
| 71 | |||
| 72 | void Tests::testXml() |
||
| 73 | { |
||
| 74 | XmlHelper xml( GlobalSettings::instance()->settings().node("model.species") ); |
||
| 75 | qDebug() << xml.value("source"); |
||
| 76 | QDomElement e = xml.node("source"); |
||
| 77 | xml.setNodeValue(e, "chief wiggum"); |
||
| 78 | qDebug() << xml.value("source"); // should print chief wiggum |
||
| 79 | xml.setNodeValue("model.species.source","species"); |
||
| 80 | qDebug() << "and back:" << xml.value("source"); |
||
| 81 | |||
| 82 | // Test of environment |
||
| 83 | Environment environment; |
||
| 84 | environment.loadFromFile(GlobalSettings::instance()->path("envtest.txt")); |
||
| 85 | try { |
||
| 86 | environment.setPosition(QPointF(34,10)); // 0/0 |
||
| 87 | environment.setPosition(QPointF(150,10)); // 1/0 |
||
| 88 | environment.setPosition(QPointF(250,10)); // 2/0 |
||
| 89 | Climate *c = environment.climate(); |
||
| 90 | SpeciesSet *s = environment.speciesSet(); |
||
| 91 | qDebug() << s->count(); |
||
| 92 | qDebug() << c->daylength_h(160); |
||
| 93 | environment.setPosition(QPointF(20,200)); // 0/2 |
||
| 94 | environment.setPosition(QPointF(-34,10)); // exception |
||
| 95 | } catch (const IException &ex) { |
||
| 575 | werner | 96 | qDebug() << ex.message(); |
| 534 | werner | 97 | } |
| 98 | } |
||
| 99 | |||
| 100 | void Tests::speedOfExpression() |
||
| 101 | { |
||
| 102 | Q_ASSERT(1==0); |
||
| 103 | // int *p; |
||
| 104 | // p = 0; |
||
| 105 | // *p = 1; --> signal handler does not really work |
||
| 106 | // (1) for each |
||
| 107 | double sum; |
||
| 108 | int count; |
||
| 109 | { |
||
| 110 | DebugTimer t("plain loop"); |
||
| 111 | for (int i=0;i<10;i++) { |
||
| 112 | sum=0; |
||
| 113 | count=0; |
||
| 114 | foreach(const ResourceUnit *ru,GlobalSettings::instance()->model()->ruList()) { |
||
| 115 | foreach(const Tree &tree, ru->constTrees()) { |
||
| 116 | sum+=tree.volume(); |
||
| 117 | count++; |
||
| 118 | } |
||
| 119 | } |
||
| 120 | |||
| 121 | } |
||
| 122 | qDebug() << "Sum of volume" << sum << "count" << count; |
||
| 123 | } |
||
| 124 | { |
||
| 125 | DebugTimer t("plain loop (iterator)"); |
||
| 126 | for (int i=0;i<10;i++) { |
||
| 127 | AllTreeIterator at(GlobalSettings::instance()->model()); |
||
| 128 | sum = 0.; |
||
| 129 | count = 0; |
||
| 130 | while (Tree *tree=at.next()) { |
||
| 781 | werner | 131 | sum += pow(tree->dbh(),2.1f); count++; |
| 534 | werner | 132 | } |
| 133 | } |
||
| 134 | qDebug() << "Sum of volume" << sum << "count" << count; |
||
| 135 | } |
||
| 136 | { |
||
| 137 | TreeWrapper tw; |
||
| 138 | Expression expr("dbh^2.1", &tw); |
||
| 139 | DebugTimer t("Expression loop"); |
||
| 140 | for (int i=0;i<10;i++) { |
||
| 141 | |||
| 142 | AllTreeIterator at(GlobalSettings::instance()->model()); |
||
| 143 | sum = 0.; |
||
| 144 | while (Tree *tree=at.next()) { |
||
| 145 | tw.setTree(tree); |
||
| 146 | sum += expr.execute(); |
||
| 147 | } |
||
| 148 | } |
||
| 149 | qDebug() << "Sum of volume" << sum; |
||
| 150 | } |
||
| 151 | } |
||
| 152 | |||
| 153 | |||
| 154 | void Tests::clearTrees() |
||
| 155 | { |
||
| 156 | |||
| 157 | ResourceUnit *ru = GlobalSettings::instance()->model()->ru(); |
||
| 158 | int tc = ru->trees().count(); |
||
| 159 | // kill n percent... |
||
| 160 | ru->trees().last().die(); |
||
| 161 | ru->trees().first().die(); |
||
| 162 | if (tc>20) { |
||
| 163 | ru->trees()[15].die(); |
||
| 164 | ru->trees()[16].die(); |
||
| 165 | ru->trees()[17].die(); |
||
| 166 | ru->trees()[18].die(); |
||
| 167 | ru->trees()[19].die(); |
||
| 168 | ru->trees()[20].die(); |
||
| 169 | } |
||
| 170 | qDebug() << "killed 8 trees"; |
||
| 171 | ru->cleanTreeList(); |
||
| 172 | ru->cleanTreeList(); |
||
| 173 | } |
||
| 174 | |||
| 175 | void Tests::killTrees() |
||
| 176 | { |
||
| 177 | AllTreeIterator at(GlobalSettings::instance()->model()); |
||
| 178 | int count=0, totalcount=0, idc=0; |
||
| 179 | while (Tree *t = at.next()) { |
||
| 180 | totalcount++; |
||
| 181 | if (drandom() < 0.20) { |
||
| 182 | t->die(); |
||
| 183 | count++; |
||
| 184 | } |
||
| 185 | } |
||
| 186 | |||
| 187 | qDebug() << "killed" << count << "of" << totalcount << "left:" << totalcount-count; |
||
| 188 | { DebugTimer t("count living"); |
||
| 189 | at.reset(); |
||
| 190 | count=0; |
||
| 191 | idc=0; |
||
| 192 | while (at.nextLiving()) { |
||
| 193 | count++; idc+=at.current()->id(); |
||
| 194 | //qDebug() << at.current()->id(); |
||
| 195 | } |
||
| 196 | } |
||
| 197 | qDebug() << count << "living trees idsum" << idc; |
||
| 198 | |||
| 199 | { |
||
| 200 | DebugTimer t("clear trees"); |
||
| 201 | foreach(ResourceUnit *ru,GlobalSettings::instance()->model()->ruList()) |
||
| 202 | ru->cleanTreeList(); |
||
| 203 | } |
||
| 204 | // count all trees |
||
| 205 | at.reset(); |
||
| 206 | count=0; |
||
| 207 | idc=0; |
||
| 208 | while (at.next()) |
||
| 209 | { |
||
| 210 | count++; idc+=at.current()->id(); |
||
| 211 | //qDebug() << (*at)->id(); |
||
| 212 | } |
||
| 213 | qDebug() << count << "trees left idsum" << idc; |
||
| 214 | } |
||
| 215 | |||
| 216 | void Tests::climate() |
||
| 217 | { |
||
| 218 | Climate clim; |
||
| 219 | try { |
||
| 220 | clim.setup(); |
||
| 221 | DebugTimer t("climate 100yrs"); |
||
| 222 | const ClimateDay *begin, *end; |
||
| 223 | int mon=0; |
||
| 224 | for (int i=0;i<100;i++) { |
||
| 225 | clim.monthRange(mon, &begin, &end); |
||
| 226 | clim.nextYear(); |
||
| 227 | mon=(mon+1)%12; |
||
| 228 | } |
||
| 229 | } catch (IException &e) { |
||
| 575 | werner | 230 | Helper::msg(e.message()); |
| 534 | werner | 231 | } |
| 232 | } |
||
| 233 | |||
| 234 | |||
| 235 | void Tests::testSun() |
||
| 236 | { |
||
| 237 | // solar radiation |
||
| 238 | Sun sun; |
||
| 239 | sun.setup(RAD(47)); |
||
| 240 | qDebug()<<sun.dump(); |
||
| 241 | sun.setup(RAD(70)); |
||
| 242 | qDebug()<<sun.dump(); |
||
| 243 | sun.setup(RAD(-20)); |
||
| 244 | qDebug()<<sun.dump(); |
||
| 245 | } |
||
| 246 | |||
| 247 | void Tests::testWater() |
||
| 248 | { |
||
| 249 | Model *model = GlobalSettings::instance()->model(); |
||
| 250 | model->createStandStatistics(); // force this! |
||
| 251 | WaterCycle wc; |
||
| 252 | wc.setup(model->ru()); |
||
| 253 | wc.run(); |
||
| 562 | werner | 254 | for (int i=0;i<12;i++) |
| 255 | qDebug() << wc.referenceEvapotranspiration()[i]; |
||
| 534 | werner | 256 | } |
| 257 | |||
| 258 | void Tests::testPheno(const Climate *clim) |
||
| 259 | { |
||
| 260 | Phenology pheno(1,clim,0.9, 4.1, 10., 11., 2. ,9.); |
||
| 261 | pheno.calculate(); |
||
| 262 | qDebug() << "Phenology is the key:"; |
||
| 263 | for (int i=0;i<12;i++) |
||
| 264 | qDebug() << i << pheno.month()[i]; |
||
| 265 | } |
||
| 266 | |||
| 267 | void Tests::climateResponse() |
||
| 268 | { |
||
| 269 | |||
| 270 | try { |
||
| 271 | |||
| 272 | DebugTimer t("climate Responses"); |
||
| 273 | |||
| 274 | // get a climate response object.... |
||
| 275 | Model *model = GlobalSettings::instance()->model(); |
||
| 276 | //model->ru()->climate()->setup(); // force setup |
||
| 277 | |||
| 278 | ResourceUnitSpecies *rus = model->ru()->ruSpecies().first(); |
||
| 279 | |||
| 280 | rus->calculate(); |
||
| 281 | const SpeciesResponse *sr = rus->speciesResponse(); |
||
| 282 | QString line; |
||
| 283 | for (int mon=0;mon<12;mon++) { |
||
| 284 | line = QString("%1;%2").arg(mon) |
||
| 285 | .arg(sr->tempResponse()[mon]); |
||
| 286 | qDebug() << line; |
||
| 287 | } |
||
| 288 | // test nitrogen response |
||
| 289 | line="Nitrogen response\n"; |
||
| 290 | for (double navailable=0.; navailable<100;navailable+=10) { |
||
| 291 | for (double cls=1.; cls<=3.; cls+=0.2) { |
||
| 292 | double res = model->ru()->speciesSet()->nitrogenResponse(navailable,cls); |
||
| 293 | line+=QString("%1;%2;%3\n").arg(navailable).arg(cls).arg(res); |
||
| 294 | } |
||
| 295 | } |
||
| 296 | qDebug() << line; |
||
| 297 | |||
| 298 | line="CO2 response\n"; |
||
| 299 | for (double nresponse=0.; nresponse<=1;nresponse+=0.1) { |
||
| 300 | for (double h2o=0.; h2o<=1; h2o+=0.1) { |
||
| 301 | for (double co2=250; co2<600; co2+=50) { |
||
| 302 | double res = model->ru()->speciesSet()->co2Response(co2, nresponse, h2o); |
||
| 303 | line+=QString("%1;%2;%3;%4\n").arg(nresponse).arg(h2o).arg(co2).arg(res); |
||
| 304 | } |
||
| 305 | } |
||
| 306 | } |
||
| 307 | qDebug() << line; |
||
| 308 | |||
| 309 | // sun |
||
| 310 | //testSun(); |
||
| 311 | testPheno(model->ru()->climate()); |
||
| 312 | |||
| 313 | // test phenology |
||
| 314 | } catch (IException &e) { |
||
| 575 | werner | 315 | Helper::msg(e.message()); |
| 534 | werner | 316 | } |
| 317 | } |
||
| 318 | |||
| 319 | void Tests::multipleLightRuns(const QString &fileName) |
||
| 320 | { |
||
| 321 | XmlHelper xml(fileName); |
||
| 322 | Model *model = GlobalSettings::instance()->model(); |
||
| 323 | if (!model || !model->ru()) |
||
| 324 | return; |
||
| 325 | QVector<Tree> &mTrees = model->ru()->trees(); |
||
| 326 | |||
| 327 | QString outPath = xml.value("outputpath"); |
||
| 328 | QString inPath = xml.value("inputpath"); |
||
| 329 | QString inFile = xml.value("stands"); |
||
| 330 | qDebug() << "standlist:" << inFile << "inpath:"<<inPath << "save to:"<<outPath; |
||
| 331 | QStringList fileList = Helper::loadTextFile(inFile).remove('\r').split('\n', QString::SkipEmptyParts); |
||
| 332 | |||
| 333 | StandLoader loader(model); |
||
| 334 | try { |
||
| 335 | foreach (QString file, fileList) { |
||
| 336 | file = inPath + "\\" + file; |
||
| 337 | qDebug() << "processing" << file; |
||
| 338 | Tree::resetStatistics(); |
||
| 339 | mTrees.clear(); |
||
| 340 | loader.loadPicusFile(file); |
||
| 341 | // do a cycle... |
||
| 342 | model->runYear(); |
||
| 343 | |||
| 344 | // create output file |
||
| 345 | QFileInfo fi(file); |
||
| 346 | QString outFileName = QString("%1\\out_%2.csv").arg(outPath, fi.baseName()); |
||
| 347 | Helper::saveToTextFile(outFileName, dumpTreeList() ); |
||
| 348 | qDebug() << mTrees.size() << "trees loaded, saved to" << outFileName; |
||
| 349 | mParent->metaObject()->invokeMethod(mParent,"repaint"); |
||
| 350 | QApplication::processEvents(); |
||
| 351 | |||
| 352 | } |
||
| 353 | } catch (IException &e) { |
||
| 575 | werner | 354 | Helper::msg(e.message()); |
| 534 | werner | 355 | } |
| 356 | } |
||
| 357 | |||
| 358 | QString Tests::dumpTreeList() |
||
| 359 | { |
||
| 360 | |||
| 361 | Model *model = GlobalSettings::instance()->model(); |
||
| 362 | |||
| 363 | AllTreeIterator at(model); |
||
| 364 | DebugList treelist; |
||
| 365 | QString line; |
||
| 366 | QStringList result; |
||
| 367 | result << "id;species;dbh;height;x;y;RU#;LRI;mWoody;mRoot;mFoliage;LA"; |
||
| 368 | |||
| 369 | while (Tree *tree = at.next()) { |
||
| 370 | treelist.clear(); |
||
| 371 | tree->dumpList(treelist); |
||
| 372 | line = ""; |
||
| 373 | foreach(QVariant value, treelist) |
||
| 374 | line+=value.toString() + ";"; |
||
| 375 | *(line.end()-1)=' '; |
||
| 376 | result << line; |
||
| 377 | } |
||
| 378 | |||
| 379 | QString resStr = result.join("\n"); |
||
| 380 | return resStr; |
||
| 381 | } |
||
| 382 | |||
| 383 | void Tests::testCSVFile() |
||
| 384 | { |
||
| 385 | CSVFile file; |
||
| 386 | file.loadFile("e:\\csvtest.txt"); |
||
| 387 | for (int row=0;row<file.rowCount(); row++) |
||
| 388 | for (int col=0;col<file.colCount(); col++) |
||
| 389 | qDebug() << "row col" << row << file.columnName(col) << "value" << file.value(row, col); |
||
| 390 | } |
||
| 391 | |||
| 707 | werner | 392 | #include "../3rdparty/MersenneTwister.h" |
| 534 | werner | 393 | void Tests::testRandom() |
| 394 | { |
||
| 1164 | werner | 395 | |
| 396 | QStringList list; |
||
| 397 | for (int i=0;i<1000;i++) |
||
| 398 | list << QString::number(irandom(0,5)); |
||
| 399 | qDebug() << "irandom test (0,5): " << list; |
||
| 400 | return; |
||
| 401 | |||
| 402 | |||
| 757 | werner | 403 | RandomGenerator::setup(RandomGenerator::ergMersenneTwister, 1); |
| 404 | // RandomCustomPDF pdf("x^2"); |
||
| 405 | // RandomCustomPDF *pdf2 = new RandomCustomPDF("x^3"); |
||
| 406 | // QStringList list; |
||
| 407 | // for (int i=0;i<1000;i++) |
||
| 408 | // list << QString::number(pdf.get()); |
||
| 409 | // qDebug() << list.join("\n"); |
||
| 410 | // delete pdf2; |
||
| 411 | // // simple random test |
||
| 412 | // list.clear(); |
||
| 413 | // for (int i=0;i<1000;i++) |
||
| 414 | // list << QString::number(irandom(0,5)); |
||
| 415 | // qDebug() << "irandom test (0,5): " << list; |
||
| 534 | werner | 416 | |
| 707 | werner | 417 | // random generator timings |
| 418 | { DebugTimer t("mersenne"); |
||
| 757 | werner | 419 | RandomGenerator::setup(RandomGenerator::ergMersenneTwister, 0); |
| 420 | int n8=0, n7=0, n6=0, n5=0, n4=0, n3=0; |
||
| 421 | for (unsigned int i=0;i<4000000000U;i++) { |
||
| 422 | double r = drandom(); |
||
| 423 | if (r< 0.00000001) ++n8; |
||
| 424 | if (r< 0.0000001) ++n7; |
||
| 425 | if (r< 0.000001) ++n6; |
||
| 426 | if (r< 0.00001) ++n5; |
||
| 427 | if (r< 0.0001) ++n4; |
||
| 428 | if (r< 0.001) ++n3; |
||
| 707 | werner | 429 | } |
| 757 | werner | 430 | qDebug() << "Mersenne" << n3 << n4 << n5 << n6 << n7 << n8; |
| 707 | werner | 431 | } |
| 757 | werner | 432 | |
| 433 | { DebugTimer t("mersenne"); |
||
| 434 | RandomGenerator::setup(RandomGenerator::ergMersenneTwister, 0); |
||
| 435 | int n8=0, n7=0, n6=0, n5=0, n4=0, n3=0; |
||
| 436 | for (int j=0;j<1000;j++) { |
||
| 437 | RandomGenerator::seed(0); // new seed |
||
| 438 | for (int i=0;i<4000000;i++) { |
||
| 439 | double r = drandom(); |
||
| 440 | if (r< 0.00000001) ++n8; |
||
| 441 | if (r< 0.0000001) ++n7; |
||
| 442 | if (r< 0.000001) ++n6; |
||
| 443 | if (r< 0.00001) ++n5; |
||
| 444 | if (r< 0.0001) ++n4; |
||
| 445 | if (r< 0.001) ++n3; |
||
| 446 | } |
||
| 447 | } |
||
| 448 | qDebug() << "Mersenne2" << n3 << n4 << n5 << n6 << n7 << n8; |
||
| 449 | } |
||
| 707 | werner | 450 | { DebugTimer t("XORShift"); |
| 757 | werner | 451 | RandomGenerator::setup(RandomGenerator::ergXORShift96, 0); |
| 452 | int n8=0, n7=0, n6=0, n5=0, n4=0, n3=0; |
||
| 453 | for (unsigned int i=0;i<4000000000U;i++) { |
||
| 454 | double r = drandom(); |
||
| 455 | if (r< 0.00000001) ++n8; |
||
| 456 | if (r< 0.0000001) ++n7; |
||
| 457 | if (r< 0.000001) ++n6; |
||
| 458 | if (r< 0.00001) ++n5; |
||
| 459 | if (r< 0.0001) ++n4; |
||
| 460 | if (r< 0.001) ++n3; |
||
| 707 | werner | 461 | } |
| 757 | werner | 462 | qDebug() << "XORShift" << n3 << n4 << n5 << n6 << n7 << n8; |
| 707 | werner | 463 | } |
| 464 | { DebugTimer t("WellRNG"); |
||
| 757 | werner | 465 | RandomGenerator::setup(RandomGenerator::ergWellRNG512, 0); |
| 466 | int n8=0, n7=0, n6=0, n5=0, n4=0, n3=0; |
||
| 467 | for (unsigned int i=0;i<4000000000U;i++) { |
||
| 468 | double r = drandom(); |
||
| 469 | if (r< 0.00000001) ++n8; |
||
| 470 | if (r< 0.0000001) ++n7; |
||
| 471 | if (r< 0.000001) ++n6; |
||
| 472 | if (r< 0.00001) ++n5; |
||
| 473 | if (r< 0.0001) ++n4; |
||
| 474 | if (r< 0.001) ++n3; |
||
| 707 | werner | 475 | } |
| 757 | werner | 476 | qDebug() << "WellRNG" << n3 << n4 << n5 << n6 << n7 << n8; |
| 707 | werner | 477 | } |
| 478 | { DebugTimer t("FastRandom"); |
||
| 757 | werner | 479 | RandomGenerator::setup(RandomGenerator::ergFastRandom, 0); |
| 480 | int n8=0, n7=0, n6=0, n5=0, n4=0, n3=0; |
||
| 481 | for (unsigned int i=0;i<4000000000U;i++) { |
||
| 482 | double r = drandom(); |
||
| 483 | if (r< 0.00000001) ++n8; |
||
| 484 | if (r< 0.0000001) ++n7; |
||
| 485 | if (r< 0.000001) ++n6; |
||
| 486 | if (r< 0.00001) ++n5; |
||
| 487 | if (r< 0.0001) ++n4; |
||
| 488 | if (r< 0.001) ++n3; |
||
| 707 | werner | 489 | } |
| 757 | werner | 490 | qDebug() << "FastRandom" << n3 << n4 << n5 << n6 << n7 << n8; |
| 707 | werner | 491 | } |
| 492 | { DebugTimer t("MersenneTwister - Reference"); |
||
| 493 | MTRand rand; |
||
| 494 | rand.seed(1); |
||
| 495 | |||
| 496 | double sum = 0; |
||
| 497 | for (int i=0;i<100000000;i++) { |
||
| 757 | werner | 498 | if (drandom()<0.0000001) ++sum; |
| 707 | werner | 499 | } |
| 500 | qDebug() << "Mersenne Reference" << sum; |
||
| 501 | } |
||
| 534 | werner | 502 | } |
| 503 | |||
| 504 | void Tests::testGridRunner() |
||
| 505 | { |
||
| 1170 | werner | 506 | Grid<int> test_grid(10,5,5); // 5x5 grid |
| 507 | test_grid.initialize(1); |
||
| 508 | GridRunner<int> r(&test_grid); |
||
| 509 | int sum = 0; |
||
| 510 | while (int *p=r.next()) |
||
| 511 | ++sum; |
||
| 512 | qDebug() << "runner all, expected: 25:" << sum; |
||
| 513 | sum=0; |
||
| 514 | for (int *p=test_grid.begin();p!=test_grid.end();++p) |
||
| 515 | ++sum; |
||
| 516 | qDebug() << "loop expected: 25:" << sum; |
||
| 517 | sum=0; |
||
| 518 | r=GridRunner<int>(test_grid,test_grid.rectangle()); |
||
| 519 | sum = 0; |
||
| 520 | while (int *p=r.next()) |
||
| 521 | ++sum; |
||
| 522 | qDebug() << "runner rectangle, expected: 25:" << sum; |
||
| 523 | |||
| 524 | r=GridRunner<int>(test_grid,test_grid.metricRect()); |
||
| 525 | sum = 0; |
||
| 526 | while (int *p=r.next()) |
||
| 527 | ++sum; |
||
| 528 | qDebug() << "runner rectangle metric, expected: 25:" << sum; |
||
| 529 | |||
| 530 | // smaller parts |
||
| 531 | QRect r1(0,0,3,3); |
||
| 532 | r=GridRunner<int>(test_grid,r1); |
||
| 533 | sum = 0; |
||
| 534 | while (int *p=r.next()) { |
||
| 535 | ++sum; |
||
| 536 | qDebug() << r.currentIndex(); |
||
| 537 | } |
||
| 538 | qDebug() << "runner QRect(0,0,3,3), expected: 4:" << sum; |
||
| 539 | |||
| 540 | |||
| 541 | r1=QRect(0,0,1,1); |
||
| 542 | r=GridRunner<int>(test_grid,r1); |
||
| 543 | sum = 0; |
||
| 544 | while (int *p=r.next()) { |
||
| 545 | ++sum; |
||
| 546 | qDebug() << r.currentIndex(); |
||
| 547 | } |
||
| 548 | qDebug() << "runner QRect(0,0,1,1), expected: 1:" << sum; |
||
| 549 | |||
| 550 | r1=QRect(0,0,2,2); |
||
| 551 | r=GridRunner<int>(test_grid,r1); |
||
| 552 | sum = 0; |
||
| 553 | while (int *p=r.next()) { |
||
| 554 | ++sum; |
||
| 555 | qDebug() << r.currentIndex(); |
||
| 556 | } |
||
| 557 | qDebug() << "runner QRect(0,0,2,2), expected: 1:" << sum; |
||
| 558 | |||
| 559 | r=GridRunner<int>(test_grid,QRect(QPoint(0,0), QPoint(3,3))); |
||
| 560 | sum = 0; |
||
| 561 | while (int *p=r.next()) { |
||
| 562 | ++sum; |
||
| 563 | qDebug() << r.currentIndex(); |
||
| 564 | } |
||
| 565 | qDebug() << "runner QRect(QPoint(0,0), QPoint(3,3)), expected: 4:" << sum; |
||
| 566 | |||
| 567 | r=GridRunner<int>(test_grid,QRect(QPoint(0,0), QPoint(0,0))); |
||
| 568 | sum = 0; |
||
| 569 | while (int *p=r.next()) { |
||
| 570 | ++sum; |
||
| 571 | qDebug() << r.currentIndex(); |
||
| 572 | } |
||
| 573 | qDebug() << "runner QRect(QPoint(0,0), QPoint(0,0)), expected: 1:" << sum; |
||
| 574 | |||
| 575 | r=GridRunner<int>(test_grid,QRect(QPoint(0,0), QPoint(1,1))); |
||
| 576 | sum = 0; |
||
| 577 | while (int *p=r.next()) { |
||
| 578 | ++sum; |
||
| 579 | qDebug() << r.currentIndex(); |
||
| 580 | } |
||
| 581 | qDebug() << "runner QRect(QPoint(0,0), QPoint(1,1)), expected: 1:" << sum; |
||
| 582 | |||
| 583 | |||
| 584 | if (!GlobalSettings::instance()->model()) |
||
| 585 | return; |
||
| 586 | |||
| 534 | werner | 587 | Grid<float> &lif = *GlobalSettings::instance()->model()->grid(); |
| 588 | //QRectF box = GlobalSettings::instance()->model()->ru(0)->boundingBox(); |
||
| 589 | QRectF box2 = QRectF(10,10,10,10); |
||
| 590 | GridRunner<float> runner(lif, box2); |
||
| 591 | int i=0; |
||
| 592 | while (float *p=runner.next()) { |
||
| 593 | QPoint point = lif.indexOf(p); |
||
| 594 | qDebug() << i++ << point.x() << point.y() << *p << p; |
||
| 595 | } |
||
| 651 | werner | 596 | QRect index_rect = QRect(lif.indexAt(QPointF(10., 10.)), lif.indexAt(QPointF(20., 20))); |
| 597 | GridRunner<float> runner2(lif, index_rect); |
||
| 598 | qDebug() << "index variables for rect" << index_rect; |
||
| 599 | i = 0; |
||
| 600 | while (float *p=runner2.next()) { |
||
| 601 | QPoint point = lif.indexOf(p); |
||
| 602 | qDebug() << i++ << point.x() << point.y() << *p << p; |
||
| 603 | } |
||
| 1170 | werner | 604 | index_rect = QRect(0,0,1,1); |
| 605 | runner2 = GridRunner<float>(lif, index_rect); |
||
| 606 | qDebug() << "test index: rect" << index_rect; |
||
| 607 | i = 0; |
||
| 608 | while (float *p=runner2.next()) { |
||
| 609 | QPoint point = lif.indexOf(p); |
||
| 610 | qDebug() << i++ << point.x() << point.y() << *p << p; |
||
| 611 | } |
||
| 612 | |||
| 613 | index_rect = QRect(0,0,3,3); |
||
| 614 | runner2 = GridRunner<float>(lif, index_rect); |
||
| 615 | qDebug() << "test index: rect" << index_rect; |
||
| 616 | i = 0; |
||
| 617 | while (float *p=runner2.next()) { |
||
| 618 | QPoint point = lif.indexOf(p); |
||
| 619 | qDebug() << i++ << point.x() << point.y() << *p << p; |
||
| 620 | } |
||
| 621 | |||
| 1157 | werner | 622 | for (int i=0;i<GlobalSettings::instance()->model()->ruList().size();++i) { |
| 623 | if (i==10) break; |
||
| 624 | ResourceUnit *ru = GlobalSettings::instance()->model()->ruList()[i]; |
||
| 625 | GridRunner<float> runner(lif, ru->boundingBox()); |
||
| 626 | int n=0; |
||
| 627 | while (runner.next()) |
||
| 628 | ++n; |
||
| 629 | qDebug() << "RU" << ru->index() << ru->boundingBox() << lif.indexOf(runner.first()) << lif.indexOf(runner.last()) << "n:" << n; |
||
| 630 | } |
||
| 534 | werner | 631 | } |
| 632 | |||
| 633 | void Tests::testSeedDispersal() |
||
| 634 | { |
||
| 635 | SeedDispersal sd; |
||
| 636 | sd.setup(); |
||
| 637 | sd.loadFromImage("e:\\temp\\test.jpg"); |
||
| 638 | QImage img = gridToImage(sd.seedMap(), true, -1., 1.); |
||
| 639 | img.save("e:\\temp\\seedmap.png"); |
||
| 640 | sd.execute(); |
||
| 641 | // sd.edgeDetection(); |
||
| 642 | // gridToImage(sd.seedMap(), true, -1., 1.).save("seedmap_edge.png"); |
||
| 643 | // sd.distribute(); |
||
| 644 | img = gridToImage(sd.seedMap(), true, -1., 1.); |
||
| 645 | img.save("e:\\temp\\seedmap_e.png"); |
||
| 646 | |||
| 647 | } |
||
| 648 | |||
| 649 | Expression tme_exp; |
||
| 650 | int tme_count; |
||
| 651 | double tme_test1(const double &x) { |
||
| 652 | double result = tme_exp.calculate(x); |
||
| 653 | if (result!=1) |
||
| 654 | tme_count++; |
||
| 655 | return result; |
||
| 656 | } |
||
| 657 | double tme_test2(const double &x) { |
||
| 658 | tme_exp.setVar("x",x); |
||
| 659 | double result = tme_exp.execute(); |
||
| 660 | if (result!=1) |
||
| 661 | tme_count++; |
||
| 662 | return result; |
||
| 663 | } |
||
| 664 | QMutex tme_mutex; |
||
| 665 | double tme_test3(const double &x) { |
||
| 666 | QMutexLocker l(&tme_mutex); |
||
| 667 | tme_exp.setVar("x",x); |
||
| 668 | double result = tme_exp.execute(); |
||
| 669 | if (result!=1) |
||
| 670 | tme_count++; |
||
| 671 | return result; |
||
| 672 | } |
||
| 1157 | werner | 673 | |
| 674 | static double testf_sum = 0.; |
||
| 675 | void testF(double *begin, double *end) { |
||
| 676 | double sum = 0.; |
||
| 677 | for (double *p=begin;p!=end;++p) { |
||
| 678 | sum += *p; |
||
| 679 | } |
||
| 680 | QMutexLocker l(&tme_mutex); |
||
| 681 | testf_sum += sum; |
||
| 682 | //qDebug() << "testf min:" << *begin << "end:" << *(end-1); |
||
| 683 | } |
||
| 684 | |||
| 534 | werner | 685 | void Tests::testMultithreadExecute() |
| 686 | { |
||
| 1157 | werner | 687 | // test of threadrunner |
| 688 | int l=12345; |
||
| 689 | double *dat=new double[l]; |
||
| 690 | for (int i=0;i<l;++i) |
||
| 691 | dat[i] = i+1; |
||
| 692 | ThreadRunner tr; |
||
| 693 | tr.runGrid(testF, dat, &dat[l], true, 100); |
||
| 694 | qDebug() << "sum" << testf_sum << l*(l+1)/2. - testf_sum; |
||
| 695 | |||
| 696 | testf_sum = 0.; |
||
| 697 | tr.runGrid(testF, dat, &dat[l], true, 17); |
||
| 698 | qDebug() << "sum" << testf_sum << l*(l+1)/2. - testf_sum; |
||
| 699 | |||
| 700 | testf_sum = 0.; |
||
| 701 | tr.runGrid(testF, dat, &dat[l], true, 170000); |
||
| 702 | qDebug() << "sum" << testf_sum << l*(l+1)/2. - testf_sum; |
||
| 703 | |||
| 704 | testf_sum = 0.; |
||
| 705 | tr.runGrid(testF, dat, &dat[l], true); |
||
| 706 | qDebug() << "sum" << testf_sum << l*(l+1)/2. - testf_sum; |
||
| 707 | return; |
||
| 708 | |||
| 534 | werner | 709 | tme_exp.setExpression("(x+x+x+x)/(sqrt(x*x)*4)"); |
| 710 | tme_exp.setStrict(false); |
||
| 711 | try { |
||
| 712 | tme_exp.parse(); |
||
| 713 | } catch(const IException &e) { |
||
| 575 | werner | 714 | QString error_msg = e.message(); |
| 534 | werner | 715 | Helper::msg(error_msg); |
| 716 | qDebug() << error_msg; |
||
| 717 | } |
||
| 718 | QVector<double> data; |
||
| 719 | for (int i=1;i<1000000;i++) |
||
| 720 | data.push_back(nrandom(1,2000)); |
||
| 721 | qDebug() << "test multiple threads, part 1 (thread safe):"; |
||
| 722 | tme_count = 0; |
||
| 723 | { DebugTimer t("threadsafe"); |
||
| 724 | QtConcurrent::blockingMap(data,(*tme_test1)); } |
||
| 725 | qDebug() << "finished: errors: " << tme_count; |
||
| 726 | |||
| 727 | tme_exp.setStrict(false); |
||
| 728 | tme_exp.addVar("x"); |
||
| 729 | try { |
||
| 730 | tme_exp.parse(); |
||
| 731 | } catch(const IException &e) { |
||
| 575 | werner | 732 | QString error_msg = e.message(); |
| 534 | werner | 733 | Helper::msg(error_msg); |
| 734 | qDebug() << error_msg; |
||
| 735 | } |
||
| 736 | qDebug() << "test multiple threads, part 2 (non thread safe):"; |
||
| 737 | tme_count = 0; |
||
| 738 | { DebugTimer t("no_thread"); |
||
| 739 | QtConcurrent::blockingMap(data,(*tme_test2)); |
||
| 740 | } |
||
| 741 | qDebug() << "finished: errors: " << tme_count; |
||
| 742 | |||
| 743 | qDebug() << "test multiple threads, part 3 (mutex locker):"; |
||
| 744 | tme_count = 0; |
||
| 745 | { DebugTimer t("mutex"); |
||
| 746 | QtConcurrent::blockingMap(data,(*tme_test3)); |
||
| 747 | } |
||
| 748 | qDebug() << "finished: errors: " << tme_count; |
||
| 749 | } |
||
| 750 | |||
| 751 | void Tests::testEstablishment() |
||
| 752 | { |
||
| 753 | Model *model = GlobalSettings::instance()->model(); |
||
| 1158 | werner | 754 | //model->saplings()->clearStats(); |
| 1111 | werner | 755 | |
| 756 | { |
||
| 757 | DebugTimer test("test establishment"); |
||
| 758 | foreach (ResourceUnit *ru, model->ruList()) |
||
| 759 | model->saplings()->establishment(ru); |
||
| 760 | } |
||
| 1158 | werner | 761 | //qDebug() << "pixel tested" << model->saplings()->pixelTested() << "saplings added" << model->saplings()->saplingsAdded(); |
| 1111 | werner | 762 | |
| 1115 | werner | 763 | { |
| 764 | DebugTimer test("test sapling growth"); |
||
| 765 | foreach (ResourceUnit *ru, model->ruList()) |
||
| 766 | model->saplings()->saplingGrowth(ru); |
||
| 767 | } |
||
| 1158 | werner | 768 | //qDebug() << "pixel tested" << model->saplings()->pixelTested() << "saplings added" << model->saplings()->saplingsAdded(); |
| 1111 | werner | 769 | |
| 534 | werner | 770 | //model->ru(0) |
| 1111 | werner | 771 | //Establishment est(model->ru(0)->climate(),model->ru(0)->ruSpecies().first()); |
| 772 | //est.calculate(); |
||
| 534 | werner | 773 | } |
| 774 | |||
| 775 | void Tests::testLinearExpressions() |
||
| 776 | { |
||
| 777 | Expression::setLinearizationEnabled(true); // enable |
||
| 778 | Expression a("40*(1-(1-(x/40)^(1/3))*exp(-0.03))^3"); // test function: sapling growth |
||
| 779 | Expression b("40*(1-(1-(x/40)^(1/3))*exp(-0.03))^3"); |
||
| 780 | |||
| 781 | b.linearize(2., 140.); |
||
| 782 | b.calculate(58.); |
||
| 783 | for (double x=5.; x<150.; x++) { |
||
| 784 | double ra = a.calculate(x); |
||
| 785 | double rb = b.calculate(x); |
||
| 786 | qDebug() << x << ra << rb << rb-ra; |
||
| 787 | } |
||
| 788 | // test b |
||
| 789 | { DebugTimer t("straight"); |
||
| 790 | for (int i=0;i<1000000;i++) { |
||
| 791 | a.calculate(3.4); |
||
| 792 | } |
||
| 793 | } |
||
| 794 | { DebugTimer t("linear"); |
||
| 795 | for (int i=0;i<1000000;i++) { |
||
| 796 | b.calculate(3.4); |
||
| 797 | } |
||
| 798 | } |
||
| 799 | { DebugTimer t("nativ"); |
||
| 800 | double x=3; |
||
| 801 | double y=0; |
||
| 802 | for (int i=0;i<1000000;i++) { |
||
| 803 | y+= 40.*pow(1.-(1.-pow((x/40.),1./3.))*exp(-0.03),3.); |
||
| 804 | } |
||
| 805 | qDebug() << y; |
||
| 806 | } |
||
| 807 | /// test matrix-case |
||
| 808 | // Expression c("x*x+exp(y*y)"); |
||
| 809 | // Expression d("x*x+exp(y*y)"); |
||
| 810 | Expression c("exp(2*ln(x)*(1-y/2))"); |
||
| 811 | Expression d("exp(2*ln(x)*(1-y/2))"); |
||
| 812 | d.linearize2d(0., 1., 0., 1.); |
||
| 813 | qDebug() << d.calculate(0, 0.8); |
||
| 814 | for (double x=0.; x<1.; x+=0.1) { |
||
| 815 | for (double y=0.; y<1.; y+=0.1) { |
||
| 816 | double ra = c.calculate(x,y); |
||
| 817 | double rb = d.calculate(x,y); |
||
| 818 | qDebug() << x << y << ra << rb << rb-ra; |
||
| 819 | } |
||
| 820 | } |
||
| 821 | qDebug() << "random test"; |
||
| 822 | for (int i=0;i<100;i++) { |
||
| 823 | double x = drandom(); |
||
| 824 | double y = drandom(); |
||
| 825 | double ra = c.calculate(x,y); |
||
| 826 | double rb = d.calculate(x,y); |
||
| 827 | qDebug() << x << y << ra << rb << rb-ra; |
||
| 828 | |||
| 829 | } |
||
| 830 | qDebug() << "performance"; |
||
| 831 | { DebugTimer t("straight"); |
||
| 832 | for (int i=0;i<1000000;i++) { |
||
| 833 | c.calculate(0.33, 0.454); |
||
| 834 | } |
||
| 835 | } |
||
| 836 | { DebugTimer t("linear"); |
||
| 837 | for (int i=0;i<1000000;i++) { |
||
| 838 | d.calculate(0.33, 0.454); |
||
| 839 | } |
||
| 840 | } |
||
| 841 | } |
||
| 842 | |||
| 843 | void Tests::testSoil() |
||
| 844 | { |
||
| 845 | // Y0l=33 # young C labile, i.e. litter |
||
| 846 | // Y0Nl=0.45 # young N labile |
||
| 847 | // Y0r=247 # young C refractory , i.e. dwd |
||
| 848 | // Y0Nr=0.99 # young N refractory |
||
| 849 | // O0=375 # old C, i.e soil organic matter (SOM) |
||
| 850 | // O0N=15 # old N |
||
| 538 | werner | 851 | // kyl=0.15 # litter decomposition rate |
| 852 | // kyr=0.0807 # downed woody debris (dwd) decomposition rate |
||
| 534 | werner | 853 | Soil soil; |
| 854 | qDebug() << "soil object has" << sizeof(soil) << "bytes."; |
||
| 855 | // values from RS R script script_snagsoil_v5.r |
||
| 538 | werner | 856 | soil.setInitialState(CNPool(33000., 450., 0.15), // young lab |
| 857 | CNPool(247000., 990., 0.0807), // young ref |
||
| 534 | werner | 858 | CNPair(375000., 15000.)); // SOM |
| 538 | werner | 859 | CNPool in_l(3040, 3040/75., 0.15); |
| 860 | CNPool in_r(11970, 11970./250., 0.0807); |
||
| 534 | werner | 861 | double re = 1.1; |
| 862 | |||
| 863 | QStringList out; |
||
| 864 | out << "year;iLabC;iLabN;ikyl;iRefC;iRefN;ikyr;RE;kyl;kyr;ylC;ylN;yrC;yrN;somC;somN;NAvailable"; |
||
| 865 | |||
| 866 | for (int i=0;i<100;i++) { |
||
| 867 | // run the soil model |
||
| 868 | soil.setClimateFactor(re); |
||
| 869 | if (i==1) { |
||
| 870 | soil.setSoilInput(in_l*5., in_r*5.); // pulse in year 2 |
||
| 871 | } else if (i>1 && i<10) { |
||
| 872 | soil.setSoilInput(CNPool(), CNPool()); // no input for years 3..10 |
||
| 873 | } else{ |
||
| 874 | soil.setSoilInput(in_l, in_r); // normal input |
||
| 875 | } |
||
| 876 | |||
| 877 | soil.calculateYear(); |
||
| 878 | QList<QVariant> list = soil.debugList(); |
||
| 879 | QString line=QString::number(i)+";"; |
||
| 880 | foreach(QVariant v, list) |
||
| 881 | line+=v.toString() + ";"; |
||
| 882 | line.chop(1); |
||
| 883 | out << line; |
||
| 884 | } |
||
| 885 | if (fabs(soil.availableNitrogen() - 77.2196456179288) < 0.000001) |
||
| 886 | qDebug() << "PASSED."; |
||
| 887 | else |
||
| 888 | qDebug() << "ERROR"; |
||
| 889 | qDebug() << "ICBM/2N run: saved to e:/soil.txt"; |
||
| 890 | //qDebug() << out; |
||
| 891 | Helper::saveToTextFile("e:/soil.txt", out.join("\r\n")); |
||
| 892 | |||
| 893 | qDebug()<< "test weighting of CNPool (C/N/rate):"; |
||
| 894 | CNPool c(1,1,1); |
||
| 895 | for (int i=0;i<10;i++) { |
||
| 896 | c.addBiomass(1,1,2); |
||
| 897 | qDebug() << c.C << c.N << c.parameter(); |
||
| 898 | } |
||
| 899 | |||
| 545 | werner | 900 | out.clear(); |
| 901 | out << "year;iLabC;iLabN;ikyl;iRefC;iRefN;ikyr;RE;kyl;kyr;ylC;ylN;yrC;yrN;somC;somN;NAvailable"; |
||
| 902 | // test soil with debug output from iLand ... |
||
| 903 | CSVFile dbg_iland(GlobalSettings::instance()->path("debug_carboncycle.csv")); |
||
| 904 | Soil *model_soil = GlobalSettings::instance()->model()->ru()->soil(); // should now have the 'right' parameters... |
||
| 905 | |||
| 906 | model_soil->setInitialState(CNPool(dbg_iland.value(0, dbg_iland.columnIndex("ylC")).toDouble()*1000., |
||
| 907 | dbg_iland.value(0, dbg_iland.columnIndex("ylN")).toDouble()*1000., |
||
| 908 | dbg_iland.value(0, dbg_iland.columnIndex("kyl")).toDouble()), // young lab |
||
| 909 | CNPool(dbg_iland.value(0, dbg_iland.columnIndex("yrC")).toDouble()*1000., |
||
| 910 | dbg_iland.value(0, dbg_iland.columnIndex("yrN")).toDouble()*1000., |
||
| 911 | dbg_iland.value(0, dbg_iland.columnIndex("kyr")).toDouble()), // young ref |
||
| 912 | CNPair(dbg_iland.value(0, dbg_iland.columnIndex("somC")).toDouble()*1000., |
||
| 913 | dbg_iland.value(0, dbg_iland.columnIndex("somN")).toDouble()*1000.)); // SOM |
||
| 914 | |||
| 915 | for (int i=1;i<dbg_iland.rowCount();i++) { |
||
| 916 | // run the soil model |
||
| 917 | model_soil->setClimateFactor(dbg_iland.value(i, dbg_iland.columnIndex("re")).toDouble()); |
||
| 918 | model_soil->setSoilInput(CNPool(dbg_iland.value(i, dbg_iland.columnIndex("iLabC")).toDouble()*1000., |
||
| 919 | dbg_iland.value(i, dbg_iland.columnIndex("iLabN")).toDouble()*1000., |
||
| 920 | dbg_iland.value(i, dbg_iland.columnIndex("iKyl")).toDouble()), |
||
| 921 | CNPool(dbg_iland.value(i, dbg_iland.columnIndex("iRefC")).toDouble()*1000., |
||
| 922 | dbg_iland.value(i, dbg_iland.columnIndex("iRefN")).toDouble()*1000., |
||
| 923 | dbg_iland.value(i, dbg_iland.columnIndex("iKyr")).toDouble())); |
||
| 924 | model_soil->calculateYear(); |
||
| 925 | QList<QVariant> list = model_soil->debugList(); |
||
| 926 | QString line=QString::number(i)+";"; |
||
| 927 | foreach(QVariant v, list) |
||
| 928 | line+=v.toString() + ";"; |
||
| 929 | line.chop(1); |
||
| 930 | out << line; |
||
| 931 | } |
||
| 932 | Helper::saveToTextFile("e:/soil2.txt", out.join("\r\n")); |
||
| 534 | werner | 933 | } |
| 539 | werner | 934 | |
| 935 | void Tests::testMap() |
||
| 936 | { |
||
| 937 | QString fileName = GlobalSettings::instance()->path("gis/montafon_grid.txt"); |
||
| 543 | werner | 938 | |
| 544 | werner | 939 | int test_id = 127; |
| 543 | werner | 940 | // test the map-grid class |
| 941 | MapGrid map(fileName); |
||
| 544 | werner | 942 | QList<ResourceUnit*> ru_list = map.resourceUnits(test_id); |
| 943 | qDebug() << "total bounding box" << map.boundingBox(test_id) << "has an area of" << map.area(test_id); |
||
| 543 | werner | 944 | foreach (const ResourceUnit *ru, ru_list) |
| 945 | qDebug() << ru->index() << ru->boundingBox(); |
||
| 946 | |||
| 803 | werner | 947 | qDebug() << "neighbors of" << test_id << "are" << map.neighborsOf(test_id); |
| 948 | return; |
||
| 539 | werner | 949 | HeightGrid *hgrid = GlobalSettings::instance()->model()->heightGrid(); |
| 543 | werner | 950 | for (int i=0;i<map.grid().count();i++) |
| 951 | hgrid->valueAtIndex(i).height = map.grid().constValueAtIndex(map.grid().indexOf(i)); |
||
| 544 | werner | 952 | |
| 953 | QList<Tree*> tree_list = map.trees(test_id); |
||
| 954 | qDebug() << "no of trees:" << tree_list.count(); |
||
| 955 | Management mgmt; |
||
| 956 | |||
| 957 | mgmt.loadFromTreeList(tree_list); |
||
| 824 | werner | 958 | mgmt.killAll(); |
| 539 | werner | 959 | } |
| 642 | werner | 960 | |
| 961 | DEM *_dem = 0; |
||
| 962 | void Tests::testDEM() |
||
| 963 | { |
||
| 964 | QString fileName = GlobalSettings::instance()->path("gis/dtm1m_clip.txt"); |
||
| 965 | |||
| 966 | if (_dem) { |
||
| 967 | int choice = QInputDialog::getInt(0, "enter type", "Type to show: 0: dem, 1: slope, 2: aspect, 3: view, 4: exit"); |
||
| 968 | switch (choice) { |
||
| 649 | werner | 969 | case 0: GlobalSettings::instance()->controller()->addGrid(_dem, "dem height", GridViewRainbow, 0, 1000); break; |
| 970 | case 1: GlobalSettings::instance()->controller()->addGrid(_dem->slopeGrid(), "slope", GridViewRainbow, 0, 3); break; |
||
| 971 | case 2: GlobalSettings::instance()->controller()->addGrid(_dem->aspectGrid(), "aspect", GridViewRainbow, 0, 360); break; |
||
| 972 | case 3: GlobalSettings::instance()->controller()->addGrid(_dem->viewGrid(), "dem", GridViewGray, 0, 1); break; |
||
| 642 | werner | 973 | default: return; |
| 974 | } |
||
| 975 | return; |
||
| 976 | } |
||
| 977 | |||
| 978 | try { |
||
| 979 | if (!_dem) |
||
| 980 | _dem = new DEM(fileName); |
||
| 981 | |||
| 982 | qDebug() << "slope grid: avg max" << _dem->slopeGrid()->avg() << _dem->slopeGrid()->max(); |
||
| 983 | qDebug() << "aspect grid: avg max" << _dem->aspectGrid()->avg() << _dem->aspectGrid()->max(); |
||
| 984 | qDebug() << "view grid: avg max" << _dem->viewGrid()->avg() << _dem->viewGrid()->max(); |
||
| 985 | |||
| 986 | // note: dem gets not released! (to be still there when painting happens) |
||
| 987 | float slope, aspect; |
||
| 988 | for (float y=0.;y<1000.;y+=100.) |
||
| 989 | for (float x=0.;x<1000.;x+=100.) { |
||
| 990 | float h = _dem->orientation(x, y, slope, aspect ); |
||
| 991 | qDebug() << "at point "<< x << y << "height"<< h <<"slope" << slope << "aspect" << aspect; |
||
| 992 | } |
||
| 993 | |||
| 994 | } catch (IException &e) { |
||
| 995 | Helper::msg(e.message()); |
||
| 996 | } |
||
| 649 | werner | 997 | GlobalSettings::instance()->controller()->addGrid(_dem, "dem height", GridViewRainbow, 0, 1000); |
| 998 | GlobalSettings::instance()->controller()->addGrid(_dem->slopeGrid(), "slope", GridViewRainbow, 0, 3); |
||
| 999 | GlobalSettings::instance()->controller()->addGrid(_dem->aspectGrid(), "aspect", GridViewRainbow, 0, 360); |
||
| 1000 | GlobalSettings::instance()->controller()->addGrid(_dem->viewGrid(), "dem", GridViewGray, 0, 1); |
||
| 642 | werner | 1001 | |
| 643 | werner | 1002 | |
| 642 | werner | 1003 | } |
| 652 | werner | 1004 | |
| 1005 | void Tests::testFire() |
||
| 1006 | { |
||
| 1007 | // get fire module |
||
| 761 | werner | 1008 | FirePlugin *plugin =dynamic_cast<FirePlugin *>(GlobalSettings::instance()->model()->modules()->module("fire")); |
| 652 | werner | 1009 | if (plugin) { |
| 1010 | FireModule *fire = plugin->fireModule(); |
||
| 655 | werner | 1011 | try { |
| 1012 | fire->testSpread(); |
||
| 1013 | } catch (const IException &e) { |
||
| 1014 | Helper::msg(e.message()); |
||
| 1015 | } |
||
| 652 | werner | 1016 | } |
| 1017 | GlobalSettings::instance()->controller()->repaint(); |
||
| 1018 | |||
| 1019 | } |
||
| 702 | werner | 1020 | |
| 1021 | void Tests::testWind() |
||
| 1022 | { |
||
| 1023 | // get fire module |
||
| 1024 | WindPlugin *plugin = dynamic_cast<WindPlugin *>(GlobalSettings::instance()->model()->modules()->module("wind")); |
||
| 1025 | if (plugin) { |
||
| 715 | werner | 1026 | static int iteration = 0; |
| 702 | werner | 1027 | WindModule *wind = plugin->windModule(); |
| 715 | werner | 1028 | if (iteration>0) { |
| 716 | werner | 1029 | if (Helper::question("next iteration (yes) or stop (no)?")) { |
| 1030 | qDebug() << ".... iteration ..." << iteration+1; |
||
| 1031 | wind->run(iteration++); |
||
| 1032 | return; |
||
| 1033 | } else { |
||
| 1034 | iteration = 0; |
||
| 1035 | } |
||
| 715 | werner | 1036 | } |
| 702 | werner | 1037 | double direction = Helper::userValue("direction (degrees), -1 for roundabout", "100").toDouble(); |
| 1038 | bool loop = false; |
||
| 1039 | if (direction==-1) |
||
| 1040 | loop = true; |
||
| 1041 | try { |
||
| 1042 | if (loop) { |
||
| 1043 | for (direction = 0.; direction<360; direction+=10.) { |
||
| 1044 | wind->run(); |
||
| 1045 | wind->testFetch(direction); |
||
| 1046 | GlobalSettings::instance()->controller()->repaint(); |
||
| 1047 | qDebug() << "looping..." << direction; |
||
| 1048 | } |
||
| 1049 | } else { |
||
| 712 | werner | 1050 | DebugTimer t; |
| 715 | werner | 1051 | |
| 713 | werner | 1052 | if (direction!=0) |
| 1053 | wind->setSimulationMode(true); |
||
| 712 | werner | 1054 | double speed = Helper::userValue("wind speed (m/2)", "300").toDouble(); |
| 1055 | wind->setWindProperties(direction * M_PI / 180., speed); |
||
| 715 | werner | 1056 | wind->run(iteration++); |
| 712 | werner | 1057 | |
| 702 | werner | 1058 | qDebug() << "fetch finished. ms:" << t.elapsed(); |
| 1059 | } |
||
| 1060 | } catch (const IException &e) { |
||
| 1061 | Helper::msg(e.message()); |
||
| 1062 | } |
||
| 1063 | |||
| 1064 | } |
||
| 1065 | GlobalSettings::instance()->controller()->repaint(); |
||
| 1066 | } |
||
| 766 | werner | 1067 | |
| 1068 | void Tests::testRumple() |
||
| 1069 | { |
||
| 1070 | RumpleIndex rumple_index; |
||
| 1071 | double result = rumple_index.value(); |
||
| 1072 | qDebug() << "rumple index test: " << result; |
||
| 1073 | |||
| 1074 | Helper::saveToTextFile("rumple_test.txt", gridToESRIRaster(rumple_index.rumpleGrid()) ); |
||
| 1075 | qDebug() << "surface area test triangle" << rumple_index.test_triangle_area(); |
||
| 1076 | |||
| 1077 | } |
||
| 808 | werner | 1078 | |
| 905 | werner | 1079 | ABE::ForestManagementEngine *fome=0; |
| 808 | werner | 1080 | void Tests::testFOMEsetup() |
| 1081 | { |
||
| 914 | werner | 1082 | |
| 1083 | fome = GlobalSettings::instance()->model()->ABEngine(); |
||
| 1084 | if (!fome) |
||
| 1085 | fome = new ABE::ForestManagementEngine(); |
||
| 873 | werner | 1086 | //fome.test(); |
| 1087 | try { |
||
| 876 | werner | 1088 | |
| 905 | werner | 1089 | ABE::FMSTP::setVerbose(true); |
| 876 | werner | 1090 | fome->setup(); |
| 873 | werner | 1091 | } catch(const IException &e) { |
| 1092 | Helper::msg(e.message()); |
||
| 1093 | } |
||
| 1094 | // todo: re-enable delete!!! |
||
| 1095 | // delete fome; |
||
| 808 | werner | 1096 | } |
| 876 | werner | 1097 | |
| 1098 | void Tests::testFOMEstep() |
||
| 1099 | { |
||
| 1100 | if (!fome) |
||
| 1101 | return; |
||
| 1102 | int n = Helper::userValue("how many years?", "1").toInt(); |
||
| 1103 | for (int i=0;i<n;++i) { |
||
| 909 | werner | 1104 | qDebug()<< "running ABE year" << i; |
| 876 | werner | 1105 | fome->run(1); |
| 1106 | } |
||
| 1107 | } |
||
| 1002 | werner | 1108 | |
| 1109 | void Tests::testDbgEstablishment() |
||
| 1110 | { |
||
| 1111 | // test speeed of debugtimer |
||
| 1112 | DebugTimer::clearAllTimers(); |
||
| 1113 | DebugTimer total("total"); |
||
| 1114 | // 4x1,000,000 timers: ca 1.5 seks... |
||
| 1115 | for (int i=0;i<100000;i++) { |
||
| 1116 | { DebugTimer a("aasdf asdfasdf"); } |
||
| 1117 | { DebugTimer a("basdflk asdf"); } |
||
| 1118 | { DebugTimer a("bas asdfasdf"); } |
||
| 1119 | { DebugTimer a("dasdflkj asdfalskfj"); } |
||
| 1120 | } |
||
| 1121 | |||
| 1122 | qDebug() << "finsihed timers"; |
||
| 1123 | |||
| 1162 | werner | 1124 | // const Grid<float> &seed_map = GlobalSettings::instance()->model()->speciesSet()->activeSpecies()[0]->seedDispersal()->seedMap(); |
| 1002 | werner | 1125 | |
| 1162 | werner | 1126 | // int n_established = 0; |
| 1127 | // int n_tested = 0, n_seed=0, n_dropped=0; |
||
| 1128 | // double mPAbiotic = 0.8; |
||
| 1129 | // DebugTimer runt("run est"); |
||
| 1130 | // int dbg_numbers = RandomGenerator::debugNRandomNumbers(); |
||
| 1131 | // // define a height map for the current resource unit on the stack |
||
| 1132 | // float sapling_map[cPxPerRU*cPxPerRU]; |
||
| 1002 | werner | 1133 | // set the map and initialize it: |
| 1134 | |||
| 1162 | werner | 1135 | // for (int i=0;i<100;i++) |
| 1136 | // foreach ( ResourceUnit *ru, GlobalSettings::instance()->model()->ruList()) { |
||
| 1002 | werner | 1137 | |
| 1162 | werner | 1138 | // ru->setSaplingHeightMap(sapling_map); |
| 1139 | // ResourceUnitSpecies *mRUS = ru->ruSpecies()[0]; |
||
| 1002 | werner | 1140 | |
| 1162 | werner | 1141 | // const QRectF &ru_rect = ru->boundingBox(); |
| 1002 | werner | 1142 | |
| 1162 | werner | 1143 | // // test the establishment code.... |
| 1144 | // GridRunner<float> seed_runner(seed_map, ru_rect); |
||
| 1145 | // Grid<float> *lif_map = GlobalSettings::instance()->model()->grid(); |
||
| 1002 | werner | 1146 | |
| 1162 | werner | 1147 | // // check every pixel inside the bounding box of the pixel with |
| 1148 | // while (float *p = seed_runner.next()) { |
||
| 1149 | // ++n_tested; |
||
| 1150 | // if (*p>0.f) { |
||
| 1151 | // ++n_seed; |
||
| 1152 | // //double p_establish = drandom(mRUS->ru()->randomGenerator()); |
||
| 1153 | // //if (p_establish > mPAbiotic) |
||
| 1154 | // // continue; |
||
| 1155 | // // pixel with seeds: now really iterate over lif pixels |
||
| 1156 | // GridRunner<float> lif_runner(lif_map, seed_map.cellRect(seed_runner.currentIndex())); |
||
| 1157 | // while (float *lif_px = lif_runner.next()) { |
||
| 1158 | // DBGMODE( |
||
| 1159 | // if (!ru_rect.contains(lif_map->cellCenterPoint(lif_map->indexOf(lif_px)))) |
||
| 1160 | // qDebug() << "(b) establish problem:" << lif_map->indexOf(lif_px) << "point: " << lif_map->cellCenterPoint(lif_map->indexOf(lif_px)) << "not in" << ru_rect; |
||
| 1161 | // ); |
||
| 1162 | // double p_establish = drandom(); |
||
| 1163 | // if (p_establish < mPAbiotic) { |
||
| 1164 | // //if (establishTree(lif_map->indexOf(lif_px), *lif_px ,*p)) |
||
| 1165 | // QPoint pos_lif = lif_map->indexOf(lif_px); |
||
| 1002 | werner | 1166 | |
| 1162 | werner | 1167 | // if (ru->saplingHeightAt(pos_lif) > 1.3f) |
| 1168 | // ++n_dropped; |
||
| 1002 | werner | 1169 | |
| 1162 | werner | 1170 | // // check if sapling of the current tree species is already established -> if so, no establishment. |
| 1171 | // if (mRUS->hasSaplingAt(pos_lif)) |
||
| 1172 | // ++n_dropped; |
||
| 1002 | werner | 1173 | |
| 1162 | werner | 1174 | // const HeightGridValue &hgv = GlobalSettings::instance()->model()->heightGrid()->constValueAtIndex(pos_lif.x()/cPxPerHeight, pos_lif.y()/cPxPerHeight); |
| 1175 | // if (hgv.count()>0) n_dropped++; |
||
| 1176 | // double h_height_grid = hgv.height; |
||
| 1177 | // double rel_height = 4. / h_height_grid; |
||
| 1002 | werner | 1178 | |
| 1162 | werner | 1179 | // double lif_corrected = mRUS->species()->speciesSet()->LRIcorrection(*lif_px, rel_height); |
| 1180 | // if (lif_corrected < drandom()) |
||
| 1181 | // ++n_dropped; |
||
| 1002 | werner | 1182 | |
| 1183 | |||
| 1162 | werner | 1184 | // n_established++; |
| 1185 | // } |
||
| 1186 | // } |
||
| 1187 | // } |
||
| 1188 | // } |
||
| 1189 | // } |
||
| 1190 | // qDebug() << "Orig tested " << n_tested << "with seed" << n_seed << "est: " << n_established << "time" << runt.elapsed() << "debug numbers used:" << RandomGenerator::debugNRandomNumbers()-dbg_numbers << "n_dropped" << n_dropped; |
||
| 1002 | werner | 1191 | |
| 1162 | werner | 1192 | // runt.start(); |
| 1193 | // n_established = 0; |
||
| 1194 | // n_tested = 0; n_seed=0; n_dropped=0; |
||
| 1195 | // for (int i=0;i<100;i++) |
||
| 1196 | // foreach ( ResourceUnit *ru, GlobalSettings::instance()->model()->ruList()) { |
||
| 1002 | werner | 1197 | |
| 1162 | werner | 1198 | // ru->setSaplingHeightMap(sapling_map); |
| 1199 | // ResourceUnitSpecies *mRUS = ru->ruSpecies()[0]; |
||
| 1200 | // const std::bitset<cPxPerRU*cPxPerRU> &pos_bitset = mRUS->sapling().presentPositions(); |
||
| 1002 | werner | 1201 | |
| 1162 | werner | 1202 | // const QRectF &ru_rect = ru->boundingBox(); |
| 1203 | // //DebugTimer estasd("establish:from_search"); |
||
| 1204 | // // 3rd step: check actual pixels in the LIF grid |
||
| 1002 | werner | 1205 | |
| 1162 | werner | 1206 | // // a large part has available seeds. simply scan the pixels... |
| 1207 | // QPoint lif_index; |
||
| 1208 | // Grid<float> *lif_map = GlobalSettings::instance()->model()->grid(); |
||
| 1209 | // GridRunner<float> lif_runner(lif_map, ru_rect); |
||
| 1210 | // float *sap_height = sapling_map; |
||
| 1211 | // size_t bit_idx=0; |
||
| 1212 | // while (float *lif_px = lif_runner.next()) { |
||
| 1213 | // ++n_tested; |
||
| 1214 | // // check for height of sapling < 1.3m (for all species |
||
| 1215 | // // and for presence of a sapling of the given species |
||
| 1216 | // if (*sap_height < 1.3 && pos_bitset[bit_idx]==false) { |
||
| 1217 | // ++n_seed; |
||
| 1218 | // lif_index = lif_map->indexOf(lif_px); |
||
| 1219 | // double sap_random_number = drandom(); |
||
| 1220 | // const float seed_map_value = seed_map.constValueAt( lif_runner.currentCoord() ); |
||
| 1002 | werner | 1221 | |
| 1162 | werner | 1222 | // if (seed_map_value*mPAbiotic > sap_random_number) |
| 1223 | // ++n_dropped; |
||
| 1002 | werner | 1224 | |
| 1225 | |||
| 1162 | werner | 1226 | //// DBGMODE( |
| 1227 | //// if (!ru_rect.contains(lif_map->cellCenterPoint(lif_index))) |
||
| 1228 | //// qDebug() << "(a) establish problem:" << lif_index << "point: " << lif_map->cellCenterPoint(lif_index) << "not in" << ru_rect; |
||
| 1229 | //// ); |
||
| 1002 | werner | 1230 | |
| 1231 | |||
| 1162 | werner | 1232 | // const HeightGridValue &hgv = GlobalSettings::instance()->model()->heightGrid()->constValueAtIndex(lif_index.x()/cPxPerHeight, lif_index.y()/cPxPerHeight); |
| 1002 | werner | 1233 | |
| 1162 | werner | 1234 | // double h_height_grid = hgv.height; |
| 1235 | // double rel_height = 4. / h_height_grid; |
||
| 1002 | werner | 1236 | |
| 1162 | werner | 1237 | // double lif_corrected = mRUS->species()->speciesSet()->LRIcorrection(*lif_px, rel_height); |
| 1238 | // if (lif_corrected < sap_random_number) |
||
| 1239 | // ++n_dropped; |
||
| 1002 | werner | 1240 | |
| 1162 | werner | 1241 | // if (seed_map_value*mPAbiotic*lif_corrected < sap_random_number) |
| 1242 | // n_established++; |
||
| 1243 | // } |
||
| 1244 | // ++sap_height; |
||
| 1245 | // ++bit_idx; |
||
| 1246 | // } |
||
| 1002 | werner | 1247 | |
| 1248 | |||
| 1162 | werner | 1249 | // } |
| 1250 | // qDebug() << "Upd tested " << n_tested << "with seed" << n_seed << "est: " << n_established << "time" << runt.elapsed() << "debug numbers used:" << RandomGenerator::debugNRandomNumbers()-dbg_numbers << "n_dropped" << n_dropped; |
||
| 1002 | werner | 1251 | |
| 1162 | werner | 1252 | // DebugTimer::printAllTimers(); |
| 1002 | werner | 1253 | } |
| 1067 | werner | 1254 | |
| 1255 | void Tests::testGridIndexHack() |
||
| 1256 | { |
||
| 1257 | Grid<float> m2, m10; |
||
| 1258 | const int n=10000; |
||
| 1259 | m2.setup(2,n,n); |
||
| 1260 | m10.setup(10, n/5, n/5 ); |
||
| 1261 | m10.wipe(); |
||
| 1262 | |||
| 1263 | for (float *f = m2.begin(), i=0.; f!=m2.end(); ++f, ++i) |
||
| 1264 | *f = i; |
||
| 1265 | |||
| 1266 | // for (int i=0;i<m2.count();++i) { |
||
| 1267 | // qDebug() << i << m2.index5(i); |
||
| 1268 | // } |
||
| 1269 | // return; |
||
| 1270 | |||
| 1271 | int el; |
||
| 1272 | { DebugTimer t("optimized"); |
||
| 1273 | for (float *f = m2.begin(); f!=m2.end(); ++f) { |
||
| 1274 | m10[ m2.index5(f-m2.begin()) ] += *f; |
||
| 1275 | } |
||
| 1276 | el = t.elapsed(); |
||
| 1277 | } |
||
| 1278 | |||
| 1279 | float s=m10.sum() / m10.count(); |
||
| 1280 | qDebug() << "test average value (new):" << s << "time" << el; |
||
| 1281 | |||
| 1282 | m10.wipe(); |
||
| 1283 | |||
| 1284 | { DebugTimer t("old"); |
||
| 1285 | for (float *f = m2.begin(); f!=m2.end(); ++f) { |
||
| 1286 | m10.valueAt(m2.cellCenterPoint(m2.indexOf(f))) += *f; |
||
| 1287 | } |
||
| 1288 | el = t.elapsed(); |
||
| 1289 | } |
||
| 1290 | s=m10.sum() / m10.count(); |
||
| 1291 | qDebug() << "test average value (old):" << s << "time" << el; |
||
| 1292 | |||
| 1293 | int errors=0; |
||
| 1294 | { DebugTimer t("compare"); |
||
| 1295 | for (float *f = m2.begin(); f!=m2.end(); ++f) { |
||
| 1296 | if ( m10.valueAtIndex( m2.index5(f-m2.begin()) ) != m10.valueAt(m2.cellCenterPoint(m2.indexOf(f)))) |
||
| 1297 | errors++; |
||
| 1298 | } |
||
| 1299 | } |
||
| 1300 | qDebug() << "test e:differenczes" << errors; |
||
| 1301 | |||
| 1302 | { DebugTimer t("optimized"); |
||
| 1303 | for (int i=0;i<m2.count();++i) { |
||
| 1304 | m10[ m2.index5(i)] += m2[i]; |
||
| 1305 | } |
||
| 1306 | el = t.elapsed(); |
||
| 1307 | } |
||
| 1308 | |||
| 1309 | s=m10.sum() / m10.count(); |
||
| 1310 | qDebug() << "test average value (square brackets):" << s << "time" << el; |
||
| 1311 | |||
| 1312 | } |