Improve isharvestable() function for ALMaSS crop model

The data/crops/almass/crop_data_general.csv file now stores a column
`phase_before_harvest` to indicate in which growing phase the plant
can be harvested.  This is `marchfirst` for most crops, and `sow` for
maize.

In theory, with the current setup there could be a problem for
peas/beans if they are planted after March, 1 (the planting dates are
15 February to 15 March).  Planted before March 1st, the last growth
phase before harvest would be `marchfirst`.  But if planted after
March 1st, the last growth phase before harvest would be `sow`.  The
GDD curves for this crop only exist for `sow`, so probably planting
before March 1st should not be allowed for this crop.

data/crops/almass/crop_data_general.csv

-name,minsowdate,maxsowdate,minharvestdate,maxharvestdate,mingrowthtemp,group,biomass_scale,is_c4_plant,sowingdensity
-spring rape,NA,NA,NA,NA,NA,grain,1.071,false,75
-winter rape,20 August,25 August,NA,NA,NA,grain,1.071,false,75
-winter wheat,15 October,31 October,NA,NA,0,grain,1.0,false,400
-spring wheat,NA,NA,NA,NA,NA,grain,1.0,false,475
-winter barley,15 September,30 September,NA,NA,0,grain,0.857,false,340
-spring barley,1 March,10 April,NA,NA,0,grain,0.857,false,360
-undersown spring barley,NA,NA,NA,NA,0,grain,0.857,false,NA
-winter rye,23 September,15 October,NA,NA,NA,grain,0.857,false,290
-triticale,25 September,10 October,NA,NA,NA,grain,1.0,false,320
-oats,NA,NA,NA,NA,NA,grain,0.857,false,385
-maize,15 April,30 April,NA,NA,8,grain,1.0,true,13
-potatoes,NA,NA,NA,NA,4,root,0.857,false,5
-carrots,NA,NA,NA,NA,NA,root,0.7857,false,160
-beet,15 March,10 May,NA,NA,NA,root,0.857,false,16
-sunflower,25 March,15 April,NA,NA,NA,other,1.0,false,8
-lucerne,NA,NA,NA,NA,NA,legumes,1.2,false,NA
-peas/beans,15 February,15 March,NA,NA,5,legumes,0.857,false,35/80
-silage clover/grass,NA,NA,NA,NA,NA,legumes,1.1,false,NA
-fodder/clover,NA,NA,NA,NA,NA,legumes,1.2,false,NA
-lawn,NA,NA,NA,NA,NA,grass,0.5,false,NA
-permanent grassland (grazed),NA,NA,NA,NA,NA,grass,1.1,false,NA
-permanent grassland (seeded),NA,NA,NA,NA,NA,grass,1.1,false,NA
-permanent grassland (low yield),NA,NA,NA,NA,NA,grass,1.0,false,NA
-permanent set-aside,NA,NA,NA,NA,NA,semi-natural,0.7857,false,NA
-natural grass,NA,NA,NA,NA,NA,semi-natural,0.567,false,NA
-no growth,NA,NA,NA,NA,NA,semi-natural,0.0,false,NA
-heath,NA,NA,NA,NA,NA,semi-natural,0.567,false,NA
+name,minsowdate,maxsowdate,minharvestdate,maxharvestdate,mingrowthtemp,group,biomass_scale,is_c4_plant,sowingdensity,phase_before_harvest
+spring rape,NA,NA,NA,NA,NA,grain,1.071,false,75,marchfirst
+winter rape,20 August,25 August,NA,NA,NA,grain,1.071,false,75,marchfirst
+winter wheat,15 October,31 October,NA,NA,0,grain,1.0,false,400,marchfirst
+spring wheat,NA,NA,NA,NA,NA,grain,1.0,false,475,marchfirst
+winter barley,15 September,30 September,NA,NA,0,grain,0.857,false,340,marchfirst
+spring barley,1 March,10 April,NA,NA,0,grain,0.857,false,360,marchfirst
+undersown spring barley,NA,NA,NA,NA,0,grain,0.857,false,NA,marchfirst
+winter rye,23 September,15 October,NA,NA,NA,grain,0.857,false,290,marchfirst
+triticale,25 September,10 October,NA,NA,NA,grain,1.0,false,320,marchfirst
+oats,NA,NA,NA,NA,NA,grain,0.857,false,385,marchfirst
+maize,15 April,30 April,NA,NA,8,grain,1.0,true,13,sow
+potatoes,NA,NA,NA,NA,4,root,0.857,false,5,marchfirst
+carrots,NA,NA,NA,NA,NA,root,0.7857,false,160,marchfirst
+beet,15 March,10 May,NA,NA,NA,root,0.857,false,16,marchfirst
+sunflower,25 March,15 April,NA,NA,NA,other,1.0,false,8,marchfirst
+lucerne,NA,NA,NA,NA,NA,legumes,1.2,false,NA,marchfirst
+peas/beans,15 February,15 March,NA,NA,5,legumes,0.857,false,35/80,marchfirst
+silage clover/grass,NA,NA,NA,NA,NA,legumes,1.1,false,NA,marchfirst
+fodder/clover,NA,NA,NA,NA,NA,legumes,1.2,false,NA,marchfirst
+lawn,NA,NA,NA,NA,NA,grass,0.5,false,NA,marchfirst
+permanent grassland (grazed),NA,NA,NA,NA,NA,grass,1.1,false,NA,NA
+permanent grassland (seeded),NA,NA,NA,NA,NA,grass,1.1,false,NA,NA
+permanent grassland (low yield),NA,NA,NA,NA,NA,grass,1.0,false,NA,NA
+permanent set-aside,NA,NA,NA,NA,NA,semi-natural,0.7857,false,NA,NA
+natural grass,NA,NA,NA,NA,NA,semi-natural,0.567,false,NA,NA
+no growth,NA,NA,NA,NA,NA,semi-natural,0.0,false,NA,NA
+heath,NA,NA,NA,NA,NA,semi-natural,0.567,false,NA,NA

src/crop/almass.jl

@@ -95,6 +95,21 @@ struct CropCurveParams
     height::Dict{GrowthPhase, Vector{typeof(1.0cm)}}
 end
 
+function Base.show(io::IO, ::MIME"text/plain", curve::CropCurveParams)
+    println(io, "CropCurveParams")
+    println(io, "  curveID       = ", curve.curveID)
+    println(io, "  highnutrients = ", curve.highnutrients)
+    println(io, "\nGDD:")
+    show(io, MIME"text/plain"(), curve.GDD)
+    println(io, "\n\nLAItotal:")
+    show(io, MIME"text/plain"(), curve.LAItotal)
+    println(io, "\n\nLAIgreen:")
+    show(io, MIME"text/plain"(), curve.LAIgreen)
+    println(io, "\n\nheight:")
+    show(io, MIME"text/plain"(), curve.height)
+    return nothing
+end
+
 """
     CropType
 
@@ -107,6 +122,7 @@ struct CropType
     name::String
     group::String
     is_c4_plant::Bool  # false means it is a C3 plant
+    phase_before_harvest::Union{Missing,GrowthPhase}  # the phase where the crop can be harvested
     minsowdate::Union{Missing,AnnualDate}
     maxsowdate::Union{Missing,AnnualDate}
     minharvestdate::Union{Missing,AnnualDate}
@@ -117,6 +133,20 @@ struct CropType
     biomass_scale::Float64
 end
 
+function Base.show(io::IO, ::MIME"text/plain", ct::CropType)
+    println(io, "CropType")
+    println(io, "=========")
+    for f in filter(x -> ! (x in (:highnutrientgrowth, :lownutrientgrowth)), fieldnames(typeof(ct)))
+        println(io, "  $f = ", string(getfield(ct, f)))
+    end
+    println(io, "\n\nlownutrientgrowth:")
+    println(io, "--------------------")
+    show(io, MIME"text/plain"(), ct.lownutrientgrowth)
+    println(io, "\n\nhighnutrientgrowth:")
+    println(io, "--------------------")
+    show(io, MIME"text/plain"(), ct.highnutrientgrowth)
+end
+
 cropname(ct::CropType) = ct.name
 
 """
@@ -161,7 +191,13 @@ cropheight(cs::CropState) = cs.veg_height
 cropcover(cs::CropState) = cs.veg_cover
 cropyield(cs::CropState) = cs.veg_biomass  # TODO: correct? units? dry or wet?
 function isharvestable(cs::CropState)
-    return false
+    phase_before_harvest = croptype(cs).phase_before_harvest
+    if ismissing(phase_before_harvest) || cs.phase != phase_before_harvest
+        return false
+    end
+    curve = growthcurve(cs)
+    gdd = curve.GDD[phase_before_harvest]
+    return cs.ddegs >= last(gdd)
 end
 
 # Constant in original ALMaSS code:
@@ -363,11 +399,11 @@ Parse a CSV file containing the required parameter values for each crop
 """
 function readcropparameters(cropdirectory::String)
     @debug "Reading crop parameters"
-    # TODO: the last column (sowingdensity) uses `String` as type.
-    # This is because the entry for "peas/beans" has a value of
+    # TODO: the second-last column (sowingdensity) uses `String` as
+    # type.  This is because the entry for "peas/beans" has a value of
     # "35/80" for the sowingdensity.
     cropdata = CSV.File(joinpath(cropdirectory, CROPFILE), missingstring="NA",
-                        types=[String,AnnualDate,AnnualDate,AnnualDate,AnnualDate,Float64,String,Float64,Bool,String])
+                        types=[String,AnnualDate,AnnualDate,AnnualDate,AnnualDate,Float64,String,Float64,Bool,String,GrowthPhase])
     growthdata = CSV.File(joinpath(cropdirectory, GROWTHFILE), missingstring="NA",
                           types=[Int,String,String,GrowthPhase,String,
                                  Float64,Float64,Float64,Float64])
@@ -379,8 +415,8 @@ function readcropparameters(cropdirectory::String)
                                     filter(x -> x.nutrient_status=="high"))
         lownuts = buildgrowthcurve(cropgrowthdata |>
                                    filter(x -> x.nutrient_status=="low"))
-        croptypes[crop.name] = CropType(crop.name, crop.group, crop.is_c4_plant, crop.minsowdate,
-                                        crop.maxsowdate, crop.minharvestdate, crop.maxharvestdate,
+        croptypes[crop.name] = CropType(crop.name, crop.group, crop.is_c4_plant, crop.phase_before_harvest,
+                                        crop.minsowdate, crop.maxsowdate, crop.minharvestdate, crop.maxharvestdate,
                                         crop.mingrowthtemp, highnuts, lownuts, crop.biomass_scale)
     end
     croptypes

test/crop_tests.jl

@@ -19,7 +19,7 @@ const TESTPARAM_SIMPLECROP = joinpath(pkgdir(Persefone), "test", "test_parameter
 end
 
 @testset "Submodule ALMaSS" begin
-    ct = Ps.ALMaSS.CropType("olive tree", "no_crop_group", false,
+    ct = Ps.ALMaSS.CropType("olive tree", "no_crop_group", false, missing,
                             missing, missing, missing, missing, missing, missing, missing, 1.0)
     id = 0
     pixels = [(0, 0)]