From e96c7ccc51442510141b165aaa47633b98e90eef Mon Sep 17 00:00:00 2001
From: Le Roux Erwan <erwan.le-roux@irstea.fr>
Date: Wed, 3 Feb 2021 19:26:55 +0100
Subject: [PATCH] [contrasting] add gap between study. and add safran 2019 from
 max files. add main_comparison_histo.py

---
 .../scm_models_data/abstract_study.py         | 34 +++++++----
 .../safran/gap_between_study.py               | 52 ++++++++++++++++
 .../scm_models_data/safran/safran.py          | 23 +++++--
 .../scm_models_data/safran/safran_2020.py     | 34 -----------
 .../safran/safran_max_snowf.py                | 46 ++++++++++++++
 .../scm_models_data/safran/safran_variable.py | 39 +++++++-----
 .../visualization/main_study_visualizer.py    |  9 +++
 .../altitudes_fit/main_altitudes_studies.py   | 40 ++++++------
 .../main_gap_altitudes_studies.py             | 61 +++++++++++++++++++
 .../one_fold_analysis/plot_total_aic.py       |  4 +-
 .../test_meteo_france_data/test_SCM_study.py  |  4 +-
 11 files changed, 256 insertions(+), 90 deletions(-)
 create mode 100644 extreme_data/meteo_france_data/scm_models_data/safran/gap_between_study.py
 delete mode 100644 extreme_data/meteo_france_data/scm_models_data/safran/safran_2020.py
 create mode 100644 extreme_data/meteo_france_data/scm_models_data/safran/safran_max_snowf.py
 create mode 100644 projects/altitude_spatial_model/altitudes_fit/main_gap_altitudes_studies.py

diff --git a/extreme_data/meteo_france_data/scm_models_data/abstract_study.py b/extreme_data/meteo_france_data/scm_models_data/abstract_study.py
index 7a0a0f7e..9d70c864 100644
--- a/extreme_data/meteo_france_data/scm_models_data/abstract_study.py
+++ b/extreme_data/meteo_france_data/scm_models_data/abstract_study.py
@@ -546,19 +546,18 @@ class AbstractStudy(object):
         # extracted for a csv file, and used only for display purposes
         df = cls.load_df_centroid()
         # Lower a bit the Mercantour massif
-        df.loc['Mercantour', 'coord_x'] += 14000 # shift to the right
-        df.loc['Mercantour', 'coord_y'] -= 7000 # shift down
+        df.loc['Mercantour', 'coord_x'] += 14000  # shift to the right
+        df.loc['Mercantour', 'coord_y'] -= 7000  # shift down
         # Lower a bit the Maurienne massif
         # df.loc['Mercantour', 'coord_x'] += 14000 # shift to the right
-        df.loc['Maurienne', 'coord_y'] -= 6000 # shift down
-        df.loc['Maurienne', 'coord_y'] -= 5000 # shift down
-        df.loc['Maurienne', 'coord_x'] += 3000 # shift down
-        df.loc['Vanoise', 'coord_y'] -= 4000 # shift down
-        df.loc['Ubaye', 'coord_y'] -= 4000 # shift down
+        df.loc['Maurienne', 'coord_y'] -= 6000  # shift down
+        df.loc['Maurienne', 'coord_y'] -= 5000  # shift down
+        df.loc['Maurienne', 'coord_x'] += 3000  # shift down
+        df.loc['Vanoise', 'coord_y'] -= 4000  # shift down
+        df.loc['Ubaye', 'coord_y'] -= 4000  # shift down
         # Filter, keep massifs present at the altitude of interest
         df = df.loc[massif_names, :]
 
-
         # Build coordinate object from df_centroid
         return AbstractSpatialCoordinates.from_df(df)
 
@@ -889,7 +888,8 @@ class AbstractStudy(object):
             self._cache_for_pointwise_fit[t] = res
             return res
 
-    def _massif_name_to_stationary_gev_params_and_confidence(self, quantile_level=None, confidence_interval_based_on_delta_method=True):
+    def _massif_name_to_stationary_gev_params_and_confidence(self, quantile_level=None,
+                                                             confidence_interval_based_on_delta_method=True):
         """ at least 90% of values must be above zero"""
         print('study computation')
         massif_name_to_stationary_gev_params = {}
@@ -899,16 +899,17 @@ class AbstractStudy(object):
             if percentage_of_non_zeros > 90:
                 start = time.time()
                 gev_params, mean_estimate, confidence = fitted_stationary_gev_with_uncertainty_interval(annual_maxima,
-                                                                                         fit_method=MarginFitMethod.extremes_fevd_mle,
-                                                                                         quantile_level=quantile_level,
-                                                                                         confidence_interval_based_on_delta_method=confidence_interval_based_on_delta_method)
+                                                                                                        fit_method=MarginFitMethod.extremes_fevd_mle,
+                                                                                                        quantile_level=quantile_level,
+                                                                                                        confidence_interval_based_on_delta_method=confidence_interval_based_on_delta_method)
                 end = time.time()
                 duration = end - start
                 print('Multiprocessing for study duration', duration)
 
                 if -0.5 <= gev_params.shape <= 0.5:
                     massif_name_to_stationary_gev_params[massif_name] = gev_params
-                    massif_name_to_confidence[massif_name] = EurocodeConfidenceIntervalFromExtremes(mean_estimate, confidence)
+                    massif_name_to_confidence[massif_name] = EurocodeConfidenceIntervalFromExtremes(mean_estimate,
+                                                                                                    confidence)
         return massif_name_to_stationary_gev_params, massif_name_to_confidence
 
     def massif_name_to_gev_param_list(self, year_min_and_max_list):
@@ -926,3 +927,10 @@ class AbstractStudy(object):
                 gev_params_list.append(fitted_stationary_gev(annual_maxima_sliced))
             massif_name_to_stationary_gev_params[massif_name] = gev_params_list
         return massif_name_to_stationary_gev_params
+
+    @property
+    def mean_annual_maxima(self):
+        annual_maxima = []
+        for maxima in self.year_to_annual_maxima.values():
+            annual_maxima.extend(maxima)
+        return np.mean(annual_maxima)
diff --git a/extreme_data/meteo_france_data/scm_models_data/safran/gap_between_study.py b/extreme_data/meteo_france_data/scm_models_data/safran/gap_between_study.py
new file mode 100644
index 00000000..a07f6be3
--- /dev/null
+++ b/extreme_data/meteo_france_data/scm_models_data/safran/gap_between_study.py
@@ -0,0 +1,52 @@
+from collections import OrderedDict
+
+from cached_property import cached_property
+
+from extreme_data.meteo_france_data.scm_models_data.safran.safran import SafranSnowfall1Day, \
+    SafranSnowfallCenterOnDay1day
+from extreme_data.meteo_france_data.scm_models_data.safran.safran_max_snowf import SafranSnowfall2020, \
+    SafranSnowfall2019
+
+
+class AbstractGapBetweenTwoStudyClass(SafranSnowfall1Day):
+
+    def __init__(self, study_class1, study_class2, **kwargs):
+        super().__init__(**kwargs)
+        self.study_1 = study_class1(**kwargs)
+        self.study_2 = study_class2(**kwargs)
+
+    @property
+    def ordered_years(self):
+        return [i for i in list(range(1959, 2020)) if self.year_min <= i <= self.year_max]
+
+    @cached_property
+    def year_to_annual_maxima(self) -> OrderedDict:
+        year_to_annual_maxima = OrderedDict()
+        for year in self.ordered_years:
+            annual_maxima_1 = self.study_1.year_to_annual_maxima[year]
+            annual_maxima_2 = self.study_2.year_to_annual_maxima[year]
+            year_to_annual_maxima[year] = annual_maxima_2 - annual_maxima_1
+        return year_to_annual_maxima
+
+
+class GapBetweenSafranSnowfall2019And2020(AbstractGapBetweenTwoStudyClass):
+
+    def __init__(self, **kwargs):
+        super().__init__(SafranSnowfall2019, SafranSnowfall2020, **kwargs)
+
+
+class GapBetweenSafranSnowfall2019AndMySafranSnowfall2019Recentered(AbstractGapBetweenTwoStudyClass):
+
+    def __init__(self, **kwargs):
+        super().__init__(SafranSnowfall2019, SafranSnowfallCenterOnDay1day, **kwargs)
+
+
+class GapBetweenSafranSnowfall2019AndMySafranSnowfall2019(AbstractGapBetweenTwoStudyClass):
+
+    def __init__(self, **kwargs):
+        super().__init__(SafranSnowfall2019, SafranSnowfall1Day, **kwargs)
+
+
+if __name__ == '__main__':
+    study = GapBetweenSafranSnowfall2019AndMySafranSnowfall2019Recentered(altitude=1800)
+    study.year_to_annual_maxima[1960]
diff --git a/extreme_data/meteo_france_data/scm_models_data/safran/safran.py b/extreme_data/meteo_france_data/scm_models_data/safran/safran.py
index b10fa0b8..29ff2124 100644
--- a/extreme_data/meteo_france_data/scm_models_data/safran/safran.py
+++ b/extreme_data/meteo_france_data/scm_models_data/safran/safran.py
@@ -10,11 +10,12 @@ from extreme_data.meteo_france_data.scm_models_data.safran.cumulated_study impor
 from extreme_data.meteo_france_data.scm_models_data.safran.safran_variable import SafranSnowfallVariable, \
     SafranRainfallVariable, SafranTemperatureVariable, SafranTotalPrecipVariable, \
     SafranNormalizedPrecipitationRateOnWetDaysVariable, SafranNormalizedPrecipitationRateVariable, \
-    SafranDateFirstSnowfallVariable
+    SafranDateFirstSnowfallVariable, SafranSnowfallVariableCenterOnDay
 
 
 class Safran(AbstractStudy):
     SAFRAN_VARIABLES = [SafranSnowfallVariable,
+                        SafranSnowfallVariableCenterOnDay,
                         SafranRainfallVariable,
                         SafranTemperatureVariable,
                         SafranTotalPrecipVariable,
@@ -43,6 +44,18 @@ class SafranSnowfall1Day(SafranSnowfall):
         super().__init__(nb_consecutive_days=1, **kwargs)
 
 
+class SafranSnowfallCenterOnDay(Safran, CumulatedStudy):
+
+    def __init__(self, **kwargs):
+        super().__init__(SafranSnowfallVariableCenterOnDay, **kwargs)
+
+
+class SafranSnowfallCenterOnDay1day(SafranSnowfallCenterOnDay):
+
+    def __init__(self, **kwargs):
+        super().__init__(nb_consecutive_days=1, **kwargs)
+
+
 class SafranDateFirstSnowfall(Safran, CumulatedStudy):
 
     def __init__(self, **kwargs):
@@ -211,9 +224,7 @@ if __name__ == '__main__':
     altitude = 900
     year_min = 1959
     year_max = 2019
-    study = SafranSnowfall(altitude=altitude, year_min=year_min, year_max=year_max)
-    print(study.study_massif_names)
-    # print(study.massif_name_to_annual_maxima)
-    # print(study.year_to_daily_time_serie_array[1959].shape)
+    # study = SafranSnowfall(altitude=altitude, year_min=year_min, year_max=year_max)
+    study = SafranSnowfallCenterOnDay1day(altitude=altitude, year_min=year_min, year_max=year_max)
+    print(study.year_to_daily_time_serie_array[1959].shape)
     # print(study.massif_name_to_daily_time_series['Vanoise'].shape)
-    study._save_excel_with_longitutde_and_latitude()
diff --git a/extreme_data/meteo_france_data/scm_models_data/safran/safran_2020.py b/extreme_data/meteo_france_data/scm_models_data/safran/safran_2020.py
deleted file mode 100644
index 78157550..00000000
--- a/extreme_data/meteo_france_data/scm_models_data/safran/safran_2020.py
+++ /dev/null
@@ -1,34 +0,0 @@
-from collections import OrderedDict
-
-import numpy as np
-from cached_property import cached_property
-from netCDF4._netCDF4 import Dataset
-
-from extreme_data.meteo_france_data.scm_models_data.safran.safran import SafranSnowfall1Day
-
-
-class SnowfallSafran2020(SafranSnowfall1Day):
-    nc_filepath = """/home/erwan/Documents/projects/spatiotemporalextremes/local/spatio_temporal_datasets/SAFRAN_montagne_CROCUS_2020/max-1day-snowf_SAFRAN.nc"""
-
-    def __init__(self, **kwargs):
-        super().__init__(**kwargs)
-
-    @property
-    def ordered_years(self):
-        return list(range(1959, 2021))
-
-    @cached_property
-    def year_to_annual_maxima(self) -> OrderedDict:
-        year_to_annual_maxima = OrderedDict()
-        dataset = Dataset(SnowfallSafran2020.nc_filepath)
-        annual_maxima = np.array(dataset.variables['max-1day-snowf'])
-        assert annual_maxima.shape[1] == len(self.column_mask)
-        annual_maxima = annual_maxima[:, self.column_mask]
-        for year, a in zip(self.ordered_years, annual_maxima):
-            year_to_annual_maxima[year] = a
-        return year_to_annual_maxima
-
-
-
-if __name__ == '__main__':
-    test_safran_2020_loader()
diff --git a/extreme_data/meteo_france_data/scm_models_data/safran/safran_max_snowf.py b/extreme_data/meteo_france_data/scm_models_data/safran/safran_max_snowf.py
new file mode 100644
index 00000000..724304dd
--- /dev/null
+++ b/extreme_data/meteo_france_data/scm_models_data/safran/safran_max_snowf.py
@@ -0,0 +1,46 @@
+from collections import OrderedDict
+import os.path as op
+
+import numpy as np
+from cached_property import cached_property
+from netCDF4._netCDF4 import Dataset
+
+from extreme_data.meteo_france_data.scm_models_data.safran.safran import SafranSnowfall1Day, Safran
+
+
+class AbstractSafranSnowfallMaxFiles(SafranSnowfall1Day):
+    path = """/home/erwan/Documents/projects/spatiotemporalextremes/local/spatio_temporal_datasets"""
+
+    def __init__(self, safran_year, **kwargs):
+        super().__init__(**kwargs)
+        self.year_max = max(safran_year, self.year_max)
+        self.nc_filepath = op.join(self.path, 'SAFRAN_montagne-CROCUS_{}'.format(safran_year),
+                                   'max-1day-snowf_SAFRAN.nc')
+        self.safran_year = safran_year
+
+    @property
+    def ordered_years(self):
+        return [i for i in list(range(1959, self.safran_year+1)) if self.year_min <= i <= self.year_max]
+
+    @cached_property
+    def year_to_annual_maxima(self) -> OrderedDict:
+        year_to_annual_maxima = OrderedDict()
+        dataset = Dataset(self.nc_filepath)
+        annual_maxima = np.array(dataset.variables['max-1day-snowf'])
+        assert annual_maxima.shape[1] == len(self.column_mask)
+        annual_maxima = annual_maxima[:, self.column_mask]
+        for year, a in zip(self.ordered_years, annual_maxima):
+            year_to_annual_maxima[year] = a
+        return year_to_annual_maxima
+
+
+class SafranSnowfall2020(AbstractSafranSnowfallMaxFiles):
+
+    def __init__(self, **kwargs):
+        super().__init__(2020, **kwargs)
+
+
+class SafranSnowfall2019(AbstractSafranSnowfallMaxFiles):
+
+    def __init__(self, **kwargs):
+        super().__init__(2019, **kwargs)
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 477f2101..5c60111d 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
@@ -14,7 +14,7 @@ class SafranSnowfallVariable(AbstractVariable):
         and multiply that by 60 x 60 which corresponds to the number of seconds in one hour
 
     -How do how I define the limit of a day ?
-        From the start, i.e. in August at 4am something like that,then if I add a 24H duration, I arrive to the next day
+        From the start, i.e. in August at 6am,then if I add a 24H duration, I arrive to the next day
 
     -How do you aggregate several days ?
         We aggregate all the N consecutive days into a value x_i, then we take the max
@@ -34,21 +34,17 @@ class SafranSnowfallVariable(AbstractVariable):
         super().__init__(variable_array)
         self.nb_consecutive_days_of_snowfall = nb_consecutive_days
         # Compute the daily snowfall in kg/m2
-        snowfall_rates = variable_array
-
-        # Compute the mean snowrate, then multiply it by 60 * 60 * 24
-        # day_duration_in_seconds = 24 * 60 * 60
-        # nb_days = len(snowfall_rates) // 24
-        # print(nb_days)
-        # daily_snowrate = [np.mean(snowfall_rates[24 * i:24 * (i + 1) + 1], axis=0) for i in range(nb_days)]
-        # self.daily_snowfall = day_duration_in_seconds * np.array(daily_snowrate)
-
-        # Compute the hourly snowfall first, then aggregate
-        mean_snowfall_rates = 0.5 * (snowfall_rates[:-1] + snowfall_rates[1:])
+        mean_snowfall_rates = self.get_snowfall_rates(variable_array)
         hourly_snowfall = 60 * 60 * mean_snowfall_rates
         # Transform the snowfall amount into a dataframe
         nb_days = len(hourly_snowfall) // 24
-        self.daily_snowfall = [sum(hourly_snowfall[24 * i:24 * (i + 1)]) for i in range(nb_days)]
+        self.daily_snowfall = self.daily_snowfall(hourly_snowfall, nb_days)
+
+    def get_snowfall_rates(self, variable_array):
+        return 0.5 * (variable_array[:-1] + variable_array[1:])
+
+    def daily_snowfall(self, hourly_snowfall, nb_days):
+        return [sum(hourly_snowfall[24 * i:24 * (i + 1)]) for i in range(nb_days)]
 
     @property
     def daily_time_serie_array(self) -> np.ndarray:
@@ -67,6 +63,20 @@ class SafranSnowfallVariable(AbstractVariable):
         return snowfall_in_consecutive_days
 
 
+class SafranSnowfallVariableCenterOnDay(SafranSnowfallVariable):
+
+    def daily_snowfall(self, hourly_snowfall, nb_days):
+        hourly_snowfall_without_first_and_last_days = hourly_snowfall[18:-6]
+        assert len(hourly_snowfall_without_first_and_last_days) % 24 == 0
+        daily_snowfall = super().daily_snowfall(hourly_snowfall[18:-5], nb_days - 2)
+        zero_array = daily_snowfall[0] * 0
+        daily_snowfall = [zero_array] + daily_snowfall + [zero_array]
+        return daily_snowfall
+
+    def get_snowfall_rates(self, variable_array):
+        return variable_array[:-1]
+
+
 class SafranDateFirstSnowfallVariable(SafranSnowfallVariable):
     NAME = 'Date First Snow'
     UNIT = 'days'
@@ -83,11 +93,12 @@ class SafranDateFirstSnowfallVariable(SafranSnowfallVariable):
                 min = np.nan
             # first_date = 1 - min / 366
             first_date = 1 - min / 366
-            first_date_repeated = np.ones(len(s))  * first_date
+            first_date_repeated = np.ones(len(s)) * first_date
             new_daily_time_series.append(first_date_repeated)
         new_daily_time_series_array = np.array(new_daily_time_series).transpose()
         return new_daily_time_series_array
 
+
 class SafranRainfallVariable(SafranSnowfallVariable):
     """Warning: this corresponds to water falling. Total precipitaiton equals Rainfall + Snowfall"""
     NAME = 'Rainfall'
diff --git a/extreme_data/meteo_france_data/scm_models_data/visualization/main_study_visualizer.py b/extreme_data/meteo_france_data/scm_models_data/visualization/main_study_visualizer.py
index 8ddf624c..f0b11500 100644
--- a/extreme_data/meteo_france_data/scm_models_data/visualization/main_study_visualizer.py
+++ b/extreme_data/meteo_france_data/scm_models_data/visualization/main_study_visualizer.py
@@ -10,11 +10,15 @@ from extreme_data.meteo_france_data.scm_models_data.crocus.crocus import CrocusD
     CrocusSnowLoadEurocode, CrocusSnowLoad5Days, CrocusSnowLoad7Days
 from extreme_data.meteo_france_data.scm_models_data.crocus.crocus_snow_density import CrocusSnowDensity
 from extreme_data.meteo_france_data.scm_models_data.crocus.crocus_variables import CrocusDensityVariable
+from extreme_data.meteo_france_data.scm_models_data.safran.gap_between_study import GapBetweenSafranSnowfall2019And2020, \
+    GapBetweenSafranSnowfall2019AndMySafranSnowfall2019Recentered, GapBetweenSafranSnowfall2019AndMySafranSnowfall2019
 from extreme_data.meteo_france_data.scm_models_data.safran.safran import SafranSnowfall, ExtendedSafranSnowfall, \
     SafranRainfall, \
     SafranTemperature, SafranPrecipitation, SafranSnowfall1Day, SafranSnowfall3Days, SafranSnowfall5Days, \
     SafranSnowfall7Days, SafranPrecipitation1Day, SafranPrecipitation3Days, SafranPrecipitation5Days, \
     SafranPrecipitation7Days, SafranDateFirstSnowfall
+from extreme_data.meteo_france_data.scm_models_data.safran.safran_max_snowf import SafranSnowfall2020, \
+    SafranSnowfall2019
 from extreme_data.meteo_france_data.scm_models_data.visualization.study_visualizer import \
     StudyVisualizer
 from projects.exceeding_snow_loads.section_discussion.crocus_study_comparison_with_eurocode import \
@@ -38,6 +42,11 @@ SCM_STUDY_NAME_TO_SCM_STUDY = dict(zip(SCM_STUDIES_NAMES, SCM_STUDIES))
 SCM_STUDY_CLASS_TO_ABBREVIATION = {
     SafranSnowfall: 'SF3',
     SafranSnowfall1Day: 'daily snowfall',
+    SafranSnowfall2020: 'daily snowfall',
+    SafranSnowfall2019: 'daily snowfall',
+    GapBetweenSafranSnowfall2019And2020: 'daily snowfall\n bias = SAFRAN 2020 minus SAFRAN 2019',
+    GapBetweenSafranSnowfall2019AndMySafranSnowfall2019Recentered: 'daily snowfall\n my SAFRAN 2019 recentered minus SAFRAN 2019',
+    GapBetweenSafranSnowfall2019AndMySafranSnowfall2019: 'daily snowfall\n my SAFRAN 2019 minus SAFRAN 2019',
     SafranSnowfall3Days: 'SF3',
     SafranSnowfall5Days: 'SF5',
     SafranSnowfall7Days: 'SF7',
diff --git a/projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py b/projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py
index 05598359..c7654f38 100644
--- a/projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py
+++ b/projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py
@@ -3,21 +3,22 @@ import time
 from typing import List
 
 import matplotlib
+matplotlib.use('Agg')
+
+from extreme_data.meteo_france_data.scm_models_data.safran.safran_max_snowf import SafranSnowfall2019
+from projects.altitude_spatial_model.altitudes_fit.plots.plot_histogram_altitude_studies import \
+    plot_shoe_plot_changes_against_altitude, plot_histogram_all_trends_against_altitudes, \
+    plot_shoe_plot_ratio_interval_size_against_altitude
+
+
 
 from extreme_fit.model.result_from_model_fit.result_from_extremes.abstract_extract_eurocode_return_level import \
     AbstractExtractEurocodeReturnLevel
 
-matplotlib.use('Agg')
-
 import matplotlib as mpl
 
-from extreme_data.meteo_france_data.scm_models_data.crocus.crocus import CrocusSnowLoadTotal, CrocusSnowLoad3Days
 from extreme_fit.model.utils import set_seed_for_test
 from projects.altitude_spatial_model.altitudes_fit.plots.plot_coherence_curves import plot_coherence_curves
-from projects.altitude_spatial_model.altitudes_fit.plots.plot_histogram_altitude_studies import \
-    plot_histogram_all_models_against_altitudes, plot_histogram_all_trends_against_altitudes, \
-    plot_shoe_plot_changes_against_altitude, plot_shoe_plot_changes_against_altitude_for_maxima_and_total, \
-    plot_shoe_plot_ratio_interval_size_against_altitude
 
 mpl.rcParams['text.usetex'] = True
 mpl.rcParams['text.latex.preamble'] = [r'\usepackage{amsmath}']
@@ -33,15 +34,17 @@ from extreme_data.meteo_france_data.scm_models_data.utils import Season
 
 
 def main():
-    study_classes = [SafranSnowfall1Day, SafranSnowfall3Days, SafranSnowfall5Days, SafranSnowfall7Days][:1]
+    study_classes = [SafranSnowfall2019,
+                     SafranSnowfall1Day, SafranSnowfall3Days,
+                     SafranSnowfall5Days, SafranSnowfall7Days][:1]
     seasons = [Season.annual, Season.winter, Season.spring, Season.automn][:1]
 
     set_seed_for_test()
 
-    fast = True
+    fast = False
     if fast is None:
-        massif_names = None
-        altitudes_list = altitudes_for_groups[2:3]
+        massif_names = ['Vanoise']
+        altitudes_list = altitudes_for_groups[:]
     elif fast:
         AbstractExtractEurocodeReturnLevel.NB_BOOTSTRAP = 10
         massif_names = ['Vanoise', 'Haute-Maurienne', 'Vercors'][2:]
@@ -56,7 +59,6 @@ def main():
     duration = str(datetime.timedelta(seconds=end - start))
     print('Total duration', duration)
 
-
 def main_loop(altitudes_list, massif_names, seasons, study_classes):
     assert isinstance(altitudes_list, List)
     assert isinstance(altitudes_list[0], List)
@@ -73,23 +75,23 @@ def main_loop(altitudes_list, massif_names, seasons, study_classes):
 
 
 def plot_visualizers(massif_names, visualizer_list):
-    # plot_histogram_all_trends_against_altitudes(massif_names, visualizer_list)
-    # plot_shoe_plot_ratio_interval_size_against_altitude(massif_names, visualizer_list)
-    # for relative in [True, False]:
-    #     plot_shoe_plot_changes_against_altitude(massif_names, visualizer_list, relative=relative)
-    # plot_coherence_curves(massif_names, visualizer_list)
+    plot_histogram_all_trends_against_altitudes(massif_names, visualizer_list)
+    plot_shoe_plot_ratio_interval_size_against_altitude(massif_names, visualizer_list)
+    for relative in [True, False]:
+        plot_shoe_plot_changes_against_altitude(massif_names, visualizer_list, relative=relative)
+    plot_coherence_curves(massif_names, visualizer_list)
     # plot_coherence_curves(['Vanoise'], visualizer_list)
     pass
 
+
 def plot_visualizer(massif_names, visualizer):
     # Plot time series
-    # visualizer.studies.plot_maxima_time_series(massif_names)
+    visualizer.studies.plot_maxima_time_series(massif_names)
     # visualizer.studies.plot_maxima_time_series(['Vanoise'])
 
     # Plot the results for the model that minimizes the individual aic
     plot_individual_aic(visualizer)
 
-
     # Plot the results for the model that minimizes the total aic
     # plot_total_aic(model_classes, visualizer)
     pass
diff --git a/projects/altitude_spatial_model/altitudes_fit/main_gap_altitudes_studies.py b/projects/altitude_spatial_model/altitudes_fit/main_gap_altitudes_studies.py
new file mode 100644
index 00000000..9801d985
--- /dev/null
+++ b/projects/altitude_spatial_model/altitudes_fit/main_gap_altitudes_studies.py
@@ -0,0 +1,61 @@
+import datetime
+import time
+from typing import List
+
+import matplotlib
+
+from extreme_data.meteo_france_data.scm_models_data.safran.gap_between_study import GapBetweenSafranSnowfall2019And2020, \
+    GapBetweenSafranSnowfall2019AndMySafranSnowfall2019Recentered
+from projects.altitude_spatial_model.altitudes_fit.altitudes_studies import AltitudesStudies
+
+matplotlib.use('Agg')
+
+from extreme_data.meteo_france_data.scm_models_data.safran.safran_max_snowf import SafranSnowfall2020
+
+import matplotlib as mpl
+
+from extreme_fit.model.utils import set_seed_for_test
+
+mpl.rcParams['text.usetex'] = True
+mpl.rcParams['text.latex.preamble'] = [r'\usepackage{amsmath}']
+
+from projects.altitude_spatial_model.altitudes_fit.one_fold_analysis.altitude_group import altitudes_for_groups
+
+from extreme_data.meteo_france_data.scm_models_data.safran.safran import SafranSnowfall1Day, SafranSnowfall3Days, \
+    SafranSnowfall5Days, SafranSnowfall7Days
+from extreme_data.meteo_france_data.scm_models_data.utils import Season
+
+
+def main():
+    study_classes = [GapBetweenSafranSnowfall2019AndMySafranSnowfall2019Recentered,
+                     GapBetweenSafranSnowfall2019And2020, SafranSnowfall2020, SafranSnowfall1Day, SafranSnowfall3Days,
+                     SafranSnowfall5Days, SafranSnowfall7Days][:1]
+    seasons = [Season.annual, Season.winter, Season.spring, Season.automn][:1]
+
+    set_seed_for_test()
+
+    fast = False
+    if fast:
+        altitudes_list = altitudes_for_groups[:1]
+    else:
+        altitudes_list = altitudes_for_groups[1:]
+
+    start = time.time()
+    main_loop_time_series_only(altitudes_list, seasons, study_classes)
+    end = time.time()
+    duration = str(datetime.timedelta(seconds=end - start))
+    print('Total duration', duration)
+
+
+def main_loop_time_series_only(altitudes_list, seasons, study_classes):
+    assert isinstance(altitudes_list, List)
+    assert isinstance(altitudes_list[0], List)
+    for season in seasons:
+        for study_class in study_classes:
+            print('Inner loop', season, study_class)
+            for altitudes in altitudes_list:
+                studies = AltitudesStudies(study_class, altitudes, season=season)
+                studies.plot_maxima_time_series()
+
+if __name__ == '__main__':
+    main()
diff --git a/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/plot_total_aic.py b/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/plot_total_aic.py
index ad0fa2a4..45120d24 100644
--- a/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/plot_total_aic.py
+++ b/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/plot_total_aic.py
@@ -15,8 +15,8 @@ from projects.exceeding_snow_loads.utils import dpi_paper1_figure
 
 def plots(visualizer: AltitudesStudiesVisualizerForNonStationaryModels):
     visualizer.plot_shape_map()
-    # visualizer.plot_moments()
-    visualizer.plot_qqplots()
+    visualizer.plot_moments()
+    # visualizer.plot_qqplots()
     # for std in [True, False]:
     #     visualizer.studies.plot_mean_maxima_against_altitude(std=std)
 
diff --git a/test/test_extreme_data/test_meteo_france_data/test_SCM_study.py b/test/test_extreme_data/test_meteo_france_data/test_SCM_study.py
index 8351ee61..6ee4eb2d 100644
--- a/test/test_extreme_data/test_meteo_france_data/test_SCM_study.py
+++ b/test/test_extreme_data/test_meteo_france_data/test_SCM_study.py
@@ -10,7 +10,7 @@ from extreme_data.meteo_france_data.scm_models_data.crocus.crocus import CrocusS
 from extreme_data.meteo_france_data.scm_models_data.safran.cumulated_study import NB_DAYS
 from extreme_data.meteo_france_data.scm_models_data.safran.safran import SafranSnowfall, SafranTemperature, \
     SafranPrecipitation, SafranSnowfall3Days, SafranRainfall3Days, SafranNormalizedPreciptationRateOnWetDays
-from extreme_data.meteo_france_data.scm_models_data.safran.safran_2020 import SnowfallSafran2020
+from extreme_data.meteo_france_data.scm_models_data.safran.safran_max_snowf import SafranSnowfall2020
 from extreme_data.meteo_france_data.scm_models_data.utils import Season
 from extreme_data.meteo_france_data.scm_models_data.visualization.main_study_visualizer import \
     study_iterator_global, SCM_STUDIES, ALL_ALTITUDES, SCM_STUDY_CLASS_TO_ABBREVIATION
@@ -166,7 +166,7 @@ class TestSafranTemperature(TestSCMStudy):
 class TestSafran2020(TestSCMStudy):
 
     def test_safran_2020_loader(self):
-        study = SnowfallSafran2020(altitude=1800)
+        study = SafranSnowfall2020(altitude=1800)
         annual_maxima = study.year_to_annual_maxima[1959]
         annual_maxima = study.year_to_annual_maxima[2020]
         annual_maxima = study.massif_name_to_annual_maxima['Vercors']
-- 
GitLab