diff --git a/CHANGELOG.md b/CHANGELOG.md index 83eb83c533fb1962623280fc730dbf303eb6796a..d78f731a36797f3806a8a8e73ca3123d350e02a3 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -24,6 +24,22 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 - added soil map section to GIS docs +- the `Weather` type now stores all weather information for the whole + simulation as a "struct-of-arrays" and with a function interface, + e.g. `sunshine(weather, date)` (previously, `Weather` was a struct + that stored the weather data for one day only, and + `AgricultureModel` had a field `weather::Dict{Date,Weather}`, making + it a "dict-of-struct" data layout) + +- `AgricultureModel` and `Weather` are now defined with `@kwdef`, + allowing for keyword arguments in their constructors + +- when reading weather data, we now throw an error when there are + missing values for the fields `min_temperature`, `max_temperature`, + `mean_temperature`, `precipitation`, and + `potential_evapotranspiration` (in the future missing values could + also be imputed) + ### Deprecated ### Removed diff --git a/src/core/simulation.jl b/src/core/simulation.jl index b19f4024a459ef835feb6e2eed0f636528a1a97e..e57f075b9b722ae78c920bd058b4a40d50413523 100644 --- a/src/core/simulation.jl +++ b/src/core/simulation.jl @@ -12,20 +12,20 @@ This is the heart of the model - a struct that holds all data and state for one simulation run. It is created by [`initialise`](@ref) and passed as input to most model functions. """ -mutable struct AgricultureModel{Tcroptype,Tcropstate} <: SimulationModel - settings::Dict{String,Any} +@kwdef mutable struct AgricultureModel{Tcroptype, Tcropstate} <: SimulationModel + settings = Dict{String, Any}() rng::AbstractRNG logger::AbstractLogger - dataoutputs::Dict{String,DataOutput} + dataoutputs = Dict{String, DataOutput}() date::Date - landscape::Matrix{Pixel} - weather::Dict{Date,Weather} - crops::Dict{String,Tcroptype} - farmers::Vector{Farmer} - farmplots::Vector{FarmPlot{Tcropstate}} - animals::Vector{Union{Animal,Nothing}} - migrants::Vector{Pair{Animal,AnnualDate}} - events::Vector{FarmEvent} + landscape = Matrix{Pixel}() + weather::Weather + crops = Dict{String, Tcroptype}() + farmers = Vector{Farmer}() + farmplots = Vector{FarmPlot{Tcropstate}}() + animals = Vector{Union{Animal, Nothing}}() + migrants = Vector{Pair{Animal, AnnualDate}}() + events = Vector{FarmEvent}() end """ @@ -36,7 +36,31 @@ Return the total number of agents in a model object. function nagents(model::AgricultureModel) length(model.animals)+length(model.farmers)+length(model.farmplots) end - + +windspeed(model::AgricultureModel, date::Date) = windspeed(model.weather, date) +windspeed(model::AgricultureModel) = windspeed(model, model.date) + +precipitation(model::AgricultureModel, date::Date) = precipitation(model.weather, date) +precipitation(model::AgricultureModel) = precipitation(model, model.date) + +sunshine(model::AgricultureModel, date::Date) = sunshine(model.weather, date) +sunshine(model::AgricultureModel) = sunshine(model, model.date) + +humidity(model::AgricultureModel, date::Date) = humidity(model.weather, date) +humidity(model::AgricultureModel) = humidity(model, model.date) + +meantemp(model::AgricultureModel, date::Date) = meantemp(model.weather, date) +meantemp(model::AgricultureModel) = meantemp(model, model.date) + +maxtemp(model::AgricultureModel, date::Date) = maxtemp(model.weather, date) +maxtemp(model::AgricultureModel) = maxtemp(model, model.date) + +mintemp(model::AgricultureModel, date::Date) = mintemp(model.weather, date) +mintemp(model::AgricultureModel) = mintemp(model, model.date) + +evapotranspiration(model::AgricultureModel, date::Date) = evapotranspiration(model.weather, date) +evapotranspiration(model::AgricultureModel) = evapotranspiration(model, model.date) + """ stepagent!(agent, model) @@ -117,22 +141,14 @@ function initmodel(settings::Dict{String, Any}) settings["core.enddate"]) crops, Tcroptype, Tcropstate = initcropmodel(settings["crop.cropmodel"], settings["crop.cropdirectory"]) - farmers = Vector{Farmer}() - farmplots = Vector{FarmPlot{Tcropstate}}() - model = AgricultureModel{Tcroptype,Tcropstate}( + model = AgricultureModel{Tcroptype, Tcropstate}(; settings, - StableRNG(settings["core.seed"]), + rng = StableRNG(settings["core.seed"]), logger, - Dict{String,DataOutput}(), - settings["core.startdate"], + date = settings["core.startdate"], landscape, weather, - crops, - farmers, - farmplots, - Vector{Union{Animal,Nothing}}(), - Vector{Pair{Animal, Date}}(), - Vector{FarmEvent}() + crops ) saveinputfiles(model) diff --git a/src/world/weather.jl b/src/world/weather.jl index 055d2e0115629438ee8477c5d03fa0a2448bcc42..12af77f2c4f4dea13e4a1af0d33b555dca320828 100644 --- a/src/world/weather.jl +++ b/src/world/weather.jl @@ -7,18 +7,81 @@ """ Weather -A single weather datum, combining the observations from one day. -""" -struct Weather - windspeed::Union{Missing,Float64} - precipitation::Union{Missing,Float64} - sunshine::Union{Missing,Float64} - cloudcover::Union{Missing,Float64} - humidity::Union{Missing,Float64} - meantemp::Union{Missing,Float64} - maxtemp::Union{Missing,Float64} - mintemp::Union{Missing,Float64} - evapotranspiration::Union{Missing,Float64} +Holds the weather information for the whole simulation period. +""" +@kwdef struct Weather + firstdate::Date + lastdate::Date + # TODO: types for some of these fields can be just `Vector{Float64}` after + # fixing missing measurement values (replacing with values + # approximated from previous/later measurements) + windspeed::Vector{Union{Missing,Float64}} + precipitation::Vector{Float64} + sunshine::Vector{Union{Missing,Float64}} + cloudcover::Vector{Union{Missing,Float64}} + humidity::Vector{Union{Missing,Float64}} + meantemp::Vector{Float64} + maxtemp::Vector{Float64} + mintemp::Vector{Float64} + evapotranspiration::Vector{Float64} +end + +function Base.length(weather::Weather) + return Dates.value(weather.lastdate - weather.firstdate) + 1 +end + +""" + daynumber(weather, date) + +Returns the number of days, counting `weather.firstdate` as day `1`. +""" +function daynumber(weather::Weather, date::Date) + @assert weather.firstdate <= date <= weather.lastdate + return Dates.value(date - weather.firstdate) + 1 +end + +""" + findspans(predicate_fn, array) -> Vector{Tuple{Int, Int}} + +Returns spans of indices in a 1-d `array` where a `predicate_fn` +returns `true`. The spans are returned as a `Vector{Tuple{Int, +Int}}`, where each tuple is of the form `(start_index, end_index)`. +""" +function findspans(predicate_fn, arr::AbstractVector) + spans = Tuple{Int, Int}[] + startidx = nothing + for i = firstindex(arr):lastindex(arr) + if predicate_fn(arr[i]) + if startidx === nothing + startidx = i + end + elseif startidx !== nothing + push!(spans, (startidx, i - 1)) + startidx = nothing + end + end + if startidx !== nothing + push!(spans, (startidx, lastindex(arr))) + end + return spans +end + +""" + check_missing_weatherdata(dataframe) + +Check the weather input data for missing values in columns where input +values are required. +""" +function check_missing_weatherdata(df::DataFrame) + colnames = ["min_temperature", "max_temperature", "mean_temperature", + "precipitation", "potential_evapotranspiration"] + for colname in colnames + spans = findspans(ismissing, getproperty(df, colname)) + # TODO: fix missing values and warn here instead of error + if ! isempty(spans) + error("Column $colname has missing values: $spans") + end + end end """ @@ -33,97 +96,103 @@ mapped to dates. """ function initweather(weatherfile::String, startdate::Date, enddate::Date) @debug "Initialising weather" - data = CSV.File(weatherfile, missingstring="NA", dateformat="yyyy-mm-dd", - types=[Date, Float64, Float64, Float64, Float64, Float64, - Float64, Float64, Float64, Float64]) - weather = Dict{Date,Weather}() - for row in data - if row.date >= startdate && row.date <= enddate - weather[row.date] = Weather(row.mean_windspeed, row.precipitation, - row.sunshine_hours, row.mean_cloud_cover, - row.mean_humidity, row.mean_temperature, - row.max_temperature, row.min_temperature, - row.potential_evapotranspiration) - end + if startdate > enddate + error("Startdate is after enddate: $startdate > $enddate.") end - if length(weather) <= Dates.value(enddate-startdate) - @warn ("There are missing days in the weather input file:" - * " expected $(Dates.value(enddate-startdate) + 1), got $(length(weather)).") + Tfloat = Union{Missing, Float64} + df = CSV.read(weatherfile, DataFrame; + missingstring="NA", dateformat="yyyy-mm-dd", + types=[Date, Tfloat, Tfloat, Tfloat, Tfloat, Tfloat, + Tfloat, Tfloat, Tfloat, Tfloat]) + needed_cols = [ + "date", "mean_windspeed", "precipitation", "sunshine_hours", + "mean_cloud_cover", "mean_humidity", "mean_temperature", + "max_temperature", "min_temperature", "potential_evapotranspiration" + ] + if ! issubset(needed_cols, names(df)) + error("Missing columns in weather data file: $(setdiff(needed_cols, names(df)))") end - weather + filter!(r -> r.date >= startdate && r.date <= enddate, df) + sort!(df, :date) + if startdate:enddate != df.date + error("There are missing days (rows) in the weather data file." + * " expected $(Dates.value(enddate - startdate) + 1), got $(nrow(df)))." + * " Missing dates are $(setdiff(startdate:enddate, df.date)).") + end + check_missing_weatherdata(df) + return Weather(; firstdate = minimum(df.date), + lastdate = maximum(df.date), + windspeed = df.mean_windspeed, + precipitation = df.precipitation, + sunshine = df.sunshine_hours, + cloudcover = df.mean_cloud_cover, + humidity = df.mean_humidity, + meantemp = df.mean_temperature, + maxtemp = df.max_temperature, + mintemp = df.min_temperature, + evapotranspiration = df.potential_evapotranspiration) end """ - windspeed(model) + windspeed(weather, date) -Return today's average windspeed in m/s. +Return the average windspeed in m/s on `date`. """ -function windspeed(model::SimulationModel) - model.weather[model.date].windspeed -end +windspeed(weather::Weather, date::Date) = + weather.windspeed[daynumber(weather, date)] """ - precipitation(model) + precipitation(weather, date) -Return today's total precipitation in mm. +Return the total precipitation in mm on `date`. """ -function precipitation(model::SimulationModel) - model.weather[model.date].precipitation -end +precipitation(weather::Weather, date::Date) = + weather.precipitation[daynumber(weather, date)] """ - sunshine(model) + sunshine(weather, date) -Return today's sunshine duration in hours. +Return the sunshine duration in hours on `date`. """ -function sunshine(model::SimulationModel) - model.weather[model.date].sunshine -end +sunshine(weather::Weather, date::Date) = + weather.sunshine[daynumber(weather, date)] """ - humidity(model) + humidity(weather, date) Return today's average vapour pressure in %. """ -function humidity(model::SimulationModel) - model.weather[model.date].humidity -end +humidity(weather::Weather, date::Date) = + weather.humidity[daynumber(weather, date)] """ - meantemp(model) + meantemp(weather, date) -Return today's mean temperature in °C. +Return the mean temperature in °C on `date`. """ -function meantemp(model::SimulationModel) - model.weather[model.date].meantemp -end +meantemp(weather::Weather, date::Date) = + weather.meantemp[daynumber(weather, date)] """ - maxtemp(model) + maxtemp(weather, date) -Return today's maximum temperature in °C. +Return the maximum temperature in °C on `date`. """ -function maxtemp(model::SimulationModel) - model.weather[model.date].maxtemp -end +maxtemp(weather::Weather, date::Date) = + weather.maxtemp[daynumber(weather, date)] """ - mintemp(model) + mintemp(weather, date) -Return today's minimum temperature in °C. +Return the minimum temperature in °C on `date`. """ -function mintemp(model::SimulationModel) - model.weather[model.date].mintemp -end +mintemp(weather::Weather, date::Date) = + weather.mintemp[daynumber(weather, date)] """ - evapotranspiration(model) + evapotranspiration(weather, date) -Return today's potential evapotranspiration (ETo), based on -the +Return today's potential evapotranspiration (ETo) on `date`. """ -function evapotranspiration(model::SimulationModel) - model.weather[model.date].evapotranspiration -end - -#TODO add functions for evapotranspiration +evapotranspiration(weather::Weather, date::Date) = + weather.evapotranspiration[daynumber(weather, date)] diff --git a/test/io_tests.jl b/test/io_tests.jl index 79d9fa103bb28e9dfb896220ef27680e325a8303..8aa70f38b4b61d5e7cd0e4a9df739541dced9e70 100644 --- a/test/io_tests.jl +++ b/test/io_tests.jl @@ -9,8 +9,8 @@ @test @param(core.startdate) == Date(2022, 2, 1) @test @param(core.loglevel) == "warn" @test @param(nature.targetspecies) == ["Wolpertinger", "Wyvern"] - @param(core.enddate) = Date(2022,1,3) - @test @param(core.enddate) == Date(2022,1,3) + @param(core.enddate) = Date(2022,2,3) + @test @param(core.enddate) == Date(2022,2,3) end @testset "Output functions" begin diff --git a/test/landscape_tests.jl b/test/landscape_tests.jl index beba14f22ef5bcb8a29e754aae8f058fa66c8ccd..7f1fe4266a15ca6651d0057510a0cf2229336c32 100644 --- a/test/landscape_tests.jl +++ b/test/landscape_tests.jl @@ -64,10 +64,10 @@ end @test Ps.distancetoedge((6,6), model) == 20m end -@testset "Weather initialisation" begin +@testset "Weather interface" begin # these tests are specific to the Jena weather file model = inittestmodel() - @test length(keys(model.weather)) == 59 + @test length(model.weather) == 59 @test ismissing(Ps.windspeed(model)) @test Ps.precipitation(model) == 1.3 @test ismissing(Ps.sunshine(model)) diff --git a/test/runtests.jl b/test/runtests.jl index f4401eb979260e76796ca17e722146e7c8ceea80..a655c01849151bd01c0fdffea176676d1c11ac1e 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -45,22 +45,14 @@ function inittestmodel(smallmap=true) TESTSETTINGS["core.enddate"]) # TODO: support other crop models besides ALMaSS crops = Ps.ALMaSS.readcropparameters(TESTSETTINGS["crop.cropdirectory"]) - farmers = Vector{Farmer}() - farmplots = Vector{Ps.FarmPlot{Ps.ALMaSS.CropState}}() - model = AgricultureModel{Ps.ALMaSS.CropType,Ps.ALMaSS.CropState}( - TESTSETTINGS, - StableRNG(TESTSETTINGS["core.seed"]), - global_logger(), - Dict{String,DataOutput}(), - TESTSETTINGS["core.startdate"], + model = AgricultureModel{Ps.ALMaSS.CropType,Ps.ALMaSS.CropState}(; + settings = copy(TESTSETTINGS), + rng = StableRNG(TESTSETTINGS["core.seed"]), + logger = global_logger(), + date = TESTSETTINGS["core.startdate"], landscape, weather, - crops, - farmers, - farmplots, - Vector{Union{Animal,Nothing}}(), - Vector{Pair{Animal, Date}}(), - Vector{FarmEvent}() + crops ) model end @@ -100,6 +92,9 @@ end include("landscape_tests.jl") include("simulation_tests.jl") end + @testset "Weather model" begin + include("weather_tests.jl") + end @testset "Nature model" begin include("nature_tests.jl") end diff --git a/test/weather_tests.jl b/test/weather_tests.jl new file mode 100644 index 0000000000000000000000000000000000000000..4f838c58b040b292f925f78d33fbd96e48f3bc85 --- /dev/null +++ b/test/weather_tests.jl @@ -0,0 +1,27 @@ + +@testset "Constructor and interface" begin + float3 = [0.0, 1.0, 2.0] + missing3 = Union{Missing,Float64}[missing, missing, missing] + weather = Weather(; firstdate = Date("2000-01-01"), + lastdate = Date("2000-01-03"), + windspeed = missing3, + precipitation = float3, + sunshine = missing3, + cloudcover = missing3, + humidity = missing3, + meantemp = float3, + maxtemp = float3, + mintemp = float3, + evapotranspiration = float3) + date = Date("2000-01-02") + @test length(weather) == 3 + @test Ps.daynumber(weather, date) == 2 + @test Ps.windspeed(weather, date) |> ismissing + @test Ps.precipitation(weather, date) == 1.0 + @test Ps.sunshine(weather, date) |> ismissing + @test Ps.humidity(weather, date) |> ismissing + @test Ps.meantemp(weather, date) == 1.0 + @test Ps.maxtemp(weather, date) == 1.0 + @test Ps.mintemp(weather, date) == 1.0 + @test Ps.evapotranspiration(weather, date) == 1.0 +end