From 35b1dbb57e37b43c94c576fcbc0c550137fc2c7d Mon Sep 17 00:00:00 2001
From: remic <remic@iri.columbia.edu>
Date: Tue, 17 Feb 2026 16:57:44 -0500
Subject: [PATCH] length of frost season

---
 pepsico/app_calc.py         | 32 +++++++++++++++++++++++++-------
 pepsico/proj_wwc/layout.py  | 12 +++++++++---
 pepsico/proj_wwc/maproom.py |  1 +
 3 files changed, 35 insertions(+), 10 deletions(-)

diff --git a/pepsico/app_calc.py b/pepsico/app_calc.py
index 1e1a86e3..169a54b5 100644
--- a/pepsico/app_calc.py
+++ b/pepsico/app_calc.py
@@ -143,17 +143,31 @@ def seasonal_wwc(
                 norm.ppf, quantile, kwargs={"loc": data_ds[var], "scale": data_std[var]},
             )
         wwc_units = "˚C"
-    # This option is also commented out as it didn't work in its present form. Didn't
-    # really expect that it would though.
     if variable == "frost_season_length":
-        data_ds = (labelled_season_data <= frost_threshold).groupby(
-            labelled_season_data["seasons_starts"]
+        # This is the duration from first to last frost days
+        # Maybe I should make a function out of this?
+        data_ds = (
+            (labelled_season_data <= frost_threshold)
+            .where(~np.isnan(labelled_season_data))
+            .groupby(labelled_season_data["seasons_starts"])
         )
         data_ds = (
-            data_ds[-1:0].idxmax()
-            - data_ds.idxmax()
+            # cumsum is just to make the last days with max value
+            data_ds.cumsum()
+            # groupby methods drop grouped coordinates so need to reassign T
+            .assign_coords(T=labelled_season_data["T"])
+            # grouping again
+            .groupby(labelled_season_data["seasons_starts"])
+            # idxmax returns label of coord where max is, so last frosting day
+            .apply(lambda x : x.idxmax(dim="T"))
+            # idxmax returns the label of first coord when multiple maxima so more 
+            # straightforward
+            - data_ds.apply(lambda x : x.idxmax(dim="T"))
             + np.timedelta64(1, "D")
         )
+        for var in data_ds.data_vars :
+            # we want to plot timedelta as days
+            data_ds[var] = data_ds[var].dt.days
         wwc_units = "days"
     if variable == "wet_days":
         data_ds = (
@@ -233,7 +247,11 @@ def seasonal_wwc(
     # Can revisit later if this code has a future
     data_ds = data_ds.rename({"seasons_starts" : "T"})
     seasons_ends = labelled_season_data["seasons_ends"].rename({"group": "T"})
-    data_ds = data_ds.drop_vars(["seasons_ends", "group"]).assign_coords({"seasons_ends" : seasons_ends})
+    data_ds = (
+        data_ds
+        .drop_vars(["seasons_ends", "group"], errors="ignore")
+        .assign_coords({"seasons_ends" : seasons_ends})
+    )
     #
     for var in data_ds.data_vars:
         data_ds[var].attrs["units"] = wwc_units
diff --git a/pepsico/proj_wwc/layout.py b/pepsico/proj_wwc/layout.py
index a3b88ecc..96ace797 100644
--- a/pepsico/proj_wwc/layout.py
+++ b/pepsico/proj_wwc/layout.py
@@ -43,8 +43,8 @@ def app_layout():
                         "mean_Tmin",
                         "dry_days",
                         "wet_days",
-                        # "heatwave_duration",
-                        #"frost_season_length",
+                        # "heatwave_duration", # needs a daily clim
+                        "frost_season_length",
                         "frost_days",
                         "Tmax_90",
                         "Tmin_10",
@@ -63,7 +63,7 @@ def app_layout():
                         "Count of Dry Days",
                         "Count of Wet Days",
                         # "Mean Heatwaves Duration",
-                        #"Frost Season Length",
+                        "Frost Season Length",
                         "Count of Frost Days",
                         "Max Temperature 90th %-ile",
                         "Min Temperature 10th %-ile",
@@ -268,6 +268,12 @@ def app_layout():
                     lesser or equal than a user-defined thredhold.
                 """
             ]),
+            html.P([
+                html.B("Frost Season Length (frost_season_length):"),"""
+                    Number of days between the first and last days where minimm 
+                    temperature is lesser or equal than a user-defined threshold.
+                """
+            ]),
         ),
         
         lou.map(GLOBAL_CONFIG["zoom"]),
diff --git a/pepsico/proj_wwc/maproom.py b/pepsico/proj_wwc/maproom.py
index b7caaa07..9bdab37c 100644
--- a/pepsico/proj_wwc/maproom.py
+++ b/pepsico/proj_wwc/maproom.py
@@ -524,6 +524,7 @@ def map_attributes(data, variable):
         if variable in [
             "frost_days", "dry_days", "longest_dry_spell", "dry_day_persistence",
             "dry_spells_mean_length", "dry_spells_median_length",
+            "frost_season_length",
         ]:
             colorscale = colorscale.reversed()
         map_amp = np.max(np.abs(data)).values