#!/usr/bin/env python # coding: utf-8 # In[2]: import cmocean.cm as cm import datetime as dt import matplotlib.pyplot as plt import matplotlib.patches as mpatches from matplotlib.colors import LogNorm import numpy as np import pandas as pd # import statfieldodels.api as sm import xarray as xr from salishsea_tools import viz_tools import netCDF4 as nc import arrow import copy import math import matplotlib.dates as mdates from scipy.io import loadmat # In[3]: mycmap = copy.copy(cm.oxy) mycmap.set_bad('darkgreen') # ERDDAP server - https://nwem.apl.washington.edu/erddap/index.html # https://nwem.apl.washington.edu/erddap/search/index.html?page=1&itemsPerPage=1000&searchFor=Oxygen - data # In[4]: with xr.open_dataset('/ocean/atall/MOAD/Obs/PugetSound/ERDDAP/L2_gridded_025_TW_20070101-20071231.nc') as TW_data: oxy_TW_obs = TW_data['mass_concentration_of_oxygen_in_sea_water']*44.64 temp_TW_obs = TW_data['sea_water_temperature'] salt_TW_obs = TW_data['sea_water_practical_salinity'] no3_TW_obs = TW_data['mole_concentration_of_nitrate_in_sea_water'] chl_TW_obs = TW_data['mass_concentration_of_chlorophyll_a_in_sea_water'] time_TW = TW_data['cast_start_time'] # In[5]: oxy_TW_obs.shape # In[10]: fig, ax = plt.subplots(5, 1, figsize=(15, 18)) temp_TW_obs.plot(ax=ax[0], y="sea_water_pressure", yincrease=False, cmap='gist_rainbow_r', vmin=3, vmax=21, ylim=(30, 0)); salt_TW_obs.plot(ax=ax[1], y="sea_water_pressure", yincrease=False, cmap= 'gist_rainbow_r', vmin=16, vmax=31, ylim=(30, 0)); oxy_TW_obs.plot(ax=ax[2], y="sea_water_pressure", yincrease=False, cmap='gist_rainbow_r', vmin=0, vmax=400, ylim=(30, 0)); chl_TW_obs.plot(ax=ax[3], y="sea_water_pressure", yincrease=False, cmap='gist_rainbow_r', vmin=0, vmax=25, ylim=(30, 0)); no3_TW_obs.plot(ax=ax[4], y="sea_water_pressure", yincrease=False, cmap='gist_rainbow_r', vmin=0, vmax=50, ylim=(30, 0)); # ORCA RAWdata - https://nwem.apl.washington.edu/prod_DataReq.shtml # In[112]: HC_HP_data = loadmat('/ocean/atall/MOAD/Obs/PugetSound/ORCA_rawdata/HC_HP_CTD_data_bin_web.mat') # Hood Canal - Hoodsport HC_TW_data = loadmat('/ocean/atall/MOAD/Obs/PugetSound/ORCA_rawdata/HC_TW_CTD_data_bin_web.mat') # Hood Canal - Lynch Cove - Twanoh HC_DB_data = loadmat('/ocean/atall/MOAD/Obs/PugetSound/ORCA_rawdata/HC_DB_CTD_data_bin_web.mat') # Hood Canal - Duckabush HC_DK_data = loadmat('/ocean/atall/MOAD/Obs/PugetSound/ORCA_rawdata/HC_DK_CTD_data_bin_web.mat') # Hood Canal - North buoy HC_NB_data = loadmat('/ocean/atall/MOAD/Obs/PugetSound/ORCA_rawdata/HC_NB_CTD_data_bin_web.mat') # Hood Canal - Hoodsport CI_data = loadmat('/ocean/atall/MOAD/Obs/PugetSound/ORCA_rawdata/CI_CTD_data_bin_web.mat') # Carr Inlet PW_data = loadmat('/ocean/atall/MOAD/Obs/PugetSound/ORCA_rawdata/PW_CTD_data_bin_web.mat') # Point Wells time_HP = HC_HP_data.get('Btime') time_TW = HC_TW_data.get('Btime') time_DB = HC_DB_data.get('Btime') time_DK = HC_DK_data.get('Btime') time_NB = HC_NB_data.get('Btime') time_CI = CI_data.get('Btime') time_PW = PW_data.get('Btime') depth_HP_obs = HC_HP_data.get('Bdepth') depth_TW_obs = HC_TW_data.get('Bdepth') depth_DB_obs = HC_DB_data.get('Bdepth') depth_DK_obs = HC_DK_data.get('Bdepth') depth_NB_obs = HC_NB_data.get('Bdepth') depth_CI_obs = CI_data.get('Bdepth') depth_PW_obs = PW_data.get('Bdepth') oxy_HP_obs = HC_HP_data.get('Boxy_umol') oxy_TW_obs = HC_TW_data.get('Boxy_umol') oxy_DB_obs = HC_DB_data.get('Boxy_umol') oxy_DK_obs = HC_DK_data.get('Boxy_umol') oxy_NB_obs = HC_NB_data.get('Boxy_umol') oxy_CI_obs = CI_data.get('Boxy_umol') oxy_PW_obs = PW_data.get('Boxy_umol') # In[70]: time_HP.shape # In[118]: depth_TW_obs = np.ma.masked_array(depth_TW_obs) oxy_TW_obs = np.ma.masked_array(oxy_TW_obs) #time_TW = np.ma.masked_invalid(time_TW) #plt.contourf(time_HP,-depth_HP_obs,oxy_HP_obs, levels=[10, 30, 500]) plt.pcolormesh(time_TW) # In[52]: con_list = [[element for element in upperElement] for upperElement in HC_HP_data['__globals__']] #con_list.shape # In[55]: # zip provides us with both the x and y in a tuple. newData = list(zip(con_list)) columns = ['obj_contour_x', 'obj_contour_y'] df = pd.DataFrame(newData) print(df) # In[51]: df.plot() # In[11]: with xr.open_dataset('/home/sallen/MEOPAR/grid/mesh_mask202108.nc') as mesh: tmask = mesh.tmask mbathy = mesh.mbathy # In[12]: # Constant and data rang syear = 2022 eyear = 2022 smonth = 8 emonth = 12 startdate = arrow.get(syear, smonth, 1) enddate = arrow.get(eyear, emonth, 31) print (startdate, enddate) # In[13]: with xr.open_dataset('/results2/SalishSea/nowcast-green.202111/01apr23/SalishSea_1d_20230401_20230401_chem_T.nc') as data_chem: oxy = data_chem['dissolved_oxygen'] with xr.open_dataset('/results2/SalishSea/nowcast-green.202111/01apr23/SalishSea_1d_20230401_20230401_biol_T.nc') as data_biol: no3 = data_biol['nitrate'] don = data_biol['dissolved_organic_nitrogen'] pon = data_biol['particulate_organic_nitrogen'] dia = data_biol['diatoms'] fla = data_biol['flagellates'] with xr.open_dataset('/results2/SalishSea/nowcast-green.202111/01apr23/SalishSea_1d_20230401_20230401_grid_T.nc') as data_grid: temp = data_grid['votemper'] salt = data_grid['vosaline'] # dens = data_grid['sigma_theta'] # almost same as salinity # In[14]: tmask2_oxy = tmask.assign_coords(z=("z", oxy.deptht.data)) tmask2_oxy = tmask2_oxy.rename({"z": "deptht"}) tmask2_no3 = tmask.assign_coords(z=("z", no3.deptht.data)) tmask2_no3 = tmask2_no3.rename({"z": "deptht"}) tmask2_don = tmask.assign_coords(z=("z", don.deptht.data)) tmask2_don = tmask2_don.rename({"z": "deptht"}) tmask2_pon = tmask.assign_coords(z=("z", pon.deptht.data)) tmask2_pon = tmask2_pon.rename({"z": "deptht"}) tmask2_dia = tmask.assign_coords(z=("z", dia.deptht.data)) tmask2_dia = tmask2_dia.rename({"z": "deptht"}) tmask2_fla = tmask.assign_coords(z=("z", fla.deptht.data)) tmask2_fla = tmask2_fla.rename({"z": "deptht"}) tmask2_temp = tmask.assign_coords(z=("z", temp.deptht.data)) tmask2_temp = tmask2_temp.rename({"z": "deptht"}) tmask2_salt = tmask.assign_coords(z=("z", salt.deptht.data)) tmask2_salt = tmask2_salt.rename({"z": "deptht"}) #tmask2_dens = tmask.assign_coords(z=("z", dens.deptht.data)) #tmask2_dens = tmask2_dens.rename({"z": "deptht"}) # In[15]: print (oxy.deptht) # In[16]: first = True i1, i2, j1, j2 = (95, 200, 68, 160) # Hood Canal all #i1, i2, j1, j2 = (95, 170, 76, 120) # Central Hood Canal #i1, i2, j1, j2 = (106, 119, 69, 79) # #i1lc, i2lc, j1lc, j2lc = (120, 150, 68, 75) # Lynch Cove i1hs, i2hs, j1hs, j2hs = (107, 110, 83, 86) # Hoodsport buoy point in Lynch Cove (Hood Canal) - 47°25'18.5"N 123°06'45.3"W i1lc, i2lc, j1lc, j2lc = (122, 125, 70, 73) # Twanoh buoy point in Lynch Cove (Hood Canal) - 47.37°N, 123.01°W for ii, day in enumerate(arrow.Arrow.range('day', startdate, enddate)): year = day.year yr2 = day.strftime("%y") month = day.month Month = day.strftime("%b").lower() day = day.day # set up filename to follow NEMO conventions filename_chem = f'/results2/SalishSea/nowcast-green.202111/{day:02}{Month}{yr2}/SalishSea_1d_{year}{month:02}{day:02}_{year}{month:02}{day:02}_chem_T.nc' with xr.open_dataset(filename_chem) as ds_chem: oxyhs = ds_chem['dissolved_oxygen'][0, :, j1hs:j2hs, i1hs:i2hs].where(tmask2_oxy[0, :, j1hs:j2hs, i1hs:i2hs] == 1).mean(axis=1).mean(axis=1) oxylc = ds_chem['dissolved_oxygen'][0, :, j1lc:j2lc, i1lc:i2lc].where(tmask2_oxy[0, :, j1lc:j2lc, i1lc:i2lc] == 1).mean(axis=1).mean(axis=1) if first: first = False oxyhs_ts = oxyhs oxylc_ts = oxylc else: oxyhs_ts = xr.concat([oxyhs_ts, oxyhs], dim='time_counter') oxylc_ts = xr.concat([oxylc_ts, oxylc], dim='time_counter') # In[17]: oxylc_ts.shape # In[18]: first = True i1, i2, j1, j2 = (95, 200, 68, 160) # Hood Canal all #i1, i2, j1, j2 = (95, 170, 76, 120) # Central Hood Canal #i1, i2, j1, j2 = (106, 119, 69, 79) # #i1lc, i2lc, j1lc, j2lc = (120, 150, 68, 75) # Lynch Cove i1hs, i2hs, j1hs, j2hs = (107, 110, 83, 86) # Hoodsport buoy point in Lynch Cove (Hood Canal) - 47°25'18.5"N 123°06'45.3"W i1lc, i2lc, j1lc, j2lc = (122, 125, 70, 73) # Twanoh buoy point in Lynch Cove (Hood Canal) - 47.37°N, 123.01°W for ii, day in enumerate(arrow.Arrow.range('day', startdate, enddate)): year = day.year yr2 = day.strftime("%y") month = day.month Month = day.strftime("%b").lower() day = day.day # set up filename to follow NEMO conventions filename_biol = f'/results2/SalishSea/nowcast-green.202111/{day:02}{Month}{yr2}/SalishSea_1d_{year}{month:02}{day:02}_{year}{month:02}{day:02}_biol_T.nc' with xr.open_dataset(filename_biol) as ds_biol: no3hs = ds_biol['nitrate'][0, :, j1hs:j2hs, i1hs:i2hs].where(tmask2_no3[0, :, j1hs:j2hs, i1hs:i2hs] == 1).mean(axis=1).mean(axis=1) no3lc = ds_biol['nitrate'][0, :, j1lc:j2lc, i1lc:i2lc].where(tmask2_no3[0, :, j1lc:j2lc, i1lc:i2lc] == 1).mean(axis=1).mean(axis=1) donlc = ds_biol['dissolved_organic_nitrogen'][0, :, j1lc:j2lc, i1lc:i2lc].where(tmask2_don[0, :, j1lc:j2lc, i1lc:i2lc] == 1).mean(axis=1).mean(axis=1) donhs = ds_biol['dissolved_organic_nitrogen'][0, :, j1hs:j2hs, i1hs:i2hs].where(tmask2_don[0, :, j1hs:j2hs, i1hs:i2hs] == 1).mean(axis=1).mean(axis=1) ponlc = ds_biol['particulate_organic_nitrogen'][0, :, j1lc:j2lc, i1lc:i2lc].where(tmask2_pon[0, :, j1lc:j2lc, i1lc:i2lc] == 1).mean(axis=1).mean(axis=1) ponhs = ds_biol['particulate_organic_nitrogen'][0, :, j1hs:j2hs, i1hs:i2hs].where(tmask2_pon[0, :, j1hs:j2hs, i1hs:i2hs] == 1).mean(axis=1).mean(axis=1) dialc = ds_biol['diatoms'][0, :, j1lc:j2lc, i1lc:i2lc].where(tmask2_dia[0, :, j1lc:j2lc, i1lc:i2lc] == 1).mean(axis=1).mean(axis=1) diahs = ds_biol['diatoms'][0, :, j1hs:j2hs, i1hs:i2hs].where(tmask2_dia[0, :, j1hs:j2hs, i1hs:i2hs] == 1).mean(axis=1).mean(axis=1) flalc = ds_biol['flagellates'][0, :, j1lc:j2lc, i1lc:i2lc].where(tmask2_fla[0, :, j1lc:j2lc, i1lc:i2lc] == 1).mean(axis=1).mean(axis=1) flahs = ds_biol['flagellates'][0, :, j1hs:j2hs, i1hs:i2hs].where(tmask2_fla[0, :, j1hs:j2hs, i1hs:i2hs] == 1).mean(axis=1).mean(axis=1) if first: first = False no3hs_ts = no3hs no3lc_ts = no3lc donlc_ts = donlc donhs_ts = donhs ponlc_ts = ponlc ponhs_ts = ponhs dialc_ts = dialc diahs_ts = diahs flalc_ts = flalc flahs_ts = flahs else: no3hs_ts = xr.concat([no3hs_ts, no3hs], dim='time_counter') no3lc_ts = xr.concat([no3lc_ts, no3lc], dim='time_counter') donlc_ts = xr.concat([donlc_ts, donlc], dim='time_counter') donhs_ts = xr.concat([donhs_ts, donhs], dim='time_counter') ponlc_ts = xr.concat([ponlc_ts, ponlc], dim='time_counter') ponhs_ts = xr.concat([ponhs_ts, ponhs], dim='time_counter') dialc_ts = xr.concat([dialc_ts, dialc], dim='time_counter') diahs_ts = xr.concat([diahs_ts, diahs], dim='time_counter') flalc_ts = xr.concat([flalc_ts, flalc], dim='time_counter') diahs_ts = xr.concat([diahs_ts, diahs], dim='time_counter') # In[19]: chllc_ts = (dialc_ts + flalc_ts) * 2 chlhs_ts = (diahs_ts + flahs_ts) * 2 # In[20]: first = True i1, i2, j1, j2 = (95, 200, 68, 160) # Hood Canal all #i1, i2, j1, j2 = (95, 170, 76, 120) # Central Hood Canal #i1, i2, j1, j2 = (106, 119, 69, 79) # #i1lc, i2lc, j1lc, j2lc = (120, 150, 68, 75) # Lynch Cove i1hs, i2hs, j1hs, j2hs = (107, 110, 83, 86) # Hoodsport buoy point in Lynch Cove (Hood Canal) - 47°25'18.5"N 123°06'45.3"W i1lc, i2lc, j1lc, j2lc = (122, 125, 70, 73) # Twanoh buoy point in Lynch Cove (Hood Canal) - 47.37°N, 123.01°W for ii, day in enumerate(arrow.Arrow.range('day', startdate, enddate)): year = day.year yr2 = day.strftime("%y") month = day.month Month = day.strftime("%b").lower() day = day.day # set up filename to follow NEMO conventions filename_grid = f'/results2/SalishSea/nowcast-green.202111/{day:02}{Month}{yr2}/SalishSea_1d_{year}{month:02}{day:02}_{year}{month:02}{day:02}_grid_T.nc' with xr.open_dataset(filename_grid) as ds_grid: temphs = ds_grid['votemper'][0, :, j1hs:j2hs, i1hs:i2hs].where(tmask2_temp[0, :, j1hs:j2hs, i1hs:i2hs] == 1).mean(axis=1).mean(axis=1) templc = ds_grid['votemper'][0, :, j1lc:j2lc, i1lc:i2lc].where(tmask2_temp[0, :, j1lc:j2lc, i1lc:i2lc] == 1).mean(axis=1).mean(axis=1) salths = ds_grid['vosaline'][0, :, j1hs:j2hs, i1hs:i2hs].where(tmask2_salt[0, :, j1hs:j2hs, i1hs:i2hs] == 1).mean(axis=1).mean(axis=1) saltlc = ds_grid['vosaline'][0, :, j1lc:j2lc, i1lc:i2lc].where(tmask2_salt[0, :, j1lc:j2lc, i1lc:i2lc] == 1).mean(axis=1).mean(axis=1) #dens = ds_grid['sigma_theta'][0, :, j1:j2, i1:i2].where(tmask2_dens[0, :, j1:j2, i1:i2] == 1).mean(axis=1).mean(axis=1) #denslc = ds_grid['sigma_theta'][0, :, j1lc:j2lc, i1lc:i2lc].where(tmask2_dens[0, :, j1lc:j2lc, i1lc:i2lc] == 1).mean(axis=1).mean(axis=1) if first: first = False temphs_ts = temphs templc_ts = templc salths_ts = salths saltlc_ts = saltlc #dens_ts = dens #denslc_ts = denslc else: temphs_ts = xr.concat([temphs_ts, temphs], dim='time_counter') templc_ts = xr.concat([templc_ts, templc], dim='time_counter') salths_ts = xr.concat([salths_ts, salths], dim='time_counter') saltlc_ts = xr.concat([saltlc_ts, saltlc], dim='time_counter') #dens_ts = xr.concat([dens_ts, dens], dim='time_counter') #denslc_ts = xr.concat([denslc_ts, denslc], dim='time_counter') # In[21]: fig, ax = plt.subplots(5, 1, figsize=(15, 18)) templc_ts[:, :24].plot(ax=ax[0], y="deptht", yincrease=False, cmap='gist_rainbow_r', vmin=3, vmax=21, ylim=(30, 0)); saltlc_ts[:, :24].plot(ax=ax[1], y="deptht", yincrease=False, cmap= 'gist_rainbow_r', vmin=16, vmax=31, ylim=(30, 0)); oxylc_ts[:, :24].plot(ax=ax[2], y="deptht", yincrease=False, cmap='gist_rainbow_r', vmin=0, vmax=400, ylim=(30, 0)); chllc_ts[:, :24].plot(ax=ax[3], y="deptht", yincrease=False, cmap='gist_rainbow_r', vmin=0, vmax=25, ylim=(30, 0)); no3lc_ts[:, :24].plot(ax=ax[4], y="deptht", yincrease=False, cmap='gist_rainbow_r', vmin=0, vmax=50, ylim=(30, 0)); ax[0].set_ylabel('Depth (m)') ax[1].set_ylabel('Depth (m)') ax[2].set_ylabel('Depth (m)') ax[3].set_ylabel('Depth (m)') ax[0].set_xlabel('') ax[1].set_xlabel('') ax[2].set_xlabel('') ax[3].set_xlabel('') ax[4].set_xlabel('Hood Canal (Twanoh - Lynch Cove)') ax[0].set_title('Temperature (°C)') ax[1].set_title('Salinity (g/kg)') ax[2].set_title('Oxygen (µmol/L)') ax[3].set_title('Chl (mg/m3)') ax[4].set_title('Nitrate (µmol/L)') # Make ticks on occurrences of each month: ax[0].xaxis.set_major_locator(mdates.MonthLocator()) ax[1].xaxis.set_major_locator(mdates.MonthLocator()) ax[2].xaxis.set_major_locator(mdates.MonthLocator()) ax[3].xaxis.set_major_locator(mdates.MonthLocator()) ax[4].xaxis.set_major_locator(mdates.MonthLocator()) # Get only the month to show in the x-axis: ax[0].xaxis.set_major_formatter(mdates.DateFormatter('%b')) ax[1].xaxis.set_major_formatter(mdates.DateFormatter('%b')) ax[2].xaxis.set_major_formatter(mdates.DateFormatter('%b')) ax[3].xaxis.set_major_formatter(mdates.DateFormatter('%b')) ax[4].xaxis.set_major_formatter(mdates.DateFormatter('%b')) # In[22]: fig, ax = plt.subplots(5, 1, figsize=(15, 18)) temphs_ts[:, :24].plot(ax=ax[0], y="deptht", yincrease=False, cmap='gist_rainbow_r', vmin=3, vmax=21, ylim=(30, 0)); salths_ts[:, :24].plot(ax=ax[1], y="deptht", yincrease=False, cmap= 'gist_rainbow_r', vmin=16, vmax=31, ylim=(30, 0)); oxyhs_ts[:, :24].plot(ax=ax[2], y="deptht", yincrease=False, cmap='gist_rainbow_r', vmin=0, vmax=400, ylim=(30, 0)); chlhs_ts[:, :24].plot(ax=ax[3], y="deptht", yincrease=False, cmap='gist_rainbow_r', vmin=0, vmax=25, ylim=(30, 0)); no3hs_ts[:, :24].plot(ax=ax[4], y="deptht", yincrease=False, cmap='gist_rainbow_r', vmin=0, vmax=50, ylim=(30, 0)); ax[0].set_ylabel('Depth (m)') ax[1].set_ylabel('Depth (m)') ax[2].set_ylabel('Depth (m)') ax[3].set_ylabel('Depth (m)') ax[0].set_xlabel('') ax[1].set_xlabel('') ax[2].set_xlabel('') ax[3].set_xlabel('') ax[4].set_xlabel('Hood Canal (Hoodsport - Central Hood Canal)') ax[0].set_title('Temperature (°C)') ax[1].set_title('Salinity (g/kg)') ax[2].set_title('Oxygen (µmol/L)') ax[3].set_title('Chl (mg/m3)') ax[4].set_title('Nitrate (µmol/L)') # Make ticks on occurrences of each month: ax[0].xaxis.set_major_locator(mdates.MonthLocator()) ax[1].xaxis.set_major_locator(mdates.MonthLocator()) ax[2].xaxis.set_major_locator(mdates.MonthLocator()) ax[3].xaxis.set_major_locator(mdates.MonthLocator()) ax[4].xaxis.set_major_locator(mdates.MonthLocator()) # Get only the month to show in the x-axis: ax[0].xaxis.set_major_formatter(mdates.DateFormatter('%b')) ax[1].xaxis.set_major_formatter(mdates.DateFormatter('%b')) ax[2].xaxis.set_major_formatter(mdates.DateFormatter('%b')) ax[3].xaxis.set_major_formatter(mdates.DateFormatter('%b')) ax[4].xaxis.set_major_formatter(mdates.DateFormatter('%b')) # In[23]: fig, ax = plt.subplots(3, 1, figsize=(18, 10)) donlc_ts[:, :22].plot(ax=ax[0], y="deptht", yincrease=False, cmap= 'gist_rainbow_r', vmin=0, vmax=5, ylim=(30, 0)); ponlc_ts[:, :22].plot(ax=ax[1], y="deptht", yincrease=False, cmap='gist_rainbow_r', vmin=0, vmax=2, ylim=(30, 0)); oxylc_ts[:, :22].plot(ax=ax[2], y="deptht", yincrease=False, cmap='gist_rainbow_r', vmin=60, vmax=300, ylim=(30, 0)); ax[2].set_xlim(startdate, enddate) ax[0].set_ylabel('Depth (m)') ax[1].set_ylabel('Depth (m)') ax[2].set_ylabel('Depth (m)') ax[0].set_xlabel('') ax[1].set_xlabel('') ax[2].set_xlabel('Hood Canal (Twanoh - Lynch Cove)') ax[0].set_title('DON (µmol/L)') ax[1].set_title('PON (µmol/L)') ax[2].set_title('Oxygen (µmol/L)') # Make ticks on occurrences of each month: ax[0].xaxis.set_major_locator(mdates.MonthLocator()) ax[1].xaxis.set_major_locator(mdates.MonthLocator()) ax[2].xaxis.set_major_locator(mdates.MonthLocator()) # Get only the month to show in the x-axis: ax[0].xaxis.set_major_formatter(mdates.DateFormatter('%b')) ax[1].xaxis.set_major_formatter(mdates.DateFormatter('%b')) ax[2].xaxis.set_major_formatter(mdates.DateFormatter('%b')) # In[24]: fig, ax = plt.subplots(3, 1, figsize=(18, 10)) donhs_ts[:, :22].plot(ax=ax[0], y="deptht", yincrease=False, cmap= 'gist_rainbow_r', vmin=0, vmax=5, ylim=(30, 0)); ponhs_ts[:, :22].plot(ax=ax[1], y="deptht", yincrease=False, cmap='gist_rainbow_r', vmin=0, vmax=2, ylim=(30, 0)); oxyhs_ts[:, :22].plot(ax=ax[2], y="deptht", yincrease=False, cmap='gist_rainbow_r', vmin=60, vmax=300, ylim=(30, 0)); ax[2].set_xlim(startdate, enddate) ax[0].set_ylabel('Depth (m)') ax[1].set_ylabel('Depth (m)') ax[2].set_ylabel('Depth (m)') ax[0].set_xlabel('') ax[1].set_xlabel('') ax[2].set_xlabel('Hood Canal (Twanoh - Lynch Cove)') ax[0].set_title('DON (µmol/L)') ax[1].set_title('PON (µmol/L)') ax[2].set_title('Oxygen (µmol/L)') # Make ticks on occurrences of each month: ax[0].xaxis.set_major_locator(mdates.MonthLocator()) ax[1].xaxis.set_major_locator(mdates.MonthLocator()) ax[2].xaxis.set_major_locator(mdates.MonthLocator()) # Get only the month to show in the x-axis: ax[0].xaxis.set_major_formatter(mdates.DateFormatter('%b')) ax[1].xaxis.set_major_formatter(mdates.DateFormatter('%b')) ax[2].xaxis.set_major_formatter(mdates.DateFormatter('%b')) # In[25]: fig, ax = plt.subplots(1, 2, figsize=(10, 9)) ax[0].plot(no3hs_ts[91:273, :22], oxyhs_ts[91:273, :22], '.'); ax[1].plot(no3lc_ts[91:273, :22], oxylc_ts[91:273, :22], '.'); ax[0].set_xlim(-1, 30) ax[1].set_xlim(-1, 30) ax[0].set_ylim(90, 380) ax[1].set_ylim(90, 380) ax[0].set_ylabel('Dissolved_oxygen') ax[0].set_xlabel('NO3') ax[1].set_xlabel('NO3') ax[0].text(10, 350, 'Hoodsport (0-30 m, Apr-Sep)') ax[1].text(10, 350, 'Lynch Cove (0-30 m, Apr-Sep)') # In[26]: fig, ax = plt.subplots(6, 1, figsize=(15, 14)) #specify x-axis locations #x_ticks = [1, 32, 60, 91, 121, 152, 182, 213, 244, 274, 305, 335] x_ticks = [1, 366, 732, 1097, 1462, 1827, 2193, 2558] #specify x-axis labels #x_labels = ['Jan', 'Feb', 'Mar', 'Apr', 'May', 'Jun', 'Jul', 'Aug', 'Sep', 'Oct', 'Nov', 'Dec'] x_labels = ['2007', '2008', '2009', '2010', '2011', '2012', '2013', '2014'] t1 = np.mean(templc_ts[:, :4], axis=1)-np.mean(templc_ts[:, :4], axis=(0, 1)) t2 = np.mean(templc_ts[:, 5:], axis=1)-np.mean(templc_ts[:, 5:], axis=(0, 1)) plt.setp(ax, xticks=x_ticks, xticklabels=x_labels) ax[0].plot(t1,'b') #ax[0].set_title('Top 5m') ax[0].set_ylabel('Temp (°C)') ax[0].plot(t2,'g') ax[0].grid() s1 = np.mean(saltlc_ts[:, :4], axis=1)-np.mean(saltlc_ts[:, :4], axis=(0, 1)) s2 = np.mean(saltlc_ts[:, 5:], axis=1)-np.mean(saltlc_ts[:, 5:], axis=(0, 1)) ax[1].plot(s1,'b') ax[1].set_ylabel('Salt (g/kg)') ax[1].plot(s2,'g') ax[1].grid() n1 = np.mean(no3lc_ts[:, :4], axis=1)-np.mean(no3lc_ts[:, :4], axis=(0, 1)) n2 = np.mean(no3lc_ts[:, 5:], axis=1)-np.mean(no3lc_ts[:, 5:], axis=(0, 1)) ax[2].plot(n1,'b') ax[2].set_ylabel('NO3 (µmol/L)') ax[2].plot(n2,'g') ax[2].grid() d1 = np.mean(donlc_ts[:, :4], axis=1)-np.mean(donlc_ts[:, :4], axis=(0, 1)) d2 = np.mean(donlc_ts[:, 5:], axis=1)-np.mean(donlc_ts[:, 5:], axis=(0, 1)) ax[3].plot(d1,'b') ax[3].set_ylabel('DON (µmol/L)') ax[3].plot(d2,'g') ax[3].grid() p1 = np.mean(ponlc_ts[:, :4], axis=1)-np.mean(ponlc_ts[:, :4], axis=(0, 1)) p2 = np.mean(ponlc_ts[:, 5:], axis=1)-np.mean(ponlc_ts[:, 5:], axis=(0, 1)) ax[4].plot(p1,'b') ax[4].set_ylabel('PON (µmol/L)') ax[4].plot(p2,'g') ax[4].grid() o1 = np.mean(oxylc_ts[:, :4], axis=1)-np.mean(oxylc_ts[:, :4], axis=(0, 1)) o2 = np.mean(oxylc_ts[:, 5:], axis=1)-np.mean(oxylc_ts[:, 5:], axis=(0, 1)) ax[5].plot(o1,'b') ax[5].set_ylabel('Oxy (µmol/L)') ax[5].plot(o2,'g') ax[5].grid() # In[27]: fig, ax = plt.subplots(6, 1, figsize=(15, 14)) #specify x-axis locations #x_ticks = [1, 32, 60, 91, 121, 152, 182, 213, 244, 274, 305, 335] x_ticks = [1, 366, 732, 1097, 1462, 1827, 2193, 2558] #specify x-axis labels #x_labels = ['Jan', 'Feb', 'Mar', 'Apr', 'May', 'Jun', 'Jul', 'Aug', 'Sep', 'Oct', 'Nov', 'Dec'] x_labels = ['2007', '2008', '2009', '2010', '2011', '2012', '2013', '2014'] t1 = np.mean(templc_ts[:, :4], axis=1)-np.mean(templc_ts[:, :4], axis=(0, 1)) t2 = np.mean(templc_ts[:, 5:], axis=1)-np.mean(templc_ts[:, 5:], axis=(0, 1)) plt.setp(ax, xticks=x_ticks, xticklabels=x_labels) ax[0].set_title('Below 5m') ax[0].set_ylabel('Temp (°C)') ax[0].plot(t2,'g') ax[0].grid() s1 = np.mean(saltlc_ts[:, :4], axis=1)-np.mean(saltlc_ts[:, :4], axis=(0, 1)) s2 = np.mean(saltlc_ts[:, 5:], axis=1)-np.mean(saltlc_ts[:, 5:], axis=(0, 1)) ax[1].set_ylabel('Salt (g/kg)') ax[1].plot(s2,'g') ax[1].grid() n1 = np.mean(no3lc_ts[:, :4], axis=1)-np.mean(no3lc_ts[:, :4], axis=(0, 1)) n2 = np.mean(no3lc_ts[:, 5:], axis=1)-np.mean(no3lc_ts[:, 5:], axis=(0, 1)) ax[2].set_ylabel('NO3 (µmol/L)') ax[2].plot(n2,'g') ax[2].grid() d1 = np.mean(donlc_ts[:, :4], axis=1)-np.mean(donlc_ts[:, :4], axis=(0, 1)) d2 = np.mean(donlc_ts[:, 5:], axis=1)-np.mean(donlc_ts[:, 5:], axis=(0, 1)) ax[3].set_ylabel('DON (µmol/L)') ax[3].plot(d2,'g') ax[3].grid() p1 = np.mean(ponlc_ts[:, :4], axis=1)-np.mean(ponlc_ts[:, :4], axis=(0, 1)) p2 = np.mean(ponlc_ts[:, 5:], axis=1)-np.mean(ponlc_ts[:, 5:], axis=(0, 1)) ax[4].set_ylabel('PON (µmol/L)') ax[4].plot(p2,'g') ax[4].grid() o1 = np.mean(oxylc_ts[:, :4], axis=1)-np.mean(oxylc_ts[:, :4], axis=(0, 1)) o2 = np.mean(oxylc_ts[:, 5:], axis=1)-np.mean(oxylc_ts[:, 5:], axis=(0, 1)) ax[5].set_ylabel('Oxy (µmol/L)') ax[5].plot(o2,'g') ax[5].grid() # In[28]: i1, i2, j1, j2 = (95, 200, 68, 160) # Hood Canal all fig, ax = plt.subplots(1, 1, figsize=(5, 9)) oxy[0, mbathy[0, :, :]-1, :, :].where(tmask[0, mbathy[0, :, :]-1, :, :] == 1).plot(ax=ax, cmap=cm.oxy, vmin=20, vmax=250) viz_tools.set_aspect(ax); # Oxygen in the Hood Canal (j=60:120, i=100:160) left, bottom, width, height = (95, 68, 105, 100) rect1=mpatches.Rectangle((left,bottom),width,height, fill=False, #alpha=0.1 color="purple", linewidth=2) #facecolor="red") left, bottom, width, height = (180, 330, 35, 60) rect2=mpatches.Rectangle((left,bottom),width,height, fill=False, #alpha=0.1 color="purple", linewidth=2) #facecolor="red") plt.gca().add_patch(rect1) plt.gca().add_patch(rect2) fig.suptitle('Bottom Oxygen - 202111-nowcast - Year 2011'); #120, 150, 68, 75 # Lynch Cove mooring #107, 110, 83, 86 # Hoodsport mooring ax.set_ylim(75, 100) ax.set_xlim(100, 120) # In[ ]: # In[ ]: