From 760640d75b0b3c3b077f7124e6426acbdff19e6a Mon Sep 17 00:00:00 2001 From: remic Date: Thu, 6 Nov 2025 11:54:15 -0500 Subject: [PATCH 01/12] app layout --- pepsico/proj_wwc/layout.py | 178 +++++++++++++++++++++++++++++++++++++ 1 file changed, 178 insertions(+) create mode 100644 pepsico/proj_wwc/layout.py diff --git a/pepsico/proj_wwc/layout.py b/pepsico/proj_wwc/layout.py new file mode 100644 index 00000000..d22424ad --- /dev/null +++ b/pepsico/proj_wwc/layout.py @@ -0,0 +1,178 @@ +from dash import html +from fieldsets import Block, Select, PickPoint, Month, Number +import layout_utilities as lou + + +from globals_ import GLOBAL_CONFIG + +IRI_BLUE = "rgb(25,57,138)" +IRI_GRAY = "rgb(113,112,116)" +LIGHT_GRAY = "#eeeeee" + + +def app_layout(): + + return lou.app_1( + + lou.navbar("CCA - WWC", + Block("Region", + Select( + id="region", + options=["Thailand", "US-CA"], + labels=["Thailand", "United States and Canada"], + init=3, + ), + ), + PickPoint(width="8em"), + Block("Submit", + Block("Scenario", Select( + id="scenario", + options=["picontrol", "ssp126", "ssp370", "ssp585"], + init=1, + )), + Block("Model", Select(id="model", options=[ + "Multi-Model-Average", "GFDL-ESM4", "IPSL-CM6A-LR", + "MPI-ESM1-2-HR", "MRI-ESM2-0", "UKESM1-0-LL", + ])), + Block("WWC statistics", Select( + id="variable", + options=[ + "killing_frost", + # "warm_nights", + # "rain_events", + "dry_days", + "mean_Tmax", + "mean_Tmin", + "Tmax_90th-ile", + "Tmin_10th-ile", + # "heatwave_duration", + "frost_season_length", + "frost_days", + "wet_days", + # "longest_dry_spell", + # "longest_wet_spell", + "wet_day_persistence", + "dry_day_persistence", + "dry_spells_mean_length", + "dry_spells_median_length", + ], + labels=[ + "Killing Frost", + # "Warm Nights", + # "Count of Rain Events", + "Count of Dry Days", + "Mean Max Temperature", + "Mean Min Temperature", + "Max Temperature 90th %-ile", + "Min Temperature 10th %-ile", + # "Mean Heatwaves Duration", + "Frost Season Length", + "Count of Frost Days", + "Count of Wet Days", + # "Longest Dry Spell", + # "Longest Wet Spell", + "Mean Wet Day Persistence", + "Mean Dry Day Persistence", + "Mean Dry Spells Length", + "Median Dry Spells Length", + ], + init=2, + )), + Block("Season", + Number( + id="start_day", + default=1, + min=1, + max=31, + width="5em", + debounce=False, + ), + Month(id="start_month", default="Jan"), + "-", + Number( + id="end_day", + default=31, + min=1, + max=31, + width="5em", + debounce=False, + ), + Month(id="end_month", default="Mar"), + ), + Block("Projected Years", + Number( + id="start_year", + default=2015, + min=2015, + max=2099, + width="5em", + debounce=False, + ), + "-", + Number( + id="end_year", + default=2019, + min=2015, + max=2099, + width="5em", + debounce=False, + ), + ), + Block("Reference Years", + Number( + id="start_year_ref", + default=1981, + min=1951, + max=2014, + width="5em", + debounce=False, + ), + "-", + Number( + id="end_year_ref", + default=2010, + min=1951, + max=2014, + width="5em", + debounce=False, + ), + ), + button_id="submit_controls", + ), + ), + + lou.description( + "Weather-within-Climate Change Analysis", + """ + This Maproom displays seasonal projected change of key + weather-within-climate agronomic and climatic variables with respect to + historical records. + """ + html.P( + """ + Use the controls in the top banner to choose other variables, models, + scenarios, seasons, projected years and reference to compare with. + """ + ), + html.P( + """ + Click the map (or enter coordinates) to show historical seasonal time + series for this variable of this model, followed by a plume of + possible projected scenarios. + """ + ), + html.P( + """ + Change is expressed as the difference between average over projected + years and average over reference historical years (in the variables + units). + """ + ), + ), + + lou.map(GLOBAL_CONFIG["zoom"]), + + lou.local_single_tabbed( + "Local History and Projections", download_button=True + ), + ) From b934dea518c8e4bc17f52fe0fcff418aeaf51acc Mon Sep 17 00:00:00 2001 From: remic Date: Fri, 7 Nov 2025 09:35:29 -0500 Subject: [PATCH 02/12] setting up maproom --- pepsico/app_calc.py | 191 ++++++++++- pepsico/proj_wwc/layout.py | 41 ++- pepsico/proj_wwc/maproom.py | 617 ++++++++++++++++++++++++++++++++++++ 3 files changed, 837 insertions(+), 12 deletions(-) create mode 100644 pepsico/proj_wwc/maproom.py diff --git a/pepsico/app_calc.py b/pepsico/app_calc.py index ff8d8323..774a289f 100644 --- a/pepsico/app_calc.py +++ b/pepsico/app_calc.py @@ -1,4 +1,5 @@ import xarray as xr +import numpy as np #This is what we should need for the app @@ -19,7 +20,9 @@ #-- then maybe some other things to beautify map tbd -def read_data(scenario, model, variable, region, unit_convert=False): +def read_data( + scenario, model, variable, region, time_res="monthly", unit_convert=False +): if region == "US-CA": xslice = slice(-154, -45) yslice = slice(60, 15) @@ -34,7 +37,7 @@ def read_data(scenario, model, variable, region, unit_convert=False): yslice = slice(28, 2) data = xr.open_zarr( f'/Data/data24/ISIMIP3b/InputData/climate/atmosphere/bias-adjusted/global' - f'/monthly/{scenario}/{model}/zarr/{variable}' + f'/{time_res}/{scenario}/{model}/zarr/{variable}' )[variable].sel(X=xslice, Y=yslice) if unit_convert : data = unit_conversion(data) @@ -92,6 +95,45 @@ def seasonal_data(monthly_data, start_month, end_month, start_year=None, end_yea ) +def seasonal_wwc( + labelled_season_data, data_var, variable, frost_threshold, wet_threshold +): + if variable == "frost_days": + data_ds = (labelled_season_data[data_var] <= frost_threshold).groupby( + labelled_season_data["seasons_starts"] + ).sum() + if variable == "dry_days": + data_ds = (labelled_season_data[data_var] <= wet_threshold).groupby( + labelled_season_data["seasons_starts"] + ).sum() + if variable in ["mean_Tmax", "mean_Tmin"]: + data_ds = labelled_season_data[data_var].groupby( + labelled_season_data["seasons_starts"] + ).mean() + if variable == "Tmax_90th-ile": + data_ds = labelled_season_data[data_var].groupby( + labelled_season_data["seasons_starts"] + ).quantile(0.9) + if variable == "Tmin_10th-ile": + data_ds = labelled_season_data[data_var].groupby( + labelled_season_data["seasons_starts"] + ).quantile(0.1) + if variable == "frost_season_length": + data_ds = (labelled_season_data[data_var] <= frost_threshold).groupby( + labelled_season_data["seasons_starts"] + ) + data_ds = ( + data_ds[-1:0].idxmax() + - data_ds.idxmax() + + np.timedelta64(1, "D") + ) + if variable == "wet_days": + data_ds = (labelled_season_data[data_var] > wet_threshold).groupby( + labelled_season_data["seasons_starts"] + ).sum() + return data_ds + + def unit_conversion(variable): #if precipitation variable, change from kg to mm per day if variable.name == 'pr': @@ -109,6 +151,151 @@ def unit_conversion(variable): # This is in enacts calc.py + + +def daily_tobegroupedby_season( + daily_data, start_day, start_month, end_day, end_month, time_dim="T" +): + """Group daily data by season. + + Returns dataset ready to be grouped by with the daily data where all days not in season of interest are dropped. + + season_starts: + An array where the non-dropped days are indexed by the first day of their season. + -- to use to groupby + seasons_ends: + An array with the dates of the end of the seasons. + Can then apply groupby on daily_data against seasons_starts, and preserving seasons_ends for the record. + + If starting day-month is 29-Feb, uses 1-Mar. + + If ending day-month is 29-Feb, uses 1-Mar and uses < rather than <= + That means that the last day included in the season will be 29-Feb in leap years and 28-Feb otherwise. + + Parameters + ----------- + daily_data : DataArray + Daily data to be grouped. + start_day : int + Day of the start date of the season. + start_month : int + Month of the start date of the season. + end_day : int + Day of the end date of the season. + end_month : int + Day of the end date of the season. + time_dim : str, optional + Time coordinate in `daily_data` (default `time_dim`="T"). + Returns + ------- + daily_tobegroupedby_season : Dataset + Daily data grouped by season using season start date. Dataset includes grouped data + and `season_starts`, `season_ends` as output variables. + See Also + -------- + Notes + ----- + Examples + -------- + """ + # Deal with leap year cases + if start_day == 29 and start_month == 2: + start_day = 1 + start_month = 3 + # Find seasons edges + start_edges = sel_day_and_month(daily_data[time_dim], start_day, start_month) + if end_day == 29 and end_month == 2: + end_edges = sel_day_and_month(daily_data[time_dim], 1 , 3, offset=-1) + else: + end_edges = sel_day_and_month(daily_data[time_dim], end_day, end_month) + # Drop dates outside very first and very last edges + # -- this ensures we get complete seasons with regards to edges, later on + daily_data = daily_data.sel(**{time_dim: slice(start_edges[0], end_edges[-1])}) + start_edges = start_edges.sel(**{time_dim: slice(start_edges[0], end_edges[-1])}) + end_edges = end_edges.sel( + **{time_dim: slice(start_edges[0], end_edges[-1])} + ).assign_coords(**{time_dim: start_edges[time_dim]}) + # Drops daily data not in seasons of interest + days_in_season = ( + daily_data[time_dim] >= start_edges.rename({time_dim: "group"}) + ) & (daily_data[time_dim] <= end_edges.rename({time_dim: "group"})) + days_in_season = days_in_season.sum(dim="group") + daily_data = daily_data.where(days_in_season == 1, drop=True) + # Creates seasons_starts that will be used for grouping + # and seasons_ends that is one of the outputs + seasons_groups = (daily_data[time_dim].dt.day == start_day) & ( + daily_data[time_dim].dt.month == start_month + ) + seasons_groups = seasons_groups.cumsum() - 1 + seasons_starts = ( + start_edges.rename({time_dim: "toto"})[seasons_groups] + .drop_vars("toto") + .rename("seasons_starts") + ) + seasons_ends = end_edges.rename({time_dim: "group"}).rename("seasons_ends") + # Dataset output + daily_tobegroupedby_season = xr.merge([daily_data, seasons_starts, seasons_ends]) + return daily_tobegroupedby_season + + +def sel_day_and_month(daily_dim, day, month, offset=0): + """Return a subset of `daily_dim` daily time dimension of corresponding + `day`/`month` + `offset` day(s) for all years. + + The returned time dimension can then be used to select daily DataArrays. + Offset is convenient to get days prior to a 1st of March, + that are not identifiable by a common `day` (28, 29). + + Parameters + ---------- + daily_dim : DataArray[datetime64[ns]] + A daily time dimension. + day : int + day of the `month`. + month : int + month of the year. + offset : int, optional + number of days to add to `day`/`month` to offset the selection + (the default is 0, which implies no offset). + + Returns + ------- + DataArray[datetime64[ns]] + a subset of `daily_dim` with all and only `day`-`month` points, offset by `offset` days. + + See Also + -------- + + Examples + -------- + >>> t = pd.date_range(start="2000-05-01", end="20002-04-30", freq="1D") + >>> values = numpy.arrange(t.size) + >>> toto = xarray.DataArray(numpy.arrange(61), dims=["T"], coords={"T": t}) + >>> sel_day_and_month(toto["T"], 6, 5) + + array(['2000-05-06T00:00:00.000000000', '2001-05-06T00:00:00.000000000',] + dtype='datetime64[ns]') + Coordinates: + * T (T) datetime64[ns] 2000-05-06 2001-05-06 + + With an offset of -1 day + + >>> t = pd.date_range(start="2000-01-01", end="20002-01-30", freq="1D") + >>> toto = xarray.DataArray(numpy.arrange(t.size), dims=["T"], coords={"T": t}) + >>> sel_day_and_month(toto["T"], 1, 3, -1) + + array(['2000-02-29T00:00:00.000000000', '2001-02-28T00:00:00.000000000',] + dtype='datetime64[ns]') + Coordinates: + * T (T) datetime64[ns] 2000-02-29 2001-02-28 + """ + return daily_dim.where( + lambda x: ((x - np.timedelta64(offset, "D")).dt.day == day) + & ((x - np.timedelta64(offset, "D")).dt.month == month), + drop=True + ) + + def strftimeb2int(strftimeb): """Convert month values to integers (1-12) from strings. diff --git a/pepsico/proj_wwc/layout.py b/pepsico/proj_wwc/layout.py index d22424ad..247c1f1a 100644 --- a/pepsico/proj_wwc/layout.py +++ b/pepsico/proj_wwc/layout.py @@ -37,7 +37,6 @@ def app_layout(): Block("WWC statistics", Select( id="variable", options=[ - "killing_frost", # "warm_nights", # "rain_events", "dry_days", @@ -51,13 +50,12 @@ def app_layout(): "wet_days", # "longest_dry_spell", # "longest_wet_spell", - "wet_day_persistence", - "dry_day_persistence", - "dry_spells_mean_length", - "dry_spells_median_length", + # "wet_day_persistence", + # "dry_day_persistence", + # "dry_spells_mean_length", + # "dry_spells_median_length", ], labels=[ - "Killing Frost", # "Warm Nights", # "Count of Rain Events", "Count of Dry Days", @@ -71,13 +69,36 @@ def app_layout(): "Count of Wet Days", # "Longest Dry Spell", # "Longest Wet Spell", - "Mean Wet Day Persistence", - "Mean Dry Day Persistence", - "Mean Dry Spells Length", - "Median Dry Spells Length", + # "Mean Wet Day Persistence", + # "Mean Dry Day Persistence", + # "Mean Dry Spells Length", + # "Median Dry Spells Length", ], init=2, )), + Block("Definitions", + "Frost <=", + Number( + id="frost", + default=-2, + min=-99, + max=0, + width="5em", + debounce=False, + ), + "˚C", + "Dry/Wet day <=/>", + Number( + id="wet", + default=0, + min=0, + max=999, + width="5em", + debounce=False, + ), + "mm", + + ) Block("Season", Number( id="start_day", diff --git a/pepsico/proj_wwc/maproom.py b/pepsico/proj_wwc/maproom.py new file mode 100644 index 00000000..637aaa8f --- /dev/null +++ b/pepsico/proj_wwc/maproom.py @@ -0,0 +1,617 @@ +import dash +import dash_bootstrap_components as dbc +from dash.dependencies import Output, Input, State +import pingrid +from pingrid import CMAPS +from . import layout +import plotly.graph_objects as pgo +import xarray as xr +import pandas as pd +from dateutil.relativedelta import * +from globals_ import FLASK, GLOBAL_CONFIG +import app_calc as ac +import maproom_utilities as mru +import numpy as np + + +def register(FLASK, config): + PFX = f"{GLOBAL_CONFIG['url_path_prefix']}/{config['core_path']}" + TILE_PFX = f"{PFX}/tile" + + # App + + APP = dash.Dash( + __name__, + server=FLASK, + external_stylesheets=[ + dbc.themes.BOOTSTRAP, + ], + url_base_pathname=f"{PFX}/", + meta_tags=[ + {"name": "description", "content": "WWC-CCA"}, + {"name": "viewport", "content": "width=device-width, initial-scale=1.0"}, + ], + ) + APP.title = "WWC-CCA" + + APP.layout = layout.app_layout() + + + @APP.callback( + Output("lat_input", "min"), + Output("lat_input", "max"), + Output("lat_input_tooltip", "children"), + Output("lng_input", "min"), + Output("lng_input", "max"), + Output("lng_input_tooltip", "children"), + Output("map", "center"), + Output("map", "zoom"), + Input("region", "value"), + Input("location", "pathname"), + ) + def initialize(region, path): + scenario = "ssp126" + model = "GFDL-ESM4" + variable = "pr" + data = ac.read_data(scenario, model, variable, region, time_res="daily") + zoom = {"US-CA": 4, "Thailand": 5} + return mru.initialize_map(data) + (zoom[region],) + + + @APP.callback( + Output("loc_marker", "position"), + Output("lat_input", "value"), + Output("lng_input", "value"), + Input("submit_lat_lng","n_clicks"), + Input("map", "click_lat_lng"), + Input("region", "value"), + State("lat_input", "value"), + State("lng_input", "value"), + ) + def pick_location(n_clicks, click_lat_lng, region, latitude, longitude): + # Reading + scenario = "ssp126" + model = "GFDL-ESM4" + variable = "pr" + data = ac.read_data(scenario, model, variable, region, time_res="daily") + initialization_cases = ["region"] + return mru.picked_location( + data, initialization_cases, click_lat_lng, latitude, longitude + ) + + + def select_var(variable): + if variable in [ + "killing_frost", + "mean_Tmin", + "Tmin_10th-ile", + "frost_season_length", + "frost_days", + ]: + data_var = "tasmin" + if variable in [ + # "warm_nights", + "mean_Tmax", + "Tmax_90th-ile", + # "heatwave_duration", + ]: + data_var = "tasmx" + if variable in [ + # "rain_events", + "dry_days", + "wet_days", + # "longest_dry_spell", + # "longest_wet_spell", + "wet_day_persistence", + "dry_day_persistence", + "dry_spells_mean_length", + "dry_spells_median_length", + ]: + data_var = "pr" + return data_var + + + def local_data( + lat, lng, region, + model, variable, + start_day, start_month, end_day, end_month, + frost_threshold, wet_threshold, + ): + model = [model] if model != "Multi-Model-Average" else [ + "GFDL-ESM4", "IPSL-CM6A-LR", "MPI-ESM1-2-HR","MRI-ESM2-0", "UKESM1-0-LL" + ] + data_var = select_var(variable) + data_ds = xr.concat([xr.Dataset({ + "histo" : ac.read_data( + "historical", m, data_var, region, + time_res="daily", unit_convert=True, + ), + "picontrol" : ac.read_data( + "picontrol", m, data_var, region, + time_res="daily", unit_convert=True, + ), + "ssp126" : ac.read_data( + "ssp126", m, data_var, region, + time_res="daily", unit_convert=True, + ), + "ssp370" : ac.read_data( + "ssp370", m, data_var, region, + time_res="daily", unit_convert=True, + ), + "ssp585" : ac.read_data( + "ssp585", m, data_var, region, + time_res="daily", unit_convert=True, + ), + }) for m in model], "M").assign_coords({"M": [m for m in model]}) + error_msg = None + missing_ds = xr.Dataset() + if any([var is None for var in data_ds.data_vars.values()]): + #This is not supposed to happen: + #would mean something happened to that data + data_ds = missing_ds + error_msg="Data missing for this model or variable" + try: + data_ds = pingrid.sel_snap(data_ds, lat, lng) + except KeyError: + data_ds = missing_ds + error_msg="Grid box out of data domain" + if error_msg == None : + data_ds = ac.seasonal_wwc( + ac.daily_tobegroupedby_season( + data_ds, start_day, start_month, end_day, end_month, + ), + data_var, variable, frost_threshold, wet_threshold, + ) + return data_ds, error_msg + + + def plot_ts(ts, name, color, start_format, units): + return pgo.Scatter( + x=ts["T"].dt.strftime(mru.STD_TIME_FORMAT), + y=ts.values, + customdata=ts["seasons_ends"].dt.strftime("%B %Y"), + hovertemplate=("%{x|"+start_format+"}%{customdata}: %{y:.2f} " + units), + name=name, + line=pgo.scatter.Line(color=color), + connectgaps=False, + ) + + + def add_period_shape( + graph, data, start_year, end_year, fill_color, line_color, annotation + ): + return graph.add_vrect( + x0=data["seasons_starts"].where( + lambda x : (x.dt.year == int(start_year)), drop=True + ).dt.strftime(mru.STD_TIME_FORMAT).values[0], + #it's hard to believe this is how it is done + x1=( + pd.to_datetime(data["seasons_ends"].where( + lambda x : (x.dt.year == int(end_year)), drop=True + ).dt.strftime(mru.STD_TIME_FORMAT).values[0] + ) + relativedelta(months=+1)).strftime(mru.STD_TIME_FORMAT), + fillcolor=fill_color, opacity=0.2, + line_color=line_color, line_width=3, + layer="below", + annotation_text=annotation, annotation_position="top left", + #editable=True, #a reminder it might be the way to interact + ) + + + @APP.callback( + Output("btn_csv", "disabled"), + Input("lat_input", "value"), + Input("lng_input", "value"), + Input("lat_input", "min"), + Input("lng_input", "min"), + Input("lat_input", "max"), + Input("lng_input", "max"), + ) + def invalid_button(lat, lng, lat_min, lng_min, lat_max, lng_max): + return ( + lat < float(lat_min) or lat > float(lat_max) + or lng < float(lng_min) or lng > float(lng_max) + ) + + + @APP.callback( + Output("download-dataframe-csv", "data"), + Input("btn_csv", "n_clicks"), + State("loc_marker", "position"), + State("region", "value"), + State("variable", "value"), + State("start_day", "value"), + State("start_month", "value"), + State("end_day", "value"), + State("end_month", "value"), + State("frost", "value"), + State("wet", "value"), + prevent_initial_call=True, + ) + def send_data_as_csv( + n_clicks, marker_pos, region, variable, + start_day, start_month, end_day, end_month, + frost_threshold, wet_threshold, + ): + lat = marker_pos[0] + lng = marker_pos[1] + start_month = ac.strftimeb2int(start_month) + end_month = ac.strftimeb2int(end_month) + model = "Multi-Model-Average" + data_ds, error_msg = local_data( + lat, lng, region, model, variable, + start_day, start_month, end_day, end_month, + frost_threshold, wet_threshold, + ) + if error_msg == None : + lng_units = "E" if (lng >= 0) else "W" + lat_units = "N" if (lat >= 0) else "S" + file_name = ( + f'{data_ds["histo"]["T"].dt.strftime("%m%d")[0].values}-' + f'{data_ds["histo"]["seasons_ends"].dt.strftime("%m%d")[0].values}' + f'_{variable}_{abs(lat)}{lat_units}_{abs(lng)}{lng_units}' + f'.csv' + ) + df = data_ds.to_dataframe() + return dash.dcc.send_data_frame(df.to_csv, file_name) + else : + return None + + + @APP.callback( + Output("local_graph", "figure"), + Input("loc_marker", "position"), + Input("region", "value"), + Input("submit_controls","n_clicks"), + State("model", "value"), + State("variable", "value"), + State("start_day", "value"), + State("end_day", "value"), + State("start_month", "value"), + State("end_month", "value"), + State("start_year", "value"), + State("end_year", "value"), + State("start_year_ref", "value"), + State("end_year_ref", "value"), + State("frost", "value"), + State("wet", "value"), + ) + def local_plots( + marker_pos, region, n_clicks, model, variable, + start_day, end_day, start_month, end_month, + start_year, end_year, start_year_ref, end_year_ref, + frost_threshold, wet_threshold, + ): + lat = marker_pos[0] + lng = marker_pos[1] + start_month = ac.strftimeb2int(start_month) + end_month = ac.strftimeb2int(end_month) + data_ds, error_msg = local_data( + lat, lng, region, model, variable, + start_day, start_month, end_day, end_month, + frost_threshold, wet_threshold, + ) + if error_msg != None : + local_graph = pingrid.error_fig(error_msg) + else : + if (end_month < start_month) : + start_format = "%d %b %Y - " + else: + start_format = "%d %b-" + local_graph = pgo.Figure() + data_color = { + "histo": "blue", "picontrol": "green", + "ssp126": "yellow", "ssp370": "orange", "ssp585": "red", + } + lng_units = "˚E" if (lng >= 0) else "˚W" + lat_units = "˚N" if (lat >= 0) else "˚S" + for var in data_ds.data_vars: + local_graph.add_trace(plot_ts( + data_ds[var].mean("M", keep_attrs=True), var, data_color[var], + start_format, data_ds[var].attrs["units"] + )) + add_period_shape( + local_graph, + data_ds, + start_year_ref, + end_year_ref, + "blue", + "RoyalBlue", + "reference period", + ) + add_period_shape( + local_graph, + data_ds, + start_year, + end_year, + "LightPink", + "Crimson", + "projected period", + ) + local_graph.update_layout( + xaxis_title="Time", + yaxis_title=( + f'{data_ds["histo"].attrs["long_name"]} ' + f'({data_ds["histo"].attrs["units"]})' + ), + title={ + "text": ( + f'{data_ds["histo"]["T"].dt.strftime("%d %b")[0].values}-' + f'{data_ds["histo"]["seasons_ends"].dt.strftime("%d %b")[0].values}' + f' {variable} seasonal average from model {model} ' + f'at ({abs(lat)}{lat_units}, {abs(lng)}{lng_units})' + ), + "font": dict(size=14), + }, + margin=dict(l=30, r=30, t=30, b=30), + ) + return local_graph + + + @APP.callback( + Output("map_description", "children"), + Input("submit_controls", "n_clicks"), + State("scenario", "value"), + State("model", "value"), + State("variable", "value"), + State("start_day", "value"), + State("end_day", "value"), + State("start_month", "value"), + State("end_month", "value"), + State("start_year", "value"), + State("end_year", "value"), + State("start_year_ref", "value"), + State("end_year_ref", "value"), + ) + def write_map_description( + n_clicks, + scenario, + model, + variable, + start_day, + end_day, + start_month, + end_month, + start_year, + end_year, + start_year_ref, + end_year_ref, + ): + return ( + f'The Map displays the change in ' + f'({start_day} {start_month} - {end_day} {end_month}) season of ' + f'{variable} from {model} model under {scenario} scenario projected for ' + f'{start_year}-{end_year} with respect to historical {start_year_ref}-' + f'{end_year_ref}' + ) + + + @APP.callback( + Output("map_title", "children"), + Input("submit_controls","n_clicks"), + State("scenario", "value"), + State("model", "value"), + State("variable", "value"), + State("start_day", "value"), + State("end_day", "value"), + State("start_month", "value"), + State("end_month", "value"), + State("start_year", "value"), + State("end_year", "value"), + State("start_year_ref", "value"), + State("end_year_ref", "value"), + ) + def write_map_title( + n_clicks, + scenario, + model, + variable, + start_day, + end_day, + start_month, + end_month, + start_year, + end_year, + start_year_ref, + end_year_ref, + ): + return ( + f'({start_day} {start_month} - {end_day} {end_month}) ' + f'{start_year}-{end_year} ' + f'{scenario} {model} {variable} change from ' + f'{start_year_ref}-{end_year_ref}' + ) + + + APP.clientside_callback( + """function(start_year, end_year, start_year_ref, end_year_ref) { + if (start_year && end_year) { + invalid_start_year = (start_year > end_year) + invalid_end_year = invalid_start_year + } else { + invalid_start_year = !start_year + invalid_end_year = !end_year + } + if (start_year_ref && end_year_ref) { + invalid_start_year_ref = (start_year_ref > end_year_ref) + invalid_end_year_ref = invalid_start_year_ref + } else { + invalid_start_year_ref = !start_year_ref + invalid_end_year_ref = !end_year_ref + } + return [ + invalid_start_year, invalid_end_year, + invalid_start_year_ref, invalid_end_year_ref, + ( + invalid_start_year || invalid_end_year + || invalid_start_year_ref || invalid_end_year_ref + ), + ] + } + """, + Output("start_year", "invalid"), + Output("end_year", "invalid"), + Output("start_year_ref", "invalid"), + Output("end_year_ref", "invalid"), + Output("submit_controls", "disabled"), + Input("start_year", "value"), + Input("end_year", "value"), + Input("start_year_ref", "value"), + Input("end_year_ref", "value"), + ) + + + def seasonal_change( + scenario, model, variable, region, + start_day, end_day, start_month, end_month, + start_year, end_year, start_year_ref, end_year_ref, + frost_threshold, wet_threshold, + ): + model = [model] if model != "Multi-Model-Average" else [ + "GFDL-ESM4", "IPSL-CM6A-LR", "MPI-ESM1-2-HR","MRI-ESM2-0", "UKESM1-0-LL" + ] + data_var = select_var(variable) + ref = xr.concat([ + ac.seasonal_wwc( + ac.daily_tobegrouped_by( + ac.read_data( + "historical", m, variable, region, + time_res="daily", unit_convert=True + ).sel(T=slice(start_year_ref, end_year_ref)), + start_day, start_month, end_day, end_month, + ), + data_var, variable, frost_threshold, wet_threshold, + ).mean(dim="T", keep_attrs=True) for m in model + ], "M").mean("M", keep_attrs=True) + data = xr.concat([ + ac.seasonal_wwc( + ac.daily_tobegrouped_by( + ac.read_data( + scenario, m, variable, region, + time_res="daily", unit_convert=True + ).sel(T=slice(start_year, end_year)), + start_day, start_month, end_day, end_month, + ), + data_var, variable, frost_threshold, wet_threshold, + ).mean(dim="T", keep_attrs=True) for m in model + ], "M").mean("M", keep_attrs=True) + #Tedious way to make a subtraction only to keep attributes + return xr.apply_ufunc( + np.subtract, data, ref, dask="allowed", keep_attrs="drop_conflicts", + ).rename({"X": "lon", "Y": "lat"}) + + + def map_attributes(data): + map_amp = np.max(np.abs(data)).values + colorscale = CMAPS["correlation"] + colorscale = colorscale.rescaled(-1*map_amp, map_amp) + return colorscale, colorscale.scale[0], colorscale.scale[-1] + + + @APP.callback( + Output("colorbar", "colorscale"), + Output("colorbar", "min"), + Output("colorbar", "max"), + Output("colorbar", "unit"), + Input("region", "value"), + Input("submit_controls","n_clicks"), + State("scenario", "value"), + State("model", "value"), + State("variable", "value"), + State("start_day", "value"), + State("end_day", "value"), + State("start_month", "value"), + State("end_month", "value"), + State("start_year", "value"), + State("end_year", "value"), + State("start_year_ref", "value"), + State("end_year_ref", "value"), + State("frost", "value"), + State("wet", "value"), + ) + def draw_colorbar( + region, n_clicks, scenario, model, variable, + start_day, end_day, start_month, end_month, + start_year, end_year, start_year_ref, end_year_ref, + frost_threshold, wet_threshold, + ): + data = seasonal_change( + scenario, model, variable, region, + start_day, end_day, start_month, end_month, + start_year, end_year, start_year_ref, end_year_ref, + frost_threshold, wet_threshold, + ) + colorbar, min, max = map_attributes(data) + return colorbar.to_dash_leaflet(), min, max, data.attrs["units"] + + + @APP.callback( + Output("layers_control", "children"), + Output("map_warning", "is_open"), + Input("region", "value"), + Input("submit_controls", "n_clicks"), + State("scenario", "value"), + State("model", "value"), + State("variable", "value"), + State("start_day", "value"), + State("end_day", "value"), + State("start_month", "value"), + State("end_month", "value"), + State("start_year", "value"), + State("end_year", "value"), + State("start_year_ref", "value"), + State("end_year_ref", "value"), + State("frost", "value"), + State("wet", "value"), + ) + def make_map( + region, n_clicks, scenario, model, variable, + start_day, end_day, start_month, end_month, + start_year, end_year, start_year_ref, end_year_ref, + frost_threshold, wet_threshold, + ): + try: + send_alarm = False + url_str = ( + f"{TILE_PFX}/{{z}}/{{x}}/{{y}}/{region}/{scenario}/{model}/{variable}/" + f"{start_day}/{end_day}/{start_month}/{end_month}/" + f"{start_year}/{end_year}/{start_year_ref}/{end_year_ref}/" + f"{frost_threshold}/{wet_threshold}" + ) + except: + url_str= "" + send_alarm = True + return mru.layers_controls( + url_str, f"wwc_change_{region}", "WWC Change", + GLOBAL_CONFIG["datasets"][f"shapes_adm_{region}"], GLOBAL_CONFIG, + adm_id_suffix=region, + ), send_alarm + + + @FLASK.route( + ( + f"{TILE_PFX}////////" + f"/////" + f"" + ), + endpoint=f"{config['core_path']}" + ) + def fcst_tiles(tz, tx, ty, + region, scenario, model, variable, + start_day, end_day, start_month, end_month, + start_year, end_year, start_year_ref, end_year_ref, + frost_threshold, wet_threshold, + ): + data = seasonal_change( + scenario, model, variable, region, + start_day, end_day, start_month, end_month, + start_year, end_year, start_year_ref, end_year_ref, + frost_threshold, wet_threshold, + ) + ( + data.attrs["colormap"], + data.attrs["scale_min"], + data.attrs["scale_max"], + ) = map_attributes(data) + resp = pingrid.tile(data, tx, ty, tz) + return resp From 3591e1c84ca4de324594411ab3b54e9c9452e2cf Mon Sep 17 00:00:00 2001 From: remic Date: Fri, 7 Nov 2025 09:38:59 -0500 Subject: [PATCH 03/12] setting up config --- pepsico/config-dev-sample.yaml | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/pepsico/config-dev-sample.yaml b/pepsico/config-dev-sample.yaml index f71905db..c361213a 100644 --- a/pepsico/config-dev-sample.yaml +++ b/pepsico/config-dev-sample.yaml @@ -15,4 +15,8 @@ maprooms: # App set up - title: Projections core_path: projections + + proj_wwc: + - title: WWC Projections + core_path: proj_wwc \ No newline at end of file From a2550cfa2c9a3425250135ad46f9907dfcf0ffcd Mon Sep 17 00:00:00 2001 From: remic Date: Fri, 7 Nov 2025 09:41:54 -0500 Subject: [PATCH 04/12] I might have too many repeats here... --- pepsico/config-defaults.yaml | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/pepsico/config-defaults.yaml b/pepsico/config-defaults.yaml index 4df5fed6..6baa20da 100644 --- a/pepsico/config-defaults.yaml +++ b/pepsico/config-defaults.yaml @@ -92,3 +92,7 @@ maprooms: # App set up - title: Projections core_path: projections + + proj_wwc: + - title: WWC Projections + core_path: proj_wwc From 20c96d11ef7d1e12e95114f5e53d030c574ec256 Mon Sep 17 00:00:00 2001 From: remic Date: Fri, 7 Nov 2025 09:55:22 -0500 Subject: [PATCH 05/12] need this to build --- pepsico/proj_wwc/__init__.py | 1 + 1 file changed, 1 insertion(+) create mode 100644 pepsico/proj_wwc/__init__.py diff --git a/pepsico/proj_wwc/__init__.py b/pepsico/proj_wwc/__init__.py new file mode 100644 index 00000000..0379666a --- /dev/null +++ b/pepsico/proj_wwc/__init__.py @@ -0,0 +1 @@ +from .maproom import register From 2bd72541c2440a56114c850b43deb86e20bbefba Mon Sep 17 00:00:00 2001 From: remic Date: Mon, 10 Nov 2025 16:57:06 -0500 Subject: [PATCH 06/12] first instance of this new app --- pepsico/app_calc.py | 56 ++++++++++++++---- pepsico/proj_wwc/layout.py | 25 ++++---- pepsico/proj_wwc/maproom.py | 115 ++++++++++++++++++++---------------- 3 files changed, 120 insertions(+), 76 deletions(-) diff --git a/pepsico/app_calc.py b/pepsico/app_calc.py index 774a289f..abba7205 100644 --- a/pepsico/app_calc.py +++ b/pepsico/app_calc.py @@ -41,6 +41,10 @@ def read_data( )[variable].sel(X=xslice, Y=yslice) if unit_convert : data = unit_conversion(data) + #Turns out that some models are centered on noon + if time_res == "daily" : + if data["T"][0].dt.hour == 12 : + data = data.assign_coords({"T" : data["T"] - np.timedelta64(12, "h")}) return data @@ -96,30 +100,42 @@ def seasonal_data(monthly_data, start_month, end_month, start_year=None, end_yea def seasonal_wwc( - labelled_season_data, data_var, variable, frost_threshold, wet_threshold + labelled_season_data, variable, frost_threshold, wet_threshold ): + # For some rason the summed groups return 0s where all NaNs (begining/end of + # projections/histo). A where ~np.isnan didn't help, neither a skipNa True. To be + # explored later if variable == "frost_days": - data_ds = (labelled_season_data[data_var] <= frost_threshold).groupby( + data_ds = (labelled_season_data <= frost_threshold).groupby( labelled_season_data["seasons_starts"] ).sum() + wwc_units = "days" if variable == "dry_days": - data_ds = (labelled_season_data[data_var] <= wet_threshold).groupby( + data_ds = (labelled_season_data <= wet_threshold).groupby( labelled_season_data["seasons_starts"] ).sum() + wwc_units = "days" if variable in ["mean_Tmax", "mean_Tmin"]: - data_ds = labelled_season_data[data_var].groupby( + data_ds = labelled_season_data.groupby( labelled_season_data["seasons_starts"] ).mean() - if variable == "Tmax_90th-ile": - data_ds = labelled_season_data[data_var].groupby( + wwc_units = "˚C" + # Percentiles options have been commented out in the layout at the calculations + # seem to have thrown off the app (got a time series once but after that... just + # spinning around) + if variable == "Tmax_90": + data_ds = labelled_season_data.groupby( labelled_season_data["seasons_starts"] ).quantile(0.9) - if variable == "Tmin_10th-ile": - data_ds = labelled_season_data[data_var].groupby( + wwc_units = "˚C" + if variable == "Tmin_10": + data_ds = labelled_season_data.groupby( labelled_season_data["seasons_starts"] ).quantile(0.1) + # 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[data_var] <= frost_threshold).groupby( + data_ds = (labelled_season_data <= frost_threshold).groupby( labelled_season_data["seasons_starts"] ) data_ds = ( @@ -127,10 +143,21 @@ def seasonal_wwc( - data_ds.idxmax() + np.timedelta64(1, "D") ) + wwc_units = "days" if variable == "wet_days": - data_ds = (labelled_season_data[data_var] > wet_threshold).groupby( + data_ds = ((labelled_season_data > wet_threshold)+0).groupby( labelled_season_data["seasons_starts"] ).sum() + wwc_units = "days" + # This is all a bit tedious but I didn't figure out another way to keep + # seasons_ends and renaming time dim T + # 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}) + # + for var in data_ds.data_vars: + data_ds[var].attrs["units"] = wwc_units return data_ds @@ -233,8 +260,13 @@ def daily_tobegroupedby_season( .rename("seasons_starts") ) seasons_ends = end_edges.rename({time_dim: "group"}).rename("seasons_ends") - # Dataset output - daily_tobegroupedby_season = xr.merge([daily_data, seasons_starts, seasons_ends]) + # I actually changed this from enacts by making both seasons_starts/ends cords + # rather than vars. This seems better following a comment Aaron made on xcdat and + # because we local data case provides an xr.Dataset of vars to groupby and we + # don't want seasons_starts/ends to be groupedby + daily_tobegroupedby_season = daily_data.assign_coords( + {"seasons_starts" : seasons_starts, "seasons_ends" : seasons_ends} + ) return daily_tobegroupedby_season diff --git a/pepsico/proj_wwc/layout.py b/pepsico/proj_wwc/layout.py index 247c1f1a..2be12442 100644 --- a/pepsico/proj_wwc/layout.py +++ b/pepsico/proj_wwc/layout.py @@ -20,7 +20,7 @@ def app_layout(): id="region", options=["Thailand", "US-CA"], labels=["Thailand", "United States and Canada"], - init=3, + init=1, ), ), PickPoint(width="8em"), @@ -42,10 +42,10 @@ def app_layout(): "dry_days", "mean_Tmax", "mean_Tmin", - "Tmax_90th-ile", - "Tmin_10th-ile", + #"Tmax_90", + #"Tmin_10", # "heatwave_duration", - "frost_season_length", + #"frost_season_length", "frost_days", "wet_days", # "longest_dry_spell", @@ -61,10 +61,10 @@ def app_layout(): "Count of Dry Days", "Mean Max Temperature", "Mean Min Temperature", - "Max Temperature 90th %-ile", - "Min Temperature 10th %-ile", + #"Max Temperature 90th %-ile", + #"Min Temperature 10th %-ile", # "Mean Heatwaves Duration", - "Frost Season Length", + #"Frost Season Length", "Count of Frost Days", "Count of Wet Days", # "Longest Dry Spell", @@ -74,7 +74,7 @@ def app_layout(): # "Mean Dry Spells Length", # "Median Dry Spells Length", ], - init=2, + init=1, )), Block("Definitions", "Frost <=", @@ -86,8 +86,8 @@ def app_layout(): width="5em", debounce=False, ), - "˚C", - "Dry/Wet day <=/>", + "˚C; ", + "Dry/Wet day <= / >", Number( id="wet", default=0, @@ -97,8 +97,7 @@ def app_layout(): debounce=False, ), "mm", - - ) + ), Block("Season", Number( id="start_day", @@ -168,7 +167,7 @@ def app_layout(): This Maproom displays seasonal projected change of key weather-within-climate agronomic and climatic variables with respect to historical records. - """ + """, html.P( """ Use the controls in the top banner to choose other variables, models, diff --git a/pepsico/proj_wwc/maproom.py b/pepsico/proj_wwc/maproom.py index 637aaa8f..17890b22 100644 --- a/pepsico/proj_wwc/maproom.py +++ b/pepsico/proj_wwc/maproom.py @@ -84,7 +84,7 @@ def select_var(variable): if variable in [ "killing_frost", "mean_Tmin", - "Tmin_10th-ile", + "Tmin_10", "frost_season_length", "frost_days", ]: @@ -92,10 +92,10 @@ def select_var(variable): if variable in [ # "warm_nights", "mean_Tmax", - "Tmax_90th-ile", + "Tmax_90", # "heatwave_duration", ]: - data_var = "tasmx" + data_var = "tasmax" if variable in [ # "rain_events", "dry_days", @@ -160,7 +160,7 @@ def local_data( ac.daily_tobegroupedby_season( data_ds, start_day, start_month, end_day, end_month, ), - data_var, variable, frost_threshold, wet_threshold, + variable, frost_threshold, wet_threshold, ) return data_ds, error_msg @@ -181,7 +181,7 @@ def add_period_shape( graph, data, start_year, end_year, fill_color, line_color, annotation ): return graph.add_vrect( - x0=data["seasons_starts"].where( + x0=data["T"].where( lambda x : (x.dt.year == int(start_year)), drop=True ).dt.strftime(mru.STD_TIME_FORMAT).values[0], #it's hard to believe this is how it is done @@ -235,8 +235,12 @@ def send_data_as_csv( ): lat = marker_pos[0] lng = marker_pos[1] + start_day = int(start_day) + end_day = int(end_day) start_month = ac.strftimeb2int(start_month) end_month = ac.strftimeb2int(end_month) + frost_threshold = float(frost_threshold) + wet_threshold = float(wet_threshold) model = "Multi-Model-Average" data_ds, error_msg = local_data( lat, lng, region, model, variable, @@ -284,8 +288,12 @@ def local_plots( ): lat = marker_pos[0] lng = marker_pos[1] + start_day = int(start_day) + end_day = int(end_day) start_month = ac.strftimeb2int(start_month) end_month = ac.strftimeb2int(end_month) + frost_threshold = float(frost_threshold) + wet_threshold = float(wet_threshold) data_ds, error_msg = local_data( lat, lng, region, model, variable, start_day, start_month, end_day, end_month, @@ -331,7 +339,10 @@ def local_plots( local_graph.update_layout( xaxis_title="Time", yaxis_title=( - f'{data_ds["histo"].attrs["long_name"]} ' + # Compared to the monthly case, I presume that groupby dropped + # the long_name, which is fine since it's indeed transformed. Can + # consider resetting it at some point in the calculation + #f'{data_ds["histo"].attrs["long_name"]} ' f'({data_ds["histo"].attrs["units"]})' ), title={ @@ -364,18 +375,9 @@ def local_plots( State("end_year_ref", "value"), ) def write_map_description( - n_clicks, - scenario, - model, - variable, - start_day, - end_day, - start_month, - end_month, - start_year, - end_year, - start_year_ref, - end_year_ref, + n_clicks, scenario, model, variable, + start_day, end_day, start_month, end_month, + start_year, end_year, start_year_ref, end_year_ref, ): return ( f'The Map displays the change in ' @@ -402,18 +404,9 @@ def write_map_description( State("end_year_ref", "value"), ) def write_map_title( - n_clicks, - scenario, - model, - variable, - start_day, - end_day, - start_month, - end_month, - start_year, - end_year, - start_year_ref, - end_year_ref, + n_clicks, scenario, model, variable, + start_day, end_day, start_month, end_month, + start_year, end_year, start_year_ref, end_year_ref, ): return ( f'({start_day} {start_month} - {end_day} {end_month}) ' @@ -473,29 +466,33 @@ def seasonal_change( data_var = select_var(variable) ref = xr.concat([ ac.seasonal_wwc( - ac.daily_tobegrouped_by( + ac.daily_tobegroupedby_season( ac.read_data( - "historical", m, variable, region, + "historical", m, data_var, region, time_res="daily", unit_convert=True - ).sel(T=slice(start_year_ref, end_year_ref)), + ).sel( + T=slice(str(start_year_ref), str(end_year_ref)) + ).to_dataset(), start_day, start_month, end_day, end_month, ), - data_var, variable, frost_threshold, wet_threshold, + variable, frost_threshold, wet_threshold, ).mean(dim="T", keep_attrs=True) for m in model ], "M").mean("M", keep_attrs=True) data = xr.concat([ ac.seasonal_wwc( - ac.daily_tobegrouped_by( + ac.daily_tobegroupedby_season( ac.read_data( - scenario, m, variable, region, + scenario, m, data_var, region, time_res="daily", unit_convert=True - ).sel(T=slice(start_year, end_year)), + ).sel( + T=slice(str(start_year), str(end_year)) + ).to_dataset(), start_day, start_month, end_day, end_month, ), - data_var, variable, frost_threshold, wet_threshold, + variable, frost_threshold, wet_threshold, ).mean(dim="T", keep_attrs=True) for m in model ], "M").mean("M", keep_attrs=True) - #Tedious way to make a subtraction only to keep attributes + #Tedious way to make a subtraction only to keep attributes (units) return xr.apply_ufunc( np.subtract, data, ref, dask="allowed", keep_attrs="drop_conflicts", ).rename({"X": "lon", "Y": "lat"}) @@ -535,14 +532,23 @@ def draw_colorbar( start_year, end_year, start_year_ref, end_year_ref, frost_threshold, wet_threshold, ): + start_day = int(start_day) + end_day = int(end_day) + start_month = ac.strftimeb2int(start_month) + end_month = ac.strftimeb2int(end_month) + frost_threshold = float(frost_threshold) + wet_threshold = float(wet_threshold) data = seasonal_change( scenario, model, variable, region, start_day, end_day, start_month, end_month, start_year, end_year, start_year_ref, end_year_ref, frost_threshold, wet_threshold, ) - colorbar, min, max = map_attributes(data) - return colorbar.to_dash_leaflet(), min, max, data.attrs["units"] + colorbar, min, max = map_attributes(data[select_var(variable)]) + return ( + colorbar.to_dash_leaflet(), min, max, + data[select_var(variable)].attrs["units"] + ) @APP.callback( @@ -573,8 +579,8 @@ def make_map( try: send_alarm = False url_str = ( - f"{TILE_PFX}/{{z}}/{{x}}/{{y}}/{region}/{scenario}/{model}/{variable}/" - f"{start_day}/{end_day}/{start_month}/{end_month}/" + f"{TILE_PFX}/{{z}}/{{x}}/{{y}}/{region}/{scenario}/{model}/" + f"{variable}/{start_day}/{end_day}/{start_month}/{end_month}/" f"{start_year}/{end_year}/{start_year_ref}/{end_year_ref}/" f"{frost_threshold}/{wet_threshold}" ) @@ -590,9 +596,10 @@ def make_map( @FLASK.route( ( - f"{TILE_PFX}////////" - f"/////" - f"" + f"{TILE_PFX}///////" + f"/////" + f"////" + f"//" ), endpoint=f"{config['core_path']}" ) @@ -602,6 +609,12 @@ def fcst_tiles(tz, tx, ty, start_year, end_year, start_year_ref, end_year_ref, frost_threshold, wet_threshold, ): + start_day = int(start_day) + end_day = int(end_day) + start_month = ac.strftimeb2int(start_month) + end_month = ac.strftimeb2int(end_month) + frost_threshold = float(frost_threshold) + wet_threshold = float(wet_threshold) data = seasonal_change( scenario, model, variable, region, start_day, end_day, start_month, end_month, @@ -609,9 +622,9 @@ def fcst_tiles(tz, tx, ty, frost_threshold, wet_threshold, ) ( - data.attrs["colormap"], - data.attrs["scale_min"], - data.attrs["scale_max"], - ) = map_attributes(data) - resp = pingrid.tile(data, tx, ty, tz) + data[select_var(variable)].attrs["colormap"], + data[select_var(variable)].attrs["scale_min"], + data[select_var(variable)].attrs["scale_max"], + ) = map_attributes(data[select_var(variable)]) + resp = pingrid.tile(data[select_var(variable)], tx, ty, tz) return resp From 99b3e5281c6cd2972f059fd3bcf286b4588f8992 Mon Sep 17 00:00:00 2001 From: remic Date: Fri, 14 Nov 2025 08:55:36 -0500 Subject: [PATCH 07/12] adding day to season end graph popup --- pepsico/proj_wwc/maproom.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pepsico/proj_wwc/maproom.py b/pepsico/proj_wwc/maproom.py index 17890b22..9ace6d9e 100644 --- a/pepsico/proj_wwc/maproom.py +++ b/pepsico/proj_wwc/maproom.py @@ -169,7 +169,7 @@ def plot_ts(ts, name, color, start_format, units): return pgo.Scatter( x=ts["T"].dt.strftime(mru.STD_TIME_FORMAT), y=ts.values, - customdata=ts["seasons_ends"].dt.strftime("%B %Y"), + customdata=ts["seasons_ends"].dt.strftime("%d %B %Y"), hovertemplate=("%{x|"+start_format+"}%{customdata}: %{y:.2f} " + units), name=name, line=pgo.scatter.Line(color=color), From 6aff108e532df5abdf701124d7c31bd49c6d77e2 Mon Sep 17 00:00:00 2001 From: remic Date: Fri, 14 Nov 2025 09:20:19 -0500 Subject: [PATCH 08/12] additional manips so that all NaNs seasons don't return 0s --- pepsico/app_calc.py | 35 +++++++++++++++++++++++------------ pepsico/proj_wwc/layout.py | 2 +- 2 files changed, 24 insertions(+), 13 deletions(-) diff --git a/pepsico/app_calc.py b/pepsico/app_calc.py index abba7205..4fc5b0c3 100644 --- a/pepsico/app_calc.py +++ b/pepsico/app_calc.py @@ -102,18 +102,26 @@ def seasonal_data(monthly_data, start_month, end_month, start_year=None, end_yea def seasonal_wwc( labelled_season_data, variable, frost_threshold, wet_threshold ): - # For some rason the summed groups return 0s where all NaNs (begining/end of - # projections/histo). A where ~np.isnan didn't help, neither a skipNa True. To be - # explored later + # Boolean variables need the additional where to return NaNs from False/0 to Nans + # and the sum parameters for entirely NaNs seasons to remain NaNs and not turn to + # 0. This is necessary because histo and projections scenarios don't cover the + # full time period and it's not creating fake results as there are no missing + # data if variable == "frost_days": - data_ds = (labelled_season_data <= frost_threshold).groupby( - labelled_season_data["seasons_starts"] - ).sum() + data_ds = ( + (labelled_season_data <= frost_threshold) + .where(~np.isnan(labelled_season_data)) + .groupby(labelled_season_data["seasons_starts"]) + .sum(skipna=True, min_count=1) + ) wwc_units = "days" if variable == "dry_days": - data_ds = (labelled_season_data <= wet_threshold).groupby( - labelled_season_data["seasons_starts"] - ).sum() + data_ds = ( + (labelled_season_data <= wet_threshold) + .where(~np.isnan(labelled_season_data)) + .groupby(labelled_season_data["seasons_starts"]) + .sum(skipna=True, min_count=1) + ) wwc_units = "days" if variable in ["mean_Tmax", "mean_Tmin"]: data_ds = labelled_season_data.groupby( @@ -145,9 +153,12 @@ def seasonal_wwc( ) wwc_units = "days" if variable == "wet_days": - data_ds = ((labelled_season_data > wet_threshold)+0).groupby( - labelled_season_data["seasons_starts"] - ).sum() + data_ds = ( + (labelled_season_data > wet_threshold) + .where(~np.isnan(labelled_season_data)) + .groupby(labelled_season_data["seasons_starts"]) + .sum(skipna=True, min_count=1) + ) wwc_units = "days" # This is all a bit tedious but I didn't figure out another way to keep # seasons_ends and renaming time dim T diff --git a/pepsico/proj_wwc/layout.py b/pepsico/proj_wwc/layout.py index 2be12442..04b8f55b 100644 --- a/pepsico/proj_wwc/layout.py +++ b/pepsico/proj_wwc/layout.py @@ -74,7 +74,7 @@ def app_layout(): # "Mean Dry Spells Length", # "Median Dry Spells Length", ], - init=1, + init=0, )), Block("Definitions", "Frost <=", From 95f6c06614b2e8669db932b8b4371eff48b268c5 Mon Sep 17 00:00:00 2001 From: remic Date: Fri, 14 Nov 2025 10:57:36 -0500 Subject: [PATCH 09/12] update on quantiles problematic comment --- pepsico/app_calc.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/pepsico/app_calc.py b/pepsico/app_calc.py index 4fc5b0c3..d01969e4 100644 --- a/pepsico/app_calc.py +++ b/pepsico/app_calc.py @@ -128,18 +128,18 @@ def seasonal_wwc( labelled_season_data["seasons_starts"] ).mean() wwc_units = "˚C" - # Percentiles options have been commented out in the layout at the calculations - # seem to have thrown off the app (got a time series once but after that... just - # spinning around) + # It takes several 10s of minutes to get the quantiles maps so the options are + # commented out in the list of options until a solution is found if variable == "Tmax_90": data_ds = labelled_season_data.groupby( labelled_season_data["seasons_starts"] - ).quantile(0.9) + ).quantile(0.9, method="closest_observation") wwc_units = "˚C" if variable == "Tmin_10": data_ds = labelled_season_data.groupby( labelled_season_data["seasons_starts"] ).quantile(0.1) + 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": From d96a7195e05e8857130626ec70c5212588fdcc60 Mon Sep 17 00:00:00 2001 From: remic Date: Fri, 14 Nov 2025 12:12:58 -0500 Subject: [PATCH 10/12] more appropriate color scales --- pepsico/proj_wwc/layout.py | 2 +- pepsico/proj_wwc/maproom.py | 14 ++++++++++---- pingrid/impl.py | 6 ++++-- 3 files changed, 15 insertions(+), 7 deletions(-) diff --git a/pepsico/proj_wwc/layout.py b/pepsico/proj_wwc/layout.py index 04b8f55b..0670f642 100644 --- a/pepsico/proj_wwc/layout.py +++ b/pepsico/proj_wwc/layout.py @@ -74,7 +74,7 @@ def app_layout(): # "Mean Dry Spells Length", # "Median Dry Spells Length", ], - init=0, + init=3, )), Block("Definitions", "Frost <=", diff --git a/pepsico/proj_wwc/maproom.py b/pepsico/proj_wwc/maproom.py index 9ace6d9e..2597851b 100644 --- a/pepsico/proj_wwc/maproom.py +++ b/pepsico/proj_wwc/maproom.py @@ -498,9 +498,15 @@ def seasonal_change( ).rename({"X": "lon", "Y": "lat"}) - def map_attributes(data): + def map_attributes(data, variable): + data = data[select_var(variable)] + if data.name == "pr": + colorscale = CMAPS["prcp_anomaly"] + if data.name in ["tasmin", "tasmax"]: + colorscale = CMAPS["temp_anomaly"] + if variable in ["frost_days", "dry_days"]: + colorscale = colorscale.reversed() map_amp = np.max(np.abs(data)).values - colorscale = CMAPS["correlation"] colorscale = colorscale.rescaled(-1*map_amp, map_amp) return colorscale, colorscale.scale[0], colorscale.scale[-1] @@ -544,7 +550,7 @@ def draw_colorbar( start_year, end_year, start_year_ref, end_year_ref, frost_threshold, wet_threshold, ) - colorbar, min, max = map_attributes(data[select_var(variable)]) + colorbar, min, max = map_attributes(data, variable) return ( colorbar.to_dash_leaflet(), min, max, data[select_var(variable)].attrs["units"] @@ -625,6 +631,6 @@ def fcst_tiles(tz, tx, ty, data[select_var(variable)].attrs["colormap"], data[select_var(variable)].attrs["scale_min"], data[select_var(variable)].attrs["scale_max"], - ) = map_attributes(data[select_var(variable)]) + ) = map_attributes(data, variable) resp = pingrid.tile(data[select_var(variable)], tx, ty, tz) return resp diff --git a/pingrid/impl.py b/pingrid/impl.py index 20af3351..5fbe5dbb 100644 --- a/pingrid/impl.py +++ b/pingrid/impl.py @@ -883,8 +883,10 @@ def produce_shape_tile( _TEMP_ANOMALY_CS = ColorScale( "temp_anomaly", - [PURPLE, CYAN, WHITE, WHITE, YELLOW, RED, FIREBRICK], - [-10, -1, -1, 1, 1, 10, 10], + # The duplicate first anchor makes the scale symmetric and thus allows to use the + # reversed method + [PURPLE, PURPLE, CYAN, WHITE, WHITE, YELLOW, RED, FIREBRICK], + [-10, -10, -1, -1, 1, 1, 10, 10], ) _RAIN_PNE_CS = _RAIN_POE_CS.reversed(name="rain_pne") From 3345b51ccf9f043f654144a14ad784ff82ab400e Mon Sep 17 00:00:00 2001 From: remic Date: Fri, 14 Nov 2025 12:32:01 -0500 Subject: [PATCH 11/12] Description of statistics definitions --- pepsico/proj_wwc/layout.py | 22 +++++++++++++++++----- 1 file changed, 17 insertions(+), 5 deletions(-) diff --git a/pepsico/proj_wwc/layout.py b/pepsico/proj_wwc/layout.py index 0670f642..a531acaf 100644 --- a/pepsico/proj_wwc/layout.py +++ b/pepsico/proj_wwc/layout.py @@ -39,15 +39,15 @@ def app_layout(): options=[ # "warm_nights", # "rain_events", - "dry_days", "mean_Tmax", "mean_Tmin", + "dry_days", + "wet_days", #"Tmax_90", #"Tmin_10", # "heatwave_duration", #"frost_season_length", "frost_days", - "wet_days", # "longest_dry_spell", # "longest_wet_spell", # "wet_day_persistence", @@ -58,15 +58,15 @@ def app_layout(): labels=[ # "Warm Nights", # "Count of Rain Events", - "Count of Dry Days", "Mean Max Temperature", "Mean Min Temperature", + "Count of Dry Days", + "Count of Wet Days", #"Max Temperature 90th %-ile", #"Min Temperature 10th %-ile", # "Mean Heatwaves Duration", #"Frost Season Length", "Count of Frost Days", - "Count of Wet Days", # "Longest Dry Spell", # "Longest Wet Spell", # "Mean Wet Day Persistence", @@ -74,7 +74,7 @@ def app_layout(): # "Mean Dry Spells Length", # "Median Dry Spells Length", ], - init=3, + init=0, )), Block("Definitions", "Frost <=", @@ -188,6 +188,18 @@ def app_layout(): units). """ ), + html.H4(["Definitions"]), + html.P([html.B("Mean Max/Min Temperature (mean_Tmax/min):"), """ + Average maximum/minimum temperature in season + """]), + html.P([html.B("Count of Dry/Wet Days (dry/wet_days):"), """ + Number of dry/wet days in the season. A dry/wet day is defined as + lesser or equal / greather than a user-defined threshold. + """]), + html.P([html.B("Count of Frost Days (frost_days):"), """ + Number of frost days in the season. A frost day is defined as + lesser or equal than a user-defined threshold. + """]), ), lou.map(GLOBAL_CONFIG["zoom"]), From 4c9e048cd905343dda77efdd4eaf16902cc88c77 Mon Sep 17 00:00:00 2001 From: remic Date: Tue, 18 Nov 2025 15:38:42 -0500 Subject: [PATCH 12/12] parametric quantiles much less costly than empirical ones --- pepsico/app_calc.py | 20 +++++++++++--------- pepsico/proj_wwc/layout.py | 14 ++++++++++---- 2 files changed, 21 insertions(+), 13 deletions(-) diff --git a/pepsico/app_calc.py b/pepsico/app_calc.py index d01969e4..a52b860f 100644 --- a/pepsico/app_calc.py +++ b/pepsico/app_calc.py @@ -1,5 +1,6 @@ import xarray as xr import numpy as np +from scipy.stats import norm #This is what we should need for the app @@ -128,17 +129,18 @@ def seasonal_wwc( labelled_season_data["seasons_starts"] ).mean() wwc_units = "˚C" - # It takes several 10s of minutes to get the quantiles maps so the options are - # commented out in the list of options until a solution is found - if variable == "Tmax_90": + # It takes several 10s of minutes to get empirical quantiles maps + if variable in ["Tmax_90", "Tmin_10"]: + quantile = 0.1 if variable == "Tmin_10" else 0.9 data_ds = labelled_season_data.groupby( labelled_season_data["seasons_starts"] - ).quantile(0.9, method="closest_observation") - wwc_units = "˚C" - if variable == "Tmin_10": - data_ds = labelled_season_data.groupby( - labelled_season_data["seasons_starts"] - ).quantile(0.1) + ) + data_std = data_ds.std() + data_ds = data_ds.mean() + for var in data_ds.data_vars: + data_ds[var].data = xr.apply_ufunc( + 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. diff --git a/pepsico/proj_wwc/layout.py b/pepsico/proj_wwc/layout.py index a531acaf..c069bcaf 100644 --- a/pepsico/proj_wwc/layout.py +++ b/pepsico/proj_wwc/layout.py @@ -43,11 +43,11 @@ def app_layout(): "mean_Tmin", "dry_days", "wet_days", - #"Tmax_90", - #"Tmin_10", # "heatwave_duration", #"frost_season_length", "frost_days", + "Tmax_90", + "Tmin_10", # "longest_dry_spell", # "longest_wet_spell", # "wet_day_persistence", @@ -62,11 +62,11 @@ def app_layout(): "Mean Min Temperature", "Count of Dry Days", "Count of Wet Days", - #"Max Temperature 90th %-ile", - #"Min Temperature 10th %-ile", # "Mean Heatwaves Duration", #"Frost Season Length", "Count of Frost Days", + "Max Temperature 90th %-ile", + "Min Temperature 10th %-ile", # "Longest Dry Spell", # "Longest Wet Spell", # "Mean Wet Day Persistence", @@ -200,6 +200,12 @@ def app_layout(): Number of frost days in the season. A frost day is defined as lesser or equal than a user-defined threshold. """]), + html.P([ + html.B("Max/Min Temperature 90/10th %-ile (Tmax/min_90/10):"),""" + Maximum/Minimum temperature 90/10th percentile in the season. + Obtained through parametric Normal distributions. + """ + ]), ), lou.map(GLOBAL_CONFIG["zoom"]),