##// END OF EJS Templates
schain-jc
RSS
Fix utc times and bugs in BLTR modules
Juan C. Espinoza -
r1085:0d7a50b6a171
parent child
Show More
Side by Side Unified
Context file:
r1085:0d7a50b6a171 Loading diff...:
Show file before | Show file after | Hide comments
@@ -1,362 +1,366
1 '''
1 '''
2 Created on Nov 9, 2016
2 Created on Nov 9, 2016
3
3
4 @author: roj- LouVD
4 @author: roj- LouVD
5 '''
5 '''
6
6
7
7
8 import os
8 import os
9 import sys
9 import sys
10 import time
10 import time
11 import glob
11 import glob
12 import datetime
12 import datetime
13
13
14 import numpy
14 import numpy
15
15
16 from schainpy.model.proc.jroproc_base import ProcessingUnit
16 from schainpy.model.proc.jroproc_base import ProcessingUnit
17 from schainpy.model.data.jrodata import Parameters
17 from schainpy.model.data.jrodata import Parameters
18 from schainpy.model.io.jroIO_base import JRODataReader, isNumber
18 from schainpy.model.io.jroIO_base import JRODataReader, isNumber
19 from schainpy.utils import log
19
20
20 FILE_HEADER_STRUCTURE = numpy.dtype([
21 FILE_HEADER_STRUCTURE = numpy.dtype([
21 ('FMN', '<u4'),
22 ('FMN', '<u4'),
22 ('nrec', '<u4'),
23 ('nrec', '<u4'),
23 ('fr_offset', '<u4'),
24 ('fr_offset', '<u4'),
24 ('id', '<u4'),
25 ('id', '<u4'),
25 ('site', 'u1', (32,))
26 ('site', 'u1', (32,))
26 ])
27 ])
27
28
28 REC_HEADER_STRUCTURE = numpy.dtype([
29 REC_HEADER_STRUCTURE = numpy.dtype([
29 ('rmn', '<u4'),
30 ('rmn', '<u4'),
30 ('rcounter', '<u4'),
31 ('rcounter', '<u4'),
31 ('nr_offset', '<u4'),
32 ('nr_offset', '<u4'),
32 ('tr_offset', '<u4'),
33 ('tr_offset', '<u4'),
33 ('time', '<u4'),
34 ('time', '<u4'),
34 ('time_msec', '<u4'),
35 ('time_msec', '<u4'),
35 ('tag', 'u1', (32,)),
36 ('tag', 'u1', (32,)),
36 ('comments', 'u1', (32,)),
37 ('comments', 'u1', (32,)),
37 ('lat', '<f4'),
38 ('lat', '<f4'),
38 ('lon', '<f4'),
39 ('lon', '<f4'),
39 ('gps_status', '<u4'),
40 ('gps_status', '<u4'),
40 ('freq', '<u4'),
41 ('freq', '<u4'),
41 ('freq0', '<u4'),
42 ('freq0', '<u4'),
42 ('nchan', '<u4'),
43 ('nchan', '<u4'),
43 ('delta_r', '<u4'),
44 ('delta_r', '<u4'),
44 ('nranges', '<u4'),
45 ('nranges', '<u4'),
45 ('r0', '<u4'),
46 ('r0', '<u4'),
46 ('prf', '<u4'),
47 ('prf', '<u4'),
47 ('ncoh', '<u4'),
48 ('ncoh', '<u4'),
48 ('npoints', '<u4'),
49 ('npoints', '<u4'),
49 ('polarization', '<i4'),
50 ('polarization', '<i4'),
50 ('rx_filter', '<u4'),
51 ('rx_filter', '<u4'),
51 ('nmodes', '<u4'),
52 ('nmodes', '<u4'),
52 ('dmode_index', '<u4'),
53 ('dmode_index', '<u4'),
53 ('dmode_rngcorr', '<u4'),
54 ('dmode_rngcorr', '<u4'),
54 ('nrxs', '<u4'),
55 ('nrxs', '<u4'),
55 ('acf_length', '<u4'),
56 ('acf_length', '<u4'),
56 ('acf_lags', '<u4'),
57 ('acf_lags', '<u4'),
57 ('sea_to_atmos', '<f4'),
58 ('sea_to_atmos', '<f4'),
58 ('sea_notch', '<u4'),
59 ('sea_notch', '<u4'),
59 ('lh_sea', '<u4'),
60 ('lh_sea', '<u4'),
60 ('hh_sea', '<u4'),
61 ('hh_sea', '<u4'),
61 ('nbins_sea', '<u4'),
62 ('nbins_sea', '<u4'),
62 ('min_snr', '<f4'),
63 ('min_snr', '<f4'),
63 ('min_cc', '<f4'),
64 ('min_cc', '<f4'),
64 ('max_time_diff', '<f4')
65 ('max_time_diff', '<f4')
65 ])
66 ])
66
67
67 DATA_STRUCTURE = numpy.dtype([
68 DATA_STRUCTURE = numpy.dtype([
68 ('range', '<u4'),
69 ('range', '<u4'),
69 ('status', '<u4'),
70 ('status', '<u4'),
70 ('zonal', '<f4'),
71 ('zonal', '<f4'),
71 ('meridional', '<f4'),
72 ('meridional', '<f4'),
72 ('vertical', '<f4'),
73 ('vertical', '<f4'),
73 ('zonal_a', '<f4'),
74 ('zonal_a', '<f4'),
74 ('meridional_a', '<f4'),
75 ('meridional_a', '<f4'),
75 ('corrected_fading', '<f4'), # seconds
76 ('corrected_fading', '<f4'), # seconds
76 ('uncorrected_fading', '<f4'), # seconds
77 ('uncorrected_fading', '<f4'), # seconds
77 ('time_diff', '<f4'),
78 ('time_diff', '<f4'),
78 ('major_axis', '<f4'),
79 ('major_axis', '<f4'),
79 ('axial_ratio', '<f4'),
80 ('axial_ratio', '<f4'),
80 ('orientation', '<f4'),
81 ('orientation', '<f4'),
81 ('sea_power', '<u4'),
82 ('sea_power', '<u4'),
82 ('sea_algorithm', '<u4')
83 ('sea_algorithm', '<u4')
83 ])
84 ])
84
85
85 class BLTRParamReader(JRODataReader, ProcessingUnit):
86 class BLTRParamReader(JRODataReader, ProcessingUnit):
86 '''
87 '''
87 Boundary Layer and Tropospheric Radar (BLTR) reader, Wind velocities and SNR from *.sswma files
88 Boundary Layer and Tropospheric Radar (BLTR) reader, Wind velocities and SNR from *.sswma files
88 '''
89 '''
89
90
90 ext = '.sswma'
91 ext = '.sswma'
91
92
92 def __init__(self, **kwargs):
93 def __init__(self, **kwargs):
93
94
94 ProcessingUnit.__init__(self , **kwargs)
95 ProcessingUnit.__init__(self , **kwargs)
95
96
96 self.dataOut = Parameters()
97 self.dataOut = Parameters()
97 self.counter_records = 0
98 self.counter_records = 0
98 self.flagNoMoreFiles = 0
99 self.flagNoMoreFiles = 0
99 self.isConfig = False
100 self.isConfig = False
100 self.filename = None
101 self.filename = None
101
102
102 def setup(self,
103 def setup(self,
103 path=None,
104 path=None,
104 startDate=None,
105 startDate=None,
105 endDate=None,
106 endDate=None,
106 ext=None,
107 ext=None,
107 startTime=datetime.time(0, 0, 0),
108 startTime=datetime.time(0, 0, 0),
108 endTime=datetime.time(23, 59, 59),
109 endTime=datetime.time(23, 59, 59),
109 timezone=0,
110 timezone=0,
110 status_value=0,
111 status_value=0,
111 **kwargs):
112 **kwargs):
112
113
113 self.path = path
114 self.path = path
115 self.startDate = startDate
116 self.endDate = endDate
114 self.startTime = startTime
117 self.startTime = startTime
115 self.endTime = endTime
118 self.endTime = endTime
116 self.status_value = status_value
119 self.status_value = status_value
120 self.datatime = datetime.datetime(1900,1,1)
117
121
118 if self.path is None:
122 if self.path is None:
119 raise ValueError, "The path is not valid"
123 raise ValueError, "The path is not valid"
120
124
121 if ext is None:
125 if ext is None:
122 ext = self.ext
126 ext = self.ext
123
127
124 self.search_files(self.path, startDate, endDate, ext)
128 self.search_files(self.path, startDate, endDate, ext)
125 self.timezone = timezone
129 self.timezone = timezone
126 self.fileIndex = 0
130 self.fileIndex = 0
127
131
128 if not self.fileList:
132 if not self.fileList:
129 raise Warning, "There is no files matching these date in the folder: %s. \n Check 'startDate' and 'endDate' "%(path)
133 raise Warning, "There is no files matching these date in the folder: %s. \n Check 'startDate' and 'endDate' "%(path)
130
134
131 self.setNextFile()
135 self.setNextFile()
132
136
133 def search_files(self, path, startDate, endDate, ext):
137 def search_files(self, path, startDate, endDate, ext):
134 '''
138 '''
135 Searching for BLTR rawdata file in path
139 Searching for BLTR rawdata file in path
136 Creating a list of file to proces included in [startDate,endDate]
140 Creating a list of file to proces included in [startDate,endDate]
137
141
138 Input:
142 Input:
139 path - Path to find BLTR rawdata files
143 path - Path to find BLTR rawdata files
140 startDate - Select file from this date
144 startDate - Select file from this date
141 enDate - Select file until this date
145 enDate - Select file until this date
142 ext - Extension of the file to read
146 ext - Extension of the file to read
143
147
144 '''
148 '''
145
149
146 print 'Searching file in %s ' % (path)
150 log.success('Searching files in {} '.format(path), 'BLTRParamReader')
147 foldercounter = 0
151 foldercounter = 0
148 fileList0 = glob.glob1(path, "*%s" % ext)
152 fileList0 = glob.glob1(path, "*%s" % ext)
149 fileList0.sort()
153 fileList0.sort()
150
154
151 self.fileList = []
155 self.fileList = []
152 self.dateFileList = []
156 self.dateFileList = []
153
157
154 for thisFile in fileList0:
158 for thisFile in fileList0:
155 year = thisFile[-14:-10]
159 year = thisFile[-14:-10]
156 if not isNumber(year):
160 if not isNumber(year):
157 continue
161 continue
158
162
159 month = thisFile[-10:-8]
163 month = thisFile[-10:-8]
160 if not isNumber(month):
164 if not isNumber(month):
161 continue
165 continue
162
166
163 day = thisFile[-8:-6]
167 day = thisFile[-8:-6]
164 if not isNumber(day):
168 if not isNumber(day):
165 continue
169 continue
166
170
167 year, month, day = int(year), int(month), int(day)
171 year, month, day = int(year), int(month), int(day)
168 dateFile = datetime.date(year, month, day)
172 dateFile = datetime.date(year, month, day)
169
173
170 if (startDate > dateFile) or (endDate < dateFile):
174 if (startDate > dateFile) or (endDate < dateFile):
171 continue
175 continue
172
176
173 self.fileList.append(thisFile)
177 self.fileList.append(thisFile)
174 self.dateFileList.append(dateFile)
178 self.dateFileList.append(dateFile)
175
179
176 return
180 return
177
181
178 def setNextFile(self):
182 def setNextFile(self):
179
183
180 file_id = self.fileIndex
184 file_id = self.fileIndex
181
185
182 if file_id == len(self.fileList):
186 if file_id == len(self.fileList):
183 print '\nNo more files in the folder'
187 log.success('No more files in the folder', 'BLTRParamReader')
184 print 'Total number of file(s) read : {}'.format(self.fileIndex + 1)
185 self.flagNoMoreFiles = 1
188 self.flagNoMoreFiles = 1
186 return 0
189 return 0
187
190
188 print '\n[Setting file] (%s) ...' % self.fileList[file_id]
191 log.success('Opening {}'.format(self.fileList[file_id]), 'BLTRParamReader')
189 filename = os.path.join(self.path, self.fileList[file_id])
192 filename = os.path.join(self.path, self.fileList[file_id])
190
193
191 dirname, name = os.path.split(filename)
194 dirname, name = os.path.split(filename)
192 self.siteFile = name.split('.')[0] # 'peru2' ---> Piura - 'peru1' ---> Huancayo or Porcuya
195 self.siteFile = name.split('.')[0] # 'peru2' ---> Piura - 'peru1' ---> Huancayo or Porcuya
193 if self.filename is not None:
196 if self.filename is not None:
194 self.fp.close()
197 self.fp.close()
195 self.filename = filename
198 self.filename = filename
196 self.fp = open(self.filename, 'rb')
199 self.fp = open(self.filename, 'rb')
197 self.header_file = numpy.fromfile(self.fp, FILE_HEADER_STRUCTURE, 1)
200 self.header_file = numpy.fromfile(self.fp, FILE_HEADER_STRUCTURE, 1)
198 self.nrecords = self.header_file['nrec'][0]
201 self.nrecords = self.header_file['nrec'][0]
199 self.sizeOfFile = os.path.getsize(self.filename)
202 self.sizeOfFile = os.path.getsize(self.filename)
200 self.counter_records = 0
203 self.counter_records = 0
201 self.flagIsNewFile = 0
204 self.flagIsNewFile = 0
202 self.fileIndex += 1
205 self.fileIndex += 1
203
206
204 return 1
207 return 1
205
208
206 def readNextBlock(self):
209 def readNextBlock(self):
207
210
208 while True:
211 while True:
209 if self.counter_records == self.nrecords:
212 if self.counter_records == self.nrecords:
210 self.flagIsNewFile = 1
213 self.flagIsNewFile = 1
211 if not self.setNextFile():
214 if not self.setNextFile():
212 return 0
215 return 0
213
216
214 self.readBlock()
217 self.readBlock()
215
218
216 if (self.datatime.time() < self.startTime) or (self.datatime.time() > self.endTime):
219 if (self.datatime < datetime.datetime.combine(self.startDate, self.startTime)) or \
217 print "[Reading] Record No. %d/%d -> %s [Skipping]" %(
220 (self.datatime > datetime.datetime.combine(self.endDate, self.endTime)):
218 self.counter_records,
221 log.warning(
219 self.nrecords,
222 'Reading Record No. {}/{} -> {} [Skipping]'.format(
220 self.datatime.ctime())
223 self.counter_records,
224 self.nrecords,
225 self.datatime.ctime()),
226 'BLTRParamReader')
221 continue
227 continue
222 break
228 break
223
229
224 print "[Reading] Record No. %d/%d -> %s" %(
230 log.log('Reading Record No. {}/{} -> {}'.format(
225 self.counter_records,
231 self.counter_records,
226 self.nrecords,
232 self.nrecords,
227 self.datatime.ctime())
233 self.datatime.ctime()), 'BLTRParamReader')
228
234
229 return 1
235 return 1
230
236
231 def readBlock(self):
237 def readBlock(self):
232
238
233 pointer = self.fp.tell()
239 pointer = self.fp.tell()
234 header_rec = numpy.fromfile(self.fp, REC_HEADER_STRUCTURE, 1)
240 header_rec = numpy.fromfile(self.fp, REC_HEADER_STRUCTURE, 1)
235 self.nchannels = header_rec['nchan'][0]/2
241 self.nchannels = header_rec['nchan'][0]/2
236 self.kchan = header_rec['nrxs'][0]
242 self.kchan = header_rec['nrxs'][0]
237 self.nmodes = header_rec['nmodes'][0]
243 self.nmodes = header_rec['nmodes'][0]
238 self.nranges = header_rec['nranges'][0]
244 self.nranges = header_rec['nranges'][0]
239 self.fp.seek(pointer)
245 self.fp.seek(pointer)
240 self.height = numpy.empty((self.nmodes, self.nranges))
246 self.height = numpy.empty((self.nmodes, self.nranges))
241 self.snr = numpy.empty((self.nmodes, self.nchannels, self.nranges))
247 self.snr = numpy.empty((self.nmodes, self.nchannels, self.nranges))
242 self.buffer = numpy.empty((self.nmodes, 3, self.nranges))
248 self.buffer = numpy.empty((self.nmodes, 3, self.nranges))
249 self.flagDiscontinuousBlock = 0
243
250
244 for mode in range(self.nmodes):
251 for mode in range(self.nmodes):
245 self.readHeader()
252 self.readHeader()
246 data = self.readData()
253 data = self.readData()
247 self.height[mode] = (data[0] - self.correction) / 1000.
254 self.height[mode] = (data[0] - self.correction) / 1000.
248 self.buffer[mode] = data[1]
255 self.buffer[mode] = data[1]
249 self.snr[mode] = data[2]
256 self.snr[mode] = data[2]
250
257
251 self.counter_records = self.counter_records + self.nmodes
258 self.counter_records = self.counter_records + self.nmodes
252
259
253 return
260 return
254
261
255 def readHeader(self):
262 def readHeader(self):
256 '''
263 '''
257 RecordHeader of BLTR rawdata file
264 RecordHeader of BLTR rawdata file
258 '''
265 '''
259
266
260 header_structure = numpy.dtype(
267 header_structure = numpy.dtype(
261 REC_HEADER_STRUCTURE.descr + [
268 REC_HEADER_STRUCTURE.descr + [
262 ('antenna_coord', 'f4', (2, self.nchannels)),
269 ('antenna_coord', 'f4', (2, self.nchannels)),
263 ('rx_gains', 'u4', (self.nchannels,)),
270 ('rx_gains', 'u4', (self.nchannels,)),
264 ('rx_analysis', 'u4', (self.nchannels,))
271 ('rx_analysis', 'u4', (self.nchannels,))
265 ]
272 ]
266 )
273 )
267
274
268 self.header_rec = numpy.fromfile(self.fp, header_structure, 1)
275 self.header_rec = numpy.fromfile(self.fp, header_structure, 1)
269 self.lat = self.header_rec['lat'][0]
276 self.lat = self.header_rec['lat'][0]
270 self.lon = self.header_rec['lon'][0]
277 self.lon = self.header_rec['lon'][0]
271 self.delta = self.header_rec['delta_r'][0]
278 self.delta = self.header_rec['delta_r'][0]
272 self.correction = self.header_rec['dmode_rngcorr'][0]
279 self.correction = self.header_rec['dmode_rngcorr'][0]
273 self.imode = self.header_rec['dmode_index'][0]
280 self.imode = self.header_rec['dmode_index'][0]
274 self.antenna = self.header_rec['antenna_coord']
281 self.antenna = self.header_rec['antenna_coord']
275 self.rx_gains = self.header_rec['rx_gains']
282 self.rx_gains = self.header_rec['rx_gains']
276 self.time = self.header_rec['time'][0]
283 self.time = self.header_rec['time'][0]
277 tseconds = self.header_rec['time'][0]
284 dt = datetime.datetime.utcfromtimestamp(self.time)
278 local_t1 = time.localtime(tseconds)
285 if dt.date()>self.datatime.date():
279 self.year = local_t1.tm_year
286 self.flagDiscontinuousBlock = 1
280 self.month = local_t1.tm_mon
287 self.datatime = dt
281 self.day = local_t1.tm_mday
282 self.t = datetime.datetime(self.year, self.month, self.day)
283 self.datatime = datetime.datetime.utcfromtimestamp(self.time)
284
288
285 def readData(self):
289 def readData(self):
286 '''
290 '''
287 Reading and filtering data block record of BLTR rawdata file, filtering is according to status_value.
291 Reading and filtering data block record of BLTR rawdata file, filtering is according to status_value.
288
292
289 Input:
293 Input:
290 status_value - Array data is set to NAN for values that are not equal to status_value
294 status_value - Array data is set to NAN for values that are not equal to status_value
291
295
292 '''
296 '''
293
297
294 data_structure = numpy.dtype(
298 data_structure = numpy.dtype(
295 DATA_STRUCTURE.descr + [
299 DATA_STRUCTURE.descr + [
296 ('rx_saturation', 'u4', (self.nchannels,)),
300 ('rx_saturation', 'u4', (self.nchannels,)),
297 ('chan_offset', 'u4', (2 * self.nchannels,)),
301 ('chan_offset', 'u4', (2 * self.nchannels,)),
298 ('rx_amp', 'u4', (self.nchannels,)),
302 ('rx_amp', 'u4', (self.nchannels,)),
299 ('rx_snr', 'f4', (self.nchannels,)),
303 ('rx_snr', 'f4', (self.nchannels,)),
300 ('cross_snr', 'f4', (self.kchan,)),
304 ('cross_snr', 'f4', (self.kchan,)),
301 ('sea_power_relative', 'f4', (self.kchan,))]
305 ('sea_power_relative', 'f4', (self.kchan,))]
302 )
306 )
303
307
304 data = numpy.fromfile(self.fp, data_structure, self.nranges)
308 data = numpy.fromfile(self.fp, data_structure, self.nranges)
305
309
306 height = data['range']
310 height = data['range']
307 winds = numpy.array((data['zonal'], data['meridional'], data['vertical']))
311 winds = numpy.array((data['zonal'], data['meridional'], data['vertical']))
308 snr = data['rx_snr'].T
312 snr = data['rx_snr'].T
309
313
310 winds[numpy.where(winds == -9999.)] = numpy.nan
314 winds[numpy.where(winds == -9999.)] = numpy.nan
311 winds[:, numpy.where(data['status'] != self.status_value)] = numpy.nan
315 winds[:, numpy.where(data['status'] != self.status_value)] = numpy.nan
312 snr[numpy.where(snr == -9999.)] = numpy.nan
316 snr[numpy.where(snr == -9999.)] = numpy.nan
313 snr[:, numpy.where(data['status'] != self.status_value)] = numpy.nan
317 snr[:, numpy.where(data['status'] != self.status_value)] = numpy.nan
314 snr = numpy.power(10, snr / 10)
318 snr = numpy.power(10, snr / 10)
315
319
316 return height, winds, snr
320 return height, winds, snr
317
321
318 def set_output(self):
322 def set_output(self):
319 '''
323 '''
320 Storing data from databuffer to dataOut object
324 Storing data from databuffer to dataOut object
321 '''
325 '''
322
326
323 self.dataOut.data_SNR = self.snr
327 self.dataOut.data_SNR = self.snr
324 self.dataOut.height = self.height
328 self.dataOut.height = self.height
325 self.dataOut.data_output = self.buffer
329 self.dataOut.data = self.buffer
326 self.dataOut.utctimeInit = self.time
330 self.dataOut.utctimeInit = self.time
327 self.dataOut.utctime = self.dataOut.utctimeInit
331 self.dataOut.utctime = self.dataOut.utctimeInit
328 self.dataOut.useLocalTime = False
332 self.dataOut.useLocalTime = False
329 self.dataOut.paramInterval = 157
333 self.dataOut.paramInterval = 157
330 self.dataOut.timezone = self.timezone
334 self.dataOut.timezone = self.timezone
331 self.dataOut.site = self.siteFile
335 self.dataOut.site = self.siteFile
332 self.dataOut.nrecords = self.nrecords/self.nmodes
336 self.dataOut.nrecords = self.nrecords/self.nmodes
333 self.dataOut.sizeOfFile = self.sizeOfFile
337 self.dataOut.sizeOfFile = self.sizeOfFile
334 self.dataOut.lat = self.lat
338 self.dataOut.lat = self.lat
335 self.dataOut.lon = self.lon
339 self.dataOut.lon = self.lon
336 self.dataOut.channelList = range(self.nchannels)
340 self.dataOut.channelList = range(self.nchannels)
337 self.dataOut.kchan = self.kchan
341 self.dataOut.kchan = self.kchan
338 # self.dataOut.nHeights = self.nranges
339 self.dataOut.delta = self.delta
342 self.dataOut.delta = self.delta
340 self.dataOut.correction = self.correction
343 self.dataOut.correction = self.correction
341 self.dataOut.nmodes = self.nmodes
344 self.dataOut.nmodes = self.nmodes
342 self.dataOut.imode = self.imode
345 self.dataOut.imode = self.imode
343 self.dataOut.antenna = self.antenna
346 self.dataOut.antenna = self.antenna
344 self.dataOut.rx_gains = self.rx_gains
347 self.dataOut.rx_gains = self.rx_gains
345 self.dataOut.flagNoData = False
348 self.dataOut.flagNoData = False
349 self.dataOut.flagDiscontinuousBlock = self.flagDiscontinuousBlock
346
350
347 def getData(self):
351 def getData(self):
348 '''
352 '''
349 Storing data from databuffer to dataOut object
353 Storing data from databuffer to dataOut object
350 '''
354 '''
351 if self.flagNoMoreFiles:
355 if self.flagNoMoreFiles:
352 self.dataOut.flagNoData = True
356 self.dataOut.flagNoData = True
353 print 'No file left to process'
357 log.success('No file left to process', 'BLTRParamReader')
354 return 0
358 return 0
355
359
356 if not self.readNextBlock():
360 if not self.readNextBlock():
357 self.dataOut.flagNoData = True
361 self.dataOut.flagNoData = True
358 return 0
362 return 0
359
363
360 self.set_output()
364 self.set_output()
361
365
362 return 1
366 return 1
Show file before | Show file after | Hide comments
@@ -1,403 +1,403
1 '''
1 '''
2 Created on Oct 24, 2016
2 Created on Oct 24, 2016
3
3
4 @author: roj- LouVD
4 @author: roj- LouVD
5 '''
5 '''
6
6
7 import numpy
7 import numpy
8 import copy
8 import copy
9 import datetime
9 import datetime
10 import time
10 import time
11 from time import gmtime
11 from time import gmtime
12
12
13 from numpy import transpose
13 from numpy import transpose
14
14
15 from jroproc_base import ProcessingUnit, Operation
15 from jroproc_base import ProcessingUnit, Operation
16 from schainpy.model.data.jrodata import Parameters
16 from schainpy.model.data.jrodata import Parameters
17
17
18
18
19 class BLTRParametersProc(ProcessingUnit):
19 class BLTRParametersProc(ProcessingUnit):
20 '''
20 '''
21 Processing unit for BLTR parameters data (winds)
21 Processing unit for BLTR parameters data (winds)
22
22
23 Inputs:
23 Inputs:
24 self.dataOut.nmodes - Number of operation modes
24 self.dataOut.nmodes - Number of operation modes
25 self.dataOut.nchannels - Number of channels
25 self.dataOut.nchannels - Number of channels
26 self.dataOut.nranges - Number of ranges
26 self.dataOut.nranges - Number of ranges
27
27
28 self.dataOut.data_SNR - SNR array
28 self.dataOut.data_SNR - SNR array
29 self.dataOut.data_output - Zonal, Vertical and Meridional velocity array
29 self.dataOut.data_output - Zonal, Vertical and Meridional velocity array
30 self.dataOut.height - Height array (km)
30 self.dataOut.height - Height array (km)
31 self.dataOut.time - Time array (seconds)
31 self.dataOut.time - Time array (seconds)
32
32
33 self.dataOut.fileIndex -Index of the file currently read
33 self.dataOut.fileIndex -Index of the file currently read
34 self.dataOut.lat - Latitude coordinate of BLTR location
34 self.dataOut.lat - Latitude coordinate of BLTR location
35
35
36 self.dataOut.doy - Experiment doy (number of the day in the current year)
36 self.dataOut.doy - Experiment doy (number of the day in the current year)
37 self.dataOut.month - Experiment month
37 self.dataOut.month - Experiment month
38 self.dataOut.day - Experiment day
38 self.dataOut.day - Experiment day
39 self.dataOut.year - Experiment year
39 self.dataOut.year - Experiment year
40 '''
40 '''
41
41
42 def __init__(self, **kwargs):
42 def __init__(self, **kwargs):
43 '''
43 '''
44 Inputs: None
44 Inputs: None
45 '''
45 '''
46 ProcessingUnit.__init__(self, **kwargs)
46 ProcessingUnit.__init__(self, **kwargs)
47 self.dataOut = Parameters()
47 self.dataOut = Parameters()
48 self.isConfig = False
48 self.isConfig = False
49
49
50 def setup(self, mode):
50 def setup(self, mode):
51 '''
51 '''
52 '''
52 '''
53 self.dataOut.mode = mode
53 self.dataOut.mode = mode
54
54
55 def run(self, mode, snr_threshold=None):
55 def run(self, mode, snr_threshold=None):
56 '''
56 '''
57 Inputs:
57 Inputs:
58 mode = High resolution (0) or Low resolution (1) data
58 mode = High resolution (0) or Low resolution (1) data
59 snr_threshold = snr filter value
59 snr_threshold = snr filter value
60 '''
60 '''
61
61
62 if not self.isConfig:
62 if not self.isConfig:
63 self.setup(mode)
63 self.setup(mode)
64 self.isConfig = True
64 self.isConfig = True
65
65
66 if self.dataIn.type == 'Parameters':
66 if self.dataIn.type == 'Parameters':
67 self.dataOut.copy(self.dataIn)
67 self.dataOut.copy(self.dataIn)
68
68
69 self.dataOut.data_output = self.dataOut.data_output[mode]
69 self.dataOut.data_param = self.dataOut.data[mode]
70 self.dataOut.heightList = self.dataOut.height[0]
70 self.dataOut.heightList = self.dataOut.height[0]
71 self.dataOut.data_SNR = self.dataOut.data_SNR[mode]
71 self.dataOut.data_SNR = self.dataOut.data_SNR[mode]
72
72
73 if snr_threshold is not None:
73 if snr_threshold is not None:
74 SNRavg = numpy.average(self.dataOut.data_SNR, axis=0)
74 SNRavg = numpy.average(self.dataOut.data_SNR, axis=0)
75 SNRavgdB = 10*numpy.log10(SNRavg)
75 SNRavgdB = 10*numpy.log10(SNRavg)
76 for i in range(3):
76 for i in range(3):
77 self.dataOut.data_output[i][SNRavgdB <= snr_threshold] = numpy.nan
77 self.dataOut.data_param[i][SNRavgdB <= snr_threshold] = numpy.nan
78
78
79 # TODO
79 # TODO
80 class OutliersFilter(Operation):
80 class OutliersFilter(Operation):
81
81
82 def __init__(self, **kwargs):
82 def __init__(self, **kwargs):
83 '''
83 '''
84 '''
84 '''
85 Operation.__init__(self, **kwargs)
85 Operation.__init__(self, **kwargs)
86
86
87 def run(self, svalue2, method, factor, filter, npoints=9):
87 def run(self, svalue2, method, factor, filter, npoints=9):
88 '''
88 '''
89 Inputs:
89 Inputs:
90 svalue - string to select array velocity
90 svalue - string to select array velocity
91 svalue2 - string to choose axis filtering
91 svalue2 - string to choose axis filtering
92 method - 0 for SMOOTH or 1 for MEDIAN
92 method - 0 for SMOOTH or 1 for MEDIAN
93 factor - number used to set threshold
93 factor - number used to set threshold
94 filter - 1 for data filtering using the standard deviation criteria else 0
94 filter - 1 for data filtering using the standard deviation criteria else 0
95 npoints - number of points for mask filter
95 npoints - number of points for mask filter
96 '''
96 '''
97
97
98 print ' Outliers Filter {} {} / threshold = {}'.format(svalue, svalue, factor)
98 print ' Outliers Filter {} {} / threshold = {}'.format(svalue, svalue, factor)
99
99
100
100
101 yaxis = self.dataOut.heightList
101 yaxis = self.dataOut.heightList
102 xaxis = numpy.array([[self.dataOut.utctime]])
102 xaxis = numpy.array([[self.dataOut.utctime]])
103
103
104 # Zonal
104 # Zonal
105 value_temp = self.dataOut.data_output[0]
105 value_temp = self.dataOut.data_output[0]
106
106
107 # Zonal
107 # Zonal
108 value_temp = self.dataOut.data_output[1]
108 value_temp = self.dataOut.data_output[1]
109
109
110 # Vertical
110 # Vertical
111 value_temp = numpy.transpose(self.dataOut.data_output[2])
111 value_temp = numpy.transpose(self.dataOut.data_output[2])
112
112
113 htemp = yaxis
113 htemp = yaxis
114 std = value_temp
114 std = value_temp
115 for h in range(len(htemp)):
115 for h in range(len(htemp)):
116 nvalues_valid = len(numpy.where(numpy.isfinite(value_temp[h]))[0])
116 nvalues_valid = len(numpy.where(numpy.isfinite(value_temp[h]))[0])
117 minvalid = npoints
117 minvalid = npoints
118
118
119 #only if valid values greater than the minimum required (10%)
119 #only if valid values greater than the minimum required (10%)
120 if nvalues_valid > minvalid:
120 if nvalues_valid > minvalid:
121
121
122 if method == 0:
122 if method == 0:
123 #SMOOTH
123 #SMOOTH
124 w = value_temp[h] - self.Smooth(input=value_temp[h], width=npoints, edge_truncate=1)
124 w = value_temp[h] - self.Smooth(input=value_temp[h], width=npoints, edge_truncate=1)
125
125
126
126
127 if method == 1:
127 if method == 1:
128 #MEDIAN
128 #MEDIAN
129 w = value_temp[h] - self.Median(input=value_temp[h], width = npoints)
129 w = value_temp[h] - self.Median(input=value_temp[h], width = npoints)
130
130
131 dw = numpy.std(w[numpy.where(numpy.isfinite(w))],ddof = 1)
131 dw = numpy.std(w[numpy.where(numpy.isfinite(w))],ddof = 1)
132
132
133 threshold = dw*factor
133 threshold = dw*factor
134 value_temp[numpy.where(w > threshold),h] = numpy.nan
134 value_temp[numpy.where(w > threshold),h] = numpy.nan
135 value_temp[numpy.where(w < -1*threshold),h] = numpy.nan
135 value_temp[numpy.where(w < -1*threshold),h] = numpy.nan
136
136
137
137
138 #At the end
138 #At the end
139 if svalue2 == 'inHeight':
139 if svalue2 == 'inHeight':
140 value_temp = numpy.transpose(value_temp)
140 value_temp = numpy.transpose(value_temp)
141 output_array[:,m] = value_temp
141 output_array[:,m] = value_temp
142
142
143 if svalue == 'zonal':
143 if svalue == 'zonal':
144 self.dataOut.data_output[0] = output_array
144 self.dataOut.data_output[0] = output_array
145
145
146 elif svalue == 'meridional':
146 elif svalue == 'meridional':
147 self.dataOut.data_output[1] = output_array
147 self.dataOut.data_output[1] = output_array
148
148
149 elif svalue == 'vertical':
149 elif svalue == 'vertical':
150 self.dataOut.data_output[2] = output_array
150 self.dataOut.data_output[2] = output_array
151
151
152 return self.dataOut.data_output
152 return self.dataOut.data_output
153
153
154
154
155 def Median(self,input,width):
155 def Median(self,input,width):
156 '''
156 '''
157 Inputs:
157 Inputs:
158 input - Velocity array
158 input - Velocity array
159 width - Number of points for mask filter
159 width - Number of points for mask filter
160
160
161 '''
161 '''
162
162
163 if numpy.mod(width,2) == 1:
163 if numpy.mod(width,2) == 1:
164 pc = int((width - 1) / 2)
164 pc = int((width - 1) / 2)
165 cont = 0
165 cont = 0
166 output = []
166 output = []
167
167
168 for i in range(len(input)):
168 for i in range(len(input)):
169 if i >= pc and i < len(input) - pc:
169 if i >= pc and i < len(input) - pc:
170 new2 = input[i-pc:i+pc+1]
170 new2 = input[i-pc:i+pc+1]
171 temp = numpy.where(numpy.isfinite(new2))
171 temp = numpy.where(numpy.isfinite(new2))
172 new = new2[temp]
172 new = new2[temp]
173 value = numpy.median(new)
173 value = numpy.median(new)
174 output.append(value)
174 output.append(value)
175
175
176 output = numpy.array(output)
176 output = numpy.array(output)
177 output = numpy.hstack((input[0:pc],output))
177 output = numpy.hstack((input[0:pc],output))
178 output = numpy.hstack((output,input[-pc:len(input)]))
178 output = numpy.hstack((output,input[-pc:len(input)]))
179
179
180 return output
180 return output
181
181
182 def Smooth(self,input,width,edge_truncate = None):
182 def Smooth(self,input,width,edge_truncate = None):
183 '''
183 '''
184 Inputs:
184 Inputs:
185 input - Velocity array
185 input - Velocity array
186 width - Number of points for mask filter
186 width - Number of points for mask filter
187 edge_truncate - 1 for truncate the convolution product else
187 edge_truncate - 1 for truncate the convolution product else
188
188
189 '''
189 '''
190
190
191 if numpy.mod(width,2) == 0:
191 if numpy.mod(width,2) == 0:
192 real_width = width + 1
192 real_width = width + 1
193 nzeros = width / 2
193 nzeros = width / 2
194 else:
194 else:
195 real_width = width
195 real_width = width
196 nzeros = (width - 1) / 2
196 nzeros = (width - 1) / 2
197
197
198 half_width = int(real_width)/2
198 half_width = int(real_width)/2
199 length = len(input)
199 length = len(input)
200
200
201 gate = numpy.ones(real_width,dtype='float')
201 gate = numpy.ones(real_width,dtype='float')
202 norm_of_gate = numpy.sum(gate)
202 norm_of_gate = numpy.sum(gate)
203
203
204 nan_process = 0
204 nan_process = 0
205 nan_id = numpy.where(numpy.isnan(input))
205 nan_id = numpy.where(numpy.isnan(input))
206 if len(nan_id[0]) > 0:
206 if len(nan_id[0]) > 0:
207 nan_process = 1
207 nan_process = 1
208 pb = numpy.zeros(len(input))
208 pb = numpy.zeros(len(input))
209 pb[nan_id] = 1.
209 pb[nan_id] = 1.
210 input[nan_id] = 0.
210 input[nan_id] = 0.
211
211
212 if edge_truncate == True:
212 if edge_truncate == True:
213 output = numpy.convolve(input/norm_of_gate,gate,mode='same')
213 output = numpy.convolve(input/norm_of_gate,gate,mode='same')
214 elif edge_truncate == False or edge_truncate == None:
214 elif edge_truncate == False or edge_truncate == None:
215 output = numpy.convolve(input/norm_of_gate,gate,mode='valid')
215 output = numpy.convolve(input/norm_of_gate,gate,mode='valid')
216 output = numpy.hstack((input[0:half_width],output))
216 output = numpy.hstack((input[0:half_width],output))
217 output = numpy.hstack((output,input[len(input)-half_width:len(input)]))
217 output = numpy.hstack((output,input[len(input)-half_width:len(input)]))
218
218
219 if nan_process:
219 if nan_process:
220 pb = numpy.convolve(pb/norm_of_gate,gate,mode='valid')
220 pb = numpy.convolve(pb/norm_of_gate,gate,mode='valid')
221 pb = numpy.hstack((numpy.zeros(half_width),pb))
221 pb = numpy.hstack((numpy.zeros(half_width),pb))
222 pb = numpy.hstack((pb,numpy.zeros(half_width)))
222 pb = numpy.hstack((pb,numpy.zeros(half_width)))
223 output[numpy.where(pb > 0.9999)] = numpy.nan
223 output[numpy.where(pb > 0.9999)] = numpy.nan
224 input[nan_id] = numpy.nan
224 input[nan_id] = numpy.nan
225 return output
225 return output
226
226
227 def Average(self,aver=0,nhaver=1):
227 def Average(self,aver=0,nhaver=1):
228 '''
228 '''
229 Inputs:
229 Inputs:
230 aver - Indicates the time period over which is averaged or consensus data
230 aver - Indicates the time period over which is averaged or consensus data
231 nhaver - Indicates the decimation factor in heights
231 nhaver - Indicates the decimation factor in heights
232
232
233 '''
233 '''
234 nhpoints = 48
234 nhpoints = 48
235
235
236 lat_piura = -5.17
236 lat_piura = -5.17
237 lat_huancayo = -12.04
237 lat_huancayo = -12.04
238 lat_porcuya = -5.8
238 lat_porcuya = -5.8
239
239
240 if '%2.2f'%self.dataOut.lat == '%2.2f'%lat_piura:
240 if '%2.2f'%self.dataOut.lat == '%2.2f'%lat_piura:
241 hcm = 3.
241 hcm = 3.
242 if self.dataOut.year == 2003 :
242 if self.dataOut.year == 2003 :
243 if self.dataOut.doy >= 25 and self.dataOut.doy < 64:
243 if self.dataOut.doy >= 25 and self.dataOut.doy < 64:
244 nhpoints = 12
244 nhpoints = 12
245
245
246 elif '%2.2f'%self.dataOut.lat == '%2.2f'%lat_huancayo:
246 elif '%2.2f'%self.dataOut.lat == '%2.2f'%lat_huancayo:
247 hcm = 3.
247 hcm = 3.
248 if self.dataOut.year == 2003 :
248 if self.dataOut.year == 2003 :
249 if self.dataOut.doy >= 25 and self.dataOut.doy < 64:
249 if self.dataOut.doy >= 25 and self.dataOut.doy < 64:
250 nhpoints = 12
250 nhpoints = 12
251
251
252
252
253 elif '%2.2f'%self.dataOut.lat == '%2.2f'%lat_porcuya:
253 elif '%2.2f'%self.dataOut.lat == '%2.2f'%lat_porcuya:
254 hcm = 5.#2
254 hcm = 5.#2
255
255
256 pdata = 0.2
256 pdata = 0.2
257 taver = [1,2,3,4,6,8,12,24]
257 taver = [1,2,3,4,6,8,12,24]
258 t0 = 0
258 t0 = 0
259 tf = 24
259 tf = 24
260 ntime =(tf-t0)/taver[aver]
260 ntime =(tf-t0)/taver[aver]
261 ti = numpy.arange(ntime)
261 ti = numpy.arange(ntime)
262 tf = numpy.arange(ntime) + taver[aver]
262 tf = numpy.arange(ntime) + taver[aver]
263
263
264
264
265 old_height = self.dataOut.heightList
265 old_height = self.dataOut.heightList
266
266
267 if nhaver > 1:
267 if nhaver > 1:
268 num_hei = len(self.dataOut.heightList)/nhaver/self.dataOut.nmodes
268 num_hei = len(self.dataOut.heightList)/nhaver/self.dataOut.nmodes
269 deltha = 0.05*nhaver
269 deltha = 0.05*nhaver
270 minhvalid = pdata*nhaver
270 minhvalid = pdata*nhaver
271 for im in range(self.dataOut.nmodes):
271 for im in range(self.dataOut.nmodes):
272 new_height = numpy.arange(num_hei)*deltha + self.dataOut.height[im,0] + deltha/2.
272 new_height = numpy.arange(num_hei)*deltha + self.dataOut.height[im,0] + deltha/2.
273
273
274
274
275 data_fHeigths_List = []
275 data_fHeigths_List = []
276 data_fZonal_List = []
276 data_fZonal_List = []
277 data_fMeridional_List = []
277 data_fMeridional_List = []
278 data_fVertical_List = []
278 data_fVertical_List = []
279 startDTList = []
279 startDTList = []
280
280
281
281
282 for i in range(ntime):
282 for i in range(ntime):
283 height = old_height
283 height = old_height
284
284
285 start = datetime.datetime(self.dataOut.year,self.dataOut.month,self.dataOut.day) + datetime.timedelta(hours = int(ti[i])) - datetime.timedelta(hours = 5)
285 start = datetime.datetime(self.dataOut.year,self.dataOut.month,self.dataOut.day) + datetime.timedelta(hours = int(ti[i])) - datetime.timedelta(hours = 5)
286 stop = datetime.datetime(self.dataOut.year,self.dataOut.month,self.dataOut.day) + datetime.timedelta(hours = int(tf[i])) - datetime.timedelta(hours = 5)
286 stop = datetime.datetime(self.dataOut.year,self.dataOut.month,self.dataOut.day) + datetime.timedelta(hours = int(tf[i])) - datetime.timedelta(hours = 5)
287
287
288
288
289 limit_sec1 = time.mktime(start.timetuple())
289 limit_sec1 = time.mktime(start.timetuple())
290 limit_sec2 = time.mktime(stop.timetuple())
290 limit_sec2 = time.mktime(stop.timetuple())
291
291
292 t1 = numpy.where(self.f_timesec >= limit_sec1)
292 t1 = numpy.where(self.f_timesec >= limit_sec1)
293 t2 = numpy.where(self.f_timesec < limit_sec2)
293 t2 = numpy.where(self.f_timesec < limit_sec2)
294 time_select = []
294 time_select = []
295 for val_sec in t1[0]:
295 for val_sec in t1[0]:
296 if val_sec in t2[0]:
296 if val_sec in t2[0]:
297 time_select.append(val_sec)
297 time_select.append(val_sec)
298
298
299
299
300 time_select = numpy.array(time_select,dtype = 'int')
300 time_select = numpy.array(time_select,dtype = 'int')
301 minvalid = numpy.ceil(pdata*nhpoints)
301 minvalid = numpy.ceil(pdata*nhpoints)
302
302
303 zon_aver = numpy.zeros([self.dataOut.nranges,self.dataOut.nmodes],dtype='f4') + numpy.nan
303 zon_aver = numpy.zeros([self.dataOut.nranges,self.dataOut.nmodes],dtype='f4') + numpy.nan
304 mer_aver = numpy.zeros([self.dataOut.nranges,self.dataOut.nmodes],dtype='f4') + numpy.nan
304 mer_aver = numpy.zeros([self.dataOut.nranges,self.dataOut.nmodes],dtype='f4') + numpy.nan
305 ver_aver = numpy.zeros([self.dataOut.nranges,self.dataOut.nmodes],dtype='f4') + numpy.nan
305 ver_aver = numpy.zeros([self.dataOut.nranges,self.dataOut.nmodes],dtype='f4') + numpy.nan
306
306
307 if nhaver > 1:
307 if nhaver > 1:
308 new_zon_aver = numpy.zeros([num_hei,self.dataOut.nmodes],dtype='f4') + numpy.nan
308 new_zon_aver = numpy.zeros([num_hei,self.dataOut.nmodes],dtype='f4') + numpy.nan
309 new_mer_aver = numpy.zeros([num_hei,self.dataOut.nmodes],dtype='f4') + numpy.nan
309 new_mer_aver = numpy.zeros([num_hei,self.dataOut.nmodes],dtype='f4') + numpy.nan
310 new_ver_aver = numpy.zeros([num_hei,self.dataOut.nmodes],dtype='f4') + numpy.nan
310 new_ver_aver = numpy.zeros([num_hei,self.dataOut.nmodes],dtype='f4') + numpy.nan
311
311
312 if len(time_select) > minvalid:
312 if len(time_select) > minvalid:
313 time_average = self.f_timesec[time_select]
313 time_average = self.f_timesec[time_select]
314
314
315 for im in range(self.dataOut.nmodes):
315 for im in range(self.dataOut.nmodes):
316
316
317 for ih in range(self.dataOut.nranges):
317 for ih in range(self.dataOut.nranges):
318 if numpy.sum(numpy.isfinite(self.f_zon[time_select,ih,im])) >= minvalid:
318 if numpy.sum(numpy.isfinite(self.f_zon[time_select,ih,im])) >= minvalid:
319 zon_aver[ih,im] = numpy.nansum(self.f_zon[time_select,ih,im]) / numpy.sum(numpy.isfinite(self.f_zon[time_select,ih,im]))
319 zon_aver[ih,im] = numpy.nansum(self.f_zon[time_select,ih,im]) / numpy.sum(numpy.isfinite(self.f_zon[time_select,ih,im]))
320
320
321 if numpy.sum(numpy.isfinite(self.f_mer[time_select,ih,im])) >= minvalid:
321 if numpy.sum(numpy.isfinite(self.f_mer[time_select,ih,im])) >= minvalid:
322 mer_aver[ih,im] = numpy.nansum(self.f_mer[time_select,ih,im]) / numpy.sum(numpy.isfinite(self.f_mer[time_select,ih,im]))
322 mer_aver[ih,im] = numpy.nansum(self.f_mer[time_select,ih,im]) / numpy.sum(numpy.isfinite(self.f_mer[time_select,ih,im]))
323
323
324 if numpy.sum(numpy.isfinite(self.f_ver[time_select,ih,im])) >= minvalid:
324 if numpy.sum(numpy.isfinite(self.f_ver[time_select,ih,im])) >= minvalid:
325 ver_aver[ih,im] = numpy.nansum(self.f_ver[time_select,ih,im]) / numpy.sum(numpy.isfinite(self.f_ver[time_select,ih,im]))
325 ver_aver[ih,im] = numpy.nansum(self.f_ver[time_select,ih,im]) / numpy.sum(numpy.isfinite(self.f_ver[time_select,ih,im]))
326
326
327 if nhaver > 1:
327 if nhaver > 1:
328 for ih in range(num_hei):
328 for ih in range(num_hei):
329 hvalid = numpy.arange(nhaver) + nhaver*ih
329 hvalid = numpy.arange(nhaver) + nhaver*ih
330
330
331 if numpy.sum(numpy.isfinite(zon_aver[hvalid,im])) >= minvalid:
331 if numpy.sum(numpy.isfinite(zon_aver[hvalid,im])) >= minvalid:
332 new_zon_aver[ih,im] = numpy.nansum(zon_aver[hvalid,im]) / numpy.sum(numpy.isfinite(zon_aver[hvalid,im]))
332 new_zon_aver[ih,im] = numpy.nansum(zon_aver[hvalid,im]) / numpy.sum(numpy.isfinite(zon_aver[hvalid,im]))
333
333
334 if numpy.sum(numpy.isfinite(mer_aver[hvalid,im])) >= minvalid:
334 if numpy.sum(numpy.isfinite(mer_aver[hvalid,im])) >= minvalid:
335 new_mer_aver[ih,im] = numpy.nansum(mer_aver[hvalid,im]) / numpy.sum(numpy.isfinite(mer_aver[hvalid,im]))
335 new_mer_aver[ih,im] = numpy.nansum(mer_aver[hvalid,im]) / numpy.sum(numpy.isfinite(mer_aver[hvalid,im]))
336
336
337 if numpy.sum(numpy.isfinite(ver_aver[hvalid,im])) >= minvalid:
337 if numpy.sum(numpy.isfinite(ver_aver[hvalid,im])) >= minvalid:
338 new_ver_aver[ih,im] = numpy.nansum(ver_aver[hvalid,im]) / numpy.sum(numpy.isfinite(ver_aver[hvalid,im]))
338 new_ver_aver[ih,im] = numpy.nansum(ver_aver[hvalid,im]) / numpy.sum(numpy.isfinite(ver_aver[hvalid,im]))
339 if nhaver > 1:
339 if nhaver > 1:
340 zon_aver = new_zon_aver
340 zon_aver = new_zon_aver
341 mer_aver = new_mer_aver
341 mer_aver = new_mer_aver
342 ver_aver = new_ver_aver
342 ver_aver = new_ver_aver
343 height = new_height
343 height = new_height
344
344
345
345
346 tstart = time_average[0]
346 tstart = time_average[0]
347 tend = time_average[-1]
347 tend = time_average[-1]
348 startTime = time.gmtime(tstart)
348 startTime = time.gmtime(tstart)
349
349
350 year = startTime.tm_year
350 year = startTime.tm_year
351 month = startTime.tm_mon
351 month = startTime.tm_mon
352 day = startTime.tm_mday
352 day = startTime.tm_mday
353 hour = startTime.tm_hour
353 hour = startTime.tm_hour
354 minute = startTime.tm_min
354 minute = startTime.tm_min
355 second = startTime.tm_sec
355 second = startTime.tm_sec
356
356
357 startDTList.append(datetime.datetime(year,month,day,hour,minute,second))
357 startDTList.append(datetime.datetime(year,month,day,hour,minute,second))
358
358
359
359
360 o_height = numpy.array([])
360 o_height = numpy.array([])
361 o_zon_aver = numpy.array([])
361 o_zon_aver = numpy.array([])
362 o_mer_aver = numpy.array([])
362 o_mer_aver = numpy.array([])
363 o_ver_aver = numpy.array([])
363 o_ver_aver = numpy.array([])
364 if self.dataOut.nmodes > 1:
364 if self.dataOut.nmodes > 1:
365 for im in range(self.dataOut.nmodes):
365 for im in range(self.dataOut.nmodes):
366
366
367 if im == 0:
367 if im == 0:
368 h_select = numpy.where(numpy.bitwise_and(height[0,:] >=0,height[0,:] <= hcm,numpy.isfinite(height[0,:])))
368 h_select = numpy.where(numpy.bitwise_and(height[0,:] >=0,height[0,:] <= hcm,numpy.isfinite(height[0,:])))
369 else:
369 else:
370 h_select = numpy.where(numpy.bitwise_and(height[1,:] > hcm,height[1,:] < 20,numpy.isfinite(height[1,:])))
370 h_select = numpy.where(numpy.bitwise_and(height[1,:] > hcm,height[1,:] < 20,numpy.isfinite(height[1,:])))
371
371
372
372
373 ht = h_select[0]
373 ht = h_select[0]
374
374
375 o_height = numpy.hstack((o_height,height[im,ht]))
375 o_height = numpy.hstack((o_height,height[im,ht]))
376 o_zon_aver = numpy.hstack((o_zon_aver,zon_aver[ht,im]))
376 o_zon_aver = numpy.hstack((o_zon_aver,zon_aver[ht,im]))
377 o_mer_aver = numpy.hstack((o_mer_aver,mer_aver[ht,im]))
377 o_mer_aver = numpy.hstack((o_mer_aver,mer_aver[ht,im]))
378 o_ver_aver = numpy.hstack((o_ver_aver,ver_aver[ht,im]))
378 o_ver_aver = numpy.hstack((o_ver_aver,ver_aver[ht,im]))
379
379
380 data_fHeigths_List.append(o_height)
380 data_fHeigths_List.append(o_height)
381 data_fZonal_List.append(o_zon_aver)
381 data_fZonal_List.append(o_zon_aver)
382 data_fMeridional_List.append(o_mer_aver)
382 data_fMeridional_List.append(o_mer_aver)
383 data_fVertical_List.append(o_ver_aver)
383 data_fVertical_List.append(o_ver_aver)
384
384
385
385
386 else:
386 else:
387 h_select = numpy.where(numpy.bitwise_and(height[0,:] <= hcm,numpy.isfinite(height[0,:])))
387 h_select = numpy.where(numpy.bitwise_and(height[0,:] <= hcm,numpy.isfinite(height[0,:])))
388 ht = h_select[0]
388 ht = h_select[0]
389 o_height = numpy.hstack((o_height,height[im,ht]))
389 o_height = numpy.hstack((o_height,height[im,ht]))
390 o_zon_aver = numpy.hstack((o_zon_aver,zon_aver[ht,im]))
390 o_zon_aver = numpy.hstack((o_zon_aver,zon_aver[ht,im]))
391 o_mer_aver = numpy.hstack((o_mer_aver,mer_aver[ht,im]))
391 o_mer_aver = numpy.hstack((o_mer_aver,mer_aver[ht,im]))
392 o_ver_aver = numpy.hstack((o_ver_aver,ver_aver[ht,im]))
392 o_ver_aver = numpy.hstack((o_ver_aver,ver_aver[ht,im]))
393
393
394 data_fHeigths_List.append(o_height)
394 data_fHeigths_List.append(o_height)
395 data_fZonal_List.append(o_zon_aver)
395 data_fZonal_List.append(o_zon_aver)
396 data_fMeridional_List.append(o_mer_aver)
396 data_fMeridional_List.append(o_mer_aver)
397 data_fVertical_List.append(o_ver_aver)
397 data_fVertical_List.append(o_ver_aver)
398
398
399
399
400 return startDTList, data_fHeigths_List, data_fZonal_List, data_fMeridional_List, data_fVertical_List
400 return startDTList, data_fHeigths_List, data_fZonal_List, data_fMeridional_List, data_fVertical_List
401
401
402
402
403 No newline at end of file
403
General Comments 0
You need to be logged in to leave comments. Login now