From 0d9a9405c1a4f95c5cfe58aa7b27638276294e41 Mon Sep 17 00:00:00 2001
From: Le Roux Erwan <erwan.le-roux@irstea.fr>
Date: Tue, 7 Apr 2020 17:36:05 +0200
Subject: [PATCH] [contrasting project] improve
 StudyVisualizerForReturnLevelChange.

---
 .../estimator/margin_estimator/utils.py       |  4 +++
 .../abstract_temporal_linear_margin_model.py  | 11 ++++++-
 .../figure1_mean_ratio_return_level_ratio.py  | 32 +++++++++++++------
 3 files changed, 36 insertions(+), 11 deletions(-)

diff --git a/extreme_fit/estimator/margin_estimator/utils.py b/extreme_fit/estimator/margin_estimator/utils.py
index df580526..6617a91f 100644
--- a/extreme_fit/estimator/margin_estimator/utils.py
+++ b/extreme_fit/estimator/margin_estimator/utils.py
@@ -1,3 +1,5 @@
+import warnings
+
 import pandas as pd
 
 from extreme_fit.distribution.gev.gev_params import GevParams
@@ -32,4 +34,6 @@ def fitted_stationary_gev(x_gev, fit_method=TemporalMarginFitMethod.is_mev_gev_f
     estimator = fitted_linear_margin_estimator(model_class, coordinates, dataset, starting_year, fit_method)
     first_coordinate = coordinates.coordinates_values()[0]
     gev_param = estimator.function_from_fit.get_gev_params(first_coordinate)
+    if not -0.5 < gev_param.shape < 0.5:
+        warnings.warn('fitted shape parameter is outside physical bounds {}'.format(gev_param.shape))
     return gev_param
diff --git a/extreme_fit/model/margin_model/linear_margin_model/abstract_temporal_linear_margin_model.py b/extreme_fit/model/margin_model/linear_margin_model/abstract_temporal_linear_margin_model.py
index 467f8906..f097aad5 100644
--- a/extreme_fit/model/margin_model/linear_margin_model/abstract_temporal_linear_margin_model.py
+++ b/extreme_fit/model/margin_model/linear_margin_model/abstract_temporal_linear_margin_model.py
@@ -24,6 +24,10 @@ class TemporalMarginFitMethod(Enum):
     extremes_fevd_l_moments = 4
 
 
+def fitmethod_to_str(fit_method):
+    return str(fit_method).split('.')[-1]
+
+
 class AbstractTemporalLinearMarginModel(LinearMarginModel):
     """Linearity only with respect to the temporal coordinates"""
 
@@ -46,7 +50,8 @@ class AbstractTemporalLinearMarginModel(LinearMarginModel):
                                   df_coordinates_temp: pd.DataFrame) -> AbstractResultFromModelFit:
         data = data[0]
         assert len(data) == len(df_coordinates_temp.values), 'len(data)={} != len(temp)={}'.format(len(data),
-                                                                                                   len(df_coordinates_temp.values))
+                                                                                                   len(
+                                                                                                       df_coordinates_temp.values))
         x = ro.FloatVector(data)
         if self.params_class is GevParams:
             if self.fit_method == TemporalMarginFitMethod.is_mev_gev_fit:
@@ -148,3 +153,7 @@ class AbstractTemporalLinearMarginModel(LinearMarginModel):
     @property
     def siglink(self):
         return r('identity')
+
+
+if __name__ == '__main__':
+    print(fitmethod_to_str(TemporalMarginFitMethod.extremes_fevd_l_moments))
diff --git a/projects/contrasting_trends_in_snow_loads/gorman_figures/figure1_mean_ratio_return_level_ratio.py b/projects/contrasting_trends_in_snow_loads/gorman_figures/figure1_mean_ratio_return_level_ratio.py
index 996c4e72..a41904db 100644
--- a/projects/contrasting_trends_in_snow_loads/gorman_figures/figure1_mean_ratio_return_level_ratio.py
+++ b/projects/contrasting_trends_in_snow_loads/gorman_figures/figure1_mean_ratio_return_level_ratio.py
@@ -3,13 +3,16 @@ from collections import OrderedDict
 import numpy as np
 
 from extreme_data.meteo_france_data.scm_models_data.abstract_study import AbstractStudy
-from extreme_data.meteo_france_data.scm_models_data.safran.safran import SafranSnowfall1Day
+from extreme_data.meteo_france_data.scm_models_data.crocus.crocus import CrocusSnowLoad1Day, CrocusSnowLoad3Days, \
+    CrocusSnowLoadTotal
+from extreme_data.meteo_france_data.scm_models_data.safran.safran import SafranSnowfall1Day, SafranSnowfall3Days
 from extreme_data.meteo_france_data.scm_models_data.visualization.create_shifted_cmap import get_shifted_map, \
     get_colors, shiftedColorMap, ticks_values_and_labels_for_percentages
+from extreme_data.meteo_france_data.scm_models_data.visualization.main_study_visualizer import ALL_ALTITUDES_WITHOUT_NAN
 from extreme_data.meteo_france_data.scm_models_data.visualization.study_visualizer import StudyVisualizer
 from extreme_fit.estimator.margin_estimator.utils import fitted_stationary_gev
 from extreme_fit.model.margin_model.linear_margin_model.abstract_temporal_linear_margin_model import \
-    TemporalMarginFitMethod
+    TemporalMarginFitMethod, fitmethod_to_str
 
 import matplotlib.pyplot as plt
 
@@ -54,11 +57,15 @@ class StudyVisualizerForReturnLevelChange(StudyVisualizer):
             massif_name_to_return_levels[massif_name] = return_level
         return massif_name_to_return_levels
 
-    def plot_return_level_percentage(self):
+    def plot_return_level_evolution(self):
         massif_name_to_return_level_past = self.massif_name_to_return_level(self.study_before)
         massif_name_to_return_level_recent = self.massif_name_to_return_level(self.study_after)
+        # massif_name_to_percentage = {
+        #     m: 100 * (v - massif_name_to_return_level_past[m]) / massif_name_to_return_level_past[m]
+        #     for m, v in
+        #     massif_name_to_return_level_recent.items()}
         massif_name_to_percentage = {
-            m: 100 * (v - massif_name_to_return_level_past[m]) / massif_name_to_return_level_past[m]
+            m:  (v - massif_name_to_return_level_past[m])
             for m, v in
             massif_name_to_return_level_recent.items()}
 
@@ -95,6 +102,7 @@ class StudyVisualizerForReturnLevelChange(StudyVisualizer):
         ax.get_xaxis().set_visible(True)
         ax.set_xticks([])
         ax.set_xlabel('Altitude = {}m'.format(self.study.altitude), fontsize=15)
+        ax.set_title('Fit method is {}'.format(fitmethod_to_str(self.fit_method)))
         self.plot_name = 'change in return levels'
         self.show_or_save_to_file(add_classic_title=False, tight_layout=True, no_title=True,
                                   dpi=500)
@@ -102,10 +110,14 @@ class StudyVisualizerForReturnLevelChange(StudyVisualizer):
 
 
 if __name__ == '__main__':
-    for altitude in [900, 1200, 1500, 1800, 2100, 2400, 2700, 3000][:]:
-        study_visualizer = StudyVisualizerForReturnLevelChange(study_class=SafranSnowfall1Day, altitude=altitude,
-                                                               return_period=30,
-                                                               save_to_file=True,
-                                                               fit_method=TemporalMarginFitMethod.extremes_fevd_l_moments)
-        study_visualizer.plot_return_level_percentage()
+    # for altitude in [900, 1200, 1500, 1800, 2100, 2400, 2700, 3000][:]:
+    for altitude in ALL_ALTITUDES_WITHOUT_NAN:
+        for study_class in [SafranSnowfall1Day, SafranSnowfall3Days, CrocusSnowLoadTotal, CrocusSnowLoad3Days, CrocusSnowLoad1Day][:1]:
+            return_period = 50
+            study_visualizer = StudyVisualizerForReturnLevelChange(study_class=study_class,
+                                                                   altitude=altitude,
+                                                                   return_period=return_period,
+                                                                   save_to_file=True,
+                                                                   fit_method=TemporalMarginFitMethod.extremes_fevd_l_moments)
+            study_visualizer.plot_return_level_evolution()
 
-- 
GitLab