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

root/proc2plot/trunk/proc2plot/crow_wq_plot.py

Revision 329 (checked in by haines, 14 years ago)

import first verison proc2plot

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