##// END OF EJS Templates
Merge branch 'v3.0-WR' of http://intranet.igp.gob.pe:8082/schain into v3.0-WR
Juan C. Espinoza -
r1480:ab3389190ff4 merge
parent child
Show More

The requested changes are too big and content was truncated. Show full diff

@@ -0,0 +1,362
1 #!python
2 '''
3 '''
4
5 import os, sys
6 import datetime
7 import time
8
9 #path = os.path.dirname(os.getcwd())
10 #path = os.path.dirname(path)
11 #sys.path.insert(0, path)
12
13 from schainpy.controller import Project
14
15 desc = "USRP_test"
16 filename = "USRP_processing.xml"
17 controllerObj = Project()
18 controllerObj.setup(id = '191', name='Test_USRP', description=desc)
19
20 ############## USED TO PLOT IQ VOLTAGE, POWER AND SPECTRA #############
21
22 #######################################################################
23 ######PATH DE LECTURA, ESCRITURA, GRAFICOS Y ENVIO WEB#################
24 #######################################################################
25 #path = '/media/data/data/vientos/57.2063km/echoes/NCO_Woodman'
26 #path = '/DATA_RM/TEST_INTEGRACION'
27 #path = '/DATA_RM/PRUEBA_USRP_RP'
28 #path = '/DATA_RM/PRUEBA_USRP_RP'
29
30 path = '/DATA_RM/TEST_2M'
31 path = '/DATA_RM/TEST_2M_UD'
32 path = '/DATA_RM/2MHZ17022022'
33 path = '/DATA_RM/10MHZTEST/'
34 path = '/DATA_RM/10MHZDRONE/'
35
36
37 path= '/home/soporte/TEST_500mVPP'
38 path= '/home/soporte/TEST_1VPP+500mVDC'
39 path = '/home/soporte/TEST_500mVPP+500mVDC'
40 path = '/home/soporte/TEST_1.5VPP'
41 path = '/home/soporte/TEST_2VPP'
42 path= '/home/soporte/TEST_1VPP'
43 path = '/home/soporte/Documents/HUANCAYO/TEST_HYO_PM@2022-05-14T11-28-19/rawdata'
44
45 #HYO_PM@2022-05-28T00-00-17
46 path = '/DATA_RM/DATA/HYO_PM@2022-05-28T00-00-17/rawdata'
47
48 #figpath = '/home/soporte/Pictures/TEST_RP_0001'
49 #figpath = '/home/soporte/Pictures/TEST_RP_6000'
50 figpath = '/home/soporte/Pictures/USRP_TEST_2M'
51 figpath = '/home/soporte/Pictures/USRP_TEST_2M_UD'
52 figpath = '/home/soporte/Pictures/10MHZDRONE'
53 figpath = '/home/soporte/Pictures/500mVPP'
54 figpath = '/home/soporte/Pictures/1VPP+500mVDC'
55 figpath = '/home/soporte/Pictures/TEST_500mVPP+500mVDC'
56 figpath = '/home/soporte/Pictures/TEST_1.5VPP'
57 figpath = '/home/soporte/Pictures/TEST_2VPP'
58 figpath = '/home/soporte/Pictures/TEST_1VPP'
59
60
61
62
63 #remotefolder = "/home/wmaster/graficos"
64 #######################################################################
65 ################# RANGO DE PLOTEO######################################
66 #######################################################################
67 dBmin = '-60'#'-20'
68 dBmax = '-5'#'-85'
69 xmin = '0'
70 xmax ='24'
71 ymin = '0'
72 ymax = '10'
73 #######################################################################
74 ########################FECHA##########################################
75 #######################################################################
76 str = datetime.date.today()
77 today = str.strftime("%Y/%m/%d")
78 str2 = str - datetime.timedelta(days=1)
79 yesterday = str2.strftime("%Y/%m/%d")
80 #######################################################################
81 ######################## UNIDAD DE LECTURA#############################
82 #######################################################################
83 readUnitConfObj = controllerObj.addReadUnit(datatype='DigitalRFReader',
84 path=path,
85 startDate="2022/05/28",#today,
86 endDate="2022/05/28",#today,
87 startTime='00:00:00',# inicio libre
88 #startTime='00:00:00',
89 endTime='23:59:59',
90 delay=0,
91 #set=0,
92 online=0,
93 walk=1,
94 ippKm = 60)
95
96 opObj11 = readUnitConfObj.addOperation(name='printInfo')
97 #opObj11 = readUnitConfObj.addOperation(name='printNumberOfBlock')
98 #######################################################################
99 ################ OPERACIONES DOMINIO DEL TIEMPO########################
100 #######################################################################
101
102 procUnitConfObjA = controllerObj.addProcUnit(datatype='VoltageProc', inputId=readUnitConfObj.getId())
103
104 op3 = procUnitConfObjA.addOperation(name='ProfileSelector', optype='other')
105 op3.addParameter(name='profileRangeList', value='1,123')
106
107
108 code=[[1]]
109
110 opObj11 = procUnitConfObjA.addOperation(name='Decoder', optype='other')
111 opObj11.addParameter(name='code', value=code)
112 opObj11.addParameter(name='nCode', value='1', format='int')
113 opObj11.addParameter(name='nBaud', value='1', format='int')
114
115
116 '''
117 op3 = procUnitConfObjA.addOperation(name='ProfileSelector', optype='other')
118 op3.addParameter(name='profileRangeList', value='122,249')
119 code8=[[1,1,1,0,1,1,0,1],[1,1,1,0,0,0,1,0]]
120
121 opObj11 = procUnitConfObjA.addOperation(name='Decoder', optype='other')
122 opObj11.addParameter(name='code', value=code8)
123 opObj11.addParameter(name='nCode', value='2', format='int')
124 opObj11.addParameter(name='nBaud', value='8', format='int')
125 '''
126 op = procUnitConfObjA.addOperation(name='CohInt', optype='other') #Minimo integrar 2 perfiles por ser codigo complementario
127 op.addParameter(name='n', value=2, format='int')
128
129
130 '''
131
132 # OJO SCOPE
133 opObj10 = procUnitConfObjA.addOperation(name='ScopePlot', optype='external')
134 opObj10.addParameter(name='id', value='10', format='int')
135 opObj10.addParameter(name='xmin', value='0', format='int')
136 opObj10.addParameter(name='xmax', value='60', format='int')
137 opObj10.addParameter(name='type', value='iq')
138 #opObj10.addParameter(name='ymin', value='-0.20000', format='int')
139 #opObj10.addParameter(name='ymax', value='0.20000', format='int')
140 opObj10.addParameter(name='save', value=figpath, format='str')
141 opObj10.addParameter(name='save_period', value=1, format='int')
142 '''
143 '''
144 opObj11 = procUnitConfObjA.addOperation(name='selectHeights')
145 opObj11.addParameter(name='minIndex', value='1', format='int')
146 # opObj11.addParameter(name='maxIndex', value='10000', format='int')
147 opObj11.addParameter(name='maxIndex', value='200', format='int')
148 '''
149 #
150 # codigo64='1,1,1,0,1,1,0,1,1,1,1,0,0,0,1,0,1,1,1,0,1,1,0,1,0,0,0,1,1,1,0,1,1,1,1,0,1,1,0,1,1,1,1,0,0,0,1,0,0,0,0,1,0,0,1,0,1,1,1,0,0,0,1,0,'+\
151 # '1,1,1,0,1,1,0,1,1,1,1,0,0,0,1,0,1,1,1,0,1,1,0,1,0,0,0,1,1,1,0,1,0,0,0,1,0,0,1,0,0,0,0,1,1,1,0,1,1,1,1,0,1,1,0,1,0,0,0,1,1,1,0,1'
152
153 #opObj11 = procUnitConfObjA.addOperation(name='setRadarFrequency')
154 #opObj11.addParameter(name='frequency', value='49920000')
155
156 '''
157 opObj11 = procUnitConfObjA.addOperation(name='PulsePair', optype='other')
158 opObj11.addParameter(name='n', value='625', format='int')#10
159 opObj11.addParameter(name='removeDC', value=1, format='int')
160 '''
161
162 # Ploteo TEST
163 '''
164 opObj11 = procUnitConfObjA.addOperation(name='PulsepairPowerPlot', optype='other')
165 opObj11 = procUnitConfObjA.addOperation(name='PulsepairSignalPlot', optype='other')
166 opObj11 = procUnitConfObjA.addOperation(name='PulsepairVelocityPlot', optype='other')
167 #opObj11.addParameter(name='xmax', value=8)
168 opObj11 = procUnitConfObjA.addOperation(name='PulsepairSpecwidthPlot', optype='other')
169 '''
170 # OJO SCOPE
171 #opObj10 = procUnitConfObjA.addOperation(name='ScopePlot', optype='external')
172 #opObj10.addParameter(name='id', value='10', format='int')
173 ##opObj10.addParameter(name='xmin', value='0', format='int')
174 ##opObj10.addParameter(name='xmax', value='50', format='int')
175 #opObj10.addParameter(name='type', value='iq')
176 ##opObj10.addParameter(name='ymin', value='-5000', format='int')
177 ##opObj10.addParameter(name='ymax', value='8500', format='int')
178 #opObj11.addParameter(name='save', value=figpath, format='str')
179 #opObj11.addParameter(name='save_period', value=10, format='int')
180
181 #opObj10 = procUnitConfObjA.addOperation(name='setH0')
182 #opObj10.addParameter(name='h0', value='-5000', format='float')
183
184 #opObj11 = procUnitConfObjA.addOperation(name='filterByHeights')
185 #opObj11.addParameter(name='window', value='1', format='int')
186
187 #codigo='1,1,-1,1,1,-1,1,-1,-1,1,-1,-1,-1,1,-1,-1,-1,1,-1,-1,-1,1,1,1,1,-1,-1,-1'
188 #opObj11 = procUnitConfObjSousy.addOperation(name='Decoder', optype='other')
189 #opObj11.addParameter(name='code', value=codigo, format='floatlist')
190 #opObj11.addParameter(name='nCode', value='1', format='int')
191 #opObj11.addParameter(name='nBaud', value='28', format='int')
192
193 #opObj11 = procUnitConfObjA.addOperation(name='CohInt', optype='other')
194 #opObj11.addParameter(name='n', value='100', format='int')
195
196 #######################################################################
197 ########## OPERACIONES ParametersProc########################
198 #######################################################################
199 ###procUnitConfObjB= controllerObj.addProcUnit(datatype='ParametersProc',inputId=procUnitConfObjA.getId())
200 '''
201
202 opObj11 = procUnitConfObjA.addOperation(name='PedestalInformation')
203 opObj11.addParameter(name='path_ped', value=path_ped)
204 opObj11.addParameter(name='path_adq', value=path_adq)
205 opObj11.addParameter(name='t_Interval_p', value='0.01', format='float')
206 opObj11.addParameter(name='n_Muestras_p', value='100', format='float')
207 opObj11.addParameter(name='blocksPerfile', value='100', format='int')
208 opObj11.addParameter(name='f_a_p', value='25', format='int')
209 opObj11.addParameter(name='online', value='0', format='int')
210
211 opObj11 = procUnitConfObjA.addOperation(name='Block360')
212 opObj11.addParameter(name='n', value='40', format='int')
213
214 opObj11= procUnitConfObjA.addOperation(name='WeatherPlot',optype='other')
215 opObj11.addParameter(name='save', value=figpath)
216 opObj11.addParameter(name='save_period', value=1)
217
218 8
219 '''
220
221
222 '''
223 opObj11 = procUnitConfObjA.addOperation(name='CohInt', optype='other')
224 opObj11.addParameter(name='n', value='250', format='int')
225 '''
226 #######################################################################
227 ########## OPERACIONES DOMINIO DE LA FRECUENCIA########################
228 #######################################################################
229 '''
230 procUnitConfObjB = controllerObj.addProcUnit(datatype='SpectraProc', inputId=procUnitConfObjA.getId())
231 procUnitConfObjB.addParameter(name='nFFTPoints', value='64', format='int')
232 procUnitConfObjB.addParameter(name='nProfiles', value='64', format='int')
233 '''
234
235 procUnitConfObjB = controllerObj.addProcUnit(datatype='SpectraProc', inputId=procUnitConfObjA.getId())
236 procUnitConfObjB.addParameter(name='nFFTPoints', value='61', format='int')
237 procUnitConfObjB.addParameter(name='nProfiles', value='61', format='int')
238
239 '''
240 procUnitConfObjC = controllerObj.addProcUnit(datatype='SpectraHeisProc', inputId=procUnitConfObjA.getId())
241 #procUnitConfObjB.addParameter(name='nFFTPoints', value='64', format='int')
242 #procUnitConfObjB.addParameter(name='nProfiles', value='64', format='int')
243 opObj11 = procUnitConfObjC.addOperation(name='IncohInt4SpectraHeis', optype='other')
244 #opObj11.addParameter(name='timeInterval', value='4', format='int')
245 opObj11.addParameter(name='n', value='100', format='int')
246
247 #procUnitConfObjB.addParameter(name='pairsList', value='(0,0),(1,1),(0,1)', format='pairsList')
248
249 #opObj13 = procUnitConfObjB.addOperation(name='removeDC')
250 #opObj13.addParameter(name='mode', value='2', format='int')
251
252 #opObj11 = procUnitConfObjB.addOperation(name='IncohInt', optype='other')
253 #opObj11.addParameter(name='n', value='8', format='float')
254 #######################################################################
255 ########## PLOTEO DOMINIO DE LA FRECUENCIA#############################
256 #######################################################################
257 #----
258 '''
259 '''
260 opObj11 = procUnitConfObjC.addOperation(name='SpectraHeisPlot')
261 opObj11.addParameter(name='id', value='10', format='int')
262 opObj11.addParameter(name='wintitle', value='Spectra_Alturas', format='str')
263 #opObj11.addParameter(name='xmin', value=-100000, format='float')
264 #opObj11.addParameter(name='xmax', value=100000, format='float')
265 opObj11.addParameter(name='oneFigure', value=False,format='bool')
266 #opObj11.addParameter(name='zmin', value=-10, format='int')
267 #opObj11.addParameter(name='zmax', value=40, format='int')
268 opObj11.addParameter(name='ymin', value=10, format='int')
269 opObj11.addParameter(name='ymax', value=55, format='int')
270 opObj11.addParameter(name='grid', value=True, format
271 [Reading] 2022-05-23 12:27:32.732775: 21333 samples <> 0.010667 sec
272 ='bool')
273 #opObj11.addParameter(name='showprofile', value='1', format='int')
274 opObj11.addParameter(name='save', value=figpath, format='str')
275 #opObj11.addParameter(name='save_period', value=10, format='int')
276 '''
277 '''
278 opObj11 = procUnitConfObjC.addOperation(name='RTIHeisPlot')
279 opObj11.addParameter(name='id', value='10', format='int')
280 opObj11.addParameter(name='wintitle', value='RTI_Alturas', format='str')
281 opObj11.addParameter(name='xmin', value=11.0, format='float')
282 opObj11.addParameter(name='xmax', value=18.0, format='float')
283 opObj11.addParameter(name='zmin', value=10, format='int')
284 opObj11.addParameter(name='zmax', value=30, format='int')
285 opObj11.addParameter(name='ymin', value=5, format='int')
286 opObj11.addParameter(name='ymax', value=28, format='int')
287 opObj11.addParameter(name='showprofile', value='1', format='int')
288 opObj11.addParameter(name='save', value=figpath, format='str')
289 opObj11.addParameter(name='save_period', value=10, format='int')
290 '''
291
292 #SpectraPlot
293
294 opObj11 = procUnitConfObjB.addOperation(name='SpectraPlot', optype='external')
295 opObj11.addParameter(name='id', value='1', format='int')
296 opObj11.addParameter(name='wintitle', value='Spectra', format='str')
297 #opObj11.addParameter(name='xmin', value=-0.01, format='float')
298 #opObj11.addParameter(name='xmax', value=0.01, format='float')
299 opObj11.addParameter(name='zmin', value=dBmin, format='int')
300 opObj11.addParameter(name='zmax', value=dBmax, format='int')
301 opObj11.addParameter(name='ymin', value=ymin, format='int')
302 opObj11.addParameter(name='ymax', value=ymax, format='int')
303 opObj11.addParameter(name='showprofile', value='1', format='int')
304 opObj11.addParameter(name='save', value=figpath, format='str')
305 opObj11.addParameter(name='save_period', value=10, format='int')
306
307
308 #RTIPLOT
309 '''
310 opObj11 = procUnitConfObjB.addOperation(name='RTIPlot', optype='external')
311 opObj11.addParameter(name='id', value='2', format='int')
312 opObj11.addParameter(name='wintitle', value='RTIPlot', format='str')
313 opObj11.addParameter(name='zmin', value=dBmin, format='int')
314 opObj11.addParameter(name='zmax', value=dBmax, format='int')
315 #opObj11.addParameter(name='ymin', value=ymin, format='int')
316 #opObj11.addParameter(name='ymax', value=ymax, format='int')
317 #opObj11.addParameter(name='xmin', value=15, format='int')
318 #opObj11.addParameter(name='xmax', value=16, format='int')
319 opObj11.addParameter(name='zmin', value=dBmin, format='int')
320 opObj11.addParameter(name='zmax', value=dBmax, format='int')
321
322 opObj11.addParameter(name='showprofile', value='1', format='int')
323 opObj11.addParameter(name='save', value=figpath, format='str')
324 opObj11.addParameter(name='save_period', value=10, format='int')
325 '''
326 '''
327 # opObj11 = procUnitConfObjB.addOperation(name='CrossSpectraPlot', optype='other')
328 # opObj11.addParameter(name='id', value='3', format='int')
329 # opObj11.addParameter(name='wintitle', value='CrossSpectraPlot', format='str')
330 # opObj11.addParameter(name='ymin', value=ymin, format='int')
331 # opObj11.addParameter(name='ymax', value=ymax, format='int')
332 # opObj11.addParameter(name='phase_cmap', value='jet', format='str')
333 # opObj11.addParameter(name='zmin', value=dBmin, format='int')
334 # opObj11.addParameter(name='zmax', value=dBmax, format='int')
335 # opObj11.addParameter(name='figpath', value=figures_path, format='str')
336 # opObj11.addParameter(name='save', value=0, format='bool')
337 # opObj11.addParameter(name='pairsList', value='(0,1)', format='pairsList')
338 # #
339 # opObj11 = procUnitConfObjB.addOperation(name='CoherenceMap', optype='other')
340 # opObj11.addParameter(name='id', value='4', format='int')
341 # opObj11.addParameter(name='wintitle', value='Coherence', format='str')
342 # opObj11.addParameter(name='phase_cmap', value='jet', format='str')
343 # opObj11.addParameter(name='xmin', value=xmin, format='float')
344 # opObj11.addParameter(name='xmax', value=xmax, format='float')
345 # opObj11.addParameter(name='figpath', value=figures_path, format='str')
346 # opObj11.addParameter(name='save', value=0, format='bool')
347 # opObj11.addParameter(name='pairsList', value='(0,1)', format='pairsList')
348 #
349 '''
350 '''
351 #######################################################################
352 ############### UNIDAD DE ESCRITURA ###################################
353 #######################################################################
354 #opObj11 = procUnitConfObjB.addOperation(name='SpectraWriter', optype='other')
355 #opObj11.addParameter(name='path', value=wr_path)
356 #opObj11.addParameter(name='blocksPerFile', value='50', format='int')
357 print ("Escribiendo el archivo XML")
358 print ("Leyendo el archivo XML")
359 '''
360
361
362 controllerObj.start()
@@ -0,0 +1,378
1
2 # SOPHY PROC script
3 import os, sys, json, argparse
4 import datetime
5 import time
6
7 PATH = '/DATA_RM/DATA'
8 # PATH = '/Users/jespinoza/workspace/data/'
9 #PATH = '/home/soporte/Documents/HUANCAYO'
10 PARAM = {
11 'P': {'name': 'dataPP_POWER', 'zmin': -45, 'zmax': -25, 'colormap': 'jet', 'label': 'Power', 'wrname': 'Pow','cb_label': 'dB', 'ch':0},
12 'V': {'name': 'dataPP_DOP', 'zmin': -20, 'zmax': 20, 'colormap': 'seismic', 'label': 'Velocity', 'wrname': 'Vel', 'cb_label': 'm/s', 'ch':0},
13 'RH': {'name': 'RhoHV_R', 'zmin': 0, 'zmax': 1, 'colormap': 'jet', 'label': 'Coef.Correlacion', 'wrname':'R', 'cb_label': '*', 'ch':0},
14 'FD': {'name': 'PhiD_P', 'zmin': -180,'zmax': 180,'colormap': 'RdBu_r', 'label': 'Fase Diferencial', 'wrname':'P' , 'cb_label': 'º', 'ch':0},
15 'ZD': {'name': 'Zdb_D', 'zmin': -20, 'zmax': 60, 'colormap': 'viridis','label': 'Reflect.Diferencial','wrname':'D' , 'cb_label': 'dBz','ch':0},
16 'Z': {'name': 'Zdb', 'zmin': -20, 'zmax': 70, 'colormap': 'gist_ncar','label': 'Reflectividad', 'wrname':'Z', 'cb_label': 'dBz','ch':1},
17 'W': {'name': 'Sigmav_W', 'zmin': 0, 'zmax':5, 'colormap': 'viridis','label': 'AnchoEspectral', 'wrname':'S', 'cb_label': 'hz', 'ch':1}
18 }
19
20 #
21 def max_index(r, sample_rate, ipp):
22
23 return int(sample_rate*ipp*1e6 * r / 60) + int(sample_rate*ipp*1e6 * 1.2 / 60)
24
25 def main(args):
26
27 experiment = args.experiment
28 fp = open(os.path.join(PATH, experiment, 'experiment.conf'))
29 conf = json.loads(fp.read())
30
31 ipp_km = conf['usrp_tx']['ipp']
32 ipp = ipp_km * 2 /300000
33 sample_rate = conf['usrp_rx']['sample_rate']
34 axis = ['0' if x=='elevation' else '1' for x in conf['pedestal']['axis']] # AZIMUTH 1 ELEVACION 0
35 speed_axis = conf['pedestal']['speed']
36 steeps = conf['pedestal']['table']
37 time_offset = args.time_offset
38 parameters = args.parameters
39 start_date = experiment.split('@')[1].split('T')[0].replace('-', '/')
40 end_date = start_date
41 start_time = experiment.split('@')[1].split('T')[1].replace('-', ':')
42 end_time = '23:59:59'
43 N = int(1/(speed_axis[0]*ipp)) # 1 GRADO DE RESOLUCION
44 path = os.path.join(PATH, experiment, 'rawdata')
45 path_ped = os.path.join(PATH, experiment, 'position')
46 path_plots = os.path.join(PATH, experiment, 'plots_ch0')
47 path_save = os.path.join(PATH, experiment, 'param')
48 RMIX = 20
49
50 from schainpy.controller import Project
51
52 project = Project()
53 project.setup(id='1', name='Sophy', description='sophy proc')
54
55 reader = project.addReadUnit(datatype='DigitalRFReader',
56 path=path,
57 startDate=start_date,
58 endDate=end_date,
59 startTime=start_time,
60 endTime=end_time,
61 delay=30,
62 channelList='0',
63 online=args.online,
64 walk=1,
65 ippKm = ipp_km,
66 getByBlock = 1,
67 nProfileBlocks = N,
68 )
69
70 if not conf['usrp_tx']['enable_2']: # One Pulse
71 voltage = project.addProcUnit(datatype='VoltageProc', inputId=reader.getId())
72
73 if conf['usrp_tx']['code_type_1']:
74 code = [c.split() for c in conf['usrp']['code_1']]
75 op = voltage.addOperation(name='Decoder', optype='other')
76 op.addParameter(name='code', value=code)
77 op.addParameter(name='nCode', value=len(code), format='int')
78 op.addParameter(name='nBaud', value=len(code[0]), format='int')
79
80 op = voltage.addOperation(name='setH0')
81 op.addParameter(name='h0', value='-1.2')
82
83 if args.range >= 0:
84 op = voltage.addOperation(name='selectHeights')
85 op.addParameter(name='minIndex', value='0', format='int')
86 op.addParameter(name='maxIndex', value=max_index(RMIX, sample_rate, ipp), format='int')
87
88 code=[[1]]
89 opObj11 = voltage.addOperation(name='Decoder', optype='other')
90 opObj11.addParameter(name='code', value=code)
91 opObj11.addParameter(name='nCode', value='1', format='int')
92 opObj11.addParameter(name='nBaud', value='1', format='int')
93
94 op = voltage2.addOperation(name='CohInt', optype='other') #Minimo integrar 2 perfiles por ser codigo complementario
95 op.addParameter(name='n', value=2*len(code), format='int')
96
97 #op = voltage.addOperation(name='PulsePair_vRF', optype='other')
98 #op.addParameter(name='n', value=int(N), format='int')
99
100 if args.range >= 0:
101 op = voltage.addOperation(name='selectHeights')
102 op.addParameter(name='minIndex', value='0', format='int')
103 op.addParameter(name='maxIndex', value=max_index(RMIX, sample_rate, ipp), format='int')
104
105
106 op = voltage.addOperation(name='PulsePair_vRF', optype='other')
107 op.addParameter(name='n', value=125, format='int')
108
109
110 proc = project.addProcUnit(datatype='ParametersProc', inputId=voltage.getId())
111 #procUnitConfObjB.addParameter(name='runNextUnit', value=True)
112
113 opObj10 = proc.addOperation(name="WeatherRadar")
114 opObj10.addParameter(name='variableList',value='Reflectividad,VelocidadRadial,AnchoEspectral')
115
116 # {"latitude": -12.0404828587, "longitude": -75.2147483647, "altitude": 3379.2147483647}
117
118 op = proc.addOperation(name='PedestalInformation')
119 op.addParameter(name='path', value=path_ped, format='str')
120 op.addParameter(name='interval', value='0.04')
121 op.addParameter(name='time_offset', value=time_offset)
122 op.addParameter(name='az_offset', value=-26.2)
123
124 for param in parameters:
125 op = proc.addOperation(name='Block360_vRF4')
126 #op.addParameter(name='axis', value=','.join(axis))
127 op.addParameter(name='attr_data', value=PARAM[param]['name'])
128 op.addParameter(name='runNextOp', value=True)
129
130 op= proc.addOperation(name='WeatherParamsPlot')
131 if args.save: op.addParameter(name='save', value=path_plots, format='str')
132 op.addParameter(name='save_period', value=-1)
133 op.addParameter(name='show', value=args.show)
134 op.addParameter(name='channels', value='0,')
135 op.addParameter(name='zmin', value=PARAM[param]['zmin'])
136 op.addParameter(name='zmax', value=PARAM[param]['zmax'])
137 op.addParameter(name='attr_data', value=PARAM[param]['name'], format='str')
138 op.addParameter(name='labels', value=[PARAM[param]['label']])
139 op.addParameter(name='save_code', value=param)
140 op.addParameter(name='cb_label', value=PARAM[param]['cb_label'])
141 op.addParameter(name='colormap', value=PARAM[param]['colormap'])
142
143 desc = {
144 'Data': {
145 PARAM[param]['name']: PARAM[param]['label'],
146 'utctime': 'time'
147 },
148 'Metadata': {
149 'heightList': 'range',
150 'data_azi': 'azimuth',
151 'data_ele': 'elevation',
152 }
153 }
154
155 if args.save:
156 opObj10 = proc.addOperation(name='HDFWriter')
157 opObj10.addParameter(name='path',value=path_save+'-{}'.format(param), format='str')
158 opObj10.addParameter(name='Reset',value=True)
159 opObj10.addParameter(name='setType',value='weather')
160 opObj10.addParameter(name='description',value='desc')
161 opObj10.addParameter(name='blocksPerFile',value='1',format='int')
162 opObj10.addParameter(name='metadataList',value='heightList,data_azi,data_ele')
163 opObj10.addParameter(name='dataList',value='{},utctime'.format(PARAM[param]['name']))
164
165 else: #Two pulses
166
167 voltage1 = project.addProcUnit(datatype='VoltageProc', inputId=reader.getId())
168
169 print("repetions",conf['usrp_tx']['repetitions_1'])
170
171 op = voltage1.addOperation(name='ProfileSelector')
172 op.addParameter(name='profileRangeList', value='0,{}'.format(conf['usrp_tx']['repetitions_1']-1))
173
174
175 #op3 = voltage1.addOperation(name='ProfileSelector', optype='other')
176 #op3.addParameter(name='profileRangeList', value='1,123')
177
178 '''
179 if conf['usrp_tx']['code_type_1'] != 'None':
180 code = [c.split() for c in conf['usrp_tx']['code_1']]
181 op = voltage1.addOperation(name='Decoder', optype='other')
182 op.addParameter(name='code', value=code)
183 op.addParameter(name='nCode', value=len(code), format='int')
184 op.addParameter(name='nBaud', value=len(code[0]), format='int')
185 '''
186
187 code=[[1]]
188
189 opObj11 = voltage1.addOperation(name='Decoder', optype='other')
190 opObj11.addParameter(name='code', value=code)
191 opObj11.addParameter(name='nCode', value='1', format='int')
192 opObj11.addParameter(name='nBaud', value='1', format='int')
193
194 op = voltage1.addOperation(name='setH0')
195 op.addParameter(name='h0', value='-1.2')
196
197 if args.range >= 0:
198 op = voltage1.addOperation(name='selectHeights')
199 op.addParameter(name='minIndex', value='0', format='int')
200 op.addParameter(name='maxIndex', value=max_index(RMIX, sample_rate, ipp), format='int')
201
202 op = voltage1.addOperation(name='CohInt', optype='other') #Minimo integrar 2 perfiles por ser codigo complementario
203 op.addParameter(name='n', value=2, format='int')
204
205 op = voltage1.addOperation(name='PulsePair_vRF', optype='other')
206 #op.addParameter(name='n', value=int(N), format='int')
207 op.addParameter(name='n', value=61, format='int')
208 #op.addParameter(name='removeDC',value=True)
209
210 '''
211 if args.range >= 0:
212 print("corto",max_index(RMIX, sample_rate, ipp))
213 op = voltage1.addOperation(name='selectHeights')
214 op.addParameter(name='minIndex', value='0', format='int')
215 op.addParameter(name='maxIndex', value=max_index(RMIX, sample_rate, ipp), format='int')
216 '''
217 proc1 = project.addProcUnit(datatype='ParametersProc', inputId=voltage1.getId())
218 proc1.addParameter(name='runNextUnit', value=True)
219
220 opObj10 = proc1.addOperation(name="WeatherRadar")
221 opObj10.addParameter(name='variableList',value='Reflectividad,VelocidadRadial,AnchoEspectral')
222 opObj10.addParameter(name='tauW',value=0.4*1e-6)
223 opObj10.addParameter(name='Pt',value=0.2)
224
225 # {"latitude": -12.0404828587, "longitude": -75.2147483647, "altitude": 3379.2147483647}
226
227 op = proc1.addOperation(name='PedestalInformation')
228 op.addParameter(name='path', value=path_ped, format='str')
229 op.addParameter(name='interval', value='0.04')
230 op.addParameter(name='time_offset', value=time_offset)
231 op.addParameter(name='az_offset', value=-26.2)
232
233 for param in parameters:
234 op = proc1.addOperation(name='Block360_vRF4')
235 op.addParameter(name='attr_data', value=PARAM[param]['name'])
236 op.addParameter(name='runNextOp', value=True)
237
238 voltage2 = project.addProcUnit(datatype='VoltageProc', inputId=reader.getId())
239
240 op = voltage2.addOperation(name='ProfileSelector')
241 op.addParameter(name='profileRangeList', value='{},{}'.format(conf['usrp_tx']['repetitions_1'], conf['usrp_tx']['repetitions_1']+conf['usrp_tx']['repetitions_2']-1))
242
243
244 if conf['usrp_tx']['code_type_2']:
245 print(conf['usrp_tx']['code_2'])
246 codes = [ c.strip() for c in conf['usrp_tx']['code_2'].split(',')]
247 code = []
248 for c in codes:
249 code.append([int(x) for x in c])
250 print(code)
251 print(code[0])
252 op = voltage2.addOperation(name='Decoder', optype='other')
253 op.addParameter(name='code', value=code)
254 op.addParameter(name='nCode', value=len(code), format='int')
255 op.addParameter(name='nBaud', value=len(code[0]), format='int')
256 import numpy
257 pwcode = numpy.sum(numpy.array(code[0])**2)
258 print("pwcode",pwcode)
259
260 op = voltage2.addOperation(name='CohInt', optype='other') #Minimo integrar 2 perfiles por ser codigo complementario
261 op.addParameter(name='n', value=len(code), format='int')
262 ncode = len(code)
263 else:
264 ncode = 1
265
266 op = voltage2.addOperation(name='setH0')
267 op.addParameter(name='h0', value='-1.2')
268
269 if args.range >= 0:
270 if args.range==0:
271 args.range= ipp_km
272 op = voltage2.addOperation(name='selectHeights')
273 op.addParameter(name='minIndex', value=max_index(RMIX, sample_rate, ipp), format='int')
274 op.addParameter(name='maxIndex', value=max_index(args.range, sample_rate, ipp), format='int')
275
276 #op = voltage2.addOperation(name='PulsePair_vRF', optype='other')
277 #op.addParameter(name='n', value=int(N)/ncode, format='int')
278 op = voltage2.addOperation(name='PulsePair_vRF', optype='other')
279 op.addParameter(name='n', value=64, format='int')
280 #op.addParameter(name='removeDC',value=True)
281
282 '''
283
284 if args.range >= 0:
285 if args.range==0:
286 args.range= ipp_km
287 op = voltage2.addOperation(name='selectHeights')
288 print("largo",max_index(RMIX, sample_rate, ipp))
289 print("largo2",max_index(args.range, sample_rate, ipp))
290
291 op.addParameter(name='minIndex', value=max_index(RMIX, sample_rate, ipp), format='int')
292 op.addParameter(name='maxIndex', value=max_index(args.range, sample_rate, ipp), format='int')
293 '''
294
295 proc2 = project.addProcUnit(datatype='ParametersProc', inputId=voltage2.getId())
296
297 opObj10 = proc2.addOperation(name="WeatherRadar")
298 opObj10.addParameter(name='variableList',value='Reflectividad,AnchoEspectral')
299 opObj10.addParameter(name='tauW',value=6.3*1e-6)
300 opObj10.addParameter(name='Pt',value=3.2)
301
302
303 # {"latitude": -12.0404828587, "longitude": -75.2147483647, "altitude": 3379.2147483647}
304
305 op = proc2.addOperation(name='PedestalInformation')
306 op.addParameter(name='path', value=path_ped, format='str')
307 op.addParameter(name='interval', value='0.04')
308 op.addParameter(name='time_offset', value=time_offset)
309 op.addParameter(name='az_offset', value=-26.2)
310
311 for param in parameters:
312 op = proc2.addOperation(name='Block360_vRF4')
313 #op.addParameter(name='axis', value=','.join(axis))
314 op.addParameter(name='attr_data', value=PARAM[param]['name'])
315 op.addParameter(name='runNextOp', value=True)
316
317 merge = project.addProcUnit(datatype='MergeProc', inputId=[proc1.getId(), proc2.getId()])
318 merge.addParameter(name='attr_data', value=PARAM[param]['name'])
319 merge.addParameter(name='mode', value='7') #RM
320
321 op= merge.addOperation(name='WeatherParamsPlot')
322 if args.save: op.addParameter(name='save', value=path_plots, format='str')
323 op.addParameter(name='save_period', value=-1)
324 op.addParameter(name='show', value=args.show)
325 op.addParameter(name='channels', value='0,')
326 op.addParameter(name='zmin', value=PARAM[param]['zmin'])
327 op.addParameter(name='zmax', value=PARAM[param]['zmax'])
328 op.addParameter(name='attr_data', value=PARAM[param]['name'], format='str')
329 op.addParameter(name='labels', value=[PARAM[param]['label']])
330 op.addParameter(name='save_code', value=param)
331 op.addParameter(name='cb_label', value=PARAM[param]['cb_label'])
332 op.addParameter(name='colormap', value=PARAM[param]['colormap'])
333
334 desc = {
335 'Data': {
336 PARAM[param]['name']: PARAM[param]['label'],
337 'utctime': 'time'
338 },
339 'Metadata': {
340 'heightList': 'range',
341 'data_azi': 'azimuth',
342 'data_ele': 'elevation',
343 }
344 }
345
346 if args.save:
347 opObj10 = merge.addOperation(name='HDFWriter')
348 opObj10.addParameter(name='path',value=path_save, format='str')
349 opObj10.addParameter(name='Reset',value=True)
350 opObj10.addParameter(name='setType',value='weather')
351 opObj10.addParameter(name='description',value='desc')
352 opObj10.addParameter(name='blocksPerFile',value='1',format='int')
353 opObj10.addParameter(name='metadataList',value='heightList,data_azi,data_ele')
354 opObj10.addParameter(name='dataList',value='{},utctime'.format(PARAM[param]['name']))
355
356 project.start()
357
358 if __name__ == '__main__':
359
360 parser = argparse.ArgumentParser(description='Script to process SOPHy data.')
361 parser.add_argument('experiment',
362 help='Experiment name')
363 parser.add_argument('--parameters', nargs='*', default=['P'],
364 help='Variables to process: P, Z, V')
365 parser.add_argument('--time_offset', default=0,
366 help='Fix time offset')
367 parser.add_argument('--range', default=0, type=float,
368 help='Max range to plot')
369 parser.add_argument('--save', action='store_true',
370 help='Create output files')
371 parser.add_argument('--show', action='store_true',
372 help='Show matplotlib plot.')
373 parser.add_argument('--online', action='store_true',
374 help='Set online mode.')
375
376 args = parser.parse_args()
377
378 main(args)
@@ -0,0 +1,376
1 # SOPHY PROC script
2 import os, sys, json, argparse
3 import datetime
4 import time
5
6 PATH = '/DATA_RM/DATA'
7 # PATH = '/Users/jespinoza/workspace/data/'
8 #PATH = '/home/soporte/Documents/HUANCAYO'
9 PARAM = {
10 'P': {'name': 'dataPP_POWER', 'zmin': -50, 'zmax': -15, 'colormap': 'jet', 'label': 'Power', 'wrname': 'Pow','cb_label': 'dB', 'ch':0},
11 'V': {'name': 'dataPP_DOP', 'zmin': -20, 'zmax': 20, 'colormap': 'seismic', 'label': 'Velocity', 'wrname': 'Vel', 'cb_label': 'm/s', 'ch':0},
12 'RH': {'name': 'RhoHV_R', 'zmin': 0, 'zmax': 1, 'colormap': 'jet', 'label': 'Coef.Correlacion', 'wrname':'R', 'cb_label': '*', 'ch':0},
13 'FD': {'name': 'PhiD_P', 'zmin': -180,'zmax': 180,'colormap': 'RdBu_r', 'label': 'Fase Diferencial', 'wrname':'P' , 'cb_label': 'º', 'ch':0},
14 'ZD': {'name': 'Zdb_D', 'zmin': -20, 'zmax': 60, 'colormap': 'viridis','label': 'Reflect.Diferencial','wrname':'D' , 'cb_label': 'dBz','ch':0},
15 'Z': {'name': 'Zdb', 'zmin': -20, 'zmax': 70, 'colormap': 'gist_ncar','label': 'Reflectividad', 'wrname':'Z', 'cb_label': 'dBz','ch':1},
16 'W': {'name': 'Sigmav_W', 'zmin': 0, 'zmax':5, 'colormap': 'viridis','label': 'AnchoEspectral', 'wrname':'S', 'cb_label': 'hz', 'ch':1}
17 }
18
19 def max_index(r, sample_rate, ipp):
20
21 return int(sample_rate*ipp*1e6 * r / 60) + int(sample_rate*ipp*1e6 * 1.2 / 60)
22
23 def main(args):
24
25 experiment = args.experiment
26 fp = open(os.path.join(PATH, experiment, 'experiment.conf'))
27 conf = json.loads(fp.read())
28
29 ipp_km = conf['usrp_tx']['ipp']
30 ipp = ipp_km * 2 /300000
31 sample_rate = conf['usrp_rx']['sample_rate']
32 axis = ['0' if x=='elevation' else '1' for x in conf['pedestal']['axis']] # AZIMUTH 1 ELEVACION 0
33 speed_axis = conf['pedestal']['speed']
34 steeps = conf['pedestal']['table']
35 time_offset = args.time_offset
36 parameters = args.parameters
37 start_date = experiment.split('@')[1].split('T')[0].replace('-', '/')
38 end_date = start_date
39 start_time = experiment.split('@')[1].split('T')[1].replace('-', ':')
40 end_time = '23:59:59'
41 N = int(1/(speed_axis[0]*ipp)) # 1 GRADO DE RESOLUCION
42 path = os.path.join(PATH, experiment, 'rawdata')
43 path_ped = os.path.join(PATH, experiment, 'position')
44 path_plots = os.path.join(PATH, experiment, 'plots_ch1')
45 path_save = os.path.join(PATH, experiment, 'param')
46 RMIX = 2
47
48 from schainpy.controller import Project
49
50 project = Project()
51 project.setup(id='1', name='Sophy', description='sophy proc')
52
53 reader = project.addReadUnit(datatype='DigitalRFReader',
54 path=path,
55 startDate=start_date,
56 endDate=end_date,
57 startTime=start_time,
58 endTime=end_time,
59 delay=30,
60 channelList='0',
61 online=args.online,
62 walk=1,
63 ippKm = ipp_km,
64 getByBlock = 1,
65 nProfileBlocks = N,
66 )
67
68 if not conf['usrp_tx']['enable_2']: # One Pulse
69 voltage = project.addProcUnit(datatype='VoltageProc', inputId=reader.getId())
70
71 if conf['usrp_tx']['code_type_1']:
72 code = [c.split() for c in conf['usrp']['code_1']]
73 op = voltage.addOperation(name='Decoder', optype='other')
74 op.addParameter(name='code', value=code)
75 op.addParameter(name='nCode', value=len(code), format='int')
76 op.addParameter(name='nBaud', value=len(code[0]), format='int')
77
78 op = voltage.addOperation(name='setH0')
79 op.addParameter(name='h0', value='-1.2')
80
81 if args.range >= 0:
82 op = voltage.addOperation(name='selectHeights')
83 op.addParameter(name='minIndex', value='0', format='int')
84 op.addParameter(name='maxIndex', value=max_index(RMIX, sample_rate, ipp), format='int')
85
86 code=[[1]]
87 opObj11 = voltage.addOperation(name='Decoder', optype='other')
88 opObj11.addParameter(name='code', value=code)
89 opObj11.addParameter(name='nCode', value='1', format='int')
90 opObj11.addParameter(name='nBaud', value='1', format='int')
91
92 op = voltage2.addOperation(name='CohInt', optype='other') #Minimo integrar 2 perfiles por ser codigo complementario
93 op.addParameter(name='n', value=2*len(code), format='int')
94
95 #op = voltage.addOperation(name='PulsePair_vRF', optype='other')
96 #op.addParameter(name='n', value=int(N), format='int')
97
98 if args.range >= 0:
99 op = voltage.addOperation(name='selectHeights')
100 op.addParameter(name='minIndex', value='0', format='int')
101 op.addParameter(name='maxIndex', value=max_index(RMIX, sample_rate, ipp), format='int')
102
103
104 op = voltage.addOperation(name='PulsePair_vRF', optype='other')
105 op.addParameter(name='n', value=125, format='int')
106
107
108 proc = project.addProcUnit(datatype='ParametersProc', inputId=voltage.getId())
109 #procUnitConfObjB.addParameter(name='runNextUnit', value=True)
110
111 opObj10 = proc.addOperation(name="WeatherRadar")
112 opObj10.addParameter(name='variableList',value='Reflectividad,VelocidadRadial,AnchoEspectral')
113
114 # {"latitude": -12.0404828587, "longitude": -75.2147483647, "altitude": 3379.2147483647}
115
116 op = proc.addOperation(name='PedestalInformation')
117 op.addParameter(name='path', value=path_ped, format='str')
118 op.addParameter(name='interval', value='0.04')
119 op.addParameter(name='time_offset', value=time_offset)
120 op.addParameter(name='az_offset', value=-26.2)
121
122 for param in parameters:
123 op = proc.addOperation(name='Block360_vRF4')
124 #op.addParameter(name='axis', value=','.join(axis))
125 op.addParameter(name='attr_data', value=PARAM[param]['name'])
126 op.addParameter(name='runNextOp', value=True)
127
128 op= proc.addOperation(name='WeatherParamsPlot')
129 if args.save: op.addParameter(name='save', value=path_plots, format='str')
130 op.addParameter(name='save_period', value=-1)
131 op.addParameter(name='show', value=args.show)
132 op.addParameter(name='channels', value='1,')
133 op.addParameter(name='zmin', value=PARAM[param]['zmin'])
134 op.addParameter(name='zmax', value=PARAM[param]['zmax'])
135 op.addParameter(name='attr_data', value=PARAM[param]['name'], format='str')
136 op.addParameter(name='labels', value=[PARAM[param]['label']])
137 op.addParameter(name='save_code', value=param)
138 op.addParameter(name='cb_label', value=PARAM[param]['cb_label'])
139 op.addParameter(name='colormap', value=PARAM[param]['colormap'])
140
141 desc = {
142 'Data': {
143 PARAM[param]['name']: PARAM[param]['label'],
144 'utctime': 'time'
145 },
146 'Metadata': {
147 'heightList': 'range',
148 'data_azi': 'azimuth',
149 'data_ele': 'elevation',
150 }
151 }
152
153 if args.save:
154 opObj10 = proc.addOperation(name='HDFWriter')
155 opObj10.addParameter(name='path',value=path_save+'-{}'.format(param), format='str')
156 opObj10.addParameter(name='Reset',value=True)
157 opObj10.addParameter(name='setType',value='weather')
158 opObj10.addParameter(name='description',value='desc')
159 opObj10.addParameter(name='blocksPerFile',value='1',format='int')
160 opObj10.addParameter(name='metadataList',value='heightList,data_azi,data_ele')
161 opObj10.addParameter(name='dataList',value='{},utctime'.format(PARAM[param]['name']))
162
163 else: #Two pulses
164
165 voltage1 = project.addProcUnit(datatype='VoltageProc', inputId=reader.getId())
166
167 print("repetions",conf['usrp_tx']['repetitions_1'])
168
169 op = voltage1.addOperation(name='ProfileSelector')
170 op.addParameter(name='profileRangeList', value='0,{}'.format(conf['usrp_tx']['repetitions_1']-1))
171
172
173 #op3 = voltage1.addOperation(name='ProfileSelector', optype='other')
174 #op3.addParameter(name='profileRangeList', value='1,123')
175
176 '''
177 if conf['usrp_tx']['code_type_1'] != 'None':
178 code = [c.split() for c in conf['usrp_tx']['code_1']]
179 op = voltage1.addOperation(name='Decoder', optype='other')
180 op.addParameter(name='code', value=code)
181 op.addParameter(name='nCode', value=len(code), format='int')
182 op.addParameter(name='nBaud', value=len(code[0]), format='int')
183 '''
184
185 code=[[1]]
186
187 opObj11 = voltage1.addOperation(name='Decoder', optype='other')
188 opObj11.addParameter(name='code', value=code)
189 opObj11.addParameter(name='nCode', value='1', format='int')
190 opObj11.addParameter(name='nBaud', value='1', format='int')
191
192 op = voltage1.addOperation(name='setH0')
193 op.addParameter(name='h0', value='-1.2')
194
195 if args.range >= 0:
196 op = voltage1.addOperation(name='selectHeights')
197 op.addParameter(name='minIndex', value='0', format='int')
198 op.addParameter(name='maxIndex', value=max_index(RMIX, sample_rate, ipp), format='int')
199
200 op = voltage1.addOperation(name='CohInt', optype='other') #Minimo integrar 2 perfiles por ser codigo complementario
201 op.addParameter(name='n', value=2, format='int')
202
203 op = voltage1.addOperation(name='PulsePair_vRF', optype='other')
204 #op.addParameter(name='n', value=int(N), format='int')
205 op.addParameter(name='n', value=61, format='int')
206 #op.addParameter(name='removeDC',value=True)
207
208 '''
209 if args.range >= 0:
210 print("corto",max_index(RMIX, sample_rate, ipp))
211 op = voltage1.addOperation(name='selectHeights')
212 op.addParameter(name='minIndex', value='0', format='int')
213 op.addParameter(name='maxIndex', value=max_index(RMIX, sample_rate, ipp), format='int')
214 '''
215 proc1 = project.addProcUnit(datatype='ParametersProc', inputId=voltage1.getId())
216 proc1.addParameter(name='runNextUnit', value=True)
217
218 opObj10 = proc1.addOperation(name="WeatherRadar")
219 opObj10.addParameter(name='variableList',value='Reflectividad,VelocidadRadial,AnchoEspectral')
220 opObj10.addParameter(name='tauW',value=0.4*1e-6)
221 opObj10.addParameter(name='Pt',value=0.2)
222
223 # {"latitude": -12.0404828587, "longitude": -75.2147483647, "altitude": 3379.2147483647}
224
225 op = proc1.addOperation(name='PedestalInformation')
226 op.addParameter(name='path', value=path_ped, format='str')
227 op.addParameter(name='interval', value='0.04')
228 op.addParameter(name='time_offset', value=time_offset)
229 op.addParameter(name='az_offset', value=-26.2)
230
231 for param in parameters:
232 op = proc1.addOperation(name='Block360_vRF4')
233 op.addParameter(name='attr_data', value=PARAM[param]['name'])
234 op.addParameter(name='runNextOp', value=True)
235
236 voltage2 = project.addProcUnit(datatype='VoltageProc', inputId=reader.getId())
237
238 op = voltage2.addOperation(name='ProfileSelector')
239 op.addParameter(name='profileRangeList', value='{},{}'.format(conf['usrp_tx']['repetitions_1'], conf['usrp_tx']['repetitions_1']+conf['usrp_tx']['repetitions_2']-1))
240
241
242 if conf['usrp_tx']['code_type_2']:
243 print(conf['usrp_tx']['code_2'])
244 codes = [ c.strip() for c in conf['usrp_tx']['code_2'].split(',')]
245 code = []
246 for c in codes:
247 code.append([int(x) for x in c])
248 print(code)
249 print(code[0])
250 op = voltage2.addOperation(name='Decoder', optype='other')
251 op.addParameter(name='code', value=code)
252 op.addParameter(name='nCode', value=len(code), format='int')
253 op.addParameter(name='nBaud', value=len(code[0]), format='int')
254 import numpy
255 pwcode = numpy.sum(numpy.array(code[0])**2)
256 print("pwcode",pwcode)
257
258 op = voltage2.addOperation(name='CohInt', optype='other') #Minimo integrar 2 perfiles por ser codigo complementario
259 op.addParameter(name='n', value=len(code), format='int')
260 ncode = len(code)
261 else:
262 ncode = 1
263
264 op = voltage2.addOperation(name='setH0')
265 op.addParameter(name='h0', value='-1.2')
266
267 if args.range >= 0:
268 if args.range==0:
269 args.range= ipp_km
270 op = voltage2.addOperation(name='selectHeights')
271 op.addParameter(name='minIndex', value=max_index(RMIX, sample_rate, ipp), format='int')
272 op.addParameter(name='maxIndex', value=max_index(args.range, sample_rate, ipp), format='int')
273
274 #op = voltage2.addOperation(name='PulsePair_vRF', optype='other')
275 #op.addParameter(name='n', value=int(N)/ncode, format='int')
276 op = voltage2.addOperation(name='PulsePair_vRF', optype='other')
277 op.addParameter(name='n', value=64, format='int')
278 #op.addParameter(name='removeDC',value=True)
279
280 '''
281
282 if args.range >= 0:
283 if args.range==0:
284 args.range= ipp_km
285 op = voltage2.addOperation(name='selectHeights')
286 print("largo",max_index(RMIX, sample_rate, ipp))
287 print("largo2",max_index(args.range, sample_rate, ipp))
288
289 op.addParameter(name='minIndex', value=max_index(RMIX, sample_rate, ipp), format='int')
290 op.addParameter(name='maxIndex', value=max_index(args.range, sample_rate, ipp), format='int')
291 '''
292
293 proc2 = project.addProcUnit(datatype='ParametersProc', inputId=voltage2.getId())
294
295 opObj10 = proc2.addOperation(name="WeatherRadar")
296 opObj10.addParameter(name='variableList',value='Reflectividad,AnchoEspectral')
297 opObj10.addParameter(name='tauW',value=3.2*1e-6)
298 opObj10.addParameter(name='Pt',value=1.6)
299
300
301 # {"latitude": -12.0404828587, "longitude": -75.2147483647, "altitude": 3379.2147483647}
302
303 op = proc2.addOperation(name='PedestalInformation')
304 op.addParameter(name='path', value=path_ped, format='str')
305 op.addParameter(name='interval', value='0.04')
306 op.addParameter(name='time_offset', value=time_offset)
307 op.addParameter(name='az_offset', value=-26.2)
308
309 for param in parameters:
310 op = proc2.addOperation(name='Block360_vRF4')
311 #op.addParameter(name='axis', value=','.join(axis))
312 op.addParameter(name='attr_data', value=PARAM[param]['name'])
313 op.addParameter(name='runNextOp', value=True)
314
315 merge = project.addProcUnit(datatype='MergeProc', inputId=[proc1.getId(), proc2.getId()])
316 merge.addParameter(name='attr_data', value=PARAM[param]['name'])
317 merge.addParameter(name='mode', value='7') #RM
318
319 op= merge.addOperation(name='WeatherParamsPlot')
320 if args.save: op.addParameter(name='save', value=path_plots, format='str')
321 op.addParameter(name='save_period', value=-1)
322 op.addParameter(name='show', value=args.show)
323 op.addParameter(name='channels', value='1,')
324 op.addParameter(name='zmin', value=PARAM[param]['zmin'])
325 op.addParameter(name='zmax', value=PARAM[param]['zmax'])
326 op.addParameter(name='attr_data', value=PARAM[param]['name'], format='str')
327 op.addParameter(name='labels', value=[PARAM[param]['label']])
328 op.addParameter(name='save_code', value=param)
329 op.addParameter(name='cb_label', value=PARAM[param]['cb_label'])
330 op.addParameter(name='colormap', value=PARAM[param]['colormap'])
331
332 desc = {
333 'Data': {
334 PARAM[param]['name']: PARAM[param]['label'],
335 'utctime': 'time'
336 },
337 'Metadata': {
338 'heightList': 'range',
339 'data_azi': 'azimuth',
340 'data_ele': 'elevation',
341 }
342 }
343
344 if args.save:
345 opObj10 = merge.addOperation(name='HDFWriter')
346 opObj10.addParameter(name='path',value=path_save, format='str')
347 opObj10.addParameter(name='Reset',value=True)
348 opObj10.addParameter(name='setType',value='weather')
349 opObj10.addParameter(name='description',value='desc')
350 opObj10.addParameter(name='blocksPerFile',value='1',format='int')
351 opObj10.addParameter(name='metadataList',value='heightList,data_azi,data_ele')
352 opObj10.addParameter(name='dataList',value='{},utctime'.format(PARAM[param]['name']))
353
354 project.start()
355
356 if __name__ == '__main__':
357
358 parser = argparse.ArgumentParser(description='Script to process SOPHy data.')
359 parser.add_argument('experiment',
360 help='Experiment name')
361 parser.add_argument('--parameters', nargs='*', default=['P'],
362 help='Variables to process: P, Z, V')
363 parser.add_argument('--time_offset', default=0,
364 help='Fix time offset')
365 parser.add_argument('--range', default=0, type=float,
366 help='Max range to plot')
367 parser.add_argument('--save', action='store_true',
368 help='Create output files')
369 parser.add_argument('--show', action='store_true',
370 help='Show matplotlib plot.')
371 parser.add_argument('--online', action='store_true',
372 help='Set online mode.')
373
374 args = parser.parse_args()
375
376 main(args)
@@ -1,716 +1,739
1 1 # Copyright (c) 2012-2020 Jicamarca Radio Observatory
2 2 # All rights reserved.
3 3 #
4 4 # Distributed under the terms of the BSD 3-clause license.
5 5 """Base class to create plot operations
6 6
7 7 """
8 8
9 9 import os
10 10 import sys
11 11 import zmq
12 12 import time
13 13 import numpy
14 14 import datetime
15 15 from collections import deque
16 16 from functools import wraps
17 17 from threading import Thread
18 import matplotlib
18 import matplotlib,re
19 19
20 20 if 'BACKEND' in os.environ:
21 21 matplotlib.use(os.environ['BACKEND'])
22 22 elif 'linux' in sys.platform:
23 matplotlib.use("TkAgg")
23 matplotlib.use("Agg")#TkAgg
24 24 elif 'darwin' in sys.platform:
25 25 matplotlib.use('MacOSX')
26 26 else:
27 27 from schainpy.utils import log
28 28 log.warning('Using default Backend="Agg"', 'INFO')
29 29 matplotlib.use('Agg')
30 30
31 31 import matplotlib.pyplot as plt
32 32 from matplotlib.patches import Polygon
33 33 from mpl_toolkits.axes_grid1 import make_axes_locatable
34 34 from matplotlib.ticker import FuncFormatter, LinearLocator, MultipleLocator
35 35
36 36 from schainpy.model.data.jrodata import PlotterData
37 37 from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation, MPDecorator
38 38 from schainpy.utils import log
39 39
40 40 jet_values = matplotlib.pyplot.get_cmap('jet', 100)(numpy.arange(100))[10:90]
41 41 blu_values = matplotlib.pyplot.get_cmap(
42 42 'seismic_r', 20)(numpy.arange(20))[10:15]
43 43 ncmap = matplotlib.colors.LinearSegmentedColormap.from_list(
44 44 'jro', numpy.vstack((blu_values, jet_values)))
45 45 matplotlib.pyplot.register_cmap(cmap=ncmap)
46 46
47 rwg=matplotlib.colors.LinearSegmentedColormap.from_list('rwg',["r", "w", "g"], N=256)
48 matplotlib.pyplot.register_cmap(cmap=rwg)
49
47 50 CMAPS = [plt.get_cmap(s) for s in ('jro', 'jet', 'viridis',
48 'plasma', 'inferno', 'Greys', 'seismic', 'bwr', 'coolwarm')]
51 'plasma', 'inferno', 'Greys', 'seismic', 'bwr', 'coolwarm','rwg')]
49 52
50 53 EARTH_RADIUS = 6.3710e3
51 54
52 55 def ll2xy(lat1, lon1, lat2, lon2):
53 56
54 57 p = 0.017453292519943295
55 58 a = 0.5 - numpy.cos((lat2 - lat1) * p)/2 + numpy.cos(lat1 * p) * \
56 59 numpy.cos(lat2 * p) * (1 - numpy.cos((lon2 - lon1) * p)) / 2
57 60 r = 12742 * numpy.arcsin(numpy.sqrt(a))
58 61 theta = numpy.arctan2(numpy.sin((lon2-lon1)*p)*numpy.cos(lat2*p), numpy.cos(lat1*p)
59 62 * numpy.sin(lat2*p)-numpy.sin(lat1*p)*numpy.cos(lat2*p)*numpy.cos((lon2-lon1)*p))
60 63 theta = -theta + numpy.pi/2
61 64 return r*numpy.cos(theta), r*numpy.sin(theta)
62 65
63 66
64 67 def km2deg(km):
65 68 '''
66 69 Convert distance in km to degrees
67 70 '''
68 71
69 72 return numpy.rad2deg(km/EARTH_RADIUS)
70 73
71 74
72 75 def figpause(interval):
73 76 backend = plt.rcParams['backend']
74 77 if backend in matplotlib.rcsetup.interactive_bk:
75 78 figManager = matplotlib._pylab_helpers.Gcf.get_active()
76 79 if figManager is not None:
77 80 canvas = figManager.canvas
78 81 if canvas.figure.stale:
79 82 canvas.draw()
80 83 try:
81 84 canvas.start_event_loop(interval)
82 85 except:
83 86 pass
84 87 return
85 88
86 89 def popup(message):
87 90 '''
88 91 '''
89 92
90 93 fig = plt.figure(figsize=(12, 8), facecolor='r')
91 94 text = '\n'.join([s.strip() for s in message.split(':')])
92 95 fig.text(0.01, 0.5, text, ha='left', va='center',
93 96 size='20', weight='heavy', color='w')
94 97 fig.show()
95 98 figpause(1000)
96 99
97 100
98 101 class Throttle(object):
99 102 '''
100 103 Decorator that prevents a function from being called more than once every
101 104 time period.
102 105 To create a function that cannot be called more than once a minute, but
103 106 will sleep until it can be called:
104 107 @Throttle(minutes=1)
105 108 def foo():
106 109 pass
107 110
108 111 for i in range(10):
109 112 foo()
110 113 print "This function has run %s times." % i
111 114 '''
112 115
113 116 def __init__(self, seconds=0, minutes=0, hours=0):
114 117 self.throttle_period = datetime.timedelta(
115 118 seconds=seconds, minutes=minutes, hours=hours
116 119 )
117 120
118 121 self.time_of_last_call = datetime.datetime.min
119 122
120 123 def __call__(self, fn):
121 124 @wraps(fn)
122 125 def wrapper(*args, **kwargs):
123 126 coerce = kwargs.pop('coerce', None)
124 127 if coerce:
125 128 self.time_of_last_call = datetime.datetime.now()
126 129 return fn(*args, **kwargs)
127 130 else:
128 131 now = datetime.datetime.now()
129 132 time_since_last_call = now - self.time_of_last_call
130 133 time_left = self.throttle_period - time_since_last_call
131 134
132 135 if time_left > datetime.timedelta(seconds=0):
133 136 return
134 137
135 138 self.time_of_last_call = datetime.datetime.now()
136 139 return fn(*args, **kwargs)
137 140
138 141 return wrapper
139 142
140 143 def apply_throttle(value):
141 144
142 145 @Throttle(seconds=value)
143 146 def fnThrottled(fn):
144 147 fn()
145 148
146 149 return fnThrottled
147 150
148 151
149 152 @MPDecorator
150 153 class Plot(Operation):
151 154 """Base class for Schain plotting operations
152 155
153 156 This class should never be use directtly you must subclass a new operation,
154 157 children classes must be defined as follow:
155 158
156 159 ExamplePlot(Plot):
157 160
158 161 CODE = 'code'
159 162 colormap = 'jet'
160 163 plot_type = 'pcolor' # options are ('pcolor', 'pcolorbuffer', 'scatter', 'scatterbuffer')
161 164
162 165 def setup(self):
163 166 pass
164 167
165 168 def plot(self):
166 169 pass
167 170
168 171 """
169 172
170 173 CODE = 'Figure'
171 174 colormap = 'jet'
172 175 bgcolor = 'white'
173 176 buffering = True
174 177 __missing = 1E30
175 178
176 179 __attrs__ = ['show', 'save', 'ymin', 'ymax', 'zmin', 'zmax', 'title',
177 180 'showprofile']
178 181
179 182 def __init__(self):
180 183
181 184 Operation.__init__(self)
182 185 self.isConfig = False
183 186 self.isPlotConfig = False
184 187 self.save_time = 0
185 188 self.sender_time = 0
186 189 self.data = None
187 190 self.firsttime = True
188 191 self.sender_queue = deque(maxlen=10)
189 192 self.plots_adjust = {'left': 0.125, 'right': 0.9, 'bottom': 0.15, 'top': 0.9, 'wspace': 0.2, 'hspace': 0.2}
190 193
191 194 def __fmtTime(self, x, pos):
192 195 '''
193 196 '''
194 197
195 198 return '{}'.format(self.getDateTime(x).strftime('%H:%M'))
196 199
197 200 def __setup(self, **kwargs):
198 201 '''
199 202 Initialize variables
200 203 '''
201 204
202 205 self.figures = []
203 206 self.axes = []
204 207 self.cb_axes = []
205 208 self.localtime = kwargs.pop('localtime', True)
206 209 self.show = kwargs.get('show', True)
207 210 self.save = kwargs.get('save', False)
208 211 self.save_period = kwargs.get('save_period', 0)
209 212 self.colormap = kwargs.get('colormap', self.colormap)
210 213 self.colormap_coh = kwargs.get('colormap_coh', 'jet')
211 214 self.colormap_phase = kwargs.get('colormap_phase', 'RdBu_r')
212 215 self.colormaps = kwargs.get('colormaps', None)
213 216 self.bgcolor = kwargs.get('bgcolor', self.bgcolor)
214 217 self.showprofile = kwargs.get('showprofile', False)
215 218 self.title = kwargs.get('wintitle', self.CODE.upper())
216 219 self.cb_label = kwargs.get('cb_label', None)
217 220 self.cb_labels = kwargs.get('cb_labels', None)
218 221 self.labels = kwargs.get('labels', None)
219 222 self.xaxis = kwargs.get('xaxis', 'frequency')
220 223 self.zmin = kwargs.get('zmin', None)
221 224 self.zmax = kwargs.get('zmax', None)
222 225 self.zlimits = kwargs.get('zlimits', None)
223 226 self.xmin = kwargs.get('xmin', None)
224 227 self.xmax = kwargs.get('xmax', None)
225 228 self.xrange = kwargs.get('xrange', 12)
226 229 self.xscale = kwargs.get('xscale', None)
227 230 self.ymin = kwargs.get('ymin', None)
228 231 self.ymax = kwargs.get('ymax', None)
229 232 self.yscale = kwargs.get('yscale', None)
230 233 self.xlabel = kwargs.get('xlabel', None)
231 234 self.attr_time = kwargs.get('attr_time', 'utctime')
232 235 self.attr_data = kwargs.get('attr_data', 'data_param')
233 236 self.decimation = kwargs.get('decimation', None)
234 237 self.oneFigure = kwargs.get('oneFigure', True)
235 238 self.width = kwargs.get('width', None)
236 239 self.height = kwargs.get('height', None)
237 240 self.colorbar = kwargs.get('colorbar', True)
238 241 self.factors = kwargs.get('factors', [1, 1, 1, 1, 1, 1, 1, 1])
239 242 self.channels = kwargs.get('channels', None)
240 243 self.titles = kwargs.get('titles', [])
241 244 self.polar = False
242 245 self.type = kwargs.get('type', 'iq')
243 246 self.grid = kwargs.get('grid', False)
244 247 self.pause = kwargs.get('pause', False)
245 248 self.save_code = kwargs.get('save_code', self.CODE)
246 249 self.throttle = kwargs.get('throttle', 0)
247 250 self.exp_code = kwargs.get('exp_code', None)
248 251 self.server = kwargs.get('server', False)
249 252 self.sender_period = kwargs.get('sender_period', 60)
250 253 self.tag = kwargs.get('tag', '')
251 254 self.height_index = kwargs.get('height_index', None)
252 255 self.__throttle_plot = apply_throttle(self.throttle)
253 256 code = self.attr_data if self.attr_data else self.CODE
254 257 self.data = PlotterData(self.CODE, self.exp_code, self.localtime)
255 258 self.ang_min = kwargs.get('ang_min', None)
256 259 self.ang_max = kwargs.get('ang_max', None)
257 260 self.mode = kwargs.get('mode', None)
258 261
259 262
260 263
261 264 if self.server:
262 265 if not self.server.startswith('tcp://'):
263 266 self.server = 'tcp://{}'.format(self.server)
264 267 log.success(
265 268 'Sending to server: {}'.format(self.server),
266 269 self.name
267 270 )
268 271
269 272 if isinstance(self.attr_data, str):
270 273 self.attr_data = [self.attr_data]
271 274
272 275 def __setup_plot(self):
273 276 '''
274 277 Common setup for all figures, here figures and axes are created
275 278 '''
276 279
277 280 self.setup()
278 281
279 282 self.time_label = 'LT' if self.localtime else 'UTC'
280 283
281 284 if self.width is None:
282 285 self.width = 8
283 286
284 287 self.figures = []
285 288 self.axes = []
286 289 self.cb_axes = []
287 290 self.pf_axes = []
288 291 self.cmaps = []
289 292
290 293 size = '15%' if self.ncols == 1 else '30%'
291 294 pad = '4%' if self.ncols == 1 else '8%'
292 295
293 296 if self.oneFigure:
294 297 if self.height is None:
295 298 self.height = 1.4 * self.nrows + 1
296 299 fig = plt.figure(figsize=(self.width, self.height),
297 300 edgecolor='k',
298 301 facecolor='w')
299 302 self.figures.append(fig)
300 303 for n in range(self.nplots):
301 304 ax = fig.add_subplot(self.nrows, self.ncols,
302 305 n + 1, polar=self.polar)
303 306 ax.tick_params(labelsize=8)
304 307 ax.firsttime = True
305 308 ax.index = 0
306 309 ax.press = None
307 310 self.axes.append(ax)
308 311 if self.showprofile:
309 312 cax = self.__add_axes(ax, size=size, pad=pad)
310 313 cax.tick_params(labelsize=8)
311 314 self.pf_axes.append(cax)
312 315 else:
313 316 if self.height is None:
314 317 self.height = 3
315 318 for n in range(self.nplots):
316 319 fig = plt.figure(figsize=(self.width, self.height),
317 320 edgecolor='k',
318 321 facecolor='w')
319 322 ax = fig.add_subplot(1, 1, 1, polar=self.polar)
320 323 ax.tick_params(labelsize=8)
321 324 ax.firsttime = True
322 325 ax.index = 0
323 326 ax.press = None
324 327 self.figures.append(fig)
325 328 self.axes.append(ax)
326 329 if self.showprofile:
327 330 cax = self.__add_axes(ax, size=size, pad=pad)
328 331 cax.tick_params(labelsize=8)
329 332 self.pf_axes.append(cax)
330 333
331 334 for n in range(self.nrows):
332 335 if self.colormaps is not None:
333 336 cmap = plt.get_cmap(self.colormaps[n])
334 337 else:
335 338 cmap = plt.get_cmap(self.colormap)
336 339 cmap.set_bad(self.bgcolor, 1.)
337 340 self.cmaps.append(cmap)
338 341
339 342 def __add_axes(self, ax, size='30%', pad='8%'):
340 343 '''
341 344 Add new axes to the given figure
342 345 '''
343 346 divider = make_axes_locatable(ax)
344 347 nax = divider.new_horizontal(size=size, pad=pad)
345 348 ax.figure.add_axes(nax)
346 349 return nax
347 350
348 351 def fill_gaps(self, x_buffer, y_buffer, z_buffer):
349 352 '''
350 353 Create a masked array for missing data
351 354 '''
352 355 if x_buffer.shape[0] < 2:
353 356 return x_buffer, y_buffer, z_buffer
354 357
355 358 deltas = x_buffer[1:] - x_buffer[0:-1]
356 359 x_median = numpy.median(deltas)
357 360
358 361 index = numpy.where(deltas > 5 * x_median)
359 362
360 363 if len(index[0]) != 0:
361 364 z_buffer[::, index[0], ::] = self.__missing
362 365 z_buffer = numpy.ma.masked_inside(z_buffer,
363 366 0.99 * self.__missing,
364 367 1.01 * self.__missing)
365 368
366 369 return x_buffer, y_buffer, z_buffer
367 370
368 371 def decimate(self):
369 372
370 373 # dx = int(len(self.x)/self.__MAXNUMX) + 1
371 374 dy = int(len(self.y) / self.decimation) + 1
372 375
373 376 # x = self.x[::dx]
374 377 x = self.x
375 378 y = self.y[::dy]
376 379 z = self.z[::, ::, ::dy]
377 380
378 381 return x, y, z
379 382
380 383 def format(self):
381 384 '''
382 385 Set min and max values, labels, ticks and titles
383 386 '''
384 387
385 388 for n, ax in enumerate(self.axes):
386 389 if ax.firsttime:
387 390 if self.xaxis != 'time':
388 391 xmin = self.xmin
389 392 xmax = self.xmax
390 393 else:
391 394 xmin = self.tmin
392 395 xmax = self.tmin + self.xrange*60*60
393 396 ax.xaxis.set_major_formatter(FuncFormatter(self.__fmtTime))
394 397 ax.xaxis.set_major_locator(LinearLocator(9))
395 398 ymin = self.ymin if self.ymin is not None else numpy.nanmin(self.y[numpy.isfinite(self.y)])
396 399 ymax = self.ymax if self.ymax is not None else numpy.nanmax(self.y[numpy.isfinite(self.y)])
397 400 ax.set_facecolor(self.bgcolor)
398 401 if self.xscale:
399 402 ax.xaxis.set_major_formatter(FuncFormatter(
400 403 lambda x, pos: '{0:g}'.format(x*self.xscale)))
401 404 if self.yscale:
402 405 ax.yaxis.set_major_formatter(FuncFormatter(
403 406 lambda x, pos: '{0:g}'.format(x*self.yscale)))
404 407 if self.xlabel is not None:
405 408 ax.set_xlabel(self.xlabel)
406 409 if self.ylabel is not None:
407 410 ax.set_ylabel(self.ylabel)
408 411 if self.showprofile:
409 412 self.pf_axes[n].set_ylim(ymin, ymax)
410 413 self.pf_axes[n].set_xlim(self.zmin, self.zmax)
411 414 self.pf_axes[n].set_xlabel('dB')
412 415 self.pf_axes[n].grid(b=True, axis='x')
413 416 [tick.set_visible(False)
414 417 for tick in self.pf_axes[n].get_yticklabels()]
415 418 if self.colorbar:
416 419 ax.cbar = plt.colorbar(
417 ax.plt, ax=ax, fraction=0.05, pad=0.02, aspect=10)
420 ax.plt, ax=ax, fraction=0.05, pad=0.06, aspect=10)
418 421 ax.cbar.ax.tick_params(labelsize=8)
419 422 ax.cbar.ax.press = None
420 423 if self.cb_label:
421 424 ax.cbar.set_label(self.cb_label, size=8)
422 425 elif self.cb_labels:
423 426 ax.cbar.set_label(self.cb_labels[n], size=8)
424 427 else:
425 428 ax.cbar = None
426 429 ax.set_xlim(xmin, xmax)
427 430 ax.set_ylim(ymin, ymax)
428 431 ax.firsttime = False
429 432 if self.grid:
430 433 ax.grid(True)
431 434 if not self.polar:
432 435 ax.set_title('{} {} {}'.format(
433 436 self.titles[n],
434 437 self.getDateTime(self.data.max_time).strftime(
435 438 '%Y-%m-%d %H:%M:%S'),
436 439 self.time_label),
437 440 size=8)
438 441 else:
439 442 #ax.set_title('{}'.format(self.titles[n]), size=8)
440 443 ax.set_title('{} {} {}'.format(
441 444 self.titles[n],
442 445 self.getDateTime(self.data.max_time).strftime(
443 446 '%Y-%m-%d %H:%M:%S'),
444 447 self.time_label),
445 448 size=8)
446 449 ax.set_ylim(0, self.ymax)
447 450 #ax.set_yticks(numpy.arange(0, self.ymax, 20))
448 ax.yaxis.labelpad = 20
451 ax.yaxis.labelpad = 28
449 452
450 453 if self.firsttime:
451 454 for n, fig in enumerate(self.figures):
452 455 fig.subplots_adjust(**self.plots_adjust)
453 456 self.firsttime = False
454 457
455 458 def clear_figures(self):
456 459 '''
457 460 Reset axes for redraw plots
458 461 '''
459 462
460 463 for ax in self.axes+self.pf_axes+self.cb_axes:
461 464 ax.clear()
462 465 ax.firsttime = True
463 466 if hasattr(ax, 'cbar') and ax.cbar:
464 467 ax.cbar.remove()
465 468
466 469 def __plot(self):
467 470 '''
468 471 Main function to plot, format and save figures
469 472 '''
470 473
471 474 self.plot()
472 475 self.format()
473 476
474 477 for n, fig in enumerate(self.figures):
475 478 if self.nrows == 0 or self.nplots == 0:
476 479 log.warning('No data', self.name)
477 480 fig.text(0.5, 0.5, 'No Data', fontsize='large', ha='center')
478 481 fig.canvas.manager.set_window_title(self.CODE)
479 482 continue
480 483
481 484 fig.canvas.manager.set_window_title('{} - {}'.format(self.title,
482 485 self.getDateTime(self.data.max_time).strftime('%Y/%m/%d')))
483 486 fig.canvas.draw()
484 487 if self.show:
485 488 fig.show()
486 489 figpause(0.01)
487 490
488 491 if self.save:
492 if self.CODE=="PPI" or self.CODE=="RHI":
493 self.save_figure(n,stitle =self.titles)
494 else:
489 495 self.save_figure(n)
490 496
491 497 if self.server:
492 498 self.send_to_server()
493 499
494 500 def __update(self, dataOut, timestamp):
495 501 '''
496 502 '''
497 503
498 504 metadata = {
499 505 'yrange': dataOut.heightList,
500 506 'interval': dataOut.timeInterval,
501 507 'channels': dataOut.channelList
502 508 }
503 509
504 510 data, meta = self.update(dataOut)
505 511 metadata.update(meta)
506 512 self.data.update(data, timestamp, metadata)
507 513
508 def save_figure(self, n):
514 def save_figure(self, n,stitle=None):
509 515 '''
510 516 '''
517 if stitle is not None:
518 s_string = re.sub(r"[^A-Z0-9.]","",str(stitle))
519 new_string=s_string[:3]+"_"+"_"+s_string[4:6]+"_"+s_string[6:]
520
511 521 if self.oneFigure:
512 522 if (self.data.max_time - self.save_time) <= self.save_period:
513 523 return
514 524
515 525 self.save_time = self.data.max_time
516 526
517 527 fig = self.figures[n]
518 528
519 529 if self.throttle == 0:
520 530 if self.oneFigure:
531 if stitle is not None:
532 figname = os.path.join(
533 self.save,
534 self.save_code,
535 '{}_{}_{}.png'.format(
536 self.save_code,
537 self.getDateTime(self.data.max_time).strftime(
538 '%Y%m%d_%H%M%S',
539 ),
540 new_string,
541 )
542 )
543 else:
521 544 figname = os.path.join(
522 545 self.save,
523 546 self.save_code,
524 547 '{}_{}.png'.format(
525 548 self.save_code,
526 549 self.getDateTime(self.data.max_time).strftime(
527 550 '%Y%m%d_%H%M%S'
528 551 ),
529 552 )
530 553 )
531 554 else:
532 555 figname = os.path.join(
533 556 self.save,
534 557 self.save_code,
535 558 '{}_ch{}_{}.png'.format(
536 559 self.save_code,n,
537 560 self.getDateTime(self.data.max_time).strftime(
538 561 '%Y%m%d_%H%M%S'
539 562 ),
540 563 )
541 564 )
542 565 log.log('Saving figure: {}'.format(figname), self.name)
543 566 if not os.path.isdir(os.path.dirname(figname)):
544 567 os.makedirs(os.path.dirname(figname))
545 568 fig.savefig(figname)
546 569
547 570 figname = os.path.join(
548 571 self.save,
549 572 '{}_{}.png'.format(
550 573 self.save_code,
551 574 self.getDateTime(self.data.min_time).strftime(
552 575 '%Y%m%d'
553 576 ),
554 577 )
555 578 )
556 579
557 580 log.log('Saving figure: {}'.format(figname), self.name)
558 581 if not os.path.isdir(os.path.dirname(figname)):
559 582 os.makedirs(os.path.dirname(figname))
560 583 fig.savefig(figname)
561 584
562 585 def send_to_server(self):
563 586 '''
564 587 '''
565 588
566 589 if self.exp_code == None:
567 590 log.warning('Missing `exp_code` skipping sending to server...')
568 591
569 592 last_time = self.data.max_time
570 593 interval = last_time - self.sender_time
571 594 if interval < self.sender_period:
572 595 return
573 596
574 597 self.sender_time = last_time
575 598
576 599 attrs = ['titles', 'zmin', 'zmax', 'tag', 'ymin', 'ymax']
577 600 for attr in attrs:
578 601 value = getattr(self, attr)
579 602 if value:
580 603 if isinstance(value, (numpy.float32, numpy.float64)):
581 604 value = round(float(value), 2)
582 605 self.data.meta[attr] = value
583 606 if self.colormap == 'jet':
584 607 self.data.meta['colormap'] = 'Jet'
585 608 elif 'RdBu' in self.colormap:
586 609 self.data.meta['colormap'] = 'RdBu'
587 610 else:
588 611 self.data.meta['colormap'] = 'Viridis'
589 612 self.data.meta['interval'] = int(interval)
590 613
591 614 self.sender_queue.append(last_time)
592 615
593 616 while True:
594 617 try:
595 618 tm = self.sender_queue.popleft()
596 619 except IndexError:
597 620 break
598 621 msg = self.data.jsonify(tm, self.save_code, self.plot_type)
599 622 self.socket.send_string(msg)
600 623 socks = dict(self.poll.poll(2000))
601 624 if socks.get(self.socket) == zmq.POLLIN:
602 625 reply = self.socket.recv_string()
603 626 if reply == 'ok':
604 627 log.log("Response from server ok", self.name)
605 628 time.sleep(0.1)
606 629 continue
607 630 else:
608 631 log.warning(
609 632 "Malformed reply from server: {}".format(reply), self.name)
610 633 else:
611 634 log.warning(
612 635 "No response from server, retrying...", self.name)
613 636 self.sender_queue.appendleft(tm)
614 637 self.socket.setsockopt(zmq.LINGER, 0)
615 638 self.socket.close()
616 639 self.poll.unregister(self.socket)
617 640 self.socket = self.context.socket(zmq.REQ)
618 641 self.socket.connect(self.server)
619 642 self.poll.register(self.socket, zmq.POLLIN)
620 643 break
621 644
622 645 def setup(self):
623 646 '''
624 647 This method should be implemented in the child class, the following
625 648 attributes should be set:
626 649
627 650 self.nrows: number of rows
628 651 self.ncols: number of cols
629 652 self.nplots: number of plots (channels or pairs)
630 653 self.ylabel: label for Y axes
631 654 self.titles: list of axes title
632 655
633 656 '''
634 657 raise NotImplementedError
635 658
636 659 def plot(self):
637 660 '''
638 661 Must be defined in the child class, the actual plotting method
639 662 '''
640 663 raise NotImplementedError
641 664
642 665 def update(self, dataOut):
643 666 '''
644 667 Must be defined in the child class, update self.data with new data
645 668 '''
646 669
647 670 data = {
648 671 self.CODE: getattr(dataOut, 'data_{}'.format(self.CODE))
649 672 }
650 673 meta = {}
651 674
652 675 return data, meta
653 676
654 677 def run(self, dataOut, **kwargs):
655 678 '''
656 679 Main plotting routine
657 680 '''
658 681
659 682 if self.isConfig is False:
660 683 self.__setup(**kwargs)
661 684
662 685 if self.localtime:
663 686 self.getDateTime = datetime.datetime.fromtimestamp
664 687 else:
665 688 self.getDateTime = datetime.datetime.utcfromtimestamp
666 689
667 690 self.data.setup()
668 691 self.isConfig = True
669 692 if self.server:
670 693 self.context = zmq.Context()
671 694 self.socket = self.context.socket(zmq.REQ)
672 695 self.socket.connect(self.server)
673 696 self.poll = zmq.Poller()
674 697 self.poll.register(self.socket, zmq.POLLIN)
675 698
676 699 tm = getattr(dataOut, self.attr_time)
677 700
678 701 if self.data and 'time' in self.xaxis and (tm - self.tmin) >= self.xrange*60*60:
679 702 self.save_time = tm
680 703 self.__plot()
681 704 self.tmin += self.xrange*60*60
682 705 self.data.setup()
683 706 self.clear_figures()
684 707
685 708 self.__update(dataOut, tm)
686 709
687 710 if self.isPlotConfig is False:
688 711 self.__setup_plot()
689 712 self.isPlotConfig = True
690 713 if self.xaxis == 'time':
691 714 dt = self.getDateTime(tm)
692 715 if self.xmin is None:
693 716 self.tmin = tm
694 717 self.xmin = dt.hour
695 718 minutes = (self.xmin-int(self.xmin)) * 60
696 719 seconds = (minutes - int(minutes)) * 60
697 720 self.tmin = (dt.replace(hour=int(self.xmin), minute=int(minutes), second=int(seconds)) -
698 721 datetime.datetime(1970, 1, 1)).total_seconds()
699 722 if self.localtime:
700 723 self.tmin += time.timezone
701 724
702 725 if self.xmin is not None and self.xmax is not None:
703 726 self.xrange = self.xmax - self.xmin
704 727
705 728 if self.throttle == 0:
706 729 self.__plot()
707 730 else:
708 731 self.__throttle_plot(self.__plot)#, coerce=coerce)
709 732
710 733 def close(self):
711 734
712 735 if self.data and not self.data.flagNoData:
713 736 self.save_time = 0
714 737 self.__plot()
715 738 if self.data and not self.data.flagNoData and self.pause:
716 739 figpause(10)
@@ -1,1936 +1,1936
1 1 import os
2 2 import datetime
3 3 import numpy
4 4 from mpl_toolkits.axisartist.grid_finder import FixedLocator, DictFormatter
5 5
6 6 from schainpy.model.graphics.jroplot_base import Plot, plt
7 7 from schainpy.model.graphics.jroplot_spectra import SpectraPlot, RTIPlot, CoherencePlot, SpectraCutPlot
8 8 from schainpy.utils import log
9 9 # libreria wradlib
10 10 #import wradlib as wrl
11 11
12 12 EARTH_RADIUS = 6.3710e3
13 13
14 14
15 15 def ll2xy(lat1, lon1, lat2, lon2):
16 16
17 17 p = 0.017453292519943295
18 18 a = 0.5 - numpy.cos((lat2 - lat1) * p)/2 + numpy.cos(lat1 * p) * \
19 19 numpy.cos(lat2 * p) * (1 - numpy.cos((lon2 - lon1) * p)) / 2
20 20 r = 12742 * numpy.arcsin(numpy.sqrt(a))
21 21 theta = numpy.arctan2(numpy.sin((lon2-lon1)*p)*numpy.cos(lat2*p), numpy.cos(lat1*p)
22 22 * numpy.sin(lat2*p)-numpy.sin(lat1*p)*numpy.cos(lat2*p)*numpy.cos((lon2-lon1)*p))
23 23 theta = -theta + numpy.pi/2
24 24 return r*numpy.cos(theta), r*numpy.sin(theta)
25 25
26 26
27 27 def km2deg(km):
28 28 '''
29 29 Convert distance in km to degrees
30 30 '''
31 31
32 32 return numpy.rad2deg(km/EARTH_RADIUS)
33 33
34 34
35 35
36 36 class SpectralMomentsPlot(SpectraPlot):
37 37 '''
38 38 Plot for Spectral Moments
39 39 '''
40 40 CODE = 'spc_moments'
41 41 # colormap = 'jet'
42 42 # plot_type = 'pcolor'
43 43
44 44 class DobleGaussianPlot(SpectraPlot):
45 45 '''
46 46 Plot for Double Gaussian Plot
47 47 '''
48 48 CODE = 'gaussian_fit'
49 49 # colormap = 'jet'
50 50 # plot_type = 'pcolor'
51 51
52 52 class DoubleGaussianSpectraCutPlot(SpectraCutPlot):
53 53 '''
54 54 Plot SpectraCut with Double Gaussian Fit
55 55 '''
56 56 CODE = 'cut_gaussian_fit'
57 57
58 58 class SnrPlot(RTIPlot):
59 59 '''
60 60 Plot for SNR Data
61 61 '''
62 62
63 63 CODE = 'snr'
64 64 colormap = 'jet'
65 65
66 66 def update(self, dataOut):
67 67
68 68 data = {
69 69 'snr': 10*numpy.log10(dataOut.data_snr)
70 70 }
71 71
72 72 return data, {}
73 73
74 74 class DopplerPlot(RTIPlot):
75 75 '''
76 76 Plot for DOPPLER Data (1st moment)
77 77 '''
78 78
79 79 CODE = 'dop'
80 80 colormap = 'jet'
81 81
82 82 def update(self, dataOut):
83 83
84 84 data = {
85 85 'dop': 10*numpy.log10(dataOut.data_dop)
86 86 }
87 87
88 88 return data, {}
89 89
90 90 class PowerPlot(RTIPlot):
91 91 '''
92 92 Plot for Power Data (0 moment)
93 93 '''
94 94
95 95 CODE = 'pow'
96 96 colormap = 'jet'
97 97
98 98 def update(self, dataOut):
99 99 data = {
100 100 'pow': 10*numpy.log10(dataOut.data_pow/dataOut.normFactor)
101 101 }
102 102 return data, {}
103 103
104 104 class SpectralWidthPlot(RTIPlot):
105 105 '''
106 106 Plot for Spectral Width Data (2nd moment)
107 107 '''
108 108
109 109 CODE = 'width'
110 110 colormap = 'jet'
111 111
112 112 def update(self, dataOut):
113 113
114 114 data = {
115 115 'width': dataOut.data_width
116 116 }
117 117
118 118 return data, {}
119 119
120 120 class SkyMapPlot(Plot):
121 121 '''
122 122 Plot for meteors detection data
123 123 '''
124 124
125 125 CODE = 'param'
126 126
127 127 def setup(self):
128 128
129 129 self.ncols = 1
130 130 self.nrows = 1
131 131 self.width = 7.2
132 132 self.height = 7.2
133 133 self.nplots = 1
134 134 self.xlabel = 'Zonal Zenith Angle (deg)'
135 135 self.ylabel = 'Meridional Zenith Angle (deg)'
136 136 self.polar = True
137 137 self.ymin = -180
138 138 self.ymax = 180
139 139 self.colorbar = False
140 140
141 141 def plot(self):
142 142
143 143 arrayParameters = numpy.concatenate(self.data['param'])
144 144 error = arrayParameters[:, -1]
145 145 indValid = numpy.where(error == 0)[0]
146 146 finalMeteor = arrayParameters[indValid, :]
147 147 finalAzimuth = finalMeteor[:, 3]
148 148 finalZenith = finalMeteor[:, 4]
149 149
150 150 x = finalAzimuth * numpy.pi / 180
151 151 y = finalZenith
152 152
153 153 ax = self.axes[0]
154 154
155 155 if ax.firsttime:
156 156 ax.plot = ax.plot(x, y, 'bo', markersize=5)[0]
157 157 else:
158 158 ax.plot.set_data(x, y)
159 159
160 160 dt1 = self.getDateTime(self.data.min_time).strftime('%y/%m/%d %H:%M:%S')
161 161 dt2 = self.getDateTime(self.data.max_time).strftime('%y/%m/%d %H:%M:%S')
162 162 title = 'Meteor Detection Sky Map\n %s - %s \n Number of events: %5.0f\n' % (dt1,
163 163 dt2,
164 164 len(x))
165 165 self.titles[0] = title
166 166
167 167
168 168 class GenericRTIPlot(Plot):
169 169 '''
170 170 Plot for data_xxxx object
171 171 '''
172 172
173 173 CODE = 'param'
174 174 colormap = 'viridis'
175 175 plot_type = 'pcolorbuffer'
176 176
177 177 def setup(self):
178 178 self.xaxis = 'time'
179 179 self.ncols = 1
180 180 self.nrows = self.data.shape('param')[0]
181 181 self.nplots = self.nrows
182 182 self.plots_adjust.update({'hspace':0.8, 'left': 0.1, 'bottom': 0.08, 'right':0.95, 'top': 0.95})
183 183
184 184 if not self.xlabel:
185 185 self.xlabel = 'Time'
186 186
187 187 self.ylabel = 'Range [km]'
188 188 if not self.titles:
189 189 self.titles = ['Param {}'.format(x) for x in range(self.nrows)]
190 190
191 191 def update(self, dataOut):
192 192
193 193 data = {
194 194 'param' : numpy.concatenate([getattr(dataOut, attr) for attr in self.attr_data], axis=0)
195 195 }
196 196
197 197 meta = {}
198 198
199 199 return data, meta
200 200
201 201 def plot(self):
202 202 # self.data.normalize_heights()
203 203 self.x = self.data.times
204 204 self.y = self.data.yrange
205 205 self.z = self.data['param']
206 206 self.z = 10*numpy.log10(self.z)
207 207 self.z = numpy.ma.masked_invalid(self.z)
208 208
209 209 if self.decimation is None:
210 210 x, y, z = self.fill_gaps(self.x, self.y, self.z)
211 211 else:
212 212 x, y, z = self.fill_gaps(*self.decimate())
213 213
214 214 for n, ax in enumerate(self.axes):
215 215
216 216 self.zmax = self.zmax if self.zmax is not None else numpy.max(
217 217 self.z[n])
218 218 self.zmin = self.zmin if self.zmin is not None else numpy.min(
219 219 self.z[n])
220 220
221 221 if ax.firsttime:
222 222 if self.zlimits is not None:
223 223 self.zmin, self.zmax = self.zlimits[n]
224 224
225 225 ax.plt = ax.pcolormesh(x, y, z[n].T * self.factors[n],
226 226 vmin=self.zmin,
227 227 vmax=self.zmax,
228 228 cmap=self.cmaps[n]
229 229 )
230 230 else:
231 231 if self.zlimits is not None:
232 232 self.zmin, self.zmax = self.zlimits[n]
233 233 ax.collections.remove(ax.collections[0])
234 234 ax.plt = ax.pcolormesh(x, y, z[n].T * self.factors[n],
235 235 vmin=self.zmin,
236 236 vmax=self.zmax,
237 237 cmap=self.cmaps[n]
238 238 )
239 239
240 240
241 241 class PolarMapPlot(Plot):
242 242 '''
243 243 Plot for weather radar
244 244 '''
245 245
246 246 CODE = 'param'
247 247 colormap = 'seismic'
248 248
249 249 def setup(self):
250 250 self.ncols = 1
251 251 self.nrows = 1
252 252 self.width = 9
253 253 self.height = 8
254 254 self.mode = self.data.meta['mode']
255 255 if self.channels is not None:
256 256 self.nplots = len(self.channels)
257 257 self.nrows = len(self.channels)
258 258 else:
259 259 self.nplots = self.data.shape(self.CODE)[0]
260 260 self.nrows = self.nplots
261 261 self.channels = list(range(self.nplots))
262 262 if self.mode == 'E':
263 263 self.xlabel = 'Longitude'
264 264 self.ylabel = 'Latitude'
265 265 else:
266 266 self.xlabel = 'Range (km)'
267 267 self.ylabel = 'Height (km)'
268 268 self.bgcolor = 'white'
269 269 self.cb_labels = self.data.meta['units']
270 270 self.lat = self.data.meta['latitude']
271 271 self.lon = self.data.meta['longitude']
272 272 self.xmin, self.xmax = float(
273 273 km2deg(self.xmin) + self.lon), float(km2deg(self.xmax) + self.lon)
274 274 self.ymin, self.ymax = float(
275 275 km2deg(self.ymin) + self.lat), float(km2deg(self.ymax) + self.lat)
276 276 # self.polar = True
277 277
278 278 def plot(self):
279 279
280 280 for n, ax in enumerate(self.axes):
281 281 data = self.data['param'][self.channels[n]]
282 282
283 283 zeniths = numpy.linspace(
284 284 0, self.data.meta['max_range'], data.shape[1])
285 285 if self.mode == 'E':
286 286 azimuths = -numpy.radians(self.data.yrange)+numpy.pi/2
287 287 r, theta = numpy.meshgrid(zeniths, azimuths)
288 288 x, y = r*numpy.cos(theta)*numpy.cos(numpy.radians(self.data.meta['elevation'])), r*numpy.sin(
289 289 theta)*numpy.cos(numpy.radians(self.data.meta['elevation']))
290 290 x = km2deg(x) + self.lon
291 291 y = km2deg(y) + self.lat
292 292 else:
293 293 azimuths = numpy.radians(self.data.yrange)
294 294 r, theta = numpy.meshgrid(zeniths, azimuths)
295 295 x, y = r*numpy.cos(theta), r*numpy.sin(theta)
296 296 self.y = zeniths
297 297
298 298 if ax.firsttime:
299 299 if self.zlimits is not None:
300 300 self.zmin, self.zmax = self.zlimits[n]
301 301 ax.plt = ax.pcolormesh( # r, theta, numpy.ma.array(data, mask=numpy.isnan(data)),
302 302 x, y, numpy.ma.array(data, mask=numpy.isnan(data)),
303 303 vmin=self.zmin,
304 304 vmax=self.zmax,
305 305 cmap=self.cmaps[n])
306 306 else:
307 307 if self.zlimits is not None:
308 308 self.zmin, self.zmax = self.zlimits[n]
309 309 ax.collections.remove(ax.collections[0])
310 310 ax.plt = ax.pcolormesh( # r, theta, numpy.ma.array(data, mask=numpy.isnan(data)),
311 311 x, y, numpy.ma.array(data, mask=numpy.isnan(data)),
312 312 vmin=self.zmin,
313 313 vmax=self.zmax,
314 314 cmap=self.cmaps[n])
315 315
316 316 if self.mode == 'A':
317 317 continue
318 318
319 319 # plot district names
320 320 f = open('/data/workspace/schain_scripts/distrito.csv')
321 321 for line in f:
322 322 label, lon, lat = [s.strip() for s in line.split(',') if s]
323 323 lat = float(lat)
324 324 lon = float(lon)
325 325 # ax.plot(lon, lat, '.b', ms=2)
326 326 ax.text(lon, lat, label.decode('utf8'), ha='center',
327 327 va='bottom', size='8', color='black')
328 328
329 329 # plot limites
330 330 limites = []
331 331 tmp = []
332 332 for line in open('/data/workspace/schain_scripts/lima.csv'):
333 333 if '#' in line:
334 334 if tmp:
335 335 limites.append(tmp)
336 336 tmp = []
337 337 continue
338 338 values = line.strip().split(',')
339 339 tmp.append((float(values[0]), float(values[1])))
340 340 for points in limites:
341 341 ax.add_patch(
342 342 Polygon(points, ec='k', fc='none', ls='--', lw=0.5))
343 343
344 344 # plot Cuencas
345 345 for cuenca in ('rimac', 'lurin', 'mala', 'chillon', 'chilca', 'chancay-huaral'):
346 346 f = open('/data/workspace/schain_scripts/{}.csv'.format(cuenca))
347 347 values = [line.strip().split(',') for line in f]
348 348 points = [(float(s[0]), float(s[1])) for s in values]
349 349 ax.add_patch(Polygon(points, ec='b', fc='none'))
350 350
351 351 # plot grid
352 352 for r in (15, 30, 45, 60):
353 353 ax.add_artist(plt.Circle((self.lon, self.lat),
354 354 km2deg(r), color='0.6', fill=False, lw=0.2))
355 355 ax.text(
356 356 self.lon + (km2deg(r))*numpy.cos(60*numpy.pi/180),
357 357 self.lat + (km2deg(r))*numpy.sin(60*numpy.pi/180),
358 358 '{}km'.format(r),
359 359 ha='center', va='bottom', size='8', color='0.6', weight='heavy')
360 360
361 361 if self.mode == 'E':
362 362 title = 'El={}$^\circ$'.format(self.data.meta['elevation'])
363 363 label = 'E{:02d}'.format(int(self.data.meta['elevation']))
364 364 else:
365 365 title = 'Az={}$^\circ$'.format(self.data.meta['azimuth'])
366 366 label = 'A{:02d}'.format(int(self.data.meta['azimuth']))
367 367
368 368 self.save_labels = ['{}-{}'.format(lbl, label) for lbl in self.labels]
369 369 self.titles = ['{} {}'.format(
370 370 self.data.parameters[x], title) for x in self.channels]
371 371
372 372 class WeatherPlot(Plot):
373 373 CODE = 'weather'
374 374 plot_name = 'weather'
375 375 plot_type = 'ppistyle'
376 376 buffering = False
377 377
378 378 def setup(self):
379 379 self.ncols = 1
380 380 self.nrows = 1
381 381 self.width =8
382 382 self.height =8
383 383 self.nplots= 1
384 384 self.ylabel= 'Range [Km]'
385 385 self.titles= ['Weather']
386 386 self.colorbar=False
387 387 self.ini =0
388 388 self.len_azi =0
389 389 self.buffer_ini = None
390 390 self.buffer_azi = None
391 391 self.plots_adjust.update({'wspace': 0.4, 'hspace':0.4, 'left': 0.1, 'right': 0.9, 'bottom': 0.08})
392 392 self.flag =0
393 393 self.indicador= 0
394 394 self.last_data_azi = None
395 395 self.val_mean = None
396 396
397 397 def update(self, dataOut):
398 398
399 399 data = {}
400 400 meta = {}
401 401 if hasattr(dataOut, 'dataPP_POWER'):
402 402 factor = 1
403 403 if hasattr(dataOut, 'nFFTPoints'):
404 404 factor = dataOut.normFactor
405 405 #print("DIME EL SHAPE PORFAVOR",dataOut.data_360.shape)
406 406 data['weather'] = 10*numpy.log10(dataOut.data_360[1]/(factor))
407 407 data['azi'] = dataOut.data_azi
408 408 data['ele'] = dataOut.data_ele
409 409 return data, meta
410 410
411 411 def get2List(self,angulos):
412 412 list1=[]
413 413 list2=[]
414 414 for i in reversed(range(len(angulos))):
415 415 diff_ = angulos[i]-angulos[i-1]
416 416 if diff_ >1.5:
417 417 list1.append(i-1)
418 418 list2.append(diff_)
419 419 return list(reversed(list1)),list(reversed(list2))
420 420
421 421 def fixData360(self,list_,ang_):
422 422 if list_[0]==-1:
423 423 vec = numpy.where(ang_<ang_[0])
424 424 ang_[vec] = ang_[vec]+360
425 425 return ang_
426 426 return ang_
427 427
428 428 def fixData360HL(self,angulos):
429 429 vec = numpy.where(angulos>=360)
430 430 angulos[vec]=angulos[vec]-360
431 431 return angulos
432 432
433 433 def search_pos(self,pos,list_):
434 434 for i in range(len(list_)):
435 435 if pos == list_[i]:
436 436 return True,i
437 437 i=None
438 438 return False,i
439 439
440 440 def fixDataComp(self,ang_,list1_,list2_):
441 441 size = len(ang_)
442 442 size2 = 0
443 443 for i in range(len(list2_)):
444 444 size2=size2+round(list2_[i])-1
445 445 new_size= size+size2
446 446 ang_new = numpy.zeros(new_size)
447 447 ang_new2 = numpy.zeros(new_size)
448 448
449 449 tmp = 0
450 450 c = 0
451 451 for i in range(len(ang_)):
452 452 ang_new[tmp +c] = ang_[i]
453 453 ang_new2[tmp+c] = ang_[i]
454 454 condition , value = self.search_pos(i,list1_)
455 455 if condition:
456 456 pos = tmp + c + 1
457 457 for k in range(round(list2_[value])-1):
458 458 ang_new[pos+k] = ang_new[pos+k-1]+1
459 459 ang_new2[pos+k] = numpy.nan
460 460 tmp = pos +k
461 461 c = 0
462 462 c=c+1
463 463 return ang_new,ang_new2
464 464
465 465 def globalCheckPED(self,angulos):
466 466 l1,l2 = self.get2List(angulos)
467 467 if len(l1)>0:
468 468 angulos2 = self.fixData360(list_=l1,ang_=angulos)
469 469 l1,l2 = self.get2List(angulos2)
470 470
471 471 ang1_,ang2_ = self.fixDataComp(ang_=angulos2,list1_=l1,list2_=l2)
472 472 ang1_ = self.fixData360HL(ang1_)
473 473 ang2_ = self.fixData360HL(ang2_)
474 474 else:
475 475 ang1_= angulos
476 476 ang2_= angulos
477 477 return ang1_,ang2_
478 478
479 479 def analizeDATA(self,data_azi):
480 480 list1 = []
481 481 list2 = []
482 482 dat = data_azi
483 483 for i in reversed(range(1,len(dat))):
484 484 if dat[i]>dat[i-1]:
485 485 diff = int(dat[i])-int(dat[i-1])
486 486 else:
487 487 diff = 360+int(dat[i])-int(dat[i-1])
488 488 if diff > 1:
489 489 list1.append(i-1)
490 490 list2.append(diff-1)
491 491 return list1,list2
492 492
493 493 def fixDATANEW(self,data_azi,data_weather):
494 494 list1,list2 = self.analizeDATA(data_azi)
495 495 if len(list1)== 0:
496 496 return data_azi,data_weather
497 497 else:
498 498 resize = 0
499 499 for i in range(len(list2)):
500 500 resize= resize + list2[i]
501 501 new_data_azi = numpy.resize(data_azi,resize)
502 502 new_data_weather= numpy.resize(date_weather,resize)
503 503
504 504 for i in range(len(list2)):
505 505 j=0
506 506 position=list1[i]+1
507 507 for j in range(list2[i]):
508 508 new_data_azi[position+j]=new_data_azi[position+j-1]+1
509 509 return new_data_azi
510 510
511 511 def fixDATA(self,data_azi):
512 512 data=data_azi
513 513 for i in range(len(data)):
514 514 if numpy.isnan(data[i]):
515 515 data[i]=data[i-1]+1
516 516 return data
517 517
518 518 def replaceNAN(self,data_weather,data_azi,val):
519 519 data= data_azi
520 520 data_T= data_weather
521 521 if data.shape[0]> data_T.shape[0]:
522 522 data_N = numpy.ones( [data.shape[0],data_T.shape[1]])
523 523 c = 0
524 524 for i in range(len(data)):
525 525 if numpy.isnan(data[i]):
526 526 data_N[i,:]=numpy.ones(data_T.shape[1])*numpy.nan
527 527 else:
528 528 data_N[i,:]=data_T[c,:]
529 529 c=c+1
530 530 return data_N
531 531 else:
532 532 for i in range(len(data)):
533 533 if numpy.isnan(data[i]):
534 534 data_T[i,:]=numpy.ones(data_T.shape[1])*numpy.nan
535 535 return data_T
536 536
537 537 def const_ploteo(self,data_weather,data_azi,step,res):
538 538 if self.ini==0:
539 539 #-------
540 540 n = (360/res)-len(data_azi)
541 541 #--------------------- new -------------------------
542 542 data_azi_new ,data_azi_old= self.globalCheckPED(data_azi)
543 543 #------------------------
544 544 start = data_azi_new[-1] + res
545 545 end = data_azi_new[0] - res
546 546 #------ new
547 547 self.last_data_azi = end
548 548 if start>end:
549 549 end = end + 360
550 550 azi_vacia = numpy.linspace(start,end,int(n))
551 551 azi_vacia = numpy.where(azi_vacia>360,azi_vacia-360,azi_vacia)
552 552 data_azi = numpy.hstack((data_azi_new,azi_vacia))
553 553 # RADAR
554 554 val_mean = numpy.mean(data_weather[:,-1])
555 555 self.val_mean = val_mean
556 556 data_weather_cmp = numpy.ones([(360-data_weather.shape[0]),data_weather.shape[1]])*val_mean
557 557 data_weather = self.replaceNAN(data_weather=data_weather,data_azi=data_azi_old,val=self.val_mean)
558 558 data_weather = numpy.vstack((data_weather,data_weather_cmp))
559 559 else:
560 560 # azimuth
561 561 flag=0
562 562 start_azi = self.res_azi[0]
563 563 #-----------new------------
564 564 data_azi ,data_azi_old= self.globalCheckPED(data_azi)
565 565 data_weather = self.replaceNAN(data_weather=data_weather,data_azi=data_azi_old,val=self.val_mean)
566 566 #--------------------------
567 567 start = data_azi[0]
568 568 end = data_azi[-1]
569 569 self.last_data_azi= end
570 570 if start< start_azi:
571 571 start = start +360
572 572 if end <start_azi:
573 573 end = end +360
574 574
575 575 pos_ini = int((start-start_azi)/res)
576 576 len_azi = len(data_azi)
577 577 if (360-pos_ini)<len_azi:
578 578 if pos_ini+1==360:
579 579 pos_ini=0
580 580 else:
581 581 flag=1
582 582 dif= 360-pos_ini
583 583 comp= len_azi-dif
584 584 #-----------------
585 585 if flag==0:
586 586 # AZIMUTH
587 587 self.res_azi[pos_ini:pos_ini+len_azi] = data_azi
588 588 # RADAR
589 589 self.res_weather[pos_ini:pos_ini+len_azi,:] = data_weather
590 590 else:
591 591 # AZIMUTH
592 592 self.res_azi[pos_ini:pos_ini+dif] = data_azi[0:dif]
593 593 self.res_azi[0:comp] = data_azi[dif:]
594 594 # RADAR
595 595 self.res_weather[pos_ini:pos_ini+dif,:] = data_weather[0:dif,:]
596 596 self.res_weather[0:comp,:] = data_weather[dif:,:]
597 597 flag=0
598 598 data_azi = self.res_azi
599 599 data_weather = self.res_weather
600 600
601 601 return data_weather,data_azi
602 602
603 603 def plot(self):
604 604 thisDatetime = datetime.datetime.utcfromtimestamp(self.data.times[-1]).strftime('%Y-%m-%d %H:%M:%S')
605 605 data = self.data[-1]
606 606 r = self.data.yrange
607 607 delta_height = r[1]-r[0]
608 608 r_mask = numpy.where(r>=0)[0]
609 609 r = numpy.arange(len(r_mask))*delta_height
610 610 self.y = 2*r
611 611 # RADAR
612 612 #data_weather = data['weather']
613 613 # PEDESTAL
614 614 #data_azi = data['azi']
615 615 res = 1
616 616 # STEP
617 617 step = (360/(res*data['weather'].shape[0]))
618 618
619 619 self.res_weather, self.res_azi = self.const_ploteo(data_weather=data['weather'][:,r_mask],data_azi=data['azi'],step=step,res=res)
620 620 self.res_ele = numpy.mean(data['ele'])
621 621 ################# PLOTEO ###################
622 622 for i,ax in enumerate(self.axes):
623 623 self.zmin = self.zmin if self.zmin else 20
624 624 self.zmax = self.zmax if self.zmax else 80
625 625 if ax.firsttime:
626 626 plt.clf()
627 627 cgax, pm = wrl.vis.plot_ppi(self.res_weather,r=r,az=self.res_azi,fig=self.figures[0], proj='cg', vmin=self.zmin, vmax=self.zmax)
628 628 else:
629 629 plt.clf()
630 630 cgax, pm = wrl.vis.plot_ppi(self.res_weather,r=r,az=self.res_azi,fig=self.figures[0], proj='cg', vmin=self.zmin, vmax=self.zmax)
631 631 caax = cgax.parasites[0]
632 632 paax = cgax.parasites[1]
633 633 cbar = plt.gcf().colorbar(pm, pad=0.075)
634 634 caax.set_xlabel('x_range [km]')
635 635 caax.set_ylabel('y_range [km]')
636 636 plt.text(1.0, 1.05, 'Azimuth '+str(thisDatetime)+" Step "+str(self.ini)+ " EL: "+str(round(self.res_ele, 1)), transform=caax.transAxes, va='bottom',ha='right')
637 637
638 638 self.ini= self.ini+1
639 639
640 640
641 641 class WeatherRHIPlot(Plot):
642 642 CODE = 'weather'
643 643 plot_name = 'weather'
644 644 plot_type = 'rhistyle'
645 645 buffering = False
646 646 data_ele_tmp = None
647 647
648 648 def setup(self):
649 649 print("********************")
650 650 print("********************")
651 651 print("********************")
652 652 print("SETUP WEATHER PLOT")
653 653 self.ncols = 1
654 654 self.nrows = 1
655 655 self.nplots= 1
656 656 self.ylabel= 'Range [Km]'
657 657 self.titles= ['Weather']
658 658 if self.channels is not None:
659 659 self.nplots = len(self.channels)
660 660 self.nrows = len(self.channels)
661 661 else:
662 662 self.nplots = self.data.shape(self.CODE)[0]
663 663 self.nrows = self.nplots
664 664 self.channels = list(range(self.nplots))
665 665 print("channels",self.channels)
666 666 print("que saldra", self.data.shape(self.CODE)[0])
667 667 self.titles = ['{} Channel {}'.format(self.CODE.upper(), x) for x in range(self.nrows)]
668 668 print("self.titles",self.titles)
669 669 self.colorbar=False
670 670 self.width =12
671 671 self.height =8
672 672 self.ini =0
673 673 self.len_azi =0
674 674 self.buffer_ini = None
675 675 self.buffer_ele = None
676 676 self.plots_adjust.update({'wspace': 0.4, 'hspace':0.4, 'left': 0.1, 'right': 0.9, 'bottom': 0.08})
677 677 self.flag =0
678 678 self.indicador= 0
679 679 self.last_data_ele = None
680 680 self.val_mean = None
681 681
682 682 def update(self, dataOut):
683 683
684 684 data = {}
685 685 meta = {}
686 686 if hasattr(dataOut, 'dataPP_POWER'):
687 687 factor = 1
688 688 if hasattr(dataOut, 'nFFTPoints'):
689 689 factor = dataOut.normFactor
690 690 print("dataOut",dataOut.data_360.shape)
691 691 #
692 692 data['weather'] = 10*numpy.log10(dataOut.data_360/(factor))
693 693 #
694 694 #data['weather'] = 10*numpy.log10(dataOut.data_360[1]/(factor))
695 695 data['azi'] = dataOut.data_azi
696 696 data['ele'] = dataOut.data_ele
697 697 #print("UPDATE")
698 698 #print("data[weather]",data['weather'].shape)
699 699 #print("data[azi]",data['azi'])
700 700 return data, meta
701 701
702 702 def get2List(self,angulos):
703 703 list1=[]
704 704 list2=[]
705 705 for i in reversed(range(len(angulos))):
706 706 if not i==0:#el caso de i=0 evalula el primero de la lista con el ultimo y no es relevante
707 707 diff_ = angulos[i]-angulos[i-1]
708 708 if abs(diff_) >1.5:
709 709 list1.append(i-1)
710 710 list2.append(diff_)
711 711 return list(reversed(list1)),list(reversed(list2))
712 712
713 713 def fixData90(self,list_,ang_):
714 714 if list_[0]==-1:
715 715 vec = numpy.where(ang_<ang_[0])
716 716 ang_[vec] = ang_[vec]+90
717 717 return ang_
718 718 return ang_
719 719
720 720 def fixData90HL(self,angulos):
721 721 vec = numpy.where(angulos>=90)
722 722 angulos[vec]=angulos[vec]-90
723 723 return angulos
724 724
725 725
726 726 def search_pos(self,pos,list_):
727 727 for i in range(len(list_)):
728 728 if pos == list_[i]:
729 729 return True,i
730 730 i=None
731 731 return False,i
732 732
733 733 def fixDataComp(self,ang_,list1_,list2_,tipo_case):
734 734 size = len(ang_)
735 735 size2 = 0
736 736 for i in range(len(list2_)):
737 737 size2=size2+round(abs(list2_[i]))-1
738 738 new_size= size+size2
739 739 ang_new = numpy.zeros(new_size)
740 740 ang_new2 = numpy.zeros(new_size)
741 741
742 742 tmp = 0
743 743 c = 0
744 744 for i in range(len(ang_)):
745 745 ang_new[tmp +c] = ang_[i]
746 746 ang_new2[tmp+c] = ang_[i]
747 747 condition , value = self.search_pos(i,list1_)
748 748 if condition:
749 749 pos = tmp + c + 1
750 750 for k in range(round(abs(list2_[value]))-1):
751 751 if tipo_case==0 or tipo_case==3:#subida
752 752 ang_new[pos+k] = ang_new[pos+k-1]+1
753 753 ang_new2[pos+k] = numpy.nan
754 754 elif tipo_case==1 or tipo_case==2:#bajada
755 755 ang_new[pos+k] = ang_new[pos+k-1]-1
756 756 ang_new2[pos+k] = numpy.nan
757 757
758 758 tmp = pos +k
759 759 c = 0
760 760 c=c+1
761 761 return ang_new,ang_new2
762 762
763 763 def globalCheckPED(self,angulos,tipo_case):
764 764 l1,l2 = self.get2List(angulos)
765 765 ##print("l1",l1)
766 766 ##print("l2",l2)
767 767 if len(l1)>0:
768 768 #angulos2 = self.fixData90(list_=l1,ang_=angulos)
769 769 #l1,l2 = self.get2List(angulos2)
770 770 ang1_,ang2_ = self.fixDataComp(ang_=angulos,list1_=l1,list2_=l2,tipo_case=tipo_case)
771 771 #ang1_ = self.fixData90HL(ang1_)
772 772 #ang2_ = self.fixData90HL(ang2_)
773 773 else:
774 774 ang1_= angulos
775 775 ang2_= angulos
776 776 return ang1_,ang2_
777 777
778 778
779 779 def replaceNAN(self,data_weather,data_ele,val):
780 780 data= data_ele
781 781 data_T= data_weather
782 782 if data.shape[0]> data_T.shape[0]:
783 783 data_N = numpy.ones( [data.shape[0],data_T.shape[1]])
784 784 c = 0
785 785 for i in range(len(data)):
786 786 if numpy.isnan(data[i]):
787 787 data_N[i,:]=numpy.ones(data_T.shape[1])*numpy.nan
788 788 else:
789 789 data_N[i,:]=data_T[c,:]
790 790 c=c+1
791 791 return data_N
792 792 else:
793 793 for i in range(len(data)):
794 794 if numpy.isnan(data[i]):
795 795 data_T[i,:]=numpy.ones(data_T.shape[1])*numpy.nan
796 796 return data_T
797 797
798 798 def check_case(self,data_ele,ang_max,ang_min):
799 799 start = data_ele[0]
800 800 end = data_ele[-1]
801 801 number = (end-start)
802 802 len_ang=len(data_ele)
803 803 print("start",start)
804 804 print("end",end)
805 805 print("number",number)
806 806
807 807 print("len_ang",len_ang)
808 808
809 809 #exit(1)
810 810
811 811 if start<end and (round(abs(number)+1)>=len_ang or (numpy.argmin(data_ele)==0)):#caso subida
812 812 return 0
813 813 #elif start>end and (round(abs(number)+1)>=len_ang or(numpy.argmax(data_ele)==0)):#caso bajada
814 814 # return 1
815 815 elif round(abs(number)+1)>=len_ang and (start>end or(numpy.argmax(data_ele)==0)):#caso bajada
816 816 return 1
817 817 elif round(abs(number)+1)<len_ang and data_ele[-2]>data_ele[-1]:# caso BAJADA CAMBIO ANG MAX
818 818 return 2
819 819 elif round(abs(number)+1)<len_ang and data_ele[-2]<data_ele[-1] :# caso SUBIDA CAMBIO ANG MIN
820 820 return 3
821 821
822 822
823 823 def const_ploteo(self,val_ch,data_weather,data_ele,step,res,ang_max,ang_min):
824 824 ang_max= ang_max
825 825 ang_min= ang_min
826 826 data_weather=data_weather
827 827 val_ch=val_ch
828 828 ##print("*********************DATA WEATHER**************************************")
829 829 ##print(data_weather)
830 830 if self.ini==0:
831 831 '''
832 832 print("**********************************************")
833 833 print("**********************************************")
834 834 print("***************ini**************")
835 835 print("**********************************************")
836 836 print("**********************************************")
837 837 '''
838 838 #print("data_ele",data_ele)
839 839 #----------------------------------------------------------
840 840 tipo_case = self.check_case(data_ele,ang_max,ang_min)
841 841 print("check_case",tipo_case)
842 842 #exit(1)
843 843 #--------------------- new -------------------------
844 844 data_ele_new ,data_ele_old= self.globalCheckPED(data_ele,tipo_case)
845 845
846 846 #-------------------------CAMBIOS RHI---------------------------------
847 847 start= ang_min
848 848 end = ang_max
849 849 n= (ang_max-ang_min)/res
850 850 #------ new
851 851 self.start_data_ele = data_ele_new[0]
852 852 self.end_data_ele = data_ele_new[-1]
853 853 if tipo_case==0 or tipo_case==3: # SUBIDA
854 854 n1= round(self.start_data_ele)- start
855 855 n2= end - round(self.end_data_ele)
856 856 print(self.start_data_ele)
857 857 print(self.end_data_ele)
858 858 if n1>0:
859 859 ele1= numpy.linspace(ang_min+1,self.start_data_ele-1,n1)
860 860 ele1_nan= numpy.ones(n1)*numpy.nan
861 861 data_ele = numpy.hstack((ele1,data_ele_new))
862 862 print("ele1_nan",ele1_nan.shape)
863 863 print("data_ele_old",data_ele_old.shape)
864 864 data_ele_old = numpy.hstack((ele1_nan,data_ele_old))
865 865 if n2>0:
866 866 ele2= numpy.linspace(self.end_data_ele+1,end,n2)
867 867 ele2_nan= numpy.ones(n2)*numpy.nan
868 868 data_ele = numpy.hstack((data_ele,ele2))
869 869 print("ele2_nan",ele2_nan.shape)
870 870 print("data_ele_old",data_ele_old.shape)
871 871 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
872 872
873 873 if tipo_case==1 or tipo_case==2: # BAJADA
874 874 data_ele_new = data_ele_new[::-1] # reversa
875 875 data_ele_old = data_ele_old[::-1]# reversa
876 876 data_weather = data_weather[::-1,:]# reversa
877 877 vec= numpy.where(data_ele_new<ang_max)
878 878 data_ele_new = data_ele_new[vec]
879 879 data_ele_old = data_ele_old[vec]
880 880 data_weather = data_weather[vec[0]]
881 881 vec2= numpy.where(0<data_ele_new)
882 882 data_ele_new = data_ele_new[vec2]
883 883 data_ele_old = data_ele_old[vec2]
884 884 data_weather = data_weather[vec2[0]]
885 885 self.start_data_ele = data_ele_new[0]
886 886 self.end_data_ele = data_ele_new[-1]
887 887
888 888 n1= round(self.start_data_ele)- start
889 889 n2= end - round(self.end_data_ele)-1
890 890 print(self.start_data_ele)
891 891 print(self.end_data_ele)
892 892 if n1>0:
893 893 ele1= numpy.linspace(ang_min+1,self.start_data_ele-1,n1)
894 894 ele1_nan= numpy.ones(n1)*numpy.nan
895 895 data_ele = numpy.hstack((ele1,data_ele_new))
896 896 data_ele_old = numpy.hstack((ele1_nan,data_ele_old))
897 897 if n2>0:
898 898 ele2= numpy.linspace(self.end_data_ele+1,end,n2)
899 899 ele2_nan= numpy.ones(n2)*numpy.nan
900 900 data_ele = numpy.hstack((data_ele,ele2))
901 901 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
902 902 # RADAR
903 903 # NOTA data_ele y data_weather es la variable que retorna
904 904 val_mean = numpy.mean(data_weather[:,-1])
905 905 self.val_mean = val_mean
906 906 data_weather = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
907 907 self.data_ele_tmp[val_ch]= data_ele_old
908 908 else:
909 909 #print("**********************************************")
910 910 #print("****************VARIABLE**********************")
911 911 #-------------------------CAMBIOS RHI---------------------------------
912 912 #---------------------------------------------------------------------
913 913 ##print("INPUT data_ele",data_ele)
914 914 flag=0
915 915 start_ele = self.res_ele[0]
916 916 tipo_case = self.check_case(data_ele,ang_max,ang_min)
917 917 #print("TIPO DE DATA",tipo_case)
918 918 #-----------new------------
919 919 data_ele ,data_ele_old = self.globalCheckPED(data_ele,tipo_case)
920 920 data_weather = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
921 921
922 922 #-------------------------------NEW RHI ITERATIVO-------------------------
923 923
924 924 if tipo_case==0 : # SUBIDA
925 925 vec = numpy.where(data_ele<ang_max)
926 926 data_ele = data_ele[vec]
927 927 data_ele_old = data_ele_old[vec]
928 928 data_weather = data_weather[vec[0]]
929 929
930 930 vec2 = numpy.where(0<data_ele)
931 931 data_ele= data_ele[vec2]
932 932 data_ele_old= data_ele_old[vec2]
933 933 ##print(data_ele_new)
934 934 data_weather= data_weather[vec2[0]]
935 935
936 936 new_i_ele = int(round(data_ele[0]))
937 937 new_f_ele = int(round(data_ele[-1]))
938 938 #print(new_i_ele)
939 939 #print(new_f_ele)
940 940 #print(data_ele,len(data_ele))
941 941 #print(data_ele_old,len(data_ele_old))
942 942 if new_i_ele< 2:
943 943 self.data_ele_tmp[val_ch] = numpy.ones(ang_max-ang_min)*numpy.nan
944 944 self.res_weather[val_ch] = self.replaceNAN(data_weather=self.res_weather[val_ch],data_ele=self.data_ele_tmp[val_ch],val=self.val_mean)
945 945 self.data_ele_tmp[val_ch][new_i_ele:new_i_ele+len(data_ele)]=data_ele_old
946 946 self.res_ele[new_i_ele:new_i_ele+len(data_ele)]= data_ele
947 947 self.res_weather[val_ch][new_i_ele:new_i_ele+len(data_ele),:]= data_weather
948 948 data_ele = self.res_ele
949 949 data_weather = self.res_weather[val_ch]
950 950
951 951 elif tipo_case==1 : #BAJADA
952 952 data_ele = data_ele[::-1] # reversa
953 953 data_ele_old = data_ele_old[::-1]# reversa
954 954 data_weather = data_weather[::-1,:]# reversa
955 955 vec= numpy.where(data_ele<ang_max)
956 956 data_ele = data_ele[vec]
957 957 data_ele_old = data_ele_old[vec]
958 958 data_weather = data_weather[vec[0]]
959 959 vec2= numpy.where(0<data_ele)
960 960 data_ele = data_ele[vec2]
961 961 data_ele_old = data_ele_old[vec2]
962 962 data_weather = data_weather[vec2[0]]
963 963
964 964
965 965 new_i_ele = int(round(data_ele[0]))
966 966 new_f_ele = int(round(data_ele[-1]))
967 967 #print(data_ele)
968 968 #print(ang_max)
969 969 #print(data_ele_old)
970 970 if new_i_ele <= 1:
971 971 new_i_ele = 1
972 972 if round(data_ele[-1])>=ang_max-1:
973 973 self.data_ele_tmp[val_ch] = numpy.ones(ang_max-ang_min)*numpy.nan
974 974 self.res_weather[val_ch] = self.replaceNAN(data_weather=self.res_weather[val_ch],data_ele=self.data_ele_tmp[val_ch],val=self.val_mean)
975 975 self.data_ele_tmp[val_ch][new_i_ele-1:new_i_ele+len(data_ele)-1]=data_ele_old
976 976 self.res_ele[new_i_ele-1:new_i_ele+len(data_ele)-1]= data_ele
977 977 self.res_weather[val_ch][new_i_ele-1:new_i_ele+len(data_ele)-1,:]= data_weather
978 978 data_ele = self.res_ele
979 979 data_weather = self.res_weather[val_ch]
980 980
981 981 elif tipo_case==2: #bajada
982 982 vec = numpy.where(data_ele<ang_max)
983 983 data_ele = data_ele[vec]
984 984 data_weather= data_weather[vec[0]]
985 985
986 986 len_vec = len(vec)
987 987 data_ele_new = data_ele[::-1] # reversa
988 988 data_weather = data_weather[::-1,:]
989 989 new_i_ele = int(data_ele_new[0])
990 990 new_f_ele = int(data_ele_new[-1])
991 991
992 992 n1= new_i_ele- ang_min
993 993 n2= ang_max - new_f_ele-1
994 994 if n1>0:
995 995 ele1= numpy.linspace(ang_min+1,new_i_ele-1,n1)
996 996 ele1_nan= numpy.ones(n1)*numpy.nan
997 997 data_ele = numpy.hstack((ele1,data_ele_new))
998 998 data_ele_old = numpy.hstack((ele1_nan,data_ele_new))
999 999 if n2>0:
1000 1000 ele2= numpy.linspace(new_f_ele+1,ang_max,n2)
1001 1001 ele2_nan= numpy.ones(n2)*numpy.nan
1002 1002 data_ele = numpy.hstack((data_ele,ele2))
1003 1003 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
1004 1004
1005 1005 self.data_ele_tmp[val_ch] = data_ele_old
1006 1006 self.res_ele = data_ele
1007 1007 self.res_weather[val_ch] = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
1008 1008 data_ele = self.res_ele
1009 1009 data_weather = self.res_weather[val_ch]
1010 1010
1011 1011 elif tipo_case==3:#subida
1012 1012 vec = numpy.where(0<data_ele)
1013 1013 data_ele= data_ele[vec]
1014 1014 data_ele_new = data_ele
1015 1015 data_ele_old= data_ele_old[vec]
1016 1016 data_weather= data_weather[vec[0]]
1017 1017 pos_ini = numpy.argmin(data_ele)
1018 1018 if pos_ini>0:
1019 1019 len_vec= len(data_ele)
1020 1020 vec3 = numpy.linspace(pos_ini,len_vec-1,len_vec-pos_ini).astype(int)
1021 1021 #print(vec3)
1022 1022 data_ele= data_ele[vec3]
1023 1023 data_ele_new = data_ele
1024 1024 data_ele_old= data_ele_old[vec3]
1025 1025 data_weather= data_weather[vec3]
1026 1026
1027 1027 new_i_ele = int(data_ele_new[0])
1028 1028 new_f_ele = int(data_ele_new[-1])
1029 1029 n1= new_i_ele- ang_min
1030 1030 n2= ang_max - new_f_ele-1
1031 1031 if n1>0:
1032 1032 ele1= numpy.linspace(ang_min+1,new_i_ele-1,n1)
1033 1033 ele1_nan= numpy.ones(n1)*numpy.nan
1034 1034 data_ele = numpy.hstack((ele1,data_ele_new))
1035 1035 data_ele_old = numpy.hstack((ele1_nan,data_ele_new))
1036 1036 if n2>0:
1037 1037 ele2= numpy.linspace(new_f_ele+1,ang_max,n2)
1038 1038 ele2_nan= numpy.ones(n2)*numpy.nan
1039 1039 data_ele = numpy.hstack((data_ele,ele2))
1040 1040 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
1041 1041
1042 1042 self.data_ele_tmp[val_ch] = data_ele_old
1043 1043 self.res_ele = data_ele
1044 1044 self.res_weather[val_ch] = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
1045 1045 data_ele = self.res_ele
1046 1046 data_weather = self.res_weather[val_ch]
1047 1047 #print("self.data_ele_tmp",self.data_ele_tmp)
1048 1048 return data_weather,data_ele
1049 1049
1050 1050
1051 1051 def plot(self):
1052 1052 thisDatetime = datetime.datetime.utcfromtimestamp(self.data.times[-1]).strftime('%Y-%m-%d %H:%M:%S')
1053 1053 data = self.data[-1]
1054 1054 r = self.data.yrange
1055 1055 delta_height = r[1]-r[0]
1056 1056 r_mask = numpy.where(r>=0)[0]
1057 1057 ##print("delta_height",delta_height)
1058 1058 #print("r_mask",r_mask,len(r_mask))
1059 1059 r = numpy.arange(len(r_mask))*delta_height
1060 1060 self.y = 2*r
1061 1061 res = 1
1062 1062 ###print("data['weather'].shape[0]",data['weather'].shape[0])
1063 1063 ang_max = self.ang_max
1064 1064 ang_min = self.ang_min
1065 1065 var_ang =ang_max - ang_min
1066 1066 step = (int(var_ang)/(res*data['weather'].shape[0]))
1067 1067 ###print("step",step)
1068 1068 #--------------------------------------------------------
1069 1069 ##print('weather',data['weather'].shape)
1070 1070 ##print('ele',data['ele'].shape)
1071 1071
1072 1072 ###self.res_weather, self.res_ele = self.const_ploteo(data_weather=data['weather'][:,r_mask],data_ele=data['ele'],step=step,res=res,ang_max=ang_max,ang_min=ang_min)
1073 1073 ###self.res_azi = numpy.mean(data['azi'])
1074 1074 ###print("self.res_ele",self.res_ele)
1075 1075 plt.clf()
1076 1076 subplots = [121, 122]
1077 1077 cg={'angular_spacing': 20.}
1078 1078 if self.ini==0:
1079 1079 self.data_ele_tmp = numpy.ones([self.nplots,int(var_ang)])*numpy.nan
1080 1080 self.res_weather= numpy.ones([self.nplots,int(var_ang),len(r_mask)])*numpy.nan
1081 1081 print("SHAPE",self.data_ele_tmp.shape)
1082 1082
1083 1083 for i,ax in enumerate(self.axes):
1084 1084 self.res_weather[i], self.res_ele = self.const_ploteo(val_ch=i, data_weather=data['weather'][i][:,r_mask],data_ele=data['ele'],step=step,res=res,ang_max=ang_max,ang_min=ang_min)
1085 1085 self.res_azi = numpy.mean(data['azi'])
1086 1086 if i==0:
1087 1087 print("*****************************************************************************to plot**************************",self.res_weather[i].shape)
1088 1088 self.zmin = self.zmin if self.zmin else 20
1089 1089 self.zmax = self.zmax if self.zmax else 80
1090 1090 if ax.firsttime:
1091 1091 #plt.clf()
1092 1092 cgax, pm = wrl.vis.plot_rhi(self.res_weather[i],r=r,th=self.res_ele,ax=subplots[i], proj=cg,vmin=self.zmin, vmax=self.zmax)
1093 1093 #fig=self.figures[0]
1094 1094 else:
1095 1095 #plt.clf()
1096 1096 if i==0:
1097 1097 print(self.res_weather[i])
1098 1098 print(self.res_ele)
1099 1099 cgax, pm = wrl.vis.plot_rhi(self.res_weather[i],r=r,th=self.res_ele,ax=subplots[i], proj=cg,vmin=self.zmin, vmax=self.zmax)
1100 1100 caax = cgax.parasites[0]
1101 1101 paax = cgax.parasites[1]
1102 1102 cbar = plt.gcf().colorbar(pm, pad=0.075)
1103 1103 caax.set_xlabel('x_range [km]')
1104 1104 caax.set_ylabel('y_range [km]')
1105 1105 plt.text(1.0, 1.05, 'Elevacion '+str(thisDatetime)+" Step "+str(self.ini)+ " Azi: "+str(round(self.res_azi,2)), transform=caax.transAxes, va='bottom',ha='right')
1106 1106 print("***************************self.ini****************************",self.ini)
1107 1107 self.ini= self.ini+1
1108 1108
1109 1109 class Weather_vRF_Plot(Plot):
1110 1110 CODE = 'PPI'
1111 1111 plot_name = 'PPI'
1112 1112 #plot_type = 'ppistyle'
1113 1113 buffering = False
1114 1114
1115 1115 def setup(self):
1116 1116
1117 1117 self.ncols = 1
1118 1118 self.nrows = 1
1119 1119 self.width =8
1120 1120 self.height =8
1121 1121 self.nplots= 1
1122 1122 self.ylabel= 'Range [Km]'
1123 1123 self.xlabel= 'Range [Km]'
1124 1124 self.titles= ['PPI']
1125 1125 self.polar = True
1126 1126 if self.channels is not None:
1127 1127 self.nplots = len(self.channels)
1128 1128 self.nrows = len(self.channels)
1129 1129 else:
1130 1130 self.nplots = self.data.shape(self.CODE)[0]
1131 1131 self.nrows = self.nplots
1132 1132 self.channels = list(range(self.nplots))
1133 1133
1134 1134 if self.CODE == 'POWER':
1135 1135 self.cb_label = r'Power (dB)'
1136 1136 elif self.CODE == 'DOPPLER':
1137 1137 self.cb_label = r'Velocity (m/s)'
1138 1138 self.colorbar=True
1139 1139 self.width = 9
1140 1140 self.height =8
1141 1141 self.ini =0
1142 1142 self.len_azi =0
1143 1143 self.buffer_ini = None
1144 1144 self.buffer_ele = None
1145 1145 self.plots_adjust.update({'wspace': 0.4, 'hspace':0.4, 'left': 0.15, 'right': 0.9, 'bottom': 0.08})
1146 1146 self.flag =0
1147 1147 self.indicador= 0
1148 1148 self.last_data_ele = None
1149 1149 self.val_mean = None
1150 1150
1151 1151 def update(self, dataOut):
1152 1152
1153 1153 data = {}
1154 1154 meta = {}
1155 1155 if hasattr(dataOut, 'dataPP_POWER'):
1156 1156 factor = 1
1157 1157 if hasattr(dataOut, 'nFFTPoints'):
1158 1158 factor = dataOut.normFactor
1159 1159
1160 1160 if 'pow' in self.attr_data[0].lower():
1161 1161 data['data'] = 10*numpy.log10(getattr(dataOut, self.attr_data[0])/(factor))
1162 1162 else:
1163 1163 data['data'] = getattr(dataOut, self.attr_data[0])/(factor)
1164 1164
1165 1165 data['azi'] = dataOut.data_azi
1166 1166 data['ele'] = dataOut.data_ele
1167 1167
1168 1168 return data, meta
1169 1169
1170 1170 def plot(self):
1171 1171 data = self.data[-1]
1172 1172 r = self.data.yrange
1173 1173 delta_height = r[1]-r[0]
1174 1174 r_mask = numpy.where(r>=0)[0]
1175 1175 self.r_mask = r_mask
1176 1176 r = numpy.arange(len(r_mask))*delta_height
1177 1177 self.y = 2*r
1178 1178
1179 1179 try:
1180 1180 z = data['data'][self.channels[0]][:,r_mask]
1181 1181
1182 1182 except:
1183 1183 z = data['data'][0][:,r_mask]
1184 1184
1185 1185 self.titles = []
1186 1186
1187 1187 self.ymax = self.ymax if self.ymax else numpy.nanmax(r)
1188 1188 self.ymin = self.ymin if self.ymin else numpy.nanmin(r)
1189 1189 self.zmax = self.zmax if self.zmax else numpy.nanmax(z)
1190 1190 self.zmin = self.zmin if self.zmin else numpy.nanmin(z)
1191 1191 self.ang_min = self.ang_min if self.ang_min else 0
1192 1192 self.ang_max = self.ang_max if self.ang_max else 360
1193 1193
1194 1194 r, theta = numpy.meshgrid(r, numpy.radians(data['azi']) )
1195 1195
1196 1196 for i,ax in enumerate(self.axes):
1197 1197
1198 1198 if ax.firsttime:
1199 1199 ax.set_xlim(numpy.radians(self.ang_min),numpy.radians(self.ang_max))
1200 1200 ax.plt = ax.pcolormesh(theta, r, z, cmap=self.colormap, vmin=self.zmin, vmax=self.zmax)
1201 1201 ax.set_theta_direction(-1)
1202 1202
1203 1203 else:
1204 1204 ax.set_xlim(numpy.radians(self.ang_min),numpy.radians(self.ang_max))
1205 1205 ax.plt = ax.pcolormesh(theta, r, z, cmap=self.colormap, vmin=self.zmin, vmax=self.zmax)
1206 1206 ax.set_theta_direction(-1)
1207 1207
1208 1208 ax.grid(True)
1209 1209
1210 1210 if len(self.channels) !=1:
1211 1211 self.titles = ['PPI {} at EL: {} Channel {}'.format(self.self.labels[x], str(round(numpy.mean(data['ele']),1)), x) for x in range(self.nrows)]
1212 1212 else:
1213 1213 self.titles = ['PPI {} at EL: {} Channel {}'.format(self.labels[0], str(round(numpy.mean(data['ele']),1)), self.channels[0])]
1214 1214
1215 1215 class WeatherRHI_vRF2_Plot(Plot):
1216 1216 CODE = 'weather'
1217 1217 plot_name = 'weather'
1218 1218 plot_type = 'rhistyle'
1219 1219 buffering = False
1220 1220 data_ele_tmp = None
1221 1221
1222 1222 def setup(self):
1223 1223 print("********************")
1224 1224 print("********************")
1225 1225 print("********************")
1226 1226 print("SETUP WEATHER PLOT")
1227 1227 self.ncols = 1
1228 1228 self.nrows = 1
1229 1229 self.nplots= 1
1230 1230 self.ylabel= 'Range [Km]'
1231 1231 self.titles= ['Weather']
1232 1232 if self.channels is not None:
1233 1233 self.nplots = len(self.channels)
1234 1234 self.nrows = len(self.channels)
1235 1235 else:
1236 1236 self.nplots = self.data.shape(self.CODE)[0]
1237 1237 self.nrows = self.nplots
1238 1238 self.channels = list(range(self.nplots))
1239 1239 print("channels",self.channels)
1240 1240 print("que saldra", self.data.shape(self.CODE)[0])
1241 1241 self.titles = ['{} Channel {}'.format(self.CODE.upper(), x) for x in range(self.nrows)]
1242 1242 print("self.titles",self.titles)
1243 1243 self.colorbar=False
1244 1244 self.width =8
1245 1245 self.height =8
1246 1246 self.ini =0
1247 1247 self.len_azi =0
1248 1248 self.buffer_ini = None
1249 1249 self.buffer_ele = None
1250 1250 self.plots_adjust.update({'wspace': 0.4, 'hspace':0.4, 'left': 0.1, 'right': 0.9, 'bottom': 0.08})
1251 1251 self.flag =0
1252 1252 self.indicador= 0
1253 1253 self.last_data_ele = None
1254 1254 self.val_mean = None
1255 1255
1256 1256 def update(self, dataOut):
1257 1257
1258 1258 data = {}
1259 1259 meta = {}
1260 1260 if hasattr(dataOut, 'dataPP_POWER'):
1261 1261 factor = 1
1262 1262 if hasattr(dataOut, 'nFFTPoints'):
1263 1263 factor = dataOut.normFactor
1264 1264 print("dataOut",dataOut.data_360.shape)
1265 1265 #
1266 1266 data['weather'] = 10*numpy.log10(dataOut.data_360/(factor))
1267 1267 #
1268 1268 #data['weather'] = 10*numpy.log10(dataOut.data_360[1]/(factor))
1269 1269 data['azi'] = dataOut.data_azi
1270 1270 data['ele'] = dataOut.data_ele
1271 1271 data['case_flag'] = dataOut.case_flag
1272 1272 #print("UPDATE")
1273 1273 #print("data[weather]",data['weather'].shape)
1274 1274 #print("data[azi]",data['azi'])
1275 1275 return data, meta
1276 1276
1277 1277 def get2List(self,angulos):
1278 1278 list1=[]
1279 1279 list2=[]
1280 1280 for i in reversed(range(len(angulos))):
1281 1281 if not i==0:#el caso de i=0 evalula el primero de la lista con el ultimo y no es relevante
1282 1282 diff_ = angulos[i]-angulos[i-1]
1283 1283 if abs(diff_) >1.5:
1284 1284 list1.append(i-1)
1285 1285 list2.append(diff_)
1286 1286 return list(reversed(list1)),list(reversed(list2))
1287 1287
1288 1288 def fixData90(self,list_,ang_):
1289 1289 if list_[0]==-1:
1290 1290 vec = numpy.where(ang_<ang_[0])
1291 1291 ang_[vec] = ang_[vec]+90
1292 1292 return ang_
1293 1293 return ang_
1294 1294
1295 1295 def fixData90HL(self,angulos):
1296 1296 vec = numpy.where(angulos>=90)
1297 1297 angulos[vec]=angulos[vec]-90
1298 1298 return angulos
1299 1299
1300 1300
1301 1301 def search_pos(self,pos,list_):
1302 1302 for i in range(len(list_)):
1303 1303 if pos == list_[i]:
1304 1304 return True,i
1305 1305 i=None
1306 1306 return False,i
1307 1307
1308 1308 def fixDataComp(self,ang_,list1_,list2_,tipo_case):
1309 1309 size = len(ang_)
1310 1310 size2 = 0
1311 1311 for i in range(len(list2_)):
1312 1312 size2=size2+round(abs(list2_[i]))-1
1313 1313 new_size= size+size2
1314 1314 ang_new = numpy.zeros(new_size)
1315 1315 ang_new2 = numpy.zeros(new_size)
1316 1316
1317 1317 tmp = 0
1318 1318 c = 0
1319 1319 for i in range(len(ang_)):
1320 1320 ang_new[tmp +c] = ang_[i]
1321 1321 ang_new2[tmp+c] = ang_[i]
1322 1322 condition , value = self.search_pos(i,list1_)
1323 1323 if condition:
1324 1324 pos = tmp + c + 1
1325 1325 for k in range(round(abs(list2_[value]))-1):
1326 1326 if tipo_case==0 or tipo_case==3:#subida
1327 1327 ang_new[pos+k] = ang_new[pos+k-1]+1
1328 1328 ang_new2[pos+k] = numpy.nan
1329 1329 elif tipo_case==1 or tipo_case==2:#bajada
1330 1330 ang_new[pos+k] = ang_new[pos+k-1]-1
1331 1331 ang_new2[pos+k] = numpy.nan
1332 1332
1333 1333 tmp = pos +k
1334 1334 c = 0
1335 1335 c=c+1
1336 1336 return ang_new,ang_new2
1337 1337
1338 1338 def globalCheckPED(self,angulos,tipo_case):
1339 1339 l1,l2 = self.get2List(angulos)
1340 1340 ##print("l1",l1)
1341 1341 ##print("l2",l2)
1342 1342 if len(l1)>0:
1343 1343 #angulos2 = self.fixData90(list_=l1,ang_=angulos)
1344 1344 #l1,l2 = self.get2List(angulos2)
1345 1345 ang1_,ang2_ = self.fixDataComp(ang_=angulos,list1_=l1,list2_=l2,tipo_case=tipo_case)
1346 1346 #ang1_ = self.fixData90HL(ang1_)
1347 1347 #ang2_ = self.fixData90HL(ang2_)
1348 1348 else:
1349 1349 ang1_= angulos
1350 1350 ang2_= angulos
1351 1351 return ang1_,ang2_
1352 1352
1353 1353
1354 1354 def replaceNAN(self,data_weather,data_ele,val):
1355 1355 data= data_ele
1356 1356 data_T= data_weather
1357 1357 if data.shape[0]> data_T.shape[0]:
1358 1358 data_N = numpy.ones( [data.shape[0],data_T.shape[1]])
1359 1359 c = 0
1360 1360 for i in range(len(data)):
1361 1361 if numpy.isnan(data[i]):
1362 1362 data_N[i,:]=numpy.ones(data_T.shape[1])*numpy.nan
1363 1363 else:
1364 1364 data_N[i,:]=data_T[c,:]
1365 1365 c=c+1
1366 1366 return data_N
1367 1367 else:
1368 1368 for i in range(len(data)):
1369 1369 if numpy.isnan(data[i]):
1370 1370 data_T[i,:]=numpy.ones(data_T.shape[1])*numpy.nan
1371 1371 return data_T
1372 1372
1373 1373 def check_case(self,data_ele,ang_max,ang_min):
1374 1374 start = data_ele[0]
1375 1375 end = data_ele[-1]
1376 1376 number = (end-start)
1377 1377 len_ang=len(data_ele)
1378 1378 print("start",start)
1379 1379 print("end",end)
1380 1380 print("number",number)
1381 1381
1382 1382 print("len_ang",len_ang)
1383 1383
1384 1384 #exit(1)
1385 1385
1386 1386 if start<end and (round(abs(number)+1)>=len_ang or (numpy.argmin(data_ele)==0)):#caso subida
1387 1387 return 0
1388 1388 #elif start>end and (round(abs(number)+1)>=len_ang or(numpy.argmax(data_ele)==0)):#caso bajada
1389 1389 # return 1
1390 1390 elif round(abs(number)+1)>=len_ang and (start>end or(numpy.argmax(data_ele)==0)):#caso bajada
1391 1391 return 1
1392 1392 elif round(abs(number)+1)<len_ang and data_ele[-2]>data_ele[-1]:# caso BAJADA CAMBIO ANG MAX
1393 1393 return 2
1394 1394 elif round(abs(number)+1)<len_ang and data_ele[-2]<data_ele[-1] :# caso SUBIDA CAMBIO ANG MIN
1395 1395 return 3
1396 1396
1397 1397
1398 1398 def const_ploteo(self,val_ch,data_weather,data_ele,step,res,ang_max,ang_min,case_flag):
1399 1399 ang_max= ang_max
1400 1400 ang_min= ang_min
1401 1401 data_weather=data_weather
1402 1402 val_ch=val_ch
1403 1403 ##print("*********************DATA WEATHER**************************************")
1404 1404 ##print(data_weather)
1405 1405 if self.ini==0:
1406 1406 '''
1407 1407 print("**********************************************")
1408 1408 print("**********************************************")
1409 1409 print("***************ini**************")
1410 1410 print("**********************************************")
1411 1411 print("**********************************************")
1412 1412 '''
1413 1413 #print("data_ele",data_ele)
1414 1414 #----------------------------------------------------------
1415 1415 tipo_case = case_flag[-1]
1416 1416 #tipo_case = self.check_case(data_ele,ang_max,ang_min)
1417 1417 print("check_case",tipo_case)
1418 1418 #exit(1)
1419 1419 #--------------------- new -------------------------
1420 1420 data_ele_new ,data_ele_old= self.globalCheckPED(data_ele,tipo_case)
1421 1421
1422 1422 #-------------------------CAMBIOS RHI---------------------------------
1423 1423 start= ang_min
1424 1424 end = ang_max
1425 1425 n= (ang_max-ang_min)/res
1426 1426 #------ new
1427 1427 self.start_data_ele = data_ele_new[0]
1428 1428 self.end_data_ele = data_ele_new[-1]
1429 1429 if tipo_case==0 or tipo_case==3: # SUBIDA
1430 1430 n1= round(self.start_data_ele)- start
1431 1431 n2= end - round(self.end_data_ele)
1432 1432 print(self.start_data_ele)
1433 1433 print(self.end_data_ele)
1434 1434 if n1>0:
1435 1435 ele1= numpy.linspace(ang_min+1,self.start_data_ele-1,n1)
1436 1436 ele1_nan= numpy.ones(n1)*numpy.nan
1437 1437 data_ele = numpy.hstack((ele1,data_ele_new))
1438 1438 print("ele1_nan",ele1_nan.shape)
1439 1439 print("data_ele_old",data_ele_old.shape)
1440 1440 data_ele_old = numpy.hstack((ele1_nan,data_ele_old))
1441 1441 if n2>0:
1442 1442 ele2= numpy.linspace(self.end_data_ele+1,end,n2)
1443 1443 ele2_nan= numpy.ones(n2)*numpy.nan
1444 1444 data_ele = numpy.hstack((data_ele,ele2))
1445 1445 print("ele2_nan",ele2_nan.shape)
1446 1446 print("data_ele_old",data_ele_old.shape)
1447 1447 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
1448 1448
1449 1449 if tipo_case==1 or tipo_case==2: # BAJADA
1450 1450 data_ele_new = data_ele_new[::-1] # reversa
1451 1451 data_ele_old = data_ele_old[::-1]# reversa
1452 1452 data_weather = data_weather[::-1,:]# reversa
1453 1453 vec= numpy.where(data_ele_new<ang_max)
1454 1454 data_ele_new = data_ele_new[vec]
1455 1455 data_ele_old = data_ele_old[vec]
1456 1456 data_weather = data_weather[vec[0]]
1457 1457 vec2= numpy.where(0<data_ele_new)
1458 1458 data_ele_new = data_ele_new[vec2]
1459 1459 data_ele_old = data_ele_old[vec2]
1460 1460 data_weather = data_weather[vec2[0]]
1461 1461 self.start_data_ele = data_ele_new[0]
1462 1462 self.end_data_ele = data_ele_new[-1]
1463 1463
1464 1464 n1= round(self.start_data_ele)- start
1465 1465 n2= end - round(self.end_data_ele)-1
1466 1466 print(self.start_data_ele)
1467 1467 print(self.end_data_ele)
1468 1468 if n1>0:
1469 1469 ele1= numpy.linspace(ang_min+1,self.start_data_ele-1,n1)
1470 1470 ele1_nan= numpy.ones(n1)*numpy.nan
1471 1471 data_ele = numpy.hstack((ele1,data_ele_new))
1472 1472 data_ele_old = numpy.hstack((ele1_nan,data_ele_old))
1473 1473 if n2>0:
1474 1474 ele2= numpy.linspace(self.end_data_ele+1,end,n2)
1475 1475 ele2_nan= numpy.ones(n2)*numpy.nan
1476 1476 data_ele = numpy.hstack((data_ele,ele2))
1477 1477 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
1478 1478 # RADAR
1479 1479 # NOTA data_ele y data_weather es la variable que retorna
1480 1480 val_mean = numpy.mean(data_weather[:,-1])
1481 1481 self.val_mean = val_mean
1482 1482 data_weather = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
1483 1483 print("eleold",data_ele_old)
1484 1484 print(self.data_ele_tmp[val_ch])
1485 1485 print(data_ele_old.shape[0])
1486 1486 print(self.data_ele_tmp[val_ch].shape[0])
1487 1487 if (data_ele_old.shape[0]==91 or self.data_ele_tmp[val_ch].shape[0]==91):
1488 1488 import sys
1489 1489 print("EXIT",self.ini)
1490 1490
1491 1491 sys.exit(1)
1492 1492 self.data_ele_tmp[val_ch]= data_ele_old
1493 1493 else:
1494 1494 #print("**********************************************")
1495 1495 #print("****************VARIABLE**********************")
1496 1496 #-------------------------CAMBIOS RHI---------------------------------
1497 1497 #---------------------------------------------------------------------
1498 1498 ##print("INPUT data_ele",data_ele)
1499 1499 flag=0
1500 1500 start_ele = self.res_ele[0]
1501 1501 #tipo_case = self.check_case(data_ele,ang_max,ang_min)
1502 1502 tipo_case = case_flag[-1]
1503 1503 #print("TIPO DE DATA",tipo_case)
1504 1504 #-----------new------------
1505 1505 data_ele ,data_ele_old = self.globalCheckPED(data_ele,tipo_case)
1506 1506 data_weather = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
1507 1507
1508 1508 #-------------------------------NEW RHI ITERATIVO-------------------------
1509 1509
1510 1510 if tipo_case==0 : # SUBIDA
1511 1511 vec = numpy.where(data_ele<ang_max)
1512 1512 data_ele = data_ele[vec]
1513 1513 data_ele_old = data_ele_old[vec]
1514 1514 data_weather = data_weather[vec[0]]
1515 1515
1516 1516 vec2 = numpy.where(0<data_ele)
1517 1517 data_ele= data_ele[vec2]
1518 1518 data_ele_old= data_ele_old[vec2]
1519 1519 ##print(data_ele_new)
1520 1520 data_weather= data_weather[vec2[0]]
1521 1521
1522 1522 new_i_ele = int(round(data_ele[0]))
1523 1523 new_f_ele = int(round(data_ele[-1]))
1524 1524 #print(new_i_ele)
1525 1525 #print(new_f_ele)
1526 1526 #print(data_ele,len(data_ele))
1527 1527 #print(data_ele_old,len(data_ele_old))
1528 1528 if new_i_ele< 2:
1529 1529 self.data_ele_tmp[val_ch] = numpy.ones(ang_max-ang_min)*numpy.nan
1530 1530 self.res_weather[val_ch] = self.replaceNAN(data_weather=self.res_weather[val_ch],data_ele=self.data_ele_tmp[val_ch],val=self.val_mean)
1531 1531 self.data_ele_tmp[val_ch][new_i_ele:new_i_ele+len(data_ele)]=data_ele_old
1532 1532 self.res_ele[new_i_ele:new_i_ele+len(data_ele)]= data_ele
1533 1533 self.res_weather[val_ch][new_i_ele:new_i_ele+len(data_ele),:]= data_weather
1534 1534 data_ele = self.res_ele
1535 1535 data_weather = self.res_weather[val_ch]
1536 1536
1537 1537 elif tipo_case==1 : #BAJADA
1538 1538 data_ele = data_ele[::-1] # reversa
1539 1539 data_ele_old = data_ele_old[::-1]# reversa
1540 1540 data_weather = data_weather[::-1,:]# reversa
1541 1541 vec= numpy.where(data_ele<ang_max)
1542 1542 data_ele = data_ele[vec]
1543 1543 data_ele_old = data_ele_old[vec]
1544 1544 data_weather = data_weather[vec[0]]
1545 1545 vec2= numpy.where(0<data_ele)
1546 1546 data_ele = data_ele[vec2]
1547 1547 data_ele_old = data_ele_old[vec2]
1548 1548 data_weather = data_weather[vec2[0]]
1549 1549
1550 1550
1551 1551 new_i_ele = int(round(data_ele[0]))
1552 1552 new_f_ele = int(round(data_ele[-1]))
1553 1553 #print(data_ele)
1554 1554 #print(ang_max)
1555 1555 #print(data_ele_old)
1556 1556 if new_i_ele <= 1:
1557 1557 new_i_ele = 1
1558 1558 if round(data_ele[-1])>=ang_max-1:
1559 1559 self.data_ele_tmp[val_ch] = numpy.ones(ang_max-ang_min)*numpy.nan
1560 1560 self.res_weather[val_ch] = self.replaceNAN(data_weather=self.res_weather[val_ch],data_ele=self.data_ele_tmp[val_ch],val=self.val_mean)
1561 1561 self.data_ele_tmp[val_ch][new_i_ele-1:new_i_ele+len(data_ele)-1]=data_ele_old
1562 1562 self.res_ele[new_i_ele-1:new_i_ele+len(data_ele)-1]= data_ele
1563 1563 self.res_weather[val_ch][new_i_ele-1:new_i_ele+len(data_ele)-1,:]= data_weather
1564 1564 data_ele = self.res_ele
1565 1565 data_weather = self.res_weather[val_ch]
1566 1566
1567 1567 elif tipo_case==2: #bajada
1568 1568 vec = numpy.where(data_ele<ang_max)
1569 1569 data_ele = data_ele[vec]
1570 1570 data_weather= data_weather[vec[0]]
1571 1571
1572 1572 len_vec = len(vec)
1573 1573 data_ele_new = data_ele[::-1] # reversa
1574 1574 data_weather = data_weather[::-1,:]
1575 1575 new_i_ele = int(data_ele_new[0])
1576 1576 new_f_ele = int(data_ele_new[-1])
1577 1577
1578 1578 n1= new_i_ele- ang_min
1579 1579 n2= ang_max - new_f_ele-1
1580 1580 if n1>0:
1581 1581 ele1= numpy.linspace(ang_min+1,new_i_ele-1,n1)
1582 1582 ele1_nan= numpy.ones(n1)*numpy.nan
1583 1583 data_ele = numpy.hstack((ele1,data_ele_new))
1584 1584 data_ele_old = numpy.hstack((ele1_nan,data_ele_new))
1585 1585 if n2>0:
1586 1586 ele2= numpy.linspace(new_f_ele+1,ang_max,n2)
1587 1587 ele2_nan= numpy.ones(n2)*numpy.nan
1588 1588 data_ele = numpy.hstack((data_ele,ele2))
1589 1589 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
1590 1590
1591 1591 self.data_ele_tmp[val_ch] = data_ele_old
1592 1592 self.res_ele = data_ele
1593 1593 self.res_weather[val_ch] = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
1594 1594 data_ele = self.res_ele
1595 1595 data_weather = self.res_weather[val_ch]
1596 1596
1597 1597 elif tipo_case==3:#subida
1598 1598 vec = numpy.where(0<data_ele)
1599 1599 data_ele= data_ele[vec]
1600 1600 data_ele_new = data_ele
1601 1601 data_ele_old= data_ele_old[vec]
1602 1602 data_weather= data_weather[vec[0]]
1603 1603 pos_ini = numpy.argmin(data_ele)
1604 1604 if pos_ini>0:
1605 1605 len_vec= len(data_ele)
1606 1606 vec3 = numpy.linspace(pos_ini,len_vec-1,len_vec-pos_ini).astype(int)
1607 1607 #print(vec3)
1608 1608 data_ele= data_ele[vec3]
1609 1609 data_ele_new = data_ele
1610 1610 data_ele_old= data_ele_old[vec3]
1611 1611 data_weather= data_weather[vec3]
1612 1612
1613 1613 new_i_ele = int(data_ele_new[0])
1614 1614 new_f_ele = int(data_ele_new[-1])
1615 1615 n1= new_i_ele- ang_min
1616 1616 n2= ang_max - new_f_ele-1
1617 1617 if n1>0:
1618 1618 ele1= numpy.linspace(ang_min+1,new_i_ele-1,n1)
1619 1619 ele1_nan= numpy.ones(n1)*numpy.nan
1620 1620 data_ele = numpy.hstack((ele1,data_ele_new))
1621 1621 data_ele_old = numpy.hstack((ele1_nan,data_ele_new))
1622 1622 if n2>0:
1623 1623 ele2= numpy.linspace(new_f_ele+1,ang_max,n2)
1624 1624 ele2_nan= numpy.ones(n2)*numpy.nan
1625 1625 data_ele = numpy.hstack((data_ele,ele2))
1626 1626 data_ele_old = numpy.hstack((data_ele_old,ele2_nan))
1627 1627
1628 1628 self.data_ele_tmp[val_ch] = data_ele_old
1629 1629 self.res_ele = data_ele
1630 1630 self.res_weather[val_ch] = self.replaceNAN(data_weather=data_weather,data_ele=data_ele_old,val=self.val_mean)
1631 1631 data_ele = self.res_ele
1632 1632 data_weather = self.res_weather[val_ch]
1633 1633 #print("self.data_ele_tmp",self.data_ele_tmp)
1634 1634 return data_weather,data_ele
1635 1635
1636 1636
1637 1637 def plot(self):
1638 1638 thisDatetime = datetime.datetime.utcfromtimestamp(self.data.times[-1]).strftime('%Y-%m-%d %H:%M:%S')
1639 1639 data = self.data[-1]
1640 1640 r = self.data.yrange
1641 1641 delta_height = r[1]-r[0]
1642 1642 r_mask = numpy.where(r>=0)[0]
1643 1643 ##print("delta_height",delta_height)
1644 1644 #print("r_mask",r_mask,len(r_mask))
1645 1645 r = numpy.arange(len(r_mask))*delta_height
1646 1646 self.y = 2*r
1647 1647 res = 1
1648 1648 ###print("data['weather'].shape[0]",data['weather'].shape[0])
1649 1649 ang_max = self.ang_max
1650 1650 ang_min = self.ang_min
1651 1651 var_ang =ang_max - ang_min
1652 1652 step = (int(var_ang)/(res*data['weather'].shape[0]))
1653 1653 ###print("step",step)
1654 1654 #--------------------------------------------------------
1655 1655 ##print('weather',data['weather'].shape)
1656 1656 ##print('ele',data['ele'].shape)
1657 1657
1658 1658 ###self.res_weather, self.res_ele = self.const_ploteo(data_weather=data['weather'][:,r_mask],data_ele=data['ele'],step=step,res=res,ang_max=ang_max,ang_min=ang_min)
1659 1659 ###self.res_azi = numpy.mean(data['azi'])
1660 1660 ###print("self.res_ele",self.res_ele)
1661 1661 plt.clf()
1662 1662 subplots = [121, 122]
1663 1663 try:
1664 1664 if self.data[-2]['ele'].max()<data['ele'].max():
1665 1665 self.ini=0
1666 1666 except:
1667 1667 pass
1668 1668 if self.ini==0:
1669 1669 self.data_ele_tmp = numpy.ones([self.nplots,int(var_ang)])*numpy.nan
1670 1670 self.res_weather= numpy.ones([self.nplots,int(var_ang),len(r_mask)])*numpy.nan
1671 1671 print("SHAPE",self.data_ele_tmp.shape)
1672 1672
1673 1673 for i,ax in enumerate(self.axes):
1674 1674 self.res_weather[i], self.res_ele = self.const_ploteo(val_ch=i, data_weather=data['weather'][i][:,r_mask],data_ele=data['ele'],step=step,res=res,ang_max=ang_max,ang_min=ang_min,case_flag=self.data['case_flag'])
1675 1675 self.res_azi = numpy.mean(data['azi'])
1676 1676
1677 1677 if ax.firsttime:
1678 1678 #plt.clf()
1679 1679 print("Frist Plot")
1680 1680 cgax, pm = wrl.vis.plot_rhi(self.res_weather[i],r=r,th=self.res_ele,ax=subplots[i], proj='cg',vmin=20, vmax=80)
1681 1681 #fig=self.figures[0]
1682 1682 else:
1683 1683 #plt.clf()
1684 1684 print("ELSE PLOT")
1685 1685 cgax, pm = wrl.vis.plot_rhi(self.res_weather[i],r=r,th=self.res_ele,ax=subplots[i], proj='cg',vmin=20, vmax=80)
1686 1686 caax = cgax.parasites[0]
1687 1687 paax = cgax.parasites[1]
1688 1688 cbar = plt.gcf().colorbar(pm, pad=0.075)
1689 1689 caax.set_xlabel('x_range [km]')
1690 1690 caax.set_ylabel('y_range [km]')
1691 1691 plt.text(1.0, 1.05, 'Elevacion '+str(thisDatetime)+" Step "+str(self.ini)+ " Azi: "+str(round(self.res_azi,2)), transform=caax.transAxes, va='bottom',ha='right')
1692 1692 print("***************************self.ini****************************",self.ini)
1693 1693 self.ini= self.ini+1
1694 1694
1695 1695
1696 1696
1697 1697
1698 1698
1699 1699 class WeatherRHI_vRF4_Plot(Plot):
1700 1700 CODE = 'RHI'
1701 1701 plot_name = 'RHI'
1702 1702 #plot_type = 'rhistyle'
1703 1703 buffering = False
1704 1704
1705 1705 def setup(self):
1706 1706
1707 1707 self.ncols = 1
1708 1708 self.nrows = 1
1709 1709 self.nplots= 1
1710 1710 self.ylabel= 'Range [Km]'
1711 1711 self.xlabel= 'Range [Km]'
1712 1712 self.titles= ['RHI']
1713 1713 self.polar = True
1714 1714 self.grid = True
1715 1715 if self.channels is not None:
1716 1716 self.nplots = len(self.channels)
1717 1717 self.nrows = len(self.channels)
1718 1718 else:
1719 1719 self.nplots = self.data.shape(self.CODE)[0]
1720 1720 self.nrows = self.nplots
1721 1721 self.channels = list(range(self.nplots))
1722 1722
1723 1723 if self.CODE == 'Power':
1724 1724 self.cb_label = r'Power (dB)'
1725 1725 elif self.CODE == 'Doppler':
1726 1726 self.cb_label = r'Velocity (m/s)'
1727 1727 self.colorbar=True
1728 1728 self.width =8
1729 1729 self.height =8
1730 1730 self.ini =0
1731 1731 self.len_azi =0
1732 1732 self.buffer_ini = None
1733 1733 self.buffer_ele = None
1734 1734 self.plots_adjust.update({'wspace': 0.4, 'hspace':0.4, 'left': 0.1, 'right': 0.9, 'bottom': 0.08})
1735 1735 self.flag =0
1736 1736 self.indicador= 0
1737 1737 self.last_data_ele = None
1738 1738 self.val_mean = None
1739 1739
1740 1740 def update(self, dataOut):
1741 1741
1742 1742 data = {}
1743 1743 meta = {}
1744 1744 if hasattr(dataOut, 'dataPP_POWER'):
1745 1745 factor = 1
1746 1746 if hasattr(dataOut, 'nFFTPoints'):
1747 1747 factor = dataOut.normFactor
1748 1748
1749 1749 if 'pow' in self.attr_data[0].lower():
1750 1750 data['data'] = 10*numpy.log10(getattr(dataOut, self.attr_data[0])/(factor))
1751 1751 else:
1752 1752 data['data'] = getattr(dataOut, self.attr_data[0])/(factor)
1753 1753
1754 1754 data['azi'] = dataOut.data_azi
1755 1755 data['ele'] = dataOut.data_ele
1756 1756
1757 1757 return data, meta
1758 1758
1759 1759 def plot(self):
1760 1760 data = self.data[-1]
1761 1761 r = self.data.yrange
1762 1762 delta_height = r[1]-r[0]
1763 1763 r_mask = numpy.where(r>=0)[0]
1764 1764 self.r_mask =r_mask
1765 1765 r = numpy.arange(len(r_mask))*delta_height
1766 1766 self.y = 2*r
1767 1767
1768 1768 try:
1769 1769 z = data['data'][self.channels[0]][:,r_mask]
1770 1770 except:
1771 1771 z = data['data'][0][:,r_mask]
1772 1772
1773 1773 self.titles = []
1774 1774
1775 1775 self.ymax = self.ymax if self.ymax else numpy.nanmax(r)
1776 1776 self.ymin = self.ymin if self.ymin else numpy.nanmin(r)
1777 1777 self.zmax = self.zmax if self.zmax else numpy.nanmax(z)
1778 1778 self.zmin = self.zmin if self.zmin else numpy.nanmin(z)
1779 1779 self.ang_min = self.ang_min if self.ang_min else 0
1780 1780 self.ang_max = self.ang_max if self.ang_max else 90
1781 1781
1782 1782 r, theta = numpy.meshgrid(r, numpy.radians(data['ele']) )
1783 1783
1784 1784 for i,ax in enumerate(self.axes):
1785 1785
1786 1786 if ax.firsttime:
1787 1787 ax.set_xlim(numpy.radians(self.ang_min),numpy.radians(self.ang_max))
1788 1788 ax.plt = ax.pcolormesh(theta, r, z, cmap=self.colormap, vmin=self.zmin, vmax=self.zmax)
1789 1789
1790 1790 else:
1791 1791 ax.set_xlim(numpy.radians(self.ang_min),numpy.radians(self.ang_max))
1792 1792 ax.plt = ax.pcolormesh(theta, r, z, cmap=self.colormap, vmin=self.zmin, vmax=self.zmax)
1793 1793 ax.grid(True)
1794 1794 if len(self.channels) !=1:
1795 1795 self.titles = ['RHI {} at AZ: {} Channel {}'.format(self.labels[x], str(round(numpy.mean(data['azi']),1)), x) for x in range(self.nrows)]
1796 1796 else:
1797 1797 self.titles = ['RHI {} at AZ: {} Channel {}'.format(self.labels[0], str(round(numpy.mean(data['azi']),1)), self.channels[0])]
1798 1798
1799 1799 class WeatherParamsPlot(Plot):
1800 1800 #CODE = 'RHI'
1801 1801 #plot_name = 'RHI'
1802 1802 #plot_type = 'rhistyle'
1803 1803 buffering = False
1804 1804
1805 1805 def setup(self):
1806 1806
1807 1807 self.ncols = 1
1808 1808 self.nrows = 1
1809 1809 self.nplots= 1
1810 1810 self.ylabel= 'Range [km]'
1811 1811 self.xlabel= 'Range [km]'
1812 1812 self.polar = True
1813 1813 self.grid = True
1814 1814 if self.channels is not None:
1815 1815 self.nplots = len(self.channels)
1816 1816 self.nrows = len(self.channels)
1817 1817 else:
1818 1818 self.nplots = self.data.shape(self.CODE)[0]
1819 1819 self.nrows = self.nplots
1820 1820 self.channels = list(range(self.nplots))
1821 1821
1822 1822 self.colorbar=True
1823 1823 self.width =8
1824 1824 self.height =8
1825 1825 self.ini =0
1826 1826 self.len_azi =0
1827 1827 self.buffer_ini = None
1828 1828 self.buffer_ele = None
1829 self.plots_adjust.update({'wspace': 0.4, 'hspace':0.4, 'left': 0.08, 'right': 0.92, 'bottom': 0.01,'top':0.99})
1829 self.plots_adjust.update({'wspace': 0.4, 'hspace':0.4, 'left': 0.1, 'right': 0.9, 'bottom': 0.08})
1830 1830 self.flag =0
1831 1831 self.indicador= 0
1832 1832 self.last_data_ele = None
1833 1833 self.val_mean = None
1834 1834
1835 1835 def update(self, dataOut):
1836 1836
1837 1837 data = {}
1838 1838 meta = {}
1839 1839 if hasattr(dataOut, 'dataPP_POWER'):
1840 1840 factor = 1
1841 1841 if hasattr(dataOut, 'nFFTPoints'):
1842 1842 factor = dataOut.normFactor
1843 1843
1844 1844 if 'pow' in self.attr_data[0].lower():
1845 1845 data['data'] = 10*numpy.log10(getattr(dataOut, self.attr_data[0])/(factor))
1846 1846 else:
1847 1847 data['data'] = getattr(dataOut, self.attr_data[0])/(factor)
1848 1848
1849 1849 if dataOut.mode_op == 'PPI':
1850 1850 self.CODE = 'PPI'
1851 1851 self.title = self.CODE
1852 1852 elif dataOut.mode_op == 'RHI':
1853 1853 self.CODE = 'RHI'
1854 1854 self.title = self.CODE
1855 1855
1856 1856 data['azi'] = dataOut.data_azi
1857 1857 data['ele'] = dataOut.data_ele
1858 1858 data['mode_op'] = dataOut.mode_op
1859 1859
1860 1860 return data, meta
1861 1861
1862 1862 def plot(self):
1863 1863 data = self.data[-1]
1864 1864 r = self.data.yrange
1865 1865 delta_height = r[1]-r[0]
1866 1866 r_mask = numpy.where(r>=0)[0]
1867 1867 self.r_mask =r_mask
1868 1868 r = numpy.arange(len(r_mask))*delta_height
1869 1869 self.y = 2*r
1870 1870
1871 1871 try:
1872 1872 z = data['data'][self.channels[0]][:,r_mask]
1873 1873 except:
1874 1874 z = data['data'][0][:,r_mask]
1875 1875
1876 1876 self.titles = []
1877 1877
1878 1878 self.ymax = self.ymax if self.ymax else numpy.nanmax(r)
1879 1879 self.ymin = self.ymin if self.ymin else numpy.nanmin(r)
1880 1880 self.zmax = self.zmax if self.zmax else numpy.nanmax(z)
1881 1881 self.zmin = self.zmin if self.zmin else numpy.nanmin(z)
1882 1882 print("mode inside plot",self.data['mode_op'],data['mode_op'])
1883 1883 if data['mode_op'] == 'RHI':
1884 1884 try:
1885 1885 if self.data['mode_op'][-2] == 'PPI':
1886 1886 self.ang_min = None
1887 1887 self.ang_max = None
1888 1888 except:
1889 1889 pass
1890 1890 self.ang_min = self.ang_min if self.ang_min else 0
1891 1891 self.ang_max = self.ang_max if self.ang_max else 90
1892 1892 r, theta = numpy.meshgrid(r, numpy.radians(data['ele']) )
1893 1893 elif data['mode_op'] == 'PPI':
1894 1894 try:
1895 1895 if self.data['mode_op'][-2] == 'RHI':
1896 1896 self.ang_min = None
1897 1897 self.ang_max = None
1898 1898 except:
1899 1899 pass
1900 1900 self.ang_min = self.ang_min if self.ang_min else 0
1901 1901 self.ang_max = self.ang_max if self.ang_max else 360
1902 1902 r, theta = numpy.meshgrid(r, numpy.radians(data['azi']) )
1903 1903
1904 1904 self.clear_figures()
1905 1905
1906 1906 for i,ax in enumerate(self.axes):
1907 1907
1908 1908 if ax.firsttime:
1909 1909 ax.set_xlim(numpy.radians(self.ang_min),numpy.radians(self.ang_max))
1910 1910 ax.plt = ax.pcolormesh(theta, r, z, cmap=self.colormap, vmin=self.zmin, vmax=self.zmax)
1911 1911 if data['mode_op'] == 'PPI':
1912 1912 ax.set_theta_direction(-1)
1913 1913 ax.set_theta_offset(numpy.pi/2)
1914 1914
1915 1915 else:
1916 1916 ax.set_xlim(numpy.radians(self.ang_min),numpy.radians(self.ang_max))
1917 1917 ax.plt = ax.pcolormesh(theta, r, z, cmap=self.colormap, vmin=self.zmin, vmax=self.zmax)
1918 1918 if data['mode_op'] == 'PPI':
1919 1919 ax.set_theta_direction(-1)
1920 1920 ax.set_theta_offset(numpy.pi/2)
1921 1921
1922 1922 ax.grid(True)
1923 1923 if data['mode_op'] == 'RHI':
1924 1924 len_aux = int(data['azi'].shape[0]/4)
1925 1925 mean = numpy.mean(data['azi'][len_aux:-len_aux])
1926 1926 if len(self.channels) !=1:
1927 1927 self.titles = ['RHI {} at AZ: {} Channel {}'.format(self.labels[x], str(round(mean,1)), x) for x in range(self.nrows)]
1928 1928 else:
1929 1929 self.titles = ['RHI {} at AZ: {} Channel {}'.format(self.labels[0], str(round(mean,1)), self.channels[0])]
1930 1930 elif data['mode_op'] == 'PPI':
1931 1931 len_aux = int(data['ele'].shape[0]/4)
1932 1932 mean = numpy.mean(data['ele'][len_aux:-len_aux])
1933 1933 if len(self.channels) !=1:
1934 1934 self.titles = ['PPI {} at EL: {} Channel {}'.format(self.self.labels[x], str(round(mean,1)), x) for x in range(self.nrows)]
1935 1935 else:
1936 1936 self.titles = ['PPI {} at EL: {} Channel {}'.format(self.labels[0], str(round(mean,1)), self.channels[0])]
1 NO CONTENT: modified file
The requested commit or file is too big and content was truncated. Show full diff
General Comments 0
You need to be logged in to leave comments. Login now