diff --git a/pyglider/plotting.py b/pyglider/plotting.py index 5f757e7..65de7e0 100644 --- a/pyglider/plotting.py +++ b/pyglider/plotting.py @@ -50,24 +50,54 @@ def timeseries_plots(fname, plottingyaml): except: pass - with xr.open_dataset(fname) as ds0: + #with xr.open_dataset(fname) as ds0: #jpnote commented + with xr.open_dataset(fname, decode_times=True) as ds0: + ds = ds0.sel(time=slice(starttime, None)) # map! - fig, axs = plt.subplots(3, 1, gridspec_kw={'height_ratios': [3, 1, 1]}) + # fig, axs = plt.subplots(3, 1, gridspec_kw={'height_ratios': [3, 1, 1]}) + fig1, axs1 = plt.subplots() + # ax = axs[0] + ax = axs1 + good = (ds.longitude < -125) ## jp added + ax.plot(ds.longitude[good], ds.latitude[good], '.',markersize=1) + + #ax.plot(ds.longitude, ds.latitude, '.') #jp commented + # ax.set_aspect(1 / np.cos(np.deg2rad(ds.latitude.mean()))) + fig1.set_size_inches(7.5, 7.5) + ax.set_ylabel('Lat [degrees north]') + ax.set_xlabel('Lon [degrees west]') + ax.grid() + fig1.savefig(config['figdir'] + '/jp1_map%s.png'%ds.attrs['deployment_name'], dpi=200) + + #locator = mdates.AutoDateLocator() + #formatter = mdates.ConciseDateFormatter(locator) + #formatter.offset_formats = ['', + # '', + # '',] + + #fig, axs = plt.subplots(3, 1, gridspec_kw={'height_ratios': [3, 1, 1]}) + fig, axs = plt.subplots(2, 1, gridspec_kw={'height_ratios': [1, 1]}) + # ax = axs[1] ax = axs[0] - ax.plot(ds.longitude, ds.latitude, '.') - ax.set_aspect(1 / np.cos(np.deg2rad(ds.latitude.mean()))) + ax.plot(ds.time[good], ds.longitude[good], '.',markersize=1) + # # ax.plot(ds.time, ds.longitude, '.') #jp commented + ax.set_ylabel('Lon [degrees west]') + ax.grid() + # ax.xaxis.set_major_locator(locator) + # ax.xaxis.set_major_formatter(formatter) + + # ax = axs[2] + ax = axs[1] + ax.plot(ds.time[good], ds.latitude[good], '.',markersize=1) + # # ax.plot(ds.time, ds.latitude, '.')#jp commented ax.set_ylabel('Lat [degrees north]') - ax.set_xlabel('Lon [degrees east]') + ax.grid() - ax = axs[1] - ax.plot(ds.time, ds.longitude, '.') - ax.set_ylabel('Lon [degrees east]') + # ax.xaxis.set_major_locator(locator) + # ax.xaxis.set_major_formatter(formatter) - ax = axs[2] - ax.plot(ds.time, ds.latitude, '.') - ax.set_ylabel('Lat [degrees north]') - fig.savefig(config['figdir'] + '/map%s.png'%ds.attrs['deployment_name'], dpi=200) + fig.savefig(config['figdir'] + '/jp2_map%s.png'%ds.attrs['deployment_name'], dpi=200) # timeseries of things.... _log.info('Plotting timeseries data') @@ -78,18 +108,125 @@ def timeseries_plots(fname, plottingyaml): sharex=True, sharey=False) axs = axs.flat for n, k in enumerate(keys): + a = [0,1,2,3] + b = [4,5] _log.debug(f'key {k}') + # print(f'key {k}') if config['timeseries'][k] == 'True': ax = axs[n] + #convert 9999 to nan + ds[k] = ds[k].where(ds[k] != 9999, np.nan) + #ValueError: zero-size array to reduction operation minimum which has no identity good = np.where(~np.isnan(ds[k]))[0] - pc = ax.plot(ds.time[good], ds[k][good], '.') - min, max = _autoclim(ds[k][good]) - ax.set_ylim(min, max) - ax.set_title(ds[k].attrs['long_name'] + ' [' + - ds[k].attrs['units'] + ']', loc='left', fontsize=9) + if good.size != 0: + + pc = ax.plot(ds.time[good], ds[k][good], '.',markersize=1) + min, max = _autoclim(ds[k][good]) + ax.set_ylim(min, max) + ax.grid() + + ax.set_ylabel(ds[k].attrs['long_name'] + ' [' + + ds[k].attrs['units'] + ']') + if n in a: + ax.xaxis.set_tick_params(labelbottom=True) + #ax.set_title(ds[k].attrs['long_name'] + ' [' + + # ds[k].attrs['units'] + ']', loc='left', fontsize=9) + # if n in b: + #jpnote: offset format fix -remove hanging date + # locator = mdates.AutoDateLocator() + # formatter = mdates.ConciseDateFormatter(locator) + # formatter.offset_formats = ['', + # '', + # '',] + # ax.xaxis.set_major_locator(locator) + # ax.xaxis.set_major_formatter(formatter) + #end offset format fix + ax.set_title(ds[k].attrs['long_name'] + ' over time', loc='left', fontsize=9) fig.savefig(config['figdir'] + '/ts_%s.png'%ds.attrs['deployment_name'], dpi=200) + _log.info('Plotting timeseries data') + keys = config['timeseries'].keys() + N = len(keys) + if 1: + fig, axs = plt.subplots(int(N / 2), 2, figsize=(7, 12), + sharex=False, sharey=False) + axs = axs.flat + for n, k in enumerate(keys): + _log.debug(f'key {k}') + print(f'key {k}') + if config['timeseries'][k] == 'True': + ax = axs[n] + #convert 9999 to nan + # ds[k] = ds[k].where(ds[k] != 9999, np.nan) + # ds[k] = ds[k].where(ds[k] != 0, np.nan) + #ValueError: zero-size array to reduction operation minimum which has no identity + # good = np.where(~np.isnan(ds[k]))[0] + + if ax == axs[1]: + ax.set_xlim(30,35) #jp hardcode? + # ds[k] = ds[k].where(ds[k] != 0, np.nan) + # good = np.where(~np.isnan(ds[k]))[0] + if ax == axs[3]: + ax.set_xlim(0,0.004) #jp hardcode + if ax == axs[4]: + ax.set_xlim(0,7.5) + if ax == axs[5]: + ax.set_xlim(0,4) + + good = np.where(~np.isnan(ds[k]))[0] + pc = ax.plot(ds[k][good],ds.depth[good], '.',markersize=1) + ax.grid() + ax.set_ylabel('DEPTH [m]') + ax.xaxis.set_tick_params(labelbottom=True) + ax.invert_yaxis() + + ax.set_xlabel(ds[k].attrs['long_name'] + ' [' + + ds[k].attrs['units'] + ']') + ax.set_title(ds[k].attrs['long_name'] + ' over depth', loc='left', fontsize=9) + + fig.savefig(config['figdir'] + '/vv_%s.png'%ds.attrs['deployment_name'], dpi=200) + + + # _log.info('Plotting timeseries data') + keys = config['timeseries'].keys() + # N = len(keys) + ##water temp v depth plot + if 1: + #fig, axs = plt.subplots(int(N / 2), 2,figsize=(7.5, 7), + # sharex=True, sharey=False) + fig, axs = plt.subplots(1, 1,figsize=(7.5, 7), + sharex=True, sharey=False) + # axs = axs.flat + for n, k in enumerate(keys): + if n == 0: + a = [0,1,2,3] + b = [4,5] + #_log.debug(f'key {k}') + #print(f'key {k}') + if config['timeseries'][k] == 'True': + ax = axs + + good = np.where(~np.isnan(ds[k]))[0] + pc = ax.scatter(ds.temperature[good], ds.depth[good], c=ds.longitude[good], cmap='gray', s=1) + + ax.grid() + ax.set_ylabel('DEPTH [m]') + x = np.random.randint(low=0, high=10, size=13) + plt.xticks(np.arange(4, len(x)+1, 1)) + + ax.set_title('water temp over depth') + + ax.set_xlabel(ds[k].attrs['long_name'] + ' [' + + ds[k].attrs['units'] + ']') + plt.gca().invert_yaxis() + + cbar = fig.colorbar(pc, ax=axs) + cbar.ax.set_yticklabels(["{:.4}".format(i) for i in cbar.get_ticks()]) + cbar.ax.set_ylabel('Longitude [degrees west]') + fig.savefig(config['figdir'] + '/jp_ts_%s.png'%ds.attrs['deployment_name'], dpi=200) + + # colorline? if 0: fig, axs = plt.subplots(int(N / 2), 2, figsize=(7.5, 7), @@ -156,6 +293,7 @@ def timeseries_plots(fname, plottingyaml): fig, axs = plt.subplots(1, 2, constrained_layout=True, sharey=True) ax = axs[0] + ds['oxygen_concentration'][ ds['oxygen_concentration']==9999]=np.nan good = ~np.isnan(ds.salinity) smin, smax = _autoclim(ds.salinity[good]) good = ~np.isnan(ds.temperature) @@ -175,13 +313,14 @@ def timeseries_plots(fname, plottingyaml): ax.set_ylabel('$T\ [^oC]$') ax.set_xlim(smin, smax) ax.set_ylim(tmin, tmax) - + ax.grid() try: ax = axs[1] ax.plot(ds['oxygen_concentration'], ds['temperature'], '.', markersize=1) ax.set_xlabel('$O^2\ [mmol/L]$') + ax.grid() except: pass add_suptitle(fig, ds) @@ -213,7 +352,8 @@ def grid_plots(fname, plottingyaml): except: pass - with xr.open_dataset(fname) as ds0: + #with xr.open_dataset(fname) as ds0: #jpnote changed + with xr.open_dataset(fname, decode_times=True) as ds0: ds = ds0.sel(time=slice(starttime, None)) _log.debug(str(ds)) @@ -223,7 +363,9 @@ def grid_plots(fname, plottingyaml): tmean = ds.temperature.mean(axis=1) indmax = np.where(~np.isnan(tmean))[0][-1] depmax = ds.depth[indmax] - fig, axs = plt.subplots(math.ceil(N/2), 2, figsize=(7.5, 7), + # fig, axs = plt.subplots(math.ceil(N/2), 2, figsize=(7.5, 7), #jpnote changed + # sharex=True, sharey=True) + fig, axs = plt.subplots(int(N / 2), 2, figsize=(7.5, 7), sharex=True, sharey=True) axs = axs.flat for n, k in enumerate(keys): @@ -231,13 +373,19 @@ def grid_plots(fname, plottingyaml): pconf = config['pcolor']['vars'][k] _log.debug(pconf) - cmap = pconf.get('cmap', 'viridis') + + cmap = pconf.get('cmap','viridis') vmin = pconf.get('vmin', None) vmax = pconf.get('vmax', None) ax = axs[n] - locator = mdates.AutoDateLocator(minticks=3, maxticks=7) + ax.yaxis.set_tick_params(labelbottom=True) + ax.xaxis.set_tick_params(labelbottom=True) + locator = mdates.AutoDateLocator() formatter = mdates.ConciseDateFormatter(locator) + formatter.offset_formats = ['', + '', + '',] ax.xaxis.set_major_locator(locator) ax.xaxis.set_major_formatter(formatter) @@ -247,13 +395,20 @@ def grid_plots(fname, plottingyaml): if vmax is not None: max = vmax # make time windows: + if ax == axs[0]: + vmax = 300 + ax.set_ylim([300, 0]) + + vmax = 300 # get good profiles. i.e. those w data ind = np.where(np.sum(np.isfinite(ds[k].values), axis=0)>10)[0] _log.debug(ind) _log.debug(len(ds.time)) if len(ind) > 1: time = ds.time[ind[1:]] + np.diff(ds.time[ind]) / 2 + # print('tim1',time) time = np.hstack((time[0] - (time[1]-time[0]) / 2, time)) + # print('tim2',time) depth = ds.depth[:-1] - np.diff(ds.depth) depth = np.hstack((depth, depth[-1] + np.diff(ds.depth)[-1])) pc = ax.pcolormesh(time, depth, ds[k][:, ind], @@ -263,9 +418,12 @@ def grid_plots(fname, plottingyaml): _log.debug(ds[k]) - fig.colorbar(pc, ax=ax, extend='both', shrink=0.6) + # fig.colorbar(pc, ax=ax, extend='both', shrink=0.6) #jpnote commented + cbar = fig.colorbar(pc, ax=ax, extend='both', shrink=0.6) + cbar.ax.set_ylabel(ds[k].attrs['long_name'] + ' [' + + ds[k].attrs['units'] + ']', fontsize=9) ax.set_title(ds[k].attrs['long_name'] + ' [' + - ds[k].attrs['units'] + ']', loc='left', fontsize=9) + ds[k].attrs['units'] + ']', loc='left', fontsize=9) #jpnote commented t0 = ds.time[0] t1 = ds.time[-1] ax.set_xlim([t0, t1]) diff --git a/pyglider/website.py b/pyglider/website.py deleted file mode 100644 index f91cb0e..0000000 --- a/pyglider/website.py +++ /dev/null @@ -1,229 +0,0 @@ -import yaml -import glob -import os -import xarray as xr -import numpy as np -from jinja2 import Environment, FileSystemLoader -import geojson -import datetime -import simplekml -import pyglider.utils as pygu - -import logging - -""" -Utilities to make smart directories for the websites. - -A driving script may look like - -``` -import pyglider.website as pyweb -import logging - -logging.basicConfig(level=logging.DEBUG) -# uses the templat in `.templates/deploymentsIndex.html` to render -# index.html in each of the subdirectories: - -pyweb.geojson_deployments('./') - -if 1: - pyweb.index_deployments('./dfo-walle652/') - pyweb.index_deployments('./dfo-bb046/') -``` - -where `dfo-walle652` will contain subdirectories, each with a deployment: - -``` -./dfo-walle652/dfo-walle652-20210903 -./dfo-walle652/dfo-walle652-20220101 -... -``` - -Output are index files in each subdirectory and a geojson of all -the deployments (if that is desired) - -``` -./deployments.geojson # from the first script.... -./dfo-walle652/index.html -./dfo-walle652/dfo-walle652-20210903/index.html -... - -Requires an html template. The ones used for cproof are in -../example_html_templates of this repo and you need -to supply - -`deploymentsIndex.html` -`deploymentsInfoNew.html` -and maybe `deploymentsInfo.html` - -These are example templates, and written in jinga templating language. In -particular, see the strings near the bottom with double braces: eg: -"{{ att['deployment_name'] }}" These are what get filled by calling -`pyweb.index_deployments`. - -Note there are some hardcoded things in this tool, so please read the source -code! - -""" - -_log = logging.getLogger(__name__) - -def index_deployments(dir, templatedir='./.templates/'): - """ - Get useful info from deployments under "dir" and add to an - index.html in "dir" - - The structure is meant to be simple: - - dir/deployment1/deployment.yml - dir/deployment2/deployment.yml - dir/deployment3/deployment.yml - - Makes a file `dir/index.html` that has table of the deployments. - """ - - file_loader = FileSystemLoader(templatedir) - env = Environment(loader=file_loader) - template = env.get_template('deploymentsIndex.html') - - - subdirs = glob.glob(dir + '/*') - atts = [] - for d in subdirs: - if os.path.isdir(d): - if 1: - _log.info(d) - nc = glob.glob(d+'/L0-timeseries/*.nc') - if len(nc) < 1: - nc = glob.glob(d+'/L1-timeseries/*.nc') - if len(nc) < 1: - continue - with xr.open_dataset(nc[0]) as ds: - att = ds.attrs - atts.append(att) - else: - pass - if len(atts) > 0: - output = template.render(atts=atts, - title=atts[-1]['glider_name'] + atts[-1]['glider_serial']) - with open(dir + '/index.html', 'w') as fout: - fout.write(output) - - # now make the individual deployment index pages... - templateOld = env.get_template('deploymentsInfo.html') - templateNew = env.get_template('deploymentsInfoNew.html') - subdirs = glob.glob(dir + '/*') - atts = [] - for d in subdirs: - if os.path.isdir(d): - _log.info(d) - if 1: - figs = glob.glob(d + '/figs/*.png') - for n, fig in enumerate(figs): - figs[n] = './figs/' + os.path.split(fig)[1] - nc = glob.glob(d+'/L0-timeseries/*.nc') - template = templateNew - - if len(nc) < 1: - # try old style - nc = glob.glob(d+'/L1-timeseries/*.nc') - template = templateOld - if len(nc) < 1: - continue - with xr.open_dataset(nc[0]) as ds: - att = ds.attrs - _log.debug(type(ds.keys())) - keys = [] - units = [] - for key in ds.keys(): - _log.debug(ds[key].attrs) - try: - unit = ds[key].attrs['units'] - except KeyError: - unit = 'no units' - keys.append(key + ' [' + unit +']') - depname = att['deployment_name'] - output = template.render(deploy_name=depname, - title=att['glider_name'] + att['glider_serial'], - figs=figs, att=att, keys=keys) - with open(d + '/index.html', 'w') as fout: - fout.write(output) - else: - pass - -def geojson_deployments(dir, outfile='cproof-deployments.geojson'): - props = ['deployment_start', 'deployment_end', 'platform_type', - 'glider_model', 'glider_name', 'glider_serial', - 'deployment_name', 'project', 'institution', 'comment'] - subdirs = glob.glob(dir + '/*') - features = [] - - kml = simplekml.Kml() - - np.random.seed(20190101) - _log.debug(f'subdirs, {subdirs}') - colornum = 0; - for d in subdirs: - _log.info(d) - if os.path.isdir(d): - subdirs2 = glob.glob(d + '/*') - for d2 in subdirs2: - _log.info(d2) - if os.path.isdir(d2): - try: - nc = glob.glob(d2+'/L0-gridfiles/*.nc') - if len(nc) < 1: - # old style - nc = glob.glob(d2+'/L2-gridfiles/*.nc') - - with xr.open_dataset(nc[0]) as ds: - _log.info(f'opened {nc[0]}') - att = ds.attrs - good = (ds.longitude < -125) - line = np.vstack((ds.longitude[good], ds.latitude[good])).T - ls = geojson.LineString(line.tolist()) - feat = geojson.Feature(geometry=ls) - for prop in props: - if prop in ds.attrs.keys(): - feat.properties[prop] = ds.attrs[prop] - else: - feat.properties[prop] = '' - - # get URL.... - feat.properties['url'] = ('' + - 'http://cproof.uvic.ca/gliderdata/deployments/' + - d2[2:]) - # get color: - cols = np.random.randint(0, 200, 3) - # cols = pygu.get_html_non_blue(colornum) - colornum += 1 - feat.properties['color'] = '#%02X%02X%02X' % (cols[0], cols[1], cols[2]) - if ds['time'][-1] > np.datetime64(datetime.datetime.now()) - np.timedelta64(2, 'D'): - feat.properties['active'] = True - else: - feat.properties['active'] = False - - features += [feat] - - # make the kml: - pnt = kml.newpoint(coords=[line[-1]]) - pnt.style.iconstyle.icon.href = 'http://cproof.uvic.ca/deployments/assets/images/slocum_glider.png' - coords = [] - for thelon, thelat in zip(ds.longitude.values, ds.latitude.values): - coords += [(thelon, thelat)] - pnt.timestamp.when = f'{ds.time.values[-1]}'[:-3] - ls = kml.newlinestring(coords=coords, - name=att['deployment_name'], - ) - ls.timespan.begin = f'{ds.time.values[0]}'[:-3] - ls.timespan.end = f'{ds.time.values[-1]}'[:-3] - ls.style.linestyle.color = 'ee' + '%02X%02X%02X' % (cols[2], cols[1], cols[0]) - ls.style.linestyle.width = 3; - kml.save(d2[2:]+'/'+att['deployment_name']+'.kml') - - except: - _log.info(f'Could not find grid file {d2}') - feature_collection = geojson.FeatureCollection(features) - with open(outfile, 'w') as fout: - s = geojson.dumps(feature_collection) - fout.write(s)