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