treatments.hpp

Go to the documentation of this file.

/*
 * PoPS model - treatments
 *
 * Copyright (C) 2015-2020 by the authors.
 *
 * Authors: Anna Petrasova <akratoc gmail com>
 *          Vaclav Petras <wenzeslaus gmail com>
 *
 * The code contained herein is licensed under the GNU General Public
 * License. You may obtain a copy of the GNU General Public License
 * Version 2 or later at the following locations:
 *
 * http://www.opensource.org/licenses/gpl-license.html
 * http://www.gnu.org/copyleft/gpl.html
 */
 
#ifndef POPS_TREATMENTS_HPP
#define POPS_TREATMENTS_HPP
 
#include "raster.hpp"
#include "date.hpp"
#include "scheduling.hpp"
#include "utils.hpp"
 
#include <map>
#include <vector>
#include <string>
#include <functional>
#include <stdexcept>
 
namespace pops {
 
enum class TreatmentApplication
{
    Ratio,  
    AllInfectedInCell  
};
 
inline TreatmentApplication treatment_app_enum_from_string(const std::string& text)
{
    std::map<std::string, TreatmentApplication> mapping{
        {"ratio_to_all", TreatmentApplication::Ratio},
        {"ratio", TreatmentApplication::Ratio},
        {"all_infected_in_cell", TreatmentApplication::AllInfectedInCell},
        {"all infected", TreatmentApplication::AllInfectedInCell}};
    try {
        return mapping.at(text);
    }
    catch (const std::out_of_range&) {
        throw std::invalid_argument(
            "treatment_application_enum_from_string:"
            " Invalid value '"
            + text + "' provided");
    }
}
 
template<typename IntegerRaster, typename FloatRaster>
class AbstractTreatment
{
public:
    virtual unsigned get_start() = 0;
    virtual unsigned get_end() = 0;
    virtual bool should_start(unsigned step) = 0;
    virtual bool should_end(unsigned step) = 0;
    virtual void apply_treatment(
        IntegerRaster& infected,
        std::vector<IntegerRaster>& exposed,
        IntegerRaster& susceptible,
        IntegerRaster& resistant,
        IntegerRaster& total_hosts,
        const std::vector<std::vector<int>>& spatial_indeices) = 0;
    virtual void end_treatment(
        IntegerRaster& susceptible,
        IntegerRaster& resistant,
        const std::vector<std::vector<int>>& spatial_indeices) = 0;
    virtual void apply_treatment_mortality(
        IntegerRaster& infected,
        const std::vector<std::vector<int>>& spatial_indeices) = 0;
    virtual ~AbstractTreatment() {}
};
 
template<typename IntegerRaster, typename FloatRaster>
class BaseTreatment : public AbstractTreatment<IntegerRaster, FloatRaster>
{
protected:
    unsigned start_step_;
    unsigned end_step_;
    FloatRaster map_;
    TreatmentApplication application_;
 
public:
    BaseTreatment(
        const FloatRaster& map,
        unsigned start,
        TreatmentApplication treatment_application)
        : start_step_(start),
          end_step_(start),
          map_(map),
          application_(treatment_application)
    {}
    unsigned get_start() override
    {
        return start_step_;
    }
    unsigned get_end() override
    {
        return end_step_;
    }
    void apply_treatment_mortality(
        IntegerRaster& infected,
        const std::vector<std::vector<int>>& suitable_cells) override
    {
        for (auto indices : suitable_cells) {
            int i = indices[0];
            int j = indices[1];
            if (application_ == TreatmentApplication::Ratio) {
                infected(i, j) = infected(i, j) - (infected(i, j) * map_(i, j));
            }
            else if (application_ == TreatmentApplication::AllInfectedInCell) {
                infected(i, j) = map_(i, j) ? 0 : infected(i, j);
            }
        }
    }
};
 
template<typename IntegerRaster, typename FloatRaster>
class SimpleTreatment : public BaseTreatment<IntegerRaster, FloatRaster>
{
public:
    SimpleTreatment(
        const FloatRaster& map,
        unsigned start,
        TreatmentApplication treatment_application)
        : BaseTreatment<IntegerRaster, FloatRaster>(map, start, treatment_application)
    {}
    bool should_start(unsigned step) override
    {
        if (this->start_step_ == step)
            return true;
        return false;
    }
    bool should_end(unsigned) override
    {
        return false;
    }
    void apply_treatment(
        IntegerRaster& infected,
        std::vector<IntegerRaster>& exposed,
        IntegerRaster& susceptible,
        IntegerRaster& resistant,
        IntegerRaster& total_hosts,
        const std::vector<std::vector<int>>& suitable_cells) override
    {
        for (auto indices : suitable_cells) {
            int i = indices[0];
            int j = indices[1];
            int new_exposed_total = 0;
            int new_infected = 0;
            int new_susceptible = 0;
            int new_exposed_individual = 0;
            if (this->application_ == TreatmentApplication::Ratio) {
                new_infected = infected(i, j) - (infected(i, j) * this->map_(i, j));
                infected(i, j) = new_infected;
            }
            else if (this->application_ == TreatmentApplication::AllInfectedInCell) {
                new_infected = this->map_(i, j) ? 0 : infected(i, j);
                infected(i, j) = new_infected;
            }
            for (auto& raster : exposed) {
                if (this->application_ == TreatmentApplication::Ratio) {
                    new_exposed_individual =
                        raster(i, j) - (raster(i, j) * this->map_(i, j));
                    raster(i, j) = new_exposed_individual;
                    new_exposed_total += new_exposed_individual;
                }
                else if (
                    this->application_ == TreatmentApplication::AllInfectedInCell) {
                    new_exposed_individual =
                        raster(i, j) - (raster(i, j) * this->map_(i, j));
                    raster(i, j) = new_exposed_individual;
                    new_exposed_total += new_exposed_individual;
                }
            }
            new_susceptible =
                susceptible(i, j) - (susceptible(i, j) * this->map_(i, j));
            susceptible(i, j) = new_susceptible;
            total_hosts(i, j) =
                new_infected + new_susceptible + new_exposed_total + resistant(i, j);
        }
    }
    void end_treatment(
        IntegerRaster&, IntegerRaster&, const std::vector<std::vector<int>>&) override
    {
        return;
    }
};
 
template<typename IntegerRaster, typename FloatRaster>
class PesticideTreatment : public BaseTreatment<IntegerRaster, FloatRaster>
{
public:
    PesticideTreatment(
        const FloatRaster& map,
        unsigned start,
        unsigned end,
        TreatmentApplication treatment_application)
        : BaseTreatment<IntegerRaster, FloatRaster>(map, start, treatment_application)
    {
        this->end_step_ = end;
    }
    bool should_start(unsigned step) override
    {
        if (this->start_step_ == step)
            return true;
        return false;
    }
    bool should_end(unsigned step) override
    {
        if (this->end_step_ == step)
            return true;
        return false;
    }
 
    void apply_treatment(
        IntegerRaster& infected,
        std::vector<IntegerRaster>& exposed_vector,
        IntegerRaster& susceptible,
        IntegerRaster& resistant,
        IntegerRaster& total_hosts,
        const std::vector<std::vector<int>>& suitable_cells) override
    {
        UNUSED(total_hosts);
        for (auto indices : suitable_cells) {
            int i = indices[0];
            int j = indices[1];
            int infected_resistant = 0;
            int exposed_resistant_sum = 0;
            int susceptible_resistant = susceptible(i, j) * this->map_(i, j);
            int current_resistant = resistant(i, j);
            if (this->application_ == TreatmentApplication::Ratio) {
                infected_resistant = infected(i, j) * this->map_(i, j);
            }
            else if (this->application_ == TreatmentApplication::AllInfectedInCell) {
                infected_resistant = this->map_(i, j) ? infected(i, j) : 0;
            }
            infected(i, j) -= infected_resistant;
            for (auto& exposed : exposed_vector) {
                int exposed_resistant = 0;
                if (this->application_ == TreatmentApplication::Ratio) {
                    exposed_resistant = exposed(i, j) * this->map_(i, j);
                }
                else if (
                    this->application_ == TreatmentApplication::AllInfectedInCell) {
                    exposed_resistant = this->map_(i, j) ? exposed(i, j) : 0;
                }
                exposed(i, j) -= exposed_resistant;
                exposed_resistant_sum += exposed_resistant;
            }
            resistant(i, j) = infected_resistant + exposed_resistant_sum
                              + susceptible_resistant + current_resistant;
            susceptible(i, j) -= susceptible_resistant;
        }
    }
    void end_treatment(
        IntegerRaster& susceptible,
        IntegerRaster& resistant,
        const std::vector<std::vector<int>>& suitable_cells) override
    {
        for (auto indices : suitable_cells) {
            int i = indices[0];
            int j = indices[1];
            if (this->map_(i, j) > 0) {
                susceptible(i, j) += resistant(i, j);
                resistant(i, j) = 0;
            }
        }
    }
};
 
template<typename IntegerRaster, typename FloatRaster>
class Treatments
{
private:
    std::vector<AbstractTreatment<IntegerRaster, FloatRaster>*> treatments;
    Scheduler scheduler_;
 
public:
    Treatments(const Scheduler& scheduler) : scheduler_(scheduler) {}
    ~Treatments()
    {
        for (auto item : treatments) {
            delete item;
        }
    }
    void add_treatment(
        const FloatRaster& map,
        const Date& start_date,
        int num_days,
        TreatmentApplication treatment_application)
    {
        unsigned start = scheduler_.schedule_action_date(start_date);
        if (num_days == 0)
            treatments.push_back(new SimpleTreatment<IntegerRaster, FloatRaster>(
                map, start, treatment_application));
        else {
            Date end_date(start_date);
            end_date.add_days(num_days);
            unsigned end = scheduler_.schedule_action_date(end_date);
            treatments.push_back(new PesticideTreatment<IntegerRaster, FloatRaster>(
                map, start, end, treatment_application));
        }
    }
    bool manage(
        unsigned current,
        IntegerRaster& infected,
        std::vector<IntegerRaster>& exposed,
        IntegerRaster& susceptible,
        IntegerRaster& resistant,
        IntegerRaster& total_hosts,
        const std::vector<std::vector<int>>& suitable_cells)
    {
        bool changed = false;
        for (unsigned i = 0; i < treatments.size(); i++) {
            if (treatments[i]->should_start(current)) {
                treatments[i]->apply_treatment(
                    infected,
                    exposed,
                    susceptible,
                    resistant,
                    total_hosts,
                    suitable_cells);
                changed = true;
            }
            else if (treatments[i]->should_end(current)) {
                treatments[i]->end_treatment(susceptible, resistant, suitable_cells);
                changed = true;
            }
        }
        return changed;
    }
    bool manage_mortality(
        unsigned current,
        IntegerRaster& infected,
        const std::vector<std::vector<int>>& suitable_cells)
    {
        bool applied = false;
        for (unsigned i = 0; i < treatments.size(); i++)
            if (treatments[i]->should_start(current)) {
                treatments[i]->apply_treatment_mortality(infected, suitable_cells);
                applied = true;
            }
        return applied;
    }
    void clear_after_step(unsigned step)
    {
        for (auto& treatment : treatments) {
            if (treatment->get_start() > step) {
                delete treatment;
                treatment = nullptr;
            }
        }
        treatments.erase(
            std::remove(treatments.begin(), treatments.end(), nullptr),
            treatments.end());
    }
};
 
}  // namespace pops
#endif  // POPS_TREATMENTS_HPP