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

root/sodarplot/trunk/sodarplot/scintec/winddist.py

Revision 80 (checked in by cbc, 13 years ago)

Refactor meanspeed.py and make most recent plot the subject of index.html for both meanspeed.py and winddist.py.

Line 
1 import os,datetime,glob
2 import multiprocessing
3 import numpy as np
4 import matplotlib as mpl
5 mpl.use("Agg")
6 from windrose import WindroseAxes
7 from matplotlib import pyplot as plt
8 import pycdf
9
10 ncDir = '/seacoos/data/nccoos/level1/billymitchell/sodar1'
11 ncFileGlob = 'billymitchell_sfas_*.nc'
12  
13 pngDir = '/seacoos/data/nccoos/level1/billymitchell/sodar1/plots/winddist'
14 if not os.path.exists(pngDir):
15     os.mkdir(pngDir)
16
17 pngExt = 'png'
18 htmlExt = 'html'
19  
20 ncFilePattern = os.path.join(ncDir,ncFileGlob)
21 files = glob.glob(ncFilePattern)
22 previous = [files[-1]] + files[:-1]
23 next = files[1:] + [files[0]]
24 files = zip(previous,files,next)
25
26 html1 = """<html>
27     <head>
28         <title>
29             Billy Mitchell Sodar :: """
30
31
32 html2 = """ :: Wind Distribution at """
33
34 html3 = """ Meters Above Sea Level
35         </title>
36     </head>
37     <body>
38         <form>
39             <a href=\""""
40
41 html4 = """\">&lt;Previous Month</a>&nbsp;&nbsp;&nbsp;
42             <a href=\""""
43
44 html5 = """\">Next Month&gt;</a>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;
45             <a href=\""""
46
47 html6 = """\">&lt;Previous Elevation</a>&nbsp;&nbsp;&nbsp;
48             <a href=\""""
49
50 html7 = """\">Next Elevation&gt;</a>
51         </form>
52         <br/>
53         <img src=\""""
54
55 html8 = """\"/>
56     </body>
57 </html>"""
58
59 def _components(fileSpec):
60     "Figure out year and month from file specification."
61
62     fileName = os.path.splitext(os.path.basename(fileSpec))[0]
63     year = fileName[-7:-3]
64     month = fileName[-2:]
65     monthName = datetime.datetime(int(year),int(month),1).strftime('%B')
66     return (fileName, year, month, monthName)
67  
68 def winddist((previous, ncFile, next), pngDir, genHtml=False, lastMonth=False):
69     print 'Processing',ncFile
70     ncFileName, year, month, monthName = _components(ncFile)
71     previous, previousYear, previousMonth, previousMonthName = _components(previous)
72     next, nextYear, nextMonth, nextMonthName = _components(next)
73
74     nc = pycdf.CDF(ncFile)
75
76     t = nc.var('time')[:]
77     z = nc.var('z')[:]
78     u = nc.var('u')[:]
79     v = nc.var('v')[:]
80     w = nc.var('w')[:]
81  
82     mu = np.ma.masked_invalid(u)
83     mv = np.ma.masked_invalid(v)
84     mw = np.ma.masked_invalid(w)
85  
86     mumean = mu.mean(axis=0)
87     mvmean = mv.mean(axis=0)
88     mwmean = mw.mean(axis=0)
89  
90     ucount = mu.count(axis=0)
91     vcount = mv.count(axis=0)
92     wcount = mw.count(axis=0)
93
94     mrho = np.sqrt((mu**2) + (mv**2))
95     mrho = mrho.T
96     mtheta = (180 * np.arctan2(mv, mu)/np.pi)
97     theta = np.piecewise(mtheta,
98                          (mtheta <= 90, mtheta > 90),
99                          (lambda x: 90 - x, lambda x: 450 - x))
100     mtheta = np.ma.array(theta,mask=mtheta.mask)
101     mtheta = mtheta.T
102
103     rho = [np.array([elem
104                      for elem,mask
105                      in zip(row.data,row.mask)
106                      if not mask])
107            for row
108            in mrho]
109     theta = [np.array([elem
110                        for elem,mask
111                        in zip(row.data,row.mask)
112                        if not mask])
113              for row
114              in mtheta]
115
116     firstElevation = True
117     indices = range(len(rho))
118     previousX = [indices[-1]] + indices[:-1]
119     nextX = indices[1:] + [indices[0]]
120     indices = zip(previousX,indices,nextX)
121     for previousX,x,nextX in indices:
122         if rho[x].any():
123             fig = plt.figure(figsize=(8, 8), dpi=80, facecolor='w', edgecolor='w')
124             spt = fig.suptitle('Billy Mitchell Sodar :: ' + monthName + " " + year)
125             rect = [0.16, 0.16, 0.68, 0.68]
126             ax = WindroseAxes(fig, rect, axisbg='w')
127             ax.set_title('Wind Distribution at %d Meters Above Sea Level\nPercentage of Wind Magnitude From Direction\n\n' % z[x],fontsize=spt.get_fontsize()*0.8)
128             fig.add_axes(ax)
129             ax.bar(theta[x], rho[x], bins = np.arange(0.0,20.0,2.5), normed=True, opening=0.8, edgecolor='white')
130             yticklabels = ['',] + \
131                           ["  %s%%" % plt.getp(textobj,'text')
132                                for textobj in plt.getp(ax,'yticklabels')[1:]]
133             plt.setp(ax,yticklabels=yticklabels)
134             l = ax.legend(axespad=-0.20,title="Magnitude (m/s)")
135             plt.setp(l.get_texts(), fontsize=8)
136             dirName = "%4u_%02u" % (int(year), int(month))
137             if not os.path.exists(os.path.join(pngDir, dirName)):
138                 os.mkdir(os.path.join(pngDir, dirName))
139             outFile = os.path.join(pngDir,dirName,('%03um' % z[x]) + os.extsep + pngExt)
140             print 'Saving', outFile
141             fig.savefig(outFile)
142             fig.clear()
143
144             if genHtml:
145                 htmlFile = os.path.join(pngDir,ncFileName + "_" + ('%dm' % z[x]) + os.extsep + htmlExt)
146                 html = html1 + monthName + " " + year
147                 html = html + html2 + ('%dm' % z[x])
148                 html = html + html3 + previous + "_" + ('%dm' % z[x]) + os.extsep + htmlExt
149                 html = html + html4 + next + "_" + ('%dm' % z[x]) + os.extsep + htmlExt
150                 html = html + html5 + ncFileName + "_" + ('%dm' % z[previousX]) + os.extsep + htmlExt
151                 html = html + html6 + ncFileName + "_" + ('%dm' % z[nextX]) + os.extsep + htmlExt
152                 html = html + html7 + os.path.join(dirName,os.path.basename(outFile)) + html8
153                 handle = open(htmlFile,'w')
154                 handle.write(html)
155                 handle.close()
156
157                 if lastMonth and firstElevation:
158                     firstElevation = False
159                     htmlFile = os.path.join(pngDir,"index" + os.extsep + htmlExt)
160                     handle = open(htmlFile,'w')
161                     handle.write(html)
162                     handle.close()
163
164 if __name__ == "__main__":
165     lastMonth = False
166     for previous,ncFile,next in files:
167         if next <= ncFile:
168             lastMonth = True
169         p = multiprocessing.Process(target=winddist, args=((previous, ncFile, next), pngDir, True, lastMonth))
170         p.start()
171         p.join()
172         if p.exitcode:
173            print "Exitcode %s from processing %s" % (p.exitcode, ncFile)
Note: See TracBrowser for help on using the browser.