NCCOOS Trac Projects: Top | Web | Platforms | Processing | Viz | Sprints | Sandbox | (Wind)

root/proc2plot/trunk/proc2plot/meet_wq_plot.py

Revision 455 (checked in by haines, 13 years ago)

Fix SVN locks and update code.

  • Property svn:executable set to
Line 
1 #!/usr/bin/env /opt/env/haines/dataproc/bin/python
2 # Last modified:  Time-stamp: <2010-04-14 11:37:45 haines>
3 """meet_wq_plot"""
4
5 # plot each month
6
7 import os
8 # import sys, glob, re
9 import datetime, time, dateutil, dateutil.tz
10 import pycdf
11 import numpy
12
13 # sys.path.append('/opt/env/haines/dataproc/raw2proc')
14 # del(sys)
15
16 os.environ["MPLCONFIGDIR"]="/home/haines/.matplotlib/"
17
18 from pylab import figure, twinx, savefig, setp, getp, cm, colorbar
19 from matplotlib.dates import DayLocator, HourLocator, MinuteLocator, DateFormatter, date2num, num2date
20 import procutil
21
22 def timeseries(pi, si, yyyy_mm, plot_type='latest'):
23     """
24     """
25
26
27     print 'meet_wq_plot ...'
28
29     #
30
31     prev_month, this_month, next_month = procutil.find_months(yyyy_mm)
32     yyyy_mm_str = this_month.strftime('%Y_%m')
33
34     ncFile1='/seacoos/data/nccoos/level1/meet/wq/meet_wq_'+prev_month.strftime('%Y_%m')+'.nc'
35     ncFile2='/seacoos/data/nccoos/level1/meet/wq/meet_wq_'+this_month.strftime('%Y_%m')+'.nc'
36
37     have_ncFile1 = os.path.exists(ncFile1)
38     have_ncFile2 = os.path.exists(ncFile2)
39
40     print ' ... loading data for graph from ...'
41     print ' ... ... ' + ncFile1 + ' ... ' + str(have_ncFile1)
42     print ' ... ... ' + ncFile2 + ' ... ' + str(have_ncFile2)
43
44     # open netcdf data
45     if have_ncFile1 and have_ncFile2:
46         nc = pycdf.CDFMF((ncFile1, ncFile2))
47     elif not have_ncFile1 and have_ncFile2:
48         nc = pycdf.CDFMF((ncFile2,))
49     elif have_ncFile1 and not have_ncFile2:
50         nc = pycdf.CDFMF((ncFile1,))
51     else:
52         print ' ... both files do not exist -- NO DATA LOADED'
53         return
54
55     # ncvars = nc.variables()
56     # print ncvars
57     es = nc.var('time')[:]
58     units = nc.var('time').units
59     dt = [procutil.es2dt(e) for e in es]
60     # set timezone info to UTC (since data from level1 should be in UTC!!)
61     dt = [e.replace(tzinfo=dateutil.tz.tzutc()) for e in dt]
62     # return new datetime based on computer local
63     dt_local = [e.astimezone(dateutil.tz.tzlocal()) for e in dt]
64     dn = date2num(dt)
65
66     z = nc.var('z')[:]
67     wtemp = nc.var('wtemp')[:]
68     cond = nc.var('cond')[:]
69     turb = nc.var('turb')[:]
70     ph = nc.var('ph')[:]
71     do_mg = nc.var('do_mg')[:]
72     do_sat = nc.var('do_sat')[:]
73     batt = nc.var('battvolts')[:]
74     nc.close()
75
76     ncFile1='/seacoos/data/nccoos/level1/meet/flow/meet_flow_'+prev_month.strftime('%Y_%m')+'.nc'
77     ncFile2='/seacoos/data/nccoos/level1/meet/flow/meet_flow_'+this_month.strftime('%Y_%m')+'.nc'
78
79     have_ncFile1 = os.path.exists(ncFile1)
80     have_ncFile2 = os.path.exists(ncFile2)
81
82     print ' ... loading data for graph from ...'
83     print ' ... ... ' + ncFile1 + ' ... ' + str(have_ncFile1)
84     print ' ... ... ' + ncFile2 + ' ... ' + str(have_ncFile2)
85
86     # open netcdf data
87     if have_ncFile1 and have_ncFile2:
88         nc = pycdf.CDFMF((ncFile1, ncFile2))
89     elif not have_ncFile1 and have_ncFile2:
90         nc = pycdf.CDFMF((ncFile2,))
91     elif have_ncFile1 and not have_ncFile2:
92         nc = pycdf.CDFMF((ncFile1,))
93     else:
94         print ' ... both files do not exist -- NO DATA LOADED'
95         return
96
97     # ncvars = nc.variables()
98     # print ncvars
99     es2 = nc.var('time')[:]
100     units = nc.var('time').units
101     dt2 = [procutil.es2dt(e) for e in es2]
102     # set timezone info to UTC (since data from level1 should be in UTC!!)
103     dt2 = [e.replace(tzinfo=dateutil.tz.tzutc()) for e in dt2]
104     dn2 = date2num(dt2)
105
106     r2 = nc.var('rain')[:] # inches of rain in past 15 min
107     have_sontek = False # don't plot sontek
108     # have_sontek = 'sontek_wl' in ncvars.keys() or 'sontek_flow' in ncvars.keys()
109     pwl2 = nc.var('press_wl')[:] # feet
110     pfl2 = nc.var('press_flow')[:] # cfs
111
112     if have_sontek:
113         swl2 = nc.var('sontek_wl')[:] # feet
114         sfl2 = nc.var('sontek_flow')[:] # cfs
115
116     nc.close()
117
118     # last dt in data for labels
119     dtu = dt[-1]
120     dtl = dt_local[-1]
121
122     diff = abs(dtu - dtl)
123     if diff.days>0:
124         last_dt_str = dtu.strftime("%H:%M %Z on %b %d, %Y") + ' (' + dtl.strftime("%H:%M %Z, %b %d") + ')'
125     else:
126         last_dt_str = dtu.strftime("%H:%M %Z") + ' (' + dtl.strftime("%H:%M %Z") + ')' \
127                       + dtl.strftime(" on %b %d, %Y")
128
129     #######################################
130     # Plot month all wq (plus flow) the timeseries data
131     #######################################
132
133     fig = figure(figsize=(10, 9))
134     fig.subplots_adjust(left=0.1, bottom=0.1, right=0.9, top=0.9, wspace=0.1, hspace=0.1)
135
136     ax = fig.add_subplot(6,1,1)
137     axs = [ax]
138
139     (x, y) = procutil.addnan(dt, wtemp)
140     ibad = y <= -6999.
141     y[ibad] = numpy.nan
142
143     # ax.plot returns a list of lines, so unpack tuple
144     l1, = ax.plot_date(x, y, fmt='b-')
145     l1.set_label('Water Temperature')
146
147     ax.set_ylabel('TEMP\n (deg C)')
148     ax.set_xlim(date2num(this_month), date2num(next_month-datetime.timedelta(seconds=1)))
149     ax.xaxis.set_major_locator( DayLocator(range(2,32,2)) )
150     ax.xaxis.set_minor_locator( HourLocator(range(0,25,12)) )
151     ax.xaxis.set_major_formatter( DateFormatter('%m/%d') )
152     ax.set_xticklabels([])
153
154     # right-hand side scale
155     ax2 = twinx(ax)
156     ax2.yaxis.tick_right()
157     # convert (lhs) deg C to (rhs) deg F
158     f = [procutil.celsius2fahrenheit(val) for val in ax.get_ylim()]
159     ax2.set_ylim(f)
160     ax2.set_ylabel('(deg F)')
161     ax2.set_xlim(date2num(this_month), date2num(next_month-datetime.timedelta(seconds=1)))
162     ax2.xaxis.set_major_locator( DayLocator(range(2,32,2)) )
163     ax2.xaxis.set_minor_locator( HourLocator(range(0,25,12)) )
164     ax2.set_xticklabels([])
165
166     # this only moves the label not the tick labels
167     ax.xaxis.set_label_position('top')
168     ax.set_xlabel('MEET Water Quality -- ' + yyyy_mm_str)
169
170     # legend
171     ls1 = l1.get_label()
172     leg = ax.legend((l1,), (ls1,), loc='upper left')
173     ltext  = leg.get_texts()  # all the text.Text instance in the legend
174     llines = leg.get_lines()  # all the lines.Line2D instance in the legend
175     frame  = leg.get_frame()  # the patch.Rectangle instance surrounding the legend
176     frame.set_facecolor('0.80')      # set the frame face color to light gray
177     frame.set_alpha(0.5)             # set alpha low to see through
178     setp(ltext, fontsize='small')    # the legend text fontsize
179     setp(llines, linewidth=1.5)      # the legend linewidth
180     # leg.draw_frame(False)           # don't draw the legend frame
181
182     #######################################
183     #
184     ax = fig.add_subplot(6,1,2)
185     axs.append(ax)
186
187     # ax.plot returns a list of lines, so unpack tuple
188     (x, y) = procutil.addnan(dt, cond)
189     l1, = ax.plot_date(x, y, fmt='b-')
190     l1.set_label('Specific Conductivity')
191
192     ax.set_ylabel('Cond\n (mS/cm)')
193     # ax.set_ylim(0,10)
194     ax.set_xlim(date2num(this_month), date2num(next_month-datetime.timedelta(seconds=1)))
195     ax.xaxis.set_major_locator( DayLocator(range(2,32,2)) )
196     ax.xaxis.set_minor_locator( HourLocator(range(0,25,12)) )
197     ax.xaxis.set_major_formatter( DateFormatter('%m/%d') )
198     ax.set_xticklabels([])
199
200     # legend
201     ls1 = l1.get_label()
202     leg = ax.legend((l1,), (ls1,), loc='upper left')
203     ltext  = leg.get_texts()  # all the text.Text instance in the legend
204     llines = leg.get_lines()  # all the lines.Line2D instance in the legend
205     frame  = leg.get_frame()  # the patch.Rectangle instance surrounding the legend
206     frame.set_facecolor('0.80')      # set the frame face color to light gray
207     frame.set_alpha(0.5)             # set alpha low to see through
208     setp(ltext, fontsize='small')    # the legend text fontsize
209     setp(llines, linewidth=1.5)      # the legend linewidth
210     # leg.draw_frame(False)           # don't draw the legend frame
211
212     #######################################
213     #
214     ax = fig.add_subplot(6,1,3)
215     axs.append(ax)
216
217     # ax.plot returns a list of lines, so unpack tuple
218     (x, y) = procutil.addnan(dt, do_mg)
219     l1, = ax.plot_date(x, y, fmt='b-')
220     l1.set_label('Dissolved Oxygen (DO)')
221
222     ax.set_ylabel('DO \n(mg/l)')
223     # ax.set_ylim(20,120)
224     ax.set_xlim(date2num(this_month), date2num(next_month-datetime.timedelta(seconds=1)))
225     ax.xaxis.set_major_locator( DayLocator(range(2,32,2)) )
226     ax.xaxis.set_minor_locator( HourLocator(range(0,25,12)) )
227     ax.xaxis.set_major_formatter( DateFormatter('%m/%d') )
228     ax.set_xticklabels([])
229
230     # right-hand side scale
231     ax2 = twinx(ax)
232     ax2.yaxis.tick_right()
233     (x, y) = procutil.addnan(dt, do_sat)
234     l2, = ax2.plot_date(x, y, fmt='c-')
235     l2.set_label('Percent Air Saturation')   
236     # convert (lhs) deg C to (rhs) deg F
237     # ax2.set_ylim(0,200)
238     ax2.set_ylabel('DO (%)')
239
240     ax2.set_xlim(date2num(this_month), date2num(next_month-datetime.timedelta(seconds=1)))
241     ax2.xaxis.set_major_locator( DayLocator(range(2,32,2)) )
242     ax2.xaxis.set_minor_locator( HourLocator(range(0,25,12)) )
243     ax2.set_xticklabels([])
244
245     # legend
246     ls1 = l1.get_label()
247     ls2 = l2.get_label()
248     leg = ax.legend((l1,l2), (ls1,ls2), loc='upper left')
249     ltext  = leg.get_texts()  # all the text.Text instance in the legend
250     llines = leg.get_lines()  # all the lines.Line2D instance in the legend
251     frame  = leg.get_frame()  # the patch.Rectangle instance surrounding the legend
252     frame.set_facecolor('0.80')      # set the frame face color to light gray
253     frame.set_alpha(0.5)             # set alpha low to see through
254     setp(ltext, fontsize='small')    # the legend text fontsize
255     setp(llines, linewidth=1.5)      # the legend linewidth
256     # leg.draw_frame(False)           # don't draw the legend frame
257     #######################################
258     #
259     ax = fig.add_subplot(6,1,4)
260     axs.append(ax)
261
262     # ax.plot returns a list of lines, so unpack tuple
263     (x, y) = procutil.addnan(dt, ph)
264     l1, = ax.plot_date(x, y, fmt='b-')
265     l1.set_label('pH')
266
267     ax.set_ylabel('pH')
268     # ax.set_ylim(6, 12)
269     ax.set_xlim(date2num(this_month), date2num(next_month-datetime.timedelta(seconds=1)))
270     ax.xaxis.set_major_locator( DayLocator(range(2,32,2)) )
271     ax.xaxis.set_minor_locator( HourLocator(range(0,25,12)) )
272     ax.xaxis.set_major_formatter( DateFormatter('%m/%d') )
273     ax.set_xticklabels([])
274
275     # legend
276     ls1 = l1.get_label()
277     leg = ax.legend((l1,), (ls1,), loc='upper left')
278     ltext  = leg.get_texts()  # all the text.Text instance in the legend
279     llines = leg.get_lines()  # all the lines.Line2D instance in the legend
280     frame  = leg.get_frame()  # the patch.Rectangle instance surrounding the legend
281     frame.set_facecolor('0.80')      # set the frame face color to light gray
282     frame.set_alpha(0.5)             # set alpha low to see through
283     setp(ltext, fontsize='small')    # the legend text fontsize
284     setp(llines, linewidth=1.5)      # the legend linewidth
285     # leg.draw_frame(False)           # don't draw the legend frame
286
287     #######################################
288     #
289     ax = fig.add_subplot(6,1,5)
290     axs.append(ax)
291
292     # ax.plot returns a list of lines, so unpack tuple
293     (x, y) = procutil.addnan(dt, turb)
294     l1, = ax.plot_date(x, y, fmt='b-')
295     l1.set_label('Turbidity')
296
297     ax.set_ylabel('Turb \n(NTU)')
298     # ax.set_ylim(0,100)
299     ax.set_xlim(date2num(this_month), date2num(next_month-datetime.timedelta(seconds=1)))
300     ax.xaxis.set_major_locator( DayLocator(range(2,32,2)) )
301     ax.xaxis.set_minor_locator( HourLocator(range(0,25,12)) )
302     ax.xaxis.set_major_formatter( DateFormatter('%m/%d') )
303     ax.set_xticklabels([])
304
305     # legend
306     ls1 = l1.get_label()
307     leg = ax.legend((l1,), (ls1,), loc='upper left')
308     ltext  = leg.get_texts()  # all the text.Text instance in the legend
309     llines = leg.get_lines()  # all the lines.Line2D instance in the legend
310     frame  = leg.get_frame()  # the patch.Rectangle instance surrounding the legend
311     frame.set_facecolor('0.80')      # set the frame face color to light gray
312     frame.set_alpha(0.5)             # set alpha low to see through
313     setp(ltext, fontsize='small')    # the legend text fontsize
314     setp(llines, linewidth=1.5)      # the legend linewidth
315     # leg.draw_frame(False)           # don't draw the legend frame
316
317     #######################################
318     #
319     ax = fig.add_subplot(6,1,6)
320     axs.append(ax)
321
322     # ax.plot returns a list of lines, so unpack tuple
323     (x, y) = procutil.addnan(dt2, pfl2)
324     l1, = ax.plot_date(x, y, fmt='b-')
325     l1.set_label('Derived from pressure')
326
327     ax.set_ylim(0,350)
328     ax.set_ylabel('Flow\n (cfs)')
329     ax.set_xlim(date2num(this_month), date2num(next_month-datetime.timedelta(seconds=1)))
330     ax.plot_date((x[1],x[-1]), (0,0), fmt='k:')
331     ax.xaxis.set_major_locator( DayLocator(range(2,32,2)) )
332     ax.xaxis.set_minor_locator( HourLocator(range(0,25,12)) )
333
334     # right-hand side scale
335     ax2 = twinx(ax)
336     ax2.yaxis.tick_right()
337     (x, y) = procutil.addnan(dt2, r2)
338     l2, = ax2.plot_date(x, y, fmt='c-')
339     l2.set_label('Precipitation')   
340
341     # ax2.set_ylim(-0.5,0.5)
342     ax2.set_ylabel('Rain\n (in)')
343     ax2.set_xlim(date2num(this_month), date2num(next_month-datetime.timedelta(seconds=1)))
344     ax2.xaxis.set_major_locator( DayLocator(range(2,32,2)) )
345     ax2.xaxis.set_minor_locator( HourLocator(range(0,25,12)) )
346     ax2.set_xticklabels([])
347
348     ax.set_xlim(date2num(this_month), date2num(next_month-datetime.timedelta(seconds=1)))
349     ax.xaxis.set_major_formatter( DateFormatter('%m/%d') )
350     ax.set_xlabel('MEET Water Quality -- ' + yyyy_mm_str)
351
352     # legend
353     ls1 = l1.get_label()
354     ls2 = l2.get_label()
355     leg = ax.legend((l1,l2), (ls1,ls2), loc='upper left')
356     ltext  = leg.get_texts()  # all the text.Text instance in the legend
357     llines = leg.get_lines()  # all the lines.Line2D instance in the legend
358     frame  = leg.get_frame()  # the patch.Rectangle instance surrounding the legend
359     frame.set_facecolor('0.80')      # set the frame face color to light gray
360     frame.set_alpha(0.5)             # set alpha low to see through
361     setp(ltext, fontsize='small')    # the legend text fontsize
362     setp(llines, linewidth=1.5)      # the legend linewidth
363     # leg.draw_frame(False)           # don't draw the legend frame
364
365     # save figure for this month
366     ofn = '/home/haines/rayleigh/img/meet/meet_wq_'+yyyy_mm_str+'.png'
367     print '... ... write: %s' % (ofn,)
368     savefig(ofn)
369
370
371     #######################################
372     # Last 30 days
373     #######################################
374     if plot_type=='latest':
375         print ' ... Last 30 days'
376         for idx, ax in enumerate(axs):
377             ax.set_xlim(date2num(dt[-1])-30, date2num(dt[-1]))
378             ax.xaxis.set_major_locator( DayLocator(range(2,32,2)) )
379             ax.xaxis.set_minor_locator( HourLocator(range(0,25,12)) )
380             ax.set_xticklabels([])
381             if idx==0:
382                 ax.set_xlabel('MEET Water Quality -- Last 30 days from ' + last_dt_str)
383             elif idx==len(axs)-1:
384                 ax.xaxis.set_major_formatter( DateFormatter('%m/%d') )
385                 ax.set_xlabel('MEET Water Quality -- Last 30 days from ' + last_dt_str)
386
387         savefig('/home/haines/rayleigh/img/meet_wq_last30days.png')
388
389
390     #######################################
391     # Last 7 days
392     #######################################
393     if plot_type=='latest':
394         print ' ... Last 7 days'
395         for idx, ax in enumerate(axs):
396             ax.set_xlim(date2num(dt[-1])-7, date2num(dt[-1]))
397             ax.xaxis.set_major_locator( DayLocator(range(0,32,1)) )
398             ax.xaxis.set_minor_locator( HourLocator(range(0,25,6)) )
399             ax.set_xticklabels([])
400             if idx==0:
401                 ax.set_xlabel('MEET Water Quality -- Last 7 days from ' + last_dt_str)
402             elif idx==len(axs)-1:
403                 ax.xaxis.set_major_formatter( DateFormatter('%m/%d') )
404                 ax.set_xlabel('MEET Water Quality -- Last 7 days from ' + last_dt_str)
405                
406         savefig('/home/haines/rayleigh/img/meet_wq_last07days.png')
407
408
409     #######################################
410     # Last 1 day (24hrs)
411     #######################################
412     if plot_type=='latest':
413         print ' ... Last 1 days'
414
415         for idx, ax in enumerate(axs):
416             ax.set_xlim(date2num(dt[-1])-1, date2num(dt[-1]))
417             ax.xaxis.set_major_locator( HourLocator(range(0,25,1)) )
418             ax.xaxis.set_minor_locator( MinuteLocator(range(0,61,30)) )
419             ax.set_xticklabels([])
420             if idx==0:
421                 ax.set_xlabel('MEET Water Quality -- Last 24 hours from ' + last_dt_str)
422             elif idx==len(axs)-1:
423                 ax.xaxis.set_major_formatter( DateFormatter('%H') )
424                 ax.set_xlabel('MEET Water Quality -- Last 24 hours from ' + last_dt_str)
425
426         savefig('/home/haines/rayleigh/img/meet_wq_last01days.png')
427
428
429
430
Note: See TracBrowser for help on using the browser.