Reorganise makecropstate function, move to crop models

src/Persefone.jl

@@ -148,6 +148,9 @@ model has to define a type `CropState <: AbstractCropState`.
 """
 abstract type AbstractCropState end
 
+function makecropstate end
+setsoiltype!(cropstate::AbstractCropState, soiltype) = nothing
+
 function stepagent! end
 
 ## include all module files (note that the order matters - if file

src/crop/almass.jl

@@ -49,6 +49,7 @@ using Persefone:
     Management,
     cm,
     SimulationModel,
+    SoilType,
     fertiliser,
     maxtemp,
     mintemp,
@@ -62,6 +63,7 @@ import Persefone:
     cropheight,
     cropcover,
     cropyield,
+    makecropstate,
     sow!,
     harvest!,
     isharvestable,
@@ -202,6 +204,12 @@ function isharvestable(cs::CropState)
     return cs.ddegs >= last(gdd)
 end
 
+function makecropstate(croptype::CropType, model::SimulationModel, soiltype::SoilType)
+    phase = month(model.date) < 3 ? janfirst : marchfirst
+    return CropState(; croptype, phase)
+end
+
+
 # Constant in original ALMaSS code:
 #     `EL_VEG_START_LAIT` from `Landscape/Elements.cpp`, line 238-239
 const VEG_START_LAIT::Float64 = 1.08

src/crop/aquacrop.jl

@@ -25,9 +25,11 @@ import Persefone:
     cropheight,
     cropcover,
     cropyield,
+    makecropstate,
     sow!,
     harvest!,
-    isharvestable
+    isharvestable,
+    setsoiltype!
 
 using Unitful: @u_str
 
@@ -167,6 +169,14 @@ cropheight(cs::CropState) = cs.height  # TODO: calculate from AquaCrop state inf
 cropcover(cs::CropState) = AquaCrop.canopycover(cs.cropstate)
 cropyield(cs::CropState) = AquaCrop.dryyield(cs.cropstate)  # TODO: there is also freshyield
 isharvestable(cs::CropState) = true  # TODO: implement this correctly
+makecropstate(crop_type::CropType, model::SimulationModel, soiltype::SoilType) =
+    CropState(crop_type, soiltype, model)
+function setsoiltype!(cs::CropState, soiltype::SoilType)
+    # TODO: this does not affect the actual soil type of the state of
+    # the aquacrop simulation
+    cs.soiltype = soiltype
+    return nothing
+end
 
 """
     stepagent!(cropstate, model)

src/crop/cropmodels.jl

@@ -32,6 +32,12 @@ end
 Initialise the farm plots in the simulation model.
 """
 function initfields!(model::SimulationModel)
+    # TODO: we take the soiltype from the first pixel in the farmplot,
+    # then later we fix the soiltypes of the farmplots and cropstates
+    # to be the most common soiltype of a farmplot.  It would be
+    # better to first collect all the data for pixels belonging
+    # together, and then create the farmplots and cropstates with the
+    # correct soiltype directly
     convertid = Dict{Int64,Int64}()
     width, height = size(model.landscape)
     for x in 1:width
@@ -46,7 +52,7 @@ function initfields!(model::SimulationModel)
                 push!(model.farmplots[objectid].pixels, (x,y))
             else
                 soiltype = model.landscape[x,y].soiltype
-                cropstate = makecropstate(model, soiltype)
+                cropstate = makecropstate(model.crops["natural grass"], model, soiltype)
                 fp = FarmPlot(length(model.farmplots) + 1, [(x, y)], -1, soiltype, cropstate)
                 push!(model.farmplots, fp)
                 model.landscape[x,y].fieldid = fp.id
@@ -58,48 +64,50 @@ function initfields!(model::SimulationModel)
     # Adjust farmplot soil type to most common soil type of its
     # pixels, and set the cropstate soil type to this as well.
     for fp in model.farmplots
-        fp.soiltype = mode(map(p -> model.landscape[p[1], p[2]].soiltype, fp.pixels))
-        set_cropstate_soiltype!(fp, model)
+        soiltype = mode(map(p -> model.landscape[p[1], p[2]].soiltype, fp.pixels))
+        setsoiltype!(fp, soiltype)
     end
 
     @info "Initialised $(length(model.farmplots)) farm plots."
 end
 
-# TODO: this function should be moved to the individual crop models,
-# and overloaded on the type of cropstate
-"""
-    set_cropstate_soiltype!(farmplot, model)
-"""
-function set_cropstate_soiltype!(farmplot::FarmPlot, model::SimulationModel)
-    if @param(crop.cropmodel) == "aquacrop"
-        farmplot.cropstate.soiltype = farmplot.soiltype
-    else
-        # do nothing for other cropmodels
-    end
-end
+# # TODO: this function should be moved to the individual crop models,
+# # and overloaded on the type of cropstate (see notes above)
+# """
+#     setsoiltype!(farmplot, soiltype, model)
+# """
+# function setsoiltype!(farmplot::FarmPlot, soiltype::SoilType, model::SimulationModel)
+#     farmplot.soiltype = soiltype
+#     if @param(crop.cropmodel) == "aquacrop"
+#         farmplot.cropstate.soiltype = soiltype
+#     else
+#         # do nothing for other cropmodels
+#     end
+# end
 
-"""
-    makecropstate(model, soiltype)
-
-An internal utility function to initialise one instance of the configured crop growth model.
-"""
-function makecropstate(model::SimulationModel, soiltype::SoilType)
-    if @param(crop.cropmodel) == "almass"
-        phase = (month(model.date) < 3 ? ALMaSS.janfirst : ALMaSS.marchfirst)
-        cs = ALMaSS.CropState(
-            croptype = model.crops["natural grass"],
-            phase = phase,
-        )
-    elseif @param(crop.cropmodel) == "simple"
-        cs = SimpleCrop.CropState(
-            model.crops["natural grass"],
-            0.0m
-        )
-    elseif @param(crop.cropmodel) == "aquacrop"
-        croptype = model.crops["natural grass"]
-        cs = AquaCropWrapper.CropState(croptype, soiltype, model)
-    else
-        Base.error("Unhandled crop model '$(@param(crop.cropmodel))' in makecropstate().")
-    end
-    return cs
-end
+# # TODO: move this function to individual crop model files
+# """
+#     makecropstate(model, soiltype)
+#
+# An internal utility function to initialise one instance of the configured crop growth model.
+# """
+# function makecropstate(model::SimulationModel, soiltype::SoilType)
+#     if @param(crop.cropmodel) == "almass"
+#         phase = (month(model.date) < 3 ? ALMaSS.janfirst : ALMaSS.marchfirst)
+#         cs = ALMaSS.CropState(
+#             croptype = model.crops["natural grass"],
+#             phase = phase,
+#         )
+#     elseif @param(crop.cropmodel) == "simple"
+#         cs = SimpleCrop.CropState(
+#             model.crops["natural grass"],
+#             0.0m
+#         )
+#     elseif @param(crop.cropmodel) == "aquacrop"
+#         croptype = model.crops["natural grass"]
+#         cs = AquaCropWrapper.CropState(croptype, soiltype, model)
+#     else
+#         Base.error("Unhandled crop model '$(@param(crop.cropmodel))' in makecropstate().")
+#     end
+#     return cs
+# end

src/crop/farmplot.jl

@@ -23,6 +23,12 @@ cropcover(f::FarmPlot) = cropcover(f.cropstate)
 cropyield(f::FarmPlot) = cropyield(f.cropstate)
 isharvestable(f::FarmPlot) = isharvestable(f.cropstate)
 
+function setsoiltype!(f::FarmPlot, soiltype::SoilType)
+    f.soiltype = soiltype
+    setsoiltype!(f.cropstate, soiltype)
+    return nothing
+end
+
 """
     stepagent!(farmplot, model)
 

src/crop/simplecrop.jl

@@ -4,9 +4,9 @@ using Persefone:
     AbstractCropState,
     AbstractCropType,
     AnnualDate,
-    FarmPlot,
     Length,
     cm,
+    SoilType,
     SimulationModel
 
 import Persefone:
@@ -16,6 +16,7 @@ import Persefone:
     cropheight,
     cropcover,
     cropyield,
+    makecropstate,
     sow!,
     harvest!,
     isharvestable
@@ -43,6 +44,9 @@ cropheight(cs::CropState) = cs.height
 cropcover(cs::CropState) = 0.0
 cropyield(cs::CropState) = 0.0  # TODO: units?
 isharvestable(cs::CropState) = true
+makecropstate(crop_type::CropType, model::SimulationModel, soiltype::SoilType) =
+    CropState(crop_type, 0.0cm)
+
 
 """
     stepagent!(cropstate, model)