PlantArchitecture.h

Go to the documentation of this file.

 
#ifndef PLANT_ARCHITECTURE
#define PLANT_ARCHITECTURE
 
#include <functional>
#include <utility>
#include "Context.h"
#include "Hungarian.h"
 
// Constants
constexpr float C_molecular_wt = 12.01; // g C mol^-1
 
// forward declarations of classes/structs
class PlantArchitecture;
class CollisionDetection;
struct Shoot;
struct Phytomer;
 
struct RandomParameter_float {
public:
 
    explicit RandomParameter_float() {
        constval = 0.f;
        distribution = "constant";
        generator = nullptr;
        sampled = false;
    }
 
 
    explicit RandomParameter_float(float val) {
        constval = val;
        distribution = "constant";
        generator = nullptr;
        sampled = false;
    }
 
 
    explicit RandomParameter_float(std::minstd_rand0 *rand_generator) {
        constval = 0.f;
        distribution = "constant";
        generator = rand_generator;
        sampled = false;
    }
 
    void initialize(float a_val, std::minstd_rand0 *rand_generator) {
        constval = a_val;
        distribution = "constant";
        generator = rand_generator;
        sampled = false;
    }
 
    void initialize(std::minstd_rand0 *rand_generator) {
        constval = 1.f;
        distribution = "constant";
        generator = rand_generator;
        sampled = false;
    }
 
    RandomParameter_float &operator=(float a) {
        this->distribution = "constant";
        this->constval = a;
        this->sampled = false;
        return *this;
    }
 
    void uniformDistribution(float minval, float maxval) {
        if (minval > maxval) {
            throw(std::runtime_error("ERROR (PlantArchitecture): RandomParameter_float::uniformDistribution() - minval must be less than or equal to maxval."));
        }
        distribution = "uniform";
        distribution_parameters = {minval, maxval};
        sampled = false;
    }
 
    void normalDistribution(float mean, float std_dev) {
        distribution = "normal";
        distribution_parameters = {mean, std_dev};
        sampled = false;
    }
 
    void weibullDistribution(float shape, float scale) {
        distribution = "weibull";
        distribution_parameters = {shape, scale};
        sampled = false;
    }
 
    float val() {
        if (!sampled) {
            constval = resample();
        }
        return constval;
    }
 
    float resample() {
        sampled = true;
        if (distribution != "constant") {
            if (generator == nullptr) {
                throw(std::runtime_error("ERROR (PlantArchitecture): Random parameter was not properly initialized with random number generator."));
            }
            if (distribution == "uniform") {
                std::uniform_real_distribution<float> unif_distribution;
                constval = distribution_parameters.at(0) + unif_distribution(*generator) * (distribution_parameters.at(1) - distribution_parameters.at(0));
            } else if (distribution == "normal") {
                std::normal_distribution<float> norm_distribution(distribution_parameters.at(0), distribution_parameters.at(1));
                constval = norm_distribution(*generator);
            } else if (distribution == "weibull") {
                std::weibull_distribution<float> wbull_distribution(distribution_parameters.at(0), distribution_parameters.at(1));
                constval = wbull_distribution(*generator);
            }
        }
        return constval;
    }
 
    std::string distribution;
    std::vector<float> distribution_parameters;
 
private:
    bool sampled;
    float constval;
 
 
    std::minstd_rand0 *generator;
};
 
struct RandomParameter_int {
public:
    explicit RandomParameter_int() {
        constval = 1;
        distribution = "constant";
        generator = nullptr;
        sampled = false;
    }
 
    void initialize(int a_val, std::minstd_rand0 *rand_generator) {
        constval = a_val;
        distribution = "constant";
        generator = rand_generator;
        sampled = false;
    }
 
    void initialize(std::minstd_rand0 *rand_generator) {
        constval = 1;
        distribution = "constant";
        generator = rand_generator;
        sampled = false;
    }
 
    RandomParameter_int &operator=(int a) {
        this->distribution = "constant";
        this->constval = a;
        this->sampled = false;
        return *this;
    }
 
    void uniformDistribution(int minval, int maxval) {
        if (minval > maxval) {
            throw(std::runtime_error("ERROR (PlantArchitecture): RandomParameter_int::uniformDistribution() - minval must be less than or equal to maxval."));
        }
        distribution = "uniform";
        distribution_parameters = {minval, maxval};
        sampled = false;
    }
 
    void discreteValues(const std::vector<int> &values) {
        distribution = "discretevalues";
        distribution_parameters = values;
        sampled = false;
    }
 
    int val() {
        if (!sampled) {
            constval = resample();
        }
        return constval;
    }
 
    int resample() {
        sampled = true;
        if (distribution != "constant") {
            if (generator == nullptr) {
                throw(std::runtime_error("ERROR (PlantArchitecture): Random parameter was not properly initialized with random number generator."));
            }
            if (distribution == "uniform") {
                std::uniform_int_distribution<> unif_distribution(distribution_parameters.at(0), distribution_parameters.at(1));
                constval = unif_distribution(*generator);
            } else if (distribution == "discretevalues") {
                std::uniform_int_distribution<> unif_distribution(0, distribution_parameters.size() - 1);
                constval = distribution_parameters.at(unif_distribution(*generator));
            }
        }
        return constval;
    }
 
    std::string distribution;
    std::vector<int> distribution_parameters;
 
private:
    bool sampled;
    int constval;
    std::minstd_rand0 *generator;
};
 
struct AxisRotation {
public:
    AxisRotation() {
        pitch = 0;
        yaw = 0;
        roll = 0;
        azimuth = 0;
        peduncle_axis = helios::make_vec3(0, 0, 1);
    }
 
    AxisRotation(float a_pitch, float a_yaw, float a_roll) {
        pitch = a_pitch;
        yaw = a_yaw;
        roll = a_roll;
        azimuth = 0;
        peduncle_axis = helios::make_vec3(0, 0, 1);
    }
 
    float pitch;
    float yaw;
    float roll;
    float azimuth; // azimuth rotation for inflorescences (aligns to peduncle orientation)
    helios::vec3 peduncle_axis; // direction vector of peduncle at attachment point (for compound yaw rotation)
 
    AxisRotation operator+(const AxisRotation &a) const;
    AxisRotation operator-(const AxisRotation &a) const;
 
    friend std::ostream &operator<<(std::ostream &os, const AxisRotation &rot) {
        return os << "AxisRotation<" << rot.pitch << ", " << rot.yaw << ", " << rot.roll << ">";
    }
};
 
inline AxisRotation make_AxisRotation(float a_pitch, float a_yaw, float a_roll) {
    return {a_pitch, a_yaw, a_roll};
}
 
inline AxisRotation AxisRotation::operator+(const AxisRotation &a) const {
    return {a.pitch + pitch, a.yaw + yaw, a.roll + roll};
}
 
inline AxisRotation AxisRotation::operator-(const AxisRotation &a) const {
    return {a.pitch - pitch, a.yaw - yaw, a.roll - roll};
}
 
enum BudState { BUD_DORMANT = 0, BUD_ACTIVE = 1, BUD_FLOWER_CLOSED = 2, BUD_FLOWER_OPEN = 3, BUD_FRUITING = 4, BUD_DEAD = 5 };
 
struct CarbohydrateParameters {
 
    // -- Stem Growth Parameters -- //
    float stem_density = 54000;
    float stem_carbon_percentage = 0.4559;
    float carbohydrate_percentage = 0.10;
    float carbohydrate_carbon_percentage = 0.4;
    float stem_structural_carbon_percentage = stem_carbon_percentage - (carbohydrate_carbon_percentage * carbohydrate_percentage);
    float maturity_age = 180;
    float initial_density_ratio = 0.2;
    float shoot_root_ratio = 3;
 
    // -- Leaf Growth Parameters -- //
    float SLA = 2.5e-2;
    float leaf_carbon_percentage = 0.4444;
 
    // -- Flower Growth Parameters -- //
    float total_flower_cost = 8.33e-4;
 
    // -- Fruit Growth Parameters -- //
    float fruit_density = 525000;
    float fruit_carbon_percentage = 0.4786;
 
    // -- Respiration Parameters -- //
    float stem_maintenance_respiration_rate = 3.5024e-05;
    float root_maintenance_respiration_rate = 3.5024e-05;
    float living_wood_fraction = 0.5;
    float dormant_respiration_fraction = 0.9;
    float growth_respiration_fraction = 0.28;
 
    // -- Organ Abortion Thresholds -- //
    float carbohydrate_abortion_threshold = 0.15;
    float carbohydrate_pruning_threshold = 0.025;
    float bud_death_threshold_days = 2;
    float branch_death_threshold_days = 5;
 
    // -- Phyllochron Adjustment Parameters -- //
    float carbohydrate_phyllochron_threshold = 0.05;
    float carbohydrate_vegetative_break_threshold = 0.15;
 
    float carbohydrate_growth_threshold = 0.15;
 
    // -- Carbon Transfer Parameters -- //
    float carbohydrate_transfer_threshold_down = 0.025;
    float carbohydrate_transfer_threshold_up = 0.04;
    float carbon_conductance_down = 0.75; //<= 1.0
    float carbon_conductance_up = carbon_conductance_down; // Conductance of carbon from parent to child shoots << conductance from child to parent
};
 
struct NitrogenParameters {
 
    // -- Leaf Nitrogen Content (Area Basis) -- //
    float target_leaf_N_area = 1.5f;
    float minimum_leaf_N_area = 0.5f;
 
    // -- Allocation Parameters -- //
    float root_allocation_fraction = 0.15f;
 
    // -- Rate Limiting Parameters -- //
    float max_N_accumulation_rate = 0.1f;
 
    // -- Remobilization Parameters -- //
    float leaf_remobilization_efficiency = 0.70f;
    float remobilization_age_threshold = 0.70f;
 
    // -- Fruit Nitrogen Parameters -- //
    float fruit_N_area = 1.0f;
};
 
 
std::vector<uint> makeTubeFromCones(uint radial_subdivisions, const std::vector<helios::vec3> &vertices, const std::vector<float> &radii, const std::vector<helios::RGBcolor> &colors, helios::Context *context_ptr);
 
struct VegetativeBud {
 
    // state of the bud
    BudState state = BUD_DORMANT;
    // label of the shoot type that will be produced if the bud breaks into a shoot
    std::string shoot_type_label;
    // ID of the shoot that the bud will produce if it breaks into a shoot
    uint shoot_ID = -1;
};
 
struct FloralBud {
 
    // state of the bud
    BudState state = BUD_DORMANT;
    // amount of time since the bud flowered (=0 if it has not yet flowered)
    float time_counter = 0;
    //=0 for axillary buds, =1 for terminal buds
    bool isterminal = false;
    // For axillary buds: index of the petiole within the internode that this floral bud originates from
    // For terminal buds: index of the phytomer within the shoot that this floral bud originates from
    uint parent_index = 0;
    // Index of the bud within the petiole that this floral bud originates from
    uint bud_index = 0;
    // Scaling factor fraction of the fruit (if present), ranging from 0 to 1
    float current_fruit_scale_factor = 1;
    float previous_fruit_scale_factor = 0;
 
 
    helios::vec3 base_position;
    AxisRotation base_rotation;
    helios::vec3 bending_axis;
 
    std::vector<helios::vec3> inflorescence_bases;
    std::vector<AxisRotation> inflorescence_rotation; // pitch, yaw, roll for each flower/fruit
    std::vector<float> inflorescence_base_scales; // individual base scale for each flower/fruit (before growth scaling)
    std::vector<uint> peduncle_objIDs;
    std::vector<uint> inflorescence_objIDs;
};
 
struct LeafPrototype {
public:
    explicit LeafPrototype(std::minstd_rand0 *generator);
 
    LeafPrototype() = default;
 
    uint (*prototype_function)(helios::Context *, LeafPrototype *prototype_parameters, int compound_leaf_index) = nullptr;
 
 
    std::string OBJ_model_file;
 
 
    std::map<int, std::string> leaf_texture_file;
 
    // Ratio of leaf width to leaf length
    RandomParameter_float leaf_aspect_ratio;
 
    RandomParameter_float midrib_fold_fraction;
 
    // Parameters for leaf curvature
    RandomParameter_float longitudinal_curvature;
    RandomParameter_float lateral_curvature;
 
    RandomParameter_float petiole_roll;
 
    // Parameters for leaf wave/wrinkles
    RandomParameter_float wave_period;
    RandomParameter_float wave_amplitude;
 
    // Parameters for leaf buckling
    RandomParameter_float leaf_buckle_length;
    RandomParameter_float leaf_buckle_angle;
 
    helios::vec3 leaf_offset;
 
    uint subdivisions = 1;
 
    uint unique_prototypes = 1;
 
    bool build_petiolule = false;
 
    uint unique_prototype_identifier = 0;
 
    void duplicate(const LeafPrototype &a) {
        this->leaf_texture_file = a.leaf_texture_file;
        this->OBJ_model_file = a.OBJ_model_file;
        this->leaf_aspect_ratio = a.leaf_aspect_ratio;
        this->midrib_fold_fraction = a.midrib_fold_fraction;
        this->longitudinal_curvature = a.longitudinal_curvature;
        this->lateral_curvature = a.lateral_curvature;
        this->petiole_roll = a.petiole_roll;
        this->wave_period = a.wave_period;
        this->wave_amplitude = a.wave_amplitude;
        this->leaf_buckle_length = a.leaf_buckle_length;
        this->leaf_buckle_angle = a.leaf_buckle_angle;
        this->leaf_offset = a.leaf_offset;
        this->subdivisions = a.subdivisions;
        this->unique_prototypes = a.unique_prototypes;
        this->unique_prototype_identifier = a.unique_prototype_identifier;
        this->build_petiolule = a.build_petiolule;
        this->prototype_function = a.prototype_function;
        this->generator = a.generator;
    }
 
    LeafPrototype &operator=(const LeafPrototype &a) {
        if (this != &a) {
            this->leaf_texture_file = a.leaf_texture_file;
            this->OBJ_model_file = a.OBJ_model_file;
            this->leaf_aspect_ratio = a.leaf_aspect_ratio;
            if (a.leaf_aspect_ratio.distribution != "constant")
                this->leaf_aspect_ratio.resample();
            this->midrib_fold_fraction = a.midrib_fold_fraction;
            if (a.midrib_fold_fraction.distribution != "constant")
                this->midrib_fold_fraction.resample();
            this->longitudinal_curvature = a.longitudinal_curvature;
            if (a.longitudinal_curvature.distribution != "constant")
                this->longitudinal_curvature.resample();
            this->lateral_curvature = a.lateral_curvature;
            if (a.lateral_curvature.distribution != "constant")
                this->lateral_curvature.resample();
            this->petiole_roll = a.petiole_roll;
            if (a.petiole_roll.distribution != "constant")
                this->petiole_roll.resample();
            this->wave_period = a.wave_period;
            if (a.wave_period.distribution != "constant")
                this->wave_period.resample();
            this->wave_amplitude = a.wave_amplitude;
            if (a.wave_amplitude.distribution != "constant")
                this->wave_amplitude.resample();
            this->leaf_buckle_length = a.leaf_buckle_length;
            if (a.leaf_buckle_length.distribution != "constant")
                this->leaf_buckle_length.resample();
            this->leaf_buckle_angle = a.leaf_buckle_angle;
            if (a.leaf_buckle_angle.distribution != "constant")
                this->leaf_buckle_angle.resample();
            this->leaf_offset = a.leaf_offset;
            this->subdivisions = a.subdivisions;
            this->unique_prototypes = a.unique_prototypes;
            this->unique_prototype_identifier = a.unique_prototype_identifier;
            this->build_petiolule = a.build_petiolule;
            this->prototype_function = a.prototype_function;
            this->generator = a.generator;
            if (this->generator != nullptr) {
                this->sampleIdentifier();
            }
        }
        return *this;
    }
 
    void sampleIdentifier() {
        assert(generator != nullptr);
        std::uniform_int_distribution<uint> unif_distribution;
        this->unique_prototype_identifier = unif_distribution(*generator);
    }
 
private:
    std::minstd_rand0 *generator{};
};
 
struct PhytomerParameters {
private:
    struct InternodeParameters {
 
        RandomParameter_float pitch;
        RandomParameter_float phyllotactic_angle;
        RandomParameter_float radius_initial;
        RandomParameter_int max_vegetative_buds_per_petiole;
        RandomParameter_int max_floral_buds_per_petiole;
        helios::RGBcolor color;
        std::string image_texture;
        uint length_segments;
        uint radial_subdivisions;
 
        InternodeParameters &operator=(const InternodeParameters &a) {
            if (this != &a) {
                this->pitch = a.pitch;
                if (a.pitch.distribution != "constant")
                    this->pitch.resample();
                this->phyllotactic_angle = a.phyllotactic_angle;
                if (a.phyllotactic_angle.distribution != "constant")
                    this->phyllotactic_angle.resample();
                this->radius_initial = a.radius_initial;
                if (a.radius_initial.distribution != "constant")
                    this->radius_initial.resample();
                this->max_vegetative_buds_per_petiole = a.max_vegetative_buds_per_petiole;
                if (a.max_vegetative_buds_per_petiole.distribution != "constant")
                    this->max_vegetative_buds_per_petiole.resample();
                this->max_floral_buds_per_petiole = a.max_floral_buds_per_petiole;
                if (a.max_floral_buds_per_petiole.distribution != "constant")
                    this->max_floral_buds_per_petiole.resample();
                this->color = a.color;
                this->image_texture = a.image_texture;
                this->length_segments = a.length_segments;
                this->radial_subdivisions = a.radial_subdivisions;
            }
            return *this;
        }
    };
 
    struct PetioleParameters {
 
        uint petioles_per_internode;
        RandomParameter_float pitch;
        RandomParameter_float radius;
        RandomParameter_float length;
        RandomParameter_float curvature;
        RandomParameter_float taper;
        helios::RGBcolor color;
        uint length_segments;
        uint radial_subdivisions;
 
        PetioleParameters &operator=(const PetioleParameters &a) {
            if (this != &a) {
                this->petioles_per_internode = a.petioles_per_internode;
                this->pitch = a.pitch;
                if (a.pitch.distribution != "constant")
                    this->pitch.resample();
                this->radius = a.radius;
                if (a.radius.distribution != "constant")
                    this->radius.resample();
                this->length = a.length;
                if (a.length.distribution != "constant")
                    this->length.resample();
                this->curvature = a.curvature;
                if (a.curvature.distribution != "constant")
                    this->curvature.resample();
                this->taper = a.taper;
                if (a.taper.distribution != "constant")
                    this->taper.resample();
                this->color = a.color;
                this->length_segments = a.length_segments;
                this->radial_subdivisions = a.radial_subdivisions;
            }
            return *this;
        }
    };
 
    struct LeafParameters {
        RandomParameter_int leaves_per_petiole;
        RandomParameter_float pitch;
        RandomParameter_float yaw;
        RandomParameter_float roll;
        RandomParameter_float leaflet_offset;
        RandomParameter_float leaflet_scale;
        RandomParameter_float prototype_scale;
        LeafPrototype prototype;
 
        LeafParameters &operator=(const LeafParameters &a) {
            if (this != &a) {
                this->leaves_per_petiole = a.leaves_per_petiole;
                if (a.leaves_per_petiole.distribution != "constant")
                    this->leaves_per_petiole.resample();
                this->pitch = a.pitch;
                if (a.pitch.distribution != "constant")
                    this->pitch.resample();
                this->yaw = a.yaw;
                if (a.yaw.distribution != "constant")
                    this->yaw.resample();
                this->roll = a.roll;
                if (a.roll.distribution != "constant")
                    this->roll.resample();
                this->leaflet_offset = a.leaflet_offset;
                if (a.leaflet_offset.distribution != "constant")
                    this->leaflet_offset.resample();
                this->leaflet_scale = a.leaflet_scale;
                if (a.leaflet_scale.distribution != "constant")
                    this->leaflet_scale.resample();
                this->prototype_scale = a.prototype_scale;
                if (a.prototype_scale.distribution != "constant")
                    this->prototype_scale.resample();
                this->prototype.duplicate(a.prototype);
            }
            return *this;
        }
    };
 
    struct PeduncleParameters {
        RandomParameter_float length;
        RandomParameter_float radius;
        RandomParameter_float pitch;
        RandomParameter_float roll;
        RandomParameter_float curvature;
        helios::RGBcolor color;
        uint length_segments;
        uint radial_subdivisions;
 
        PeduncleParameters &operator=(const PeduncleParameters &a) {
            if (this != &a) {
                this->length = a.length;
                if (a.length.distribution != "constant")
                    this->length.resample();
                this->radius = a.radius;
                if (a.radius.distribution != "constant")
                    this->radius.resample();
                this->pitch = a.pitch;
                if (a.pitch.distribution != "constant")
                    this->pitch.resample();
                this->roll = a.roll;
                if (a.roll.distribution != "constant")
                    this->roll.resample();
                this->curvature = a.curvature;
                if (a.curvature.distribution != "constant")
                    this->curvature.resample();
                this->color = a.color;
                this->length_segments = a.length_segments;
                this->radial_subdivisions = a.radial_subdivisions;
            }
            return *this;
        }
    };
 
    struct InflorescenceParameters {
 
        RandomParameter_int flowers_per_peduncle;
        RandomParameter_float flower_offset;
        RandomParameter_float pitch;
        RandomParameter_float roll;
        RandomParameter_float flower_prototype_scale;
        uint (*flower_prototype_function)(helios::Context *, uint subdivisions, bool flower_is_open) = nullptr;
        RandomParameter_float fruit_prototype_scale;
        uint (*fruit_prototype_function)(helios::Context *, uint subdivisions) = nullptr;
        RandomParameter_float fruit_gravity_factor_fraction;
        uint unique_prototypes;
 
        InflorescenceParameters &operator=(const InflorescenceParameters &a) {
            if (this != &a) {
                this->flowers_per_peduncle = a.flowers_per_peduncle;
                this->flowers_per_peduncle.resample();
                this->flower_offset = a.flower_offset;
                this->flower_offset.resample();
                this->pitch = a.pitch;
                this->pitch.resample();
                this->roll = a.roll;
                this->roll.resample();
                this->flower_prototype_scale = a.flower_prototype_scale;
                this->flower_prototype_scale.resample();
                this->flower_prototype_function = a.flower_prototype_function;
                this->fruit_prototype_scale = a.fruit_prototype_scale;
                this->fruit_prototype_scale.resample();
                this->fruit_prototype_function = a.fruit_prototype_function;
                this->fruit_gravity_factor_fraction = a.fruit_gravity_factor_fraction;
                this->fruit_gravity_factor_fraction.resample();
                this->unique_prototypes = a.unique_prototypes;
            }
            return *this;
        }
    };
 
public:
    InternodeParameters internode;
 
    PetioleParameters petiole;
 
    LeafParameters leaf;
 
    PeduncleParameters peduncle;
 
    InflorescenceParameters inflorescence;
 
    // Custom user-defined function that is called when a phytomer is created
    void (*phytomer_creation_function)(std::shared_ptr<Phytomer> phytomer_ptr, uint shoot_node_index, uint parent_shoot_node_index, uint shoot_max_nodes, float plant_age) = nullptr;
 
    // Custom user-defined function that is called for each phytomer on every time step
    void (*phytomer_callback_function)(std::shared_ptr<Phytomer> phytomer_ptr) = nullptr;
 
 
    PhytomerParameters();
 
    explicit PhytomerParameters(std::minstd_rand0 *generator);
 
    friend class PlantArchitecture;
    friend struct Phytomer;
    friend struct Shoot;
};
 
struct ShootParameters {
 
    ShootParameters();
 
    explicit ShootParameters(std::minstd_rand0 *generator);
 
    PhytomerParameters phytomer_parameters;
 
    // ---- Geometric Parameters ---- //
 
    RandomParameter_int max_nodes;
    RandomParameter_int max_nodes_per_season;
    RandomParameter_float girth_area_factor;
    RandomParameter_float insertion_angle_tip;
    RandomParameter_float insertion_angle_decay_rate;
    RandomParameter_float internode_length_max;
    RandomParameter_float internode_length_min;
    RandomParameter_float internode_length_decay_rate;
    RandomParameter_float base_roll;
    RandomParameter_float base_yaw;
    RandomParameter_float gravitropic_curvature;
    RandomParameter_float tortuosity;
 
    // --- Growth Parameters --- //
 
    RandomParameter_float phyllochron_min;
    RandomParameter_float elongation_rate_max;
    RandomParameter_float vegetative_bud_break_probability_min;
    RandomParameter_float vegetative_bud_break_probability_max;
 
    RandomParameter_float vegetative_bud_break_probability_decay_rate;
    RandomParameter_int max_terminal_floral_buds;
    RandomParameter_float flower_bud_break_probability;
    RandomParameter_float fruit_set_probability;
    RandomParameter_float vegetative_bud_break_time;
    bool flowers_require_dormancy;
    bool growth_requires_dormancy;
    bool determinate_shoot_growth;
 
    // ---- Custom Functions ---- //
 
    void defineChildShootTypes(const std::vector<std::string> &child_shoot_type_labels, const std::vector<float> &child_shoot_type_probabilities);
 
    ShootParameters &operator=(const ShootParameters &a) {
        this->phytomer_parameters = a.phytomer_parameters;
        this->max_nodes = a.max_nodes;
        this->max_nodes.resample();
        this->max_nodes_per_season = a.max_nodes_per_season;
        this->max_nodes_per_season.resample();
        this->phyllochron_min = a.phyllochron_min;
        this->phyllochron_min.resample();
        this->girth_area_factor = a.girth_area_factor;
        this->girth_area_factor.resample();
        this->vegetative_bud_break_probability_min = a.vegetative_bud_break_probability_min;
        this->vegetative_bud_break_probability_min.resample();
        this->vegetative_bud_break_probability_max = a.vegetative_bud_break_probability_max;
        this->vegetative_bud_break_probability_max.resample();
        this->flower_bud_break_probability = a.flower_bud_break_probability;
        this->flower_bud_break_probability.resample();
        this->fruit_set_probability = a.fruit_set_probability;
        this->fruit_set_probability.resample();
        this->gravitropic_curvature = a.gravitropic_curvature;
        this->gravitropic_curvature.resample();
        this->tortuosity = a.tortuosity;
        this->tortuosity.resample();
        this->vegetative_bud_break_probability_min = a.vegetative_bud_break_probability_min;
        this->vegetative_bud_break_probability_min.resample();
        this->vegetative_bud_break_probability_decay_rate = a.vegetative_bud_break_probability_decay_rate;
        this->vegetative_bud_break_probability_decay_rate.resample();
        this->max_terminal_floral_buds = a.max_terminal_floral_buds;
        this->max_terminal_floral_buds.resample();
        this->flower_bud_break_probability = a.flower_bud_break_probability;
        this->flower_bud_break_probability.resample();
        this->fruit_set_probability = a.fruit_set_probability;
        this->fruit_set_probability.resample();
        this->vegetative_bud_break_time = a.vegetative_bud_break_time;
        this->vegetative_bud_break_time.resample();
        this->insertion_angle_tip = a.insertion_angle_tip;
        this->insertion_angle_tip.resample();
        this->insertion_angle_decay_rate = a.insertion_angle_decay_rate;
        this->insertion_angle_decay_rate.resample();
        this->internode_length_max = a.internode_length_max;
        this->internode_length_max.resample();
        this->internode_length_min = a.internode_length_min;
        this->internode_length_min.resample();
        this->internode_length_decay_rate = a.internode_length_decay_rate;
        this->internode_length_decay_rate.resample();
        this->base_roll = a.base_roll;
        this->base_roll.resample();
        this->base_yaw = a.base_yaw;
        this->base_yaw.resample();
        this->flowers_require_dormancy = a.flowers_require_dormancy;
        this->growth_requires_dormancy = a.growth_requires_dormancy;
        this->child_shoot_type_labels = a.child_shoot_type_labels;
        this->child_shoot_type_probabilities = a.child_shoot_type_probabilities;
        this->determinate_shoot_growth = a.determinate_shoot_growth;
        this->child_shoot_type_labels = a.child_shoot_type_labels;
        this->child_shoot_type_probabilities = a.child_shoot_type_probabilities;
        return *this;
    }
 
    friend class PlantArchitecture;
    friend struct Shoot;
 
protected:
    std::vector<std::string> child_shoot_type_labels;
    std::vector<float> child_shoot_type_probabilities;
};
 
struct Phytomer {
public:
    Phytomer(const PhytomerParameters &params, Shoot *parent_shoot, uint phytomer_index, const helios::vec3 &parent_internode_axis, const helios::vec3 &parent_petiole_axis, helios::vec3 internode_base_origin, const AxisRotation &shoot_base_rotation,
             float internode_radius, float internode_length_max, float internode_length_scale_factor_fraction, float leaf_scale_factor_fraction, uint rank, PlantArchitecture *plantarchitecture_ptr, helios::Context *context_ptr);
 
    // ---- query info about the phytomer ---- //
 
    [[nodiscard]] std::vector<helios::vec3> getInternodeNodePositions() const;
 
    [[nodiscard]] std::vector<float> getInternodeNodeRadii() const;
 
    [[nodiscard]] float calculatePhytomerVolume(uint node_number) const;
 
    [[nodiscard]] helios::vec3 getInternodeAxisVector(float stem_fraction) const;
 
    [[nodiscard]] helios::vec3 getPetioleAxisVector(float stem_fraction, uint petiole_index) const;
 
    [[nodiscard]] helios::vec3 getPeduncleAxisVector(float stem_fraction, uint petiole_index, uint bud_index) const;
 
    [[nodiscard]] static helios::vec3 getAxisVector(float stem_fraction, const std::vector<helios::vec3> &axis_vertices);
 
    [[nodiscard]] float getInternodeLength() const;
 
    [[nodiscard]] float getInternodeRadius() const;
 
    [[nodiscard]] float getPetioleLength() const;
 
    [[nodiscard]] float getInternodeRadius(float stem_fraction) const;
 
    [[nodiscard]] float getLeafArea() const;
 
    [[nodiscard]] helios::vec3 getLeafBasePosition(uint petiole_index, uint leaf_index) const;
 
    [[nodiscard]] bool hasLeaf() const;
 
    [[nodiscard]] float calculateDownstreamLeafArea() const;
 
    // ---- modify the phytomer ---- //
 
 
    void setInternodeLengthScaleFraction(float internode_scale_factor_fraction, bool update_context_geometry);
 
 
    void scaleInternodeMaxLength(float scale_factor);
 
 
    void setInternodeMaxLength(float internode_length_max_new);
 
 
    void setInternodeMaxRadius(float internode_radius_max_new);
 
 
    void setLeafScaleFraction(uint petiole_index, float leaf_scale_factor_fraction);
 
 
    void setLeafScaleFraction(float leaf_scale_factor_fraction);
 
 
    void setLeafPrototypeScale(uint petiole_index, float leaf_prototype_scale);
 
 
    void setLeafPrototypeScale(float leaf_prototype_scale);
 
 
    void scaleLeafPrototypeScale(uint petiole_index, float scale_factor);
 
 
    void scaleLeafPrototypeScale(float scale_factor);
 
 
    void scalePetioleGeometry(uint petiole_index, float target_length, float target_base_radius);
 
    void setInflorescenceScaleFraction(FloralBud &fbud, float inflorescence_scale_factor_fraction) const;
 
    void setPetioleBase(const helios::vec3 &base_position);
 
    void rotateLeaf(uint petiole_index, uint leaf_index, const AxisRotation &rotation);
 
    void setVegetativeBudState(BudState state);
 
    void setVegetativeBudState(BudState state, uint petiole_index, uint bud_index);
 
    static void setVegetativeBudState(BudState state, VegetativeBud &vbud);
 
    void setFloralBudState(BudState state);
 
    void setFloralBudState(BudState state, uint petiole_index, uint bud_index);
 
    void setFloralBudState(BudState state, FloralBud &fbud);
 
    void removeLeaf();
 
    void deletePhytomer();
 
private:
    helios::vec3 calculateCollisionAvoidanceDirection(const helios::vec3 &internode_base_origin, const helios::vec3 &internode_axis, bool &collision_detection_active) const;
 
    bool applySolidObstacleAvoidance(const helios::vec3 &current_position, helios::vec3 &internode_axis) const;
 
    helios::vec3 calculateAttractionPointDirection(const helios::vec3 &internode_base_origin, const helios::vec3 &internode_axis, bool &attraction_active) const;
 
    helios::vec3 calculatePetioleCollisionAvoidanceDirection(const helios::vec3 &petiole_base_origin, const helios::vec3 &proposed_petiole_axis, bool &collision_detection_active) const;
 
 
    helios::vec3 calculateFruitCollisionAvoidanceDirection(const helios::vec3 &fruit_base_origin, const helios::vec3 &proposed_fruit_axis, bool &collision_detection_active) const;
 
 
public:
    // ---- phytomer data ---- //
 
    std::vector<std::vector<helios::vec3>> petiole_vertices; // first index is petiole within internode, second index is tube segment within petiole tube
    std::vector<std::vector<helios::vec3>> leaf_bases; // first index is petiole within internode, second index is leaf within petiole
    std::vector<std::vector<std::vector<helios::vec3>>> peduncle_vertices; // first index is petiole within internode, second index is floral bud within petiole, third index is tube segment within peduncle
    std::vector<std::vector<std::vector<float>>> peduncle_radii; // first index is petiole within internode, second index is floral bud within petiole, third index is tube segment within peduncle
    std::vector<std::vector<float>> peduncle_length; // actual sampled length for each peduncle - first index is petiole, second is bud
    std::vector<std::vector<float>> peduncle_radius; // actual sampled radius for each peduncle - first index is petiole, second is bud
    std::vector<std::vector<float>> peduncle_pitch; // actual sampled pitch for each peduncle - first index is petiole, second is bud
    std::vector<std::vector<float>> peduncle_curvature; // actual sampled curvature for each peduncle - first index is petiole, second is bud
    float internode_pitch, internode_phyllotactic_angle;
 
    std::vector<std::vector<float>> petiole_radii; // first index is petiole within internode, second index is segment within petiole tube
    std::vector<float> petiole_length; // index is petiole within internode
    std::vector<float> petiole_pitch; // index is petiole within internode
    std::vector<float> petiole_curvature; // index is petiole within internode
    std::vector<float> petiole_taper; // taper value for each petiole (tip_radius/base_radius ratio)
    std::vector<helios::vec3> petiole_axis_initial; // initial petiole axis (before curvature) for each petiole
    std::vector<helios::vec3> petiole_rotation_axis; // rotation axis for curvature application for each petiole
    std::vector<std::vector<float>> leaf_size_max; // first index is petiole within internode, second index is leaf within petiole
    std::vector<std::vector<AxisRotation>> leaf_rotation; // first index is petiole within internode, second index is leaf within petiole
 
    std::vector<helios::RGBcolor> internode_colors; // index is segment within internode tube
    std::vector<helios::RGBcolor> petiole_colors; // index is segment within petiole tube
 
    std::vector<std::vector<uint>> petiole_objIDs; // first index is petiole within internode, second index is segment within petiole tube
    std::vector<std::vector<uint>> leaf_objIDs; // first index is petiole within internode, second index is leaf within petiole tube
 
    PhytomerParameters phytomer_parameters;
 
    uint rank;
    helios::int3 shoot_index;
 
    uint plantID;
    uint parent_shoot_ID;
    Shoot *parent_shoot_ptr;
 
    float age = 0;
    bool isdormant = false;
 
    float current_internode_scale_factor = 1;
    std::vector<float> current_leaf_scale_factor; // index is petiole within internode
 
    float old_phytomer_volume = 0;
 
    float downstream_leaf_area = 0;
 
    std::vector<std::vector<VegetativeBud>> axillary_vegetative_buds; // first index is petiole within internode, second index is bud within petiole
    std::vector<std::vector<FloralBud>> floral_buds; // first index is petiole within internode, second index is bud within petiole
 
    float internode_radius_initial;
    float internode_radius_max;
    float internode_length_max;
 
    std::vector<float> internode_curvature_perturbations; 
    std::vector<float> internode_yaw_perturbations; 
 
    bool build_context_geometry_petiole = true;
    bool build_context_geometry_peduncle = true;
 
protected:
    helios::vec3 inflorescence_bending_axis;
 
    helios::Context *context_ptr;
 
    PlantArchitecture *plantarchitecture_ptr;
 
    void updateInflorescence(FloralBud &fbud);
 
    void createInflorescenceGeometry(FloralBud &fbud, const helios::vec3 &fruit_base, const helios::vec3 &peduncle_axis, float pitch, float roll, float azimuth, float yaw_compound, float scale_factor, bool is_open_flower);
 
    [[nodiscard]] float calculatePhytomerConstructionCosts() const;
 
    [[nodiscard]] float calculateFlowerConstructionCosts(const FloralBud &fbud) const;
 
    [[nodiscard]] float calculateFruitConstructionCosts(const FloralBud &fbud) const;
 
    friend struct Shoot;
    friend class PlantArchitecture;
};
 
struct Shoot {
 
    Shoot(uint plant_ID, int shoot_ID, int parent_shoot_ID, uint parent_node, uint parent_petiole_index, uint rank, const helios::vec3 &shoot_base_position, const AxisRotation &shoot_base_rotation, uint current_node_number,
          float internode_length_shoot_initial, ShootParameters &shoot_params, std::string shoot_type_label, PlantArchitecture *plant_architecture_ptr);
 
 
    void buildShootPhytomers(float internode_radius, float internode_length, float internode_length_scale_factor_fraction, float leaf_scale_factor_fraction, float radius_taper);
 
 
    int appendPhytomer(float internode_radius, float internode_length_max, float internode_length_scale_factor_fraction, float leaf_scale_factor_fraction, const PhytomerParameters &phytomer_parameters);
 
 
    [[nodiscard]] std::string sampleChildShootType() const;
 
 
    [[nodiscard]] bool sampleVegetativeBudBreak(uint node_index) const;
 
    [[nodiscard]] bool sampleVegetativeBudBreak_carb(uint node_index) const;
 
 
    uint sampleEpicormicShoot(float dt, std::vector<float> &epicormic_positions_fraction) const;
 
    void terminateApicalBud();
 
    void terminateAxillaryVegetativeBuds();
 
    void addTerminalFloralBud();
 
    [[nodiscard]] float calculateShootInternodeVolume() const;
 
    [[nodiscard]] float calculateShootLength() const;
 
    [[nodiscard]] helios::vec3 getShootAxisVector(float shoot_fraction) const;
 
    [[nodiscard]] float sumShootLeafArea(uint start_node_index = 0) const;
 
    [[nodiscard]] float sumChildVolume(uint start_node_index = 0) const;
 
    void propagateDownstreamLeafArea(const Shoot *shoot, uint node_index, float leaf_area);
 
    void updateShootNodes(bool update_context_geometry = true);
 
    uint current_node_number = 0;
    uint nodes_this_season = 0;
 
    helios::vec3 base_position;
    AxisRotation base_rotation;
    helios::vec3 radial_outward_axis;
 
    const int ID;
    const int parent_shoot_ID;
    const uint plantID;
    const uint parent_node_index;
    const uint rank;
    const uint parent_petiole_index;
 
    float carbohydrate_pool_molC = 0; // mol C
 
    std::map<uint, float> leaf_nitrogen_gN_m2;
 
    float old_shoot_volume = 0;
 
    float phyllochron_increase = 5;
    float phyllochron_recovery = phyllochron_increase;
 
    float days_with_negative_carbon_balance = 0;
 
    void breakDormancy();
    void makeDormant();
 
    bool isdormant;
    uint dormancy_cycles = 0;
 
    bool meristem_is_alive = true;
 
    float phyllochron_counter = 0;
    float phyllochron_min = 6.f;
    float elongation_max = 0.25;
 
    float curvature_perturbation = 0;
    float yaw_perturbation = 0;
 
    float gravitropic_curvature = 0;
 
    const float internode_length_max_shoot_initial;
 
    uint internode_tube_objID = 4294967294;
 
    std::vector<std::vector<helios::vec3>> shoot_internode_vertices; // first index is phytomer within shoot, second index is segment within phytomer internode tube
    std::vector<std::vector<float>> shoot_internode_radii; // first index is phytomer within shoot, second index is segment within phytomer internode tube
 
    bool build_context_geometry_internode = true;
 
    // map of node number (key) to IDs of shoot children (value)
    std::map<int, std::vector<int>> childIDs;
 
    ShootParameters shoot_parameters;
 
    std::string shoot_type_label;
 
    float phyllochron_instantaneous;
    float elongation_rate_instantaneous;
 
    std::vector<std::shared_ptr<Phytomer>> phytomers;
 
    PlantArchitecture *plantarchitecture_ptr;
 
    helios::Context *context_ptr;
};
 
struct PlantInstance {
 
    PlantInstance(const helios::vec3 &a_base_position, float a_current_age, const std::string &a_plant_name, helios::Context *a_context_ptr) :
        base_position(a_base_position), current_age(a_current_age), plant_name(a_plant_name), context_ptr(a_context_ptr) {
    }
    std::vector<std::shared_ptr<Shoot>> shoot_tree;
    helios::vec3 base_position;
    float current_age;
    float time_since_dormancy = 0;
    helios::Context *context_ptr;
    std::string plant_name;
    std::pair<std::string, float> epicormic_shoot_probability_perlength_per_day; //.first is the epicormic shoot label string, .second is the probability
 
    // Phenological thresholds
    float dd_to_dormancy_break = 0;
    float dd_to_flower_initiation = 0;
    float dd_to_flower_opening = 0;
    float dd_to_fruit_set = 0;
    float dd_to_fruit_maturity = 0;
    float dd_to_dormancy = 0;
    float max_leaf_lifespan = 1e6;
    bool is_evergreen = false;
 
    float max_age = 999;
 
    CarbohydrateParameters carb_parameters;
 
    // --- Nitrogen Model --- //
 
    NitrogenParameters nitrogen_parameters;
 
    float root_nitrogen_pool_gN = 0;
    float available_nitrogen_pool_gN = 0;
    float cumulative_N_uptake_gN = 0;
 
    std::map<std::string, ShootParameters> shoot_types_snapshot;
 
    // --- Per-plant Attraction Points --- //
 
    bool attraction_points_enabled = false;
 
    std::vector<helios::vec3> attraction_points;
 
    float attraction_cone_half_angle_rad = 80.f * M_PI / 180.f;
 
    float attraction_cone_height = 0.1f;
 
    float attraction_weight = 0.6f;
 
    float attraction_obstacle_reduction_factor = 0.5f;
};
 
class PlantArchitecture {
public:
 
    explicit PlantArchitecture(helios::Context *context_ptr);
 
    ~PlantArchitecture();
 
    static int selfTest(int argc, char **argv);
 
 
    void optionalOutputObjectData(const std::string &object_data_label);
 
 
    void optionalOutputObjectData(const std::vector<std::string> &object_data_labels);
 
    // ********* Methods for Building Plants from Existing Library ********* //
 
 
    void loadPlantModelFromLibrary(const std::string &plant_label);
 
    [[nodiscard]] std::vector<std::string> getAvailablePlantModels() const;
 
 
    uint buildPlantInstanceFromLibrary(const helios::vec3 &base_position, float age, const std::map<std::string, float> &build_parameters = {});
 
 
    std::vector<uint> buildPlantCanopyFromLibrary(const helios::vec3 &canopy_center_position, const helios::vec2 &plant_spacing_xy, const helios::int2 &plant_count_xy, float age, float germination_rate = 1.f,
                                                  const std::map<std::string, float> &build_parameters = {});
 
 
    std::vector<uint> buildPlantCanopyFromLibrary(const helios::vec3 &canopy_center_position, const helios::vec2 &canopy_extent_xy, uint plant_count, float age, const std::map<std::string, float> &build_parameters = {});
 
 
    ShootParameters getCurrentShootParameters(const std::string &shoot_type_label);
 
 
    std::map<std::string, ShootParameters> getCurrentShootParameters();
 
 
    std::map<std::string, PhytomerParameters> getCurrentPhytomerParameters();
 
 
    [[nodiscard]] std::vector<std::string> listShootTypeLabels() const;
 
 
    [[nodiscard]] std::vector<std::string> listShootTypeLabels(const std::string &plant_model_name);
 
 
    [[nodiscard]] std::vector<std::string> listShootTypeLabels(uint plantID) const;
 
 
    void updateCurrentShootParameters(const std::string &shoot_type_label, const ShootParameters &params);
 
 
    void updateCurrentShootParameters(const std::map<std::string, ShootParameters> &params);
 
    // ********* Methods for Building Custom Plant Geometry from Scratch ********* //
 
 
    uint addPlantInstance(const helios::vec3 &base_position, float current_age);
 
 
    uint duplicatePlantInstance(uint plantID, const helios::vec3 &base_position, const AxisRotation &base_rotation, float current_age);
 
 
    void deletePlantInstance(uint plantID);
 
 
    void deletePlantInstance(const std::vector<uint> &plantIDs);
 
 
    void setPlantPhenologicalThresholds(uint plantID, float time_to_dormancy_break, float time_to_flower_initiation, float time_to_flower_opening, float time_to_fruit_set, float time_to_fruit_maturity, float time_to_dormancy,
                                        float max_leaf_lifespan = 1e6, bool is_evergreen = false);
 
 
    void setPlantCarbohydrateModelParameters(uint plantID, const CarbohydrateParameters &carb_parameters);
 
 
    void setPlantCarbohydrateModelParameters(const std::vector<uint> &plantIDs, const CarbohydrateParameters &carb_parameters);
 
    void disablePlantPhenology(uint plantID);
 
 
    void advanceTime(float time_step_days);
 
 
    void advanceTime(int time_step_years, float time_step_days);
 
 
    void advanceTime(uint plantID, float time_step_days);
 
 
    void advanceTime(const std::vector<uint> &plantIDs, float time_step_days);
 
 
    void adjustFruitForObstacleCollision();
 
 
    void accumulateHourlyLeafPhotosynthesis() const;
 
    // -- plant building methods -- //
 
 
    void defineShootType(const std::string &shoot_type_label, const ShootParameters &shoot_params);
 
 
    uint addBaseStemShoot(uint plantID, uint current_node_number, const AxisRotation &base_rotation, float internode_radius, float internode_length_max, float internode_length_scale_factor_fraction, float leaf_scale_factor_fraction,
                          float radius_taper, const std::string &shoot_type_label);
 
 
    uint appendShoot(uint plantID, int parent_shoot_ID, uint current_node_number, const AxisRotation &base_rotation, float internode_radius, float internode_length_max, float internode_length_scale_factor_fraction, float leaf_scale_factor_fraction,
                     float radius_taper, const std::string &shoot_type_label);
 
 
    uint addChildShoot(uint plantID, int parent_shoot_ID, uint parent_node_index, uint current_node_number, const AxisRotation &shoot_base_rotation, float internode_radius, float internode_length_max, float internode_length_scale_factor_fraction,
                       float leaf_scale_factor_fraction, float radius_taper, const std::string &shoot_type_label, uint petiole_index = 0);
 
 
    uint addEpicormicShoot(uint plantID, int parent_shoot_ID, float parent_position_fraction, uint current_node_number, float zenith_perturbation_degrees, float internode_radius, float internode_length_max,
                           float internode_length_scale_factor_fraction, float leaf_scale_factor_fraction, float radius_taper, const std::string &shoot_type_label);
 
 
    int appendPhytomerToShoot(uint plantID, uint shootID, const PhytomerParameters &phytomer_parameters, float internode_radius, float internode_length_max, float internode_length_scale_factor_fraction, float leaf_scale_factor_fraction);
 
 
    void enableEpicormicChildShoots(uint plantID, const std::string &epicormic_shoot_type_label, float epicormic_probability_perlength_perday);
 
    void disableInternodeContextBuild();
 
    void disablePetioleContextBuild();
 
    void disablePeduncleContextBuild();
 
 
    void enableGroundClipping(float ground_height = 0.f);
 
    // -- collision detection methods -- //
 
 
    void enableSoftCollisionAvoidance(const std::vector<uint> &target_object_UUIDs = {}, const std::vector<uint> &target_object_IDs = {}, bool enable_petiole_collision = false, bool enable_fruit_collision = false);
 
    void disableCollisionDetection();
 
 
    void setSoftCollisionAvoidanceParameters(float view_half_angle_deg, float look_ahead_distance, int sample_count, float inertia_weight);
 
 
    void setStaticObstacles(const std::vector<uint> &target_UUIDs);
 
 
    [[nodiscard]] CollisionDetection *getCollisionDetection() const;
 
 
    void setCollisionRelevantOrgans(bool include_internodes, bool include_leaves, bool include_petioles, bool include_flowers, bool include_fruit);
 
 
    void enableSolidObstacleAvoidance(const std::vector<uint> &obstacle_UUIDs, float avoidance_distance = 0.5f, bool enable_fruit_adjustment = false, bool enable_obstacle_pruning = false);
 
 
    void setGeometryUpdateScheduling(int update_frequency = 3, bool force_update_on_collision = true);
 
    // -- attraction points methods -- //
 
 
    void enableAttractionPoints(const std::vector<helios::vec3> &attraction_points, float view_half_angle_deg = 45.0f, float look_ahead_distance = 0.1f, float attraction_weight = 0.6f);
 
    void disableAttractionPoints();
 
 
    void updateAttractionPoints(const std::vector<helios::vec3> &attraction_points);
 
 
    void appendAttractionPoints(const std::vector<helios::vec3> &attraction_points);
 
 
    void setAttractionParameters(float view_half_angle_deg, float look_ahead_distance, float attraction_weight, float obstacle_reduction_factor = 0.75f);
 
 
    void enableAttractionPoints(uint plantID, const std::vector<helios::vec3> &attraction_points, float view_half_angle_deg = 80.0f, float look_ahead_distance = 0.1f, float attraction_weight = 0.6f);
 
 
    void disableAttractionPoints(uint plantID);
 
 
    void updateAttractionPoints(uint plantID, const std::vector<helios::vec3> &attraction_points);
 
 
    void appendAttractionPoints(uint plantID, const std::vector<helios::vec3> &attraction_points);
 
 
    void setAttractionParameters(uint plantID, float view_half_angle_deg, float look_ahead_distance, float attraction_weight, float obstacle_reduction_factor = 0.75f);
 
    bool detectAttractionPointsInCone(const helios::vec3 &vertex, const helios::vec3 &look_direction, float look_ahead_distance, float half_angle_degrees, helios::vec3 &direction_to_closest) const;
 
    bool detectAttractionPointsInCone(const std::vector<helios::vec3> &attraction_points, const helios::vec3 &vertex, const helios::vec3 &look_direction, float look_ahead_distance, float half_angle_degrees, helios::vec3 &direction_to_closest) const;
 
    // -- methods for modifying the current plant state -- //
 
    void initializeCarbohydratePool(float carbohydrate_concentration_molC_m3) const;
 
    void initializePlantCarbohydratePool(uint plantID, float carbohydrate_concentration_molC_m3);
 
    void initializeShootCarbohydratePool(uint plantID, uint shootID, float carbohydrate_concentration_molC_m3);
 
    void setPhytomerLeafScale(uint plantID, uint shootID, uint node_number, float leaf_scale_factor_fraction);
 
    void setPlantBasePosition(uint plantID, const helios::vec3 &base_position);
 
    void setPlantLeafElevationAngleDistribution(uint plantID, float Beta_mu_inclination, float Beta_nu_inclination) const;
 
    void setPlantLeafElevationAngleDistribution(const std::vector<uint> &plantIDs, float Beta_mu_inclination, float Beta_nu_inclination) const;
 
    void setPlantLeafAzimuthAngleDistribution(uint plantID, float eccentricity, float ellipse_rotation_degrees) const;
 
    void setPlantLeafAzimuthAngleDistribution(const std::vector<uint> &plantIDs, float eccentricity, float ellipse_rotation_degrees) const;
 
    void setPlantLeafAngleDistribution(uint plantID, float Beta_mu_inclination, float Beta_nu_inclination, float eccentricity, float ellipse_rotation_degrees) const;
 
    void setPlantLeafAngleDistribution(const std::vector<uint> &plantIDs, float Beta_mu_inclination, float Beta_nu_inclination, float eccentricity, float ellipse_rotation_degrees) const;
 
    void setPlantAge(uint plantID, float current_age);
 
    void harvestPlant(uint plantID);
 
    void removeShootLeaves(uint plantID, uint shootID);
 
    void removeShootVegetativeBuds(uint plantID, uint shootID);
 
    void removeShootFloralBuds(uint plantID, uint shootID);
 
    void removePlantLeaves(uint plantID);
 
    void makePlantDormant(uint plantID);
 
    void breakPlantDormancy(uint plantID);
 
    void pruneBranch(uint plantID, uint shootID, uint node_index);
 
    // -- methods for querying information about the plant -- //
 
    [[nodiscard]] std::string getPlantName(uint plantID) const;
 
    [[nodiscard]] float getPlantAge(uint plantID) const;
 
    [[nodiscard]] uint getShootNodeCount(uint plantID, uint shootID) const;
 
    [[nodiscard]] float getShootTaper(uint plantID, uint shootID) const;
 
    [[nodiscard]] helios::vec3 getPlantBasePosition(uint plantID) const;
 
    [[nodiscard]] std::vector<helios::vec3> getPlantBasePosition(const std::vector<uint> &plantIDs) const;
 
 
    [[nodiscard]] float sumPlantLeafArea(uint plantID) const;
 
 
    [[nodiscard]] float getPlantStemHeight(uint plantID) const;
 
 
    [[nodiscard]] float getPlantHeight(uint plantID) const;
 
 
    [[nodiscard]] std::vector<float> getPlantLeafInclinationAngleDistribution(uint plantID, uint Nbins, bool normalize = true) const;
 
 
    [[nodiscard]] std::vector<float> getPlantLeafInclinationAngleDistribution(const std::vector<uint> &plantIDs, uint Nbins, bool normalize = true) const;
 
 
    [[nodiscard]] std::vector<float> getPlantLeafAzimuthAngleDistribution(uint plantID, uint Nbins, bool normalize = true) const;
 
 
    [[nodiscard]] std::vector<float> getPlantLeafAzimuthAngleDistribution(const std::vector<uint> &plantIDs, uint Nbins, bool normalize = true) const;
 
 
    [[nodiscard]] uint getPlantLeafCount(uint plantID) const;
 
 
    [[nodiscard]] std::vector<helios::vec3> getPlantLeafBases(uint plantID) const;
 
 
    [[nodiscard]] std::vector<helios::vec3> getPlantLeafBases(const std::vector<uint> &plantIDs) const;
 
 
    [[nodiscard]] bool isPlantDormant(uint plantID) const;
 
 
    [[nodiscard]] std::string determinePhenologyStage(uint plantID) const;
 
 
    void writePlantMeshVertices(uint plantID, const std::string &filename) const;
 
 
    [[nodiscard]] std::vector<uint> getAllPlantIDs() const;
 
 
    [[nodiscard]] std::vector<uint> getAllPlantObjectIDs(uint plantID) const;
 
 
    [[nodiscard]] std::vector<uint> getAllPlantUUIDs(uint plantID) const;
 
 
    [[nodiscard]] std::vector<uint> getPlantInternodeObjectIDs(uint plantID) const;
 
 
    [[nodiscard]] std::vector<uint> getPlantInternodeObjectIDs(uint plantID, const std::string &shoot_type_label) const;
 
 
    [[nodiscard]] std::vector<uint> getPlantPetioleObjectIDs(uint plantID) const;
 
 
    [[nodiscard]] std::vector<uint> getPlantLeafObjectIDs(uint plantID) const;
 
 
    [[nodiscard]] std::vector<uint> getPlantLeafObjectIDs(const std::vector<uint> &plantIDs) const;
 
 
    [[nodiscard]] std::vector<uint> getPlantPeduncleObjectIDs(uint plantID) const;
 
 
    [[nodiscard]] std::vector<uint> getPlantFlowerObjectIDs(uint plantID) const;
 
 
    [[nodiscard]] std::vector<uint> getPlantFruitObjectIDs(uint plantID) const;
 
 
    void updateShootFruitCounts(uint plantID) const;
 
 
 
    [[nodiscard]] std::vector<uint> getShootInternodeObjectIDs(uint plantID) const;
 
 
    [[nodiscard]] std::vector<uint> getPlantCollisionRelevantObjectIDs(uint plantID) const;
 
 
    [[nodiscard]] std::vector<uint> getAllUUIDs() const;
 
 
    [[nodiscard]] std::vector<uint> getAllLeafUUIDs() const;
 
 
    [[nodiscard]] std::vector<uint> getAllInternodeUUIDs() const;
 
 
    [[nodiscard]] std::vector<uint> getAllPetioleUUIDs() const;
 
 
    [[nodiscard]] std::vector<uint> getAllPeduncleUUIDs() const;
 
 
    [[nodiscard]] std::vector<uint> getAllFlowerUUIDs() const;
 
 
    [[nodiscard]] std::vector<uint> getAllFruitUUIDs() const;
 
 
    [[nodiscard]] std::vector<uint> getAllObjectIDs() const;
 
    // -- carbohydrate model -- //
 
    void enableCarbohydrateModel();
 
    void disableCarbohydrateModel();
 
    // -- Nitrogen Model Methods -- //
 
    void enableNitrogenModel();
 
    void disableNitrogenModel();
 
    [[nodiscard]] bool isNitrogenModelEnabled() const;
 
    void setPlantNitrogenParameters(uint plantID, const NitrogenParameters &params);
 
    void setPlantNitrogenParameters(const std::vector<uint> &plantIDs, const NitrogenParameters &params);
 
    void initializeNitrogenPools(float initial_leaf_N_concentration);
 
    void initializePlantNitrogenPools(uint plantID, float initial_leaf_N_concentration);
 
    void addPlantNitrogen(uint plantID, float amount_gN);
 
    void addPlantNitrogen(const std::vector<uint> &plantIDs, float amount_gN);
 
    // -- manual plant generation from input string -- //
 
    [[nodiscard]] std::string getPlantString(uint plantID) const;
 
    uint generatePlantFromString(const std::string &generation_string, const PhytomerParameters &phytomer_parameters);
 
    uint generatePlantFromString(const std::string &generation_string, const std::map<std::string, PhytomerParameters> &phytomer_parameters);
 
    void writePlantStructureXML(uint plantID, const std::string &filename) const;
 
 
    void writeQSMCylinderFile(uint plantID, const std::string &filename) const;
 
    std::vector<uint> readPlantStructureXML(const std::string &filename, bool quiet = false);
 
    void disableMessages();
 
    void enableMessages();
 
 
    static std::string resolveTextureFile(const std::string &texture_file);
 
    friend struct Phytomer;
    friend struct Shoot;
 
private:
 
    float getParameterValue(const std::map<std::string, float> &build_parameters, const std::string &parameter_name, float default_value, float min_value, float max_value, const std::string &parameter_description) const;
 
    void clearBVHCache() const;
 
    void rebuildBVHForTimestep();
 
 
    void setPlantAttractionPoints(uint plantID, const std::vector<helios::vec3> &attraction_points, float view_half_angle_deg = 80.0f, float look_ahead_distance = 0.1f, float attraction_weight = 0.6f, float obstacle_reduction_factor = 0.75f);
 
 
    void ensureInflorescencePrototypesInitialized(const PhytomerParameters &params);
 
protected:
    helios::Context *context_ptr;
 
    std::minstd_rand0 *generator = nullptr;
 
    uint plant_count = 0;
 
    std::string current_plant_model;
 
    // Current build parameters for plant construction (set before calling builder functions)
    std::map<std::string, float> current_build_parameters;
 
    // Function pointer maps for plant model registration
    std::map<std::string, std::function<void()>> shoot_initializers;
    std::map<std::string, std::function<uint(const helios::vec3 &)>> plant_builders;
    std::map<std::string, std::string> plant_type_map;
 
    std::map<uint, PlantInstance> plant_instances;
 
    [[nodiscard]] std::string makeShootString(const std::string &current_string, const std::shared_ptr<Shoot> &shoot, const std::vector<std::shared_ptr<Shoot>> &shoot_tree) const;
 
    std::map<std::string, ShootParameters> shoot_types;
 
    // Key is the prototype function pointer; value first index is the unique leaf prototype, second index is the leaflet along a compound leaf (if applicable)
    // std::map<uint(*)(helios::Context* context_ptr, LeafPrototype* prototype_parameters, int compound_leaf_index),std::vector<std::vector<uint>> > unique_leaf_prototype_objIDs;
    std::map<uint, std::vector<std::vector<uint>>> unique_leaf_prototype_objIDs;
 
    // Key is the prototype function pointer; value index is the unique flower prototype
    std::map<uint (*)(helios::Context *context_ptr, uint subdivisions, bool flower_is_open), std::vector<uint>> unique_open_flower_prototype_objIDs;
    // Key is the prototype function pointer; value index is the unique flower prototype
    std::map<uint (*)(helios::Context *context_ptr, uint subdivisions, bool flower_is_open), std::vector<uint>> unique_closed_flower_prototype_objIDs;
    // Key is the prototype function pointer; value index is the unique fruit prototype
    std::map<uint (*)(helios::Context *context_ptr, uint subdivisions), std::vector<uint>> unique_fruit_prototype_objIDs;
 
    bool build_context_geometry_internode = true;
    bool build_context_geometry_petiole = true;
    bool build_context_geometry_peduncle = true;
 
    float ground_clipping_height = -99999;
 
    void validateShootTypes(ShootParameters &shoot_parameters, const std::map<std::string, ShootParameters> &shoot_types_ref) const;
 
    void registerPlantModel(const std::string &name, std::function<void()> shoot_init, std::function<uint(const helios::vec3 &)> plant_build, const std::string &plant_type = "herbaceous");
 
    void initializePlantModelRegistrations();
 
    void parseStringShoot(const std::string &LString_shoot, uint plantID, int parentID, uint parent_node, const std::map<std::string, PhytomerParameters> &phytomer_parameters, ShootParameters &shoot_parameters);
 
    void parseShootArgument(const std::string &shoot_argument, const std::map<std::string, PhytomerParameters> &phytomer_parameters, ShootParameters &shoot_parameters, AxisRotation &base_rotation, std::string &phytomer_label);
 
    void parseInternodeArgument(const std::string &internode_argument, float &internode_radius, float &internode_length, PhytomerParameters &phytomer_parameters);
 
    void parsePetioleArgument(const std::string &petiole_argument, PhytomerParameters &phytomer_parameters);
 
    void parseLeafArgument(const std::string &leaf_argument, PhytomerParameters &phytomer_parameters);
 
    void initializeDefaultShoots(const std::string &plant_label);
 
    [[nodiscard]] bool detectGroundCollision(uint objID);
 
    [[nodiscard]] bool detectGroundCollision(const std::vector<uint> &objID) const;
 
    void setPlantLeafAngleDistribution_private(const std::vector<uint> &plantIDs, float Beta_mu_inclination, float Beta_nu_inclination, float eccentricity_azimuth, float ellipse_rotation_azimuth_degrees, bool set_elevation, bool set_azimuth) const;
 
    static float interpolateTube(const std::vector<float> &P, float frac);
 
    static helios::vec3 interpolateTube(const std::vector<helios::vec3> &P, float frac);
 
    std::map<std::string, bool> output_object_data;
 
    // --- Plant Growth --- //
 
    void incrementPhytomerInternodeGirth(uint plantID, uint shootID, uint node_number, float dt, bool update_context_geometry);
    void incrementPhytomerInternodeGirth_carb(uint plantID, uint shootID, uint node_number, float dt, bool update_context_geometry);
 
    void pruneGroundCollisions(uint plantID);
 
    void pruneSolidBoundaryCollisions();
 
    // --- Carbohydrate Model --- //
 
    void accumulateShootPhotosynthesis() const;
 
    void subtractShootMaintenanceCarbon(float dt) const;
    void subtractShootGrowthCarbon();
 
    void checkCarbonPool_abortOrgans(float dt);
    void checkCarbonPool_adjustPhyllochron(float dt);
    void checkCarbonPool_transferCarbon(float dt);
 
    bool carbon_model_enabled = false;
 
    // --- Nitrogen Model --- //
 
    void accumulateLeafNitrogen(float dt);
    void remobilizeNitrogen(float dt);
    void updateNitrogenStressFactor();
    void removeFruitNitrogen();
 
    bool nitrogen_model_enabled = false;
 
    // --- Collision Detection --- //
 
    CollisionDetection *collision_detection_ptr = nullptr;
 
    bool owns_collision_detection = false;
 
    bool collision_detection_enabled = false;
 
    std::vector<uint> collision_target_UUIDs;
 
    std::vector<uint> collision_target_object_IDs;
 
    float collision_cone_half_angle_rad = 80.f * M_PI / 180.f;
 
    float collision_cone_height = 0.1f;
 
    int collision_sample_count = 256;
 
    float collision_inertia_weight = 0.4f;
 
    int geometry_update_frequency = 3;
 
    bool force_update_on_collision = true;
 
    bool collision_include_internodes = false;
    bool collision_include_leaves = true;
    bool collision_include_petioles = false;
    bool collision_include_flowers = false;
    bool collision_include_fruit = false;
 
    bool petiole_collision_detection_enabled = false;
 
    bool fruit_collision_detection_enabled = false;
 
    int geometry_update_counter = 0;
 
    mutable bool collision_avoidance_applied = false;
 
    bool spatial_filtering_enabled = false;
 
    float spatial_max_distance = 5.0f;
 
    mutable bool bvh_cached_for_current_growth = false;
    mutable std::vector<uint> cached_target_geometry;
    mutable std::vector<uint> cached_filtered_geometry;
 
    bool solid_obstacle_avoidance_enabled = false;
    std::vector<uint> solid_obstacle_UUIDs;
    float solid_obstacle_avoidance_distance = 0.5f;
    float solid_obstacle_minimum_distance = 0.05f;
    bool solid_obstacle_fruit_adjustment_enabled = false;
    bool solid_obstacle_pruning_enabled = false;
 
    // --- Attraction Points --- //
 
    bool attraction_points_enabled = false;
 
    std::vector<helios::vec3> attraction_points;
 
    float attraction_cone_half_angle_rad = 80.f * M_PI / 180.f;
 
    float attraction_cone_height = 0.1f;
 
    float attraction_weight = 0.6f;
 
    float attraction_obstacle_reduction_factor = 0.5f;
 
    bool printmessages = true;
 
    // --- Plant Library --- //
 
    void initializeAlmondTreeShoots();
 
    uint buildAlmondTree(const helios::vec3 &base_position);
 
    void initializeAlmondTreeAldrichShoots();
 
    uint buildAlmondTreeAldrich(const helios::vec3 &base_position);
 
    void initializeAlmondTreeWoodColonyShoots();
 
    uint buildAlmondTreeWoodColony(const helios::vec3 &base_position);
 
    void initializeAppleTreeShoots();
 
    uint buildAppleTree(const helios::vec3 &base_position);
 
    void initializeAppleFruitingWallShoots();
 
    uint buildAppleFruitingWall(const helios::vec3 &base_position);
 
    void initializeAsparagusShoots();
 
    uint buildAsparagusPlant(const helios::vec3 &base_position);
 
    void initializeBindweedShoots();
 
    uint buildBindweedPlant(const helios::vec3 &base_position);
 
    void initializeBeanShoots();
 
    uint buildBeanPlant(const helios::vec3 &base_position);
 
    void initializeBougainvilleaShoots();
 
    uint buildBougainvilleaPlant(const helios::vec3 &base_position);
 
    void initializeCapsicumShoots();
 
    uint buildCapsicumPlant(const helios::vec3 &base_position);
 
    void initializeCheeseweedShoots();
 
    uint buildCheeseweedPlant(const helios::vec3 &base_position);
 
    void initializeCowpeaShoots();
 
    uint buildCowpeaPlant(const helios::vec3 &base_position);
 
    void initializeGrapevineVSPShoots();
 
    uint buildGrapevineVSP(const helios::vec3 &base_position);
 
    void initializeGrapevineWyeShoots();
 
    uint buildGrapevineWye(const helios::vec3 &base_position);
 
    void initializeGroundCherryWeedShoots();
 
    uint buildGroundCherryWeedPlant(const helios::vec3 &base_position);
 
    void initializeMaizeShoots();
 
    uint buildMaizePlant(const helios::vec3 &base_position);
 
    void initializeOliveTreeShoots();
 
    uint buildOliveTree(const helios::vec3 &base_position);
 
    void initializePistachioTreeShoots();
 
    uint buildPistachioTree(const helios::vec3 &base_position);
 
    void initializePuncturevineShoots();
 
    uint buildPuncturevinePlant(const helios::vec3 &base_position);
 
    void initializeEasternRedbudShoots();
 
    uint buildEasternRedbudPlant(const helios::vec3 &base_position);
 
    void initializeRiceShoots();
 
    uint buildRicePlant(const helios::vec3 &base_position);
 
    void initializeButterLettuceShoots();
 
    uint buildButterLettucePlant(const helios::vec3 &base_position);
 
    void initializeSoybeanShoots();
 
    uint buildSoybeanPlant(const helios::vec3 &base_position);
 
    void initializeSorghumShoots();
 
    uint buildSorghumPlant(const helios::vec3 &base_position);
 
    void initializeStrawberryShoots();
 
    uint buildStrawberryPlant(const helios::vec3 &base_position);
 
    void initializeSugarbeetShoots();
 
    uint buildSugarbeetPlant(const helios::vec3 &base_position);
 
    void initializeTomatoShoots();
 
    uint buildTomatoPlant(const helios::vec3 &base_position);
 
    void initializeCherryTomatoShoots();
 
    uint buildCherryTomatoPlant(const helios::vec3 &base_position);
 
    void initializeWalnutTreeShoots();
 
    uint buildWalnutTree(const helios::vec3 &base_position);
 
    void initializeWheatShoots();
 
    uint buildWheatPlant(const helios::vec3 &base_position);
};
 
#include "Assets.h"
 
#endif // PLANT_ARCHITECTURE