diff --git a/extreme_data/meteo_france_data/adamont_data/abstract_adamont_study.py b/extreme_data/meteo_france_data/adamont_data/abstract_adamont_study.py index 2cd805a76af047e60b38eeccc417131dcad2851f..63fc61f977c6428c27a7b35fbc7e65c18bfedda3 100644 --- a/extreme_data/meteo_france_data/adamont_data/abstract_adamont_study.py +++ b/extreme_data/meteo_france_data/adamont_data/abstract_adamont_study.py @@ -21,6 +21,8 @@ from cached_property import cached_property from extreme_data.meteo_france_data.scm_models_data.abstract_study import AbstractStudy from extreme_data.meteo_france_data.scm_models_data.utils import Season, FrenchRegion, french_region_to_str +class WrongYearMinOrYearMax(Exception): + pass class AbstractAdamontStudy(AbstractStudy): YEAR_MIN = 1950 @@ -33,8 +35,13 @@ class AbstractAdamontStudy(AbstractStudy): scenario=AdamontScenario.histo, gcm_rcm_couple=('CNRM-CM5', 'ALADIN53')): # Load the default year_min & year_max for the scenario if not specified year_min_scenario, year_max_scenario = get_year_min_and_year_max_from_scenario(scenario, gcm_rcm_couple) + # Raise exception if year_min is None: year_min = year_min_scenario + else: + if year_min < year_min_scenario: + raise WrongYearMinOrYearMax('year min is {} and should be larger than {}'.format(year_min, year_min_scenario)) + if year_max is None: year_max = year_max_scenario super().__init__(variable_class=variable_class, altitude=altitude, year_min=year_min, year_max=year_max, 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 152971ced77200f9bcbbf2f6e38ca5304c6b0013..8ddf624c61ed13eeac8a095b8bd5007c568eb112 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 @@ -32,6 +32,9 @@ SLEurocode = 'SL from max HS with ' + eurocode_snow_density SCM_STUDIES = [SafranSnowfall, CrocusSweTotal, CrocusDepth, CrocusSwe3Days] SCM_STUDIES_NAMES = [get_display_name_from_object_type(k) for k in SCM_STUDIES] SCM_STUDY_NAME_TO_SCM_STUDY = dict(zip(SCM_STUDIES_NAMES, SCM_STUDIES)) + + +# I keep the scm study separated from the adamont study (for the tests) SCM_STUDY_CLASS_TO_ABBREVIATION = { SafranSnowfall: 'SF3', SafranSnowfall1Day: 'daily snowfall', @@ -62,6 +65,7 @@ SCM_STUDY_CLASS_TO_ABBREVIATION = { ADAMONT_STUDY_CLASS_TO_ABBREVIATION = { AdamontSnowfall: 'daily snowfall', } +STUDY_CLASS_TO_ABBREVIATION = {**ADAMONT_STUDY_CLASS_TO_ABBREVIATION, **SCM_STUDY_CLASS_TO_ABBREVIATION} altitude_massif_name_and_study_class_for_poster = [ (900, 'Chartreuse', CrocusSweTotal), diff --git a/projects/altitude_spatial_model/altitudes_fit/altitudes_studies.py b/projects/altitude_spatial_model/altitudes_fit/altitudes_studies.py index 524f58aa6422afd2a71b0440708843509467b793..2cdcdba4b63197749dbadfe3311d81f079e374f3 100644 --- a/projects/altitude_spatial_model/altitudes_fit/altitudes_studies.py +++ b/projects/altitude_spatial_model/altitudes_fit/altitudes_studies.py @@ -6,7 +6,7 @@ from cached_property import cached_property from extreme_data.meteo_france_data.scm_models_data.abstract_study import AbstractStudy from extreme_data.meteo_france_data.scm_models_data.visualization.main_study_visualizer import \ - SCM_STUDY_CLASS_TO_ABBREVIATION + SCM_STUDY_CLASS_TO_ABBREVIATION, STUDY_CLASS_TO_ABBREVIATION from extreme_data.meteo_france_data.scm_models_data.visualization.plot_utils import plot_against_altitude from extreme_data.meteo_france_data.scm_models_data.visualization.study_visualizer import StudyVisualizer from spatio_temporal_dataset.coordinates.abstract_coordinates import AbstractCoordinates @@ -134,7 +134,8 @@ class AltitudesStudies(object): y = study.massif_name_to_annual_maxima[massif_name] label = '{} m'.format(altitude) ax.plot(x, y, linewidth=2, label=label) - ax.xaxis.set_ticks(x[11::20]) + ax.xaxis.set_ticks([e for e in x if e % 10 == 0][::2]) + ax.set_xlim((x[0], x[-1])) # Plot for the paper 2 if massif_name == "Vanoise": @@ -144,10 +145,10 @@ class AltitudesStudies(object): ax.tick_params(axis='both', which='major', labelsize=20) handles, labels = ax.get_legend_handles_labels() ax.legend(handles[::-1], labels[::-1], prop={'size': 20}) - plot_name = 'Annual maxima of {} in {}'.format(SCM_STUDY_CLASS_TO_ABBREVIATION[self.study_class], + plot_name = 'Annual maxima of {} in {}'.format(STUDY_CLASS_TO_ABBREVIATION[self.study_class], massif_name.replace('_', ' ')) - # ax.set_ylabel('{} ({})'.format(plot_name, self.study.variable_unit), fontsize=15) - # ax.set_xlabel('years', fontsize=15) + ax.set_ylabel('{} ({})'.format(plot_name, self.study.variable_unit), fontsize=15) + ax.set_xlabel('years', fontsize=15) plot_name = 'time series/' + plot_name self.show_or_save_to_file(plot_name=plot_name, show=show, no_title=True, tight_layout=True) ax.clear() diff --git a/projects/altitude_spatial_model/altitudes_fit/utils_altitude_studies_visualizer.py b/projects/altitude_spatial_model/altitudes_fit/utils_altitude_studies_visualizer.py index 343de5879e157f1c3124c477e9b2281b15aefa10..3eaf91fece20944a2fda1dd21d8fb637bbc8eb65 100644 --- a/projects/altitude_spatial_model/altitudes_fit/utils_altitude_studies_visualizer.py +++ b/projects/altitude_spatial_model/altitudes_fit/utils_altitude_studies_visualizer.py @@ -5,24 +5,17 @@ from projects.altitude_spatial_model.altitudes_fit.one_fold_analysis.altitudes_s AltitudesStudiesVisualizerForNonStationaryModels -def load_visualizer_list(season, study_class, altitudes_list, massif_names): +def load_visualizer_list(season, study_class, altitudes_list, massif_names, **kwargs_study): model_classes = ALTITUDINAL_GEV_MODELS_BASED_ON_POINTWISE_ANALYSIS visualizer_list = [] # Load all studies for altitudes in altitudes_list: - # if issubclass(study_class, SimulationStudy): - # for ensemble_idx in list(range(14))[:1]: - # studies = AltitudesStudies(study_class, altitudes, season=season, - # ensemble_idx=ensemble_idx) - # plot_studies(massif_names, season, studies, study_class) - # else: - studies = AltitudesStudies(study_class, altitudes, season=season) + studies = AltitudesStudies(study_class, altitudes, season=season, **kwargs_study) visualizer = AltitudesStudiesVisualizerForNonStationaryModels(studies=studies, model_classes=model_classes, massif_names=massif_names, show=False, temporal_covariate_for_fit=None, - # temporal_covariate_for_fit=MeanAlpsTemperatureCovariate, ) visualizer_list.append(visualizer) compute_and_assign_max_abs(visualizer_list) diff --git a/projects/projected_snowfall/trends/main_projected_trends_separately.py b/projects/projected_snowfall/trends/main_projected_trends_separately.py new file mode 100644 index 0000000000000000000000000000000000000000..8ed7dcc680526d0b8ea5f7de0b4418ad1225088b --- /dev/null +++ b/projects/projected_snowfall/trends/main_projected_trends_separately.py @@ -0,0 +1,95 @@ +import time +from typing import List + +import matplotlib as mpl +import numpy as np + +from extreme_data.meteo_france_data.adamont_data.abstract_adamont_study import WrongYearMinOrYearMax +from extreme_data.meteo_france_data.adamont_data.adamont.adamont_snowfall import AdamontSnowfall +from extreme_data.meteo_france_data.adamont_data.adamont_scenario import AdamontScenario +from extreme_data.meteo_france_data.scm_models_data.crocus.crocus import CrocusSnowLoadTotal, CrocusSnowLoad3Days +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 + +mpl.rcParams['text.usetex'] = True +mpl.rcParams['text.latex.preamble'] = [r'\usepackage{amsmath}'] + +from projects.altitude_spatial_model.altitudes_fit.utils_altitude_studies_visualizer import load_visualizer_list + +from projects.altitude_spatial_model.altitudes_fit.one_fold_analysis.altitude_group import altitudes_for_groups +from projects.altitude_spatial_model.altitudes_fit.one_fold_analysis.plot_total_aic import plot_individual_aic + +from extreme_data.meteo_france_data.scm_models_data.safran.safran import SafranSnowfall1Day, SafranSnowfall3Days, \ + SafranSnowfall5Days, SafranSnowfall7Days, SafranDateFirstSnowfall, SafranPrecipitation1Day, SafranPrecipitation3Days +from extreme_data.meteo_france_data.scm_models_data.utils import Season + + +def main(): + study_class = AdamontSnowfall + scenario = AdamontScenario.rcp85_extended + gcm_rcm_couples = [('CNRM-CM5', 'CCLM4-8-17'), ('CNRM-CM5', 'RCA4')][1:] + + fast = True + if fast is None: + massif_names = None + altitudes_list = altitudes_for_groups[1:2] + elif fast: + massif_names = ['Vanoise', 'Haute-Maurienne', 'Vercors'][:1] + altitudes_list = altitudes_for_groups[1:3] + else: + massif_names = None + altitudes_list = altitudes_for_groups[:] + + # One plot for each massif name + for massif_name in massif_names: + print(massif_name) + main_loop(altitudes_list, [massif_name], gcm_rcm_couples, study_class, scenario) + + +def main_loop(altitudes_list, massif_names, gcm_rcm_couples, study_class, scenario): + assert isinstance(altitudes_list, List) + assert isinstance(altitudes_list[0], List) + altitudes_for_return_levels = [altitudes[-2] for altitudes in altitudes_list] + print(altitudes_for_return_levels) + season = Season.annual + first_year_min, first_year_max = 1951, 2010 + nb_years = first_year_max - first_year_min + 1 + temporal_windows = [(first_year_min + i, first_year_max + i) for i in [30 * j for j in range(4)]] + all_changes_in_return_levels = [] + for gcm_rcm_couple in gcm_rcm_couples: + print('Inner', gcm_rcm_couple) + changes_in_return_levels = [] + for temporal_window in temporal_windows: + year_min, year_max = temporal_window + try: + visualizer_list = load_visualizer_list(season, study_class, altitudes_list, massif_names, + scenario=scenario, gcm_rcm_couple=gcm_rcm_couple, + year_min=year_min, year_max=year_max) + + for visualizer in visualizer_list: + visualizer.studies.plot_maxima_time_series(massif_names) + massif_name = massif_names[0] + changes_for_temporal_window = [ + v.massif_name_to_one_fold_fit[massif_name].changes_of_moment(altitudes=[a], + year=year_max, + nb_years=nb_years, + order=None + )[0] + for a, v in zip(altitudes_for_return_levels, visualizer_list)] + + except WrongYearMinOrYearMax: + changes_for_temporal_window = [np.nan for _ in altitudes_for_return_levels] + + print(changes_for_temporal_window) + changes_in_return_levels.append(changes_for_temporal_window) + changes_in_return_levels = np.array(changes_in_return_levels) + all_changes_in_return_levels.append(changes_in_return_levels) + all_changes_in_return_levels = np.array(all_changes_in_return_levels) + + return all_changes_in_return_levels, temporal_windows, altitudes_for_return_levels, nb_years + + +if __name__ == '__main__': + main()