models package¶
simpest.models
special
¶
Model modules for SIMPLACE and FraNchEstYN workflows.
fr_crop_model
¶
Daily crop growth step and disease damage mechanisms.
This module computes one day of crop growth together with the four disease damage mechanisms that couple the epidemiological state to crop performance:
- light stealers reduce intercepted radiation,
- RUE reducers lower radiation use efficiency,
- assimilate sappers drain fixed carbon, and
- senescence accelerators shorten the green-canopy duration.
Crop growth is computed through one of two branches:
- Internal growth model. When no external crop-model series is supplied, the canopy, biomass, and yield are simulated from thermal time using logistic light-interception curves and a radiation-use-efficiency biomass model.
- External growth model. When a daily crop-model series is supplied, the attainable light interception, biomass, and yield are taken from that series, and the damage mechanisms are applied to the daily increments to obtain the actual (disease-limited) trajectories.
In the external branch, day_after_sowing is incremented each simulated day
so that the runner's maturity/safety stop and the calibration objective's
"is-planted" logic operate correctly, and growing_degree_days is carried
from the external series. These fields therefore reflect the simulated crop
calendar rather than being left unset.
run(input_, parameters, output, output1)
¶
Compute one daily crop growth step.
The four disease damage mechanisms are derived first from the current
severity, then crop growth is advanced through either the internal growth
model (logistic light interception plus a radiation-use-efficiency biomass
model) or the external crop-model branch, depending on whether a daily
crop-model series is attached to input_. In both branches the attainable
(disease-free) and actual (disease-limited) light interception, biomass, and
yield are written to output1.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
input_ |
InputsDaily |
Today's daily inputs, including the optional external crop-model series. |
required |
parameters |
Parameters |
Model parameters, including the crop and disease groups. |
required |
output |
Outputs |
Previous day's output state, used as the starting state. |
required |
output1 |
Outputs |
Today's output state, updated in place. |
required |
Returns:
| Type | Description |
|---|---|
None |
This function updates |
Source code in simpest/models/fr_crop_model.py
def run(input_: InputsDaily, parameters: Parameters,
output: Outputs, output1: Outputs) -> None:
"""Compute one daily crop growth step.
The four disease damage mechanisms are derived first from the current
severity, then crop growth is advanced through either the internal growth
model (logistic light interception plus a radiation-use-efficiency biomass
model) or the external crop-model branch, depending on whether a daily
crop-model series is attached to ``input_``. In both branches the attainable
(disease-free) and actual (disease-limited) light interception, biomass, and
yield are written to ``output1``.
Args:
input_ (InputsDaily): Today's daily inputs, including the optional
external crop-model series.
parameters (Parameters): Model parameters, including the crop and disease
groups.
output (Outputs): Previous day's output state, used as the starting state.
output1 (Outputs): Today's output state, updated in place.
Returns:
None: This function updates ``output1`` in place.
"""
pc = parameters.par_crop
pd = parameters.par_disease
# -----------------------------------------------------------------------
# Damage mechanisms – always computed before selecting model branch
# -----------------------------------------------------------------------
severity = output.disease.disease_severity
output1.disease.damage_mechanisms.light_stealers = (
1.0 - _light_stealers_fn(severity, pd.light_stealer_damage)
)
output1.disease.damage_mechanisms.rue_reducers = (
1.0 - _rue_reduction_fn(severity, pd.rue_reducer_damage)
)
output1.disease.damage_mechanisms.assimilate_sappers = (
_assimilate_sappers_fn(severity, pd.assimilate_sappers_damage)
)
output1.disease.damage_mechanisms.senescence_accelerators = (
_senescence_accelerator_fn(severity, pd.senescence_accelerator_damage)
)
dm = output1.disease.damage_mechanisms
cmd = input_.crop_model_data
# -----------------------------------------------------------------------
# Branch A: internal crop growth model (no external crop-model series)
# -----------------------------------------------------------------------
if cmd is None or len(cmd.f_int) == 0:
t_ave = (input_.tmin + input_.tmax) / 2.0
t_func = t_response(t_ave, pc.tbase_crop, pc.topt_crop, pc.tmax_crop)
# Accumulate GDD
output1.crop.growing_degree_days = (
output.crop.growing_degree_days
+ t_func * (pc.topt_crop - pc.tbase_crop)
)
if output.crop.cycle_completion_percentage <= 100.0:
gdd = output1.crop.growing_degree_days
# Phenological code (1 = vegetative, 2 = reproductive)
output1.crop.pheno_code = _pheno_code(
gdd, pc.flowering_start / 100.0 * pc.cycle_length
)
# Attainable light interception
f_int_att, senescence_started = _f_int_compute(
pc.cycle_length, pc.slope_growth, pc.half_int_growth,
pc.slope_senescence, pc.half_int_senescence, gdd
)
output1.crop.light_interception_attainable = f_int_att
output1.crop.senescence_started = senescence_started
# Save peak fInt at the onset of senescence
if senescence_started and output1.crop.f_int_peak == 0.0:
output1.crop.f_int_peak = output.crop.light_interception_attainable
# Shift senescence half-point due to senescence accelerators
half_int_sen_shifted = (
pc.half_int_senescence - dm.senescence_accelerators * 100.0
)
# Actual light interception (with disease pressure)
f_int_act, _ = _f_int_compute(
pc.cycle_length, pc.slope_growth, pc.half_int_growth,
pc.slope_senescence, half_int_sen_shifted, gdd
)
output1.crop.light_interception_actual = f_int_act
# Apply light stealers
output1.crop.light_interception_actual -= (
output1.crop.light_interception_actual * dm.light_stealers
)
if output1.crop.light_interception_actual < 0.0:
output1.crop.light_interception_actual = 0.0
# Potential biomass accumulation
carbon_rate_pot = _carbon_rate(
pc.radiation_use_efficiency, input_.rad, t_func,
output1.crop.light_interception_attainable
) * 10.0 # g m⁻² → kg ha⁻¹
output1.crop.agb_attainable = output.crop.agb_attainable + carbon_rate_pot
# Actual biomass accumulation
carbon_rate_act = _carbon_rate(
pc.radiation_use_efficiency
- pc.radiation_use_efficiency * dm.rue_reducers,
input_.rad, t_func,
output1.crop.light_interception_actual
) * 10.0
carbon_rate_act -= dm.assimilate_sappers
output1.crop.agb_actual = output.crop.agb_actual + carbon_rate_act
if output1.crop.agb_actual < 0.0:
output1.crop.agb_actual = 0.0
# Potential yield
output1.crop.yield_attainable = output.crop.yield_attainable + _yield_rate(
output1.crop.pheno_code, carbon_rate_pot, pc.partitioning_maximum
)
# Actual yield
output1.crop.yield_actual = output.crop.yield_actual + _yield_rate(
output1.crop.pheno_code, carbon_rate_act, pc.partitioning_maximum
)
if output1.crop.yield_actual < 0.0:
output1.crop.yield_actual = 0.0
# Days after sowing
output1.crop.day_after_sowing = output.crop.day_after_sowing + 1
# Cycle completion (clamped at 100 %)
output1.crop.cycle_completion_percentage = min(
gdd / pc.cycle_length * 100.0, 100.0
)
else:
# Crop cycle complete → reset
output1.crop = CropOutputs()
output.crop = CropOutputs()
# -----------------------------------------------------------------------
# Branch B: External crop model supplies daily f_int, AGB, yield, GDD
# -----------------------------------------------------------------------
else:
# Normalise date key (InputsDaily.date is datetime; dict keys are date)
today: date = (
input_.date.date() if isinstance(input_.date, datetime) else input_.date
)
if today in cmd.f_int:
# --- Attainable values from external model ---
output1.crop.light_interception_attainable = cmd.f_int[today]
output1.crop.agb_attainable = cmd.agb[today]
output1.crop.yield_attainable = cmd.yield_.get(today, 0.0)
# --- Actual light interception ---
f_int_att = output1.crop.light_interception_attainable
f_int_act = (f_int_att - f_int_att * dm.light_stealers - dm.senescence_accelerators)
if f_int_act < 0.0:
f_int_act = 0.0
output1.crop.light_interception_actual = f_int_act
# Senescence flag: senescence starts when yield appears
if cmd.yield_.get(today, 0.0) > 0.0:
output1.crop.senescence_started = True
# Peak fInt across entire season
output1.crop.f_int_peak = max(cmd.f_int.values())
# --- Daily rates (attainable) ---
prev_day = today - timedelta(days=1)
if prev_day in cmd.agb:
pot_agb_rate = cmd.agb[today] - cmd.agb[prev_day]
pot_yield_rate = (cmd.yield_.get(today, 0.0) - cmd.yield_.get(prev_day, 0.0))
else: # first day of season
pot_agb_rate = cmd.agb[today]
pot_yield_rate = cmd.yield_.get(today, 0.0)
# --- Actual rates ---
if f_int_att > 0.0:
damage_ratio = (f_int_att - f_int_act) / f_int_att
act_agb_rate = pot_agb_rate - pot_agb_rate * damage_ratio
act_yield_rate = pot_yield_rate - pot_yield_rate * damage_ratio
else:
# First day, no light interception yet
act_agb_rate = output.crop.agb_actual
act_yield_rate = output.crop.yield_actual
# Clamp negative rates before applying further reductions
if act_agb_rate < 0.0:
act_agb_rate = 0.0
if act_yield_rate < 0.0:
act_yield_rate = 0.0
# Apply RUE reducers and assimilate sappers
act_agb_rate = (
act_agb_rate
- act_agb_rate * dm.rue_reducers
- dm.assimilate_sappers
)
act_yield_rate = (
act_yield_rate
- act_yield_rate * dm.rue_reducers
- dm.assimilate_sappers
)
# Update state variables
output1.crop.agb_actual = output.crop.agb_actual + act_agb_rate
output1.crop.yield_actual = output.crop.yield_actual + act_yield_rate
if output1.crop.agb_actual < 0.0:
output1.crop.agb_actual = 0.0
if output1.crop.yield_actual < 0.0:
output1.crop.yield_actual = 0.0
# Cycle completion from external model
output1.crop.cycle_completion_percentage = cmd.cycle_percentage.get(today, 0.0)
# Track the crop calendar: increment days-after-sowing so the
# maturity/safety stop and the calibration objective's is-planted
# logic operate correctly (see the module docstring).
output1.crop.day_after_sowing = output.crop.day_after_sowing + 1
# Carry thermal time from the external series. This is reporting-only
# when use_gdd=False and drives cycle completion when use_gdd=True.
output1.crop.growing_degree_days = cmd.gdd.get(today, 0.0)
else:
# Date not covered by external model → crop harvested, reset
output1.crop = CropOutputs()
output.crop = CropOutputs()
fr_data
¶
Input, output, and parameter data structures for the FraNchEstYN model.
This module defines the typed containers that carry information through the simulation: hourly and daily weather inputs, external crop-model series, the crop, disease, and fungicide state written on each simulated day, the parameter groups that configure each sub-model, and the per-site simulation unit that binds sowing schedule and reference observations together.
All quantities follow the model's canonical units, which are noted inline beside each field (temperatures in °C, biomass and yield in kg ha⁻¹, radiation in MJ m⁻², light interception as a 0–1 fraction, and so on).
CropModelData
dataclass
¶
Output from an external crop model, used as input to damage mechanisms.
Keys are datetime.date objects (year, month, day).
Source code in simpest/models/fr_data.py
@dataclass
class CropModelData:
"""Output from an external crop model, used as input to damage mechanisms.
Keys are datetime.date objects (year, month, day).
"""
f_int: Dict[date, float] = field(default_factory=dict) # fraction 0-1: light interception
yield_: Dict[date, float] = field(default_factory=dict) # kg ha⁻¹: dynamic crop yield
agb: Dict[date, float] = field(default_factory=dict) # kg ha⁻¹: above-ground biomass
cycle_percentage: Dict[date, float] = field(default_factory=dict) # %: crop cycle completion
# Growing degree days (thermal time) carried alongside the crop series. It is
# populated from the external crop model's thermal-time column and used to
# derive cycle progress when GDD-based phenology is enabled.
gdd: Dict[date, float] = field(default_factory=dict) # °C·day: growing degree days
CropOutputs
dataclass
¶
Crop growth state variables output each day.
Source code in simpest/models/fr_data.py
@dataclass
class CropOutputs:
"""Crop growth state variables output each day."""
growing_season: int = 0 # season identifier
day_after_sowing: int = 0 # days since sowing
pheno_code: float = 0.0 # phenological code (1 = vegetative, 2 = reproductive)
growing_degree_days: float = 0.0 # °C·day accumulated GDD
agb_attainable: float = 0.0 # kg ha⁻¹ attainable AGB
agb_actual: float = 0.0 # kg ha⁻¹ actual AGB
yield_attainable: float = 0.0 # kg ha⁻¹ attainable yield
yield_actual: float = 0.0 # kg ha⁻¹ actual yield
light_interception_attainable: float = 0.0 # fraction
light_interception_actual: float = 0.0 # fraction
cycle_completion_percentage: float = 0.0 # %
senescence_started: bool = False
f_int_peak: float = 0.0
DamageMechanisms
dataclass
¶
Damage outputs (one per day) that feed into crop loss calculations.
Source code in simpest/models/fr_data.py
@dataclass
class DamageMechanisms:
"""Damage outputs (one per day) that feed into crop loss calculations."""
light_stealers: float = 0.0 # fraction: reduction in light interception
senescence_accelerators: float = 0.0 # fraction: acceleration of tissue aging
rue_reducers: float = 0.0 # fraction: reduction in radiation use efficiency
assimilate_sappers: float = 0.0 # kg ha⁻¹: assimilate loss from disease
DiseaseOutputs
dataclass
¶
All disease state variables output each day.
Source code in simpest/models/fr_data.py
@dataclass
class DiseaseOutputs:
"""All disease state variables output each day."""
counter_dry: float = 0.0 # hours since last wet period
temp_function: float = 0.0 # temperature suitability (0-1)
rh_function: float = 0.0 # RH suitability (0-1)
hydro_thermal_time_state: float = 0.0 # accumulated hydro-thermal time (hydro-degree days)
hydro_thermal_time_rate: float = 0.0 # daily hydro-thermal time rate
hydro_thermal_time_infection: float = 0.0 # infection accumulation
sporulation_efficiency: float = 0.0 # %
tmax_daily: float = 0.0 # °C
tmin_daily: float = 0.0 # °C
rain_daily: float = 0.0 # mm
rhmax_daily: float = 0.0 # %
rhmin_daily: float = 0.0 # %
lw_daily: float = 0.0 # % (leaf wetness fraction of day)
# Tissue tracking is handled as a list of TissueState in disease_model.py
# Summary aggregates written here each day:
tissue_state: TissueState = field(default_factory=TissueState)
latent_sum: float = 0.0 # % cumulative latent
sporulating_sum: float = 0.0 # % cumulative sporulating
dead_sum: float = 0.0 # % cumulative dead
affected_sum: float = 0.0 # % cumulative affected
disease_severity: float = 0.0 # % disease severity (sporulating + dead)
susceptible_fraction: float = 1.0 # fraction of healthy susceptible tissue
damage_mechanisms: DamageMechanisms = field(default_factory=DamageMechanisms)
first_seasonal_infection: datetime = field(default_factory=lambda: datetime(1900, 1, 1))
cycle_percentage_first_infection: float = 0.0
is_primary_inoculum_started: bool = False
outer_inoculum: float = 0.0
FungicideOutputs
dataclass
¶
Fungicide state variables output each day.
Source code in simpest/models/fr_data.py
@dataclass
class FungicideOutputs:
"""Fungicide state variables output each day."""
concentration_factor: float = 0.0 # fraction - fungicide concentration
tenacity_function: float = 0.0 # unitless - fungicide tenacity function
tenacity: float = 0.0 # fraction - fungicide tenacity (potential efficacy after wash-off)
actual_degradation: float = 0.0 # fraction (1 = not degraded)
potential_degradation: float = 0.0 # fraction
efficacy: float = 0.0 # fraction
FungicideTreatmentSchedule
dataclass
¶
Holds all scheduled fungicide treatment dates for a simulation unit.
Source code in simpest/models/fr_data.py
@dataclass
class FungicideTreatmentSchedule:
"""Holds all scheduled fungicide treatment dates for a simulation unit."""
treatments: List[datetime] = field(default_factory=list)
def add_treatment(self, treatment_date: datetime) -> None:
self.treatments.append(treatment_date)
InputsDaily
dataclass
¶
Daily input record (assembled from 24 hourly records each day).
Source code in simpest/models/fr_data.py
@dataclass
class InputsDaily:
"""Daily input record (assembled from 24 hourly records each day)."""
date: datetime = datetime.min
tmax: float = 0.0 # °C
tmin: float = 0.0 # °C
rad: float = 0.0 # MJ m⁻² d⁻¹
rhx: float = 0.0 # % (max RH)
rhn: float = 0.0 # % (min RH)
precipitation: float = 0.0 # mm
leaf_wetness: float = 0.0 # hours
dew_point: float = 0.0 # °C
latitude: float = 0.0 # degrees
date_treatment_last: datetime = datetime.min
crop_model_data: "CropModelData" = None # populated from external crop model CSV
InputsHourly
dataclass
¶
Hourly input record for the disease model.
Source code in simpest/models/fr_data.py
@dataclass
class InputsHourly:
"""Hourly input record for the disease model."""
date: datetime = datetime.min
air_temperature: float = 0.0 # °C
precipitation: float = 0.0 # mm
relative_humidity: float = 0.0 # %
leaf_wetness: float = 0.0 # 0 or 1
rad: float = 0.0 # MJ m⁻² h⁻¹
dis_ideotype_potential_rate: float = 0.0 # unitless
latitude: float = 0.0 # degrees
date_treatment_last: datetime = datetime.min
Outputs
dataclass
¶
Top-level daily output bundle.
Source code in simpest/models/fr_data.py
@dataclass
class Outputs:
"""Top-level daily output bundle."""
disease: DiseaseOutputs = field(default_factory=DiseaseOutputs)
crop: CropOutputs = field(default_factory=CropOutputs)
inputs_daily: InputsDaily = field(default_factory=InputsDaily)
fungicide: FungicideOutputs = field(default_factory=FungicideOutputs)
ParCrop
dataclass
¶
Crop model parameters.
Source code in simpest/models/fr_data.py
@dataclass
class ParCrop:
"""Crop model parameters."""
tbase_crop: float = 0.0 # °C – base temperature for crop growth
topt_crop: float = 0.0 # °C – optimum temperature
tmax_crop: float = 0.0 # °C – maximum temperature
cycle_length: float = 0.0 # degree-days – total crop cycle length
flowering_start: float = 0.0 # % – crop cycle % when flowering starts
half_int_growth: float = 0.0 # % – cycle % when light interception = 0.5 (growth phase)
half_int_senescence: float = 0.0 # % – cycle % when light interception = 0.5 (senescence phase)
slope_growth: float = 0.0 # unitless – slope of growth light interception curve
slope_senescence: float = 0.0 # unitless – slope of senescence curve
radiation_use_efficiency: float = 0.0 # g MJ⁻¹
partitioning_maximum: float = 0.0 # unitless – max yield partitioning fraction
varietal_resistance: float = 0.0 # 0–1 – varietal resistance to pathogen
ParDisease
dataclass
¶
Disease model parameters.
Source code in simpest/models/fr_data.py
@dataclass
class ParDisease:
"""Disease model parameters."""
outer_inoculum_max: float = 0.0 # unitless – initial inoculum level
outer_inoculum_shape_release: float = 0.0 # 0/1/2 – release shape
outer_inoculum_shape_parameter: float = 0.0 # empirical parameter for shape 1 or 2
pathogen_spread: float = 0.0 # dispersal potential
wetness_duration_optimum: float = 0.0 # hours – optimal leaf wetness
wetness_duration_minimum: float = 0.0 # hours – minimum required wetness
dry_critical_interruption: float = 0.0 # hours – dry spell that disrupts infection
rain50_detachment: float = 0.0 # mm – rain needed for spore detachment
cycle_percentage_onset: float = 0.0 # % – crop cycle % when disease can start
tmin: float = 0.0 # °C – minimum temperature for pathogen
topt: float = 0.0 # °C – optimum temperature
tmax: float = 0.0 # °C – maximum temperature
relative_humidity_critical: float = 0.0 # % – minimum RH for infection
relative_humidity_not_limiting: float = 0.0 # % – RH above which it is not limiting
hydro_thermal_time_onset: float = 0.0 # hydro-degree days – onset threshold
latency_duration: float = 0.0 # days – latency period duration
sporulation_duration: float = 0.0 # days – sporulation duration
light_stealer_damage: float = 0.0 # unitless – lesion size representation
rue_reducer_damage: float = 0.0 # fraction – RUE reduction
senescence_accelerator_damage: float = 0.0 # fraction – senescence acceleration rate
assimilate_sappers_damage: float = 0.0 # kg ha⁻¹ – max assimilate loss
is_splash_borne: bool = False # True if spore release requires rain
ParFungicide
dataclass
¶
Fungicide model parameters.
Source code in simpest/models/fr_data.py
@dataclass
class ParFungicide:
"""Fungicide model parameters."""
a_shape_parameter: float = 0.0 # empirical degradation shape A
b_shape_parameter: float = 0.0 # empirical degradation shape B
degradation_rate: float = 0.0 # first-order degradation rate
initial_dose: float = 0.0 # fraction – initial dose (1 = max dose)
initial_efficacy: float = 0.0 # fraction – initial efficacy (1 = max)
tenacity_factor: float = 0.0 # rain-driven degradation empirical factor
Parameter
dataclass
¶
A single parameter definition with calibration bounds.
Holds the numeric value, the lower and upper bounds used as the calibration search space, and metadata describing the parameter's class and whether it is boolean. Instances are produced by the parameter readers and consumed by both the model runner and the calibration optimizer.
Source code in simpest/models/fr_data.py
@dataclass
class Parameter:
"""A single parameter definition with calibration bounds.
Holds the numeric value, the lower and upper bounds used as the calibration
search space, and metadata describing the parameter's class and whether it is
boolean. Instances are produced by the parameter readers and consumed by both
the model runner and the calibration optimizer.
"""
minimum: float = 0.0 # Lower bound (e.g., calibration search space)
maximum: float = 0.0 # Upper bound (e.g., calibration search space)
value: float = 0.0 # Current numeric value (e.g., for simulation or optimization)
calibration: str = "" # calibration subset tag (e.g., "included", "excluded")
param_class: str = "" # class/category of the parameter (e.g., "disease", "crop")
value_bool: bool = False # For boolean parameters, the value is stored here instead of 'value'
is_boolean: bool = False # True if this parameter is boolean (e.g., "IsSplashBorne"), in which case 'value_bool' is used instead of 'value'
Parameters
dataclass
¶
Container for all three parameter groups.
Source code in simpest/models/fr_data.py
@dataclass
class Parameters:
"""Container for all three parameter groups."""
par_disease: ParDisease = field(default_factory=ParDisease)
par_crop: ParCrop = field(default_factory=ParCrop)
par_fungicide: ParFungicide = field(default_factory=ParFungicide)
ReferenceData
dataclass
¶
Observed / reference data for validation (field measurements).
Source code in simpest/models/fr_data.py
@dataclass
class ReferenceData:
"""Observed / reference data for validation (field measurements)."""
date_fint: Dict[date, float] = field(default_factory=dict)
date_agb: Dict[date, float] = field(default_factory=dict)
date_yield_actual: Dict[date, float] = field(default_factory=dict)
date_yield_attainable: Dict[date, float] = field(default_factory=dict)
# disease_date_disease_sev[pathogen][date] = severity
disease_date_disease_sev: Dict[str, Dict[date, float]] = field(default_factory=dict)
SimulationUnit
dataclass
¶
Everything needed to run one site × variety combination.
Source code in simpest/models/fr_data.py
@dataclass
class SimulationUnit:
"""Everything needed to run one site × variety combination."""
site: str = "" # site identifier (e.g., location name)
crop: str = "" # crop identifier (e.g., crop name)
latitude: float = 0.0 # degrees - site latitude
longitude: float = 0.0 # degrees - site longitude
variety: str = "" # variety identifier (e.g., variety name)
pathogen: str = "" # pathogen identifier
variety_resistance: float = 0.0 # 0–1 – varietal resistance to pathogen
year_sowing_doy: Dict[int, int] = field(default_factory=dict) # Maps simulation year → sowing day-of-year
reference_data: ReferenceData = field(default_factory=ReferenceData) # observed or reference data for validation
fungicide_treatment_schedule: FungicideTreatmentSchedule = field( # fungicide treatment scheduling
default_factory=FungicideTreatmentSchedule
)
TissueState
dataclass
¶
State of a single tissue cohort (SEIR disease model).
Source code in simpest/models/fr_data.py
@dataclass
class TissueState:
"""State of a single tissue cohort (SEIR disease model)."""
latent_state: float = 0.0 # fraction currently in latent phase
sporulating_state: float = 0.0 # fraction currently sporulating
dead_state: float = 0.0 # fraction dead
latent_counter: float = 0.0 # progress counter in latent phase (→ 1 = transitions to sporulating)
sporulating_counter: float = 0.0 # progress counter in sporulating phase (→ 1 = transitions to dead)
fr_disease_model
¶
Hourly infection accumulation and daily SEIR tissue progression.
This module implements the epidemiological core of the model. Infection pressure is accumulated at an hourly time step from temperature and leaf-wetness suitability; once per day the accumulated pressure drives new infections and the progression of existing lesions through a compartmental SEIR scheme (susceptible → latent → sporulating → dead).
The :class:DiseaseModel is stateful and must be instantiated once per growing
season (the runner creates a fresh instance at sowing). Within a day,
:meth:DiseaseModel.run_hourly is called for every hour; at the end of the day
the runner advances the daily state and calls :meth:DiseaseModel.run_daily
with the aggregated weather.
DiseaseModel
¶
Stateful SEIR disease model.
Lifecycle¶
Create a fresh instance at the start of each growing season (sowing date). For each hour: call run_hourly(). After hour 23 (handled by the runner): call run_daily().
Source code in simpest/models/fr_disease_model.py
class DiseaseModel:
"""Stateful SEIR disease model.
Lifecycle
---------
Create a fresh instance at the start of each growing season (sowing date).
For each hour: call run_hourly().
After hour 23 (handled by the runner): call run_daily().
"""
def __init__(self) -> None:
# Hourly accumulators – cleared at end of each day by run_hourly()
self._temp: List[float] = []
self._rh: List[float] = []
self._rain: List[float] = []
self._lw: List[float] = []
# Tissue cohort list – persists across the entire season
self.tissue_tracking: List[TissueState] = []
# -----------------------------------------------------------------------
# Hourly step
# -----------------------------------------------------------------------
def run_hourly(
self,
input_: InputsHourly,
parameters: Parameters,
output: Outputs,
output1: Outputs,
) -> None:
"""Process one hour of weather data.
During hours 0–22 this accumulates the hydro-thermal time rate (the
product of temperature and relative-humidity suitability) and tracks the
dry-spell counter. At hour 23 it finalises the daily infection and
sporulation efficiencies, normalises the accumulated rate, records the
daily weather summary, and clears the hourly buffers.
The same ``output1`` object is shared across all 24 hours of a day; the
daily state advance happens only after the hour-23 call returns.
Args:
input_ (InputsHourly): The current hour's weather record.
parameters (Parameters): Model parameters, including the disease group.
output (Outputs): Previous day's output state, used to carry forward
the accumulated hydro-thermal time state.
output1 (Outputs): Current day's output state, updated in place.
Returns:
None: The function mutates ``output1`` in place.
"""
par = parameters.par_disease
# Accumulate hourly observations
self._temp.append(input_.air_temperature)
self._rh.append(input_.relative_humidity)
self._rain.append(input_.precipitation)
self._lw.append(input_.leaf_wetness)
# Dry-spell counter (resets on any wet hour, otherwise increments)
if input_.leaf_wetness == 1:
output1.disease.counter_dry = 0.0
else:
output1.disease.counter_dry += 1.0
# Temperature suitability (0–1)
output1.disease.temp_function = t_response(
input_.air_temperature, par.tmin, par.topt, par.tmax
)
hour = input_.date.hour
if hour == 23:
# ----------------------------------------------------------------
# End-of-day: finalise infection metrics and daily weather summary
# ----------------------------------------------------------------
rate = output1.disease.hydro_thermal_time_rate # accumulated hours 0-22
# Sporulation efficiency = average hourly HTT per hour of the day
output1.disease.sporulation_efficiency = rate / 24.0
# Infection efficiency (before normalisation against WetnessDurationOptimum)
if rate > par.wetness_duration_minimum:
htt_inf = rate / par.wetness_duration_optimum
if htt_inf >= 1.0:
htt_inf = 1.0
output1.disease.hydro_thermal_time_infection = htt_inf
else:
output1.disease.hydro_thermal_time_infection = 0.0
# Normalise daily rate (0–1)
rate_norm = rate / par.wetness_duration_optimum
if rate_norm >= 1.0:
rate_norm = 1.0
output1.disease.hydro_thermal_time_rate = rate_norm
# Cumulate state (carried forward day to day via output)
output1.disease.hydro_thermal_time_state = (
output.disease.hydro_thermal_time_state + rate_norm
)
# Dry-spell interruption overrides infection for this day
if output1.disease.counter_dry > par.dry_critical_interruption:
output1.disease.hydro_thermal_time_infection = 0.0
# Daily weather summary (for output only)
output1.disease.tmax_daily = max(self._temp)
output1.disease.tmin_daily = min(self._temp)
output1.disease.rhmax_daily = max(self._rh)
output1.disease.rhmin_daily = min(self._rh)
output1.disease.rain_daily = sum(self._rain)
output1.disease.lw_daily = sum(self._lw)
# Clear buffers for the next day
self._temp.clear()
self._rh.clear()
self._rain.clear()
self._lw.clear()
else:
# ----------------------------------------------------------------
# Within-day: accumulate HTT rate when within wet-period threshold
# ----------------------------------------------------------------
if output1.disease.counter_dry <= par.dry_critical_interruption:
rh_fn = _rh_function(
input_.relative_humidity,
input_.leaf_wetness,
par.relative_humidity_not_limiting,
par.relative_humidity_critical,
)
output1.disease.rh_function = rh_fn
output1.disease.hydro_thermal_time_rate += (
output1.disease.temp_function * rh_fn
)
else:
output1.disease.rh_function = 0.0
# -----------------------------------------------------------------------
# Daily step (called AFTER the runner swaps output ↔ output1 at hour 23)
# -----------------------------------------------------------------------
def run_daily(
self,
input_: InputsDaily,
parameters: Parameters,
output: Outputs,
output1: Outputs,
) -> None:
"""Run the daily SEIR tissue progression.
Given the day's accumulated infection pressure, this step (1) creates a
new latent tissue cohort when onset conditions are met, combining
external (primary) and internal (secondary) infection sources scaled by
susceptible tissue, host resistance, and fungicide efficacy; (2) advances
every existing cohort through the latent → sporulating → dead transitions
according to thermal time; and (3) aggregates the compartment fractions
into the daily latent, sporulating, dead, affected, and severity outputs,
and updates the susceptible fraction for the following day.
At entry, ``output`` holds the hourly-accumulated state for the day and
``output1`` is a fresh daily output whose season-persistent fields
(growing season, peak interception, and first-infection markers) have
already been carried forward by the runner.
Args:
input_ (InputsDaily): The day's aggregated inputs.
parameters (Parameters): Model parameters, including the disease and
crop groups.
output (Outputs): Previous (hourly-accumulated) state for the day.
output1 (Outputs): Current day's output state, updated in place.
Returns:
None: The function mutates ``output1`` and the cohort list in place.
"""
par = parameters.par_disease
pc = parameters.par_crop
# ----------------------------------------------------------------
# New infection – only once onset conditions are satisfied
# ----------------------------------------------------------------
if (output.disease.hydro_thermal_time_state >= par.hydro_thermal_time_onset
and output.crop.cycle_completion_percentage >= par.cycle_percentage_onset):
# Record first infection date (once per season)
if not output1.disease.is_primary_inoculum_started:
output1.disease.first_seasonal_infection = input_.date
output1.disease.cycle_percentage_first_infection = (
output.crop.cycle_completion_percentage
)
output1.disease.is_primary_inoculum_started = True
# Rescale existing tissue fractions for expanding green area
# (growth phase only; skip during senescence)
light_today = output1.crop.light_interception_attainable
light_yesterday = output.crop.light_interception_attainable
if light_today > 0.0:
if (output.disease.affected_sum <= 1.0
and not output1.crop.senescence_started):
for tissue in self.tissue_tracking:
ratio = light_yesterday / light_today
tissue.latent_state *= ratio
tissue.sporulating_state *= ratio
tissue.dead_state *= ratio
# Primary inoculum (reduced by fungicide efficacy)
outer_ino = _inoculum_model(input_, parameters, output, output1) * (
1.0 - output1.fungicide.efficacy
)
output1.disease.outer_inoculum = outer_ino
# Splash-borne: spore detachment efficiency
spo_detach = 0.0
if par.is_splash_borne:
spo_detach = 1.0 - rain_detachment(
input_.precipitation,
par.rain50_detachment,
output1.crop.light_interception_attainable,
)
# External infection (primary inoculum × HTT infection efficiency)
ext_inf = output.disease.hydro_thermal_time_infection * outer_ino
# Internal infection (secondary spread from sporulating tissue)
int_inf = (
output.disease.sporulating_sum
* output.disease.sporulation_efficiency
* par.pathogen_spread
* (1.0 - spo_detach)
)
# New infection cohort size
latent_value = (
(ext_inf + int_inf)
* output.disease.susceptible_fraction
* output1.crop.light_interception_attainable
* (1.0 - pc.varietal_resistance)
* (1.0 - output1.fungicide.efficacy)
)
if latent_value > 0.0:
self.tissue_tracking.append(TissueState(latent_state=latent_value))
# ----------------------------------------------------------------
# SEIR progression for all existing cohorts
# ----------------------------------------------------------------
t_ave = (input_.tmin + input_.tmax) / 2.0
gdd_disease = t_response(t_ave, par.tmin, par.topt, par.tmax)
lat_progress = (gdd_disease / par.latency_duration
if par.latency_duration > 0.0 else 0.0)
spo_progress = (gdd_disease / par.sporulation_duration
if par.sporulation_duration > 0.0 else 0.0)
efficacy = output1.fungicide.efficacy
for tissue in self.tissue_tracking:
if tissue.dead_state == 1.0:
continue
# Latent → sporulating
tissue.latent_counter += lat_progress * (1.0 - efficacy)
if tissue.latent_counter >= 1.0 and tissue.latent_state > 0.0:
tissue.sporulating_state = tissue.latent_state
tissue.latent_state = 0.0
# Sporulating → dead
if tissue.sporulating_state > 0.0:
tissue.sporulating_counter += spo_progress * (1.0 - efficacy)
if tissue.sporulating_counter >= 1.0 and tissue.sporulating_state > 0.0:
tissue.dead_state = tissue.sporulating_state
tissue.sporulating_state = 0.0
# ----------------------------------------------------------------
# Aggregate tissue states (cap each fraction at 1)
# ----------------------------------------------------------------
latent_sum = min(sum(t.latent_state for t in self.tissue_tracking), 1.0)
spor_sum = min(sum(t.sporulating_state for t in self.tissue_tracking), 1.0)
dead_sum = min(sum(t.dead_state for t in self.tissue_tracking), 1.0)
affected_sum = min(latent_sum + spor_sum + dead_sum, 1.0)
disease_severity = min(spor_sum + dead_sum, 1.0)
output1.disease.latent_sum = latent_sum
output1.disease.sporulating_sum = spor_sum
output1.disease.dead_sum = dead_sum
output1.disease.affected_sum = affected_sum
output1.disease.disease_severity = disease_severity
# ----------------------------------------------------------------
# Susceptible tissue fraction
# ----------------------------------------------------------------
tissue_availability = 1.0 - affected_sum
susceptible = max(0.0, min(tissue_availability, 1.0))
light_att_today = output1.crop.light_interception_attainable
light_att_prev = output.crop.light_interception_attainable
if (light_att_today >= light_att_prev
and not output1.crop.senescence_started):
# Growth phase: susceptible fraction = unaffected tissue
output1.disease.susceptible_fraction = susceptible
else:
# Senescence phase: also subtract senesced green area
f_int_peak = output1.crop.f_int_peak
if f_int_peak > 0.0:
sen_loss = (f_int_peak - light_att_today) / f_int_peak
else:
sen_loss = 0.0
output1.disease.susceptible_fraction = susceptible - sen_loss
output1.disease.susceptible_fraction = max(
0.0, min(output1.disease.susceptible_fraction, 1.0)
)
# ----------------------------------------------------------------
# Reset aggregates when no green tissue remains
# ----------------------------------------------------------------
if light_att_today == 0.0:
output1.disease.latent_sum = 0.0
output1.disease.sporulating_sum = 0.0
output1.disease.dead_sum = 0.0
output1.disease.affected_sum = 0.0
run_daily(self, input_, parameters, output, output1)
¶
Run the daily SEIR tissue progression.
Given the day's accumulated infection pressure, this step (1) creates a new latent tissue cohort when onset conditions are met, combining external (primary) and internal (secondary) infection sources scaled by susceptible tissue, host resistance, and fungicide efficacy; (2) advances every existing cohort through the latent → sporulating → dead transitions according to thermal time; and (3) aggregates the compartment fractions into the daily latent, sporulating, dead, affected, and severity outputs, and updates the susceptible fraction for the following day.
At entry, output holds the hourly-accumulated state for the day and
output1 is a fresh daily output whose season-persistent fields
(growing season, peak interception, and first-infection markers) have
already been carried forward by the runner.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
input_ |
InputsDaily |
The day's aggregated inputs. |
required |
parameters |
Parameters |
Model parameters, including the disease and crop groups. |
required |
output |
Outputs |
Previous (hourly-accumulated) state for the day. |
required |
output1 |
Outputs |
Current day's output state, updated in place. |
required |
Returns:
| Type | Description |
|---|---|
None |
The function mutates |
Source code in simpest/models/fr_disease_model.py
def run_daily(
self,
input_: InputsDaily,
parameters: Parameters,
output: Outputs,
output1: Outputs,
) -> None:
"""Run the daily SEIR tissue progression.
Given the day's accumulated infection pressure, this step (1) creates a
new latent tissue cohort when onset conditions are met, combining
external (primary) and internal (secondary) infection sources scaled by
susceptible tissue, host resistance, and fungicide efficacy; (2) advances
every existing cohort through the latent → sporulating → dead transitions
according to thermal time; and (3) aggregates the compartment fractions
into the daily latent, sporulating, dead, affected, and severity outputs,
and updates the susceptible fraction for the following day.
At entry, ``output`` holds the hourly-accumulated state for the day and
``output1`` is a fresh daily output whose season-persistent fields
(growing season, peak interception, and first-infection markers) have
already been carried forward by the runner.
Args:
input_ (InputsDaily): The day's aggregated inputs.
parameters (Parameters): Model parameters, including the disease and
crop groups.
output (Outputs): Previous (hourly-accumulated) state for the day.
output1 (Outputs): Current day's output state, updated in place.
Returns:
None: The function mutates ``output1`` and the cohort list in place.
"""
par = parameters.par_disease
pc = parameters.par_crop
# ----------------------------------------------------------------
# New infection – only once onset conditions are satisfied
# ----------------------------------------------------------------
if (output.disease.hydro_thermal_time_state >= par.hydro_thermal_time_onset
and output.crop.cycle_completion_percentage >= par.cycle_percentage_onset):
# Record first infection date (once per season)
if not output1.disease.is_primary_inoculum_started:
output1.disease.first_seasonal_infection = input_.date
output1.disease.cycle_percentage_first_infection = (
output.crop.cycle_completion_percentage
)
output1.disease.is_primary_inoculum_started = True
# Rescale existing tissue fractions for expanding green area
# (growth phase only; skip during senescence)
light_today = output1.crop.light_interception_attainable
light_yesterday = output.crop.light_interception_attainable
if light_today > 0.0:
if (output.disease.affected_sum <= 1.0
and not output1.crop.senescence_started):
for tissue in self.tissue_tracking:
ratio = light_yesterday / light_today
tissue.latent_state *= ratio
tissue.sporulating_state *= ratio
tissue.dead_state *= ratio
# Primary inoculum (reduced by fungicide efficacy)
outer_ino = _inoculum_model(input_, parameters, output, output1) * (
1.0 - output1.fungicide.efficacy
)
output1.disease.outer_inoculum = outer_ino
# Splash-borne: spore detachment efficiency
spo_detach = 0.0
if par.is_splash_borne:
spo_detach = 1.0 - rain_detachment(
input_.precipitation,
par.rain50_detachment,
output1.crop.light_interception_attainable,
)
# External infection (primary inoculum × HTT infection efficiency)
ext_inf = output.disease.hydro_thermal_time_infection * outer_ino
# Internal infection (secondary spread from sporulating tissue)
int_inf = (
output.disease.sporulating_sum
* output.disease.sporulation_efficiency
* par.pathogen_spread
* (1.0 - spo_detach)
)
# New infection cohort size
latent_value = (
(ext_inf + int_inf)
* output.disease.susceptible_fraction
* output1.crop.light_interception_attainable
* (1.0 - pc.varietal_resistance)
* (1.0 - output1.fungicide.efficacy)
)
if latent_value > 0.0:
self.tissue_tracking.append(TissueState(latent_state=latent_value))
# ----------------------------------------------------------------
# SEIR progression for all existing cohorts
# ----------------------------------------------------------------
t_ave = (input_.tmin + input_.tmax) / 2.0
gdd_disease = t_response(t_ave, par.tmin, par.topt, par.tmax)
lat_progress = (gdd_disease / par.latency_duration
if par.latency_duration > 0.0 else 0.0)
spo_progress = (gdd_disease / par.sporulation_duration
if par.sporulation_duration > 0.0 else 0.0)
efficacy = output1.fungicide.efficacy
for tissue in self.tissue_tracking:
if tissue.dead_state == 1.0:
continue
# Latent → sporulating
tissue.latent_counter += lat_progress * (1.0 - efficacy)
if tissue.latent_counter >= 1.0 and tissue.latent_state > 0.0:
tissue.sporulating_state = tissue.latent_state
tissue.latent_state = 0.0
# Sporulating → dead
if tissue.sporulating_state > 0.0:
tissue.sporulating_counter += spo_progress * (1.0 - efficacy)
if tissue.sporulating_counter >= 1.0 and tissue.sporulating_state > 0.0:
tissue.dead_state = tissue.sporulating_state
tissue.sporulating_state = 0.0
# ----------------------------------------------------------------
# Aggregate tissue states (cap each fraction at 1)
# ----------------------------------------------------------------
latent_sum = min(sum(t.latent_state for t in self.tissue_tracking), 1.0)
spor_sum = min(sum(t.sporulating_state for t in self.tissue_tracking), 1.0)
dead_sum = min(sum(t.dead_state for t in self.tissue_tracking), 1.0)
affected_sum = min(latent_sum + spor_sum + dead_sum, 1.0)
disease_severity = min(spor_sum + dead_sum, 1.0)
output1.disease.latent_sum = latent_sum
output1.disease.sporulating_sum = spor_sum
output1.disease.dead_sum = dead_sum
output1.disease.affected_sum = affected_sum
output1.disease.disease_severity = disease_severity
# ----------------------------------------------------------------
# Susceptible tissue fraction
# ----------------------------------------------------------------
tissue_availability = 1.0 - affected_sum
susceptible = max(0.0, min(tissue_availability, 1.0))
light_att_today = output1.crop.light_interception_attainable
light_att_prev = output.crop.light_interception_attainable
if (light_att_today >= light_att_prev
and not output1.crop.senescence_started):
# Growth phase: susceptible fraction = unaffected tissue
output1.disease.susceptible_fraction = susceptible
else:
# Senescence phase: also subtract senesced green area
f_int_peak = output1.crop.f_int_peak
if f_int_peak > 0.0:
sen_loss = (f_int_peak - light_att_today) / f_int_peak
else:
sen_loss = 0.0
output1.disease.susceptible_fraction = susceptible - sen_loss
output1.disease.susceptible_fraction = max(
0.0, min(output1.disease.susceptible_fraction, 1.0)
)
# ----------------------------------------------------------------
# Reset aggregates when no green tissue remains
# ----------------------------------------------------------------
if light_att_today == 0.0:
output1.disease.latent_sum = 0.0
output1.disease.sporulating_sum = 0.0
output1.disease.dead_sum = 0.0
output1.disease.affected_sum = 0.0
run_hourly(self, input_, parameters, output, output1)
¶
Process one hour of weather data.
During hours 0–22 this accumulates the hydro-thermal time rate (the product of temperature and relative-humidity suitability) and tracks the dry-spell counter. At hour 23 it finalises the daily infection and sporulation efficiencies, normalises the accumulated rate, records the daily weather summary, and clears the hourly buffers.
The same output1 object is shared across all 24 hours of a day; the
daily state advance happens only after the hour-23 call returns.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
input_ |
InputsHourly |
The current hour's weather record. |
required |
parameters |
Parameters |
Model parameters, including the disease group. |
required |
output |
Outputs |
Previous day's output state, used to carry forward the accumulated hydro-thermal time state. |
required |
output1 |
Outputs |
Current day's output state, updated in place. |
required |
Returns:
| Type | Description |
|---|---|
None |
The function mutates |
Source code in simpest/models/fr_disease_model.py
def run_hourly(
self,
input_: InputsHourly,
parameters: Parameters,
output: Outputs,
output1: Outputs,
) -> None:
"""Process one hour of weather data.
During hours 0–22 this accumulates the hydro-thermal time rate (the
product of temperature and relative-humidity suitability) and tracks the
dry-spell counter. At hour 23 it finalises the daily infection and
sporulation efficiencies, normalises the accumulated rate, records the
daily weather summary, and clears the hourly buffers.
The same ``output1`` object is shared across all 24 hours of a day; the
daily state advance happens only after the hour-23 call returns.
Args:
input_ (InputsHourly): The current hour's weather record.
parameters (Parameters): Model parameters, including the disease group.
output (Outputs): Previous day's output state, used to carry forward
the accumulated hydro-thermal time state.
output1 (Outputs): Current day's output state, updated in place.
Returns:
None: The function mutates ``output1`` in place.
"""
par = parameters.par_disease
# Accumulate hourly observations
self._temp.append(input_.air_temperature)
self._rh.append(input_.relative_humidity)
self._rain.append(input_.precipitation)
self._lw.append(input_.leaf_wetness)
# Dry-spell counter (resets on any wet hour, otherwise increments)
if input_.leaf_wetness == 1:
output1.disease.counter_dry = 0.0
else:
output1.disease.counter_dry += 1.0
# Temperature suitability (0–1)
output1.disease.temp_function = t_response(
input_.air_temperature, par.tmin, par.topt, par.tmax
)
hour = input_.date.hour
if hour == 23:
# ----------------------------------------------------------------
# End-of-day: finalise infection metrics and daily weather summary
# ----------------------------------------------------------------
rate = output1.disease.hydro_thermal_time_rate # accumulated hours 0-22
# Sporulation efficiency = average hourly HTT per hour of the day
output1.disease.sporulation_efficiency = rate / 24.0
# Infection efficiency (before normalisation against WetnessDurationOptimum)
if rate > par.wetness_duration_minimum:
htt_inf = rate / par.wetness_duration_optimum
if htt_inf >= 1.0:
htt_inf = 1.0
output1.disease.hydro_thermal_time_infection = htt_inf
else:
output1.disease.hydro_thermal_time_infection = 0.0
# Normalise daily rate (0–1)
rate_norm = rate / par.wetness_duration_optimum
if rate_norm >= 1.0:
rate_norm = 1.0
output1.disease.hydro_thermal_time_rate = rate_norm
# Cumulate state (carried forward day to day via output)
output1.disease.hydro_thermal_time_state = (
output.disease.hydro_thermal_time_state + rate_norm
)
# Dry-spell interruption overrides infection for this day
if output1.disease.counter_dry > par.dry_critical_interruption:
output1.disease.hydro_thermal_time_infection = 0.0
# Daily weather summary (for output only)
output1.disease.tmax_daily = max(self._temp)
output1.disease.tmin_daily = min(self._temp)
output1.disease.rhmax_daily = max(self._rh)
output1.disease.rhmin_daily = min(self._rh)
output1.disease.rain_daily = sum(self._rain)
output1.disease.lw_daily = sum(self._lw)
# Clear buffers for the next day
self._temp.clear()
self._rh.clear()
self._rain.clear()
self._lw.clear()
else:
# ----------------------------------------------------------------
# Within-day: accumulate HTT rate when within wet-period threshold
# ----------------------------------------------------------------
if output1.disease.counter_dry <= par.dry_critical_interruption:
rh_fn = _rh_function(
input_.relative_humidity,
input_.leaf_wetness,
par.relative_humidity_not_limiting,
par.relative_humidity_critical,
)
output1.disease.rh_function = rh_fn
output1.disease.hydro_thermal_time_rate += (
output1.disease.temp_function * rh_fn
)
else:
output1.disease.rh_function = 0.0
fr_fungicide_model
¶
Daily fungicide degradation, tenacity, and efficacy.
This module advances the fungicide state by one day. After an application, the active ingredient decays through two mechanisms: first-order chemical degradation over time and rainfall-driven wash-off (tenacity). The resulting effective dose is mapped to a protective efficacy through a logistic response, and the treatment is retired after a fixed persistence window.
run(input_, parameters, output, output1)
¶
Advance the fungicide state by one day.
Updates the concentration factor, tenacity, degradation, and protective efficacy for the current day. The calculation proceeds in stages:
- Concentration. On the application day the concentration factor is reset to 1; thereafter it decays exponentially with the previous day's effective degradation.
- Potential degradation. First-order decay of the initial dose,
initial_dose · exp(-degradation_rate · days). - Tenacity. Rainfall wash-off reduces the retained fraction by
exp(-tenacity_factor · concentration · √precipitation), applied cumulatively across days. - Effective degradation. The product of tenacity and potential degradation, i.e. the dose still active on the canopy.
- Efficacy. A logistic response of the effective degradation, scaled by the initial efficacy.
If no treatment has yet been applied (date_treatment_last.year <= 1) the
call is a no-op. All state is cleared once the persistence window of
_MAX_DAYS days has elapsed.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
input_ |
InputsDaily |
Daily inputs, including date, precipitation, and the last treatment date. |
required |
parameters |
Parameters |
Model parameters, including the fungicide group. |
required |
output |
Outputs |
Previous day's output state. |
required |
output1 |
Outputs |
Current day's output state, updated in place. |
required |
Returns:
| Type | Description |
|---|---|
None |
The function mutates |
Source code in simpest/models/fr_fungicide_model.py
def run(
input_: InputsDaily,
parameters: Parameters,
output: Outputs,
output1: Outputs,
) -> None:
"""Advance the fungicide state by one day.
Updates the concentration factor, tenacity, degradation, and protective
efficacy for the current day. The calculation proceeds in stages:
1. **Concentration.** On the application day the concentration factor is
reset to 1; thereafter it decays exponentially with the previous day's
effective degradation.
2. **Potential degradation.** First-order decay of the initial dose,
``initial_dose · exp(-degradation_rate · days)``.
3. **Tenacity.** Rainfall wash-off reduces the retained fraction by
``exp(-tenacity_factor · concentration · √precipitation)``, applied
cumulatively across days.
4. **Effective degradation.** The product of tenacity and potential
degradation, i.e. the dose still active on the canopy.
5. **Efficacy.** A logistic response of the effective degradation, scaled by
the initial efficacy.
If no treatment has yet been applied (``date_treatment_last.year <= 1``) the
call is a no-op. All state is cleared once the persistence window of
``_MAX_DAYS`` days has elapsed.
Args:
input_ (InputsDaily): Daily inputs, including date, precipitation, and
the last treatment date.
parameters (Parameters): Model parameters, including the fungicide group.
output (Outputs): Previous day's output state.
output1 (Outputs): Current day's output state, updated in place.
Returns:
None: The function mutates ``output1`` in place.
"""
pf = parameters.par_fungicide
last = input_.date_treatment_last
# Sentinel: DateTreatmentLast.year == 1 means no treatment yet
if last.year <= 1:
return
# Days since last application
days = (input_.date - last).total_seconds() / 86_400.0
# -----------------------------------------------------------------------
# Concentration factor
# -----------------------------------------------------------------------
if input_.date == last:
# Application day: baseline values
output1.fungicide.concentration_factor = 1.0
output1.fungicide.tenacity_function = 1.0
output.fungicide.tenacity = 1.0
else:
# After application: exponential decay based on yesterday's actual degradation
output1.fungicide.concentration_factor = math.exp(
(output.fungicide.actual_degradation - 1.0) * 3.0
)
# -----------------------------------------------------------------------
# Potential degradation (first-order decay of initial dose)
# -----------------------------------------------------------------------
output1.fungicide.potential_degradation = pf.initial_dose * math.exp(
-pf.degradation_rate * days
)
# -----------------------------------------------------------------------
# Tenacity (rainfall-driven wash-off)
# -----------------------------------------------------------------------
output.fungicide.tenacity_function = math.exp(
-pf.tenacity_factor
* output1.fungicide.concentration_factor
* math.sqrt(max(input_.precipitation, 0.0))
)
output1.fungicide.tenacity = (
output.fungicide.tenacity * output.fungicide.tenacity_function
)
# -----------------------------------------------------------------------
# Actual degradation
# -----------------------------------------------------------------------
output1.fungicide.actual_degradation = (
output1.fungicide.tenacity * output1.fungicide.potential_degradation
)
# -----------------------------------------------------------------------
# Efficacy (logistic response to actual degradation)
# -----------------------------------------------------------------------
output1.fungicide.efficacy = pf.initial_efficacy / (
1.0 + math.exp(
pf.a_shape_parameter
- pf.b_shape_parameter * output1.fungicide.actual_degradation
)
)
# -----------------------------------------------------------------------
# Hard stop: zero everything after 30 days
# -----------------------------------------------------------------------
if days >= _MAX_DAYS:
output1.fungicide.efficacy = 0.0
output1.fungicide.actual_degradation = 0.0
output1.fungicide.concentration_factor = 0.0
output1.fungicide.potential_degradation = 0.0
output1.fungicide.tenacity_function = 0.0
fr_optimizer
¶
Multi-start Nelder–Mead calibration for the FraNchEstYN model.
This module estimates model parameters by minimising the run's root-mean-square error against reference observations. It uses a self-contained, pure-Python multi-start Nelder–Mead (downhill simplex) search: several simplexes are each optimised from a random starting configuration and the best result across all restarts is returned. Multiple restarts reduce the chance of settling in a poor local minimum of the objective surface.
The objective is the runner's :meth:~simpest.models.fr_runner.FranchestynRunner.compute_rmse,
and candidate parameter sets that fall outside their bounds receive a large
penalty so the search remains in the feasible region. The search is controlled
by three knobs: n_restarts (number of independent simplexes), max_iter
(iterations per simplex), and ftol (objective-spread convergence tolerance).
Because the initial simplex is drawn at random, results vary from run to run
unless a fixed seed is supplied. For deterministic, fixed-parameter
validation, run the model directly rather than through calibration.
FranchestynOptimizer
¶
Multi-start Nelder–Mead calibration for the FraNchEstYN model.
Wraps a configured runner and searches for the parameter set that minimises the run's RMSE against reference data. Only parameters flagged for calibration (and not explicitly disabled) within the requested scope are optimised; all others are held at their default values.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
runner |
FranchestynRunner |
Fully configured runner instance. |
required |
calibration_variable |
str |
Calibration target scope: |
'all' |
n_restarts |
int |
Number of independent simplexes (random restarts). |
5 |
max_iter |
int |
Maximum iterations per simplex. |
1000 |
ftol |
float |
Convergence tolerance on the objective spread across the simplex vertices. |
1e-12 |
disabled_by_class |
Optional[Dict[str, Set[str]]] |
Parameter names to exclude from calibration, keyed by parameter class. |
None |
seed |
Optional[int] |
Seed for the random number generator, giving
reproducible restarts. |
None |
Source code in simpest/models/fr_optimizer.py
class FranchestynOptimizer:
"""Multi-start Nelder–Mead calibration for the FraNchEstYN model.
Wraps a configured runner and searches for the parameter set that minimises
the run's RMSE against reference data. Only parameters flagged for
calibration (and not explicitly disabled) within the requested scope are
optimised; all others are held at their default values.
Args:
runner (FranchestynRunner): Fully configured runner instance.
calibration_variable (str): Calibration target scope: ``"crop"``,
``"disease"``, or ``"all"``.
n_restarts (int): Number of independent simplexes (random restarts).
max_iter (int): Maximum iterations per simplex.
ftol (float): Convergence tolerance on the objective spread across the
simplex vertices.
disabled_by_class (Optional[Dict[str, Set[str]]]): Parameter names to
exclude from calibration, keyed by parameter class.
seed (Optional[int]): Seed for the random number generator, giving
reproducible restarts. ``None`` (default) is non-deterministic.
"""
def __init__(
self,
runner: FranchestynRunner,
calibration_variable: str = "all",
n_restarts: int = 5,
max_iter: int = 1000,
ftol: float = 1e-12,
disabled_by_class: Optional[Dict[str, Set[str]]] = None,
seed: Optional[int] = None,
) -> None:
self.runner = runner
self.calibration_variable = calibration_variable.lower()
self.n_restarts = n_restarts
self.max_iter = max_iter
self.ftol = ftol
self.disabled_by_class = disabled_by_class or {}
# Optional RNG seed making the multi-start search reproducible.
self.seed = seed
# Select calibration parameters and record their bounds
self.calib_keys, self.bounds = self._select_calib_params()
self._n_eval = 0
self._current_restart = 0
self._iter_in_restart = 0
self._last_rmse = math.inf
# -----------------------------------------------------------------------
# Public API
# -----------------------------------------------------------------------
def calibrate(self) -> Dict[str, float]:
"""
Run multi-start simplex and return the best parameter set.
Returns:
Dict[str, float]: Best-fit parameter values keyed by
``class_ParamName``.
"""
if not self.calib_keys:
print("No calibration parameters found — returning defaults.")
return {}
best_rmse = math.inf
best_params: Dict[str, float] = {}
rng = np.random.default_rng(self.seed)
print(f"- Calibrating {len(self.calib_keys)} using multi-start simplex method.\n"f"- Parameters:\n{self.calib_keys}")
for restart in range(self.n_restarts):
self._current_restart = restart + 1
self._iter_in_restart = 0
simplex, fvals, _nit = self._nelder_mead_single_restart(rng)
best_idx = int(np.argmin(fvals))
rmse = float(fvals[best_idx])
if rmse < best_rmse:
best_rmse = rmse
best_params = dict(zip(self.calib_keys, simplex[best_idx]))
print(f"\nBest RMSE: {best_rmse:.4f}")
return best_params
# -----------------------------------------------------------------------
# Private helpers
# -----------------------------------------------------------------------
def _objective(self, x: np.ndarray) -> float:
"""Evaluate the calibration objective at parameter vector ``x``.
Runs the model with the candidate parameters and returns the resulting
RMSE. Candidates that violate the parameter bounds, or that raise during
the run, return a large penalty (``1e300``) so the search avoids the
infeasible region.
"""
self._n_eval += 1
# Penalise out-of-bounds candidates with a large finite objective value
for val, (lo, hi) in zip(x, self.bounds):
if val <= lo or val > hi:
return 1e300
param_values = dict(zip(self.calib_keys, x))
try:
date_outputs = self.runner.run(param_values)
except Exception:
return 1e300
include_crop = self.calibration_variable in ("crop", "all")
include_disease = self.calibration_variable in ("disease", "all")
rmse = self.runner.compute_rmse(
date_outputs,
include_crop=include_crop,
include_disease=include_disease,
)
self._last_rmse = rmse
return rmse
def _on_iteration(self, _xk: np.ndarray) -> None:
"""Progress callback for each simplex iteration."""
self._iter_in_restart += 1
sys.stdout.write(
f"\rRun {self._current_restart}/{self.n_restarts} Iteration {self._iter_in_restart}/{self.max_iter} CURR RMSE={self._last_rmse:.4f}"
)
sys.stdout.flush()
def _select_calib_params(self) -> Tuple[List[str], List[Tuple[float, float]]]:
"""Return (keys, bounds) for parameters flagged for calibration."""
calib_keys: List[str] = []
bounds: List[Tuple[float, float]] = []
for key, p in self.runner.name_param.items():
# Skip non-calibrated params
if not p.calibration.strip():
continue
# Restrict to the requested calibration variable
param_class = key.split("_", 1)[0].lower()
if self.calibration_variable not in ("all", param_class):
continue
# Explicitly exclude user-deactivated parameters by name.
param_name = key.split("_", 1)[1] if "_" in key else key
if param_name in self.disabled_by_class.get(param_class, set()):
continue
# Skip boolean parameters
if p.is_boolean:
continue
calib_keys.append(key)
bounds.append((p.minimum, p.maximum))
return calib_keys, bounds
def _random_simplex(self, rng: np.random.Generator) -> np.ndarray:
"""Create a random initial simplex fully contained in parameter bounds."""
dim = len(self.bounds)
simplex = np.empty((dim + 1, dim), dtype=float)
# Anchor vertex
simplex[0] = np.array([rng.uniform(lo, hi) for lo, hi in self.bounds], dtype=float)
# One perturbed vertex per dimension
for i in range(dim):
v = simplex[0].copy()
lo, hi = self.bounds[i]
span = hi - lo
delta = rng.uniform(0.05 * span, 0.25 * span)
direction = -1.0 if rng.random() < 0.5 else 1.0
v[i] = np.clip(v[i] + direction * delta, lo, hi)
simplex[i + 1] = v
return simplex
def _nelder_mead_single_restart(
self, rng: np.random.Generator
) -> Tuple[np.ndarray, np.ndarray, int]:
"""Run one Nelder-Mead restart using standard coefficients."""
# Standard Nelder-Mead coefficients
alpha = 1.0 # reflection
gamma = 2.0 # expansion
rho = 0.5 # contraction
sigma = 0.5 # shrink
simplex = self._random_simplex(rng)
fvals = np.array([self._objective(v) for v in simplex], dtype=float)
nit = 0
while nit < self.max_iter:
order = np.argsort(fvals)
simplex = simplex[order]
fvals = fvals[order]
# Convergence check driven by the objective spread across vertices
if np.max(np.abs(fvals - fvals[0])) <= self.ftol:
break
centroid = np.mean(simplex[:-1], axis=0)
worst = simplex[-1]
# Reflection
xr = centroid + alpha * (centroid - worst)
fr = self._objective(xr)
if fvals[0] <= fr < fvals[-2]:
simplex[-1] = xr
fvals[-1] = fr
elif fr < fvals[0]:
# Expansion
xe = centroid + gamma * (xr - centroid)
fe = self._objective(xe)
if fe < fr:
simplex[-1] = xe
fvals[-1] = fe
else:
simplex[-1] = xr
fvals[-1] = fr
else:
# Contraction
if fr < fvals[-1]:
# Outside contraction
xc = centroid + rho * (xr - centroid)
fc = self._objective(xc)
if fc <= fr:
simplex[-1] = xc
fvals[-1] = fc
else:
# Shrink
best = simplex[0].copy()
for i in range(1, len(simplex)):
simplex[i] = best + sigma * (simplex[i] - best)
fvals[i] = self._objective(simplex[i])
else:
# Inside contraction
xc = centroid - rho * (centroid - worst)
fc = self._objective(xc)
if fc < fvals[-1]:
simplex[-1] = xc
fvals[-1] = fc
else:
# Shrink
best = simplex[0].copy()
for i in range(1, len(simplex)):
simplex[i] = best + sigma * (simplex[i] - best)
fvals[i] = self._objective(simplex[i])
nit += 1
self._on_iteration(simplex[0])
return simplex, fvals, nit
calibrate(self)
¶
Run multi-start simplex and return the best parameter set.
Returns:
| Type | Description |
|---|---|
Dict[str, float] |
Best-fit parameter values keyed by
|
Source code in simpest/models/fr_optimizer.py
def calibrate(self) -> Dict[str, float]:
"""
Run multi-start simplex and return the best parameter set.
Returns:
Dict[str, float]: Best-fit parameter values keyed by
``class_ParamName``.
"""
if not self.calib_keys:
print("No calibration parameters found — returning defaults.")
return {}
best_rmse = math.inf
best_params: Dict[str, float] = {}
rng = np.random.default_rng(self.seed)
print(f"- Calibrating {len(self.calib_keys)} using multi-start simplex method.\n"f"- Parameters:\n{self.calib_keys}")
for restart in range(self.n_restarts):
self._current_restart = restart + 1
self._iter_in_restart = 0
simplex, fvals, _nit = self._nelder_mead_single_restart(rng)
best_idx = int(np.argmin(fvals))
rmse = float(fvals[best_idx])
if rmse < best_rmse:
best_rmse = rmse
best_params = dict(zip(self.calib_keys, simplex[best_idx]))
print(f"\nBest RMSE: {best_rmse:.4f}")
return best_params
fr_param_reader
¶
Readers for model parameter definitions.
Parameters can be supplied either as a flat CSV table or as modular JSON files
organised by crop, disease, and fungicide type. Each reader returns a dictionary
keyed by "Class_ParamName" (for example "crop_TbaseCrop") whose values
are :class:~simpest.models.fr_data.Parameter objects carrying the numeric
value together with the calibration bounds.
Parameter definition CSV columns:
| Column | Field | Notes |
|---|---|---|
| 0 | Name | Parameter class |
| 1 | Class | Parameter name |
| 2 | Description | Ignored |
| 3 | Unit | Ignored |
| 4 | Min | Lower calibration bound |
| 5 | Max | Upper calibration bound |
| 6 | Value | Default value |
| 7 | CalibrationSubset | Calibration inclusion tag |
Calibrated-output CSV columns are Name, Class, and Value.
calibrated_read(file)
¶
Read a calibrated output CSV and return a dict keyed 'Name_Class'.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
file |
str | Path |
Path to the calibrated parameters CSV. If the file does not exist, an empty dict is returned. |
required |
Returns:
| Type | Description |
|---|---|
Dict[str, float] |
Dictionary mapping 'ParamName_ClassName' → calibrated float value. |
Source code in simpest/models/fr_param_reader.py
def calibrated_read(file: str | Path) -> Dict[str, float]:
"""Read a calibrated output CSV and return a dict keyed 'Name_Class'.
Args:
file: Path to the calibrated parameters CSV. If the file does not
exist, an empty dict is returned.
Returns:
Dictionary mapping 'ParamName_ClassName' → calibrated float value.
"""
result: Dict[str, float] = {}
path = Path(file)
if not path.exists():
return result
with path.open(newline="", encoding="utf-8-sig") as fh:
reader = csv.reader(fh)
next(reader, None) # skip header
for row in reader:
if not row or len(row) < 3:
continue
name = row[0].strip()
cls = row[1].strip()
try:
value = float(row[2].strip())
except ValueError:
continue
result[f"{name}_{cls}"] = value
return result
read(file, calibration_variable='')
¶
Read a parameter definition CSV and return a dict keyed 'Name_Class'.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
file |
str | Path |
Path to the parameter CSV file. |
required |
calibration_variable |
str |
Name of the variable being calibrated (currently reserved for future use; not applied as a filter here). |
'' |
Returns:
| Type | Description |
|---|---|
Dict[str, Parameter] |
Dictionary mapping 'ParamName_ClassName' → Parameter. |
Source code in simpest/models/fr_param_reader.py
def read(file: str | Path, calibration_variable: str = "") -> Dict[str, Parameter]:
"""Read a parameter definition CSV and return a dict keyed 'Name_Class'.
Args:
file: Path to the parameter CSV file.
calibration_variable: Name of the variable being calibrated (currently
reserved for future use; not applied as a filter here).
Returns:
Dictionary mapping 'ParamName_ClassName' → Parameter.
"""
result: Dict[str, Parameter] = {}
path = Path(file)
with path.open(newline="", encoding="utf-8-sig") as fh:
reader = csv.reader(fh)
next(reader, None) # skip header
for row in reader:
if not row or all(cell.strip() == "" for cell in row):
continue
# Ensure we have at least 8 columns
while len(row) < 8:
row.append("")
model_class = row[0].strip() # col 0: Model (e.g., "crop", "disease")
param_name = row[1].strip() # col 1: Parameter (e.g., "TbaseCrop")
raw_value = row[6].strip()
param = Parameter(param_class=model_class)
# IsSplashBorne is the only boolean parameter (values are "0" or "1")
is_bool_param = param_name.lower() == "issplashborne" and raw_value.lower() in (
"true", "false", "1", "0"
)
if is_bool_param:
param.value_bool = raw_value.lower() in ("true", "1")
param.is_boolean = True
else:
try:
param.value = float(raw_value)
param.minimum = float(row[4].strip())
param.maximum = float(row[5].strip())
except ValueError:
# Skip malformed rows
continue
param.calibration = row[7].strip() if len(row) > 7 else ""
key = f"{model_class}_{param_name}"
result[key] = param
return result
read_by_crop(file, crop_type='wheat')
¶
Read parameters from JSON file organized by crop type.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
file |
str | Path |
Path to the parameters_by_crop.json file. |
required |
crop_type |
str |
Crop type key (e.g., "wheat", "rice"). Defaults to "wheat". |
'wheat' |
Returns:
| Type | Description |
|---|---|
Dict[str, Parameter] |
Dictionary mapping 'ParamName_ClassName' → Parameter for the specified crop. Returns empty dict if crop_type not found. |
Source code in simpest/models/fr_param_reader.py
def read_by_crop(file: str | Path, crop_type: str = "wheat") -> Dict[str, Parameter]:
"""Read parameters from JSON file organized by crop type.
Args:
file: Path to the parameters_by_crop.json file.
crop_type: Crop type key (e.g., "wheat", "rice"). Defaults to "wheat".
Returns:
Dictionary mapping 'ParamName_ClassName' → Parameter for the specified crop.
Returns empty dict if crop_type not found.
"""
result: Dict[str, Parameter] = {}
path = Path(file)
if not path.exists():
raise FileNotFoundError(f"Parameter file not found: {path}")
with path.open(encoding="utf-8") as fh:
data = json.load(fh)
# Get the crop section
if crop_type not in data:
raise ValueError(f"Crop type '{crop_type}' not found in parameter file. Available: {list(data.keys())}")
crop_params = data[crop_type]
# Iterate through model classes (crop, disease, fungicide)
for model_class, params in crop_params.items():
for param_name, param_dict in params.items():
param = Parameter(param_class=model_class)
# IsSplashBorne is the only boolean parameter
if param_name.lower() == "issplashborne":
param.value_bool = bool(param_dict.get("value", 0))
param.is_boolean = True
else:
param.value = float(param_dict.get("value", 0.0))
param.minimum = float(param_dict.get("min", 0.0))
param.maximum = float(param_dict.get("max", 1.0))
param.calibration = "x" if param_dict.get("calibration", False) else ""
key = f"{model_class}_{param_name}"
result[key] = param
return result
read_crop_parameters(file, crop_type='wheat')
¶
Read crop parameters from crop_parameters.json.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
file |
str | Path |
Path to the crop_parameters.json file. |
required |
crop_type |
str |
Crop type key (e.g., "wheat", "rice"). Defaults to "wheat". |
'wheat' |
Returns:
| Type | Description |
|---|---|
Dict[str, Parameter] |
Dictionary mapping 'crop_ParamName' → Parameter for the specified crop type. Raises ValueError if crop_type not found. |
Source code in simpest/models/fr_param_reader.py
def read_crop_parameters(file: str | Path, crop_type: str = "wheat") -> Dict[str, Parameter]:
"""Read crop parameters from crop_parameters.json.
Args:
file: Path to the crop_parameters.json file.
crop_type: Crop type key (e.g., "wheat", "rice"). Defaults to "wheat".
Returns:
Dictionary mapping 'crop_ParamName' → Parameter for the specified crop type.
Raises ValueError if crop_type not found.
"""
result: Dict[str, Parameter] = {}
path = Path(file)
if not path.exists():
raise FileNotFoundError(f"Crop parameter file not found: {path}")
with path.open(encoding="utf-8") as fh:
data = json.load(fh)
if crop_type not in data:
raise ValueError(f"Crop type '{crop_type}' not found. Available: {list(data.keys())}")
crop_params = data[crop_type]
for param_name, param_dict in crop_params.items():
param = Parameter(param_class="crop")
param.value = float(param_dict.get("value", 0.0))
param.minimum = float(param_dict.get("min", 0.0))
param.maximum = float(param_dict.get("max", 1.0))
param.calibration = "x" if param_dict.get("calibration", False) else ""
key = f"crop_{param_name}"
result[key] = param
return result
read_disease_parameters(file, disease_type)
¶
Read disease parameters from disease_parameters.json.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
file |
str | Path |
Path to the disease_parameters.json file. |
required |
disease_type |
str |
Disease type key (e.g., "septoria", "brown_rust", "black_rust", "yellow_rust", "wheat_blast"). |
required |
Returns:
| Type | Description |
|---|---|
Dict[str, Parameter] |
Dictionary mapping 'disease_ParamName' → Parameter for the specified disease type. Raises ValueError if disease_type not found. |
Source code in simpest/models/fr_param_reader.py
def read_disease_parameters(file: str | Path, disease_type: str) -> Dict[str, Parameter]:
"""Read disease parameters from disease_parameters.json.
Args:
file: Path to the disease_parameters.json file.
disease_type: Disease type key (e.g., "septoria", "brown_rust", "black_rust",
"yellow_rust", "wheat_blast").
Returns:
Dictionary mapping 'disease_ParamName' → Parameter for the specified disease type.
Raises ValueError if disease_type not found.
"""
result: Dict[str, Parameter] = {}
path = Path(file)
if not path.exists():
raise FileNotFoundError(f"Disease parameter file not found: {path}")
with path.open(encoding="utf-8") as fh:
data = json.load(fh)
if disease_type not in data:
raise ValueError(f"Disease type '{disease_type}' not found. Available: {list(data.keys())}")
disease_params = data[disease_type]
for param_name, param_dict in disease_params.items():
param = Parameter(param_class="disease")
# IsSplashBorne is a boolean parameter
if param_name.lower() == "issplashborne":
param.value_bool = bool(param_dict.get("value", 0))
param.is_boolean = True
else:
param.value = float(param_dict.get("value", 0.0))
param.minimum = float(param_dict.get("min", 0.0))
param.maximum = float(param_dict.get("max", 1.0))
param.calibration = "x" if param_dict.get("calibration", False) else ""
key = f"disease_{param_name}"
result[key] = param
return result
read_fungicide_parameters(file, fungicide_type='protectant')
¶
Read fungicide parameters from fungicide_parameters.json.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
file |
str | Path |
Path to the fungicide_parameters.json file. |
required |
fungicide_type |
str |
Fungicide type key (e.g., "protectant"). Defaults to "protectant". |
'protectant' |
Returns:
| Type | Description |
|---|---|
Dict[str, Parameter] |
Dictionary mapping 'fungicide_ParamName' → Parameter for the specified type. Raises ValueError if fungicide_type not found. |
Source code in simpest/models/fr_param_reader.py
def read_fungicide_parameters(file: str | Path, fungicide_type: str = "protectant") -> Dict[str, Parameter]:
"""Read fungicide parameters from fungicide_parameters.json.
Args:
file: Path to the fungicide_parameters.json file.
fungicide_type: Fungicide type key (e.g., "protectant"). Defaults to "protectant".
Returns:
Dictionary mapping 'fungicide_ParamName' → Parameter for the specified type.
Raises ValueError if fungicide_type not found.
"""
result: Dict[str, Parameter] = {}
path = Path(file)
if not path.exists():
raise FileNotFoundError(f"Fungicide parameter file not found: {path}")
with path.open(encoding="utf-8") as fh:
data = json.load(fh)
if fungicide_type not in data:
raise ValueError(f"Fungicide type '{fungicide_type}' not found. Available: {list(data.keys())}")
fungicide_params = data[fungicide_type]
for param_name, param_dict in fungicide_params.items():
param = Parameter(param_class="fungicide")
param.value = float(param_dict.get("value", 0.0))
param.minimum = float(param_dict.get("min", 0.0))
param.maximum = float(param_dict.get("max", 1.0))
param.calibration = "x" if param_dict.get("calibration", False) else ""
key = f"fungicide_{param_name}"
result[key] = param
return result
fr_reference_reader
¶
Readers for sowing schedules, reference observations, and crop-model series.
This module loads the non-weather inputs that drive and constrain a run:
- the per-site sowing schedule and fungicide treatment dates,
- field reference observations used to score the model during calibration, and
- a daily crop-model series (light interception, biomass, yield, and optional thermal time) supplied by an external crop growth model.
Each reader populates the appropriate fields of a
:class:~simpest.models.fr_data.SimulationUnit or returns a
:class:~simpest.models.fr_data.CropModelData container.
read_crop_model_data(crop_model_file, use_gdd=False)
¶
Read an external crop-model series and compute cycle progress.
Loads the daily crop-model output (light interception, biomass, yield, and optionally thermal time) and segments it into growing cycles, computing the cycle-completion percentage for each day.
Expected CSV columns (matched case-insensitively against common aliases):
year, doy, fint, agb, yield, and optionally gdd.
A new cycle is detected when the day-of-year steps backwards without a year boundary (a new sowing) or when yield resets from above to at or below the 100 kg ha⁻¹ harvest threshold.
Cycle completion is computed as follows:
- When
use_gddisTrueand thermal time is available, progress within a cycle isgdd[date] / max_gdd_in_cycle * 100, which ties phenology to accumulated thermal time. - Otherwise, progress is interpolated linearly over calendar days between the first and last day of the cycle.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
crop_model_file |
str | Path |
Path to the directory containing |
required |
use_gdd |
bool |
Whether to derive cycle completion from thermal time rather than calendar days. |
False |
Returns:
| Type | Description |
|---|---|
Populated |
Source code in simpest/models/fr_reference_reader.py
def read_crop_model_data(crop_model_file: str | Path, use_gdd: bool = False) -> CropModelData:
"""Read an external crop-model series and compute cycle progress.
Loads the daily crop-model output (light interception, biomass, yield, and
optionally thermal time) and segments it into growing cycles, computing the
cycle-completion percentage for each day.
Expected CSV columns (matched case-insensitively against common aliases):
``year``, ``doy``, ``fint``, ``agb``, ``yield``, and optionally ``gdd``.
A new cycle is detected when the day-of-year steps backwards without a
year boundary (a new sowing) or when yield resets from above to at or below
the 100 kg ha⁻¹ harvest threshold.
Cycle completion is computed as follows:
- When ``use_gdd`` is ``True`` and thermal time is available, progress within
a cycle is ``gdd[date] / max_gdd_in_cycle * 100``, which ties phenology to
accumulated thermal time.
- Otherwise, progress is interpolated linearly over calendar days between the
first and last day of the cycle.
Args:
crop_model_file: Path to the directory containing ``cropModelData.csv``,
or the path to the CSV file itself.
use_gdd: Whether to derive cycle completion from thermal time rather than
calendar days.
Returns:
Populated :class:`~simpest.models.fr_data.CropModelData`.
"""
cmd = CropModelData()
path = Path(crop_model_file)
if path.is_dir():
path = path / "cropModelData.csv"
with path.open(newline="", encoding="utf-8-sig") as fh:
reader = csv.reader(fh)
raw_headers = next(reader, [])
header_map = {_norm(h): i for i, h in enumerate(raw_headers)}
fint_col = _get_col(header_map, "fint", "f_int", "lightinterception", "lightint")
agb_col = _get_col(header_map, "agb", "abovegroundbiomass", "biomass", "wtop")
yield_col = _get_col(header_map, "yield", "yieldattainable", "yieldunlimited", "yieldpotential", "wgrn", "grainyieldpotential")
year_col = _get_col(header_map, "year", "yr")
doy_col = _get_col(header_map, "doy", "dayofyear", "dy", "d")
gdd_col = _get_col(header_map, "gdd", "growingdegreedays", "tsum", "thermaltime")
date_order: List[datetime] = []
for row in reader:
if not row or all(c.strip() == "" for c in row):
continue
if year_col < 0 or doy_col < 0:
continue
try:
y = int(row[year_col].strip())
doy = int(row[doy_col].strip())
except (ValueError, IndexError):
continue
try:
dt = datetime(y, 1, 1) + timedelta(days=doy - 1)
except ValueError:
continue
d = dt.date()
if fint_col >= 0 and fint_col < len(row):
v = _pf(row[fint_col])
if v is not None:
cmd.f_int[d] = v
if agb_col >= 0 and agb_col < len(row):
v = _pf(row[agb_col])
if v is not None:
cmd.agb[d] = v
if yield_col >= 0 and yield_col < len(row):
v = _pf(row[yield_col])
if v is not None:
cmd.yield_[d] = v
if use_gdd and gdd_col >= 0 and gdd_col < len(row):
v = _pf(row[gdd_col])
if v is not None:
cmd.gdd[d] = v
if d in cmd.f_int or d in cmd.agb or d in cmd.yield_:
date_order.append(dt)
if not date_order:
return cmd
# --- Cycle detection ---
date_order_sorted = sorted(set(date_order))
dates_d = [dt.date() for dt in date_order_sorted]
cycles: List[Tuple[date, date]] = []
cycle_start = dates_d[0]
for i in range(1, len(dates_d)):
prev = dates_d[i - 1]
curr = dates_d[i]
doy_backwards = curr.timetuple().tm_yday < prev.timetuple().tm_yday
year_wrap = prev.month == 12 and curr.month == 1
sowing_jump = doy_backwards and not year_wrap
y_prev = cmd.yield_.get(prev, 0.0)
y_curr = cmd.yield_.get(curr, 0.0)
harvest_reset = (y_curr <= 100.0 and y_prev > 100.0)
if sowing_jump or harvest_reset:
cycles.append((cycle_start, prev))
cycle_start = curr
cycles.append((cycle_start, dates_d[-1]))
# --- Compute cycle percentage ---
for start, end in cycles:
cycle_dates = [d for d in dates_d if start <= d <= end]
if not cycle_dates:
continue
gdd_max = 0.0
if use_gdd:
gdd_values = [cmd.gdd.get(d, 0.0) for d in cycle_dates]
gdd_max = max(gdd_values)
if use_gdd and gdd_max > 0.0:
# Thermal-time progression: scale by the cycle's maximum GDD
for d in cycle_dates:
gdd_d = cmd.gdd.get(d, 0.0)
cmd.cycle_percentage[d] = min(100.0, gdd_d / gdd_max * 100.0)
else:
# Calendar-day interpolation between cycle start and end
total_days = (end - start).days
if total_days <= 0:
continue
for d in cycle_dates:
frac = (d - start).days / total_days
cmd.cycle_percentage[d] = frac * 100.0
return cmd
read_reference(ref_dir, sowing_file, site, variety, start_year, end_year, sim_unit=None, disease='thisDisease')
¶
Read referenceData.csv and populate the SimulationUnit's reference_data.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
ref_dir |
str | Path |
Directory containing referenceData.csv. |
required |
sowing_file |
str | Path |
Path to sowing.csv (loaded first if sim_unit is empty). |
required |
site |
str |
Site identifier. |
required |
variety |
str |
Variety identifier. |
required |
start_year |
int |
First year. |
required |
end_year |
int |
Last year. |
required |
sim_unit |
Optional[SimulationUnit] |
Existing SimulationUnit to populate; created if None. |
None |
disease |
str |
Disease column name to look for (e.g., "thisDisease", "stripe_rust"). |
'thisDisease' |
Returns:
| Type | Description |
|---|---|
SimulationUnit |
Updated SimulationUnit with reference_data populated. |
Source code in simpest/models/fr_reference_reader.py
def read_reference(
ref_dir: str | Path,
sowing_file: str | Path,
site: str,
variety: str,
start_year: int,
end_year: int,
sim_unit: Optional[SimulationUnit] = None,
disease: str = "thisDisease",
) -> SimulationUnit:
"""Read referenceData.csv and populate the SimulationUnit's reference_data.
Args:
ref_dir: Directory containing referenceData.csv.
sowing_file: Path to sowing.csv (loaded first if sim_unit is empty).
site: Site identifier.
variety: Variety identifier.
start_year: First year.
end_year: Last year.
sim_unit: Existing SimulationUnit to populate; created if None.
disease: Disease column name to look for (e.g., "thisDisease", "stripe_rust").
Returns:
Updated SimulationUnit with reference_data populated.
"""
if sim_unit is None or not sim_unit.year_sowing_doy:
sim_unit = read_sowing(sowing_file, site, variety, start_year, end_year)
sim_unit.site = site
_rp = Path(ref_dir)
path = _rp if _rp.is_file() else _rp / "referenceData.csv"
if not path.exists():
raise FileNotFoundError(f"referenceData.csv not found at '{path}'")
# return sim_unit # no reference data available
with path.open(newline="", encoding="utf-8-sig") as fh:
reader = csv.reader(fh)
raw_headers = next(reader, [])
header_map = {_norm(h): i for i, h in enumerate(raw_headers)}
fint_col = _get_col(header_map, "fint", "f_int", "lightinterception", "lightint")
agb_col = _get_col(header_map, "agb", "abovegroundbiomass", "biomass", "wtop")
year_col = _get_col(header_map, "year", "yr")
doy_col = _get_col(header_map, "doy", "dayofyear", "dy", "d")
variety_col = _get_col(header_map, "variety", "cultivar", "cv")
dis_col = _get_col(header_map, disease, f"{disease}sev", f"{disease}severity")
yield_att_col = _get_col(header_map, "yieldattainable", "yieldunlimited", "yieldpotential",
"yield", "wgrn", "grainyieldpotential")
yield_act_col = _get_col(header_map, "yieldactual", "yielddiseased", "yieldact",
"yieldlimited", "grainyieldlimited")
if dis_col < 0:
warnings.warn(
f"[read_reference] Column for disease '{disease}' not found in {path}. "
f"Available columns: {', '.join(raw_headers)}",
stacklevel=2,
)
for row in reader:
if not row or all(c.strip() == "" for c in row):
continue
if variety_col >= 0:
row_var = row[variety_col].strip().strip('"').lower() if variety_col < len(row) else ""
if row_var != variety.lower():
continue
# Parse date from year + doy
obs_date = None
if year_col >= 0 and doy_col >= 0 and year_col < len(row) and doy_col < len(row):
try:
y = int(row[year_col].strip())
doy = int(row[doy_col].strip())
if 1 <= doy <= 366:
obs_date = (datetime(y, 1, 1) + timedelta(days=doy - 1)).date()
except ValueError:
pass
if obs_date is None:
obs_date = datetime.min.date() # sentinel date for undated observations
# FINT
if fint_col >= 0:
v = _pf(row[fint_col]) if fint_col < len(row) else None
if v is not None:
sim_unit.reference_data.date_fint[obs_date] = v
# AGB
if agb_col >= 0:
v = _pf(row[agb_col]) if agb_col < len(row) else None
if v is not None:
sim_unit.reference_data.date_agb[obs_date] = v
# Yield attainable
if yield_att_col >= 0:
v = _pf(row[yield_att_col]) if yield_att_col < len(row) else None
if v is not None:
sim_unit.reference_data.date_yield_attainable[obs_date] = v
# Yield actual
if yield_act_col >= 0:
v = _pf(row[yield_act_col]) if yield_act_col < len(row) else None
if v is not None:
sim_unit.reference_data.date_yield_actual[obs_date] = v
# Disease severity
if dis_col >= 0:
v = _pf(row[dis_col]) if dis_col < len(row) else None
if v is not None:
sim_unit.reference_data.disease_date_disease_sev.setdefault(disease, {})[obs_date] = v
return sim_unit
read_sowing(sowing_file, site, variety, start_year, end_year, all_row_includes_end_year=False)
¶
Read sowing.csv and build the SimulationUnit for one site × variety.
Expected CSV columns (case-insensitive): site, crop, variety, sowingDOY, year (+ optional treatment1, treatment2, ... columns for fungicide dates)
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
sowing_file |
str | Path |
Path to the sowing CSV. |
required |
site |
str |
Site identifier to filter on (e.g., "indiana"). |
required |
variety |
str |
Variety identifier to filter on (e.g., "Generic"). |
required |
start_year |
int |
First simulation year (inclusive). |
required |
end_year |
int |
Last simulation year (inclusive). |
required |
all_row_includes_end_year |
bool |
How far an |
False |
Returns:
| Type | Description |
|---|---|
SimulationUnit |
Populated SimulationUnit (year_sowing_doy, fungicide_treatment_schedule). |
Source code in simpest/models/fr_reference_reader.py
def read_sowing(
sowing_file: str | Path,
site: str,
variety: str,
start_year: int,
end_year: int,
all_row_includes_end_year: bool = False,
) -> SimulationUnit:
"""Read sowing.csv and build the SimulationUnit for one site × variety.
Expected CSV columns (case-insensitive):
site, crop, variety, sowingDOY, year
(+ optional treatment1, treatment2, ... columns for fungicide dates)
Args:
sowing_file: Path to the sowing CSV.
site: Site identifier to filter on (e.g., "indiana").
variety: Variety identifier to filter on (e.g., "Generic").
start_year: First simulation year (inclusive).
end_year: Last simulation year (inclusive).
all_row_includes_end_year: How far an ``"All"`` sowing row is applied.
When ``False`` (default) the ``"All"`` row is applied to the half-open
range ``[start_year, end_year - 1]`` and therefore omits the final
year; when ``True`` it is applied through ``end_year`` inclusive. This
setting affects only ``"All"`` rows; explicit per-year rows are always
honoured for every year, so it has no effect when sowing is specified
per year (for example the CSV written by
:func:`simpest.models.simplace.build_management`).
Returns:
Populated SimulationUnit (year_sowing_doy, fungicide_treatment_schedule).
"""
sim = SimulationUnit()
path = Path(sowing_file)
# print(f'sowing_file: {sowing_file}')
with path.open(newline="", encoding="utf-8-sig") as fh:
reader = csv.reader(fh)
raw_headers = next(reader, [])
headers = [h.strip().strip('"').lower().replace(" ", "").replace("_", "") for h in raw_headers]
header_map = {h: i for i, h in enumerate(headers)}
site_idx = _get_col(header_map, "site")
crop_idx = _get_col(header_map, "crop")
variety_idx = _get_col(header_map, "variety", "cultivar", "cv")
sowing_idx = _get_col(header_map, "sowingdoy", "sowingday", "doy")
year_idx = _get_col(header_map, "year", "yr")
# Collect treatment columns (any column starting with "treatment")
fung_indices = [i for i, h in enumerate(headers) if h.startswith("treatment")]
# Store "All" row separately; per-year rows in a dict
all_row: Optional[Tuple[int, List[int]]] = None
per_year: Dict[int, Tuple[int, List[int]]] = {}
for row in reader:
if not row or all(c.strip() == "" for c in row):
continue
if site_idx >= 0:
row_site = row[site_idx].strip().strip('"').lower()
if row_site != site.lower():
continue
if variety_idx >= 0:
row_var = row[variety_idx].strip().strip('"').lower()
if row_var != variety.lower():
continue
if crop_idx >= 0 and sim.crop == "":
sim.crop = row[crop_idx].strip().strip('"')
if variety_idx >= 0 and sim.variety == "":
sim.variety = row[variety_idx].strip().strip('"')
try:
sow_doy = int(row[sowing_idx].strip())
except (ValueError, IndexError):
continue
fung_doys = []
for fi in fung_indices:
if fi < len(row):
try:
d = int(row[fi].strip())
if d > 0:
fung_doys.append(d)
except ValueError:
pass
fung_doys = sorted(set(fung_doys))
year_cell = row[year_idx].strip().strip('"') if year_idx >= 0 else ""
print(f'Parsed row: site={row_site}, year={year_cell}, variety={row_var}, sow_doy={sow_doy}, fung_doys={fung_doys}')
if year_cell.lower() == "all":
all_row = (sow_doy, fung_doys)
else:
try:
y = int(year_cell)
per_year[y] = (sow_doy, fung_doys)
except ValueError:
pass
# Apply the "All" row to every year in range. By default the final year is
# excluded (half-open range); see ``all_row_includes_end_year``.
all_end = end_year + 1 if all_row_includes_end_year else end_year
if all_row is not None:
for y in range(start_year, all_end):
_apply_row(sim, y, all_row[0], all_row[1])
# Override with per-year entries
for y, (sow_doy, fung_doys) in per_year.items():
_apply_row(sim, y, sow_doy, fung_doys)
return sim
fr_runner
¶
Main simulation runner for the FraNchEstYN model.
The runner orchestrates a full simulation: it loads parameters, weather, the sowing schedule, reference observations, and any external crop-model series, then drives the hourly–daily model loop across the configured years. Each simulated day couples the crop, disease, and fungicide sub-models and records a complete output bundle. The runner also exposes the calibration objective used to score a run against reference data.
Examples:
1 2 3 4 | |
FranchestynRunner
¶
End-to-end FraNchEstYN simulation runner.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
weather_dir |
str |
Directory containing weather data. |
required |
param_file |
str |
Path to parameter file. |
required |
sowing_file |
str |
Path to sowing or management CSV. |
required |
ref_dir |
str |
Directory or file path for reference data. |
required |
crop_model_dir |
Optional[str] |
Directory or file path for external
crop model data. If |
required |
site |
str |
Site identifier. |
required |
variety |
str |
Variety identifier. |
required |
disease |
str |
Disease column/key name used in reference data. |
required |
start_year |
int |
First simulation year (inclusive). |
required |
end_year |
int |
Last simulation year (inclusive). |
required |
weather_time_step |
str |
Weather frequency, |
'daily' |
calibration_variable |
str |
Calibration target scope: |
'all' |
is_calibration |
bool |
Whether the run is used for calibration. |
False |
latitude |
float |
Site latitude in decimal degrees. |
0.0 |
crop_type |
Optional[str] |
Crop type for modular parameter loading. |
None |
crop_param_file |
Optional[str] |
Crop parameter JSON path. |
None |
disease_param_file |
Optional[str] |
Disease parameter JSON path. |
None |
disease_type |
Optional[str] |
Disease type key for modular loading. Also
acts as the on/off switch for the disease model: when |
None |
fungicide_param_file |
Optional[str] |
Fungicide parameter JSON path. |
None |
fungicide_type |
Optional[str] |
Fungicide type key for modular loading.
Also the on/off switch for the fungicide model: when |
None |
use_gdd |
bool |
How crop cycle completion is derived. Defaults to
|
False |
use_prev_day_alignment |
bool |
Day alignment used by the calibration
objective. Defaults to |
True |
all_row_includes_end_year |
bool |
How far an |
False |
Source code in simpest/models/fr_runner.py
class FranchestynRunner:
"""
End-to-end FraNchEstYN simulation runner.
Args:
weather_dir (str): Directory containing weather data.
param_file (str): Path to parameter file.
sowing_file (str): Path to sowing or management CSV.
ref_dir (str): Directory or file path for reference data.
crop_model_dir (Optional[str]): Directory or file path for external
crop model data. If ``None``, the internal crop model is used.
site (str): Site identifier.
variety (str): Variety identifier.
disease (str): Disease column/key name used in reference data.
start_year (int): First simulation year (inclusive).
end_year (int): Last simulation year (inclusive).
weather_time_step (str): Weather frequency, ``daily`` or ``hourly``.
calibration_variable (str): Calibration target scope: ``crop``,
``disease``, or ``all``.
is_calibration (bool): Whether the run is used for calibration.
latitude (float): Site latitude in decimal degrees.
crop_type (Optional[str]): Crop type for modular parameter loading.
crop_param_file (Optional[str]): Crop parameter JSON path.
disease_param_file (Optional[str]): Disease parameter JSON path.
disease_type (Optional[str]): Disease type key for modular loading. Also
acts as the on/off switch for the disease model: when ``None`` the
SEIR disease step (hourly and daily) is skipped entirely, so disease
severity stays 0 and there is no disease impact (a crop-only run).
Provide a ``disease_type`` to enable the epidemiological model.
fungicide_param_file (Optional[str]): Fungicide parameter JSON path.
fungicide_type (Optional[str]): Fungicide type key for modular loading.
Also the on/off switch for the fungicide model: when ``None`` the
fungicide step is skipped. The step is a no-op unless a treatment has
been scheduled, so it only needs to be enabled when fungicide
treatments are present (otherwise they would be ignored).
use_gdd (bool): How crop cycle completion is derived. Defaults to
``False`` (calendar-day interpolation of the external crop-model
series); set ``True`` to use the thermal-time (GDD) based cycle
percentage, which ties phenology more directly to accumulated heat.
use_prev_day_alignment (bool): Day alignment used by the calibration
objective. Defaults to ``True``, which compares the previous day's
simulated output to the current day's reference observation
(sim[d-1] vs ref[d]). Set ``False`` for same-day alignment
(sim[d] vs ref[d]), which matches the daily output table.
all_row_includes_end_year (bool): How far an ``"All"`` sowing row is
applied. Defaults to ``False`` (the final year is omitted); set
``True`` to apply it through ``end_year`` inclusive. Has no effect
for per-year sowing CSVs.
"""
def __init__(
self,
weather_dir: str,
param_file: str,
sowing_file: str,
ref_dir: str,
crop_model_dir: Optional[str],
site: str,
variety: str,
disease: str,
start_year: int,
end_year: int,
weather_time_step: str = "daily",
calibration_variable: str = "all",
is_calibration: bool = False,
latitude: float = 0.0,
crop_type: Optional[str] = None,
crop_param_file: Optional[str] = None,
disease_param_file: Optional[str] = None,
disease_type: Optional[str] = None,
fungicide_param_file: Optional[str] = None,
fungicide_type: Optional[str] = None,
use_gdd: bool = False,
use_prev_day_alignment: bool = True,
all_row_includes_end_year: bool = False,
) -> None:
self.weather_dir = weather_dir
self.param_file = param_file
self.sowing_file = sowing_file
self.ref_dir = ref_dir
self.crop_model_dir = crop_model_dir
self.site = site
self.variety = variety
self.disease = disease
self.start_year = start_year
self.end_year = end_year
self.weather_time_step = weather_time_step.lower()
self.calibration_variable = calibration_variable
self.is_calibration = is_calibration
self.latitude = latitude
self.disease_type = disease_type
self.fungicide_type = fungicide_type
self.use_gdd = use_gdd
self.use_prev_day_alignment = use_prev_day_alignment
self.all_row_includes_end_year = all_row_includes_end_year
# Read parameter definitions (bounds, defaults)
# Three modular loading scenarios:
# 1. All three files provided → modular (crop + optional disease + optional fungicide)
# 2. Only crop_type provided → legacy multi-crop JSON
# 3. Only param_file → legacy CSV
self.name_param: Dict[str, Parameter] = {}
if crop_param_file and crop_type:
# Modular loading: always load crop
self.name_param.update(read_crop_parameters(crop_param_file, crop_type))
# Optionally load disease if both file and type provided
if disease_param_file and disease_type:
self.name_param.update(read_disease_parameters(disease_param_file, disease_type))
# Optionally load fungicide if both file and type provided
if fungicide_param_file and fungicide_type:
self.name_param.update(read_fungicide_parameters(fungicide_param_file, fungicide_type))
elif crop_type:
# Legacy multi-crop JSON (parameters_by_crop.json)
self.name_param = read_by_crop(param_file, crop_type)
else:
# Legacy CSV loader
self.name_param = param_read(param_file, calibration_variable)
# Read calibrated values (override file, empty if not present)
self.param_out_calibration: Dict[str, float] = {}
# Read sowing + reference data
self.sim_unit: SimulationUnit = read_sowing(
sowing_file, site, variety, start_year, end_year,
all_row_includes_end_year=all_row_includes_end_year,
)
self.sim_unit = read_reference(
ref_dir, sowing_file, site, variety,
start_year, end_year, self.sim_unit, disease
)
# Read external crop model data (None → internal model used)
self.crop_model_data: Optional[CropModelData] = None
if crop_model_dir:
self.crop_model_data = read_crop_model_data(crop_model_dir, use_gdd=self.use_gdd)
# -----------------------------------------------------------------------
# Public API
# -----------------------------------------------------------------------
def run(
self,
param_values: Optional[Dict[str, float]] = None,
) -> Dict[datetime, Outputs]:
"""
Run the model for all configured years and return daily outputs.
Args:
param_values (Optional[Dict[str, float]]): Optional override values
keyed by ``class_ParamName``.
Returns:
Dict[datetime, Outputs]: Mapping of end-of-day timestamp to model
outputs.
"""
parameters = self._build_parameters(param_values or {})
weather_data = self._load_weather()
date_outputs: Dict[datetime, Outputs] = {}
# Per-hour accumulators aggregated into the daily input each day
hourly_temps: List[float] = []
hourly_rads: List[float] = []
hourly_precip: List[float] = []
hourly_rhs: List[float] = []
hourly_lw: List[float] = []
output = Outputs()
output_t1 = Outputs()
disease_model = DiseaseModel()
is_planted = False
last_treatment_date = datetime(1900, 1, 1)
for hour_dt in sorted(weather_data.keys()):
year = hour_dt.year
if year < self.start_year or year > self.end_year:
output = Outputs()
output_t1 = Outputs()
continue
hourly_rec = weather_data[hour_dt]
hourly_rec.date = hour_dt
hourly_rec.latitude = self.latitude
# Sowing: reset on the correct DOY
sow_doy = self.sim_unit.year_sowing_doy.get(year)
if sow_doy and hour_dt.timetuple().tm_yday == sow_doy:
is_planted = True
output = Outputs()
output_t1 = Outputs()
disease_model = DiseaseModel()
output_t1.crop.growing_season = year
if is_planted:
# Track fungicide treatments
sched = self.sim_unit.fungicide_treatment_schedule
treatment_dates = {t.date() for t in sched.treatments}
if hour_dt.date() in treatment_dates:
last_treatment_date = datetime.combine(
hour_dt.date(), datetime.min.time()
)
hourly_rec.date_treatment_last = last_treatment_date
# Accumulate hourly observations
hourly_temps.append(hourly_rec.air_temperature)
hourly_rads.append(hourly_rec.rad)
hourly_precip.append(hourly_rec.precipitation)
hourly_rhs.append(hourly_rec.relative_humidity)
hourly_lw.append(hourly_rec.leaf_wetness)
# Skip the hourly disease step when disease is disabled, and also
# for crop-only calibration runs where disease plays no role.
if self.disease_type and (self.calibration_variable != "crop" or not self.is_calibration):
disease_model.run_hourly(hourly_rec, parameters, output, output_t1)
if hour_dt.hour == 23:
# ----------------------------------------------------------
# End of day: swap, build daily input, run daily models
# ----------------------------------------------------------
output = output_t1
# Fresh daily output with season-persistent fields carried over
output_t1 = Outputs()
output_t1.crop.growing_season = output.crop.growing_season
output_t1.crop.f_int_peak = output.crop.f_int_peak
output_t1.disease.is_primary_inoculum_started = (
output.disease.is_primary_inoculum_started
)
output_t1.disease.first_seasonal_infection = (
output.disease.first_seasonal_infection
)
output_t1.disease.cycle_percentage_first_infection = (
output.disease.cycle_percentage_first_infection
)
# Assemble daily input from hourly lists
input_daily = InputsDaily()
input_daily.tmax = max(hourly_temps)
input_daily.tmin = min(hourly_temps)
input_daily.rad = sum(hourly_rads)
input_daily.precipitation = sum(hourly_precip)
input_daily.rhx = max(hourly_rhs)
input_daily.rhn = min(hourly_rhs)
input_daily.leaf_wetness = sum(hourly_lw)
# Date of this day = hour 23 minus 23 hours
input_daily.date = hour_dt - timedelta(hours=23)
input_daily.date_treatment_last = hourly_rec.date_treatment_last
input_daily.crop_model_data = self.crop_model_data
# Run the daily sub-models. The crop step always runs; the
# fungicide and disease steps are optional and gated on their
# respective *_type (None skips the step, giving a crop-only
# or no-fungicide run).
crop_run(input_daily, parameters, output, output_t1)
if self.fungicide_type:
fungicide_run(input_daily, parameters, output, output_t1)
if self.disease_type:
disease_model.run_daily(input_daily, parameters, output, output_t1)
# Store daily output
output_t1.inputs_daily = input_daily
date_outputs[hour_dt] = output_t1
# Clear hourly accumulators
hourly_temps.clear()
hourly_rads.clear()
hourly_precip.clear()
hourly_rhs.clear()
hourly_lw.clear()
# Safety stop or maturity check
if (output_t1.crop.cycle_completion_percentage >= 100
or output_t1.crop.day_after_sowing >= _MAX_DAS):
is_planted = False
else:
is_planted = False
output = Outputs()
output_t1 = Outputs()
return date_outputs
def compute_rmse(
self,
date_outputs: Dict[datetime, Outputs],
include_crop: bool = True,
include_disease: bool = True,
) -> float:
"""Compute the root-mean-square error against reference data.
This is the calibration objective. For each simulated day that has a
matching reference observation, it accumulates the squared, scaled
errors of the selected variables — above-ground biomass, attainable and
actual yield, light interception, and disease severity — and returns the
root mean of those per-day error sums. Two penalty rules sharpen the
objective: a zero simulated yield at maturity where the reference is
positive, and a near-zero simulated light interception during the
growing season, are both heavily penalised.
Args:
date_outputs (Dict[datetime, Outputs]): Daily outputs produced by
:meth:`run`.
include_crop (bool): Include the crop-related error terms (biomass,
attainable yield, light interception).
include_disease (bool): Include the disease-related error terms
(severity and actual yield).
Returns:
float: The root-mean-square error, rounded to three decimals; ``0.0``
when no reference observations overlap the simulated days.
"""
errors: List[float] = []
ref = self.sim_unit.reference_data
from datetime import date as date_type
for hour_dt, out in date_outputs.items():
if hour_dt.hour != 23:
continue
# Reference keys are datetime.date; convert for lookup.
sim_day: date_type = hour_dt.date()
# Day alignment between simulated output and reference data.
# Previous-day alignment (default): compare this simulated day to
# the next day's reference, i.e. sim[d] vs ref[d+1], which is
# equivalent to sim[d-1] vs ref[d].
# Same-day alignment: sim[d] vs ref[d], matching the daily table.
ref_key: date_type = (
sim_day + timedelta(days=1)
if self.use_prev_day_alignment else sim_day
)
total_err = 0.0
has_ref = False
# Infer crop state for the penalty multipliers below
is_matured = out.crop.cycle_completion_percentage >= 100.0
is_planted = out.crop.day_after_sowing > 0 and not is_matured
if include_crop:
agb_err = 0.0
if ref_key in ref.date_agb:
sim_agb = out.crop.agb_attainable
agb_err = ((ref.date_agb[ref_key] - sim_agb) / 200.0) ** 2
has_ref = True
yield_err = 0.0
if ref_key in ref.date_yield_attainable:
sim_y = out.crop.yield_attainable
yield_ref = ref.date_yield_attainable[ref_key]
yield_err = ((yield_ref - sim_y) / 100.0) ** 2
# Penalty: heavily penalise zero yield at maturity when ref > 0
if sim_y == 0.0 and yield_ref > 0.0 and is_matured:
yield_err *= 1000.0
has_ref = True
fint_err = 0.0
if ref_key in ref.date_fint:
sim_fi = out.crop.light_interception_attainable * 100.0
fint_err = (ref.date_fint[ref_key] * 100.0 - sim_fi) ** 2
# Penalty: heavily penalise near-zero interception during the season
if sim_fi < 0.1 and is_planted:
fint_err *= 1000.0
has_ref = True
total_err += agb_err + yield_err + fint_err
if include_disease:
dis_err = 0.0
by_date = ref.disease_date_disease_sev.get(self.disease, {})
if ref_key in by_date:
sim_ds = out.disease.disease_severity * 100.0
dis_err = (by_date[ref_key] - sim_ds) ** 2
has_ref = True
yield_err2 = 0.0
if ref_key in ref.date_yield_actual:
sim_ya = out.crop.yield_actual
yield_ref2 = ref.date_yield_actual[ref_key]
yield_err2 = ((yield_ref2 - sim_ya) / 100.0) ** 2
# Penalty: heavily penalise zero actual yield at maturity when ref > 0
if sim_ya == 0.0 and yield_ref2 > 0.0 and is_matured:
yield_err2 *= 1000.0
has_ref = True
total_err += dis_err + yield_err2
if has_ref:
errors.append(total_err)
if not errors:
return 0.0
import math
return round(math.sqrt(sum(errors) / len(errors)), 3)
# -----------------------------------------------------------------------
# Private helpers
# -----------------------------------------------------------------------
def _build_parameters(
self, param_values: Dict[str, float]
) -> Parameters:
"""Assemble a Parameters object from CSV defaults + overrides."""
parameters = Parameters()
for key, p in self.name_param.items():
# key = "class_ParamName", e.g. "crop_TbaseCrop"
parts = key.split("_", 1)
if len(parts) != 2:
continue
param_class, param_name = parts
# Use override value if supplied, else CSV default
if key in param_values:
value = param_values[key]
elif key in self.param_out_calibration:
value = self.param_out_calibration[key]
else:
value = p.value_bool if p.is_boolean else p.value
_set_param(parameters, param_class, param_name, value)
return parameters
def _load_weather(self) -> Dict[datetime, InputsHourly]:
"""Load weather data for the configured site and time step.
If ``self.weather_dir`` is itself an existing CSV file, it is used
directly (i.e. the caller already resolved the full path).
"""
wd = self.weather_dir
if os.path.isfile(wd):
weather_file = wd
else:
weather_file = os.path.join(
wd, self.weather_time_step, f"{self.site}.csv"
)
if self.weather_time_step == "hourly":
return read_hourly(
weather_file, self.start_year, self.end_year,
site=self.site, latitude=self.latitude
)
else:
return read_daily(
weather_file, self.start_year, self.end_year,
latitude=self.latitude
)
compute_rmse(self, date_outputs, include_crop=True, include_disease=True)
¶
Compute the root-mean-square error against reference data.
This is the calibration objective. For each simulated day that has a matching reference observation, it accumulates the squared, scaled errors of the selected variables — above-ground biomass, attainable and actual yield, light interception, and disease severity — and returns the root mean of those per-day error sums. Two penalty rules sharpen the objective: a zero simulated yield at maturity where the reference is positive, and a near-zero simulated light interception during the growing season, are both heavily penalised.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
date_outputs |
Dict[datetime, Outputs] |
Daily outputs produced by
:meth: |
required |
include_crop |
bool |
Include the crop-related error terms (biomass, attainable yield, light interception). |
True |
include_disease |
bool |
Include the disease-related error terms (severity and actual yield). |
True |
Returns:
| Type | Description |
|---|---|
float |
The root-mean-square error, rounded to three decimals; |
Source code in simpest/models/fr_runner.py
def compute_rmse(
self,
date_outputs: Dict[datetime, Outputs],
include_crop: bool = True,
include_disease: bool = True,
) -> float:
"""Compute the root-mean-square error against reference data.
This is the calibration objective. For each simulated day that has a
matching reference observation, it accumulates the squared, scaled
errors of the selected variables — above-ground biomass, attainable and
actual yield, light interception, and disease severity — and returns the
root mean of those per-day error sums. Two penalty rules sharpen the
objective: a zero simulated yield at maturity where the reference is
positive, and a near-zero simulated light interception during the
growing season, are both heavily penalised.
Args:
date_outputs (Dict[datetime, Outputs]): Daily outputs produced by
:meth:`run`.
include_crop (bool): Include the crop-related error terms (biomass,
attainable yield, light interception).
include_disease (bool): Include the disease-related error terms
(severity and actual yield).
Returns:
float: The root-mean-square error, rounded to three decimals; ``0.0``
when no reference observations overlap the simulated days.
"""
errors: List[float] = []
ref = self.sim_unit.reference_data
from datetime import date as date_type
for hour_dt, out in date_outputs.items():
if hour_dt.hour != 23:
continue
# Reference keys are datetime.date; convert for lookup.
sim_day: date_type = hour_dt.date()
# Day alignment between simulated output and reference data.
# Previous-day alignment (default): compare this simulated day to
# the next day's reference, i.e. sim[d] vs ref[d+1], which is
# equivalent to sim[d-1] vs ref[d].
# Same-day alignment: sim[d] vs ref[d], matching the daily table.
ref_key: date_type = (
sim_day + timedelta(days=1)
if self.use_prev_day_alignment else sim_day
)
total_err = 0.0
has_ref = False
# Infer crop state for the penalty multipliers below
is_matured = out.crop.cycle_completion_percentage >= 100.0
is_planted = out.crop.day_after_sowing > 0 and not is_matured
if include_crop:
agb_err = 0.0
if ref_key in ref.date_agb:
sim_agb = out.crop.agb_attainable
agb_err = ((ref.date_agb[ref_key] - sim_agb) / 200.0) ** 2
has_ref = True
yield_err = 0.0
if ref_key in ref.date_yield_attainable:
sim_y = out.crop.yield_attainable
yield_ref = ref.date_yield_attainable[ref_key]
yield_err = ((yield_ref - sim_y) / 100.0) ** 2
# Penalty: heavily penalise zero yield at maturity when ref > 0
if sim_y == 0.0 and yield_ref > 0.0 and is_matured:
yield_err *= 1000.0
has_ref = True
fint_err = 0.0
if ref_key in ref.date_fint:
sim_fi = out.crop.light_interception_attainable * 100.0
fint_err = (ref.date_fint[ref_key] * 100.0 - sim_fi) ** 2
# Penalty: heavily penalise near-zero interception during the season
if sim_fi < 0.1 and is_planted:
fint_err *= 1000.0
has_ref = True
total_err += agb_err + yield_err + fint_err
if include_disease:
dis_err = 0.0
by_date = ref.disease_date_disease_sev.get(self.disease, {})
if ref_key in by_date:
sim_ds = out.disease.disease_severity * 100.0
dis_err = (by_date[ref_key] - sim_ds) ** 2
has_ref = True
yield_err2 = 0.0
if ref_key in ref.date_yield_actual:
sim_ya = out.crop.yield_actual
yield_ref2 = ref.date_yield_actual[ref_key]
yield_err2 = ((yield_ref2 - sim_ya) / 100.0) ** 2
# Penalty: heavily penalise zero actual yield at maturity when ref > 0
if sim_ya == 0.0 and yield_ref2 > 0.0 and is_matured:
yield_err2 *= 1000.0
has_ref = True
total_err += dis_err + yield_err2
if has_ref:
errors.append(total_err)
if not errors:
return 0.0
import math
return round(math.sqrt(sum(errors) / len(errors)), 3)
run(self, param_values=None)
¶
Run the model for all configured years and return daily outputs.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
param_values |
Optional[Dict[str, float]] |
Optional override values
keyed by |
None |
Returns:
| Type | Description |
|---|---|
Dict[datetime, Outputs] |
Mapping of end-of-day timestamp to model outputs. |
Source code in simpest/models/fr_runner.py
def run(
self,
param_values: Optional[Dict[str, float]] = None,
) -> Dict[datetime, Outputs]:
"""
Run the model for all configured years and return daily outputs.
Args:
param_values (Optional[Dict[str, float]]): Optional override values
keyed by ``class_ParamName``.
Returns:
Dict[datetime, Outputs]: Mapping of end-of-day timestamp to model
outputs.
"""
parameters = self._build_parameters(param_values or {})
weather_data = self._load_weather()
date_outputs: Dict[datetime, Outputs] = {}
# Per-hour accumulators aggregated into the daily input each day
hourly_temps: List[float] = []
hourly_rads: List[float] = []
hourly_precip: List[float] = []
hourly_rhs: List[float] = []
hourly_lw: List[float] = []
output = Outputs()
output_t1 = Outputs()
disease_model = DiseaseModel()
is_planted = False
last_treatment_date = datetime(1900, 1, 1)
for hour_dt in sorted(weather_data.keys()):
year = hour_dt.year
if year < self.start_year or year > self.end_year:
output = Outputs()
output_t1 = Outputs()
continue
hourly_rec = weather_data[hour_dt]
hourly_rec.date = hour_dt
hourly_rec.latitude = self.latitude
# Sowing: reset on the correct DOY
sow_doy = self.sim_unit.year_sowing_doy.get(year)
if sow_doy and hour_dt.timetuple().tm_yday == sow_doy:
is_planted = True
output = Outputs()
output_t1 = Outputs()
disease_model = DiseaseModel()
output_t1.crop.growing_season = year
if is_planted:
# Track fungicide treatments
sched = self.sim_unit.fungicide_treatment_schedule
treatment_dates = {t.date() for t in sched.treatments}
if hour_dt.date() in treatment_dates:
last_treatment_date = datetime.combine(
hour_dt.date(), datetime.min.time()
)
hourly_rec.date_treatment_last = last_treatment_date
# Accumulate hourly observations
hourly_temps.append(hourly_rec.air_temperature)
hourly_rads.append(hourly_rec.rad)
hourly_precip.append(hourly_rec.precipitation)
hourly_rhs.append(hourly_rec.relative_humidity)
hourly_lw.append(hourly_rec.leaf_wetness)
# Skip the hourly disease step when disease is disabled, and also
# for crop-only calibration runs where disease plays no role.
if self.disease_type and (self.calibration_variable != "crop" or not self.is_calibration):
disease_model.run_hourly(hourly_rec, parameters, output, output_t1)
if hour_dt.hour == 23:
# ----------------------------------------------------------
# End of day: swap, build daily input, run daily models
# ----------------------------------------------------------
output = output_t1
# Fresh daily output with season-persistent fields carried over
output_t1 = Outputs()
output_t1.crop.growing_season = output.crop.growing_season
output_t1.crop.f_int_peak = output.crop.f_int_peak
output_t1.disease.is_primary_inoculum_started = (
output.disease.is_primary_inoculum_started
)
output_t1.disease.first_seasonal_infection = (
output.disease.first_seasonal_infection
)
output_t1.disease.cycle_percentage_first_infection = (
output.disease.cycle_percentage_first_infection
)
# Assemble daily input from hourly lists
input_daily = InputsDaily()
input_daily.tmax = max(hourly_temps)
input_daily.tmin = min(hourly_temps)
input_daily.rad = sum(hourly_rads)
input_daily.precipitation = sum(hourly_precip)
input_daily.rhx = max(hourly_rhs)
input_daily.rhn = min(hourly_rhs)
input_daily.leaf_wetness = sum(hourly_lw)
# Date of this day = hour 23 minus 23 hours
input_daily.date = hour_dt - timedelta(hours=23)
input_daily.date_treatment_last = hourly_rec.date_treatment_last
input_daily.crop_model_data = self.crop_model_data
# Run the daily sub-models. The crop step always runs; the
# fungicide and disease steps are optional and gated on their
# respective *_type (None skips the step, giving a crop-only
# or no-fungicide run).
crop_run(input_daily, parameters, output, output_t1)
if self.fungicide_type:
fungicide_run(input_daily, parameters, output, output_t1)
if self.disease_type:
disease_model.run_daily(input_daily, parameters, output, output_t1)
# Store daily output
output_t1.inputs_daily = input_daily
date_outputs[hour_dt] = output_t1
# Clear hourly accumulators
hourly_temps.clear()
hourly_rads.clear()
hourly_precip.clear()
hourly_rhs.clear()
hourly_lw.clear()
# Safety stop or maturity check
if (output_t1.crop.cycle_completion_percentage >= 100
or output_t1.crop.day_after_sowing >= _MAX_DAS):
is_planted = False
else:
is_planted = False
output = Outputs()
output_t1 = Outputs()
return date_outputs
fr_utilities
¶
Biophysical helper functions shared across the model.
This module provides small, stateless functions that encode reusable biophysical relationships, such as the cardinal-temperature response used by both the crop and disease sub-models and the rainfall-driven spore detachment index used by the disease model.
rain_detachment(rainfall, rain50, f_int)
¶
Rain-driven spore detachment index.
Returns a dimensionless index in [0, 1] that saturates as rainfall
increases relative to the canopy-scaled half-saturation term:
\[
\text{detachment} = \frac{\text{rainfall}}
{\text{rain50} \cdot f_{int} + \text{rainfall}}
\]
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
rainfall |
float |
Precipitation over the time step (mm). |
required |
rain50 |
float |
Half-saturation parameter (mm); the rainfall that yields
an index of about 0.5 at full canopy ( |
required |
f_int |
float |
Light interception fraction in |
required |
Returns:
| Type | Description |
|---|---|
float |
Detachment index in |
Source code in simpest/models/fr_utilities.py
def rain_detachment(rainfall: float, rain50: float, f_int: float) -> float:
"""Rain-driven spore detachment index.
Returns a dimensionless index in ``[0, 1]`` that saturates as rainfall
increases relative to the canopy-scaled half-saturation term:
$$
\\text{detachment} = \\frac{\\text{rainfall}}
{\\text{rain50} \\cdot f_{int} + \\text{rainfall}}
$$
Args:
rainfall (float): Precipitation over the time step (mm).
rain50 (float): Half-saturation parameter (mm); the rainfall that yields
an index of about 0.5 at full canopy (``f_int = 1``).
f_int (float): Light interception fraction in ``[0, 1]``.
Returns:
float: Detachment index in ``[0, 1]``; ``0.0`` when rainfall is zero or
the denominator is zero.
"""
denominator = (rain50 * f_int) + rainfall
if denominator == 0.0:
return 0.0
return rainfall / denominator
t_response(t_ave, t_base, t_opt, t_max)
¶
Beta-shaped cardinal-temperature response function.
Returns a dimensionless growth or development efficiency in the range
[0, 1] that peaks at the optimum temperature and falls to zero at the
base and maximum temperatures. For t_base < t_ave < t_max the response is
\[
f = \frac{t_{max} - t_{ave}}{t_{max} - t_{opt}}
\left(\frac{t_{ave} - t_{base}}{t_{opt} - t_{base}}\right)
^{\frac{t_{opt} - t_{base}}{t_{max} - t_{opt}}}
\]
and 0 outside the interval [t_base, t_max].
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
t_ave |
float |
Average temperature for the step (°C). |
required |
t_base |
float |
Base (minimum) cardinal temperature (°C). |
required |
t_opt |
float |
Optimum cardinal temperature (°C). |
required |
t_max |
float |
Maximum cardinal temperature (°C). |
required |
Returns:
| Type | Description |
|---|---|
float |
Temperature response factor in |
Source code in simpest/models/fr_utilities.py
def t_response(t_ave: float, t_base: float, t_opt: float, t_max: float) -> float:
"""Beta-shaped cardinal-temperature response function.
Returns a dimensionless growth or development efficiency in the range
``[0, 1]`` that peaks at the optimum temperature and falls to zero at the
base and maximum temperatures. For ``t_base < t_ave < t_max`` the response is
$$
f = \\frac{t_{max} - t_{ave}}{t_{max} - t_{opt}}
\\left(\\frac{t_{ave} - t_{base}}{t_{opt} - t_{base}}\\right)
^{\\frac{t_{opt} - t_{base}}{t_{max} - t_{opt}}}
$$
and ``0`` outside the interval ``[t_base, t_max]``.
Args:
t_ave (float): Average temperature for the step (°C).
t_base (float): Base (minimum) cardinal temperature (°C).
t_opt (float): Optimum cardinal temperature (°C).
t_max (float): Maximum cardinal temperature (°C).
Returns:
float: Temperature response factor in ``[0, 1]``.
"""
if t_ave <= t_base or t_ave >= t_max:
return 0.0
first_term = (t_max - t_ave) / (t_max - t_opt)
second_term = (t_ave - t_base) / (t_opt - t_base)
exponent = (t_opt - t_base) / (t_max - t_opt)
return first_term * math.pow(second_term, exponent)
fr_weather_reader
¶
Weather readers that build an hourly forcing series for the model.
The disease sub-model operates on an hourly time step, so these readers ingest either daily or hourly weather CSV files and return a complete hourly series. When the input is daily, a physically based diurnal synthesis fills in the 24 hourly records for each day.
Key behaviour:
- Header parsing. Column names are matched case-insensitively against a set of common aliases, so files from different sources are accepted without renaming.
- Date columns. A single
Date/Datetimecolumn or separateYear/Month/Day(plusHourfor hourly files) are both accepted. - Radiation. Measured radiation is used when present; otherwise daily global solar radiation is estimated with the Hargreaves–Samani relationship from the daily temperature range and distributed over the day by clear-sky (extraterrestrial) fractions.
- Humidity. Relative humidity is reconstructed from daily
RHx/RHnextrema with a cosine diurnal curve, or estimated from a dew-point relation when extrema are absent.
read_daily(file, start_year, end_year, latitude=0.0)
¶
Read a daily weather CSV and return synthesized hourly records.
Supports year/month/day columns OR a single date/datetime column. Radiation is estimated when absent (Hargreaves-Samani) if latitude is known.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
file |
str | Path |
Path to the daily weather CSV. |
required |
start_year |
int |
First year to include (inclusive). |
required |
end_year |
int |
Last year to include (inclusive). |
required |
latitude |
float |
Fallback latitude in decimal degrees when not in CSV. |
0.0 |
Returns:
| Type | Description |
|---|---|
Dict[datetime, InputsHourly] |
Dict keyed by |
Source code in simpest/models/fr_weather_reader.py
def read_daily(
file: str | Path,
start_year: int,
end_year: int,
latitude: float = 0.0,
) -> Dict[datetime, InputsHourly]:
"""Read a daily weather CSV and return synthesized hourly records.
Supports year/month/day columns OR a single date/datetime column.
Radiation is estimated when absent (Hargreaves-Samani) if latitude is known.
Args:
file: Path to the daily weather CSV.
start_year: First year to include (inclusive).
end_year: Last year to include (inclusive).
latitude: Fallback latitude in decimal degrees when not in CSV.
Returns:
Dict keyed by ``datetime(year, month, day, hour)`` → ``InputsHourly``.
"""
result: Dict[datetime, InputsHourly] = {}
path = Path(file)
with path.open(newline="", encoding="utf-8-sig") as fh:
sample = fh.read(2048)
fh.seek(0)
delimiter = "\t" if "\t" in sample.split("\n")[0] else ","
reader = csv.reader(fh, delimiter=delimiter)
headers = [_clean(h) for h in next(reader)]
col = {name: i for i, name in enumerate(headers)}
# Date columns
date_idx = _idx(col, ["date", "datetime", "timestamp"], optional=True)
year_idx = _idx(col, ["year"], optional=True)
month_idx = _idx(col, ["month"], optional=True)
day_idx = _idx(col, ["day"], optional=True)
has_ymd = (year_idx >= 0 and month_idx >= 0 and day_idx >= 0)
has_date = date_idx >= 0
if not has_ymd and not has_date:
raise ValueError(
"Weather file must contain (year, month, day) OR a single date/datetime column."
)
# Meteorological columns
rad_idx = _idx(col, ["rad", "radiation", "solar", "solarrad", "srad"], optional=True)
lat_idx = _idx(col, ["lat", "latitude", "sitelat", "phi"], optional=True)
tmax_idx = _idx(col, ["tmax", "t2mmax", "maxtemp", "tx"])
tmin_idx = _idx(col, ["tmin", "t2mmin", "mintemp", "tn"])
prec_idx = _idx(col, ["prec", "precip", "rain", "rainfall", "precipitation", "p"], optional=True)
rhx_idx = _idx(col, ["rhmax", "humiditymax", "relativehumiditymax", "hummax", "rhx"], optional=True)
rhn_idx = _idx(col, ["rhmin", "humiditymin", "relativehumiditymin", "hummin", "rhn"], optional=True)
if rad_idx < 0 and lat_idx < 0 and latitude == 0.0:
raise ValueError(
"Weather file must contain a radiation column OR a latitude column "
"(or pass latitude= to read_daily)."
)
for row in reader:
if not row or all(c.strip() == "" for c in row):
continue
# --- Date ---
try:
if has_ymd:
yr = int(row[year_idx])
mo = int(row[month_idx])
dy = int(row[day_idx])
day_date = date(yr, mo, dy)
else:
raw_date = row[date_idx].strip().strip('"')
day_date = datetime.fromisoformat(raw_date).date()
except (ValueError, IndexError):
continue
if not (start_year <= day_date.year <= end_year):
continue
# --- Temperatures ---
try:
tmax = _pf(row, tmax_idx)
tmin = _pf(row, tmin_idx)
except (ValueError, IndexError):
continue
# --- Radiation ---
rad = _pf(row, rad_idx) if rad_idx >= 0 else float("nan")
rad_missing = not math.isfinite(rad) or rad <= 0.0
# --- Latitude (for radiation estimation) ---
lat = _pf(row, lat_idx) if lat_idx >= 0 else latitude
if not math.isfinite(lat):
lat = latitude
if rad_missing:
if lat == 0.0:
continue # cannot estimate without latitude
rd = _day_length(day_date, lat, tmax, tmin)
rad = rd["gsr"]
if not math.isfinite(rad) or rad <= 0.0:
continue
# --- Build InputsDaily ---
prec = _pf(row, prec_idx) if prec_idx >= 0 else 0.0
if not math.isfinite(prec):
prec = 0.0
id_ = InputsDaily(
date=datetime(day_date.year, day_date.month, day_date.day),
tmax=tmax,
tmin=tmin,
rad=rad,
precipitation=prec,
latitude=lat,
)
id_.dew_point = _dew_point(id_.tmax, id_.tmin)
if rhx_idx >= 0:
v = _pf(row, rhx_idx)
if math.isfinite(v):
id_.rhx = v
if rhn_idx >= 0:
v = _pf(row, rhn_idx)
if math.isfinite(v):
id_.rhn = v
# Calculate leaf wetness from RH and precipitation
# High humidity (>80%) or significant precipitation (>0.5mm) indicate wet leaves
avg_rh = (id_.rhx + id_.rhn) / 2.0 if id_.rhx > 0 else 0.0
if prec >= 0.5 or avg_rh >= 80.0:
# Estimate wetness hours: scale from 60% RH (minimal) to 100% RH (24 hours)
if avg_rh >= 80.0:
id_.leaf_wetness = min(24.0, (avg_rh - 80.0) * 1.2) # 80% -> 0h, 100% -> 24h
else:
id_.leaf_wetness = 0.0
# Add hours if precipitation present
if prec >= 0.5:
id_.leaf_wetness = min(24.0, id_.leaf_wetness + prec * 2.0) # Add ~2h per mm rainfall
# --- Synthesize 24 hourly records ---
hourly = _estimate_hourly(id_, day_date, lat)
result.update(hourly)
return result
read_hourly(file, start_year, end_year, site='', latitude=0.0)
¶
Read an hourly weather CSV and return the hourly records.
Missing hourly radiation is filled per-day using the Hargreaves estimate distributed via clear-sky (ETR) fractions.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
file |
str | Path |
Path to the hourly weather CSV. |
required |
start_year |
int |
First year to include (inclusive). |
required |
end_year |
int |
Last year to include (inclusive). |
required |
site |
str |
Optional site filter (column 'site' or 'station'). |
'' |
latitude |
float |
Fallback latitude when not in the file. |
0.0 |
Returns:
| Type | Description |
|---|---|
Dict[datetime, InputsHourly] |
Dict keyed by |
Source code in simpest/models/fr_weather_reader.py
def read_hourly(
file: str | Path,
start_year: int,
end_year: int,
site: str = "",
latitude: float = 0.0,
) -> Dict[datetime, InputsHourly]:
"""Read an hourly weather CSV and return the hourly records.
Missing hourly radiation is filled per-day using the Hargreaves estimate
distributed via clear-sky (ETR) fractions.
Args:
file: Path to the hourly weather CSV.
start_year: First year to include (inclusive).
end_year: Last year to include (inclusive).
site: Optional site filter (column 'site' or 'station').
latitude: Fallback latitude when not in the file.
Returns:
Dict keyed by ``datetime(year, month, day, hour)`` → ``InputsHourly``.
"""
result: Dict[datetime, InputsHourly] = {}
daily_buckets: Dict[date, List[InputsHourly]] = {}
path = Path(file)
with path.open(newline="", encoding="utf-8-sig") as fh:
sample = fh.read(2048)
fh.seek(0)
delimiter = "\t" if "\t" in sample.split("\n")[0] else ","
reader = csv.reader(fh, delimiter=delimiter)
headers = [_clean(h) for h in next(reader)]
col = {name: i for i, name in enumerate(headers)}
year_idx = _idx(col, ["year"], optional=True)
month_idx = _idx(col, ["month", "mo"], optional=True)
day_idx = _idx(col, ["day", "dd", "dy"], optional=True)
hour_idx = _idx(col, ["hour", "hr", "h"], optional=True)
date_idx = _idx(col, ["date", "datetime", "timestamp"], optional=True)
has_ymdh = (year_idx >= 0 and month_idx >= 0 and day_idx >= 0 and hour_idx >= 0)
has_dateh = (date_idx >= 0 and hour_idx >= 0)
if not has_ymdh and not has_dateh:
raise ValueError(
"Hourly file must contain (year,month,day,hour) OR (date,hour)."
)
tmax_idx = _idx(col, ["tmax", "t2mmax", "maxtemp"], optional=True)
tmin_idx = _idx(col, ["tmin", "t2mmin", "mintemp"], optional=True)
temp_idx = _idx(col, ["temp", "temperature", "t2m"], optional=True)
prec_idx = _idx(col, ["prec", "precip", "precipitation", "prectotcorr", "rain", "rainfall"], optional=True)
rh_idx = _idx(col, ["rh", "humidity", "relhumidity", "relativehumidity"], optional=True)
rad_idx = _idx(col, ["rad", "radiation", "solar", "solarrad"], optional=True)
lat_idx = _idx(col, ["latitude", "lat"], optional=True)
if rad_idx < 0 and lat_idx < 0 and latitude == 0.0:
raise ValueError(
"Hourly weather file must contain a radiation column OR a latitude column."
)
for row in reader:
if not row or all(c.strip() == "" for c in row):
continue
# --- Timestamp ---
try:
if has_ymdh:
yr = int(row[year_idx])
mo = int(row[month_idx])
dy = int(row[day_idx])
hr = int(row[hour_idx])
ts = datetime(yr, mo, dy, hr)
else:
raw = row[date_idx].strip().strip('"')
ts = datetime.fromisoformat(raw)
ts = ts.replace(minute=0, second=0, microsecond=0)
hr = int(row[hour_idx])
ts = ts.replace(hour=hr)
except (ValueError, IndexError):
continue
if not (start_year <= ts.year <= end_year):
continue
gw = InputsHourly(date=ts)
# Temperature
if temp_idx >= 0:
v = _pf(row, temp_idx)
if math.isfinite(v):
gw.air_temperature = v
elif tmax_idx >= 0 and tmin_idx >= 0:
tx = _pf(row, tmax_idx)
tn = _pf(row, tmin_idx)
if math.isfinite(tx) and math.isfinite(tn):
gw.air_temperature = 0.5 * (tx + tn)
# Precipitation
if prec_idx >= 0:
v = _pf(row, prec_idx)
if math.isfinite(v) and v >= 0:
gw.precipitation = v
# Relative humidity
if rh_idx >= 0:
v = _pf(row, rh_idx)
if math.isfinite(v):
gw.relative_humidity = max(0.0, min(100.0, v))
else:
t = gw.air_temperature
dp = _dew_point(t, t)
es = 0.61121 * math.exp((17.502 * t) / (240.97 + t))
ea = 0.61121 * math.exp((17.502 * dp) / (240.97 + dp))
gw.relative_humidity = max(0.0, min(100.0, ea / es * 100.0)) if es > 0 else 0.0
# Radiation
if rad_idx >= 0:
v = _pf(row, rad_idx)
if math.isfinite(v) and v >= 0:
gw.rad = v
# Latitude
if lat_idx >= 0:
v = _pf(row, lat_idx)
if math.isfinite(v):
gw.latitude = v
elif latitude != 0.0:
gw.latitude = latitude
# Leaf wetness
gw.leaf_wetness = 1.0 if (gw.relative_humidity > 90.0 or gw.precipitation >= 0.2) else 0.0
day_key = ts.date()
daily_buckets.setdefault(day_key, []).append(gw)
# --- Post-process: fill missing hourly radiation per day ---
for day_key, records in daily_buckets.items():
if not records:
continue
lat = records[0].latitude if records[0].latitude else latitude
temps = [r.air_temperature for r in records]
tmax_d = max(temps)
tmin_d = min(temps)
all_rad_zero = all(r.rad <= 0.0 for r in records)
rd = None
if all_rad_zero and lat:
rd = _day_length(day_key, lat, tmax_d, tmin_d)
for rec in records:
hr = rec.date.hour
if rec.rad <= 0.0 and rd is not None:
rec.rad = rd["gsr_hourly"][hr]
rec.leaf_wetness = 1.0 if (rec.relative_humidity > 90.0 or rec.precipitation >= 0.2) else 0.0
result[rec.date] = rec
return result
franchestyn
¶
FranchestynConfig
dataclass
¶
Configuration for FraNchEstYN model runs.
Attributes:
| Name | Type | Description |
|---|---|---|
param_file |
str |
Path to the main parameter file. |
crop_param_file |
str |
Path to crop parameter file. |
disease_param_file |
str |
Path to disease parameter file. |
fungicide_param_file |
str |
Path to fungicide parameter file. |
reference_path |
str |
Path to reference CSV. |
crop_type |
str |
Crop type (e.g., 'wheat'). |
disease_type |
str |
Disease type (e.g., 'septoria'). Also the disease
on/off switch — |
fungicide_type |
str|None |
Fungicide type (e.g., 'protectant'). Also the
fungicide on/off switch — |
site |
str |
Site name. |
variety |
str |
Variety name. |
disease |
str |
Disease name. |
is_calibration |
bool |
Whether to run calibration. |
calibration_variable |
str |
Calibration variable ('all', 'crop', 'disease'). |
use_gdd |
bool |
Method for crop cycle completion. |
use_prev_day_alignment |
bool |
Day alignment used by the calibration
objective. |
all_row_includes_end_year |
bool |
Whether an |
n_restarts |
int |
Number of calibration restarts. |
max_iter |
int |
Maximum calibration iterations. |
crop_disabled_params |
frozenset[str] |
Crop parameter names to hard-exclude from calibration. |
disease_disabled_params |
frozenset[str] |
Disease parameter names to hard-exclude from calibration. |
Source code in simpest/models/franchestyn.py
@dataclass(frozen=True)
class FranchestynConfig:
"""
Configuration for FraNchEstYN model runs.
Attributes:
param_file (str): Path to the main parameter file.
crop_param_file (str): Path to crop parameter file.
disease_param_file (str): Path to disease parameter file.
fungicide_param_file (str): Path to fungicide parameter file.
reference_path (str): Path to reference CSV.
crop_type (str): Crop type (e.g., 'wheat').
disease_type (str): Disease type (e.g., 'septoria'). Also the disease
on/off switch — ``None`` skips the disease model entirely, giving a
crop-only run.
fungicide_type (str|None): Fungicide type (e.g., 'protectant'). Also the
fungicide on/off switch — ``None`` skips the fungicide model. Set it
whenever treatments are scheduled, or they will be ignored.
site (str): Site name.
variety (str): Variety name.
disease (str): Disease name.
is_calibration (bool): Whether to run calibration.
calibration_variable (str): Calibration variable ('all', 'crop', 'disease').
use_gdd (bool): Method for crop cycle completion. ``False`` (default)
uses calendar-day interpolation of the external crop-model series;
``True`` uses the thermal-time (GDD) based cycle percentage.
use_prev_day_alignment (bool): Day alignment used by the calibration
objective. ``True`` (default) compares the previous simulated day to
the current reference observation (sim[d-1] vs ref[d]); ``False``
uses same-day alignment (sim[d] vs ref[d]).
all_row_includes_end_year (bool): Whether an ``"All"`` sowing row covers
the final year. ``False`` (default) omits ``end_year``; ``True``
includes it. Has no effect for per-year sowing CSVs.
n_restarts (int): Number of calibration restarts.
max_iter (int): Maximum calibration iterations.
crop_disabled_params (frozenset[str]): Crop parameter names to hard-exclude
from calibration.
disease_disabled_params (frozenset[str]): Disease parameter names to
hard-exclude from calibration.
"""
param_file: str = ""
crop_param_file: str = field(default_factory=lambda: _resolve_local_model_file("fr_crop_parameters.json"))
disease_param_file: str = field(default_factory=lambda: _resolve_local_model_file("fr_disease_parameters.json"))
fungicide_param_file: str = field(default_factory=lambda: _resolve_local_model_file("fr_fungicide_parameters.json"))
reference_path: str = _resolve_default_reference()
crop_type: str = "wheat"
disease_type: str = "septoria"
fungicide_type: str | None = "protectant"
site: str = "indiana"
variety: str = "Generic"
disease: str = "thisDisease"
is_calibration: bool = True
calibration_variable: str = "all"
use_gdd: bool = False
use_prev_day_alignment: bool = True
all_row_includes_end_year: bool = False
n_restarts: int = 1
max_iter: int = 100
crop_disabled_params: frozenset[str] = frozenset()
disease_disabled_params: frozenset[str] = frozenset()
build_season_summary(df, site, variety)
¶
Aggregate daily simulation results into a per-season summary.
Rows are grouped by growing season (sowing year) and, for each season, the function reports the epidemic and yield outcomes:
- AUDPC (area under the disease progress curve) is the time integral of disease severity (expressed as a percentage) over the season, computed by the trapezoidal rule against the calendar date.
- Yield loss is the gap between attainable and actual yield, reported
both in absolute terms and as a percentage,
(attainable − actual) / attainable × 100. - Peak attainable and actual yield and above-ground biomass, peak disease severity, and season-level weather aggregates (mean temperatures and humidity, total precipitation, radiation, and leaf wetness) are also included where the corresponding columns are present.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
df |
pd.DataFrame |
Daily simulation results. |
required |
site |
str |
Site name written to each summary row. |
required |
variety |
str |
Variety name written to each summary row. |
required |
Returns:
| Type | Description |
|---|---|
pd.DataFrame |
One row per growing season; empty if the input has no valid post-sowing days. |
Source code in simpest/models/franchestyn.py
def build_season_summary(df: pd.DataFrame, site: str, variety: str) -> pd.DataFrame:
"""Aggregate daily simulation results into a per-season summary.
Rows are grouped by growing season (sowing year) and, for each season, the
function reports the epidemic and yield outcomes:
- **AUDPC** (area under the disease progress curve) is the time integral of
disease severity (expressed as a percentage) over the season, computed by
the trapezoidal rule against the calendar date.
- **Yield loss** is the gap between attainable and actual yield, reported
both in absolute terms and as a percentage,
``(attainable − actual) / attainable × 100``.
- Peak attainable and actual yield and above-ground biomass, peak disease
severity, and season-level weather aggregates (mean temperatures and
humidity, total precipitation, radiation, and leaf wetness) are also
included where the corresponding columns are present.
Args:
df (pd.DataFrame): Daily simulation results.
site (str): Site name written to each summary row.
variety (str): Variety name written to each summary row.
Returns:
pd.DataFrame: One row per growing season; empty if the input has no
valid post-sowing days.
"""
if df.empty:
return pd.DataFrame()
d = df.copy()
d["Date"] = pd.to_datetime(d["Date"], dayfirst=True, errors="coerce")
d = d.dropna(subset=["Date"])
if "DaysAfterSowing" in d.columns:
d = d[d["DaysAfterSowing"] > 0]
if d.empty:
return pd.DataFrame()
# Group by the crop's growing season (= sowing year). Fall back to the
# calendar year only where the season label is missing or unset (0), so
# seasons that cross 1 January are not split across two summary rows.
if "GrowingSeason" not in d.columns:
d["GrowingSeason"] = d["Date"].dt.year
else:
season = pd.to_numeric(d["GrowingSeason"], errors="coerce")
d["GrowingSeason"] = season.where(season > 0, d["Date"].dt.year)
rows = []
for season, g in d.groupby("GrowingSeason"):
g = g.sort_values("Date")
y = g["DiseaseSeverity"].fillna(0.0).to_numpy(dtype=float) * 100.0
x = g["Date"].map(pd.Timestamp.toordinal).to_numpy(dtype=float)
audpc = float(np.trapezoid(y, x)) if len(g) >= 2 else 0.0
yield_att = float(g["YieldAttainable"].max()) if "YieldAttainable" in g.columns else 0.0
yield_act = float(g["YieldActual"].max()) if "YieldActual" in g.columns else 0.0
agb_att = float(g["AGBattainable"].max()) if "AGBattainable" in g.columns else 0.0
agb_act = float(g["AGBactual"].max()) if "AGBactual" in g.columns else 0.0
dis_sev = float(g["DiseaseSeverity"].max()) if "DiseaseSeverity" in g.columns else 0.0
loss_raw = yield_att - yield_act
loss_perc = (loss_raw / yield_att * 100.0) if yield_att > 0 else 0.0
row = {
"GrowingSeason": int(season),
"Site": site,
"Variety": variety,
"AUDPC": audpc,
"DiseaseSeverity": dis_sev,
"YieldAttainable": yield_att,
"YieldActual": yield_act,
"YieldLossRaw": loss_raw,
"YieldLossPerc": loss_perc,
"AGBattainable": agb_att,
"AGBactual": agb_act,
}
if "Tmax" in g.columns:
row["AveTx"] = float(g["Tmax"].mean())
if "Tmin" in g.columns:
row["AveTn"] = float(g["Tmin"].mean())
if "RHx" in g.columns:
row["AveRHx"] = float(g["RHx"].mean())
if "RHn" in g.columns:
row["AveRHn"] = float(g["RHn"].mean())
if "TotalPrec" in g.columns:
row["TotalPrec"] = float(g["TotalPrec"].sum())
if "TotalRad" in g.columns:
row["TotalRad"] = float(g["TotalRad"].sum())
if "TotalLW" in g.columns:
row["TotalLW"] = float(g["TotalLW"].sum())
rows.append(row)
return pd.DataFrame(rows).sort_values("GrowingSeason").reset_index(drop=True)
deactivate_calibration(params, disable_list)
¶
Disable calibration for selected parameter names in a parameter section.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
params |
dict |
Mapping of parameter name to parameter specification dict. |
required |
disable_list |
Iterable[str] |
Parameter names for which calibration should be forced to False. |
required |
Source code in simpest/models/franchestyn.py
def deactivate_calibration(params: dict, disable_list) -> None:
"""
Disable calibration for selected parameter names in a parameter section.
Args:
params (dict): Mapping of parameter name to parameter specification dict.
disable_list (Iterable[str]): Parameter names for which calibration
should be forced to False.
"""
for param_name in disable_list:
if param_name in params:
params[param_name]["calibration"] = False
return params
run_franchestyn(weather_path, management_path, start_year, end_year, config, cropmodel_path=None, crop_param_file=None, disease_param_file=None, fungicide_param_file=None)
¶
Run the FraNchEstYN model with the given configuration and input files.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
weather_path |
str |
Path to weather input file. |
required |
management_path |
str |
Path to management input file. |
required |
start_year |
int |
Start year for simulation. |
required |
end_year |
int |
End year for simulation. |
required |
config |
FranchestynConfig |
FraNchEstYN configuration object. |
required |
cropmodel_path |
str|None |
Path to crop model data file. |
None |
crop_param_file |
str|None |
Path to crop parameter file. |
None |
disease_param_file |
str|None |
Path to disease parameter file. |
None |
fungicide_param_file |
str|None |
Path to fungicide parameter file. |
None |
Returns:
| Type | Description |
|---|---|
dict |
Dictionary with simulation outputs and summary. |
Source code in simpest/models/franchestyn.py
def run_franchestyn(
weather_path: str,
management_path: str,
start_year: int,
end_year: int,
config: FranchestynConfig,
cropmodel_path: str | None = None,
crop_param_file: str | None = None,
disease_param_file: str | None = None,
fungicide_param_file: str | None = None,
) -> dict:
"""
Run the FraNchEstYN model with the given configuration and input files.
Args:
weather_path (str): Path to weather input file.
management_path (str): Path to management input file.
start_year (int): Start year for simulation.
end_year (int): End year for simulation.
config (FranchestynConfig): FraNchEstYN configuration object.
cropmodel_path (str|None, optional): Path to crop model data file.
crop_param_file (str|None, optional): Path to crop parameter file.
disease_param_file (str|None, optional): Path to disease parameter file.
fungicide_param_file (str|None, optional): Path to fungicide parameter file.
Returns:
dict: Dictionary with simulation outputs and summary.
"""
crop_param_file = crop_param_file or config.crop_param_file
disease_param_file = disease_param_file or config.disease_param_file
fungicide_param_file = fungicide_param_file or config.fungicide_param_file
runner = FranchestynRunner(
weather_dir=weather_path,
param_file=config.param_file,
sowing_file=management_path,
ref_dir=config.reference_path,
crop_model_dir=cropmodel_path,
site=config.site,
variety=config.variety,
disease=config.disease,
start_year=start_year,
end_year=end_year,
weather_time_step="daily",
calibration_variable=config.calibration_variable,
is_calibration=config.is_calibration,
use_gdd=config.use_gdd,
use_prev_day_alignment=config.use_prev_day_alignment,
all_row_includes_end_year=config.all_row_includes_end_year,
crop_type=config.crop_type,
crop_param_file=crop_param_file,
disease_param_file=disease_param_file,
disease_type=config.disease_type,
fungicide_param_file=fungicide_param_file,
fungicide_type=config.fungicide_type,
)
best_params = {}
if config.is_calibration:
disabled_by_class = {
"crop": set(config.crop_disabled_params),
"disease": set(config.disease_disabled_params),
}
optimizer = FranchestynOptimizer(
runner=runner,
calibration_variable=config.calibration_variable,
n_restarts=config.n_restarts,
max_iter=config.max_iter,
disabled_by_class=disabled_by_class,
)
best_params = optimizer.calibrate()
date_outputs = runner.run(param_values=best_params)
else:
date_outputs = runner.run()
include_crop = config.calibration_variable in ("crop", "all")
include_disease = config.calibration_variable in ("disease", "all")
rmse = runner.compute_rmse(
date_outputs,
include_crop=include_crop,
include_disease=include_disease,
)
records = _outputs_to_records(date_outputs)
return {
"outputs": {
"simulation": records,
"summary": {
"rmse": rmse,
"is_calibration": config.is_calibration,
"calibration_variable": config.calibration_variable,
"best_params": best_params,
},
}
}
save_calibrated_parameters_csv(best_params, output_root, site, variety, filename=None, config=None, r_like=False)
¶
Save calibration parameters to a CSV file.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
best_params |
dict |
Best parameter values from calibration. |
required |
output_root |
Path |
Output root directory. |
required |
site |
str |
Site name. |
required |
variety |
str |
Variety name. |
required |
filename |
str|None |
Output filename. If None, uses a default pattern. |
None |
config |
FranchestynConfig|None |
Configuration used to load
parameter metadata when |
None |
r_like |
bool |
If True, export an extended diagnostics table that includes each parameter's default, bounds, unit, calibration flag, and calibrated value, in addition to the calibrated values. |
False |
Returns:
| Type | Description |
|---|---|
Path|None |
Path to the saved CSV file, or None if best_params is empty. |
Source code in simpest/models/franchestyn.py
def save_calibrated_parameters_csv(
best_params: dict,
output_root: Path,
site: str,
variety: str,
filename: str | None = None,
config: FranchestynConfig | None = None,
r_like: bool = False,
) -> Path | None:
"""
Save calibration parameters to a CSV file.
Args:
best_params (dict): Best parameter values from calibration.
output_root (Path): Output root directory.
site (str): Site name.
variety (str): Variety name.
filename (str|None, optional): Output filename. If None, uses a default pattern.
config (FranchestynConfig|None, optional): Configuration used to load
parameter metadata when ``r_like=True``.
r_like (bool, optional): If True, export an extended diagnostics table
that includes each parameter's default, bounds, unit, calibration
flag, and calibrated value, in addition to the calibrated values.
Returns:
Path|None: Path to the saved CSV file, or None if best_params is empty.
"""
out_dir = output_root / "SimulationExperimentTemplate" / "calibratedParameters"
out_dir.mkdir(parents=True, exist_ok=True)
if filename is None:
filename = f"calibratedParameters_{site}_{variety}.csv"
output_file = out_dir / filename
if r_like:
cfg = config or FranchestynConfig(site=site, variety=variety)
rows = []
parameter_file = f"parameters_{site}_{variety}.csv"
def _load_param_defs(json_path: str, model: str, kind: str):
path = Path(json_path)
if not path.exists():
return
data = json.loads(path.read_text(encoding="utf-8"))
section = data.get(kind, {})
for param, spec in section.items():
default_val = spec.get("value", "NA")
min_val = spec.get("min", "NA")
max_val = spec.get("max", "NA")
unit = spec.get("unit", "unitless")
calib_enabled = bool(spec.get("calibration", False))
key = f"{model}_{param}"
calibrated_val = best_params.get(key, default_val)
value_col = calibrated_val if calib_enabled else "NA"
file_col = parameter_file if calib_enabled else "NA"
rows.append(
{
"Model": model,
"Parameter": param,
"unit": unit,
"min": min_val,
"max": max_val,
"default": default_val,
"calibration": "TRUE" if calib_enabled else "FALSE",
"calibrated": calibrated_val,
"value": value_col,
"file": file_col,
"facet_label": f"{param} ({unit})",
}
)
_load_param_defs(cfg.crop_param_file, "crop", cfg.crop_type)
_load_param_defs(cfg.disease_param_file, "disease", cfg.disease_type)
with output_file.open("w", newline="", encoding="utf-8") as f_out:
fieldnames = [
"Model",
"Parameter",
"unit",
"min",
"max",
"default",
"calibration",
"calibrated",
"value",
"file",
"facet_label",
]
writer = csv.DictWriter(f_out, fieldnames=fieldnames, delimiter=",")
writer.writeheader()
for row in rows:
writer.writerow(row)
return output_file
if not best_params:
return None
with output_file.open("w", newline="", encoding="utf-8") as f_out:
writer = csv.DictWriter(f_out, fieldnames=["model", "param", "value"], delimiter=",")
writer.writeheader()
for key, value in sorted(best_params.items()):
if "_" in key:
model, param = key.split("_", 1)
else:
model, param = "", key
writer.writerow({"model": model, "param": param, "value": value})
return output_file
save_season_summary_csv(summary_df, output_root, filename='franchestyn_season_summary.csv')
¶
Save season summary DataFrame to a CSV file.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
summary_df |
pd.DataFrame |
Season summary DataFrame. |
required |
output_root |
Path |
Output root directory. |
required |
filename |
str |
Output filename. Defaults to 'franchestyn_season_summary.csv'. |
'franchestyn_season_summary.csv' |
Returns:
| Type | Description |
|---|---|
Path|None |
Path to the saved CSV file, or None if summary_df is empty. |
Source code in simpest/models/franchestyn.py
def save_season_summary_csv(summary_df: pd.DataFrame, output_root: Path, filename: str = "franchestyn_season_summary.csv") -> Path | None:
"""
Save season summary DataFrame to a CSV file.
Args:
summary_df (pd.DataFrame): Season summary DataFrame.
output_root (Path): Output root directory.
filename (str, optional): Output filename. Defaults to 'franchestyn_season_summary.csv'.
Returns:
Path|None: Path to the saved CSV file, or None if summary_df is empty.
"""
if summary_df.empty:
return None
output_file = output_root / "SimulationExperimentTemplate" / filename
summary_df.to_csv(output_file, index=False)
return output_file
save_simulation_results_csv(res_ot_simulation, output_root, filename='franchestyn_simulation_results.csv')
¶
Save FraNchEstYN simulation results to a CSV file.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
res_ot_simulation |
list of dict |
Simulation output records. |
required |
output_root |
Path |
Output root directory. |
required |
filename |
str |
Output filename. Defaults to 'franchestyn_simulation_results.csv'. |
'franchestyn_simulation_results.csv' |
Returns:
| Type | Description |
|---|---|
Path |
Path to the saved CSV file. |
Source code in simpest/models/franchestyn.py
def save_simulation_results_csv(res_ot_simulation, output_root: Path, filename: str = "franchestyn_simulation_results.csv") -> Path:
"""
Save FraNchEstYN simulation results to a CSV file.
Args:
res_ot_simulation (list of dict): Simulation output records.
output_root (Path): Output root directory.
filename (str, optional): Output filename. Defaults to 'franchestyn_simulation_results.csv'.
Returns:
Path: Path to the saved CSV file.
"""
output_file = output_root / "SimulationExperimentTemplate" / filename
df = pd.DataFrame(res_ot_simulation)
df.to_csv(output_file, index=False)
return output_file
simplace
¶
SimplaceConfig
dataclass
¶
Configuration for Simplace runs.
Attributes:
| Name | Type | Description |
|---|---|---|
install_dir |
str |
Path to Simplace installation directory. |
work_dir |
str |
Path to Simplace working directory. |
output_dir |
str |
Path to Simplace output directory. |
solution_path |
str |
Path to Simplace solution XML file. |
project_path |
str |
Path to Simplace project XML file. |
Source code in simpest/models/simplace.py
@dataclass(frozen=True)
class SimplaceConfig:
"""
Configuration for Simplace runs.
Attributes:
install_dir (str): Path to Simplace installation directory.
work_dir (str): Path to Simplace working directory.
output_dir (str): Path to Simplace output directory.
solution_path (str): Path to Simplace solution XML file.
project_path (str): Path to Simplace project XML file.
"""
install_dir: str = "C:/ParamVC/Research/simplace+/simplace_portable/workspace/"
work_dir: str = "C:/ParamVC/Research/simplace+/simplace_portable/workspace/simplace_run/simulation/"
output_dir: str = "C:/ParamVC/Research/simplace+/simplace_out/"
solution_path: str = "SimulationExperimentTemplate/solution/Lintul5_indiana.sol.xml"
project_path: str = "SimulationExperimentTemplate/project/Lintul5All_indiana.proj.xml"
build_management(output_root, project_row, crop='wheat', variety='generic', yearly_sowing_doy=None)
¶
Build and save a management CSV file for FraNchEstYN.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
output_root |
Path |
Output root directory. |
required |
project_row |
dict |
Project row dictionary. |
required |
crop |
str |
Crop name. Defaults to "wheat". |
'wheat' |
variety |
str |
Variety name. Defaults to "generic". |
'generic' |
yearly_sowing_doy |
dict[int, int] | None |
Mapping of simulation year to sowing DOY. If omitted, a single "All" row is written to preserve backward compatibility. |
None |
Returns:
| Type | Description |
|---|---|
Path |
Path to the management CSV file. |
Source code in simpest/models/simplace.py
def build_management(
output_root: Path,
project_row: dict,
crop: str = "wheat",
variety: str = "generic",
yearly_sowing_doy: Optional[dict[int, int]] = None,
) -> Path:
"""
Build and save a management CSV file for FraNchEstYN.
Args:
output_root (Path): Output root directory.
project_row (dict): Project row dictionary.
crop (str, optional): Crop name. Defaults to "wheat".
variety (str, optional): Variety name. Defaults to "generic".
yearly_sowing_doy (dict[int, int] | None, optional): Mapping of
simulation year to sowing DOY. If omitted, a single "All" row
is written to preserve backward compatibility.
Returns:
Path: Path to the management CSV file.
"""
sowing_doy = max(1, project_row["idem"] - 7)
if yearly_sowing_doy is None:
candidate = project_row.get("yearly_sowing_doy") if isinstance(project_row, dict) else None
if isinstance(candidate, dict) and candidate:
yearly_sowing_doy = {
int(year): int(doy)
for year, doy in candidate.items()
}
management_dst = output_root / "SimulationExperimentTemplate" / "management_franchestyn.csv"
with management_dst.open("w", newline="", encoding="utf-8") as f_out:
writer = csv.DictWriter(
f_out,
fieldnames=["site", "crop", "variety", "year", "sowingDOY"],
delimiter=",",
)
writer.writeheader()
if yearly_sowing_doy:
for year in sorted(yearly_sowing_doy):
writer.writerow(
{
"site": project_row["location"],
"crop": crop,
"variety": variety,
"year": year,
"sowingDOY": max(1, int(yearly_sowing_doy[year])),
}
)
else:
writer.writerow(
{
"site": project_row["location"],
"crop": crop,
"variety": variety,
"year": "All",
"sowingDOY": sowing_doy,
}
)
return management_dst
convert_weather(work_root, output_root, location)
¶
Convert Simplace weather file to FraNchEstYN-compatible CSV format.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
work_root |
Path |
Simplace workspace root directory. |
required |
output_root |
Path |
Output root directory. |
required |
location |
str |
Location name (used for weather file). |
required |
Returns:
| Type | Description |
|---|---|
Path |
Path to the converted weather CSV file. |
Source code in simpest/models/simplace.py
def convert_weather(work_root: Path, output_root: Path, location: str) -> Path:
"""
Convert Simplace weather file to FraNchEstYN-compatible CSV format.
Args:
work_root (Path): Simplace workspace root directory.
output_root (Path): Output root directory.
location (str): Location name (used for weather file).
Returns:
Path: Path to the converted weather CSV file.
"""
weather_src = work_root / "SimulationExperimentTemplate" / "data" / "weather" / f"{location}.txt"
weather_dst = output_root / "SimulationExperimentTemplate" / "weather_franchestyn.csv"
with weather_src.open(newline="", encoding="utf-8-sig") as f_in, weather_dst.open(
"w", newline="", encoding="utf-8"
) as f_out:
reader = csv.DictReader(f_in, delimiter="\t")
reader.fieldnames = [h.lstrip("\ufeff").strip() for h in reader.fieldnames]
writer = csv.DictWriter(
f_out,
fieldnames=["site", "year", "month", "day", "tx", "tn", "p", "rad", "vp", "rhx", "rhn"],
delimiter=",",
)
writer.writeheader()
for row in reader:
d = datetime.strptime(row["CURRENTDAY"], "%d.%m.%Y")
tmax = float(row["Tmax"])
tmin = float(row["Tmin"])
vp = float(row["VapourPressure"])
es_tmax = _saturation_vapor_pressure(tmax)
es_tmin = _saturation_vapor_pressure(tmin)
rhx = min(100.0, max(0.0, 100.0 * vp / es_tmin)) if es_tmin > 0 else 0.0
rhn = min(100.0, max(0.0, 100.0 * vp / es_tmax)) if es_tmax > 0 else 0.0
writer.writerow(
{
"site": location,
"year": d.year,
"month": d.month,
"day": d.day,
"tx": tmax,
"tn": tmin,
"p": row["Precipitation"],
"rad": float(row["Irradiation"]) / 1000.0,
"vp": vp,
"rhx": rhx,
"rhn": rhn,
}
)
return weather_dst
export_crop_model_data(output_root, project_row)
¶
Export crop model data from Simplace daily output to a CSV file for FraNchEstYN.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
output_root |
Path |
Output root directory. |
required |
project_row |
dict |
Project row dictionary. |
required |
Returns:
| Type | Description |
|---|---|
Path |
Path to the exported crop model CSV file. |
Source code in simpest/models/simplace.py
def export_crop_model_data(output_root: Path, project_row: dict) -> Path:
"""
Export crop model data from Simplace daily output to a CSV file for FraNchEstYN.
Args:
output_root (Path): Output root directory.
project_row (dict): Project row dictionary.
Returns:
Path: Path to the exported crop model CSV file.
"""
src = output_root / "SimulationExperimentTemplate" / f"{project_row['location']}{project_row['iopt']}_daily.csv"
dst = output_root / "SimulationExperimentTemplate" / "cropModel_data.csv"
with src.open(newline="", encoding="utf-8") as f_in, dst.open("w", newline="", encoding="utf-8") as f_out:
reader = csv.DictReader(f_in, delimiter=";")
writer = csv.DictWriter(
f_out,
fieldnames=["year", "doy", "agb", "yield", "fint", "lai", "gdd"],
delimiter=",",
)
writer.writeheader()
for row in reader:
if all(float(row[name]) == 0.0 for name in ("TAGB", "WSO", "FINT", "LAI")):
continue
d = datetime.strptime(row["CURRENT.DATE"], "%d.%m.%Y")
writer.writerow(
{
"year": d.year,
"doy": d.timetuple().tm_yday,
"agb": float(row["TAGB"]) * 10,
"yield": float(row["WSO"]) * 10,
"fint": row["FINT"],
"lai": row["LAI"],
"gdd": row["TSUM"],
}
)
return dst
extract_yearly_sowing_doy(work_root, project_row, start_year, end_year)
¶
Extract per-year sowing DOY from projectdata CSV files.
The function prefers explicit ISOW values when present. If no per-year
values are available, it falls back to a constant sowing DOY derived from
IDEM - 7 to preserve current behavior.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
work_root |
Path |
Simplace workspace root. |
required |
project_row |
dict |
Selected project row metadata. |
required |
start_year |
int |
Start year (inclusive). |
required |
end_year |
int |
End year (inclusive). |
required |
Returns:
| Type | Description |
|---|---|
dict[int, int] |
Mapping year -> sowing DOY. |
Source code in simpest/models/simplace.py
def extract_yearly_sowing_doy(
work_root: Path,
project_row: dict,
start_year: int,
end_year: int,
) -> dict[int, int]:
"""
Extract per-year sowing DOY from projectdata CSV files.
The function prefers explicit `ISOW` values when present. If no per-year
values are available, it falls back to a constant sowing DOY derived from
`IDEM - 7` to preserve current behavior.
Args:
work_root (Path): Simplace workspace root.
project_row (dict): Selected project row metadata.
start_year (int): Start year (inclusive).
end_year (int): End year (inclusive).
Returns:
dict[int, int]: Mapping year -> sowing DOY.
"""
projectdata_dir = work_root / "SimulationExperimentTemplate" / "data" / "projectdata"
if not projectdata_dir.exists():
return {}
location = str(project_row.get("location", "")).strip().lower()
iopt = _to_int(project_row.get("iopt"))
idem = _to_int(project_row.get("idem"))
default_sowing_doy = max(1, (idem or 8) - 7)
rows: list[dict] = []
for project_csv in sorted(projectdata_dir.glob("*.csv")):
with project_csv.open(newline="", encoding="utf-8-sig") as f_proj:
reader = csv.DictReader(f_proj, delimiter=";")
if not reader.fieldnames:
continue
name_to_col = {name.strip().lower(): name for name in reader.fieldnames}
col_location = name_to_col.get("location")
col_iopt = name_to_col.get("iopt")
col_start = name_to_col.get("startdate")
col_end = name_to_col.get("enddate")
col_isow = name_to_col.get("isow")
col_idem = name_to_col.get("idem")
for raw in reader:
row_location = str(raw.get(col_location, "")).strip().lower() if col_location else ""
if location and row_location and row_location != location:
continue
row_iopt = _to_int(raw.get(col_iopt)) if col_iopt else None
if iopt is not None and row_iopt is not None and row_iopt != iopt:
continue
start_dt = _parse_project_date(raw.get(col_start, "")) if col_start else None
end_dt = _parse_project_date(raw.get(col_end, "")) if col_end else None
if start_dt is None and end_dt is None:
continue
year = start_dt.year if start_dt is not None else end_dt.year
if year < start_year or year > end_year:
continue
sow_doy = _to_int(raw.get(col_isow)) if col_isow else None
if sow_doy is None:
row_idem = _to_int(raw.get(col_idem)) if col_idem else idem
sow_doy = max(1, (row_idem or 8) - 7)
rows.append(
{
"year": year,
"sowing_doy": sow_doy,
"span_days": ((end_dt - start_dt).days if start_dt and end_dt else 999999),
}
)
# Prefer the most specific row per year (shortest start-end span).
yearly_sowing: dict[int, int] = {}
for row in sorted(rows, key=lambda x: (x["year"], x["span_days"])):
yearly_sowing[row["year"]] = row["sowing_doy"]
if yearly_sowing:
return yearly_sowing
# Backward-compatible fallback: constant sowing day across all years.
return {year: default_sowing_doy for year in range(start_year, end_year + 1)}
get_project_row(work_root, selected_line)
¶
Get a row from the project CSV file by line index.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
work_root |
Path |
Root directory for Simplace workspace. |
required |
selected_line |
int |
Line number (1-based) to select from the project CSV. |
required |
Returns:
| Type | Description |
|---|---|
dict |
Dictionary with project row data. |
Source code in simpest/models/simplace.py
def get_project_row(work_root: Path, selected_line: int) -> dict:
"""
Get a row from the project CSV file by line index.
Args:
work_root (Path): Root directory for Simplace workspace.
selected_line (int): Line number (1-based) to select from the project CSV.
Returns:
dict: Dictionary with project row data.
"""
project_csv = (
work_root
/ "SimulationExperimentTemplate"
/ "data"
/ "projectdata"
/ "lintul5all_indiana.csv"
)
with project_csv.open(newline="", encoding="utf-8-sig") as f_proj:
proj_reader = csv.DictReader(f_proj, delimiter=";")
project_rows = list(proj_reader)
row = project_rows[selected_line - 1]
key_map = {k.strip().lower(): k for k in row.keys()}
def field(*aliases: str, default: str = "") -> str:
for alias in aliases:
key = key_map.get(alias.lower())
if key is not None:
return str(row.get(key, default)).strip()
return default
project_row = {
"projectid": field("projectid"),
"simulationid": field("simulationid"),
"startdate": field("startdate"),
"enddate": field("enddate"),
"location": field("location"),
"iopt": field("iopt"),
"idem": _to_int(field("idem")) or 0,
"irri": field("irri"),
}
start_dt = _parse_project_date(project_row["startdate"])
end_dt = _parse_project_date(project_row["enddate"])
if start_dt and end_dt:
project_row["yearly_sowing_doy"] = extract_yearly_sowing_doy(
work_root,
project_row,
start_dt.year,
end_dt.year,
)
return project_row
get_simplace_directories(shell)
¶
Wrapper for simplace.getSimplaceDirectories to avoid exposing simplace in notebooks.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
shell |
Simplace shell instance. |
required |
Returns:
| Type | Description |
|---|---|
dict |
Simplace directories information. |
Source code in simpest/models/simplace.py
def get_simplace_directories(shell):
"""
Wrapper for simplace.getSimplaceDirectories to avoid exposing simplace in notebooks.
Args:
shell: Simplace shell instance.
Returns:
dict: Simplace directories information.
"""
return simplace.getSimplaceDirectories(shell)
init_simplace(config)
¶
Initialize and return a Simplace instance using the given configuration.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
config |
SimplaceConfig |
Simplace configuration object. |
required |
Returns:
| Type | Description |
|---|---|
Simplace instance. |
Source code in simpest/models/simplace.py
def init_simplace(config: SimplaceConfig):
"""
Initialize and return a Simplace instance using the given configuration.
Args:
config (SimplaceConfig): Simplace configuration object.
Returns:
Simplace instance.
"""
global _SIMPLACE_INSTANCE
if _SIMPLACE_INSTANCE is not None and jpype.isJVMStarted():
return _SIMPLACE_INSTANCE
if jpype.isJVMStarted():
Wrapper = jpype.JClass('net.simplace.sim.wrapper.SimplaceWrapper')
_SIMPLACE_INSTANCE = Wrapper(config.work_dir, config.output_dir, None, None)
return _SIMPLACE_INSTANCE
_SIMPLACE_INSTANCE = simplace.initSimplace(
installDir=config.install_dir,
workDir=config.work_dir,
outputDir=config.output_dir,
)
return _SIMPLACE_INSTANCE
merge_simplace_and_franchestyn(output_root, project_row, franchestyn_df, out_name='merged_simulation_data.csv')
¶
Merge Simplace and FraNchEstYN daily outputs into a single CSV file.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
output_root |
Path |
Output root directory. |
required |
project_row |
dict |
Project row dictionary. |
required |
franchestyn_df |
pd.DataFrame |
FraNchEstYN daily output DataFrame. |
required |
out_name |
str |
Output filename. Defaults to "merged_simulation_data.csv". |
'merged_simulation_data.csv' |
Returns:
| Type | Description |
|---|---|
Path |
Path to the merged CSV file. |
Source code in simpest/models/simplace.py
def merge_simplace_and_franchestyn(
output_root: Path,
project_row: dict,
franchestyn_df: pd.DataFrame,
out_name: str = "merged_simulation_data.csv",
) -> Path:
"""
Merge Simplace and FraNchEstYN daily outputs into a single CSV file.
Args:
output_root (Path): Output root directory.
project_row (dict): Project row dictionary.
franchestyn_df (pd.DataFrame): FraNchEstYN daily output DataFrame.
out_name (str, optional): Output filename. Defaults to "merged_simulation_data.csv".
Returns:
Path: Path to the merged CSV file.
"""
simplace_daily_path = output_root / "SimulationExperimentTemplate" / f"{project_row['location']}{project_row['iopt']}_daily.csv"
simplace_df = pd.read_csv(simplace_daily_path, sep=";")
simplace_df.columns = simplace_df.columns.str.strip()
franchestyn_df.columns = franchestyn_df.columns.str.strip()
simplace_df["date"] = pd.to_datetime(simplace_df["CURRENT.DATE"], format="%d.%m.%Y")
franchestyn_df["date"] = pd.to_datetime(franchestyn_df["Date"], dayfirst=True)
sim_df = simplace_df.rename(columns={c: f"{c}_S" for c in simplace_df.columns if c != "date"})
fran_df = franchestyn_df.rename(columns={c: f"{c}_F" for c in franchestyn_df.columns if c != "date"})
merged = pd.merge(sim_df, fran_df, on="date", how="left")
merged = merged.rename(
columns={
"WSO_S": "WSO_S_g_m2",
"YieldAttainable_F": "YieldAttainable_F_kg_ha",
"YieldActual_F": "YieldActual_F_kg_ha",
}
)
out_path = output_root / "SimulationExperimentTemplate" / out_name
merged.to_csv(out_path, index=False)
return out_path
run_simplace(shell, config, project_lines)
¶
Run a Simplace project for the specified project lines.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
shell |
Simplace shell instance. |
required | |
config |
SimplaceConfig |
Simplace configuration object. |
required |
project_lines |
list[int] |
List of project line indices to run. |
required |
Source code in simpest/models/simplace.py
def run_simplace(shell, config: SimplaceConfig, project_lines: list[int]):
"""
Run a Simplace project for the specified project lines.
Args:
shell: Simplace shell instance.
config (SimplaceConfig): Simplace configuration object.
project_lines (list[int]): List of project line indices to run.
"""
simplace.openProject(shell, config.solution_path, config.project_path)
simplace.setProjectLines(shell, project_lines)
simplace.runProject(shell)