You cannot select more than 25 topics
Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
199 lines
9.8 KiB
Python
199 lines
9.8 KiB
Python
# -*- coding: utf-8 -*-
|
|
#####################################----------------------------------
|
|
#Last Updated - March 2018
|
|
#@author: z5025317 Valentin Heimhuber
|
|
#code for creating climate prioritization plots for NARCLIM variables.
|
|
#Inputs: Uses CSV files that contain all 12 NARCLIM model runs time series for 1 grid cell created with: P1_NARCliM_NC_to_CSV_CCRC_SS.py
|
|
#####################################----------------------------------
|
|
#Load packages
|
|
#####################################----------------------------------
|
|
import numpy as np
|
|
import os
|
|
import pandas as pd
|
|
import glob
|
|
import matplotlib
|
|
import matplotlib.pyplot as plt
|
|
from datetime import datetime
|
|
from datetime import timedelta
|
|
from matplotlib.backends.backend_pdf import PdfPages
|
|
matplotlib.style.use('ggplot')
|
|
#
|
|
# Set working direcotry (where postprocessed NARClIM data is located)
|
|
os.chdir('C:/Users/z5025317/WRL_Postdoc/Projects/Paper#1/NARCLIM/')
|
|
#
|
|
#
|
|
#####################################----------------------------------
|
|
#set input parameters
|
|
Base_period_start = '1990-01-01'
|
|
Base_period_end = '2080-01-01' #use last day that's not included in period as < is used for subsetting
|
|
#Clim_var_type = 'tasmean' will create pdf for all variables in folder
|
|
#####################################----------------------------------
|
|
#set directory path for output files
|
|
output_directory = 'Output'
|
|
#output_directory = 'J:/Project wrl2016032/NARCLIM_Raw_Data/Extracted'
|
|
if not os.path.exists(output_directory):
|
|
os.makedirs(output_directory)
|
|
print("output directory folder didn't exist and was generated")
|
|
Clim_Var_CSVs = glob.glob('./Site_CSVs/*')
|
|
#Clim_Var_CSV = glob.glob('./Site_CSVs/' + Clim_var_type + '*' )
|
|
#read CSV file
|
|
for clim_var_csv_path in Clim_Var_CSVs:
|
|
Filename = os.path.basename(os.path.normpath(clim_var_csv_path))
|
|
Clim_var_type = Filename.split('_', 1)[0]
|
|
print(clim_var_csv_path)
|
|
Full_df = pd.read_csv(clim_var_csv_path, index_col=0, parse_dates = True)
|
|
#pandas datestamp index to period (we don't need the 12 pm info in the index (daily periods are enough))
|
|
Full_df.index = Full_df.index.to_period('D')
|
|
#check data types of columns
|
|
#Full_df.dtypes
|
|
|
|
#substract a constant from all values (temp)
|
|
if Clim_var_type == 'tasmean':
|
|
Full_df = Full_df.iloc[:,0:26]-273.15
|
|
|
|
#index the narclim periods if needed
|
|
#Full_df.loc[(Full_df.index >= '1990-01-01') & (Full_df.index < '2010-01-01'), 'period']= '1990-2009'
|
|
#Full_df.loc[(Full_df.index >= '2020-01-01') & (Full_df.index < '2040-01-01'), 'period']= '2020-2039'
|
|
#Full_df.loc[(Full_df.index >= '2060-01-01') & (Full_df.index < '2080-01-01'), 'period']= '2060-2079'
|
|
|
|
#Subset the data to the minimum base period and above (used to set the lenght of the present day climate period)
|
|
Fdf_1900_2080 = Full_df.loc[(Full_df.index >= Base_period_start) & (Full_df.index < Base_period_end)]
|
|
#Aggregate daily df to annual time series
|
|
if (Clim_var_type == 'pracc' or Clim_var_type == 'evspsblmean' or Clim_var_type == 'potevpmean'
|
|
or Clim_var_type == 'pr1Hmaxtstep' or Clim_var_type == 'wss1Hmaxtstep'):
|
|
Fdf_1900_2080_annual = Fdf_1900_2080.resample('A').sum()
|
|
Fdf_1900_2080_annual = Fdf_1900_2080_annual.replace(0, np.nan)
|
|
Fdf_1900_2080_monthly = Fdf_1900_2080.resample('M').sum()
|
|
Fdf_1900_2080_monthly = Fdf_1900_2080_monthly.replace(0, np.nan)
|
|
Fdf_1900_2080_weekly = Fdf_1900_2080.resample('W').sum()
|
|
Fdf_1900_2080_weekly = Fdf_1900_2080_weekly.replace(0, np.nan)
|
|
Fdf_Seas_means = Fdf_1900_2080.resample('Q-NOV').sum() #seasonal means
|
|
Fdf_Seas_means = Fdf_Seas_means.replace(0, np.nan)
|
|
else:
|
|
Fdf_1900_2080_annual = Fdf_1900_2080.resample('A').mean()
|
|
Fdf_1900_2080_monthly = Fdf_1900_2080.resample('M').mean()
|
|
Fdf_1900_2080_weekly = Fdf_1900_2080.resample('W').mean()
|
|
Fdf_Seas_means = Fdf_1900_2080.resample('Q-NOV').mean() #seasonal means
|
|
#plot the mean of all model runs
|
|
print('-------------------------------------------')
|
|
print('mean of all models for climate variable: ' + Clim_var_type)
|
|
Fdf_1900_2080_means = Fdf_1900_2080.mean()
|
|
Fdf_1900_2080_means.plot(kind='bar', ylim=(16,22)).figure
|
|
print('-------------------------------------------')
|
|
|
|
#time series plots:
|
|
#ranks = Fdf_1900_2080_means[3:].rank(axis=0)
|
|
df3 = Fdf_1900_2080_annual
|
|
Fdf_annual_sorted = df3.reindex_axis(df3.mean().sort_values().index, axis=1)
|
|
Fdf_annual_sorted_subset = Fdf_annual_sorted.iloc[:,[0,1,2,3,5,7,13,15,18, 22, 24,25]]
|
|
|
|
#write the key plots to a single pdf document
|
|
pdf_out_file_name = Clim_var_type + '_start_' + Base_period_start + '_NARCliM_summary4.pdf'
|
|
pdf_out_path = output_directory +'/' + pdf_out_file_name
|
|
#open pdf and add the plots
|
|
|
|
#developement
|
|
#I want to create a simple density plot for all values in present day, near future and far future
|
|
#this should be an easy indicator of change
|
|
#subset present day simulations
|
|
Full_current_df = Fdf_1900_2080.iloc[:,[0,1,2]]
|
|
Full_current_df = Full_current_df.stack()
|
|
#nearfuture
|
|
Full_nearfuture_df = Fdf_1900_2080.iloc[:,range(3,15)]
|
|
Full_nearfuture_df = Full_nearfuture_df.stack()
|
|
#farfuture
|
|
Full_farfuture_df = Fdf_1900_2080.iloc[:,range(15,27)]
|
|
Full_farfuture_df = Full_farfuture_df.stack()
|
|
Summarized_df = pd.concat([Full_current_df, Full_nearfuture_df], axis=1, ignore_index=True)
|
|
Summarized_df = pd.concat([Summarized_df, Full_farfuture_df], axis=1, ignore_index=True)
|
|
Summarized_df.columns = ['presentday', 'nearfuture', 'farfuture']
|
|
#end of developement
|
|
with PdfPages(pdf_out_path) as pdf:
|
|
#barplot of model means
|
|
plt.title(Clim_var_type + ' - model means - full period')
|
|
ymin = min(Fdf_1900_2080_means)
|
|
ymax = max(Fdf_1900_2080_means)
|
|
Fdf_1900_2080_means.plot(kind='bar', ylim=(ymin,ymax))
|
|
Fdf_1900_2080_means.plot(kind='bar')
|
|
pdf.savefig(bbox_inches='tight', pad_inches=0.4)
|
|
plt.close()
|
|
#full period density comparison
|
|
plt.title(Clim_var_type + ' - density comparison - full period')
|
|
Summarized_df.plot.kde()
|
|
pdf.savefig(bbox_inches='tight', pad_inches=0.4)
|
|
plt.close()
|
|
#annual box
|
|
plt.title(Clim_var_type + ' - Annual means')
|
|
Fdf_1900_2080_annual.boxplot(rot=90)
|
|
pdf.savefig(bbox_inches='tight', pad_inches=0.4)
|
|
plt.close()
|
|
#monthly box
|
|
plt.title(Clim_var_type + ' - Monthly means')
|
|
Fdf_1900_2080_monthly.boxplot(rot=90)
|
|
pdf.savefig(bbox_inches='tight', pad_inches=0.4)
|
|
plt.close()
|
|
#weekly box
|
|
plt.title(Clim_var_type + ' - Weekly means')
|
|
Fdf_1900_2080_weekly.boxplot(rot=90)
|
|
pdf.savefig(bbox_inches='tight', pad_inches=0.4)
|
|
plt.close()
|
|
#daily box
|
|
plt.title(Clim_var_type + ' - Daily means')
|
|
Fdf_1900_2080.boxplot(rot=90)
|
|
pdf.savefig(bbox_inches='tight', pad_inches=0.4)
|
|
plt.close()
|
|
# time series plot annual ALL models
|
|
plt.title(Clim_var_type + ' - Time series - all models')
|
|
Fdf_1900_2080_annual.plot(legend=False)
|
|
pdf.savefig(bbox_inches='tight', pad_inches=0.4)
|
|
plt.close()
|
|
#
|
|
plt.title(Clim_var_type + ' - Time series - representative models')
|
|
Fdf_annual_sorted_subset.plot(legend=False)
|
|
pdf.savefig(bbox_inches='tight', pad_inches=0.4)
|
|
plt.close()
|
|
# seasonal mean boxplots
|
|
ymin = min(Fdf_Seas_means[Fdf_Seas_means.index.quarter==1].mean())
|
|
ymax = max(Fdf_Seas_means[Fdf_Seas_means.index.quarter==1].mean())
|
|
plt.title(Clim_var_type + ' - DJF Summer means')
|
|
Fdf_Seas_means[Fdf_Seas_means.index.quarter==1].mean().plot(kind='bar', ylim=(ymin,ymax))
|
|
pdf.savefig(bbox_inches='tight', pad_inches=0.4)
|
|
plt.close()
|
|
plt.title(Clim_var_type + ' - DJF Summer means')
|
|
Fdf_Seas_means[Fdf_Seas_means.index.quarter==1].boxplot(rot=90)
|
|
pdf.savefig(bbox_inches='tight', pad_inches=0.4)
|
|
plt.close()
|
|
ymin = min(Fdf_Seas_means[Fdf_Seas_means.index.quarter==2].mean())
|
|
ymax = max(Fdf_Seas_means[Fdf_Seas_means.index.quarter==2].mean())
|
|
plt.title(Clim_var_type + ' - MAM Autumn means')
|
|
Fdf_Seas_means[Fdf_Seas_means.index.quarter==2].mean().plot(kind='bar', ylim=(ymin,ymax))
|
|
pdf.savefig(bbox_inches='tight', pad_inches=0.4)
|
|
plt.close()
|
|
plt.title(Clim_var_type + ' - MAM Autumn means')
|
|
Fdf_Seas_means[Fdf_Seas_means.index.quarter==2].boxplot(rot=90)
|
|
pdf.savefig(bbox_inches='tight', pad_inches=0.4)
|
|
plt.close()
|
|
ymin = min(Fdf_Seas_means[Fdf_Seas_means.index.quarter==3].mean())
|
|
ymax = max(Fdf_Seas_means[Fdf_Seas_means.index.quarter==3].mean())
|
|
plt.title(Clim_var_type + ' - JJA Winter means')
|
|
Fdf_Seas_means[Fdf_Seas_means.index.quarter==3].mean().plot(kind='bar', ylim=(ymin,ymax))
|
|
pdf.savefig(bbox_inches='tight', pad_inches=0.4)
|
|
plt.close()
|
|
plt.title(Clim_var_type + ' - JJA Winter means')
|
|
Fdf_Seas_means[Fdf_Seas_means.index.quarter==3].boxplot(rot=90)
|
|
pdf.savefig(bbox_inches='tight', pad_inches=0.4)
|
|
plt.close()
|
|
ymin = min(Fdf_Seas_means[Fdf_Seas_means.index.quarter==4].mean())
|
|
ymax = max(Fdf_Seas_means[Fdf_Seas_means.index.quarter==4].mean())
|
|
plt.title(Clim_var_type + ' - SON Spring means')
|
|
Fdf_Seas_means[Fdf_Seas_means.index.quarter==4].mean().plot(kind='bar', ylim=(ymin,ymax))
|
|
pdf.savefig(bbox_inches='tight', pad_inches=0.4)
|
|
plt.close()
|
|
plt.title(Clim_var_type + ' - SON Spring means')
|
|
Fdf_Seas_means[Fdf_Seas_means.index.quarter==4].boxplot(rot=90)
|
|
pdf.savefig(bbox_inches='tight', pad_inches=0.4)
|
|
plt.close()
|
|
|
|
#plots not used
|
|
#Fdf_annual_sorted_subset.plot(legend=False, subplots=True)
|