[contrasting] harmonize coherence plot & preliminary plot. add verification...

[contrasting] harmonize coherence plot & preliminary plot. add verification script for return levels.

projects/altitude_spatial_model/altitudes_fit/plots/plot_coherence_curves.py

@@ -40,14 +40,19 @@ def plot_coherence_curves(massif_names, visualizer_list: List[AltitudesStudiesVi
                                                                                     is_altitudinal,
                                                                                     year)
                 label = labels[i]
+                label = label.capitalize()
 
                 # Set labels
                 fontsize_label = 15
+
+                for a in [ax, ax2, ax3]:
+                    a.tick_params(labelsize=fontsize_label)
+
                 if elevation_as_xaxis:
                     ax3.set_xlabel('Elevation (m)', fontsize=fontsize_label)
                     ax2.set_ylabel(label, fontsize=fontsize_label)
                 else:
-                    ax2.set_ylabel('Elevation(m)', fontsize=fontsize_label)
+                    ax2.set_ylabel('Elevation (m)', fontsize=fontsize_label)
                     if i == 3:
                         ax.set_xlabel(label, fontsize=fontsize_label)
                     else:

projects/altitude_spatial_model/preliminary_analysis.py

@@ -45,11 +45,13 @@ class PointwiseGevStudyVisualizer(AltitudesStudies):
             print(param_name, np.mean([c for c in massif_name_to_linear_coef.values()]))
 
             # Display x label
-            elevation_ticks = [1000 * i for i in range(1, 4)]
+            elevation_ticks = [500 * i for i in range(1, 8)]
             if elevation_as_xaxis:
                 ax.set_xticks(elevation_ticks)
             else:
                 ax.set_yticks(elevation_ticks)
+                if j == 2:
+                    ax.set_xlim(right=0.45)
             fontsize_label = 15
             ax.tick_params(labelsize=fontsize_label)
 

projects/altitude_spatial_model/verification_return_level.py

+import numpy as np
+
+from extreme_data.meteo_france_data.scm_models_data.safran.safran import SafranSnowfall1Day
+
+
+def compute_return_level():
+    count_all = 0
+    count_exceedance = 0
+    diff_list = []
+    for altitude in [600, 900, 1200, 1500, 1800, 2100, 2400, 2700, 3000, 3300, 3600][:]:
+        study = SafranSnowfall1Day(altitude=altitude, year_max=2008)
+        for massif_name in study.study_massif_names:
+            if massif_name in study.massif_name_to_stationary_gev_params:
+
+                gev_params = study.massif_name_to_stationary_gev_params[massif_name]
+                return_level100 = gev_params.return_level(return_period=100)
+                annual_maxima = study.massif_name_to_annual_maxima[massif_name]
+                annual_maxima = sorted(annual_maxima)
+                max_annual_maxima = max(annual_maxima)
+
+                count_all += 1
+                relative_diff =100 *  (return_level100 - max_annual_maxima) / max_annual_maxima
+                diff_list.append(relative_diff)
+
+                if return_level100 < max_annual_maxima:
+                    count_exceedance += 1
+                    print(altitude, massif_name, annual_maxima[-2:], return_level100, gev_params.shape)
+    percent_of_exceedance = 100 * count_exceedance / count_all
+    print(percent_of_exceedance)
+    mean_relatif_diff = np.mean(diff_list)
+    print(mean_relatif_diff)
+
+if __name__ == '__main__':
+    compute_return_level()
\ No newline at end of file