[PAPER 1] fix the code for taking the argmax of some value. add some new...

[PAPER 1] fix the code for taking the argmax of some value. add some new classes for making comparison between various snow variables

experiment/meteo_france_data/scm_models_data/visualization/study_visualization/main_study_visualizer.py

@@ -3,8 +3,10 @@ from typing import List
 
 from experiment.meteo_france_data.scm_models_data.visualization.study_visualization.study_visualizer import \
     StudyVisualizer
-from experiment.paper_past_snow_loads.result_data_comparison_with_eurocode.crocus_study_comparison_with_eurocode import CrocusDifferenceSnowLoad, \
-    CrocusSnowDensityAtMaxofSwe, CrocusDifferenceSnowLoadRescaledAndEurocodeToSeeSynchronization
+from experiment.paper_past_snow_loads.result_data_comparison_with_eurocode.crocus_study_comparison_with_eurocode import \
+    CrocusDifferenceSnowLoad, \
+    CrocusSnowDensityAtMaxofSwe, CrocusDifferenceSnowLoadRescaledAndEurocodeToSeeSynchronization, \
+    CrocusSnowDepthDifference, CrocusSnowDepthAtMaxofSwe
 from experiment.trend_analysis.abstract_score import MannKendall
 from experiment.meteo_france_data.scm_models_data.abstract_study import AbstractStudy
 from experiment.meteo_france_data.scm_models_data.crocus.crocus import CrocusDepth, CrocusSweTotal, ExtendedCrocusDepth, \
@@ -32,13 +34,15 @@ SCM_STUDY_CLASS_TO_ABBREVIATION = {
     SafranSnowfall: 'SF3',
     CrocusSweTotal: 'SWE',
     CrocusSwe3Days: 'SWE3',
-    CrocusSnowLoadEurocode: 'SL_Eurocode',
+    CrocusSnowLoadEurocode: 'GSL_Eurocode',
     CrocusDepth: 'SD',
-    CrocusSnowLoadTotal: 'max SL',
-    CrocusSnowLoad3Days: 'SL3',
-    CrocusSnowDensityAtMaxofSwe: '{} when the max of SL \nis reached'.format(snow_density_str),
-    CrocusDifferenceSnowLoadRescaledAndEurocodeToSeeSynchronization: 'max SL rescaled - SL from max HS \nboth with {}'.format(eurocode_snow_density),
-    CrocusDifferenceSnowLoad: ('max SL - SL from max HS \n with {}'.format(eurocode_snow_density))
+    CrocusSnowLoadTotal: 'max GSL',
+    CrocusSnowLoad3Days: 'GSL3',
+    CrocusSnowDensityAtMaxofSwe: '{} when the max of GSL \nis reached'.format(snow_density_str),
+    CrocusDifferenceSnowLoadRescaledAndEurocodeToSeeSynchronization: 'max GSL rescaled - GSL from max HS \nboth with {}'.format(eurocode_snow_density),
+    CrocusDifferenceSnowLoad: ('max GSL - GSL from max HS \n with {}'.format(eurocode_snow_density)),
+    CrocusSnowDepthDifference: 'max HS - HS at max of GSL',
+    CrocusSnowDepthAtMaxofSwe: 'HS at max of GSL',
 }
 
 altitude_massif_name_and_study_class_for_poster = [

experiment/paper_past_snow_loads/result_data_comparison_with_eurocode/crocus_study_comparison_with_eurocode.py

 from collections import OrderedDict
 
 import numpy as np
+from cached_property import cached_property
 
 from experiment.meteo_france_data.scm_models_data.crocus.crocus import Crocus, CrocusSweTotal, CrocusSnowLoadTotal, \
-    CrocusSnowLoadEurocode
+    CrocusSnowLoadEurocode, CrocusDepth
 from experiment.meteo_france_data.scm_models_data.crocus.crocus_variables import TotalSnowLoadVariable, \
-    CrocusDensityVariable
+    CrocusDensityVariable, CrocusDepthVariable
+
+
+class CrocusSnowDepthAtMaxofSwe(Crocus):
+
+    def __init__(self, *args, **kwargs):
+        super().__init__(CrocusDepthVariable, *args, **kwargs)
+        self.year_to_snow_depth_at_max_of_swe = OrderedDict()
+        study_swe = CrocusSweTotal(*args, **kwargs)
+        for year in study_swe.ordered_years:
+            argmax_swe = study_swe.year_to_annual_maxima_index[year]
+            snow_depth = self.year_to_daily_time_serie_array[year]
+            snow_depth_at_max = [snow_depth[argmax, i] for i, argmax in enumerate(argmax_swe)]
+            self.year_to_snow_depth_at_max_of_swe[year] = np.array(snow_depth_at_max)
+
+    @cached_property
+    def year_to_annual_maxima(self) -> OrderedDict:
+        return self.year_to_snow_depth_at_max_of_swe
+
+
+class CrocusSnowDepthDifference(Crocus):
+
+    def __init__(self, *args, **kwargs):
+        super().__init__(CrocusDepthVariable, *args, **kwargs)
+        snow_depth_at_max_swe = CrocusSnowDepthAtMaxofSwe(*args, **kwargs)
+        snow_depth_max = CrocusDepth(*args, **kwargs)
+        self.year_to_diff = OrderedDict()
+        for year in snow_depth_max.ordered_years:
+            self.year_to_diff[year] = snow_depth_max.year_to_annual_maxima[year] \
+                                      - snow_depth_at_max_swe.year_to_annual_maxima[year]
+
+    @property
+    def year_to_annual_maxima(self) -> OrderedDict:
+        return self.year_to_diff
 
 
 class CrocusSnowDensityAtMaxofSwe(Crocus):
@@ -18,7 +52,8 @@ class CrocusSnowDensityAtMaxofSwe(Crocus):
             max_swe = study_swe.year_to_annual_maxima[year]
             argmax_swe = study_swe.year_to_annual_maxima_index[year]
             snow_depth = self.year_to_daily_time_serie_array[year]
-            snow_depth_at_max = np.take(np.transpose(snow_depth), argmax_swe)
+            snow_depth_at_max = np.array([snow_depth[argmax, i] for i, argmax in enumerate(argmax_swe)])
+            # todo: should i take into acocunt the density of the water
             self.year_to_snow_density_at_max_of_swe[year] = max_swe / snow_depth_at_max
 
     @property

experiment/paper_past_snow_loads/result_data_comparison_with_eurocode/main_comparison_with_eurocode.py

+from experiment.meteo_france_data.scm_models_data.crocus.crocus import CrocusDepth
+from experiment.meteo_france_data.scm_models_data.crocus.crocus_variables import CrocusDepthVariable
 from experiment.meteo_france_data.scm_models_data.visualization.study_visualization.main_study_visualizer import \
     study_iterator_global, SCM_STUDY_CLASS_TO_ABBREVIATION, snow_density_str
 from experiment.meteo_france_data.scm_models_data.visualization.study_visualization.study_visualizer import \
     StudyVisualizer
 import matplotlib.pyplot as plt
 
-from experiment.paper_past_snow_loads.result_data_comparison_with_eurocode.crocus_study_comparison_with_eurocode import CrocusDifferenceSnowLoad, \
-    CrocusSnowDensityAtMaxofSwe, CrocusDifferenceSnowLoadRescaledAndEurocodeToSeeSynchronization
+from experiment.paper_past_snow_loads.result_data_comparison_with_eurocode.crocus_study_comparison_with_eurocode import \
+    CrocusDifferenceSnowLoad, \
+    CrocusSnowDensityAtMaxofSwe, CrocusDifferenceSnowLoadRescaledAndEurocodeToSeeSynchronization, \
+    CrocusSnowDepthAtMaxofSwe, CrocusSnowDepthDifference
 
 
 def max_graph_annual_maxima_comparison():
@@ -15,14 +19,20 @@ def max_graph_annual_maxima_comparison():
     :return:
     """
     save_to_file = True
-    study_classes = [CrocusSnowDensityAtMaxofSwe, CrocusDifferenceSnowLoadRescaledAndEurocodeToSeeSynchronization, CrocusDifferenceSnowLoad][:]
+    study_classes = [CrocusSnowDensityAtMaxofSwe,
+                     # CrocusDifferenceSnowLoadRescaledAndEurocodeToSeeSynchronization,
+                     CrocusDifferenceSnowLoad,
+                     # CrocusDepth,
+                     # CrocusSnowDepthAtMaxofSwe,
+                     CrocusSnowDepthDifference,
+                     ][:]
     for study_class in study_classes:
 
         marker_altitude_massif_name = [
             ('magenta', 900, 'Ubaye'),
             ('darkmagenta', 1800, 'Vercors'),
             ('mediumpurple', 2700, 'Beaufortain'),
-        ]
+        ][:]
         ax = plt.gca()
         for color, altitude, massif_name in marker_altitude_massif_name:
             for study in study_iterator_global([study_class], altitudes=[altitude]):
@@ -36,7 +46,8 @@ def max_graph_annual_maxima_comparison():
                 last_plot = altitude == 2700
                 if last_plot:
                     constant = 150 if study_class == CrocusSnowDensityAtMaxofSwe else 0
-                    label = '{} Eurocode'.format(snow_density_str) if study_class == CrocusSnowDensityAtMaxofSwe else None
+                    label = '{} Eurocode'.format(
+                        snow_density_str) if study_class == CrocusSnowDensityAtMaxofSwe else None
                     snow_density_eurocode = [constant for _ in study.ordered_years]
                     ax.plot(study.ordered_years, snow_density_eurocode, color='k', label=label)
                     ax.legend()
@@ -45,9 +56,5 @@ def max_graph_annual_maxima_comparison():
                     ax.clear()
 
 
-
-
-
-
 if __name__ == '__main__':
-    max_graph_annual_maxima_comparison()
+    max_graph_annual_maxima_comparison()
\ No newline at end of file

experiment/paper_past_snow_loads/result_trends_and_return_levels/eurocode_visualizer.py

@@ -3,7 +3,7 @@ from typing import Dict, List, Tuple
 import matplotlib.pyplot as plt
 import numpy as np
 
-from experiment.eurocode_data.utils import EUROCODE_RETURN_LEVEL_STR
+from experiment.eurocode_data.utils import EUROCODE_RETURN_LEVEL_STR, EUROCODE_ALTITUDES
 from experiment.paper_past_snow_loads.result_trends_and_return_levels.study_visualizer_for_non_stationary_trends import \
     StudyVisualizerForNonStationaryTrends
 from extreme_fit.model.result_from_model_fit.result_from_extremes.abstract_extract_eurocode_return_level import \
@@ -19,11 +19,12 @@ def plot_uncertainty_massifs(altitude_to_visualizer: Dict[int, StudyVisualizerFo
     """ Plot several uncertainty plots
     :return:
     """
+    altitude_to_visualizer = {a:v for a,v in altitude_to_visualizer.items() if a in EUROCODE_ALTITUDES}
     visualizer = list(altitude_to_visualizer.values())[-1]
     # Subdivide massif names in group of 3
     m = 1
     uncertainty_massif_names = visualizer.uncertainty_massif_names
-    n = (len(uncertainty_massif_names) // m) + 1
+    n = (len(uncertainty_massif_names) // m)
     print('total nb of massif', n)
     for i in list(range(n))[:]:
         massif_names = uncertainty_massif_names[m * i: m * (i + 1)]

experiment/paper_past_snow_loads/result_trends_and_return_levels/main_result_trends_and_return_levels.py

@@ -9,7 +9,8 @@ from experiment.meteo_france_data.scm_models_data.visualization.study_visualizat
     StudyVisualizer
 from experiment.paper_past_snow_loads.result_trends_and_return_levels.eurocode_visualizer import \
     plot_uncertainty_massifs
-from experiment.meteo_france_data.scm_models_data.crocus.crocus import CrocusSnowLoadTotal
+from experiment.meteo_france_data.scm_models_data.crocus.crocus import CrocusSnowLoadTotal, CrocusSnowLoadEurocode, \
+    CrocusSnowLoad3Days
 from experiment.paper_past_snow_loads.result_trends_and_return_levels.study_visualizer_for_non_stationary_trends import \
     StudyVisualizerForNonStationaryTrends
 from extreme_fit.model.result_from_model_fit.result_from_extremes.confidence_interval_method import \
@@ -57,17 +58,18 @@ def intermediate_result(altitudes, massif_names=None,
 
 
 def major_result():
-    altitudes = [[1500, 1800]][0]
+    altitudes = [300, 600, 900, 1200, 1500, 1800, 2100, 2400, 2700]
     uncertainty_methods = [ConfidenceIntervalMethodFromExtremes.my_bayes,
-                           ConfidenceIntervalMethodFromExtremes.ci_mle][1:]
-    massif_names = ['Chartreuse']
-    non_stationary_models_for_uncertainty = [False, True][:1]
-    #
-    # altitudes
-    # study_class
+                           ConfidenceIntervalMethodFromExtremes.ci_mle][:]
+    massif_names = None
+    non_stationary_uncertainty = [False, True][:]
+    study_classes = [CrocusSnowLoadTotal, CrocusSnowLoadEurocode, CrocusSnowLoad3Days]
+    for study_class in study_classes[2:]:
+        intermediate_result(altitudes, massif_names, non_stationary_uncertainty, uncertainty_methods, study_class)
 
 
 if __name__ == '__main__':
+    major_result()
     # minor_result(altitude=1800)
     # intermediate_result(altitudes=[1500, 1800], massif_names=['Chartreuse'],
     #                     uncertainty_methods=[ConfidenceIntervalMethodFromExtremes.ci_mle,
@@ -79,10 +81,10 @@ if __name__ == '__main__':
     # intermediate_result(altitudes=[300, 600, 900, 1200, 1500, 1800], massif_names=None,
     #                     uncertainty_methods=[ConfidenceIntervalMethodFromExtremes.ci_mle],
     #                     non_stationary_uncertainty=[False])
-    # intermediate_result(altitudes=[300, 600, 900, 1200, 1500, 1800], massif_names=None,
+    # intermediate_result(altitudes=[300, 600, 900, 1200, 1500, 1800, 2100, 2400, 2700], massif_names=None,
     #                     uncertainty_methods=[ConfidenceIntervalMethodFromExtremes.ci_mle,
     #                                          ConfidenceIntervalMethodFromExtremes.ci_bayes],
     #                     non_stationary_uncertainty=[False, True])
-    intermediate_result(altitudes=[300, 600, 900], massif_names=None,
-                        uncertainty_methods=[ConfidenceIntervalMethodFromExtremes.ci_mle],
-                        non_stationary_uncertainty=[False, True])
+    # intermediate_result(altitudes=[300, 600, 900], massif_names=None,
+    #                     uncertainty_methods=[ConfidenceIntervalMethodFromExtremes.ci_mle],
+    #                     non_stationary_uncertainty=[False, True])