windmodule.h

Go to the documentation of this file.

/********************************************************************************************
**    iLand - an individual based forest landscape and disturbance model
**    http://iland-model.org
**    Copyright (C) 2009-  Werner Rammer, Rupert Seidl
**
**    This program is free software: you can redistribute it and/or modify
**    it under the terms of the GNU General Public License as published by
**    the Free Software Foundation, either version 3 of the License, or
**    (at your option) any later version.
**
**    This program is distributed in the hope that it will be useful,
**    but WITHOUT ANY WARRANTY; without even the implied warranty of
**    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
**    GNU General Public License for more details.
**
**    You should have received a copy of the GNU General Public License
**    along with this program.  If not, see <http://www.gnu.org/licenses/>.
********************************************************************************************/
#ifndef WINDMODULE_H
#define WINDMODULE_H
#include "grid.h"
#include "layeredgrid.h"
#include "expression.h"
#include <QHash>
 
class Tree; // forward
class Species; // forward
class ResourceUnit; // forward
 
class WindCell {
public:
    WindCell() { topex = 0; n_affected=0; sum_volume_killed=0.; edge_age=0; clear(); }
    void clear() {height = edge = 0.f; n_trees=0.f; tree=0; n_killed = 0; basal_area_killed = 0.f; cws_uproot = 0.; cws_break= crown_windspeed= 0.; n_iteration = 0;}
    bool isValid() const { return height<9999.f; } 
    float topex; 
    float height; 
    float n_trees; 
    const Tree *tree; 
    float edge; 
    // statistics
    int n_iteration; 
    int n_killed; 
    int edge_age; 
    double basal_area_killed; 
    double cws_uproot; 
    double cws_break; 
    double crown_windspeed; 
    int n_affected; 
    double sum_volume_killed; 
};
// data structure for a resource unit
class WindRUCell {
public:
    WindRUCell(): flag(0), soilIsFrozen(false) {}
    int flag; // flag to indicate that the trees of the current resource are already processed
    bool soilIsFrozen;
};
 
class WindLayers: public LayeredGrid<WindCell> {
  public:
    void setGrid(const Grid<WindCell> &grid) { mGrid = &grid; }
    double value(const WindCell& data, const int index) const;
    const QVector<LayerElement> &names();
    // specifics for wind layers
    void setRUGrid(const Grid<WindRUCell> *grid) { mRUGrid = grid; }
private:
    QVector<LayerElement> mNames;
    double ruValueAt(const WindCell *cell, const int what) const; // get resource-unit level factor at cell "cell"
    const Grid<WindRUCell> *mRUGrid;
};
 
struct WindSpeciesParameters
{
    WindSpeciesParameters():  crown_area_factor(0.5), crown_length(0.5), Creg(111), MOR(30.6), wet_biomass_factor(1.86) {}
    double crown_area_factor; // empirical factor related to the crown shape (fraction of crown shape compared to rectangle)
    double crown_length; // crown length of the tree (fraction of tree height)
    double Creg; // Nm/kg, critical turning coefficient from tree pulling
    double MOR; // MPa, modulus of rupture
    double wet_biomass_factor; // conversion factor between dry and wet biomass (wet = dry*factor)
};
 
class WindModule
{
public:
    WindModule();
    static double cellsize() { return 10.; }
    void setup();
    void setupResourceUnit(const ResourceUnit* ru);
    WindLayers &layers() { return mWindLayers; }
 
    void run(const int iteration=-1, const bool execute_from_script=false);
 
    // test functions
    void setWindProperties(const double direction_rad, const double speed_ms) { mWindDirection = direction_rad; mWindSpeed = speed_ms; }
    void setSimulationMode(const bool mode) { mSimulationMode = mode; }
    bool simulationMode() const { return mSimulationMode; }
    void setMaximumIterations(const double maxit) { mMaxIteration = maxit; }
 
    void testFetch(double degree_direction);
    void testEffect();
private:
    // main functions
    bool eventTriggered(); 
    void initWindGrid(); 
    void detectEdges(bool at_startup=false); 
    void calculateFetch(); 
    int calculateWindImpact(); 
 
    // details
    bool checkFetch(const int startx, const int starty, const double direction, const double max_distance, const double threshold) ;
    bool windImpactOnPixel(const QPoint position, WindCell *cell);
    double calculateCrownWindSpeed(const Tree *tree, const WindSpeciesParameters &params, const double n_trees, const double wind_speed_10);
    double calculateCrititalWindSpeed(const Tree *tree, const WindSpeciesParameters &params, const double gap_length, double &rCWS_uproot, double &rCWS_break);
    bool isSoilFrozen(const ResourceUnit *ru, const int day_of_year) const;
 
    // helping functions
    void scanResourceUnitTrees(const QPoint &position);
    void loadSpeciesParameter(const QString &table_name);
    const WindSpeciesParameters &speciesParameter(const Species *s);
    void initializeEdgeAge(const int years);
    void incrementEdgeAge();
 
    // variables
    bool mTopexFromGrid; 
    double mWindDirection; 
    double mWindDirectionVariation; 
    double mWindSpeed; 
    double mEdgeDetectionThreshold; 
    double mFactorEdge; 
    int mWindDayOfYear; 
    bool mSimulationMode; 
    int mCurrentIteration; 
    int mMaxIteration; 
    double mGustModifier; 
    double mCurrentGustFactor; 
    enum ETopexFactorModificationType {gfMultiply, gfAdd} mTopexFactorModificationType; 
    double mIterationsPerMinute; 
    int mEdgeAgeBaseValue; 
    Expression mEdgeProbability; 
    double mEdgeBackgroundProbability; 
    // some statistics
    int mPixelAffected; 
    int mTreesKilled; 
    double mTotalKilledBasalArea; 
    double mTotalKilledVolume; 
    enum ESoilFreezeMode {esfFrozen, esfNotFrozen, esfAuto, esfInvalid} mSoilFreezeMode; 
    Grid<WindCell> mGrid; 
    Grid<WindRUCell> mRUGrid; 
    WindLayers mWindLayers; 
    // species parameters for the wind module
    QHash<const Species*, WindSpeciesParameters> mSpeciesParameters;
    Expression mLRITransferFunction;
 
    QString mAfterExecEvent;
 
    friend class WindScript;
    friend class WindOut;
    friend void nc_calculateFetch(WindCell *begin, WindCell *end);
    friend void nc_calculateWindImpact(ResourceUnit *unit);
 
};
 
 
 
#endif // WINDMODULE_H