# -*- coding: utf-8 -*- """ Created on Mon Nov 25 10:56:21 2019 @author: z5025317 """ #==========================================================# #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 import seaborn matplotlib.style.use('ggplot') #========================================================# #========================================================# #set input parameters and paths #========================================================# Gauges = [ 'Stockton', 'Raymond', 'Oakhampton', 'Paterson'] #name of all gauges or IDs in a given river River = 'Clarence' #name of plot scenario/experiment . This is a unique identifier for this particular scenario comparison Plot_scenario = 'Clarence_gauge_analysis_V1' #time frames for analysis tideanalysis_day = '2006 06 07' # the time series plots will be centred around this day. year = 2006 #set directories Folder = 'C:/Users/z5025317/OneDrive - UNSW/WRL_Postdoc_Manual_Backup/WRL_Postdoc/Projects/Paper#2_HD_modeling/01_Data/Tides/CSV/' #needs to be adjusted - should be a folder with all the MHL data in individual CSV files Out_Folder = 'C:/Users/z5025317/OneDrive - UNSW/WRL_Postdoc_Manual_Backup/WRL_Postdoc/Projects/Paper#2_HD_modeling/02_Analysis/Tides/' #needs to be adjusted #========================================================# #========================================================# #laod MHL data and store in a single dataframe #========================================================# WL_df = pd.DataFrame() for gauge in Gauges: # gauge = 'Oyster' print 'loading data for ' + gauge CSVpath = glob.glob(Folder + '*' + gauge + '*' ) # will find the csv file corresponding to the gauge name or ID print(CSVpath) Full_df = pd.read_csv(CSVpath ,header= 0, skiprows= 30, parse_dates = True) #needs to be adjusted to fit the MHL format #Full_df.index = pd.to_datetime(Full_df.Date, format = '%d/%m/%Y') + pd.to_timedelta(Full_df.Time, unit='h') #time series index may have to be adjusted, depending on the outcome of parse_dates Full_df.index = pd.to_datetime(Full_df.Date, format = '%d/%m/%Y') + pd.to_timedelta(Full_df.Time, unit='h') Full_df = Full_df.iloc[:,[-2]] Full_df.columns = ['WL[AHD]'] Full_df = Full_df.loc[(Full_df.index >= year + '-01-01 00:00:00') & (Full_df.index <= str(np.int(year)+1) + '-01-01 00:00:00')] Full_df = Full_df.loc[~Full_df.index.duplicated(keep='first')] #get rid of possible duplicate datetimes WL_df = pd.concat([WL_df, Full_df], axis=1) WL_df.columns = Gauges #subset the data in time if desired (e.g. to a given year) WL_df_1 = WL_df.loc[(WL_df.index >= '2017-01-10 00:00:00') & (WL_df.index <= '2017-01-18 00:00:00')] #better #WL_df_1 = WL_df.loc[(WL_df.index >= datetime.strptime(startdate , '%Y %m %d').date()) & (WL_df.index <= datetime.strptime(enddate , '%Y %m %d').date())] #may need to be adjusted #export dataframe to CSV if desired #export the pandas data frame as a CSV file within the output directory is desired #WL_df.to_csv(Folder + 'data_combined.csv') #========================================================# #========================================================# # Plot water level histogram at all gauges on a singel plot #========================================================# #full period density comparison all four gauges png_out_name = Out_Folder + River + '_water level histograms.pdf' fig = plt.figure(figsize=(40,18)) ax=plt.subplot(1,1,1) plt.title('Hunter water level density comparison - Jan 2016- Feb 2020') WL_df.iloc[:,[0,1,2,3]].plot.kde(xlim=(-1.5,1.5), ax=ax) fig.tight_layout() fig.savefig(png_out_name) plt.close() fig.tight_layout() fig.savefig(png_out_name) plt.close() #========================================================# #========================================================# # Plot water level histogram at all gauges #========================================================# ylims = [-0.9,2.2] #set plot parameters ALPHA_figs = 1 font = {'family' : 'sans-serif', 'weight' : 'normal', 'size' : 6} matplotlib.rc('font', **font) #plot i=1 png_out_name = Out_Folder + '/' + Plot_scenario + '/' + River + 'B_water level histogram and exceedance planes '+ startdate + '_' + enddate + '3.pdf' fig = plt.figure(figsize=(len(Gauges)*4,8)) for Gauge in Gauges: #Gauge = 'Yamba' xlocation = 0.6 #full period density comparison for individual gauges overlaid by tidal planes Wollow_13ExPr = np.nanpercentile(WL_df[Gauge],87) Wollow_6p4ExPr = np.nanpercentile(WL_df[Gauge],93.6) Wollow_1_3ExPr = np.nanpercentile(WL_df[Gauge],98.7) Wollow_0p01ExPr = np.nanpercentile(WL_df[Gauge],99.99) Wollow_50ExPr = np.nanpercentile(WL_df[Gauge],50) ax=plt.subplot(4,1,i) plt.title( Gauge + ' water level density comparison - Jan 2016- Feb 2020') #WL_df[Gauge].plot.kde(xlim=(-0.5,1),vertical=True, lw=3, ax=ax) seaborn.kdeplot(WL_df[Gauge], vertical=True, color='royalblue', lw=3, ax=ax) plt.axhline(y=Wollow_13ExPr, color='firebrick', linestyle='--', lw=3, alpha=1) plt.text( xlocation , (Wollow_13ExPr), '13% exceedance', ha='left', va='bottom') plt.axhline(y=Wollow_6p4ExPr, color='orange', linestyle='--', lw=3, alpha=1) plt.text(xlocation , (Wollow_6p4ExPr), '6.4% exceedance', ha='left', va='bottom') plt.axhline(y=Wollow_1_3ExPr, color='yellow', linestyle='--', lw=3, alpha=1) plt.text( xlocation ,Wollow_1_3ExPr , '1.3% exceedance', ha='left', va='bottom') plt.axhline(y=Wollow_0p01ExPr, color='limegreen', linestyle='--', lw=3, alpha=1) plt.text( xlocation , Wollow_0p01ExPr , '0.01% exceedance', ha='left', va='bottom') plt.axhline(y=Wollow_50ExPr, color='cyan', linestyle='--', lw=3, alpha=1) plt.text( xlocation , Wollow_50ExPr , '50% exceedance', ha='left', va='bottom') plt.ylim(-1,1.7) plt.xlabel("Probability/density") plt.ylabel("Water level [m AHD]") fig.tight_layout() fig.savefig(png_out_name) plt.close() fig.tight_layout() fig.savefig(png_out_name) i=i+1 plt.close() #========================================================# #========================================================# # Save exceedance planes to csv #========================================================# Summary={} Summary['13% exceedance'] = Wollow_13ExPr Summary['6.4% exceedance'] = Wollow_6p4ExPr Summary['1.3% exceedance'] = Wollow_1_3ExPr Summary['0.01% exceedance'] = Wollow_0p01ExPr Summary['50% exceedance'] = Wollow_50ExPr pdf1 = pd.DataFrame(Summary, index=[0]) pdf1.to_csv(os.path.join( Out_Folder + Gauge + '_entrance_stats.csv')) #========================================================# #========================================================# # Plot water level time series at control points #========================================================# #set plot parameters ALPHA_figs = 1 font = {'family' : 'sans-serif', 'weight' : 'normal', 'size' : 5} matplotlib.rc('font', **font) halfwindows = [1,2,5] #Days to add before and after central date chosen for analysis for halfwindow in halfwindows: png_out_name = Out_Folder + '/' + Plot_scenario + '/' +River+ 'water level time series '+ startdate + '_' + enddate + ' _window_of_' + str(halfwindow*2) +'_days.pdf' fig = plt.figure(figsize=(10,len(Gauges)*3)) i=1 for Gauge in Gauges: Gauge_only_df = WL_df.filter(regex=Gauge) #Gauge_only_df_1 = Gauge_only_df.loc[(Gauge_only_df.index >= '2008-05-09 12:00:00') & (Gauge_only_df.index <= '2008-05-20 12:00:00')] Gauge_only_df_1 = Gauge_only_df[datetime.strptime(tideanalysis_day , '%Y %m %d').date()- timedelta(days=halfwindow) : datetime.strptime(tideanalysis_day , '%Y %m %d').date() + timedelta(days=halfwindow) ] ax=plt.subplot(len(Gauges),1,i) plt.title('Tidal water levels at 4 locations in the ' + River +' Estuary') ax.set_prop_cycle(color=['steelblue', 'coral', 'steelblue', 'deeppink'],linestyle=[ '-', '-','--','-'], lw=[1,1,1,2]) Gauge_only_df_1.iloc[:,[0,1,2]].plot( legend=False, ax=ax) plt.axhline(y=0, color='slateblue', linestyle='--', lw=1, alpha=0.5) #plt.axhline(y=0.9, color='tomato', linestyle='--', lw=1, alpha=0.5) plt.ylim(ylims) plt.ylabel("Water level [m AHD]") plt.xlabel("Time") ax.legend() fig.tight_layout() i=i+1 fig.savefig(png_out_name) plt.close() #========================================================# # Plot water level histograms at control points #========================================================#