Map values are now read into a single landscape matrix

docs/gis.md

@@ -68,5 +68,6 @@ created above - this ensures the two layers match) as the output extent. Make su
 the "fixed value to burn" is "Not set". Then choose "Georeferenced units" as the 
 "Out raster size units" and set horizontal and vertical resolution to 10.0. In 
 the advanced parameters, set the output data type to `UInt32`. Finally, enter an 
-output file name and run.
+output file name and run. The resulting TIF file can be passed to Persephone
+as the `farmfieldmap` parameter.
 

src/Persephone.jl

@@ -38,6 +38,7 @@ const PARAMFILE = "src/parameters.toml"
 ## (if file b references something from file a, it must be included later)
 include("core/input.jl")
 include("core/output.jl")
+include("core/landscape.jl")
 include("farm/farm.jl")
 include("crop/crops.jl")
 include("nature/nature.jl")

src/core/landscape.jl

+### Persephone - a socio-economic-ecological model of European agricultural landscapes.
+###
+### This file manages the landscape maps that underlie the model.
+###
+
+## The land cover classes encoded in the Mundialis Sentinel data.
+## Do not change the order of this enum, or initlandscape() will break!
+@enum LandCover nodata forest grass water builtup soil agriculture
+
+"""
+    Pixel
+
+A pixel is a simple data structure to combine land use and ownership information
+in a single object. The model landscape consists of a matrix of pixels.
+(Note: further landscape information may be added here in future.)
+"""
+struct Pixel
+    landcover::LandCover
+    fieldid::Union{Missing, UInt32}
+end
+
+"""
+    initlandscape()
+
+Initialise the model landscape based on the map files specified in the
+configuration. Returns a matrix of pixels.
+"""
+function initlandscape()
+    @debug "Initialising landscape"
+    landcover = GeoArrays.read(param("core.landcovermap"))
+    farmfields = GeoArrays.read(param("core.farmfieldsmap"))
+    (size(landcover) != size(farmfields)) && Base.error("Input map sizes don't match.")
+    width, height = size(landcover)[1:2]
+    landscape = Matrix{Pixel}(undef, width, height)
+    for x in 1:width
+        for y in 1:height
+            # convert the numeric landcover value to LandCover, then create the Pixel object
+            lc = landcover[x,y][1]
+            (ismissing(lc)) && (lc = 0)
+            lcv = LandCover(Int(lc/10))
+            landscape[x,y] = Pixel(lcv, farmfields[x,y][1])
+        end
+    end
+    return landscape
+end
+
+
+"""
+    landcover(model, position)
+
+Return the land cover class at this position (utility wrapper).
+"""
+function landcover(model::AgentBasedModel, pos::Tuple{Int64,Int64})
+    model.landscape[pos...].landcover
+end
+
+"""
+    fieldid(model, position)
+
+Return the UID of the field at this position (utility wrapper).
+"""
+function fieldid(model::AgentBasedModel, pos::Tuple{Int64,Int64})
+    model.landscape[pos...].fieldid
+end

src/core/simulation.jl

@@ -15,13 +15,10 @@ function initialise(config::String=PARAMFILE)
     Random.seed!(param("core.seed"))
     setupdatadir()
     # initialise world-level properties
-    landcover = GeoArrays.read(param("core.landcovermap"))
-    farmfields = GeoArrays.read(param("core.farmfieldsmap"))
-    space = GridSpace(size(landcover)[1:2], periodic=false)
-    properties = Dict{Symbol,Any}(:age=>0,
-                                  :date=>param("core.startdate"),
-                                  :landcover=>landcover,
-                                  :farmfields=>farmfields)
+    landscape = initlandscape()
+    space = GridSpace(size(landscape), periodic=false)
+    properties = Dict{Symbol,Any}(:date=>param("core.startdate"),
+                                  :landscape=>landscape)
     model = AgentBasedModel(Union{Farmer,Animal,CropPlot}, space, properties=properties,
                             rng=Random.Xoshiro(param("core.seed")))
     # initialise submodels
@@ -38,13 +35,12 @@ end
 Execute one update of the model.
 """
 function stepsimulation!(model::AgentBasedModel)
-    model.age += 1
-    model.date += Day(1)
-    @info "Simulating day $(model.date)."
+    @info "Simulating day $(model.date+Day(1))."
     for a in Schedulers.ByType((Farmer,Animal,CropPlot), true)(model)
         stepagent!(getindex(model, a), model)
     end
     savepopulationdata(model)
+    model.date += Day(1)
 end
 
 """
@@ -53,6 +49,7 @@ end
 Wrap up the simulation. Output all remaining data and exit.
 """
 function finalise(model::AgentBasedModel)
+    @info "Simulated $(model.date-param("core.startdate"))."
     @info "Simulation completed at $(Dates.now()),\nwrote output to $(param("core.outdir"))."
     #TODO
     genocide!(model)
@@ -65,7 +62,7 @@ Carry out a complete simulation run.
 """
 function simulate(config::String=PARAMFILE)
     model = initialise(config)
-    runtime = Dates.value(param("core.enddate")-param("core.startdate"))
+    runtime = Dates.value(param("core.enddate")-param("core.startdate"))+1
     step!(model, dummystep, stepsimulation!, runtime)
     finalise(model)
 end

src/nature/nature.jl

@@ -6,10 +6,6 @@
 ## An enum used to assign a sex to each animal
 @enum Sex hermaphrodite male female
 
-## The land cover classes encoded in the Mundialis Sentinel data.
-## Do not change the order of this enum, or landcover() will break!
-@enum LandCover nodata forest grass water builtup soil agriculture
-
 """
     Species
 
@@ -84,15 +80,3 @@ function initnature!(model::AgentBasedModel)
         getspecies(s).initpop!(model)
     end
 end
-
-"""
-    landcover(model, position)
-
-Return the land cover class at this position.
-"""
-function landcover(model::AgentBasedModel, pos::Tuple{Int64,Int64})
-    #XXX It's probably worth converting the entire GeoArray to a LandCover array
-    lc = model.landcover[pos...][1]
-    (ismissing(lc)) && (return(nodata))
-    return LandCover(Int(lc/10))
-end

src/nature/wolpertinger.jl

@@ -12,7 +12,7 @@ Initialise a population of Wolpertingers in random locations around the landscap
 """
 function initwolpertinger!(model::AgentBasedModel)
     species = getspecies("Wolpertinger")
-    x, y = size(model.landcover)
+    x, y = size(model.landscape)
     popsize = Int(round((x*y)/10000))
     for i in 1:popsize
         add_agent!(Animal, model, species, hermaphrodite, 0, 100)