Subversion Repositories public iLand

/********************************************************************************************

**    iLand - an individual based forest landscape and disturbance model

**    http://iland.boku.ac.at

**    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/>.

********************************************************************************************/

#include "abe_global.h"

#include "global.h"

#include "fmunit.h"

#include "forestmanagementengine.h"

#include "fmstand.h"

#include "scheduler.h"

#include "agent.h"

#include "agenttype.h"

#include "fomescript.h"

#include "fmtreelist.h"

namespace ABE {

/** @class FMUnit

    @ingroup abe

    The FMUnit class encapsulates a forest management unit, comprised of forest stands. Units are the base level at which

    the scheduling works.

  */

void FMUnit::aggregate()

{

    // loop over all stands

    // collect some data....

    double age=0.;

    double volume = 0.;

    double harvest = 0.;

    double totalarea = 0.;

    const QMultiMap<FMUnit*, FMStand*> &stands = ForestManagementEngine::instance()->stands();

    QMultiMap<FMUnit*, FMStand*>::const_iterator it = stands.constFind(this);

    while (it != stands.constEnd() && it.key()==this) {

        const FMStand *s = it.value();

        age += s->age() * s->area();

        volume += s->volume() * s->area();

        totalarea += s->area();

        ++it;

    }

    if (totalarea>0.) {

        age /= totalarea;

        volume /= totalarea;

        harvest /= totalarea;

    }

    qCDebug(abe) << "unit" << id() << "volume (m3/ha)" << volume << "age" << age << "planned harvest: todo";

}

QStringList FMUnit::info() const

{

    return QStringList() << QString("(accumulated) harvest: %1").arg(mRealizedHarvest)

                         << QString("MAI: %1").arg(mMAI)

                         << QString("HDZ: %1").arg(mHDZ)

                         << QString("average age: %1").arg(mMeanAge)

                         << QString("decadal plan: %1").arg(mAnnualHarvestTarget)

                         << QString("current plan: %1").arg(constScheduler()!=0?constScheduler()->harvestTarget():0.);

}

double FMUnit::annualThinningHarvest() const

{

    const QMultiMap<FMUnit*, FMStand*> &stands = ForestManagementEngine::instance()->stands();

    QMultiMap<FMUnit*, FMStand*>::const_iterator it = stands.constFind(const_cast<FMUnit*>(this));

    double harvested=0.;

    while (it != stands.constEnd() && it.key()==this) {

        FMStand *stand = it.value();

        harvested += stand->totalThinningHarvest();

        ++it;

    }

    return harvested;

}

FMUnit::FMUnit(const Agent *agent)

{

    mAgent = agent;

    mScheduler = 0;

    mAnnualHarvestTarget = -1.;

    mRealizedHarvest = 0.;

    mMAI = 0.; mHDZ = 0.; mMeanAge = 0.;

    mTotalArea = 0.; mTotalPlanDeviation = 0.;

    mTotalVolume = 0.;

    mAnnualHarvest = 0.;

    mNumberOfStands = 0;

    mU = 100, mThinningIntensityClass = 2, mSpeciesCompositionIndex = 0;

    mAverageMAI = 0.;

    //if (agent->type()->schedulerOptions().useScheduler)

    // explicit scheduler only for stands/units that include more than one stand

    mScheduler = new Scheduler(this);

}

FMUnit::~FMUnit()

{

    if (mScheduler)

        delete mScheduler;

}

void FMUnit::setId(const QString &id)

{

    mId = id;

}

void FMUnit::resetHarvestCounter()

{

    if (scheduler())

        scheduler()->resetHarvestCounter();

}

void FMUnit::managementPlanUpdate()

{

    const double period_length = 10.;

    // calculate the planned harvest in the next planning period (i.e., 10yrs).

    // this is the sum of planned operations that are already in the scheduler.

    double plan_thinning, plan_final;

    mScheduler->plannedHarvests(plan_final, plan_thinning);

    // the actual harvests of the last planning period

    //double realized = mRealizedHarvest;

    mRealizedHarvest = 0.; // reset

    mRealizedHarvestLastYear = 0.;

    // preparations:

    // MAI-calculation for all stands:

    double total_area = 0.;

    double age = 0.;

    double mai = 0.;

    double hdz = 0.;

    double volume = 0.;

    const QMultiMap<FMUnit*, FMStand*> &stands = ForestManagementEngine::instance()->stands();

    QMultiMap<FMUnit*, FMStand*>::const_iterator it = stands.constFind(this);

    while (it != stands.constEnd() && it.key()==this) {

        FMStand *stand = it.value();

        stand->reload();

        stand->calculateMAI();

        // calculate sustainable total harvest (following Breymann)

        double area = stand->area();

        mai += stand->meanAnnualIncrementTotal() * area; // m3/yr

        age += stand->absoluteAge() * area;

        volume += stand->volume() * area;

        // HDZ: "haubarer" average increment: timber that is ready for final harvest

        if (stand->readyForFinalHarvest())

            hdz += stand->volume() / stand->absoluteAge() * area; //(0.1* stand->U()) * area; // note: changed!!!! was: volume/age * area

        total_area += area;

        ++it;

    }

    // reset

    ForestManagementEngine::instance()->scriptBridge()->treesObj()->setStand(0);

    mTotalArea = total_area;

    if (total_area==0.)

        return;

    mai /= total_area; // m3/ha*yr area weighted average of annual increment

    age /= total_area; // area weighted mean age

    hdz /= total_area; // =sum(Vol/age * share)

    mMAI = mai;

    //mMAI = mMAI * 1.15; // 15% increase, hack WR

    mHDZ = hdz;

    mMeanAge = age;

    mTotalVolume = volume;

    double rotation_length = U();

    double h_tot = mai * 2.*age / rotation_length;  //

    double h_reg = hdz * 2.*age / rotation_length;

    h_reg = h_tot * 0.85; // hack!

    h_reg *= agent()->schedulerOptions().harvestIntensity;

    h_tot *= agent()->schedulerOptions().harvestIntensity;

    double h_thi = qMax(h_tot - h_reg, 0.);

    qCDebug(abe) << "plan-update for unit" << id() << ": h-tot:" << h_tot << "h_reg:" << h_reg << "h_thi:" << h_thi << "of total volume:" << volume;

    double sf = mAgent->useSustainableHarvest();

    // we do not calculate sustainable harvest levels.

    // do a pure bottom up calculation

    double bottom_up_harvest = (plan_final / period_length) / total_area; // m3/ha*yr

    // the sustainable harvest yield is the current yield and some carry over from the last period

    double sustainable_harvest = h_reg;

//    if (mAnnualHarvestTarget>0.) {

//        double delta = realized/(total_area*period_length) - mAnnualHarvestTarget;

//        // if delta > 0: timber removal was too high -> plan less for the current period, and vice versa.

//        sustainable_harvest -= delta;

//    }

    mAnnualHarvestTarget = sustainable_harvest * sf + bottom_up_harvest * (1.-sf);

    mAnnualHarvestTarget = qMax(mAnnualHarvestTarget, 0.);

    mAnnualThinningTarget = (plan_thinning / period_length) / total_area; // m3/ha*yr

    if (scheduler())

        scheduler()->setHarvestTarget(mAnnualHarvestTarget, mAnnualThinningTarget);

}

QMutex _protect_agent_exec;

void FMUnit::runAgent()

{

    QMutexLocker m(&_protect_agent_exec); // avoid parallel execution of agent-code....

    // we need to set an execution context

    FMStand *stand = *ForestManagementEngine::instance()->stands().find(this);

    if (!stand)

        throw IException("Invalid stand in FMUnit::runAgent");

    FomeScript::setExecutionContext(stand, true); // true: add also agent as 'agent'

    QJSValue val;

    QJSValue agent_type = agent()->type()->jsObject();

    if (agent_type.property("run").isCallable()) {

        val = agent_type.property("run").callWithInstance(agent_type);

        qCDebug(abe) << "running agent-function 'run' for unit" <<  id() << ":" << val.toString();

    } else {

       qCDebug(abe) << "function 'run' is not a valid function of agent-type" << agent()->type()->name();

    }

}

void FMUnit::updatePlanOfCurrentYear()

{

    if (!scheduler())

        return;

    if (mTotalArea==0.)

        throw IException("FMUnit:updatePlan: unit area = 0???");

    // compare the harvests of the last year to the plan:

    double harvests = mRealizedHarvest - mRealizedHarvestLastYear;

    mRealizedHarvestLastYear = mRealizedHarvest;

    mAnnualHarvest = harvests;

    // difference in m3/ha

    double delta = harvests/mTotalArea - mAnnualHarvestTarget;

    mTotalPlanDeviation += delta;

    // apply decay function for deviation

    mTotalPlanDeviation *= mAgent->schedulerOptions().deviationDecayRate;

    // relative deviation: >0: too many harvests

    double rel_deviation = mAnnualHarvestTarget? mTotalPlanDeviation / mAnnualHarvestTarget : 0;

    // the current deviation is reduced to 50% in rebounce_yrs years.

    double rebounce_yrs = mAgent->schedulerOptions().scheduleRebounceDuration;

    double new_harvest = mAnnualHarvestTarget * (1. - rel_deviation/rebounce_yrs);

    // limit to minimum/maximum parameter

    new_harvest = qMax(new_harvest, mAgent->schedulerOptions().minScheduleHarvest);

    new_harvest = qMin(new_harvest, mAgent->schedulerOptions().maxScheduleHarvest);

    scheduler()->setHarvestTarget(new_harvest, mAnnualThinningTarget);

}

} // namesapce