diff --git a/src/Persefone.jl b/src/Persefone.jl index 0517e45808ceb5ad8ea46d48704378426b9fd7d8..40ffdef208e10219b35cf40d510a1f056f8ebf10 100644 --- a/src/Persefone.jl +++ b/src/Persefone.jl @@ -100,7 +100,10 @@ export populationtrends, visualiseoutput, savemodelobject, - loadmodelobject + loadmodelobject, + croptype, + cropname, + cropheight ## Import and define units and dimensions import Unitful: cm, m, km, ha, mg, g, kg, Length, Area, Mass @@ -138,9 +141,11 @@ include("analysis/makieplots.jl") include("world/landscape.jl") include("world/weather.jl") +include("crop/croptypes.jl") include("crop/farmplot.jl") include("crop/almass.jl") include("crop/simplecrop.jl") + include("farm/farm.jl") include("nature/insects.jl") diff --git a/src/core/simulation.jl b/src/core/simulation.jl index ed8b2dafd9e8296001085b95a01efb358ca849c2..8ed33f1db1923e160a4bae9455ddb12141d84e93 100644 --- a/src/core/simulation.jl +++ b/src/core/simulation.jl @@ -21,7 +21,7 @@ mutable struct AgricultureModel <: SimulationModel date::Date landscape::Matrix{Pixel} weather::Dict{Date,Weather} - crops::Dict{String,ALMaSS.CropType} + crops::Dict{String,AbstractCropType} farmers::Vector{Farmer} farmplots::Vector{AbstractFarmPlot} animals::Vector{Union{Animal,Nothing}} diff --git a/src/crop/almass.jl b/src/crop/almass.jl index 4c8f708c5663b3b259b251ab4dead07170563a81..5bac605353eab1612cdb85cbec947778c448eef2 100644 --- a/src/crop/almass.jl +++ b/src/crop/almass.jl @@ -8,14 +8,20 @@ module ALMaSS using Persefone: - AbstractFarmPlot, + AbstractCropType, + AbstractCropState, EventType, SimulationModel, fertiliser, initfields_fill_with!, maxtemp, mintemp -import Persefone: stepagent! +import Persefone: + Persefone, + stepagent!, + croptype, + cropname, + cropheight using Dates: Date, month, monthday using CSV: CSV @@ -48,7 +54,7 @@ end The type struct for all crops. Currently follows the crop growth model as implemented in ALMaSS. """ -struct CropType +struct CropType <: AbstractCropType #TODO make this into an abstract type and create subtypes for different # crop submodels (#70) name::String @@ -62,6 +68,32 @@ struct CropType #issowable::Union{Function,Bool} end +cropname(ct::CropType) = ct.name + +""" + CropState + +The state data for an ALMaSS vegetation point calculation as used in +FarmPlot. +""" +mutable struct CropState <: AbstractCropState + #TODO add Unitful + croptype::CropType + phase::GrowthPhase + growingdegreedays::Float64 + height::Float64 + LAItotal::Float64 + LAIgreen::Float64 + #biomass::Float64 #XXX I need to figure out how to calculate this + events::Vector{EventType} +end + +croptype(cs::CropState) = cs.croptype +cropname(cs::CropState) = cropname(croptype(cs)) +cropheight(cs::CropState) = cs.height + +const FarmPlot = Persefone.FarmPlot{CropState} + """ Base.tryparse(type, str) @@ -134,29 +166,6 @@ function readcropparameters(generalcropfile::String, growthfile::String) croptypes end - -#XXX not sure whether it makes sense to have this as an agent type, -# or perhaps better a grid property? -""" - FarmPlot - -This represents one field, i.e. a collection of pixels with the same management. -This is the spatial unit with which the crop growth model and the farm model work. -""" -mutable struct FarmPlot <: AbstractFarmPlot - #TODO add Unitful - const id::Int64 - pixels::Vector{Tuple{Int64, Int64}} - croptype::ALMaSS.CropType - phase::ALMaSS.GrowthPhase - growingdegreedays::Float64 - height::Float64 - LAItotal::Float64 - LAIgreen::Float64 - #biomass::Float64 #XXX I need to figure out how to calculate this - events::Vector{EventType} -end - """ initfields!(model) @@ -165,11 +174,11 @@ Initialise the model with its farm plots. function initfields!(model::SimulationModel) initfields_fill_with!(model) do model, x, y month(model.date) < 3 ? phase = ALMaSS.janfirst : phase = ALMaSS.marchfirst - FarmPlot(length(model.farmplots) + 1, - [(x,y)], - model.crops["natural grass"], - phase, - 0.0, 0.0, 0.0, 0.0, Vector{EventType}()) + FarmPlot(length(model.farmplots) + 1, [(x,y)], + CropState(model.crops["natural grass"], + phase, + 0.0, 0.0, 0.0, 0.0, Vector{EventType}()) + ) end end @@ -182,13 +191,14 @@ function stepagent!(farmplot::FarmPlot, model::SimulationModel) # update growing degree days # if no crop-specific base temperature is given, default to 5°C # (https://www.eea.europa.eu/publications/europes-changing-climate-hazards-1/heat-and-cold/heat-and-cold-2014-mean) - basetemp = farmplot.croptype.mingrowthtemp + cs = farmplot.crop_state + basetemp = cs.croptype.mingrowthtemp ismissing(basetemp) && (basetemp = 5.0) gdd = (maxtemp(model)+mintemp(model))/2 - basetemp - gdd > 0 && (farmplot.growingdegreedays += gdd) + gdd > 0 && (cs.growingdegreedays += gdd) # update the phase on key dates - monthday(model.date) == (1,1) && (farmplot.phase = ALMaSS.janfirst) - monthday(model.date) == (3,1) && (farmplot.phase = ALMaSS.marchfirst) + monthday(model.date) == (1,1) && (cs.phase = ALMaSS.janfirst) + monthday(model.date) == (3,1) && (cs.phase = ALMaSS.marchfirst) # update crop growth growcrop!(farmplot, model) end @@ -202,9 +212,10 @@ end Sow the specified crop on this farmplot. """ function sow!(cropname::String, farmplot::FarmPlot, model::SimulationModel) + cs = farmplot.crop_state createevent!(model, farmplot.pixels, sowing) - farmplot.croptype = model.crops[cropname] - farmplot.phase = ALMaSS.sow + cs.croptype = model.crops[cropname] + cs.phase = ALMaSS.sow #XXX test if the crop is sowable? end @@ -214,10 +225,11 @@ end Harvest the crop on this farmplot. """ function harvest!(farmplot::FarmPlot, model::SimulationModel) + cs = farmplot.crop_state createevent!(model, farmplot.pixels, harvesting) - farmplot.phase in [ALMaSS.harvest1, ALMaSS.harvest2] ? - farmplot.phase = ALMaSS.harvest2 : - farmplot.phase = ALMaSS.harvest1 + cs.phase in [ALMaSS.harvest1, ALMaSS.harvest2] ? + cs.phase = ALMaSS.harvest2 : + cs.phase = ALMaSS.harvest1 # height & LAI will be automatically adjusted by the growth function #TODO calculate and return yield end @@ -233,25 +245,26 @@ Apply the relevant crop growth model to update the plants on this farm plot. Currently only supports the ALMaSS crop growth model by Topping et al. """ function growcrop!(farmplot::FarmPlot, model::SimulationModel) - fertiliser in farmplot.events ? - curve = farmplot.croptype.lownutrientgrowth : - curve = farmplot.croptype.highnutrientgrowth - points = curve.GDD[farmplot.phase] + cs = farmplot.crop_state + fertiliser in cs.events ? + curve = cs.croptype.lownutrientgrowth : + curve = cs.croptype.highnutrientgrowth + points = curve.GDD[cs.phase] for p in 1:length(points) if points[p] == 99999 return # the marker that there is no further growth this phase elseif points[p] == -1 # the marker to set all variables to specified values - farmplot.height = curve.height[farmplot.phase][p] - farmplot.LAItotal = curve.LAItotal[farmplot.phase][p] - farmplot.LAIgreen = curve.LAIgreen[farmplot.phase][p] + cs.height = curve.height[cs.phase][p] + cs.LAItotal = curve.LAItotal[cs.phase][p] + cs.LAIgreen = curve.LAIgreen[cs.phase][p] return else - gdd = farmplot.growingdegreedays + gdd = cs.growingdegreedays # figure out which is the correct slope value to use for growth if p == length(points) || gdd < points[p+1] - farmplot.height += curve.height[farmplot.phase][p] - farmplot.LAItotal += curve.LAItotal[farmplot.phase][p] - farmplot.LAIgreen += curve.LAIgreen[farmplot.phase][p] + cs.height += curve.height[cs.phase][p] + cs.LAItotal += curve.LAItotal[cs.phase][p] + cs.LAIgreen += curve.LAIgreen[cs.phase][p] return end #XXX To be precise, we ought to figure out if one or more inflection @@ -263,53 +276,4 @@ function growcrop!(farmplot::FarmPlot, model::SimulationModel) end end - -## UTILITY FUNCTIONS - -""" - averagefieldsize(model) - -Calculate the average field size in hectares for the model landscape. -""" -function averagefieldsize(model::SimulationModel) - conversionfactor = 100 #our pixels are currently 10x10m, so 100 pixels per hectare - sizes::Vector{Float64} = [] - for fp in model.farmplots - push!(sizes, size(fp.pixels)[1]/conversionfactor) - end - round(sum(sizes)/size(sizes)[1], digits=2) -end - -""" - croptype(model, position) - -Return the crop at this position, or nothing if there is no crop here (utility wrapper). -""" -function croptype(pos::Tuple{Int64,Int64}, model::SimulationModel) - ismissing(model.landscape[pos...].fieldid) ? nothing : - model.farmplots[model.landscape[pos...].fieldid].croptype -end - -""" - cropname(model, position) - -Return the name of the crop at this position, or an empty string if there is no crop here -(utility wrapper). -""" -function cropname(pos::Tuple{Int64,Int64}, model::SimulationModel) - field = model.landscape[pos...].fieldid - ismissing(field) ? "" : model.farmplots[field].croptype.name -end - -""" - cropheight(model, position) - -Return the height of the crop at this position, or nothing if there is no crop here -(utility wrapper). -""" -function cropheight(pos::Tuple{Int64,Int64}, model::SimulationModel) - ismissing(model.landscape[pos...].fieldid) ? nothing : - model.farmplots[model.landscape[pos...].fieldid].height -end - end # module ALMaSS diff --git a/src/crop/croptypes.jl b/src/crop/croptypes.jl new file mode 100644 index 0000000000000000000000000000000000000000..3665322a90d17a9c7da6afd2aa855f2dcedbf53e --- /dev/null +++ b/src/crop/croptypes.jl @@ -0,0 +1,13 @@ +### Persefone.jl - a model of agricultural landscapes and ecosystems in Europe. +### +### This file contains abstract types that crop models must derive from. +### + +abstract type AbstractCropType end +cropname(t::AbstractCropType) = t.name + +abstract type AbstractCropState end +croptype(s::AbstractCropState) = s.croptype +cropname(s::AbstractCropState) = cropname(cropytype(s)) +# cropheight(s::AbstractCropState) = s.cropheight +# cropcover(s::AbstractCropState) = s.cropcover diff --git a/src/crop/farmplot.jl b/src/crop/farmplot.jl index 4272e05971a1d5d37fe0d25ac2d6aae476ab4a30..b570df0c3271a7c7d746561cfa550939289767aa 100644 --- a/src/crop/farmplot.jl +++ b/src/crop/farmplot.jl @@ -214,6 +214,16 @@ end +mutable struct FarmPlot{T <: AbstractCropState} <: AbstractFarmPlot + const id::Int64 + pixels::Vector{Tuple{Int64, Int64}} + crop_state :: T +end + +croptype(f::FarmPlot{T}) where {T} = croptype(f.crop_state) +cropname(f::FarmPlot{T}) where {T} = cropname(croptype(f)) +cropheight(f::FarmPlot{T}) where {T} = cropheight(f.crop_state) + """ initfields_fill_with!(make_farmplot_fn, model) @@ -246,3 +256,52 @@ function initfields_fill_with!(make_farmplot_fn::Function, model::SimulationMode end @info "Initialised $n farm plots." end + + +## UTILITY FUNCTIONS + +""" + averagefieldsize(model) + +Calculate the average field size in hectares for the model landscape. +""" +function averagefieldsize(model::SimulationModel) + conversionfactor = 100 #our pixels are currently 10x10m, so 100 pixels per hectare + sizes::Vector{Float64} = [] + for fp in model.farmplots + push!(sizes, size(fp.pixels)[1]/conversionfactor) + end + round(sum(sizes)/size(sizes)[1], digits=2) +end + +""" + croptype(model, position) + +Return the crop at this position, or nothing if there is no crop here (utility wrapper). +""" +function croptype(pos::Tuple{Int64,Int64}, model::SimulationModel) + ismissing(model.landscape[pos...].fieldid) ? nothing : + croptype(model.farmplots[model.landscape[pos...].fieldid]) +end + +""" + cropname(model, position) + +Return the name of the crop at this position, or an empty string if there is no crop here +(utility wrapper). +""" +function cropname(pos::Tuple{Int64,Int64}, model::SimulationModel) + field = model.landscape[pos...].fieldid + ismissing(field) ? "" : cropname(model.farmplots[field]) +end + +""" + cropheight(model, position) + +Return the height of the crop at this position, or nothing if there is no crop here +(utility wrapper). +""" +function cropheight(pos::Tuple{Int64,Int64}, model::SimulationModel) + ismissing(model.landscape[pos...].fieldid) ? nothing : + cropheight(model.farmplots[model.landscape[pos...].fieldid]) +end diff --git a/src/crop/simplecrop.jl b/src/crop/simplecrop.jl index d96cf303e3eb750d6c9f77ccddb32341abffe9fb..9259324b5451eb697520252276fc65f4808d0872 100644 --- a/src/crop/simplecrop.jl +++ b/src/crop/simplecrop.jl @@ -1,13 +1,34 @@ module SimpleCrop -using Persefone: AbstractFarmPlot, SimulationModel, initfields_fill_with! -import Persefone: stepagent! +using Persefone: + AbstractCropType, + AbstractCropState, + SimulationModel, + initfields_fill_with! +import Persefone: + Persefone, + stepagent!, + croptype, + cropname, + cropheight -mutable struct FarmPlot <: AbstractFarmPlot - const id::Int64 - pixels::Vector{Tuple{Int64, Int64}} +struct CropType <: AbstractCropType + name::String end +cropname(ct::CropType) = ct.name + +mutable struct CropState <: AbstractCropState + croptype::CropType + cropheight::Float64 +end + +croptype(cs::CropState) = cs.croptype +cropname(cs::CropState) = cropname(croptype(cs)) +cropheight(cs::CropState) = cs.height + +const FarmPlot = Persefone.FarmPlot{CropState} + """ stepagent!(farmplot, model) @@ -24,7 +45,7 @@ Initialise the model with its farm plots. """ function initfields!(model::SimulationModel) initfields_fill_with!(model) do model, x, y - FarmPlot(length(model.farmplots) + 1, [(x,y)]) + FarmPlot(length(model.farmplots) + 1, [(x,y)], CropState(CropType("nothing"), 0.0)) end end diff --git a/test/landscape_tests.jl b/test/landscape_tests.jl index 9fcc4ec7fd7f8e884e59da231a2b6006235b6984..fc6ae84b6b0f5ca8adaa411b0cbfd9031dbd8f19 100644 --- a/test/landscape_tests.jl +++ b/test/landscape_tests.jl @@ -14,7 +14,7 @@ @test Ps.landcover((400,400), model) == Ps.grass @test Ps.landcover((800,800), model) == Ps.agriculture @test Ps.landcover((1100,1100), model) == Ps.builtup - @test Ps.ALMaSS.averagefieldsize(model) == 5.37 + @test Ps.averagefieldsize(model) == 5.37 @test count(f -> ismissing(f.fieldid), model.landscape) == 1685573 @test length(Ps.farmplot((800,800), model).pixels) == 4049 end diff --git a/test/nature_tests.jl b/test/nature_tests.jl index 027fed680d4fef444553bc319cbbaecb000a6acb..1e277d80ef60b601ba0eb6865fa12a72dbc0beb1 100644 --- a/test/nature_tests.jl +++ b/test/nature_tests.jl @@ -47,10 +47,15 @@ end) # end eval @testset "Habitat macros" begin # set up the testing landscape model = inittestmodel() - model.landscape[6,6] = Pixel(Ps.agriculture, 1) - push!(model.farmplots, - FarmPlot(1, [(6,6)], model.crops["winter wheat"], Ps.ALMaSS.janfirst, - 0.0, 0.0, 0.0, 0.0, Vector{Ps.Management}())) + model.landscape[6,6] = Pixel(Ps.agriculture, 1, [], []) + fp = Ps.FarmPlot( + 1, [(6,6)], + Ps.ALMaSS.CropState( + model.crops["winter wheat"], Ps.ALMaSS.janfirst, + 0.0, 0.0, 0.0, 0.0, Vector{Ps.Management}() + ) + ) + push!(model.farmplots, fp) push!(model.animals, Ps.Mermaid(1, Ps.male, (-1,-1), (3,3), Ps.life), Ps.Mermaid(2, Ps.female, (-1,-1), (4,4), Ps.life))