##// END OF EJS Templates
Merge branch 'v3.0-WR' of http://intranet.igp.gob.pe:8082/schain into v3.0-WR
Juan C. Espinoza -
r1531:a8e64958fb6a merge
parent child
Show More
@@ -1,33 +1,51
1 1 # Importra librerías
2 2 import os
3 3 #import imageio
4 4 import imageio.v2 as imageio
5 5 # Ubicación de la base de datos
6 '''
6 7 #path = '/home/soporte/Downloads/plots_C0/S/6/6_evento0/'
7 8 path= '/home/soporte/Documents/EVENTO/HYO_PM@2022-06-09T15-05-12/plotsC0_PL_R60.0km/V_PPI_EL_8.0CH0/'
8
9 9 path= '/home/soporte/Documents/EVENTO/HYO_PM@2022-05-31T12-00-17/plotsC0_PL_R60.0km_removeDC/V_PPI_EL_8.0CH0/'
10
11 10 path= '/home/soporte/Documents/EVENTO/HYO_PM@2022-06-09T15-05-12/plotsC0_FD_PL_R15.0km_removeDC/Z_PPI_EL_2.0CH0/'
12 11 path= '/home/soporte/Documents/EVENTO/HYO_PM@2022-06-09T15-05-12/plotsC0_FD_PL_R15.0km/Z_PPI_EL_8.0CH0/'
13
14
15 12 path= '/home/soporte/Documents/EVENTO/HYO_PM@2022-06-09T15-05-12/paramC0N36.0/2022-06-09T18-00-00/W_E.2CH0/'
13 '''
14 #--------------------------------------------------------------------------------------------------------
15 #pulse pair
16 #path = '/home/soporte/Documents/EVENTO/HYO_PM@2022-06-09T15-05-12/plotsC0_PL_R28.0km/Z_PPI_EL_2.0CH0/'
17 #path = '/home/soporte/Documents/EVENTO/HYO_PM@2022-06-09T15-05-12/plotsC0_PL_R28.0km/W_PPI_WE_L8.0CH0/'
18 #---------------------------------------------------------------------------------------------------------
19 #pulse pair REMOVEDC
20 #path = '/home/soporte/Documents/EVENTO/HYO_PM@2022-06-09T15-05-12/plotsC0_PL_R28.0km_removeDC/V_PPI_EL_2.0CH0/'
21 #path = '/home/soporte/Documents/EVENTO/HYO_PM@2022-06-09T15-05-12/plotsC0_PL_R28.0km_removeDC/W_PPI_WE_L8.0CH0/'
22 #----------------------------------------------------------------------------------------------------------
23 #Spectro
24 #path = '/home/soporte/Documents/EVENTO/HYO_PM@2022-06-09T15-05-12/plotsC0_FD_PL_R28.0km/V_PPI_EL_2.0CH0/'
25 path = '/home/soporte/Documents/EVENTO/HYO_PM@2022-06-09T15-05-12/plotsC0_FD_PL_R28.0km/W_PPI_WE_L8.0CH0/'
26
27 #---------------------------------------------------------------------------------------------------------
28 #Spectro REMOVEDC
29 path = '/home/soporte/Documents/EVENTO/HYO_PM@2022-06-09T15-05-12/plotsC0_FD_PL_R28.0km_removeDC/V_PPI_EL_8.0CH0/'
30 #path = '/home/soporte/Documents/EVENTO/HYO_PM@2022-06-09T15-05-12/plotsC0_FD_PL_R28.0km_removeDC/W_PPI_WE_L2.0CH0/'
31
16 32
17 33
18 34 archivos = sorted(os.listdir(path))
19 35 img_array = []
20 36
21 37 #Leer todos los archivos formato imagen desde path
22 38 for x in range(0, len(archivos)):
23 39 nomArchivo = archivos[x]
24 40 dirArchivo = path + str(nomArchivo)
25 41
26 42 #Asignar a variable leer_imagen, el nombre de cada imagen
27 43 leer_imagen = imageio.imread(dirArchivo)
28 44
29 45 # añadir imágenes al arreglo img_array
30 46 img_array.append(leer_imagen)
31 47
32 48 #Guardar Gif
33 imageio.mimwrite('evento_PL_PP_2_28KM_W', img_array, 'GIF', duration=0.5)
49 #imageio.mimwrite('/home/soporte/Documents/giff_28km/evento_PL_PP_8_28KM_W', img_array, 'GIF', duration=0.5)
50 #imageio.mimwrite('/home/soporte/Documents/giff_28km/evento_PL_PP_8_28KM_W_reDC', img_array, 'GIF', duration=0.5)
51 imageio.mimwrite('/home/soporte/Documents/giff_28kmFD/evento_PL_FD_8_28KM_V_reDC', img_array, 'GIF', duration=0.5)
@@ -1,375 +1,375
1 1 # SOPHY PROC script
2 2 import os, sys, json, argparse
3 3 import datetime
4 4 import time
5 5
6 6 PATH = '/DATA_RM/DATA'
7 7 PATH = '/media/jespinoza/Elements'
8 8 PATH = '/media/jespinoza/data/SOPHY'
9 9 PATH = '/home/soporte/Documents/EVENTO'
10 10
11 11 PARAM = {
12 12 'S': {'zmin': -45, 'zmax': -25, 'colormap': 'jet', 'label': 'Power', 'wrname': 'power','cb_label': 'dBm', 'ch':0},
13 13 'SNR': {'zmin': -40, 'zmax': -20, 'colormap': 'jet', 'label': 'SNR', 'wrname': 'snr','cb_label': 'dB', 'ch':0},
14 14 'V': {'zmin': -12, 'zmax': 12, 'colormap': 'sophy_v', 'label': 'Velocity', 'wrname': 'velocity', 'cb_label': 'm/s', 'ch':0},
15 15 'R': {'zmin': 0, 'zmax': 1, 'colormap': 'jet', 'label': 'RhoHV', 'wrname':'rhoHV', 'cb_label': '*', 'ch':0},
16 16 'P': {'zmin': -180,'zmax': 180,'colormap': 'RdBu_r', 'label': 'PhiDP', 'wrname':'phiDP' , 'cb_label': 'º', 'ch':0},
17 17 'D': {'zmin': -30, 'zmax': 80, 'colormap': 'sophy_r','label': 'ZDR','wrname':'differential_reflectivity' , 'cb_label': 'dBz','ch':0},
18 18 'Z': {'zmin': -30, 'zmax': 80, 'colormap': 'sophy_r','label': 'Reflectivity ', 'wrname':'reflectivity', 'cb_label': 'dBz','ch':0},
19 'W': {'zmin': 0, 'zmax': 15, 'colormap': 'sophy_w','label': 'Spectral Width', 'wrname':'spectral_width', 'cb_label': 'hz', 'ch':0}
19 'W': {'zmin': 0, 'zmax': 15, 'colormap': 'sophy_w','label': 'Spectral Width', 'wrname':'spectral_width', 'cb_label': 'm/s', 'ch':0}
20 20 }
21 21
22 22 def max_index(r, sample_rate, ipp):
23 23
24 24 return int(sample_rate*ipp*1e6 * r / 60) + int(sample_rate*ipp*1e6 * 1.2 / 60)
25 25
26 26 def main(args):
27 27
28 28 experiment = args.experiment
29 29 fp = open(os.path.join(PATH, experiment, 'experiment.conf'))
30 30 conf = json.loads(fp.read())
31 31
32 32 ipp_km = conf['usrp_tx']['ipp']
33 33 ipp = ipp_km * 2 /300000
34 34 sample_rate = conf['usrp_rx']['sample_rate']
35 35 axis = ['0' if x=='elevation' else '1' for x in conf['pedestal']['axis']] # AZIMUTH 1 ELEVACION 0
36 36 speed_axis = conf['pedestal']['speed']
37 37 steps = conf['pedestal']['table']
38 38 time_offset = args.time_offset
39 39 parameters = args.parameters
40 40 start_date = experiment.split('@')[1].split('T')[0].replace('-', '/')
41 41 end_date = start_date
42 42 if args.start_time:
43 43 start_time = args.start_time
44 44 else:
45 45 start_time = experiment.split('@')[1].split('T')[1].replace('-', ':')
46 46 end_time = '23:59:59'
47 47 N = int(1/(speed_axis[0]*ipp)) # 1 GRADO DE RESOLUCION
48 48 path = os.path.join(PATH, experiment, 'rawdata')
49 49 path_ped = os.path.join(PATH, experiment, 'position')
50 50 path_plots = os.path.join(PATH, experiment, 'plotsC0N'+str(args.range))
51 51 path_save = os.path.join(PATH, experiment, 'paramC0N'+str(args.range))
52 52 RMIX = 1.62
53 53 H0 = -1.68
54 54 MASK = 0.3
55 55
56 56 from schainpy.controller import Project
57 57
58 58 project = Project()
59 59 project.setup(id='1', name='Sophy', description='sophy proc')
60 60
61 61 reader = project.addReadUnit(datatype='DigitalRFReader',
62 62 path=path,
63 63 startDate=start_date,
64 64 endDate=end_date,
65 65 startTime=start_time,
66 66 endTime=end_time,
67 67 delay=30,
68 68 online=args.online,
69 69 walk=1,
70 70 ippKm = ipp_km,
71 71 getByBlock = 1,
72 72 nProfileBlocks = N,
73 73 )
74 74
75 75 if not conf['usrp_tx']['enable_2']: # One Pulse
76 76 voltage = project.addProcUnit(datatype='VoltageProc', inputId=reader.getId())
77 77
78 78 if conf['usrp_tx']['code_type_1'] != 'None':
79 79 codes = [ c.strip() for c in conf['usrp_tx']['code_1'].split(',')]
80 80 code = []
81 81 for c in codes:
82 82 code.append([int(x) for x in c])
83 83 op = voltage.addOperation(name='Decoder', optype='other')
84 84 op.addParameter(name='code', value=code)
85 85 op.addParameter(name='nCode', value=len(code), format='int')
86 86 op.addParameter(name='nBaud', value=len(code[0]), format='int')
87 87
88 88 op = voltage.addOperation(name='CohInt', optype='other') #Minimo integrar 2 perfiles por ser codigo complementario
89 89 op.addParameter(name='n', value=len(code), format='int')
90 90 ncode = len(code)
91 91 else:
92 92 ncode = 1
93 93 code = ['0']
94 94
95 95 op = voltage.addOperation(name='setH0')
96 96 op.addParameter(name='h0', value=H0)
97 97
98 98 if args.range > 0:
99 99 op = voltage.addOperation(name='selectHeights')
100 100 op.addParameter(name='minIndex', value='0', format='int')
101 101 op.addParameter(name='maxIndex', value=max_index(args.range, sample_rate, ipp), format='int')
102 102
103 103 op = voltage.addOperation(name='PulsePair_vRF', optype='other')
104 104 op.addParameter(name='n', value=int(N)/ncode, format='int')
105 105 #op.addParameter(name='removeDC', value=1, format='int')
106 106
107 107
108 108 proc = project.addProcUnit(datatype='ParametersProc', inputId=voltage.getId())
109 109
110 110 opObj10 = proc.addOperation(name="WeatherRadar")
111 111 opObj10.addParameter(name='tauW',value=(1e-6/sample_rate)*len(code[0]))
112 112 opObj10.addParameter(name='Pt',value=((1e-6/sample_rate)*len(code[0])/ipp)*200)
113 113
114 114 op = proc.addOperation(name='PedestalInformation')
115 115 op.addParameter(name='path', value=path_ped, format='str')
116 116 op.addParameter(name='interval', value='0.04')
117 117 op.addParameter(name='time_offset', value=time_offset)
118 118 op.addParameter(name='mode', value='PPI')
119 119
120 120 for param in parameters:
121 121 op = proc.addOperation(name='Block360')
122 122 op.addParameter(name='runNextOp', value=True)
123 123
124 124 op= proc.addOperation(name='WeatherParamsPlot')
125 125 if args.save: op.addParameter(name='save', value=path_plots, format='str')
126 126 op.addParameter(name='save_period', value=-1)
127 127 op.addParameter(name='show', value=args.show)
128 128 op.addParameter(name='channels', value='1,')
129 129 op.addParameter(name='zmin', value=PARAM[param]['zmin'])
130 130 op.addParameter(name='zmax', value=PARAM[param]['zmax'])
131 131 op.addParameter(name='attr_data', value=param, format='str')
132 132 op.addParameter(name='labels', value=[PARAM[param]['label']])
133 133 op.addParameter(name='save_code', value=param)
134 134 op.addParameter(name='cb_label', value=PARAM[param]['cb_label'])
135 135 op.addParameter(name='colormap', value=PARAM[param]['colormap'])
136 136 op.addParameter(name='bgcolor', value='black')
137 137 if MASK: op.addParameter(name='mask', value=MASK, format='float')
138 138 if args.server:
139 139 op.addParameter(name='server', value='0.0.0.0:4444')
140 140 op.addParameter(name='exp_code', value='400')
141 141
142 142 desc = {
143 143 'Data': {
144 144 param: PARAM[param]['wrname'],
145 145 'utctime': 'time'
146 146 },
147 147 'Metadata': {
148 148 'heightList': 'range',
149 149 'data_azi': 'azimuth',
150 150 'data_ele': 'elevation',
151 151 'mode_op': 'scan_type',
152 152 'h0': 'range_correction',
153 153 }
154 154 }
155 155
156 156 if args.save:
157 157 opObj10 = proc.addOperation(name='HDFWriter')
158 158 writer.addParameter(name='path', value=path_save, format='str')
159 159 writer.addParameter(name='Reset', value=True)
160 160 writer.addParameter(name='setType', value='weather')
161 161 writer.addParameter(name='description', value=json.dumps(desc))
162 162 writer.addParameter(name='blocksPerFile', value='1',format='int')
163 163 writer.addParameter(name='metadataList', value='heightList,data_azi,data_ele,mode_op,latitude,longitude,altitude,heading,radar_name,institution,contact,h0,range_unit')
164 164 writer.addParameter(name='dataList', value='{},utctime'.format(param))
165 165 writer.addParameter(name='mask', value=MASK, format='float')
166 166 # meta
167 167 writer.addParameter(name='latitude', value='-12.040436')
168 168 writer.addParameter(name='longitude', value='-75.295893')
169 169 writer.addParameter(name='altitude', value='3379.2147')
170 170 writer.addParameter(name='heading', value='0')
171 171 writer.addParameter(name='radar_name', value='SOPHy')
172 172 writer.addParameter(name='institution', value='IGP')
173 173 writer.addParameter(name='contact', value='dscipion@igp.gob.pe')
174 174 writer.addParameter(name='created_by', value='Signal Chain (https://pypi.org/project/schainpy/)')
175 175 writer.addParameter(name='range_unit', value='km')
176 176
177 177 else: #Two pulses
178 178
179 179 voltage1 = project.addProcUnit(datatype='VoltageProc', inputId=reader.getId())
180 180
181 181 op = voltage1.addOperation(name='ProfileSelector')
182 182 op.addParameter(name='profileRangeList', value='0,{}'.format(conf['usrp_tx']['repetitions_1']-1))
183 183
184 184 if conf['usrp_tx']['code_type_1'] != 'None':
185 185 codes = [ c.strip() for c in conf['usrp_tx']['code_1'].split(',')]
186 186 code = []
187 187 for c in codes:
188 188 code.append([int(x) for x in c])
189 189 op = voltage1.addOperation(name='Decoder', optype='other')
190 190 op.addParameter(name='code', value=code)
191 191 op.addParameter(name='nCode', value=len(code), format='int')
192 192 op.addParameter(name='nBaud', value=len(code[0]), format='int')
193 193 else:
194 194 code = ['0']
195 195
196 196 op = voltage1.addOperation(name='CohInt', optype='other') #Minimo integrar 2 perfiles por ser codigo complementario
197 197 op.addParameter(name='n', value=2, format='int')
198 198
199 199 if args.range > 0:
200 200 op = voltage1.addOperation(name='selectHeights')
201 201 op.addParameter(name='minIndex', value='0', format='int')
202 202 op.addParameter(name='maxIndex', value=max_index(RMIX, sample_rate, ipp), format='int')
203 203
204 204 op = voltage1.addOperation(name='setH0')
205 205 op.addParameter(name='h0', value=H0, format='float')
206 206
207 207 op = voltage1.addOperation(name='PulsePair_vRF', optype='other')
208 208 op.addParameter(name='n', value=int(conf['usrp_tx']['repetitions_1'])/2, format='int')
209 209 #op.addParameter(name='removeDC', value=1, format='int')
210 210
211 211
212 212 proc1 = project.addProcUnit(datatype='ParametersProc', inputId=voltage1.getId())
213 213 proc1.addParameter(name='runNextUnit', value=True)
214 214
215 215 opObj10 = proc1.addOperation(name="WeatherRadar")
216 216 opObj10.addParameter(name='tauW',value=(1e-6/sample_rate)*len(code[0]))
217 217 opObj10.addParameter(name='Pt',value=((1e-6/sample_rate)*len(code[0])/ipp)*200)
218 218
219 219 op = proc1.addOperation(name='PedestalInformation')
220 220 op.addParameter(name='path', value=path_ped, format='str')
221 221 op.addParameter(name='interval', value='0.04')
222 222 op.addParameter(name='time_offset', value=time_offset)
223 223 op.addParameter(name='mode', value='PPI')
224 224
225 225 op = proc1.addOperation(name='Block360')
226 226 op.addParameter(name='attr_data', value='data_param')
227 227 op.addParameter(name='runNextOp', value=True)
228 228
229 229
230 230 voltage2 = project.addProcUnit(datatype='VoltageProc', inputId=reader.getId())
231 231
232 232 op = voltage2.addOperation(name='ProfileSelector')
233 233 op.addParameter(name='profileRangeList', value='{},{}'.format(conf['usrp_tx']['repetitions_1'], conf['usrp_tx']['repetitions_1']+conf['usrp_tx']['repetitions_2']-1))
234 234
235 235 if conf['usrp_tx']['code_type_2']:
236 236 codes = [ c.strip() for c in conf['usrp_tx']['code_2'].split(',')]
237 237 code = []
238 238 for c in codes:
239 239 code.append([int(x) for x in c])
240 240 op = voltage2.addOperation(name='Decoder', optype='other')
241 241 op.addParameter(name='code', value=code)
242 242 op.addParameter(name='nCode', value=len(code), format='int')
243 243 op.addParameter(name='nBaud', value=len(code[0]), format='int')
244 244
245 245 op = voltage2.addOperation(name='CohInt', optype='other') #Minimo integrar 2 perfiles por ser codigo complementario
246 246 op.addParameter(name='n', value=len(code), format='int')
247 247 ncode = len(code)
248 248 else:
249 249 ncode = 1
250 250
251 251 if args.range > 0:
252 252 op = voltage2.addOperation(name='selectHeights')
253 253 op.addParameter(name='minIndex', value=max_index(RMIX, sample_rate, ipp), format='int')
254 254 op.addParameter(name='maxIndex', value=max_index(args.range, sample_rate, ipp), format='int')
255 255
256 256 op = voltage2.addOperation(name='setH0')
257 257 op.addParameter(name='h0', value=H0, format='float')
258 258
259 259 op = voltage2.addOperation(name='PulsePair_vRF', optype='other')
260 260 op.addParameter(name='n', value=int(conf['usrp_tx']['repetitions_2'])/ncode, format='int')
261 261 #op.addParameter(name='removeDC', value=1, format='int')
262 262
263 263
264 264 proc2 = project.addProcUnit(datatype='ParametersProc', inputId=voltage2.getId())
265 265 proc2.addParameter(name='runNextUnit', value=True)
266 266
267 267 opObj10 = proc2.addOperation(name="WeatherRadar")
268 268 opObj10.addParameter(name='tauW',value=(1e-6/sample_rate)*len(code[0]))
269 269 opObj10.addParameter(name='Pt',value=((1e-6/sample_rate)*len(code[0])/ipp)*200)
270 270
271 271 op = proc2.addOperation(name='PedestalInformation')
272 272 op.addParameter(name='path', value=path_ped, format='str')
273 273 op.addParameter(name='interval', value='0.04')
274 274 op.addParameter(name='time_offset', value=time_offset)
275 275 op.addParameter(name='mode', value='PPI')
276 276
277 277 op = proc2.addOperation(name='Block360')
278 278 op.addParameter(name='attr_data', value='data_param')
279 279 op.addParameter(name='runNextOp', value=True)
280 280
281 281 merge = project.addProcUnit(datatype='MergeProc', inputId=[proc1.getId(), proc2.getId()])
282 282 merge.addParameter(name='attr_data', value='data_param')
283 283 merge.addParameter(name='mode', value='7') #RM
284 284
285 285
286 286 for param in parameters:
287 287
288 288 if args.plot:
289 289 op= merge.addOperation(name='WeatherParamsPlot')
290 290 if args.save:
291 291 op.addParameter(name='save', value=path_plots, format='str')
292 292 op.addParameter(name='save_period', value=-1)
293 293 op.addParameter(name='show', value=args.show)
294 294 op.addParameter(name='channels', value='0,')
295 295 op.addParameter(name='zmin', value=PARAM[param]['zmin'], format='int')
296 296 op.addParameter(name='zmax', value=PARAM[param]['zmax'], format='int')
297 297 op.addParameter(name='attr_data', value=param, format='str')
298 298 op.addParameter(name='labels', value=[PARAM[param]['label']])
299 299 op.addParameter(name='save_code', value=param)
300 300 op.addParameter(name='cb_label', value=PARAM[param]['cb_label'])
301 301 op.addParameter(name='colormap', value=PARAM[param]['colormap'])
302 302 op.addParameter(name='bgcolor', value='black')
303 303 if MASK: op.addParameter(name='mask', value=MASK, format='float')
304 304 if args.server:
305 305 op.addParameter(name='server', value='0.0.0.0:4444')
306 306 op.addParameter(name='exp_code', value='400')
307 307
308 308 desc = {
309 309 'Data': {
310 310 'data_param': {PARAM[param]['wrname']: ['H', 'V']},
311 311 'utctime': 'time'
312 312 },
313 313 'Metadata': {
314 314 'heightList': 'range',
315 315 'data_azi': 'azimuth',
316 316 'data_ele': 'elevation',
317 317 'mode_op': 'scan_type',
318 318 'h0': 'range_correction',
319 319 }
320 320 }
321 321
322 322 if args.save:
323 323 writer = merge.addOperation(name='HDFWriter')
324 324 writer.addParameter(name='path', value=path_save, format='str')
325 325 writer.addParameter(name='Reset', value=True)
326 326 writer.addParameter(name='setType', value='weather')
327 327 writer.addParameter(name='description', value=json.dumps(desc))
328 328 writer.addParameter(name='blocksPerFile', value='1',format='int')
329 329 writer.addParameter(name='metadataList', value='heightList,data_azi,data_ele,mode_op,latitude,longitude,altitude,heading,radar_name,institution,contact,h0,range_unit')
330 330 writer.addParameter(name='dataList', value='data_param,utctime')
331 331 writer.addParameter(name='weather_var', value=param)
332 332 writer.addParameter(name='mask', value=MASK, format='float')
333 333 # meta
334 334 writer.addParameter(name='latitude', value='-12.040436')
335 335 writer.addParameter(name='longitude', value='-75.295893')
336 336 writer.addParameter(name='altitude', value='3379.2147')
337 337 writer.addParameter(name='heading', value='0')
338 338 writer.addParameter(name='radar_name', value='SOPHy')
339 339 writer.addParameter(name='institution', value='IGP')
340 340 writer.addParameter(name='contact', value='dscipion@igp.gob.pe')
341 341 writer.addParameter(name='created_by', value='Signal Chain (https://pypi.org/project/schainpy/)')
342 342 writer.addParameter(name='range_unit', value='km')
343 343
344 344 project.start()
345 345
346 346 if __name__ == '__main__':
347 347
348 348 parser = argparse.ArgumentParser(description='Script to process SOPHy data.')
349 349 parser.add_argument('experiment',
350 350 help='Experiment name')
351 351 parser.add_argument('--parameters', nargs='*', default=['S'],
352 352 help='Variables to process: P, Z, V')
353 353 parser.add_argument('--time_offset', default=0,
354 354 help='Fix time offset')
355 355 parser.add_argument('--range', default=0, type=float,
356 356 help='Max range to plot')
357 357 parser.add_argument('--save', action='store_true',
358 358 help='Create output files')
359 359 parser.add_argument('--plot', action='store_true',
360 360 help='Create plot files')
361 361 parser.add_argument('--show', action='store_true',
362 362 help='Show matplotlib plot.')
363 363 parser.add_argument('--online', action='store_true',
364 364 help='Set online mode.')
365 365 parser.add_argument('--server', action='store_true',
366 366 help='Send to realtime')
367 367 parser.add_argument('--start_time', default='',
368 368 help='Set start time.')
369 369
370 370
371 371 args = parser.parse_args()
372 372
373 373 main(args)
374 374
375 # python sophy_proc.py HYO_PM@2022-06-09T15-05-12 --parameters V --plot --save --show --range 36
375 # python sophy_proc.py HYO_PM@2022-06-09T15-05-12 --parameters V --plot --save --show --range 36S
@@ -1,198 +1,225
1 1 # SOPHY PROC script
2 2 import os, sys, json, argparse
3 3 import datetime
4 4 import time
5 5
6 6 PATH = '/DATA_RM/DATA'
7 7 # PATH = '/Users/jespinoza/workspace/data/'
8 8 #PATH = '/home/roberto/DATA/data_WR_RHI/RHI'
9 9 PATH = '/home/soporte/Downloads/data_WR_RHI'
10 10 PATH = '/home/soporte/Documents/EVENTO/'
11 11
12 12
13 13 PARAM = {
14 14 'S': {'zmin': -45, 'zmax': -15, 'colormap': 'jet', 'label': 'Power', 'wrname': 'power','cb_label': 'dBm', 'ch':0},
15 'SNR': {'zmin': -40, 'zmax': -20, 'colormap': 'jet', 'label': 'SNR', 'wrname': 'snr','cb_label': 'dB', 'ch':0},
15 16 #'V': {'name': 'dataPP_DOP', 'zmin': -10, 'zmax': 10, 'colormap': 'sophy_v', 'label': 'Velocity', 'wrname': 'velocity', 'cb_label': 'm/s', 'ch':0},
16 17 'V': {'zmin': -10, 'zmax': 10, 'colormap': 'sophy_v', 'label': 'Velocity', 'wrname': 'velocity', 'cb_label': 'm/s', 'ch':0},
17 18 'R': {'zmin': 0, 'zmax': 1, 'colormap': 'jet', 'label': 'RhoHV', 'wrname':'rhoHV', 'cb_label': '*', 'ch':0},
18 19 'P': {'zmin': -180,'zmax': 180,'colormap': 'RdBu_r', 'label': 'PhiDP', 'wrname':'phiDP' , 'cb_label': 'º', 'ch':0},
19 'D': {'zmin': -20, 'zmax': 80, 'colormap': 'gist_ncar','label': 'ZDR','wrname':'differential_reflectivity' , 'cb_label': 'dBz','ch':0},
20 'D': {'zmin': -30, 'zmax': 80, 'colormap': 'gist_ncar','label': 'ZDR','wrname':'differential_reflectivity' , 'cb_label': 'dBz','ch':0},
20 21 'Z': {'zmin': -30, 'zmax': 80, 'colormap': 'sophy_r','label': 'Reflectivity', 'wrname':'reflectivity', 'cb_label': 'dBz','ch':1},
21 'W': {'zmin': 0, 'zmax': 12, 'colormap': 'sophy_w','label': 'Spectral Width', 'wrname':'spectral_width', 'cb_label': 'hz', 'ch':1}
22 'W': {'zmin': 0, 'zmax': 12, 'colormap': 'sophy_w','label': 'Spectral Width', 'wrname':'spectral_width', 'cb_label': 'm/s', 'ch':1}
22 23 }
23 24
24 25 def max_index(r, sample_rate, ipp):
25 26
26 27 return int(sample_rate*ipp*1e6 * r / 60) + int(sample_rate*ipp*1e6 * 1.2 / 60)
27 28
28 29 def main(args):
29 30
30 31 experiment = args.experiment
31 32 fp = open(os.path.join(PATH, experiment, 'experiment.conf'))
32 33 conf = json.loads(fp.read())
33 34
34 35 ipp_km = conf['usrp_tx']['ipp']
35 36 ipp = ipp_km * 2 /300000
36 37 sample_rate = conf['usrp_rx']['sample_rate']
37 38 axis = ['0' if x=='elevation' else '1' for x in conf['pedestal']['axis']] # AZIMUTH 1 ELEVACION 0
38 39 speed_axis = conf['pedestal']['speed']
39 40 steps = conf['pedestal']['table']
40 41 time_offset = args.time_offset
41 42 parameters = args.parameters
42 43 start_date = experiment.split('@')[1].split('T')[0].replace('-', '/')
43 44 end_date = start_date
44 45 if args.start_time:
45 46 start_time = args.start_time
46 47 else:
47 48 start_time = experiment.split('@')[1].split('T')[1].replace('-', ':')
48 49 #start_time = '16:15:00'
49 50 end_time = '23:59:59'
50 51 N = int(1/(speed_axis[0]*ipp)) # 1 GRADO DE RESOLUCION
51 52 path = os.path.join(PATH, experiment, 'rawdata')
52 53 path_ped = os.path.join(PATH, experiment, 'position')
53 54 path_plots = os.path.join(PATH, experiment, 'plotsC0_PL_R'+str(args.range)+'km_removeDC')
54 55 path_save = os.path.join(PATH, experiment, 'paramC0_PL_R'+str(args.range)+'km_removeDC')
55 56 #path_plots = os.path.join(PATH, experiment, 'plotsC0_PL_R'+str(args.range)+'km')
56 57 #path_save = os.path.join(PATH, experiment, 'paramC0_PL_R'+str(args.range)+'km')
57 58 RMIX = 1.62
59 H0=-1.68
60 MASK =0.3
58 61
59 62 from schainpy.controller import Project
60 63
61 64 project = Project()
62 65 project.setup(id='1', name='Sophy', description='sophy proc')
63 66
64 67 reader = project.addReadUnit(datatype='DigitalRFReader',
65 68 path=path,
66 69 startDate=start_date,
67 70 endDate=end_date,
68 71 startTime=start_time,
69 72 endTime=end_time,
70 73 delay=30,
71 74 online=args.online,
72 75 walk=1,
73 76 ippKm = ipp_km,
74 77 getByBlock = 1,
75 78 nProfileBlocks = N,
76 79 )
77 80 voltage = project.addProcUnit(datatype='VoltageProc', inputId=reader.getId())
78 81
79 82 op = voltage.addOperation(name='ProfileSelector')
80 83 op.addParameter(name='profileRangeList', value='{},{}'.format(conf['usrp_tx']['repetitions_1'], conf['usrp_tx']['repetitions_1']+conf['usrp_tx']['repetitions_2']-1))
81 84
82 85
83
84
85 86 if conf['usrp_tx']['code_type_2']:
86 87 codes = [ c.strip() for c in conf['usrp_tx']['code_2'].split(',')]
87 88 code = []
88 89 for c in codes:
89 90 code.append([int(x) for x in c])
90 91 op = voltage.addOperation(name='Decoder', optype='other')
91 92 op.addParameter(name='code', value=code)
92 93 op.addParameter(name='nCode', value=len(code), format='int')
93 94 op.addParameter(name='nBaud', value=len(code[0]), format='int')
94 95
95 96 op = voltage.addOperation(name='CohInt', optype='other') #Minimo integrar 2 perfiles por ser codigo complementario
96 97 op.addParameter(name='n', value=len(code), format='int')
97 98 ncode = len(code)
99 print("codes",codes)
100 print("nCode",len(code))
98 101 else:
99 102 ncode = 1
100 103
101 op = voltage.addOperation(name='setH0')
102 op.addParameter(name='h0', value='-1.68')
103 104
104 105 if args.range > 0:
105 106 op = voltage.addOperation(name='selectHeights')
106 107 op.addParameter(name='minIndex', value='0', format='int')
107 108 op.addParameter(name='maxIndex', value=max_index(args.range, sample_rate, ipp), format='int')
108 109
110
111 op = voltage.addOperation(name='setH0')
112 op.addParameter(name='h0', value='-1.68' , format='float')
113
114
109 115 op = voltage.addOperation(name='PulsePair_vRF', optype='other')
110 116 op.addParameter(name='n', value=int(conf['usrp_tx']['repetitions_2'])/2, format='int')
111 #op.addParameter(name='removeDC', value=1, format='int')
117 op.addParameter(name='removeDC', value=1, format='int')
112 118
113 119
114 120 proc = project.addProcUnit(datatype='ParametersProc', inputId=voltage.getId())
121 proc.addParameter(name='runNextUnit', value=True)
115 122
116 123 opObj10 = proc.addOperation(name="WeatherRadar")
117 opObj10.addParameter(name='variableList',value='Reflectividad,VelocidadRadial,AnchoEspectral')
124 #opObj10.addParameter(name='variableList',value='Reflectividad,VelocidadRadial,AnchoEspectral')
118 125 opObj10.addParameter(name='tauW',value=(1e-6/sample_rate)*len(code[0]))
119 126 opObj10.addParameter(name='Pt',value=((1e-6/sample_rate)*len(code[0])/ipp)*200)
120 127
121 128 # {"latitude": -12.0404828587, "longitude": -75.2147483647, "altitude": 3379.2147483647}
122 129
123 130 op = proc.addOperation(name='PedestalInformation')
124 131 op.addParameter(name='path', value=path_ped, format='str')
125 132 op.addParameter(name='interval', value='0.04')
126 133 op.addParameter(name='time_offset', value=time_offset)
127 134 #op.addParameter(name='az_offset', value=-26.2)
128 135 op.addParameter(name='mode', value='PPI')
129 136
130 137 for param in parameters:
131 138 op = proc.addOperation(name='Block360')
132 139 op.addParameter(name='attr_data', value='data_param')
133 140 op.addParameter(name='runNextOp', value=True)
134 141
135 142 op= proc.addOperation(name='WeatherParamsPlot')
136 143 if args.save: op.addParameter(name='save', value=path_plots, format='str')
137 144 op.addParameter(name='save_period', value=-1)
138 145 op.addParameter(name='show', value=args.show)
139 146 op.addParameter(name='channels', value='0,')
140 147 op.addParameter(name='zmin', value=PARAM[param]['zmin'])
141 148 op.addParameter(name='zmax', value=PARAM[param]['zmax'])
142 149 op.addParameter(name='attr_data', value=param, format='str')
143 150 op.addParameter(name='labels', value=[PARAM[param]['label']])
144 151 op.addParameter(name='save_code', value=param)
145 152 op.addParameter(name='cb_label', value=PARAM[param]['cb_label'])
146 153 op.addParameter(name='colormap', value=PARAM[param]['colormap'])
147 154 op.addParameter(name='bgcolor',value='black')
148 op.addParameter(name='snr_threshold',value=0.4)
155 if MASK: op.addParameter(name='mask', value=MASK, format='float')
156 if args.server:
157 op.addParameter(name='server', value='0.0.0.0:4444')
158 op.addParameter(name='exp_code', value='400')
149 159
150 160 desc = {
151 161 'Data': {
152 'data_param': PARAM[param]['wrname'],
162 'data_param': {PARAM[param]['wrname']: ['H', 'V']},
153 163 'utctime': 'time'
154 164 },
155 165 'Metadata': {
156 166 'heightList': 'range',
157 167 'data_azi': 'azimuth',
158 168 'data_ele': 'elevation',
169 'mode_op': 'scan_type',
170 'h0': 'range_correction',
159 171 }
160 172 }
161 173
162 174 if args.save:
163 opObj10 = proc.addOperation(name='HDFWriter')
164 opObj10.addParameter(name='path', value=path_save, format='str')
165 opObj10.addParameter(name='Reset', value=True)
166 opObj10.addParameter(name='setType', value='weather')
167 opObj10.addParameter(name=' description', value=json.dumps(desc))
168 opObj10.addParameter(name='blocksPerFile', value='1',format='int')
169 opObj10.addParameter(name='metadataList', value='heightList,data_azi,data_ele')
170 opObj10.addParameter(name='dataList', value='{},utctime'.format(PARAM[param]['name']))
175 writer = proc.addOperation(name='HDFWriter')
176 writer.addParameter(name='path', value=path_save, format='str')
177 writer.addParameter(name='Reset', value=True)
178 writer.addParameter(name='setType', value='weather')
179 writer.addParameter(name='description', value=json.dumps(desc))
180 writer.addParameter(name='blocksPerFile', value='1',format='int')
181 writer.addParameter(name='metadataList', value='heightList,data_azi,data_ele,mode_op,latitude,longitude,altitude,heading,radar_name,institution,contact,h0,range_unit')
182 writer.addParameter(name='dataList', value='data_param,utctime')
183 writer.addParameter(name='weather_var', value=param)
184 writer.addParameter(name='mask', value=MASK, format='float')
185 # meta
186 writer.addParameter(name='latitude', value='-12.040436')
187 writer.addParameter(name='longitude', value='-75.295893')
188 writer.addParameter(name='altitude', value='3379.2147')
189 writer.addParameter(name='heading', value='0')
190 writer.addParameter(name='radar_name', value='SOPHy')
191 writer.addParameter(name='institution', value='IGP')
192 writer.addParameter(name='contact', value='dscipion@igp.gob.pe')
193 writer.addParameter(name='created_by', value='Signal Chain (https://pypi.org/project/schainpy/)')
194 writer.addParameter(name='range_unit', value='km')
195
171 196 project.start()
172 197
173 198
174 199
175 200 if __name__ == '__main__':
176 201
177 202 parser = argparse.ArgumentParser(description='Script to process SOPHy data.')
178 203 parser.add_argument('experiment',
179 204 help='Experiment name')
180 205 parser.add_argument('--parameters', nargs='*', default=['S'],
181 206 help='Variables to process: P, Z, V')
182 207 parser.add_argument('--time_offset', default=0,
183 208 help='Fix time offset')
184 209 parser.add_argument('--range', default=0, type=float,
185 210 help='Max range to plot')
186 211 parser.add_argument('--save', action='store_true',
187 212 help='Create output files')
188 213 parser.add_argument('--show', action='store_true',
189 214 help='Show matplotlib plot.')
190 215 parser.add_argument('--online', action='store_true',
191 216 help='Set online mode.')
217 parser.add_argument('--server', action='store_true',
218 help='Send to realtime')
192 219 parser.add_argument('--start_time', default='',
193 220 help='Set start time.')
194 221
195 222
196 223 args = parser.parse_args()
197 224
198 225 main(args)
@@ -1,207 +1,225
1 1 # SOPHY PROC script
2 2 import os, sys, json, argparse
3 3 import datetime
4 4 import time
5 5
6 6 PATH = '/DATA_RM/DATA'
7 7 # PATH = '/Users/jespinoza/workspace/data/'
8 8 #PATH = '/home/roberto/DATA/data_WR_RHI/RHI'
9 9 PATH = '/home/soporte/Downloads/data_WR_RHI'
10 10 PATH = '/home/soporte/Documents/EVENTO/'
11 11
12 12
13 13 PARAM = {
14 14 'S': {'name': 'dataPP_POWER', 'zmin': -45, 'zmax': -15, 'colormap': 'jet', 'label': 'Power', 'wrname': 'power','cb_label': 'dBm', 'ch':0},
15 15 'S2': {'name': 'data_pow', 'zmin': -45, 'zmax': -15, 'colormap': 'jet', 'label': 'Power', 'wrname': 'power','cb_label': 'dBm', 'ch':0},
16 16 #'V': {'name': 'dataPP_DOP', 'zmin': -10, 'zmax': 10, 'colormap': 'sophy_v', 'label': 'Velocity', 'wrname': 'velocity', 'cb_label': 'm/s', 'ch':0},
17 17 'V': {'name': 'velRadial_V', 'zmin': -10, 'zmax': 10, 'colormap': 'sophy_v', 'label': 'Velocity', 'wrname': 'velocity', 'cb_label': 'm/s', 'ch':0},
18 18 'R': {'name': 'RhoHV_R', 'zmin': 0, 'zmax': 1, 'colormap': 'jet', 'label': 'RhoHV', 'wrname':'rhoHV', 'cb_label': '*', 'ch':0},
19 19 'P': {'name': 'PhiD_P', 'zmin': -180,'zmax': 180,'colormap': 'RdBu_r', 'label': 'PhiDP', 'wrname':'phiDP' , 'cb_label': 'º', 'ch':0},
20 20 'D': {'name': 'Zdb_D', 'zmin': -20, 'zmax': 80, 'colormap': 'gist_ncar','label': 'ZDR','wrname':'differential_reflectivity' , 'cb_label': 'dBz','ch':0},
21 21 'Z': {'name': 'Zdb', 'zmin': -30, 'zmax': 80, 'colormap': 'sophy_r','label': 'Reflectivity', 'wrname':'reflectivity', 'cb_label': 'dBz','ch':1},
22 'W': {'name': 'Sigmav_W', 'zmin': 0, 'zmax': 12, 'colormap': 'sophy_w','label': 'Spectral Width', 'wrname':'spectral_width', 'cb_label': 'hz', 'ch':1}
22 'W': {'name': 'Sigmav_W', 'zmin': 0, 'zmax': 12, 'colormap': 'sophy_w','label': 'Spectral Width', 'wrname':'spectral_width', 'cb_label': 'm/s', 'ch':1}
23 23 }
24 24
25 25 def max_index(r, sample_rate, ipp):
26 26
27 27 return int(sample_rate*ipp*1e6 * r / 60) + int(sample_rate*ipp*1e6 * 1.2 / 60)
28 28
29 29 def main(args):
30 30
31 31 experiment = args.experiment
32 32 fp = open(os.path.join(PATH, experiment, 'experiment.conf'))
33 33 conf = json.loads(fp.read())
34 34
35 35 ipp_km = conf['usrp_tx']['ipp']
36 36 ipp = ipp_km * 2 /300000
37 37 sample_rate = conf['usrp_rx']['sample_rate']
38 38 axis = ['0' if x=='elevation' else '1' for x in conf['pedestal']['axis']] # AZIMUTH 1 ELEVACION 0
39 39 speed_axis = conf['pedestal']['speed']
40 40 steps = conf['pedestal']['table']
41 41 time_offset = args.time_offset
42 42 parameters = args.parameters
43 43 start_date = experiment.split('@')[1].split('T')[0].replace('-', '/')
44 44 end_date = start_date
45 45 if args.start_time:
46 46 start_time = args.start_time
47 47 else:
48 48 start_time = experiment.split('@')[1].split('T')[1].replace('-', ':')
49 49 #start_time = '16:15:00'
50 50 end_time = '23:59:59'
51 51 N = int(1/(speed_axis[0]*ipp)) # 1 GRADO DE RESOLUCION
52 52 path = os.path.join(PATH, experiment, 'rawdata')
53 53 path_ped = os.path.join(PATH, experiment, 'position')
54 #path_plots = os.path.join(PATH, experiment, 'plotsC0_FD_PL_R'+str(args.range)+'km_removeDC')
55 #path_save = os.path.join(PATH, experiment, 'paramC0_FD_PL_R'+str(args.range)+'km_removeDC')
56 path_plots = os.path.join(PATH, experiment, 'plotsC0_FD_PL_R'+str(args.range)+'km')
57 path_save = os.path.join(PATH, experiment, 'paramC0_FD_PL_R'+str(args.range)+'km')
54 path_plots = os.path.join(PATH, experiment, 'plotsC0_FD_PL_R'+str(args.range)+'km_removeDC')
55 path_save = os.path.join(PATH, experiment, 'paramC0_FD_PL_R'+str(args.range)+'km_removeDC')
56 #path_plots = os.path.join(PATH, experiment, 'plotsC0_FD_PL_R'+str(args.range)+'km')
57 #path_save = os.path.join(PATH, experiment, 'paramC0_FD_PL_R'+str(args.range)+'km')
58 58 RMIX = 1.62
59 MASK =0.3
59 60
60 61 from schainpy.controller import Project
61 62
62 63 project = Project()
63 64 project.setup(id='1', name='Sophy', description='sophy proc')
64 65
65 66 reader = project.addReadUnit(datatype='DigitalRFReader',
66 67 path=path,
67 68 startDate=start_date,
68 69 endDate=end_date,
69 70 startTime=start_time,
70 71 endTime=end_time,
71 72 delay=30,
72 73 online=args.online,
73 74 walk=1,
74 75 ippKm = ipp_km,
75 76 getByBlock = 1,
76 77 nProfileBlocks = N,
77 78 )
78 79 voltage = project.addProcUnit(datatype='VoltageProc', inputId=reader.getId())
79 80
80 81 op = voltage.addOperation(name='setAttribute')
81 82 op.addParameter(name='frequency', value='9345.e6', format='float')
82 83
83 84 op = voltage.addOperation(name='ProfileSelector')
84 85 op.addParameter(name='profileRangeList', value='{},{}'.format(conf['usrp_tx']['repetitions_1'], conf['usrp_tx']['repetitions_1']+conf['usrp_tx']['repetitions_2']-1))
85 86
86 87 if conf['usrp_tx']['code_type_2']:
87 88 codes = [ c.strip() for c in conf['usrp_tx']['code_2'].split(',')]
88 89 code = []
89 90 for c in codes:
90 91 code.append([int(x) for x in c])
91 92 op = voltage.addOperation(name='Decoder', optype='other')
92 93 op.addParameter(name='code', value=code)
93 94 op.addParameter(name='nCode', value=len(code), format='int')
94 95 op.addParameter(name='nBaud', value=len(code[0]), format='int')
95 96
96 97 op = voltage.addOperation(name='CohInt', optype='other') #Minimo integrar 2 perfiles por ser codigo complementario
97 98 op.addParameter(name='n', value=len(code), format='int')
98 99 ncode = len(code)
99 100 else:
100 101 ncode = 1
101 102
102 103 op = voltage.addOperation(name='setH0')
103 104 op.addParameter(name='h0', value='-1.68')
104 105
105 106 if args.range > 0:
106 107 op = voltage.addOperation(name='selectHeights')
107 108 op.addParameter(name='minIndex', value='0', format='int')
108 109 op.addParameter(name='maxIndex', value=max_index(args.range, sample_rate, ipp), format='int')
109 110
110 111 #---------------------------------------NEW PROCESSING -----------------------------------------------------
111 112 procB = project.addProcUnit(datatype='SpectraProc', inputId=voltage.getId())
112 113 procB.addParameter(name='nFFTPoints', value=int(conf['usrp_tx']['repetitions_2'])/2, format='int')
113 114 procB.addParameter(name='nProfiles', value=int(conf['usrp_tx']['repetitions_2'])/2, format='int')
114 115
115 #opObj11 = procB.addOperation(name='removeDC')
116 #opObj11.addParameter(name='mode', value=2)
116 opObj11 = procB.addOperation(name='removeDC')
117 opObj11.addParameter(name='mode', value=2)
117 118
118 119 proc= project.addProcUnit(datatype='ParametersProc',inputId=procB.getId())
119 120
120 121 opObj10 = proc.addOperation(name='SpectralMoments')
121 122 opObj10.addParameter(name='wradar',value=True)
122 123
123 124 #---------------------------------------NEW PROCESSING -----------------------------------------------------
124 125
125 126 opObj10 = proc.addOperation(name="WeatherRadar")
126 opObj10.addParameter(name='variableList',value='Reflectividad,VelocidadRadial,AnchoEspectral')
127 #opObj10.addParameter(name='variableList',value='Reflectividad,VelocidadRadial,AnchoEspectral')
127 128 opObj10.addParameter(name='tauW',value=(1e-6/sample_rate)*len(code[0]))
128 129 opObj10.addParameter(name='Pt',value=((1e-6/sample_rate)*len(code[0])/ipp)*200)
129 130
130 131 # {"latitude": -12.0404828587, "longitude": -75.2147483647, "altitude": 3379.2147483647}
131 132
132 133 op = proc.addOperation(name='PedestalInformation')
133 134 op.addParameter(name='path', value=path_ped, format='str')
134 135 op.addParameter(name='interval', value='0.04')
135 136 op.addParameter(name='time_offset', value=time_offset)
136 137 #op.addParameter(name='az_offset', value=-26.2)
137 138 op.addParameter(name='mode', value='PPI')
138 139
139 140 for param in parameters:
140 141 op = proc.addOperation(name='Block360')
141 142 op.addParameter(name='attr_data', value='data_param')
142 143 op.addParameter(name='runNextOp', value=True)
143 144
144 145 op= proc.addOperation(name='WeatherParamsPlot')
145 146 if args.save: op.addParameter(name='save', value=path_plots, format='str')
146 147 op.addParameter(name='save_period', value=-1)
147 148 op.addParameter(name='show', value=args.show)
148 149 op.addParameter(name='channels', value='0,')
149 150 op.addParameter(name='zmin', value=PARAM[param]['zmin'])
150 151 op.addParameter(name='zmax', value=PARAM[param]['zmax'])
151 152 op.addParameter(name='attr_data', value=param, format='str')
152 153 op.addParameter(name='labels', value=[PARAM[param]['label']])
153 154 op.addParameter(name='save_code', value=param)
154 155 op.addParameter(name='cb_label', value=PARAM[param]['cb_label'])
155 156 op.addParameter(name='colormap', value=PARAM[param]['colormap'])
156 157 op.addParameter(name='bgcolor',value='black')
157 158 op.addParameter(name='snr_threshold',value=0.4)
158 159
160 if MASK: op.addParameter(name='mask', value=MASK, format='float')
161
159 162 desc = {
160 163 'Data': {
161 'data_param': PARAM[param]['wrname'],
164 'data_param': {PARAM[param]['wrname']: ['H', 'V']},
162 165 'utctime': 'time'
163 166 },
164 167 'Metadata': {
165 168 'heightList': 'range',
166 169 'data_azi': 'azimuth',
167 170 'data_ele': 'elevation',
171 'mode_op': 'scan_type',
172 'h0': 'range_correction',
168 173 }
169 174 }
170 175
171 176 if args.save:
172 opObj10 = proc.addOperation(name='HDFWriter')
173 opObj10.addParameter(name='path', value=path_save, format='str')
174 opObj10.addParameter(name='Reset', value=True)
175 opObj10.addParameter(name='setType', value='weather')
176 opObj10.addParameter(name=' description', value=json.dumps(desc))
177 opObj10.addParameter(name='blocksPerFile', value='1',format='int')
178 opObj10.addParameter(name='metadataList', value='heightList,data_azi,data_ele')
179 opObj10.addParameter(name='dataList', value='{},utctime'.format(PARAM[param]['name']))
177 writer = proc.addOperation(name='HDFWriter')
178 writer.addParameter(name='path', value=path_save, format='str')
179 writer.addParameter(name='Reset', value=True)
180 writer.addParameter(name='setType', value='weather')
181 writer.addParameter(name='description', value=json.dumps(desc))
182 writer.addParameter(name='blocksPerFile', value='1',format='int')
183 writer.addParameter(name='metadataList', value='heightList,data_azi,data_ele,mode_op,latitude,longitude,altitude,heading,radar_name,institution,contact,h0,range_unit')
184 writer.addParameter(name='dataList', value='data_param,utctime')
185 writer.addParameter(name='weather_var', value=param)
186 writer.addParameter(name='mask', value=MASK, format='float')
187 # meta
188 writer.addParameter(name='latitude', value='-12.040436')
189 writer.addParameter(name='longitude', value='-75.295893')
190 writer.addParameter(name='altitude', value='3379.2147')
191 writer.addParameter(name='heading', value='0')
192 writer.addParameter(name='radar_name', value='SOPHy')
193 writer.addParameter(name='institution', value='IGP')
194 writer.addParameter(name='contact', value='dscipion@igp.gob.pe')
195 writer.addParameter(name='created_by', value='Signal Chain (https://pypi.org/project/schainpy/)')
196 writer.addParameter(name='range_unit', value='km')
197
180 198 project.start()
181 199
182 200
183 201
184 202 if __name__ == '__main__':
185 203
186 204 parser = argparse.ArgumentParser(description='Script to process SOPHy data.')
187 205 parser.add_argument('experiment',
188 206 help='Experiment name')
189 207 parser.add_argument('--parameters', nargs='*', default=['S'],
190 208 help='Variables to process: P, Z, V')
191 209 parser.add_argument('--time_offset', default=0,
192 210 help='Fix time offset')
193 211 parser.add_argument('--range', default=0, type=float,
194 212 help='Max range to plot')
195 213 parser.add_argument('--save', action='store_true',
196 214 help='Create output files')
197 215 parser.add_argument('--show', action='store_true',
198 216 help='Show matplotlib plot.')
199 217 parser.add_argument('--online', action='store_true',
200 218 help='Set online mode.')
201 219 parser.add_argument('--start_time', default='',
202 220 help='Set start time.')
203 221
204 222
205 223 args = parser.parse_args()
206 224
207 main(args)
225 main(args) # Operator#ñ8
@@ -1,229 +1,229
1 1 #!/usr/bin/env python
2 2 # Ing. AVP
3 3 # 22/06/2022
4 4 # ARCHIVO DE LECTURA y PLOT
5 5 import matplotlib.pyplot as pl
6 6 import matplotlib
7 7 import wradlib
8 8 import numpy
9 9 import warnings
10 10 import argparse
11 11 from wradlib.io import read_generic_hdf5
12 12 from wradlib.util import get_wradlib_data_file
13 13 from plotting_codes import sophy_cb_tables
14 14 import os,time
15 15
16 16 for name, cb_table in sophy_cb_tables:
17 17 ncmap = matplotlib.colors.ListedColormap(cb_table, name=name)
18 18 matplotlib.pyplot.register_cmap(cmap=ncmap)
19 19 #LINUX bash: export WRADLIB_DATA=/path/to/wradlib-data
20 20 #example
21 21 #export WRADLIB_DATA="/home/soporte/Documents/EVENTO/HYO_PM@2022-06-09T15-05-12/paramC0N36.0/2022-06-09T18-00-00/"
22 22 warnings.filterwarnings('ignore')
23 23 PARAM = {
24 24 'S': {'var': 'power','vmin': -45, 'vmax': -15, 'cmap': 'jet', 'label': 'Power','unit': 'dBm'},
25 25 'V': {'var': 'velocity', 'vmin': -10, 'vmax': 10 , 'cmap': 'sophy_v', 'label': 'Velocity','unit': 'm/s'},
26 26 'Z': {'var': 'reflectivity','vmin': -30, 'vmax': 80 , 'cmap': 'sophy_r','label': 'Reflectivity','unit': 'dBZ'},
27 'W': {'var': 'spectral_width', 'vmin': 0 , 'vmax': 12 , 'cmap': 'sophy_w','label': 'Spectral Width','unit': 'hz'}
27 'W': {'var': 'spectral_width', 'vmin': 0 , 'vmax': 12 , 'cmap': 'sophy_w','label': 'Spectral Width','unit': 'm/s'}
28 28 }
29 29 class Readsophy():
30 30 def __init__(self):
31 31 self.list_file = None
32 32 self.grado = None
33 33 self.variable = None
34 34 self.save = None
35 35 self.range = None
36 36
37 37 def read_files(self,path_file,grado=None, variable=None):
38 38 filter= "_E"+str(grado)+".0_"+variable
39 39 validFilelist = []
40 40 fileList= os.listdir(path_file)
41 41 for thisFile in fileList:
42 42 if (os.path.splitext(thisFile)[0][-7:] != filter):
43 43 #print("s_:",os.path.splitext(thisFile)[0][-7:])
44 44 continue
45 45 validFilelist.append(thisFile)
46 46 validFilelist.sort()
47 47 return validFilelist
48 48
49 49 def setup(self, path_file,grado,range,variable,save):
50 50 self.path_file = path_file
51 51 self.range = range
52 52 self.grado = grado
53 53 self.variable = variable
54 54 self.save = save
55 55 self.list_file = self.read_files(path_file=self.path_file,grado=self.grado, variable=self.variable)
56 56
57 57 def selectHeights(self,heightList,minHei,maxHei):
58 58
59 59 if minHei and maxHei:
60 60 if (minHei < heightList[0]):
61 61 minHei = heightList[0]
62 62 if (maxHei > heightList[-1]):
63 63 maxHei = heightList[-1]
64 64 minIndex = 0
65 65 maxIndex = 0
66 66 heights = heightList
67 67
68 68 inda = numpy.where(heights >= minHei)
69 69 indb = numpy.where(heights <= maxHei)
70 70
71 71 try:
72 72 minIndex = inda[0][0]
73 73 except:
74 74 minIndex = 0
75 75
76 76 try:
77 77 maxIndex = indb[0][-1]
78 78 except:
79 79 maxIndex = len(heights)
80 80
81 81 new_heightList= self.selectHeightsByIndex(heightList=heightList,minIndex=minIndex, maxIndex=maxIndex)
82 82
83 83 return new_heightList, minIndex,maxIndex
84 84
85 85 def selectHeightsByIndex(self,heightList,minIndex, maxIndex):
86 86
87 87 if (minIndex < 0) or (minIndex > maxIndex):
88 88 raise ValueError("Height index range (%d,%d) is not valid" % (minIndex, maxIndex))
89 89
90 90 if (maxIndex >= len(heightList)):
91 91 maxIndex = len(heightList)
92 92
93 93 new_h = heightList[minIndex:maxIndex]
94 94 return new_h
95 95
96 96 def readAttributes(self,obj,variable):
97 97 var = PARAM[variable]['var']
98 98 unit = PARAM[variable]['unit']
99 99 cmap = PARAM[variable]['cmap']
100 100 vmin = PARAM[variable]['vmin']
101 101 vmax = PARAM[variable]['vmax']
102 102 label = PARAM[variable]['label']
103 103 var_ = 'Data/'+var+'/H'
104 104 data_arr = numpy.array(obj[var_]['data']) # data
105 105 utc_time = numpy.array(obj['Data/time']['data'])
106 106 data_azi = numpy.array(obj['Metadata/azimuth']['data']) # th
107 107 data_ele = numpy.array(obj["Metadata/elevation"]['data'])
108 108 heightList = numpy.array(obj["Metadata/range"]['data']) # r
109 109
110 110 return data_arr, utc_time, data_azi,data_ele, heightList,unit,cmap,vmin,vmax,label
111 111
112 112 def run(self):
113 113 count= 0
114 114 len_files = len(self.list_file)
115 115
116 116 for thisFile in self.list_file:
117 117 count = count +1
118 118 fullpathfile = self.path_file + thisFile
119 119 filename = get_wradlib_data_file(fullpathfile)
120 120 test_hdf5 = read_generic_hdf5(filename)
121 121
122 122 # LECTURA
123 123 data_arr, utc_time, data_azi,data_ele, heightList,unit,cmap,vmin,vmax,label = self.readAttributes(obj= test_hdf5,variable=self.variable)
124 124
125 125 if self.range==0:
126 126 self.range == heightList[-1]
127 127 new_heightList,minIndex,maxIndex = self.selectHeights(heightList,0.06,self.range)
128 128
129 129 # TIEMPO
130 130 my_time = time.strftime('%Y-%m-%d %H:%M:%S', time.localtime(utc_time[0]))
131 131 time_save = time.strftime('%Y%m%d_%H%M%S',time.localtime(utc_time[0]))
132 132
133 133 # PLOT DATA WITH ANNOTATION
134 134 if count ==1:
135 135 fig = pl.figure(figsize=(10,8))
136 136 cgax, pm = wradlib.vis.plot_ppi(data_arr[:,minIndex:maxIndex],r=new_heightList,az=data_azi,rf=1,fig=fig, ax=111,proj='cg',cmap=cmap,vmin=vmin, vmax=vmax)
137 137 caax = cgax.parasites[0]
138 138 title = 'Simple PPI'+"-"+ my_time+" E."+self.grado
139 139 t = pl.title(title, fontsize=12,y=1.05)
140 140 cbar = pl.gcf().colorbar(pm, pad=0.075)
141 141 pl.text(1.0, 1.05, 'azimuth', transform=caax.transAxes, va='bottom',ha='right')
142 142 cbar.set_label(label+'[' + unit + ']')
143 143 gh = cgax.get_grid_helper()
144 144 else:
145 145 cgax, pm = wradlib.vis.plot_ppi(data_arr[:,minIndex:maxIndex],r=new_heightList,az=data_azi,rf=1,fig=fig, ax=111,proj='cg',cmap=cmap,vmin=vmin, vmax=vmax)
146 146 caax = cgax.parasites[0]
147 147 title = 'Simple PPI'+"-"+my_time+" E."+self.grado
148 148 t = pl.title(title, fontsize=12,y=1.05)
149 149 cbar = pl.gcf().colorbar(pm, pad=0.075)
150 150 pl.text(1.0, 1.05, 'azimuth', transform=caax.transAxes, va='bottom',ha='right')
151 151 cbar.set_label(label+'[' + unit + ']')
152 152 gh = cgax.get_grid_helper()
153 153 if self.save == 1:
154 154 if count ==1:
155 155 filename = "SOPHY"+"_"+time_save+"_"+"E."+self.grado+"_"+self.variable+".png"
156 156 dir =self.variable+"_"+"E."+self.grado+"CH0/"
157 157 filesavepath = os.path.join(self.path_file,dir)
158 158 try:
159 159 os.mkdir(filesavepath)
160 160 except:
161 161 pass
162 162 else:
163 163 filename = "SOPHY"+"_"+time_save+"_"+"E."+self.grado+"_"+self.variable+".png"
164 164 pl.savefig(filesavepath+filename)
165 165
166 166 pl.pause(1)
167 167 pl.clf()
168 168 if count==len_files:
169 169 pl.close()
170 170 pl.show()
171 171
172 172 PATH = "/home/soporte/Documents/EVENTO/HYO_PM@2022-06-09T15-05-12/paramC0N36.0/2022-06-09T18-00-00/"
173 173 #PATH = "/home/soporte/Documents/EVENTO/HYO_PM@2022-06-09T15-05-12/paramC0N36.0/2022-06-09T19-00-00/"
174 174
175 175 #PATH = "/home/soporte/Documents/EVENTO/HYO_PM@2022-05-31T12-00-17/paramC0N36.0/2022-05-31T16-00-00/"
176 176
177 177 def main(args):
178 178 grado = args.grado
179 179 parameters = args.parameters
180 180 save = args.save
181 181 range = args.range
182 182 obj = Readsophy()
183 183 for param in parameters:
184 184 obj.setup(path_file = PATH,grado = grado,range=range, variable=param,save=int(save))
185 185 obj.run()
186 186
187 187 if __name__ == '__main__':
188 188
189 189 parser = argparse.ArgumentParser(description='Script to process SOPHy data.')
190 190 parser.add_argument('--parameters', nargs='*', default=['S'],
191 191 help='Variables to process: P, Z, V ,W')
192 192 parser.add_argument('--grado', default=2,
193 193 help='Angle in Elev to plot')
194 194 parser.add_argument('--save', default=0,
195 195 help='Save plot')
196 196 parser.add_argument('--range', default=0, type=float,
197 197 help='Max range to plot')
198 198 args = parser.parse_args()
199 199
200 200 main(args)
201 201
202 202 #python sophy_proc_rev006.py --parameters Z --grado 8 --save 1 --range 28
203 203 '''
204 204 def read_and_overview(filename):
205 205 """Read HDF5 using read_generic_hdf5 and print upper level dictionary keys
206 206 """
207 207 test = read_generic_hdf5(filename)
208 208 print("\nPrint keys for file %s" % os.path.basename(filename))
209 209 for key in test.keys():
210 210 print("\t%s" % key)
211 211 return test
212 212
213 213 file__ = "/home/soporte/Documents/EVENTO/HYO_PM@2022-06-09T15-05-12/paramC0_FD_PL_R15.0km/2022-06-09T18-00-00/SOPHY_20220609_180229_E2.0_Z.hdf5"
214 214
215 215
216 216 filename = get_wradlib_data_file(file__)
217 217
218 218 print("filename:\n",filename )
219 219
220 220 test= read_and_overview(filename)
221 221 # ANADIR INFORMACION
222 222 # informacion de los pulsos de TX
223 223 # informacion de los ruidos
224 224 # informacion de los SNR ¿?
225 225 # Aumentar la amplitud de la USRP
226 226 LAST_UPDATE
227 227 ---- Noise
228 228 ---- Mapas
229 229 '''
General Comments 0
You need to be logged in to leave comments. Login now