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.
222 lines
8.2 KiB
Python
222 lines
8.2 KiB
Python
5 years ago
|
# -*- coding: utf-8 -*-
|
||
|
|
#==========================================================#
|
||
|
|
#Last Updated - June 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
|
||
|
|
from ggplot import *
|
||
|
|
matplotlib.style.use('ggplot')
|
||
|
|
|
||
|
|
#==========================================================#
|
||
|
|
#Input parameters
|
||
|
|
#==========================================================#
|
||
|
|
# Set working direcotry (where postprocessed NARClIM data is located)
|
||
|
|
#set beginning and end years and corresponding scenario code
|
||
|
|
fs=['Hwq003', 'Hwq005']
|
||
|
|
|
||
|
|
comparison_code = '003v005'
|
||
|
|
startyear=1999
|
||
|
|
endyear=1999
|
||
|
|
year=range(startyear, endyear+1)
|
||
|
|
|
||
|
|
#set directory path for in and output files
|
||
|
|
output_directory = 'C:/Users/z5025317/OneDrive - UNSW/Hunter_CC_Modeling/08_Results/Output/Postprocessed' + comparison_code +'/'
|
||
|
|
nodes_csv = 'H:/WRL_Projects/Hunter_CC_Modeling/Module_6/03_Results/Chainages/Hunter_nodes.csv'
|
||
|
|
Flow_csv = 'H:/WRL_Projects/Hunter_CC_Modeling/Module_6/01_Input/BCGeneration/Scenarios/Calibration/Calibration_Greta.csv'
|
||
|
|
png_output_directory = output_directory + '/Figures/'
|
||
|
|
variables = ['Sal'] #Tem or Sal
|
||
|
|
|
||
|
|
#read csv file with nodes to extract data from
|
||
|
|
#node=[837, 666, 981, 59,527, 34, 1, 391]
|
||
|
|
node = pd.read_csv(nodes_csv)['Hunter'].values
|
||
|
|
chainages = pd.read_csv(nodes_csv)['x_km'].values
|
||
|
|
|
||
|
|
|
||
|
|
#set plot parameters
|
||
|
|
ALPHA_figs = 1
|
||
|
|
font = {'family' : 'sans-serif',
|
||
|
|
'weight' : 'normal',
|
||
|
|
'size' : 14}
|
||
|
|
matplotlib.rc('font', **font)
|
||
|
|
#==========================================================#
|
||
|
|
|
||
|
|
#==========================================================#
|
||
|
|
#100% automated part of the code doing the data extraction
|
||
|
|
#==========================================================
|
||
|
|
#output_directory = 'J:/Project wrl2016032/NARCLIM_Raw_Data/Extracted'
|
||
|
|
if not os.path.exists(output_directory):
|
||
|
|
os.makedirs(output_directory)
|
||
|
|
#output_directory = 'J:/Project wrl2016032/NARCLIM_Raw_Data/Extracted'
|
||
|
|
if not os.path.exists(png_output_directory):
|
||
|
|
os.makedirs(png_output_directory)
|
||
|
|
|
||
|
|
|
||
|
|
Summary_df = pd.DataFrame()
|
||
|
|
for f in fs:
|
||
|
|
for NODE in node:
|
||
|
|
NODE = str(NODE)
|
||
|
|
#set input and output directories
|
||
|
|
input_directory = 'C:/Users/z5025317/OneDrive - UNSW/Hunter_CC_Modeling/08_Results/Output/'+ f
|
||
|
|
#set directory path for individual png files
|
||
|
|
png_output_directory = output_directory + '/Figures/'
|
||
|
|
# Set working direcotry (where postprocessed NARClIM data is located)
|
||
|
|
os.chdir(input_directory)
|
||
|
|
#==========================================================#
|
||
|
|
#Load data file
|
||
|
|
Clim_Var_CSVs = glob.glob('*'+ NODE + '*WQ*')
|
||
|
|
clim_var_csv_path = Clim_Var_CSVs[0]
|
||
|
|
df = pd.read_csv(clim_var_csv_path, index_col=False, sep=' ')
|
||
|
|
df.index = pd.to_datetime(df.Year, format = '%Y') + pd.to_timedelta(df.Hour, unit='h')
|
||
|
|
df= df.drop(columns=['Year', 'Hour'])
|
||
|
|
df.columns = [NODE+'_Sal', NODE+'_Tem']
|
||
|
|
#append to summary data frame
|
||
|
|
if f in ['HCC010', 'HCC042']:
|
||
|
|
scname = 'Present'
|
||
|
|
if f == 'HCC011':
|
||
|
|
scname = 'Near_Future'
|
||
|
|
df.index = df.index - (datetime.strptime('2025 01 01', '%Y %m %d').date() - datetime.strptime('1995 01 01', '%Y %m %d').date())
|
||
|
|
if f == 'HCC012':
|
||
|
|
scname = 'Far_Future'
|
||
|
|
df.index = df.index - (datetime.strptime('2065 01 01', '%Y %m %d').date() - datetime.strptime('1995 01 01', '%Y %m %d').date())
|
||
|
|
if f in ['HCC040', 'HCC041']:
|
||
|
|
scname = 'Present_Resto'
|
||
|
|
df.columns = df.columns+'_'+ scname
|
||
|
|
Summary_df = pd.concat([Summary_df, df], axis=1, join='outer')
|
||
|
|
|
||
|
|
#cut down the summary df to common years
|
||
|
|
Summary_df = Summary_df[datetime.strptime(str(startyear) + ' 01 05', '%Y %m %d').date():datetime.strptime(str(endyear) + ' 12 31', '%Y %m %d').date()]
|
||
|
|
Summary_df.tail()
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
for variable in variables:
|
||
|
|
ii=1
|
||
|
|
for NODE in node:
|
||
|
|
NODE = str(NODE)
|
||
|
|
#=========#
|
||
|
|
#comparison time series plot at node
|
||
|
|
#=========#
|
||
|
|
#out name
|
||
|
|
png_out_file_name = comparison_code + '_' + variable + '_' + NODE + '_time_series.png'
|
||
|
|
png_out_path = png_output_directory + '/' + png_out_file_name
|
||
|
|
#
|
||
|
|
plotin_df = Summary_df.filter(regex=variable)
|
||
|
|
plotin_df = plotin_df.filter(regex=NODE)
|
||
|
|
#prepare data frame for plot
|
||
|
|
fig = plt.figure(figsize=(14,18))
|
||
|
|
ax=plt.subplot(4,1,1)
|
||
|
|
#start and end dat of plots
|
||
|
|
plt.title(variable + '_' + NODE + ' - time series')
|
||
|
|
start_date = datetime(int('1996'), int('02'), int('01'))
|
||
|
|
end_date = datetime(int('1996'), int('02'), int('5'))
|
||
|
|
plotin_df1=plotin_df[start_date:end_date]
|
||
|
|
ymin = min(plotin_df1.min()) - 0.01 *min(plotin_df1.min())
|
||
|
|
ymax = max(plotin_df1.max()) + 0.01 *max(plotin_df1.max())
|
||
|
|
plotin_df1.plot(ylim=(ymin,ymax) ,legend=False, ax=ax)
|
||
|
|
#plt.legend( loc=1)
|
||
|
|
plt.xticks([])
|
||
|
|
|
||
|
|
ax=plt.subplot(4,1,2)
|
||
|
|
start_date = datetime(int('1996'), int('02'), int('01'))
|
||
|
|
end_date = datetime(int('1996'), int('02'), int('3'))
|
||
|
|
plotin_df2=plotin_df[start_date:end_date]
|
||
|
|
ymin = min(plotin_df2.min()) - 0.1 *min(plotin_df2.min())
|
||
|
|
ymax = max(plotin_df2.max()) + 0.1 *max(plotin_df2.max())
|
||
|
|
plotin_df2.plot(ylim=(ymin,ymax),legend=False, ax=ax)
|
||
|
|
#plt.legend( loc=1)
|
||
|
|
ax.xaxis.grid(False)
|
||
|
|
|
||
|
|
ax=plt.subplot(4,1,3)
|
||
|
|
plt.title(variable + '_' + NODE + ' - monthly maximum')
|
||
|
|
Monthly_max = plotin_df.resample('M').max()
|
||
|
|
Monthly_max.plot(ax=ax,legend=False)
|
||
|
|
plt.legend( loc=1)
|
||
|
|
ax.xaxis.grid(False)
|
||
|
|
|
||
|
|
ax=plt.subplot(4,1,4)
|
||
|
|
plt.title(variable + '_' + NODE + ' - monthly mean')
|
||
|
|
Monthly_mean = plotin_df.resample('M').mean()
|
||
|
|
Monthly_mean.plot(ax=ax)
|
||
|
|
#plt.legend(bbox_to_anchor=(0., 1.02, 1., .102), loc=3,
|
||
|
|
# ncol=2, mode="expand", borderaxespad=0.)
|
||
|
|
#plt.legend( loc=1)
|
||
|
|
ax.xaxis.grid(False)
|
||
|
|
|
||
|
|
#ax.patch.set_alpha(ALPHA_figs)
|
||
|
|
fig.patch.set_alpha(ALPHA_figs)
|
||
|
|
fig.tight_layout()
|
||
|
|
fig.savefig(png_out_path)
|
||
|
|
plt.close()
|
||
|
|
#=========#
|
||
|
|
|
||
|
|
png_out_file_name = comparison_code + '_' + variable + '_' + NODE + '_exceedance.png'
|
||
|
|
png_out_path = png_output_directory + '/' + png_out_file_name
|
||
|
|
fig = plt.figure(figsize=(10,5))
|
||
|
|
ax=plt.subplot(1,1,1)
|
||
|
|
start_date = datetime(int('1995'), int('02'), int('01'))
|
||
|
|
end_date = datetime(int('1998'), int('12'), int('31'))
|
||
|
|
|
||
|
|
plotin_df3=plotin_df[start_date:end_date]
|
||
|
|
Qhist = plotin_df3.filter(regex='Present')
|
||
|
|
col = Qhist.columns[0]
|
||
|
|
Qhist = Qhist.sort_values(by=[col],ascending=False).values
|
||
|
|
Qnearf = plotin_df3.filter(regex='Near_Future')
|
||
|
|
col = Qnearf.columns[0]
|
||
|
|
Qnearf = Qnearf.sort_values(by=[col],ascending=False).values
|
||
|
|
Qmodern = plotin_df3.filter(regex='Far_Future')
|
||
|
|
col = Qmodern.columns[0]
|
||
|
|
Qmodern = Qmodern.sort_values(by=[col],ascending=False).values
|
||
|
|
Nh = len(Qhist) # number of data points
|
||
|
|
Nm = len(Qmodern)
|
||
|
|
# here's one place types matter: must divide by float not int
|
||
|
|
plt.plot(np.arange(Nh)/float(Nh)*100, Qhist)
|
||
|
|
plt.plot(np.arange(Nm)/float(Nm)*100, Qnearf)
|
||
|
|
plt.plot(np.arange(Nm)/float(Nm)*100, Qmodern)
|
||
|
|
#plt.plot(np.sort(Qmodern)[::-1])
|
||
|
|
plt.xlabel('% Exceedance')
|
||
|
|
if variable == 'Sal':
|
||
|
|
plt.ylabel('Salinity [PSU]')
|
||
|
|
if variable == 'Tem':
|
||
|
|
plt.ylabel('Temperature [C]')
|
||
|
|
plt.legend(['Present', 'Near Future', 'Far Future'])
|
||
|
|
ax.xaxis.grid(False)
|
||
|
|
#ax.patch.set_alpha(ALPHA_figs)
|
||
|
|
fig.patch.set_alpha(ALPHA_figs)
|
||
|
|
fig.tight_layout()
|
||
|
|
fig.savefig(png_out_path)
|
||
|
|
plt.close()
|
||
|
|
|