From 94e7a5e517b73ebfdd9562d9b669ca2ca2bf293d Mon Sep 17 00:00:00 2001
From: Le Roux Erwan <erwan.le-roux@irstea.fr>
Date: Tue, 3 Dec 2019 17:31:18 +0100
Subject: [PATCH] [PAPER 1] fix some display issues. add dpi to plots. improve
 the quality of plots.

---
 .../main_gelman_convergence_test.py           |  5 +--
 .../data/main_example_swe_total_plot.py       |  3 +-
 .../main_comparison_with_eurocode.py          |  5 +--
 .../paper_past_snow_loads/paper_main_utils.py | 14 ++++++++
 .../paper_past_snow_loads/paper_utils.py      | 15 ++-------
 .../main_result_trends_and_return_levels.py   | 32 +++++++++++++------
 .../plot_uncertainty_curves.py                | 31 +++++++++++-------
 .../plot_uncertainty_histogram.py             |  3 +-
 ...dy_visualizer_for_non_stationary_trends.py |  4 ++-
 9 files changed, 71 insertions(+), 41 deletions(-)
 create mode 100644 experiment/paper_past_snow_loads/paper_main_utils.py

diff --git a/experiment/paper_past_snow_loads/check_mcmc_convergence_for_return_levels/main_gelman_convergence_test.py b/experiment/paper_past_snow_loads/check_mcmc_convergence_for_return_levels/main_gelman_convergence_test.py
index b6815268..53466706 100644
--- a/experiment/paper_past_snow_loads/check_mcmc_convergence_for_return_levels/main_gelman_convergence_test.py
+++ b/experiment/paper_past_snow_loads/check_mcmc_convergence_for_return_levels/main_gelman_convergence_test.py
@@ -1,6 +1,7 @@
 import pandas as pd
-from experiment.paper_past_snow_loads.paper_utils import paper_altitudes, paper_study_classes, \
-    load_altitude_to_visualizer
+
+from experiment.paper_past_snow_loads.paper_main_utils import load_altitude_to_visualizer
+from experiment.paper_past_snow_loads.paper_utils import paper_altitudes, paper_study_classes
 from extreme_fit.model.margin_model.linear_margin_model.temporal_linear_margin_models import StationaryTemporalModel
 from root_utils import get_display_name_from_object_type
 
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 f31f388f..0e3c1a38 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
@@ -5,6 +5,7 @@ from experiment.meteo_france_data.scm_models_data.visualization.study_visualizat
     study_iterator_global, SCM_STUDY_CLASS_TO_ABBREVIATION
 from experiment.meteo_france_data.scm_models_data.visualization.study_visualization.study_visualizer import \
     StudyVisualizer
+from experiment.paper_past_snow_loads.paper_utils import dpi_paper1_figure
 
 
 def max_graph_annual_maxima_poster():
@@ -34,7 +35,7 @@ def max_graph_annual_maxima_poster():
             tight_pad = {'h_pad': 0.2}
             study_visualizer.visualize_max_graphs_poster(massif_name, altitude, snow_abbreviation, color, label,
                                                          last_plot, ax, tight_pad=tight_pad,
-                                                         dpi=1000)
+                                                         dpi=dpi_paper1_figure)
 
 
 if __name__ == '__main__':
diff --git a/experiment/paper_past_snow_loads/discussion_data_comparison_with_eurocode/main_comparison_with_eurocode.py b/experiment/paper_past_snow_loads/discussion_data_comparison_with_eurocode/main_comparison_with_eurocode.py
index 623d750c..b3d04911 100644
--- a/experiment/paper_past_snow_loads/discussion_data_comparison_with_eurocode/main_comparison_with_eurocode.py
+++ b/experiment/paper_past_snow_loads/discussion_data_comparison_with_eurocode/main_comparison_with_eurocode.py
@@ -10,6 +10,7 @@ from experiment.paper_past_snow_loads.discussion_data_comparison_with_eurocode.c
     CrocusDifferenceSnowLoad, \
     CrocusSnowDensityAtMaxofSwe, CrocusDifferenceSnowLoadRescaledAndEurocodeToSeeSynchronization, \
     CrocusSnowDepthAtMaxofSwe, CrocusSnowDepthDifference
+from experiment.paper_past_snow_loads.paper_utils import dpi_paper1_figure
 
 
 def max_graph_annual_maxima_comparison():
@@ -52,7 +53,7 @@ def max_graph_annual_maxima_comparison():
                 last_plot = altitude == 2700
                 if last_plot:
                     if study_class == CrocusSnowDensityAtMaxofSwe:
-                        label = '{} Eurocode'.format(snow_density_str)
+                        label = '{} for French standards'.format(snow_density_str)
                         snow_density_eurocode = [150 for _ in study.ordered_years]
                         ax.plot(study.ordered_years, snow_density_eurocode, color='k', label=label)
                     ax.legend()
@@ -61,7 +62,7 @@ def max_graph_annual_maxima_comparison():
                     ax.set_xlim([1957, 2018])
                     ax.yaxis.set_ticks(yticks)
                     study_visualizer.show_or_save_to_file(no_title=True, tight_layout=True,
-                                                          tight_pad=tight_pad, dpi=1000)
+                                                          tight_pad=tight_pad, dpi=dpi_paper1_figure)
                     ax.clear()
 
 
diff --git a/experiment/paper_past_snow_loads/paper_main_utils.py b/experiment/paper_past_snow_loads/paper_main_utils.py
new file mode 100644
index 00000000..94493f2f
--- /dev/null
+++ b/experiment/paper_past_snow_loads/paper_main_utils.py
@@ -0,0 +1,14 @@
+from collections import OrderedDict
+
+from experiment.paper_past_snow_loads.study_visualizer_for_non_stationary_trends import \
+    StudyVisualizerForNonStationaryTrends
+
+
+def load_altitude_to_visualizer(altitudes, massif_names, non_stationary_uncertainty, study_class, uncertainty_methods):
+    altitude_to_visualizer = OrderedDict()
+    for altitude in altitudes:
+        altitude_to_visualizer[altitude] = StudyVisualizerForNonStationaryTrends(
+            study=study_class(altitude=altitude), multiprocessing=True, save_to_file=True,
+            uncertainty_massif_names=massif_names, uncertainty_methods=uncertainty_methods,
+            non_stationary_contexts=non_stationary_uncertainty)
+    return altitude_to_visualizer
diff --git a/experiment/paper_past_snow_loads/paper_utils.py b/experiment/paper_past_snow_loads/paper_utils.py
index 7af16183..bfbeae19 100644
--- a/experiment/paper_past_snow_loads/paper_utils.py
+++ b/experiment/paper_past_snow_loads/paper_utils.py
@@ -1,19 +1,8 @@
-from collections import OrderedDict
-
 from experiment.meteo_france_data.scm_models_data.crocus.crocus import CrocusSnowLoadTotal, CrocusSnowLoadEurocode, \
     CrocusSnowLoad3Days
-from experiment.paper_past_snow_loads.study_visualizer_for_non_stationary_trends import \
-    StudyVisualizerForNonStationaryTrends
 
 paper_altitudes = [300, 600, 900, 1200, 1500, 1800, 2100, 2400, 2700]
 paper_study_classes = [CrocusSnowLoadTotal, CrocusSnowLoadEurocode, CrocusSnowLoad3Days]
+# dpi_paper1_figure = 700
+dpi_paper1_figure = None
 
-
-def load_altitude_to_visualizer(altitudes, massif_names, non_stationary_uncertainty, study_class, uncertainty_methods):
-    altitude_to_visualizer = OrderedDict()
-    for altitude in altitudes:
-        altitude_to_visualizer[altitude] = StudyVisualizerForNonStationaryTrends(
-            study=study_class(altitude=altitude), multiprocessing=True, save_to_file=True,
-            uncertainty_massif_names=massif_names, uncertainty_methods=uncertainty_methods,
-            non_stationary_contexts=non_stationary_uncertainty)
-    return altitude_to_visualizer
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 4e611d73..f0b21ab0 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
@@ -1,21 +1,19 @@
-from time import sleep
+from multiprocessing.pool import Pool
 
 import matplotlib as mpl
-import matplotlib.pyplot as plt
 
-from experiment.meteo_france_data.scm_models_data.visualization.study_visualization.main_study_visualizer import \
-    ALL_ALTITUDES_WITHOUT_NAN
-from experiment.paper_past_snow_loads.paper_utils import paper_study_classes, paper_altitudes, \
-    load_altitude_to_visualizer
+from experiment.meteo_france_data.scm_models_data.crocus.crocus import CrocusSnowLoadTotal, CrocusSnowLoadEurocode
+from experiment.paper_past_snow_loads.paper_main_utils import load_altitude_to_visualizer
+from experiment.paper_past_snow_loads.paper_utils import paper_study_classes, paper_altitudes
 from experiment.paper_past_snow_loads.result_trends_and_return_levels.plot_uncertainty_curves import \
     plot_uncertainty_massifs
-from experiment.meteo_france_data.scm_models_data.crocus.crocus import CrocusSnowLoadTotal, CrocusSnowLoadEurocode
 from experiment.paper_past_snow_loads.result_trends_and_return_levels.plot_uncertainty_histogram import \
     plot_uncertainty_histogram
 from experiment.paper_past_snow_loads.study_visualizer_for_non_stationary_trends import \
     StudyVisualizerForNonStationaryTrends
 from extreme_fit.model.result_from_model_fit.result_from_extremes.confidence_interval_method import \
     ConfidenceIntervalMethodFromExtremes
+from root_utils import NB_CORES
 
 mpl.rcParams['text.usetex'] = True
 mpl.rcParams['text.latex.preamble'] = [r'\usepackage{amsmath}']
@@ -28,10 +26,14 @@ def minor_result(altitude):
     visualizer.plot_trends()
     # plt.show()
 
+def compute_minimized_aic(visualizer):
+    _ = visualizer.massif_name_to_minimized_aic_non_stationary_trend_test
+    return True
 
 def intermediate_result(altitudes, massif_names=None,
                         non_stationary_uncertainty=None, uncertainty_methods=None,
-                        study_class=CrocusSnowLoadTotal):
+                        study_class=CrocusSnowLoadTotal,
+                        multiprocessing=False):
     """
     Plot all the trends for all altitudes
     And enable to plot uncertainty plot for some specific massif_names, uncertainty methods to be fast
@@ -46,9 +48,18 @@ def intermediate_result(altitudes, massif_names=None,
     altitude_to_visualizer = load_altitude_to_visualizer(altitudes, massif_names, non_stationary_uncertainty,
                                                          study_class, uncertainty_methods)
     # Plot trends
-    max_abs_tdrl = max([visualizer.max_abs_tdrl for visualizer in altitude_to_visualizer.values()])
+    visualizers = list(altitude_to_visualizer.values())
+    if multiprocessing:
+        with Pool(NB_CORES) as p:
+            _ = p.map(compute_minimized_aic, visualizers)
+    else:
+        for visualizer in visualizers:
+            _ = compute_minimized_aic(visualizer)
+    # Compute common max value for the colorbar
+    max_abs_tdrl = max([visualizer.max_abs_tdrl for visualizer in visualizers])
     for visualizer in altitude_to_visualizer.values():
         visualizer.plot_trends(max_abs_tdrl)
+
     # Plot graph
     plot_uncertainty_massifs(altitude_to_visualizer)
     # Plot histogram
@@ -72,7 +83,8 @@ if __name__ == '__main__':
     intermediate_result(altitudes=paper_altitudes, massif_names=['Maurienne'],
                         uncertainty_methods=[ConfidenceIntervalMethodFromExtremes.my_bayes,
                            ConfidenceIntervalMethodFromExtremes.ci_mle][:],
-                        non_stationary_uncertainty=[False, True][:])
+                        non_stationary_uncertainty=[False, True][:],
+                        multiprocessing=False)
     # intermediate_result(altitudes=[900, 1200], massif_names=None)
     # intermediate_result(ALL_ALTITUDES_WITHOUT_NAN)
     # intermediate_result(paper_altitudes)
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 0e34087f..9ba39fd5 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
@@ -7,6 +7,7 @@ from experiment.eurocode_data.utils import EUROCODE_RETURN_LEVEL_STR, EUROCODE_A
 from experiment.meteo_france_data.scm_models_data.abstract_study import AbstractStudy
 from experiment.meteo_france_data.scm_models_data.visualization.study_visualization.main_study_visualizer import \
     SCM_STUDY_CLASS_TO_ABBREVIATION
+from experiment.paper_past_snow_loads.paper_utils import dpi_paper1_figure
 from experiment.paper_past_snow_loads.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 \
@@ -61,7 +62,7 @@ def plot_single_uncertainty_massif(altitude_to_visualizer: Dict[int, StudyVisual
         massif_names_str = massif_name
         model_names_str = 'NonStationarity={}'.format(non_stationary_context)
         visualizer.plot_name = model_names_str + '_' + massif_names_str
-        visualizer.show_or_save_to_file(no_title=True, dpi=1000)
+        visualizer.show_or_save_to_file(no_title=True, dpi=dpi_paper1_figure)
         plt.close()
 
 
@@ -83,6 +84,8 @@ def plot_single_uncertainty_massif_and_non_stationary_context(ax, massif_name, n
     altitudes = list(altitude_to_visualizer.keys())
     visualizer = list(altitude_to_visualizer.values())[0]
     alpha = 0.2
+    legend_size = 20
+    fontsize_label = 20
     # Display the EUROCODE return level
     eurocode_region = massif_name_to_eurocode_region[massif_name]()
 
@@ -100,12 +103,23 @@ def plot_single_uncertainty_massif_and_non_stationary_context(ax, massif_name, n
 
         # Plot bars of TDRL only in the non stationary case
         if j == 0 and non_stationary_context:
-            plot_tdrl_bars(altitude_to_visualizer, ax, massif_name, valid_altitudes)
+            plot_tdrl_bars(altitude_to_visualizer, ax, massif_name, valid_altitudes, legend_size, fontsize_label)
 
-    ax.legend(loc=2)
+    ax.legend(loc=2, prop={'size': legend_size})
     # ax.set_ylim([-1, 16])
+    ax.set_xlim([200, 1900])
     if massif_name == 'Maurienne':
         ax.set_ylim([-1, 13])
+    # add_title(ax, eurocode_region, massif_name, non_stationary_context)
+    ax.set_xticks(altitudes)
+    ax.tick_params(labelsize=fontsize_label)
+    ylabel = EUROCODE_RETURN_LEVEL_STR.replace('GSL', SCM_STUDY_CLASS_TO_ABBREVIATION[type(visualizer.study)])
+    ax.set_ylabel(ylabel, fontsize=fontsize_label)
+    ax.set_xlabel('Altitude (m)', fontsize=fontsize_label)
+    ax.grid()
+
+
+def add_title(ax, eurocode_region, massif_name, non_stationary_context):
     massif_name_str = massif_name.replace('_', ' ')
     eurocode_region_str = get_display_name_from_object_type(type(eurocode_region))
     is_non_stationary_model = non_stationary_context if isinstance(non_stationary_context,
@@ -116,14 +130,9 @@ def plot_single_uncertainty_massif_and_non_stationary_context(ax, massif_name, n
         non_stationary_context = 'the stationary model'
     title = '{} massif with {}'.format(massif_name_str, non_stationary_context)
     ax.set_title(title)
-    ax.set_xticks(altitudes)
-    ylabel = EUROCODE_RETURN_LEVEL_STR.replace('GSL', SCM_STUDY_CLASS_TO_ABBREVIATION[type(visualizer.study)])
-    ax.set_ylabel(ylabel)
-    ax.set_xlabel('Altitude (m)')
-    ax.grid()
 
 
-def plot_tdrl_bars(altitude_to_visualizer, ax, massif_name, valid_altitudes):
+def plot_tdrl_bars(altitude_to_visualizer, ax, massif_name, valid_altitudes, legend_size, fontsize):
     visualizers = [v for a, v in altitude_to_visualizer.items() if
                    a in valid_altitudes and massif_name in v.uncertainty_massif_names]
     if len(visualizers) > 0:
@@ -145,8 +154,8 @@ def plot_tdrl_bars(altitude_to_visualizer, ax, massif_name, valid_altitudes):
     ax2 = ax.twinx()
     # ax2.legend(handles=legend_elements, bbox_to_anchor=(0.93, 0.7), loc='upper right')
     # ax2.annotate("Filled symbol = significant trend ", xy=(0.85, 0.5), xycoords='axes fraction', fontsize=7)
-    ax2.legend(handles=legend_elements, loc='upper right')
-    ax2.annotate("Filled symbol = significant trend ", xy=(0.75, 0.97), xycoords='axes fraction', fontsize=10)
+    ax2.legend(handles=legend_elements, loc='upper right', prop={'size': legend_size})
+    ax2.annotate("Filled symbol = significant trend ", xy=(0.5, 0.93), xycoords='axes fraction', fontsize=fontsize)
     ax2.set_yticks([])
 
 
diff --git a/experiment/paper_past_snow_loads/result_trends_and_return_levels/plot_uncertainty_histogram.py b/experiment/paper_past_snow_loads/result_trends_and_return_levels/plot_uncertainty_histogram.py
index 74b0cc00..c1d5f99b 100644
--- a/experiment/paper_past_snow_loads/result_trends_and_return_levels/plot_uncertainty_histogram.py
+++ b/experiment/paper_past_snow_loads/result_trends_and_return_levels/plot_uncertainty_histogram.py
@@ -3,6 +3,7 @@ import matplotlib.pyplot as plt
 import numpy as np
 
 from experiment.eurocode_data.utils import EUROCODE_RETURN_LEVEL_STR, EUROCODE_ALTITUDES
+from experiment.paper_past_snow_loads.paper_utils import dpi_paper1_figure
 from experiment.paper_past_snow_loads.study_visualizer_for_non_stationary_trends import \
     StudyVisualizerForNonStationaryTrends
 from extreme_fit.model.result_from_model_fit.result_from_extremes.confidence_interval_method import ci_method_to_color, \
@@ -45,7 +46,7 @@ def plot_histogram(altitude_to_visualizer, non_stationary_context):
     ax.set_yticks([10 * i for i in range(11)])
     visualizer.plot_name = 'Percentages of exceedance with non_stationary={}'.format(non_stationary_context)
     # visualizer.show = True
-    visualizer.show_or_save_to_file(no_title=True, dpi=1000)
+    visualizer.show_or_save_to_file(no_title=True, dpi=dpi_paper1_figure)
     ax.clear()
 
 
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 35ee348e..da62927e 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
@@ -14,6 +14,7 @@ from experiment.meteo_france_data.scm_models_data.visualization.study_visualizat
     StudyVisualizer
 from experiment.paper_past_snow_loads.check_mcmc_convergence_for_return_levels.gelman_convergence_test import \
     compute_gelman_convergence_value
+from experiment.paper_past_snow_loads.paper_utils import dpi_paper1_figure
 from experiment.trend_analysis.abstract_score import MeanScore
 from experiment.trend_analysis.univariate_test.abstract_gev_trend_test import AbstractGevTrendTest
 from experiment.trend_analysis.univariate_test.gev_trend_test_one_parameter import GevScaleTrendTest, \
@@ -104,6 +105,7 @@ class StudyVisualizerForNonStationaryTrends(StudyVisualizer):
             massif_name_to_trend_test_that_minimized_aic[massif_name] = trend_test_that_minimized_aic
         return massif_name_to_trend_test_that_minimized_aic
 
+
     # Part 1 - Trends
 
     @property
@@ -133,7 +135,7 @@ class StudyVisualizerForNonStationaryTrends(StudyVisualizer):
 
         self.plot_name = 'tdlr_trends'
         self.show_or_save_to_file(add_classic_title=False, tight_layout=True, no_title=True,
-                                  dpi=1000)
+                                  dpi=dpi_paper1_figure)
         plt.close()
 
     @property
-- 
GitLab