diff --git a/experiment/paper_past_snow_loads/check_mle_convergence_for_trends/shape/main_shape_repartition.py b/experiment/paper_past_snow_loads/check_mle_convergence_for_trends/shape/main_shape_repartition.py
index f9cb04c4c6217e3e14e201c975cb4f6aa216751f..ca663e11b455e9eded51cc1525a248f359542f56 100644
--- a/experiment/paper_past_snow_loads/check_mle_convergence_for_trends/shape/main_shape_repartition.py
+++ b/experiment/paper_past_snow_loads/check_mle_convergence_for_trends/shape/main_shape_repartition.py
@@ -28,6 +28,6 @@ if __name__ == '__main__':
     # main_shape_repartition([900], save_to_file=False)
     # main_shape_repartition([900, 1800, 2700])
     # main_shape_repartition([300, 600, 900, 1200, 1500, 1800, 2700])
-    # main_shape_repartition([900], study_visualizer_class=StudyVisualizerForShape, save_to_file=False)
-    main_shape_repartition([300, 600, 900, 1200, 1500, 1800, 2100, 2400, 2700, 3000, 3300, 3600, 3900, 4200],
-                           study_visualizer_class=StudyVisualizerForShape, save_to_file=True)
+    main_shape_repartition([900, 1800, 2700], study_visualizer_class=StudyVisualizerForShape, save_to_file=True)
+    # main_shape_repartition([300, 600, 900, 1200, 1500, 1800, 2100, 2400, 2700, 3000, 3300, 3600, 3900, 4200],
+    #                        study_visualizer_class=StudyVisualizerForShape, save_to_file=True)
diff --git a/experiment/paper_past_snow_loads/data/main_example_swe_total_plot.py b/experiment/paper_past_snow_loads/data/main_example_swe_total_plot.py
index 4d7313fde7241f58673e6544ed1397b5832830f9..d912e4c8b765d6ed549681d27d022da2e89008bb 100644
--- a/experiment/paper_past_snow_loads/data/main_example_swe_total_plot.py
+++ b/experiment/paper_past_snow_loads/data/main_example_swe_total_plot.py
@@ -7,12 +7,6 @@ from experiment.meteo_france_data.scm_models_data.visualization.study_visualizat
     StudyVisualizer
 from experiment.paper_past_snow_loads.paper_utils import dpi_paper1_figure
 
-marker_altitude_massif_name_for_paper1 = [
-    ('magenta', 900, 'Ubaye'),
-    ('darkmagenta', 1800, 'Vercors'),
-    ('mediumpurple', 2700, 'Beaufortain'),
-]
-
 
 def max_graph_annual_maxima_poster():
     """
@@ -22,17 +16,31 @@ def max_graph_annual_maxima_poster():
     """
     save_to_file = True
     study_class = CrocusSnowLoadTotal
+    examples_for_the_paper = True
 
     ax = plt.gca()
-    ax.set_ylim([0, 20])
-    ax.set_yticks(list(range(0, 21, 2)))
-    for color, altitude, massif_name in marker_altitude_massif_name_for_paper1:
+
+    if examples_for_the_paper:
+        ax.set_ylim([0, 20])
+        ax.set_yticks(list(range(0, 21, 2)))
+        marker_altitude_massif_name_for_paper1 = [
+            ('magenta', 900, 'Ubaye'),
+            ('darkmagenta', 1800, 'Vercors'),
+            ('mediumpurple', 2700, 'Beaufortain'),
+        ]
+    else:
+        marker_altitude_massif_name_for_paper1 = [
+            ('yellow', 600, 'Ubaye'),
+            ('purple', 600, 'Parpaillon'),
+        ]
+
+    for color, altitude, massif_name in marker_altitude_massif_name_for_paper1[::-1]:
         for study in study_iterator_global([study_class], altitudes=[altitude]):
             study_visualizer = StudyVisualizer(study, save_to_file=save_to_file,
                                                verbose=True,
                                                multiprocessing=True)
             snow_abbreviation = SCM_STUDY_CLASS_TO_ABBREVIATION[study_class]
-            last_plot = altitude == 2700
+            last_plot = massif_name == "Ubaye"
             label = '{} massif at {}m'.format(massif_name, altitude)
             tight_pad = {'h_pad': 0.2}
             study_visualizer.visualize_max_graphs_poster(massif_name, altitude, snow_abbreviation, color, label,
diff --git a/experiment/paper_past_snow_loads/paper_main_utils.py b/experiment/paper_past_snow_loads/paper_main_utils.py
index 8ba05340dfb5c582d86881a8c5342ddef7d8373a..78b7dd19655a15de38ecf15302bf375439dd49bf 100644
--- a/experiment/paper_past_snow_loads/paper_main_utils.py
+++ b/experiment/paper_past_snow_loads/paper_main_utils.py
@@ -13,7 +13,8 @@ def load_altitude_to_visualizer(altitudes, massif_names, model_subsets_for_uncer
                                 save_to_file=True):
     fit_method = TemporalMarginFitMethod.extremes_fevd_mle
     select_only_acceptable_shape_parameter = True
-    print('Fit method: {}, Select shape parameter: {}'.format(fit_method, select_only_acceptable_shape_parameter))
+    print('Fit method: {}, Select only acceptable shape parameter: {}'
+          .format(fit_method, select_only_acceptable_shape_parameter))
     altitude_to_visualizer = OrderedDict()
     for altitude in altitudes:
         altitude_to_visualizer[altitude] = study_visualizer_class(
diff --git a/experiment/paper_past_snow_loads/paper_utils.py b/experiment/paper_past_snow_loads/paper_utils.py
index 66dd821a91924de8628bd8a922a9d0b7e38ac1f9..ffcee10d22f69ea4fe09c4e985a51f9d3afd17f7 100644
--- a/experiment/paper_past_snow_loads/paper_utils.py
+++ b/experiment/paper_past_snow_loads/paper_utils.py
@@ -2,9 +2,12 @@ from enum import Enum
 
 from experiment.meteo_france_data.scm_models_data.crocus.crocus import CrocusSnowLoadTotal, CrocusSnowLoadEurocode, \
     CrocusSnowLoad3Days
+from experiment.meteo_france_data.scm_models_data.visualization.study_visualization.main_study_visualizer import \
+    ALL_ALTITUDES_WITHOUT_NAN
 from root_utils import get_display_name_from_object_type
 
-paper_altitudes = [300, 600, 900, 1200, 1500, 1800, 2100, 2400, 2700]
+# paper_altitudes = [300, 600, 900, 1200, 1500, 1800, 2100, 2400, 2700]
+paper_altitudes = ALL_ALTITUDES_WITHOUT_NAN
 paper_study_classes = [CrocusSnowLoadTotal, CrocusSnowLoadEurocode, CrocusSnowLoad3Days]
 # dpi_paper1_figure = 700
 dpi_paper1_figure = None
diff --git a/experiment/paper_past_snow_loads/result_trends_and_return_levels/main_result_trends_and_return_levels.py b/experiment/paper_past_snow_loads/result_trends_and_return_levels/main_result_trends_and_return_levels.py
index 2f81491c503d4606ca54e78d8ee5c038fd2cf24c..06d426b589d1a314432f65b399006ffd68c9982f 100644
--- a/experiment/paper_past_snow_loads/result_trends_and_return_levels/main_result_trends_and_return_levels.py
+++ b/experiment/paper_past_snow_loads/result_trends_and_return_levels/main_result_trends_and_return_levels.py
@@ -63,15 +63,17 @@ def intermediate_result(altitudes, massif_names=None,
         for visualizer in visualizers:
             _ = compute_minimized_aic(visualizer)
     # Compute common max value for the colorbar
-    altitudes_for_plot_trend = [900, 1800, 2700]
-    altitudes_for_plot_trend = altitudes
-    visualizers_for_altitudes = [visualizer
-                                 for altitude, visualizer in altitude_to_visualizer.items()
-                                 if altitude in altitudes_for_plot_trend]
-    max_abs_tdrl = max([visualizer.max_abs_change for visualizer in visualizers_for_altitudes])
-    for visualizer in visualizers_for_altitudes:
-        # visualizer.plot_trends(max_abs_tdrl, add_colorbar=visualizer.study.altitude == 2700)
-        visualizer.plot_trends(None, add_colorbar=True)
+    max_abs_tdrl = max([visualizer.max_abs_change for altitude, visualizer in altitude_to_visualizer.items()
+                        if altitude >= 900])
+    for altitude, visualizer in altitude_to_visualizer.items():
+        if 900 <= altitude <= 4200:
+            add_color = (visualizer.study.altitude - 1500) % 900 == 0
+            visualizer.plot_trends(max_abs_tdrl, add_colorbar=add_color)
+            # Plot 2700 also with a colorbar
+            if altitude == 2700:
+                visualizer.plot_trends(max_abs_tdrl, add_colorbar=True)
+        else:
+            visualizer.plot_trends(None, add_colorbar=True)
 
     # Plot graph
     plot_uncertainty_massifs(altitude_to_visualizer)
@@ -84,11 +86,11 @@ def major_result():
                            ConfidenceIntervalMethodFromExtremes.ci_mle][1:]
     massif_names = None
     study_classes = paper_study_classes[:2]
-    model_subsets_for_uncertainty = [ModelSubsetForUncertainty.stationary_gumbel,
-                                     ModelSubsetForUncertainty.stationary_gumbel_and_gev,
-                                     ModelSubsetForUncertainty.non_stationary_gumbel,
-                                     ModelSubsetForUncertainty.non_stationary_gumbel_and_gev]
-    # model_subsets_for_uncertainty = None
+    # model_subsets_for_uncertainty = [ModelSubsetForUncertainty.stationary_gumbel,
+    #                                  ModelSubsetForUncertainty.stationary_gumbel_and_gev,
+    #                                  ModelSubsetForUncertainty.non_stationary_gumbel,
+    #                                  ModelSubsetForUncertainty.non_stationary_gumbel_and_gev]
+    model_subsets_for_uncertainty = None
     # study_classes = [CrocusSnowLoad3Days, CrocusSnowLoad5Days, CrocusSnowLoad7Days][::-1]
     for study_class in study_classes:
         intermediate_result(paper_altitudes, massif_names, model_subsets_for_uncertainty,
@@ -100,7 +102,6 @@ if __name__ == '__main__':
     # intermediate_result(altitudes=[900, 1200], massif_names=['Vercors'],
     #                     uncertainty_methods=[ConfidenceIntervalMethodFromExtremes.my_bayes,
     #                                          ConfidenceIntervalMethodFromExtremes.ci_mle][1:],
-    #                     non_stationary_uncertainty=[False, True][1:],
     #                     multiprocessing=True)
     # intermediate_result(altitudes=[900, 1200], massif_names=['Maurienne'],
     #                     uncertainty_methods=[ConfidenceIntervalMethodFromExtremes.my_bayes,
diff --git a/experiment/paper_past_snow_loads/result_trends_and_return_levels/plot_uncertainty_curves.py b/experiment/paper_past_snow_loads/result_trends_and_return_levels/plot_uncertainty_curves.py
index 3f274fdb691f9dd23811a688f808cc8fa8336907..1948b79daf61f9e1476886fdaac900649e95e02e 100644
--- a/experiment/paper_past_snow_loads/result_trends_and_return_levels/plot_uncertainty_curves.py
+++ b/experiment/paper_past_snow_loads/result_trends_and_return_levels/plot_uncertainty_curves.py
@@ -203,4 +203,7 @@ def plot_valid_return_level_uncertainties(alpha, altitude_to_visualizer, altitud
     confidence_interval_str += '\% confidence interval'
     ax.fill_between(valid_altitudes, lower_bound, upper_bound, color=color, alpha=alpha,
                     label=label_name + confidence_interval_str)
+    # Plot error bars
+    yerr = np.array([[d[1] - d[0], d[2] - d[1]] for d in zip(lower_bound, mean, upper_bound)]).transpose()
+    ax.bar(valid_altitudes, mean,  ecolor='black', capsize=5, yerr=yerr)
     return valid_altitudes
diff --git a/experiment/paper_past_snow_loads/study_visualizer_for_non_stationary_trends.py b/experiment/paper_past_snow_loads/study_visualizer_for_non_stationary_trends.py
index bce1e2fcd907e4e19a68cf85d2487f07b4dd1d6a..fb37f57b56dd0dac12684329ad2dadd2ef56fdb1 100644
--- a/experiment/paper_past_snow_loads/study_visualizer_for_non_stationary_trends.py
+++ b/experiment/paper_past_snow_loads/study_visualizer_for_non_stationary_trends.py
@@ -190,7 +190,7 @@ class StudyVisualizerForNonStationaryTrends(StudyVisualizer):
         ax.get_xaxis().set_visible(True)
         ax.set_xticks([])
         ax.set_xlabel('Altitude = {}m'.format(self.study.altitude), fontsize=15)
-        self.plot_name = 'tdlr_trends'
+        self.plot_name = 'tdlr_trends_w' + 'o' if not add_colorbar else '' + '_colorbar'
         self.show_or_save_to_file(add_classic_title=False, tight_layout=True, no_title=True,
                                   dpi=500)
         plt.close()