diff --git a/extreme_data/meteo_france_data/scm_models_data/safran/safran_variable.py b/extreme_data/meteo_france_data/scm_models_data/safran/safran_variable.py
index 55d8d7b5982ac3f06b608cbe12285f017928b6ca..b1966f5d8503e2884d2a6c3ef731cb5889a8fcc5 100644
--- a/extreme_data/meteo_france_data/scm_models_data/safran/safran_variable.py
+++ b/extreme_data/meteo_france_data/scm_models_data/safran/safran_variable.py
@@ -126,7 +126,7 @@ class SafranNormalizedPrecipitationRateOnWetDaysVariable(SafranNormalizedPrecipi
class SafranTemperatureVariable(AbstractVariable):
NAME = 'Temperature'
- UNIT = 'Celsius Degrees'
+ UNIT = '$^oC$'
@classmethod
def keyword(cls):
diff --git a/extreme_data/meteo_france_data/scm_models_data/utils.py b/extreme_data/meteo_france_data/scm_models_data/utils.py
index 3f1c57e9b35f0fc9177eabfd0eb56471dd72ba53..97c8c2ea54e8831438e18385d65875e2ec108762 100644
--- a/extreme_data/meteo_france_data/scm_models_data/utils.py
+++ b/extreme_data/meteo_france_data/scm_models_data/utils.py
@@ -11,6 +11,7 @@ WP_PATTERN_MAX_YEAR = 2008
alps_massif_order = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 30]
pyrenees_massif_order = [64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91]
+
def date_to_str(date: datetime) -> str:
return str(date).split()[0]
@@ -21,6 +22,12 @@ class SeasonForTheMaxima(Enum):
# i could add the classical seasons if needed
+def season_to_str(season):
+ season_name = str(season).split('.')[-1]
+ season_name = ' '.join(season_name.split('_'))
+ return season_name
+
+
class FrenchRegion(Enum):
alps = 1
pyrenees = 2
diff --git a/projects/contrasting_trends_in_snow_loads/gorman_figures/figure 3/figure3_daily snowfall fraction.py b/projects/contrasting_trends_in_snow_loads/gorman_figures/daily_snowfall_fraction.py
similarity index 80%
rename from projects/contrasting_trends_in_snow_loads/gorman_figures/figure 3/figure3_daily snowfall fraction.py
rename to projects/contrasting_trends_in_snow_loads/gorman_figures/daily_snowfall_fraction.py
index 1243fc450ada9b103f94d9858ce7325066ac72ca..f1e40e47dcaaeff3cd8d1915f863df5f10a89562 100644
--- a/projects/contrasting_trends_in_snow_loads/gorman_figures/figure 3/figure3_daily snowfall fraction.py
+++ b/projects/contrasting_trends_in_snow_loads/gorman_figures/daily_snowfall_fraction.py
@@ -8,9 +8,10 @@ from extreme_data.meteo_france_data.scm_models_data.safran.safran import SafranR
Altitudes_groups = [[300, 600, 900], [1200, 1500, 1800], [2100, 2400, 2700], [3000, 3300, 3600]]
-def compute_smoother_snowfall_fraction(altitudes, year_min, year_max, lim, massif_names=None):
- all_time_series = np.concatenate(
- [get_time_series(altitude, year_min, year_max, massif_names) for altitude in altitudes], axis=1)
+def compute_smoother_snowfall_fraction(altitudes, year_min, year_max, lim, massif_names=None, all_time_series=None):
+ if all_time_series is None:
+ all_time_series = np.concatenate(
+ [get_time_series(altitude, year_min, year_max, massif_names) for altitude in altitudes], axis=1)
temperature_array, snowfall_array, rainfall_array = all_time_series
# nb_bins_for_one_celsius = 4
nb_bins_for_one_celsius = 4
@@ -25,10 +26,14 @@ def compute_smoother_snowfall_fraction(altitudes, year_min, year_max, lim, massi
return np.array(snowfall_fractions), center_bins
-def daily_snowfall_fraction(ax, altitudes, year_min=1959, year_max=2019, massif_names=None, sigma=1.0):
+def compute_daily_snowfall_fraction(ax=None, altitudes=None, year_min=1959, year_max=2019, massif_names=None, sigma=1.0, plot=True, all_time_series=None):
+ if plot:
+ assert ax is not None
+ assert altitudes is not None
+
lim = 6
snowfall_fractions, center_bins = compute_smoother_snowfall_fraction(altitudes, year_min, year_max, lim=lim,
- massif_names=massif_names)
+ massif_names=massif_names, all_time_series=all_time_series)
# Smooth the data afterwards
snowfall_fractions = 100 * gaussian_filter1d(snowfall_fractions, sigma=sigma)
@@ -39,15 +44,17 @@ def daily_snowfall_fraction(ax, altitudes, year_min=1959, year_max=2019, massif_
h = np.exp(0.06 * center_bins_interpolation) * snowfall_fractions_interpolated
temperature_optimal = center_bins_interpolation[np.argmax(h)]
+ if plot:
+ # Plot results
+ label = '{}-{} at {}'.format(year_min, year_max, ' & '.join(['{} m'.format(a) for a in altitudes]))
+ ax.set_xlim([-lim, lim])
+ p = ax.plot(center_bins, snowfall_fractions, label=label, linewidth=5)
+ color = p[0].get_color()
+ ax.plot([temperature_optimal], [0], color=color, marker='x', markersize=10)
+ ax.set_ylabel('Snowfall fraction (%)')
+ end_plot(ax, sigma)
- # Plot results
- label = '{}-{} at {}'.format(year_min, year_max, ' & '.join(['{} m'.format(a) for a in altitudes]))
- ax.set_xlim([-lim, lim])
- p = ax.plot(center_bins, snowfall_fractions, label=label, linewidth=5)
- color = p[0].get_color()
- ax.plot([temperature_optimal], [0], color=color, marker='x', markersize=10)
- ax.set_ylabel('Snowfall fraction (%)')
- end_plot(ax, sigma)
+ return snowfall_fractions_interpolated, center_bins_interpolation, temperature_optimal
def snowfall_fraction_difference_between_periods(ax, altitudes, year_min=1959, year_middle=1989, year_max=2019):
@@ -101,7 +108,7 @@ def daily_snowfall_fraction_wrt_altitude(fast=True):
massif_names = None
# massif_names = ['Vercors']
for altitudes in groups:
- daily_snowfall_fraction(ax, altitudes, massif_names=massif_names)
+ compute_daily_snowfall_fraction(ax, altitudes, massif_names=massif_names)
plt.show()
@@ -116,7 +123,7 @@ def daily_snowfall_fraction_wrt_time():
ax = plt.gca()
for altitudes in Altitudes_groups:
for year_min, year_max in [(1959, 1989), (1990, 2019)]:
- daily_snowfall_fraction(ax, altitudes, year_min, year_max)
+ compute_daily_snowfall_fraction(ax, altitudes, year_min, year_max)
plt.show()
diff --git a/projects/contrasting_trends_in_snow_loads/gorman_figures/temperature_of_snowfall_maxima.py b/projects/contrasting_trends_in_snow_loads/gorman_figures/temperature_of_snowfall_maxima.py
new file mode 100644
index 0000000000000000000000000000000000000000..a47ec2ee6e5dd40dff7e8dd3cb44ef5480034c65
--- /dev/null
+++ b/projects/contrasting_trends_in_snow_loads/gorman_figures/temperature_of_snowfall_maxima.py
@@ -0,0 +1,69 @@
+import matplotlib as mpl
+mpl.rcParams['text.usetex'] = True
+mpl.rcParams['text.latex.preamble'] = [r'\usepackage{amsmath}']
+
+from extreme_data.meteo_france_data.scm_models_data.utils import SeasonForTheMaxima, season_to_str
+
+from collections import OrderedDict
+
+import matplotlib.pyplot as plt
+import numpy as np
+
+from extreme_data.meteo_france_data.scm_models_data.safran.safran import SafranSnowfall1Day, SafranTemperature, \
+ SafranPrecipitation1Day
+from projects.contrasting_trends_in_snow_loads.gorman_figures.daily_snowfall_fraction import \
+ compute_daily_snowfall_fraction
+
+
+def distribution_temperature_of_maxima_of_snowfall(altitudes, temperature_at_maxima=True):
+ season = SeasonForTheMaxima.annual if temperature_at_maxima else SeasonForTheMaxima.winter_extended
+ # Load the temperature corresponding to the maxima of snowfall
+ altitude_to_maxima_temperature = OrderedDict()
+ altitude_to_optimal_temperature = OrderedDict()
+ for altitude in altitudes:
+ print(altitude)
+ snowfall_study = SafranSnowfall1Day(altitude=altitude, season=season)
+ temperature_study = SafranTemperature(altitude=altitude, season=season)
+ precipitation_study = SafranPrecipitation1Day(altitude=altitude, season=season)
+ # Compute optimal temperature
+ all_time_series = [temperature_study.all_daily_series, snowfall_study.all_daily_series,
+ precipitation_study.all_daily_series]
+ *_, optimal_temperature = compute_daily_snowfall_fraction(plot=False, all_time_series=all_time_series)
+ altitude_to_optimal_temperature[altitude] = optimal_temperature
+ # Compute temperature corresponding to maxima
+ maxima_temperatures = []
+ for year in snowfall_study.ordered_years[:2]:
+ a = temperature_study.year_to_daily_time_serie_array[year]
+ if temperature_at_maxima:
+ annual_maxima_index = snowfall_study.year_to_annual_maxima_tuple_indices_for_daily_time_series[year]
+ temp = [a[tuple_idx] for tuple_idx in annual_maxima_index]
+ else:
+ temp = a.flatten()
+ maxima_temperatures.append(temp)
+ altitude_to_maxima_temperature[altitude] = np.concatenate(maxima_temperatures)
+
+ ax = plt.gca()
+ # Plot the optimal temperature
+ ax.plot(altitudes, [altitude_to_optimal_temperature[a] for a in altitudes], marker='o', linestyle='--', label='Optimal temperature for snowfall')
+
+ # Plot the box plot
+ width = 150
+ ax.boxplot([altitude_to_maxima_temperature[a] for a in altitudes], positions=altitudes, widths=width)
+ ax.set_xlim([min(altitudes) - width, max(altitudes) + width])
+ suffix = 'at maxima of snowfall' if temperature_at_maxima else ''
+ ylabel = 'Daily {} temperature {} ({})'.format(season_to_str(season), suffix, temperature_study.variable_class.UNIT)
+ print(ylabel)
+ ax.set_ylabel(ylabel)
+ ax.set_xlabel('Altitude (m)')
+ ax.legend()
+ ax.grid()
+ ax.set_ylim([-8, 5])
+
+ plt.show()
+
+
+if __name__ == '__main__':
+ # altitudes = [900, 1200]
+ altitudes = [900, 1200, 1500, 1800, 2100, 2400, 2700, 3000]
+ distribution_temperature_of_maxima_of_snowfall(altitudes=altitudes,
+ temperature_at_maxima=False)
diff --git a/projects/contrasting_trends_in_snow_loads/temperature_of_snowfall_maxima.py b/projects/contrasting_trends_in_snow_loads/temperature_of_snowfall_maxima.py
deleted file mode 100644
index 6de908ab65930a8ed64fb0fa927d12ab5a061926..0000000000000000000000000000000000000000
--- a/projects/contrasting_trends_in_snow_loads/temperature_of_snowfall_maxima.py
+++ /dev/null
@@ -1,35 +0,0 @@
-from extreme_data.meteo_france_data.scm_models_data.safran.safran import SafranSnowfall1Day, SafranTemperature
-import numpy as np
-import matplotlib.pyplot as plt
-
-
-def distribution_temperature_of_maxima_of_snowfall(altitudes):
- # Load the temperature corresponding to the maxima of snowfall
- altitude_to_maxima_temperature = {}
- for altitude in altitudes:
- print(altitude)
- snowfall_study = SafranSnowfall1Day(altitude=altitude)
- temperature_study = SafranTemperature(altitude=altitude)
- # temperature_study = SafranTemperature(altitude=altitude)
- maxima_temperatures = []
- for year in snowfall_study.ordered_years[:2]:
- annual_maxima_index = snowfall_study.year_to_annual_maxima_tuple_indices_for_daily_time_series[year]
- a = temperature_study.year_to_daily_time_serie_array[year]
- temp = [a[tuple_idx] for tuple_idx in annual_maxima_index]
- maxima_temperatures.append(temp)
- altitude_to_maxima_temperature[altitude] = np.concatenate(maxima_temperatures)
-
- # Plot the box plot
- ax = plt.gca()
- width = 150
- ax.boxplot(altitude_to_maxima_temperature.values(), positions=altitudes, widths=width)
- ax.set_xlim([min(altitudes) - width, max(altitudes) + width])
- ax.set_ylabel('Temperature for maxima of snowfall ({})'.format(temperature_study.variable_class.UNIT))
- ax.set_xlabel('Altitude (m)')
- ax.grid()
-
- plt.show()
-
-
-if __name__ == '__main__':
- distribution_temperature_of_maxima_of_snowfall(altitudes=[900, 1200, 1500, 1800, 2100, 2400, 2700, 3000])