Improve ALMaSS crop growth model

- use a `setphase!` function to change vegetative growth phases

- winter crop varieties still are not being reset in size after
  harvest, the other tested crops now look reasonable

- the logic in the usage of the growth curves still is slightly wrong
  and will be improved in the next patches

src/crop/almass.jl

@@ -206,6 +206,40 @@ function readcropparameters(generalcropfile::String, growthfile::String)
     croptypes
 end
 
+function setphase!(cs::CropState, phase::GrowthPhase)
+    cs.phase = phase
+    # TODO: should cs.growingdegreedays be set here?
+    # cs.growingdegreedays = 0.0
+
+    fertiliser in cs.events ?
+        curve = cs.croptype.lownutrientgrowth :
+        curve = cs.croptype.highnutrientgrowth
+    if ismissing(curve)
+        error("Missing curve = $curve")
+    end
+    GDD = curve.GDD[cs.phase]
+    height = curve.height[cs.phase]
+    LAItotal = curve.LAItotal[cs.phase]
+    LAIgreen = curve.LAIgreen[cs.phase]
+
+    if (length(GDD) == 0
+        || length(height) == 0
+        || length(LAItotal) == 0
+        || length(LAIgreen) == 0)
+        error("One of these has length zero:\n  GDD = $GDD\n  height = $height\n  LAItotal = $LAItotal\n  LAIgreen = $LAIgreen")
+    end
+    if GDD[1] == -1.0
+        # set height, LAItotal, LAIgreen
+        # TODO: what about when
+        #           GDD[1] == curve.GDD[cs.phase][1] == 99999.0
+        #       ???
+        # TODO: cropname(cs) == "no growth" has curve.GDD[harvest1] == 99999.0
+        cs.height = height[1]
+        cs.LAItotal = LAItotal[1]
+        cs.LAIgreen = LAIgreen[1]
+    end
+end
+
 """
     stepagent!(cropstate, model)
 
@@ -214,10 +248,10 @@ Update a farm plot by one day.
 function stepagent!(cs::CropState, model::SimulationModel)
     # update the phase on key dates
     if monthday(model.date) == (1, 1)
-        cs.phase = ALMaSS.janfirst
+        setphase!(cs, ALMaSS.janfirst)
         cs.growingdegreedays = 0.0
     elseif monthday(model.date) == (3, 1)
-        cs.phase = ALMaSS.marchfirst
+        setphase!(cs, ALMaSS.marchfirst)
         cs.growingdegreedays = 0.0
     end
     # update growing degree days
@@ -243,13 +277,17 @@ function sow!(cs::CropState, model::SimulationModel, cropname::String)
     !ismissing(cs.croptype.minsowdate) && model.date < cs.croptype.minsowdate &&
         @warn "$(model.date) is earlier than the minimum sowing date for $(cropname)."
     cs.croptype = model.crops[cropname]
-    cs.phase = ALMaSS.sow
+
+    setphase!(cs, ALMaSS.sow)
     cs.growingdegreedays = 0.0
-    # TODO: set from crop curve params
-    #       cs.height = cs.croptype.lownutrientgrowth.height[sow]  # or highnutrientgrowth
-    cs.height = 0.0cm
-    cs.LAItotal = 0.0
-    cs.LAIgreen = 0.0
+    # cs.phase = ALMaSS.sow
+    # cs.growingdegreedays = 0.0
+    # # TODO: set from crop curve params
+    # #       cs.height = cs.croptype.lownutrientgrowth.height[sow]  # or highnutrientgrowth
+    # cs.height = 0.0cm
+    # cs.LAItotal = 0.0
+    # cs.LAIgreen = 0.0
+
     cs.mature = false
     cs.events = Vector{Management}()
 end
@@ -260,31 +298,12 @@ end
 Harvest the crop of this cropstate.
 """
 function harvest!(cs::CropState, model::SimulationModel)
-    cs.phase in [ALMaSS.harvest1, ALMaSS.harvest2] ?
-        cs.phase = ALMaSS.harvest2 :
-        cs.phase = ALMaSS.harvest1
-    cs.growingdegreedays = 0.0
-    cs.mature = false
-
-    # adjust height, LAItotal, LAIgreen
-    fertiliser in cs.events ?
-        curve = cs.croptype.lownutrientgrowth :
-        curve = cs.croptype.highnutrientgrowth
-    if ismissing(curve)
-        error("Missing curve = $curve")
+    if cs.phase in (ALMaSS.harvest1, ALMaSS.harvest2)
+        setphase!(cs, ALMaSS.harvest2)
+    else
+        setphase!(cs, ALMaSS.harvest1)
     end
-    # TODO: cropname == "no growth" has curve.GDD[harvest1] == 99999.0
-    if length(curve.GDD[cs.phase]) == 0 || !(curve.GDD[cs.phase][1] == -1.0 || curve.GDD[cs.phase][1] == 99999.0)
-        error("Bad crop curve GDD\n  for phase == $(cs.phase), croptype = $(cs.croptype)):\n  $(curve.GDD)")
-    end
-    if (length(curve.height[cs.phase]) == 0
-        || length(curve.LAItotal[cs.phase]) == 0
-        || length(curve.LAIgreen[cs.phase]) == 0)
-        error("One of these has length zero:\n  $curve.height\n  $curve.LAItotal\n  $curve.LAIgreen")
-    end
-    cs.height = curve.height[cs.phase][1]
-    cs.LAItotal = curve.LAItotal[cs.phase][1]
-    cs.LAIgreen = curve.LAIgreen[cs.phase][1]
+    cs.mature = false
 
     #TODO calculate and return yield
 end
@@ -325,14 +344,14 @@ function growcrop!(cs::CropState, delta_gdd::Real, model::SimulationModel)
     idx = findfirst(x -> x > total_gdd, points)
     if isnothing(idx)
         idx = lastindex(points)
-    elseif idx != firstindex(points)
-        idx -= 1  # we need the last index that is not > total_gdd
+    # elseif idx != firstindex(points)
+    #     idx -= 1  # we need the last index that is not > total_gdd
     end
     cs.height += curve.height[cs.phase][idx] * delta_gdd
     cs.LAItotal += curve.LAItotal[cs.phase][idx] * delta_gdd
     cs.LAIgreen += curve.LAIgreen[cs.phase][idx] * delta_gdd
 
-    if cs.phase == janfirst && length(points
+    #if cs.phase == janfirst && length(points
 end
 
 end # module ALMaSS