1 |
#!/usr/bin/env python |
---|
2 |
# Last modified: Time-stamp: <2009-10-28 17:27:56 haines> |
---|
3 |
""" |
---|
4 |
how to parse data, and assert what data and info goes into |
---|
5 |
creating and updating monthly netcdf files |
---|
6 |
|
---|
7 |
Nortek/AWAC processed adcp 2-D power spectrum (wds) as function of |
---|
8 |
frequency and direction |
---|
9 |
|
---|
10 |
parser : sample date and time and pressure from .wap, |
---|
11 |
energy spectrum m^2/Hz from .was, |
---|
12 |
normalized energy/deg from .wds |
---|
13 |
|
---|
14 |
based on George Voulgaris' matlab script (version 8, Feb 14, 2005, |
---|
15 |
polar_waves_cur_rdi.m) and additional parameters. |
---|
16 |
creator : lat, lon, z, time, freq, dir, Sxx(time, freq, dir), Sf(time, freq), |
---|
17 |
Stheta(time, dir), Stheta_swell(time, dir), Stheta_wind(time, dir), |
---|
18 |
Hs, Hs_swell, Hs_wind, |
---|
19 |
Tp, Tp_swell, Tp_wind, Tm, Tm_swell, Tm_wind, |
---|
20 |
Dp, Dp_swell, Dp_wind, Dm, Dm_swell, Dm_wind, |
---|
21 |
|
---|
22 |
updater : time, Sxx(time, freq, dir), Sf(time, freq), |
---|
23 |
Stheta(time, dir), Stheta_swell(time, dir), Stheta_wind(time, dir), |
---|
24 |
Hs, Hs_swell, Hs_wind, |
---|
25 |
Tp, Tp_swell, Tp_wind, Tm, Tm_swell, Tm_wind, |
---|
26 |
Dp, Dp_swell, Dp_wind, Dm, Dm_swell, Dm_wind, |
---|
27 |
|
---|
28 |
check that freq and dir have not changed from what is in current |
---|
29 |
NetCDF file |
---|
30 |
|
---|
31 |
Examples |
---|
32 |
-------- |
---|
33 |
|
---|
34 |
>> (parse, create, update) = load_processors(module_name_without_dot_py) |
---|
35 |
For example, |
---|
36 |
>> (parse, create, update) = load_processors('proc_rdi_logdata_adcp') |
---|
37 |
or |
---|
38 |
>> si = get_config(cn+'.sensor_info') |
---|
39 |
>> (parse, create, update) = load_processors(si['adcp']['proc_module']) |
---|
40 |
|
---|
41 |
Then use the generic name of processor to parse data, create or update |
---|
42 |
monthly output file |
---|
43 |
|
---|
44 |
>> lines = load_data(filename) |
---|
45 |
>> data = parse(platform_info, sensor_info, lines) |
---|
46 |
>> create(platform_info, sensor_info, data) |
---|
47 |
or |
---|
48 |
>> update(platform_info, sensor_info, data) |
---|
49 |
|
---|
50 |
""" |
---|
51 |
|
---|
52 |
from raw2proc import * |
---|
53 |
from procutil import * |
---|
54 |
from ncutil import * |
---|
55 |
|
---|
56 |
now_dt = datetime.utcnow() |
---|
57 |
now_dt.replace(microsecond=0) |
---|
58 |
|
---|
59 |
def parser(platform_info, sensor_info, lines): |
---|
60 |
""" |
---|
61 |
parse and assign wave spectra data from RDI ADCP Dspec |
---|
62 |
and compute wave statistics and parameters |
---|
63 |
|
---|
64 |
Notes |
---|
65 |
----- |
---|
66 |
1. adapted from polar_waves_cur_rdi.m (Version 8 - February 14, 2005) |
---|
67 |
by George Voulgaris |
---|
68 |
Coastal Processes & Sediment Dynamics Lab |
---|
69 |
Department of Geological Sciences |
---|
70 |
University of South Carolina, Columbia, SC 29205 |
---|
71 |
Email: gvoulgaris@geol.sc.edu |
---|
72 |
2. This parser requires date/time be parsed from .wap for each |
---|
73 |
spectum sample in .wds, strip .wds on input filename and load and |
---|
74 |
parse .wap here. If .wap with same name not available, then use sample |
---|
75 |
per hour starting at time parsed from filename. |
---|
76 |
|
---|
77 |
3. The .wds contains several bursts of full directional wave |
---|
78 |
spectrum. One for each hour in .wap. Each directional burst is |
---|
79 |
formatted in .wds file as each row is a one frequency, default |
---|
80 |
[0.02:0.01:0.99] or [0.02:0.01:0.49]. Each column is descretized |
---|
81 |
by 4 degrees 0:4:356 as degrees |
---|
82 |
|
---|
83 |
(faxed doc also states that freq's could be reported to .was file, |
---|
84 |
but I didn't find this so to be true) |
---|
85 |
|
---|
86 |
""" |
---|
87 |
|
---|
88 |
import numpy |
---|
89 |
from datetime import datetime |
---|
90 |
from time import strptime |
---|
91 |
|
---|
92 |
# get sample datetime from filename |
---|
93 |
fn = sensor_info['fn'] |
---|
94 |
sample_dt_start = filt_datetime(fn)[0] |
---|
95 |
|
---|
96 |
# try getting sample date/times from .wap |
---|
97 |
wap_fn = os.path.splitext(fn)[0] + ".wap" |
---|
98 |
if os.path.exists(wap_fn): |
---|
99 |
wap_lines = load_data(wap_fn) |
---|
100 |
|
---|
101 |
data = { |
---|
102 |
'dt' : numpy.array(numpy.ones((len(wap_lines),), dtype=object)*numpy.nan), |
---|
103 |
'time' : numpy.array(numpy.ones((len(wap_lines),), dtype=long)*numpy.nan), |
---|
104 |
'press' : numpy.array(numpy.ones((len(wap_lines),), dtype=float)*numpy.nan), |
---|
105 |
} |
---|
106 |
i=0 |
---|
107 |
|
---|
108 |
for line in wap_lines: |
---|
109 |
# split line and parse float and integers |
---|
110 |
wap = [] |
---|
111 |
sw = re.split(' ', line) |
---|
112 |
for s in sw: |
---|
113 |
m = re.search(REAL_RE_STR, s) |
---|
114 |
if m: |
---|
115 |
wap.append(float(m.groups()[0])) |
---|
116 |
|
---|
117 |
# get sample datetime from data |
---|
118 |
sample_str = '%02d-%02d-%4d %02d:%02d:%02d' % tuple(wap[0:6]) |
---|
119 |
if sensor_info['utc_offset']: |
---|
120 |
sample_dt = scanf_datetime(sample_str, fmt='%m-%d-%Y %H:%M:%S') + \ |
---|
121 |
timedelta(hours=sensor_info['utc_offset']) |
---|
122 |
else: |
---|
123 |
sample_dt = scanf_datetime(sample_str, fmt='%m-%d-%Y %H:%M:%S') |
---|
124 |
|
---|
125 |
# these items can also be teased out of raw adcp but for now get from config file |
---|
126 |
# th = sensor_info['transducer_ht'] # Transducer height above bottom (meters) |
---|
127 |
|
---|
128 |
# pressure (dbar) converted to water depth |
---|
129 |
pressure = wap[17] # pressure (dbar) at tranducer height (?) |
---|
130 |
# water_depth = th + sw_dpth(pressure, lat) |
---|
131 |
|
---|
132 |
data['dt'][i] = sample_dt |
---|
133 |
data['time'][i] = dt2es(sample_dt) |
---|
134 |
data['press'][i] = pressure # dbar |
---|
135 |
i=i+1 |
---|
136 |
|
---|
137 |
else: |
---|
138 |
print "error: No corresponding .wap file" |
---|
139 |
print " .... skipping %s" % (fn,) |
---|
140 |
return data |
---|
141 |
|
---|
142 |
# assign specific fields |
---|
143 |
nbursts = len(data['dt']) |
---|
144 |
Df = 0.01 # (Hz) |
---|
145 |
f = numpy.arange(0.02, 0.99, Df) |
---|
146 |
nfreq = len(f) # Number of frequencies (no units) |
---|
147 |
|
---|
148 |
# Did we get the number of data rows that we expected? Should equal nfreq |
---|
149 |
n = int(len(lines)/nfreq) |
---|
150 |
if n != nbursts: |
---|
151 |
print "Number of data rows %d does not match expected number %d" % (n, nbursts) |
---|
152 |
print " .... skipping %s" % (fn,) |
---|
153 |
return data |
---|
154 |
|
---|
155 |
Dtheta = 1.0 # degrees |
---|
156 |
D = numpy.arange(0.0, 360.0, 4) |
---|
157 |
D = numpy.mod(D,360) |
---|
158 |
ndir = len(D) # Number of directions (no units) |
---|
159 |
|
---|
160 |
# now get power spectra from .was |
---|
161 |
was_fn = os.path.splitext(fn)[0] + ".was" |
---|
162 |
was_Sf = numpy.array(numpy.ones((nbursts,nfreq), dtype=float)*numpy.nan) |
---|
163 |
if os.path.exists(was_fn): |
---|
164 |
was_lines = load_data(was_fn) |
---|
165 |
|
---|
166 |
i=0 |
---|
167 |
# first line is freq label for each column start at [1] |
---|
168 |
for line in was_lines[1:]: |
---|
169 |
# split line and parse float and integers |
---|
170 |
was = [] |
---|
171 |
sw = re.split(' ', line) |
---|
172 |
for s in sw: |
---|
173 |
m = re.search(REAL_RE_STR, s) |
---|
174 |
if m: |
---|
175 |
was.append(float(m.groups()[0])) |
---|
176 |
|
---|
177 |
# just the frequencies we have in directional spectra [0:nfreq] |
---|
178 |
was_Sf[i] = was[0:nfreq] # (m^2/Hz) non-directional power spectrum for each sample time |
---|
179 |
|
---|
180 |
i=i+1 |
---|
181 |
|
---|
182 |
else: |
---|
183 |
print "error: No corresponding .was file" |
---|
184 |
print " .... skipping %s" % (fn,) |
---|
185 |
return data |
---|
186 |
|
---|
187 |
# add these keys, value pairs to dictionary "data" already setup after reading .wap |
---|
188 |
data['dirs'] = numpy.array(numpy.ones((ndir,), dtype=float)*numpy.nan) |
---|
189 |
data['freqs'] = numpy.array(numpy.ones((nfreq,), dtype=float)*numpy.nan) |
---|
190 |
data['Sxx'] = numpy.array(numpy.ones((nbursts,nfreq,ndir), dtype=float)*numpy.nan) |
---|
191 |
data['Sf'] = numpy.array(numpy.ones((nbursts,nfreq), dtype=float)*numpy.nan) |
---|
192 |
data['Stheta'] = numpy.array(numpy.ones((nbursts,ndir), dtype=float)*numpy.nan) |
---|
193 |
data['Stheta_swell'] = numpy.array(numpy.ones((nbursts,ndir), dtype=float)*numpy.nan) |
---|
194 |
data['Stheta_wind'] = numpy.array(numpy.ones((nbursts,ndir), dtype=float)*numpy.nan) |
---|
195 |
data['Hs'] = numpy.array(numpy.ones((nbursts,), dtype=float)*numpy.nan) |
---|
196 |
data['Hs_swell'] = numpy.array(numpy.ones((nbursts,), dtype=float)*numpy.nan) |
---|
197 |
data['Hs_wind'] = numpy.array(numpy.ones((nbursts,), dtype=float)*numpy.nan) |
---|
198 |
data['Tm'] = numpy.array(numpy.ones((nbursts,), dtype=float)*numpy.nan) |
---|
199 |
data['Tm_swell'] = numpy.array(numpy.ones((nbursts,), dtype=float)*numpy.nan) |
---|
200 |
data['Tm_wind'] = numpy.array(numpy.ones((nbursts,), dtype=float)*numpy.nan) |
---|
201 |
data['Tp'] = numpy.array(numpy.ones((nbursts,), dtype=float)*numpy.nan) |
---|
202 |
data['Tp_swell'] = numpy.array(numpy.ones((nbursts,), dtype=float)*numpy.nan) |
---|
203 |
data['Tp_wind'] = numpy.array(numpy.ones((nbursts,), dtype=float)*numpy.nan) |
---|
204 |
data['Dm'] = numpy.array(numpy.ones((nbursts,), dtype=float)*numpy.nan) |
---|
205 |
data['Dm_swell'] = numpy.array(numpy.ones((nbursts,), dtype=float)*numpy.nan) |
---|
206 |
data['Dm_wind'] = numpy.array(numpy.ones((nbursts,), dtype=float)*numpy.nan) |
---|
207 |
data['Dp'] = numpy.array(numpy.ones((nbursts,), dtype=float)*numpy.nan) |
---|
208 |
data['Dp_swell'] = numpy.array(numpy.ones((nbursts,), dtype=float)*numpy.nan) |
---|
209 |
data['Dp_wind'] = numpy.array(numpy.ones((nbursts,), dtype=float)*numpy.nan) |
---|
210 |
|
---|
211 |
# for each burst read nfreq lines |
---|
212 |
for j in range(nbursts): |
---|
213 |
|
---|
214 |
i = 0 |
---|
215 |
Sxx = numpy.array(numpy.ones((nfreq,ndir), dtype=float)*numpy.nan) |
---|
216 |
# each line is a freq, each column is a direction |
---|
217 |
for line in lines[j*nfreq:nfreq*(j+1)]: |
---|
218 |
wds = [] |
---|
219 |
# split line and parse float and integers |
---|
220 |
sw = re.split(' ', line) |
---|
221 |
for s in sw: |
---|
222 |
m = re.search(REAL_RE_STR, s) |
---|
223 |
if m: |
---|
224 |
wds.append(float(m.groups()[0])) |
---|
225 |
|
---|
226 |
# wds[] is in units of Normalized-Energy/degree from .wds |
---|
227 |
# use power (m^2/Hz) at same time and freq from .was to get units of m^2/Hz/deg |
---|
228 |
if len(wds) == ndir: |
---|
229 |
Sxx[i,:] = numpy.array(wds[:])*was_Sf[j,i] # cross spectrum as m^2/Hz/deg |
---|
230 |
i = i+1 |
---|
231 |
|
---|
232 |
# Did we get the number of data rows that we expected? Should equal nfreq |
---|
233 |
if i != nfreq: |
---|
234 |
print "Number of data rows %d does not match expected number %d" % (i, nfreq) |
---|
235 |
|
---|
236 |
|
---|
237 |
# NOTE make fupper location dependent?? (add to config_files??) |
---|
238 |
fupper = 0.65 # upper freq limit 0.65 Hz or wave periods less than T~1.538s |
---|
239 |
iswell = f<=1/10. # swell band for T>10s |
---|
240 |
iwind = (f>1/10.) * (f<=fupper) # wind band 1/fupper<T<10s |
---|
241 |
# NOTE about python boolean overloaded operator '*' == and == bitwise_and() |
---|
242 |
iall = f<=fupper # all wave freq upper limit |
---|
243 |
|
---|
244 |
# compute non-directional spectrum by integrating over all angles |
---|
245 |
# Sxx(freq, dir) sum axis=1 is along direction |
---|
246 |
Sf = Sxx.sum(axis=1)*Dtheta |
---|
247 |
# Sxx(freq, dir) axis=0 is along freq |
---|
248 |
Stheta = Sxx[iall].sum(axis=0)*Df |
---|
249 |
Stheta_s = Sxx[iswell].sum(axis=0)*Df |
---|
250 |
Stheta_w = Sxx[iwind].sum(axis=0)*Df |
---|
251 |
|
---|
252 |
# compute zeroth-, first- and second-moment from the non-directional spectrum |
---|
253 |
# all frequency ranges |
---|
254 |
m0 = Sf[iall].sum()*Df |
---|
255 |
m1 = (f[iall]*Sf[iall]).sum()*Df |
---|
256 |
m2 = ((f[iall]**2)*Sf[iall]).sum()*Df |
---|
257 |
# swell band |
---|
258 |
m0s = Sf[iswell].sum()*Df |
---|
259 |
m1s = (f[iswell]*Sf[iswell]).sum()*Df |
---|
260 |
m2s = ((f[iswell]**2)*Sf[iswell]).sum()*Df |
---|
261 |
# wind band |
---|
262 |
m0w = Sf[iwind].sum()*Df |
---|
263 |
m1w = (f[iwind]*Sf[iwind]).sum()*Df |
---|
264 |
m2w = ((f[iwind]**2)*Sf[iwind]).sum()*Df |
---|
265 |
|
---|
266 |
# Significant Wave Height (Hs) |
---|
267 |
Hs = 4*numpy.sqrt(m0) |
---|
268 |
Hss = 4*numpy.sqrt(m0s) |
---|
269 |
Hsw = 4*numpy.sqrt(m0w) |
---|
270 |
|
---|
271 |
# Mean Wave Period (Tm) |
---|
272 |
Tm = m0/m1 |
---|
273 |
Tms = m0s/m1s |
---|
274 |
Tmw = m0w/m1w |
---|
275 |
|
---|
276 |
# Peak Wave Period (Tp) |
---|
277 |
# imax = Sf[iall]==Sf[iall].max() |
---|
278 |
# Fp = f(imax) |
---|
279 |
# Tp = 1/Fp[0] |
---|
280 |
# This wave parameters added by SH (not in GV's matlab script) |
---|
281 |
# one-liner version of above |
---|
282 |
# Tp = 1/(f[Sf[iall]==Sf[iall].max()][0]) |
---|
283 |
# Tps = 1/(f[Sf[iswell]==Sf[iswell].max()][0]) |
---|
284 |
# Tpw = 1/(f[Sf[iwind]==Sf[iwind].max()][0]) |
---|
285 |
imax = Sf[iall]==Sf[iall].max() |
---|
286 |
Tp = 1/(f[imax][0]) |
---|
287 |
imax = Sf[iswell]==Sf[iswell].max() |
---|
288 |
Tps = 1/(f[imax][0]) |
---|
289 |
|
---|
290 |
imax = Sf[iwind]==Sf[iwind].max() |
---|
291 |
# account for offset of iwind by iswell in finding peak wind freq |
---|
292 |
nswell = len(f[iswell]) |
---|
293 |
false_swell = numpy.array([False for i in range(nswell)]) |
---|
294 |
imax = numpy.concatenate((false_swell,imax)) |
---|
295 |
Tpw = 1/(f[imax][0]) |
---|
296 |
|
---|
297 |
# mean direction of wave approach used by Kuik et al (1989) |
---|
298 |
# Mean wave direction as a function of frequency |
---|
299 |
# for all freq, wind and swell bands as adapted from GV's code |
---|
300 |
# (polar_waves_cur_wds.m, version 8) |
---|
301 |
pi = numpy.pi |
---|
302 |
ac1 = numpy.cos(D*pi/180) |
---|
303 |
as1 = numpy.sin(D*pi/180) |
---|
304 |
|
---|
305 |
ch0 = (ac1*Stheta*Dtheta).sum() |
---|
306 |
sh0 = (as1*Stheta*Dtheta).sum() |
---|
307 |
Dm = numpy.arctan2(sh0,ch0)*180/pi |
---|
308 |
if Dm<0: Dm = Dm+360. |
---|
309 |
|
---|
310 |
ch0s = (ac1*Stheta_s*Dtheta).sum() |
---|
311 |
sh0s = (as1*Stheta_s*Dtheta).sum() |
---|
312 |
Dms = numpy.arctan2(sh0s,ch0s)*180/pi |
---|
313 |
if Dms<0: Dms = Dms+360. |
---|
314 |
|
---|
315 |
ch0w = (ac1*Stheta_w*Dtheta).sum() |
---|
316 |
sh0w = (as1*Stheta_w*Dtheta).sum() |
---|
317 |
Dmw = numpy.arctan2(sh0w,ch0w)*180/pi |
---|
318 |
if Dmw<0: Dmw = Dmw+360. |
---|
319 |
|
---|
320 |
# Peak Wave Direction (Dp) defined as the direction which |
---|
321 |
# corresponds to the "Peak frequency", or Fp. Peak frequency is the |
---|
322 |
# frequency at which the "Spectral density function" is at a |
---|
323 |
# maximum. The spectral density function gives the dependence |
---|
324 |
# with frequency of the energy of the waves considered. also |
---|
325 |
# known as the one-dimensional spectrum or energy spectrum. |
---|
326 |
# Definitions from Metocean Glossary |
---|
327 |
# http://www.ifremer.fr/web-com/glossary |
---|
328 |
# |
---|
329 |
# This wave parameter added by SH (not in GV's matlab script) |
---|
330 |
imax = Sf[iall]==Sf[iall].max() |
---|
331 |
idir = numpy.squeeze(Sxx[imax,:]==Sxx[imax,:].max()) |
---|
332 |
if len(idir.shape)==2: idir = idir[0] |
---|
333 |
if idir.any(): Dp = D[idir][0] |
---|
334 |
else: Dp = numpy.nan |
---|
335 |
|
---|
336 |
imax = Sf[iswell]==Sf[iswell].max() |
---|
337 |
idir = numpy.squeeze(Sxx[imax,:]==Sxx[imax,:].max()) |
---|
338 |
if len(idir.shape)==2: idir = idir[0] |
---|
339 |
if idir.any(): Dps = D[idir][0] |
---|
340 |
else: Dps = numpy.nan |
---|
341 |
|
---|
342 |
imax = Sf[iwind]==Sf[iwind].max() |
---|
343 |
idir = numpy.squeeze(Sxx[imax,:]==Sxx[imax,:].max()) |
---|
344 |
if len(idir.shape)==2: idir = idir[0] |
---|
345 |
if idir.any(): Dpw = D[idir][0] |
---|
346 |
else: Dpw = numpy.nan |
---|
347 |
|
---|
348 |
# --------------------------------------------------------------- |
---|
349 |
# data['dt'][j] = sample_dt # already have |
---|
350 |
data['dirs'] = D |
---|
351 |
data['freqs'] = f |
---|
352 |
|
---|
353 |
data['Sxx'][j] = Sxx # full directional spectrum (m^2/Hz/deg) |
---|
354 |
data['Sf'][j] = Sf # non-directional spectrum (m^2/Hz) |
---|
355 |
data['Stheta'][j] = Stheta # Energy from all freq from each direction |
---|
356 |
data['Stheta_swell'][j] = Stheta_s |
---|
357 |
data['Stheta_wind'][j] = Stheta_w |
---|
358 |
|
---|
359 |
data['Hs'][j] = Hs |
---|
360 |
data['Hs_swell'][j] = Hss |
---|
361 |
data['Hs_wind'][j] = Hsw |
---|
362 |
|
---|
363 |
data['Tm'][j] = Tm |
---|
364 |
data['Tm_swell'][j] = Tms |
---|
365 |
data['Tm_wind'][j] = Tmw |
---|
366 |
|
---|
367 |
data['Tp'][j] = Tp |
---|
368 |
data['Tp_swell'][j] = Tps |
---|
369 |
data['Tp_wind'][j] = Tpw |
---|
370 |
|
---|
371 |
data['Dm'][j] = Dm |
---|
372 |
data['Dm_swell'][j] = Dms |
---|
373 |
data['Dm_wind'][j] = Dmw |
---|
374 |
|
---|
375 |
data['Dp'][j] = Dp |
---|
376 |
data['Dp_swell'][j] = Dps |
---|
377 |
data['Dp_wind'][j] = Dpw |
---|
378 |
|
---|
379 |
# if j==0: |
---|
380 |
# print " Hs(m)\t Tp(s)\t Tm(s)\t Dp(N)\t Dm(N)" |
---|
381 |
# else: |
---|
382 |
# print "%.2g\t %.2g\t %.2g\t %g\t %g" % (Hs, Tp, Tm, Dp, Dm) |
---|
383 |
# print "%.2g\t %.2g\t %.2g\t %g\t %g" % (Hss, Tps, Tms, Dps, Dms) |
---|
384 |
# print "%.2g\t %.2g\t %.2g\t %g\t %g" % (Hsw, Tpw, Tmw, Dpw, Dmw) |
---|
385 |
|
---|
386 |
# print " Waves: All / Swell / Wind -- burst %d, start %d, end %d" % (j, j*nfreq, nfreq*(j+1)) |
---|
387 |
# print " Hs (m): %g /%g /%g" % (Hs, Hss, Hsw) |
---|
388 |
# print " Tp (s): %g /%g /%g" % (Tp, Tps, Tpw) |
---|
389 |
# print " Tm (s): %g /%g /%g" % (Tm, Tms, Tmw) |
---|
390 |
# print " Dp (N): %g /%g /%g" % (Dp, Dps, Dpw) |
---|
391 |
# print " Dm (N): %g /%g /%g" % (Dm, Dms, Dmw) |
---|
392 |
|
---|
393 |
del Sxx, Sf, Stheta, Stheta_w, Stheta_s |
---|
394 |
# for each burst |
---|
395 |
|
---|
396 |
return data |
---|
397 |
|
---|
398 |
def creator(platform_info, sensor_info, data): |
---|
399 |
# |
---|
400 |
# |
---|
401 |
title_str = sensor_info['description']+' at '+ platform_info['location'] |
---|
402 |
global_atts = { |
---|
403 |
'title' : title_str, |
---|
404 |
'institution' : 'University of North Carolina at Chapel Hill (UNC-CH)', |
---|
405 |
'institution_url' : 'http://nccoos.unc.edu', |
---|
406 |
'institution_dods_url' : 'http://nccoos.unc.edu', |
---|
407 |
'metadata_url' : 'http://nccoos.unc.edu', |
---|
408 |
'references' : 'http://nccoos.unc.edu', |
---|
409 |
'contact' : 'Sara Haines (haines@email.unc.edu)', |
---|
410 |
# |
---|
411 |
'source' : 'directional wave (acoustic doppler) observation', |
---|
412 |
'history' : 'raw2proc using ' + sensor_info['process_module'], |
---|
413 |
'comment' : 'File created using pycdf'+pycdfVersion()+' and numpy '+pycdfArrayPkg(), |
---|
414 |
# conventions |
---|
415 |
'Conventions' : 'CF-1.0; SEACOOS-CDL-v2.0', |
---|
416 |
# SEACOOS CDL codes |
---|
417 |
'format_category_code' : 'directional waves', |
---|
418 |
'institution_code' : platform_info['institution'], |
---|
419 |
'platform_code' : platform_info['id'], |
---|
420 |
'package_code' : sensor_info['id'], |
---|
421 |
# institution specific |
---|
422 |
'project' : 'North Carolina Coastal Ocean Observing System (NCCOOS)', |
---|
423 |
'project_url' : 'http://nccoos.unc.edu', |
---|
424 |
# timeframe of data contained in file yyyy-mm-dd HH:MM:SS |
---|
425 |
'start_date' : data['dt'][0].strftime("%Y-%m-%d %H:%M:%S"), |
---|
426 |
'end_date' : data['dt'][-1].strftime("%Y-%m-%d %H:%M:%S"), |
---|
427 |
'release_date' : now_dt.strftime("%Y-%m-%d %H:%M:%S"), |
---|
428 |
# |
---|
429 |
'creation_date' : now_dt.strftime("%Y-%m-%d %H:%M:%S"), |
---|
430 |
'modification_date' : now_dt.strftime("%Y-%m-%d %H:%M:%S"), |
---|
431 |
'process_level' : 'level1', |
---|
432 |
# |
---|
433 |
# must type match to data (e.g. fillvalue is real if data is real) |
---|
434 |
'_FillValue' : -99999., |
---|
435 |
} |
---|
436 |
|
---|
437 |
var_atts = { |
---|
438 |
# coordinate variables |
---|
439 |
'time' : {'short_name': 'time', |
---|
440 |
'long_name': 'Time', |
---|
441 |
'standard_name': 'time', |
---|
442 |
'units': 'seconds since 1970-1-1 00:00:00 -0', # UTC |
---|
443 |
'axis': 'T', |
---|
444 |
}, |
---|
445 |
'lat' : {'short_name': 'lat', |
---|
446 |
'long_name': 'Latitude', |
---|
447 |
'standard_name': 'latitude', |
---|
448 |
'reference':'geographic coordinates', |
---|
449 |
'units': 'degrees_north', |
---|
450 |
'valid_range':(-90.,90.), |
---|
451 |
'axis': 'Y', |
---|
452 |
}, |
---|
453 |
'lon' : {'short_name': 'lon', |
---|
454 |
'long_name': 'Longitude', |
---|
455 |
'standard_name': 'longitude', |
---|
456 |
'reference':'geographic coordinates', |
---|
457 |
'units': 'degrees_east', |
---|
458 |
'valid_range':(-180.,180.), |
---|
459 |
'axis': 'Y', |
---|
460 |
}, |
---|
461 |
'z' : {'short_name': 'z', |
---|
462 |
'long_name': 'Height', |
---|
463 |
'standard_name': 'height', |
---|
464 |
'reference':'zero at sea-surface', |
---|
465 |
'units': 'm', |
---|
466 |
'axis': 'Z', |
---|
467 |
}, |
---|
468 |
'f' : {'short_name': 'f', |
---|
469 |
'long_name': 'Frequency', |
---|
470 |
'standard_name': 'frequency', |
---|
471 |
'units': 'Hz', |
---|
472 |
}, |
---|
473 |
'd' : {'short_name': 'd', |
---|
474 |
'long_name': 'Direction', |
---|
475 |
'standard_name': 'direction', |
---|
476 |
'reference':'clock-wise from True North', |
---|
477 |
'units': 'deg', |
---|
478 |
}, |
---|
479 |
# data variables |
---|
480 |
'Sxx' : {'short_name': 'Sxx', |
---|
481 |
'long_name': 'Directional Spectral Density Function', |
---|
482 |
'definition': 'Distribution of the wave energy with both frequency and direction', |
---|
483 |
'standard_name': 'wave_directional_spectral_density', |
---|
484 |
'units': 'm2 Hz-1 deg-1', |
---|
485 |
}, |
---|
486 |
'Sf' : {'short_name': 'Sf', |
---|
487 |
'long_name': 'Spectral Density Function', |
---|
488 |
'definition': 'Distribution of the wave energy with frequency from all directions', |
---|
489 |
'standard_name': 'wave_spectral_density', |
---|
490 |
'units': 'm2 Hz-1', |
---|
491 |
}, |
---|
492 |
'Stheta' : {'short_name': 'St', |
---|
493 |
'long_name': 'Spectral Density Function', |
---|
494 |
'definition': 'Distribution of the wave energy with direction from all frequencies', |
---|
495 |
'standard_name': 'wave_directional_density', |
---|
496 |
'units': 'm2 deg-1', |
---|
497 |
}, |
---|
498 |
'Stheta_swell' : {'short_name': 'Sts', |
---|
499 |
'long_name': 'Swell Spectral Density Function', |
---|
500 |
'definition': 'Distribution of the wave energy with direction from all swell frequencies', |
---|
501 |
'standard_name': 'swell_wave_directional_density', |
---|
502 |
'units': 'm2 deg-1', |
---|
503 |
}, |
---|
504 |
'Stheta_wind' : {'short_name': 'Stw', |
---|
505 |
'long_name': 'Wind Spectral Density Function', |
---|
506 |
'definition': 'Distribution of the wave energy with direction from all Wind frequencies', |
---|
507 |
'standard_name': 'wind_wave_directional_density', |
---|
508 |
'units': 'm2 deg-1', |
---|
509 |
}, |
---|
510 |
'Hs' : {'short_name': 'Hs', |
---|
511 |
'long_name': 'Significant Wave Height', |
---|
512 |
'definition': 'Four times the square root of the first moment of the wave spectrum (4*sqrt(m0))', |
---|
513 |
'standard_name': 'significant_wave_height', |
---|
514 |
'units': 'm', |
---|
515 |
}, |
---|
516 |
'Hs_swell' : {'short_name': 'Hss', |
---|
517 |
'long_name': 'Significant Swell Wave Height', |
---|
518 |
'definition': 'Four times the square root of the first moment of the swell wave spectrum (4*sqrt(m0s))', |
---|
519 |
'standard_name': 'significant_swell_wave_height', |
---|
520 |
'units': 'm', |
---|
521 |
}, |
---|
522 |
'Hs_wind' : {'short_name': 'Hsw', |
---|
523 |
'long_name': 'Significant Wind Wave Height', |
---|
524 |
'definition': 'Four times the square root of the first moment of the wind wave spectrum (4*sqrt(m0w))', |
---|
525 |
'standard_name': 'significant_wind_wave_height', |
---|
526 |
'units': 'm', |
---|
527 |
}, |
---|
528 |
'Tp' : {'short_name': 'Tp', |
---|
529 |
'long_name': 'Peak Wave Period', |
---|
530 |
'definition': 'Period of strongest wave (Sf maximum)', |
---|
531 |
'standard_name': 'peak_wave_period', |
---|
532 |
'units': 'sec', |
---|
533 |
}, |
---|
534 |
'Tp_swell' : {'short_name': 'Tps', |
---|
535 |
'long_name': 'Peak Swell Wave Period', |
---|
536 |
'definition': 'Period of strongest swell (Sfs energy maximum)', |
---|
537 |
'standard_name': 'peak_swell_wave_period', |
---|
538 |
'units': 'sec', |
---|
539 |
}, |
---|
540 |
'Tp_wind' : {'short_name': 'Tpw', |
---|
541 |
'long_name': 'Peak Wind Wave Period', |
---|
542 |
'definition': 'Period of strongest wind wave (Sfw energy maximum)', |
---|
543 |
'standard_name': 'peak_wind_wave_period', |
---|
544 |
'units': 'sec', |
---|
545 |
}, |
---|
546 |
'Tm' : {'short_name': 'Tm', |
---|
547 |
'long_name': 'Mean Wave Period', |
---|
548 |
'definition': 'Zero-moment of the non-directional spectrum divided by the first-moment (m0/m1)', |
---|
549 |
'standard_name': 'mean_wave_period', |
---|
550 |
'units': 'sec', |
---|
551 |
}, |
---|
552 |
'Tm_swell' : {'short_name': 'Tms', |
---|
553 |
'long_name': 'Mean Swell Wave Period', |
---|
554 |
'definition': 'Zero-moment of the non-directional spectrum divided by the first-moment (m0s/m1s)', |
---|
555 |
'standard_name': 'mean_swell_wave_period', |
---|
556 |
'units': 'sec', |
---|
557 |
}, |
---|
558 |
'Tm_wind' : {'short_name': 'Tmw', |
---|
559 |
'long_name': 'Mean Wind Wave Period', |
---|
560 |
'definition': 'Zero-moment of the non-directional spectrum divided by the first-moment (m0w/m1w)', |
---|
561 |
'standard_name': 'mean_wind_wave_period', |
---|
562 |
'units': 'sec', |
---|
563 |
}, |
---|
564 |
'Dp' : {'short_name': 'Dp', |
---|
565 |
'long_name': 'Peak Wave Direction', |
---|
566 |
'definition': 'Direction from which strongest waves (wave energy) are coming (dir of max(S(Tp,dir)', |
---|
567 |
'standard_name': 'peak_wave_from_direction', |
---|
568 |
'units': 'deg from N', |
---|
569 |
'reference': 'clockwise from True North', |
---|
570 |
}, |
---|
571 |
'Dp_swell' : {'short_name': 'Dps', |
---|
572 |
'long_name': 'Peak Swell Wave Direction', |
---|
573 |
'definition': 'Direction from which strongest waves (swell energy) are coming (dir of max(S(Tps,dir)', |
---|
574 |
'standard_name': 'peak_wave_from_direction', |
---|
575 |
'units': 'deg from N', |
---|
576 |
'reference': 'clockwise from True North', |
---|
577 |
}, |
---|
578 |
'Dp_wind' : {'short_name': 'Dpw', |
---|
579 |
'long_name': 'Peak Wind Wave Direction', |
---|
580 |
'definition': 'Direction from which strongest waves (wind wave energy) are coming (dir of max(S(Tpw,dir)', |
---|
581 |
'standard_name': 'peak_wave_from_direction', |
---|
582 |
'units': 'deg from N', |
---|
583 |
'reference': 'clockwise from True North', |
---|
584 |
}, |
---|
585 |
'Dm' : {'short_name': 'Dm', |
---|
586 |
'long_name': 'Mean Wave Direction', |
---|
587 |
'definition': 'Mean direction from which strongest waves (wave energy max) are coming for all frequencies', |
---|
588 |
'standard_name': 'mean_wave_from_direction', |
---|
589 |
'units': 'deg from N', |
---|
590 |
'reference': 'clockwise from True North', |
---|
591 |
}, |
---|
592 |
'Dm_swell' : {'short_name': 'Dms', |
---|
593 |
'long_name': 'Mean Swell Wave Direction', |
---|
594 |
'definition': 'Mean direction from which strongest waves (wave energy max) are coming for swell frequencies', |
---|
595 |
'standard_name': 'mean_swell_wave_from_direction', |
---|
596 |
'units': 'deg from N', |
---|
597 |
'reference': 'clockwise from True North', |
---|
598 |
}, |
---|
599 |
'Dm_wind' : {'short_name': 'Dmw', |
---|
600 |
'long_name': 'Mean Wind Wave Direction', |
---|
601 |
'definition': 'Mean direction from which strongest waves (wave energy max) are coming for wind wave frequencies', |
---|
602 |
'standard_name': 'mean_wind_wave_from_direction', |
---|
603 |
'units': 'deg from N', |
---|
604 |
'reference': 'clockwise from True North', |
---|
605 |
}, |
---|
606 |
} |
---|
607 |
|
---|
608 |
|
---|
609 |
# dimension names use tuple so order of initialization is maintained |
---|
610 |
dim_inits = ( |
---|
611 |
('ntime', NC.UNLIMITED), |
---|
612 |
('nlat', 1), |
---|
613 |
('nlon', 1), |
---|
614 |
('nz', 1), |
---|
615 |
('nfreq', sensor_info['nfreq']), |
---|
616 |
('ndir', sensor_info['ndir']), |
---|
617 |
) |
---|
618 |
|
---|
619 |
# using tuple of tuples so order of initialization is maintained |
---|
620 |
# using dict for attributes order of init not important |
---|
621 |
# use dimension names not values |
---|
622 |
# (varName, varType, (dimName1, [dimName2], ...)) |
---|
623 |
var_inits = ( |
---|
624 |
# coordinate variables |
---|
625 |
('time', NC.INT, ('ntime',)), |
---|
626 |
('lat', NC.FLOAT, ('nlat',)), |
---|
627 |
('lon', NC.FLOAT, ('nlon',)), |
---|
628 |
('z', NC.FLOAT, ('nz',)), |
---|
629 |
('f', NC.FLOAT, ('nfreq',)), |
---|
630 |
('d', NC.FLOAT, ('ndir',)), |
---|
631 |
# data variables |
---|
632 |
('Sxx', NC.FLOAT, ('ntime','nfreq','ndir')), |
---|
633 |
('Sf', NC.FLOAT, ('ntime','nfreq')), |
---|
634 |
('Stheta', NC.FLOAT, ('ntime','ndir')), |
---|
635 |
('Stheta_swell', NC.FLOAT, ('ntime','ndir')), |
---|
636 |
('Stheta_wind', NC.FLOAT, ('ntime','ndir')), |
---|
637 |
('Hs', NC.FLOAT, ('ntime',)), |
---|
638 |
('Hs_swell', NC.FLOAT, ('ntime',)), |
---|
639 |
('Hs_wind', NC.FLOAT, ('ntime',)), |
---|
640 |
('Tp', NC.FLOAT, ('ntime',)), |
---|
641 |
('Tp_swell', NC.FLOAT, ('ntime',)), |
---|
642 |
('Tp_wind', NC.FLOAT, ('ntime',)), |
---|
643 |
('Tm', NC.FLOAT, ('ntime',)), |
---|
644 |
('Tm_swell', NC.FLOAT, ('ntime',)), |
---|
645 |
('Tm_wind', NC.FLOAT, ('ntime',)), |
---|
646 |
('Dp', NC.FLOAT, ('ntime',)), |
---|
647 |
('Dp_swell', NC.FLOAT, ('ntime',)), |
---|
648 |
('Dp_wind', NC.FLOAT, ('ntime',)), |
---|
649 |
('Dm', NC.FLOAT, ('ntime',)), |
---|
650 |
('Dm_swell', NC.FLOAT, ('ntime',)), |
---|
651 |
('Dm_wind', NC.FLOAT, ('ntime',)), |
---|
652 |
) |
---|
653 |
|
---|
654 |
# subset data only to month being processed (see raw2proc.process()) |
---|
655 |
i = data['in'] |
---|
656 |
|
---|
657 |
# var data |
---|
658 |
var_data = ( |
---|
659 |
('lat', platform_info['lat']), |
---|
660 |
('lon', platform_info['lon']), |
---|
661 |
('z', 0), |
---|
662 |
('f', data['freqs']), |
---|
663 |
('d', data['dirs']), |
---|
664 |
# |
---|
665 |
('time', data['time'][i]), |
---|
666 |
('Sxx', data['Sxx'][i]), |
---|
667 |
('Sf', data['Sf'][i]), |
---|
668 |
('Stheta', data['Stheta'][i]), |
---|
669 |
('Stheta_swell', data['Stheta_swell'][i]), |
---|
670 |
('Stheta_wind', data['Stheta_wind'][i]), |
---|
671 |
('Hs', data['Hs'][i]), |
---|
672 |
('Hs_swell', data['Hs_swell'][i]), |
---|
673 |
('Hs_wind', data['Hs_wind'][i]), |
---|
674 |
('Tp', data['Tp'][i]), |
---|
675 |
('Tp_swell', data['Tp_swell'][i]), |
---|
676 |
('Tp_wind', data['Tp_wind'][i]), |
---|
677 |
('Tm', data['Tm'][i]), |
---|
678 |
('Tm_swell', data['Tm_swell'][i]), |
---|
679 |
('Tm_wind', data['Tm_wind'][i]), |
---|
680 |
('Dp', data['Dp'][i]), |
---|
681 |
('Dp_swell', data['Dp_swell'][i]), |
---|
682 |
('Dp_wind', data['Dp_wind'][i]), |
---|
683 |
('Dm', data['Dm'][i]), |
---|
684 |
('Dm_swell', data['Tm_swell'][i]), |
---|
685 |
('Dm_wind', data['Tm_wind'][i]), |
---|
686 |
) |
---|
687 |
|
---|
688 |
return (global_atts, var_atts, dim_inits, var_inits, var_data) |
---|
689 |
|
---|
690 |
def updater(platform_info, sensor_info, data): |
---|
691 |
# |
---|
692 |
global_atts = { |
---|
693 |
# update times of data contained in file (yyyy-mm-dd HH:MM:SS) |
---|
694 |
# last date in monthly file |
---|
695 |
'end_date' : data['dt'][-1].strftime("%Y-%m-%d %H:%M:%S"), |
---|
696 |
'release_date' : now_dt.strftime("%Y-%m-%d %H:%M:%S"), |
---|
697 |
# |
---|
698 |
'modification_date' : now_dt.strftime("%Y-%m-%d %H:%M:%S"), |
---|
699 |
} |
---|
700 |
|
---|
701 |
# data variables |
---|
702 |
# update any variable attributes like range, min, max |
---|
703 |
var_atts = {} |
---|
704 |
# var_atts = { |
---|
705 |
# 'u': {'max': max(data.u), |
---|
706 |
# 'min': min(data.v), |
---|
707 |
# }, |
---|
708 |
# 'v': {'max': max(data.u), |
---|
709 |
# 'min': min(data.v), |
---|
710 |
# }, |
---|
711 |
# } |
---|
712 |
|
---|
713 |
# subset data only to month being processed (see raw2proc.process()) |
---|
714 |
i = data['in'] |
---|
715 |
|
---|
716 |
# data |
---|
717 |
var_data = ( |
---|
718 |
('time', data['time'][i]), |
---|
719 |
('Sxx', data['Sxx'][i]), |
---|
720 |
('Sf', data['Sf'][i]), |
---|
721 |
('Stheta', data['Stheta'][i]), |
---|
722 |
('Stheta_swell', data['Stheta_swell'][i]), |
---|
723 |
('Stheta_wind', data['Stheta_wind'][i]), |
---|
724 |
('Hs', data['Hs'][i]), |
---|
725 |
('Hs_swell', data['Hs_swell'][i]), |
---|
726 |
('Hs_wind', data['Hs_wind'][i]), |
---|
727 |
('Tp', data['Tp'][i]), |
---|
728 |
('Tp_swell', data['Tp_swell'][i]), |
---|
729 |
('Tp_wind', data['Tp_wind'][i]), |
---|
730 |
('Tm', data['Tm'][i]), |
---|
731 |
('Tm_swell', data['Tm_swell'][i]), |
---|
732 |
('Tm_wind', data['Tm_wind'][i]), |
---|
733 |
('Dp', data['Dp'][i]), |
---|
734 |
('Dp_swell', data['Dp_swell'][i]), |
---|
735 |
('Dp_wind', data['Dp_wind'][i]), |
---|
736 |
('Dm', data['Dm'][i]), |
---|
737 |
('Dm_swell', data['Dm_swell'][i]), |
---|
738 |
('Dm_wind', data['Dm_wind'][i]), |
---|
739 |
) |
---|
740 |
|
---|
741 |
return (global_atts, var_atts, var_data) |
---|
742 |
|
---|
743 |
# |
---|