##// END OF EJS Templates
schain
RSS
ax.plt.remove(), se añaden 3 lineas de codigo
avaldez -
r1785:6f1a660b8ab2
parent child
Show More
Side by Side Unified
Context file:
r1785:6f1a660b8ab2 Loading diff...:
Show file before | Show file after | Hide comments
@@ -1,494 +1,499
1 import os
1 import os
2 import datetime
2 import datetime
3 import numpy
3 import numpy
4
4
5 from schainpy.model.graphics.jroplot_base import Plot, plt
5 from schainpy.model.graphics.jroplot_base import Plot, plt
6 from schainpy.model.graphics.jroplot_spectra import SpectraPlot, RTIPlot, CoherencePlot, SpectraCutPlot
6 from schainpy.model.graphics.jroplot_spectra import SpectraPlot, RTIPlot, CoherencePlot, SpectraCutPlot
7 from schainpy.utils import log
7 from schainpy.utils import log
8
8
9 EARTH_RADIUS = 6.3710e3
9 EARTH_RADIUS = 6.3710e3
10
10
11
11
12 def ll2xy(lat1, lon1, lat2, lon2):
12 def ll2xy(lat1, lon1, lat2, lon2):
13
13
14 p = 0.017453292519943295
14 p = 0.017453292519943295
15 a = 0.5 - numpy.cos((lat2 - lat1) * p)/2 + numpy.cos(lat1 * p) * \
15 a = 0.5 - numpy.cos((lat2 - lat1) * p)/2 + numpy.cos(lat1 * p) * \
16 numpy.cos(lat2 * p) * (1 - numpy.cos((lon2 - lon1) * p)) / 2
16 numpy.cos(lat2 * p) * (1 - numpy.cos((lon2 - lon1) * p)) / 2
17 r = 12742 * numpy.arcsin(numpy.sqrt(a))
17 r = 12742 * numpy.arcsin(numpy.sqrt(a))
18 theta = numpy.arctan2(numpy.sin((lon2-lon1)*p)*numpy.cos(lat2*p), numpy.cos(lat1*p)
18 theta = numpy.arctan2(numpy.sin((lon2-lon1)*p)*numpy.cos(lat2*p), numpy.cos(lat1*p)
19 * numpy.sin(lat2*p)-numpy.sin(lat1*p)*numpy.cos(lat2*p)*numpy.cos((lon2-lon1)*p))
19 * numpy.sin(lat2*p)-numpy.sin(lat1*p)*numpy.cos(lat2*p)*numpy.cos((lon2-lon1)*p))
20 theta = -theta + numpy.pi/2
20 theta = -theta + numpy.pi/2
21 return r*numpy.cos(theta), r*numpy.sin(theta)
21 return r*numpy.cos(theta), r*numpy.sin(theta)
22
22
23
23
24 def km2deg(km):
24 def km2deg(km):
25 '''
25 '''
26 Convert distance in km to degrees
26 Convert distance in km to degrees
27 '''
27 '''
28
28
29 return numpy.rad2deg(km/EARTH_RADIUS)
29 return numpy.rad2deg(km/EARTH_RADIUS)
30
30
31
31
32
32
33 class SpectralMomentsPlot(SpectraPlot):
33 class SpectralMomentsPlot(SpectraPlot):
34 '''
34 '''
35 Plot for Spectral Moments
35 Plot for Spectral Moments
36 '''
36 '''
37 CODE = 'spc_moments'
37 CODE = 'spc_moments'
38 # colormap = 'jet'
38 # colormap = 'jet'
39 # plot_type = 'pcolor'
39 # plot_type = 'pcolor'
40
40
41 class DobleGaussianPlot(SpectraPlot):
41 class DobleGaussianPlot(SpectraPlot):
42 '''
42 '''
43 Plot for Double Gaussian Plot
43 Plot for Double Gaussian Plot
44 '''
44 '''
45 CODE = 'gaussian_fit'
45 CODE = 'gaussian_fit'
46 # colormap = 'jet'
46 # colormap = 'jet'
47 # plot_type = 'pcolor'
47 # plot_type = 'pcolor'
48
48
49
49
50 class DoubleGaussianSpectraCutPlot(SpectraCutPlot):
50 class DoubleGaussianSpectraCutPlot(SpectraCutPlot):
51 '''
51 '''
52 Plot SpectraCut with Double Gaussian Fit
52 Plot SpectraCut with Double Gaussian Fit
53 '''
53 '''
54 CODE = 'cut_gaussian_fit'
54 CODE = 'cut_gaussian_fit'
55
55
56
56
57 class SpectralFitObliquePlot(SpectraPlot):
57 class SpectralFitObliquePlot(SpectraPlot):
58 '''
58 '''
59 Plot for Spectral Oblique
59 Plot for Spectral Oblique
60 '''
60 '''
61 CODE = 'spc_moments'
61 CODE = 'spc_moments'
62 colormap = 'jet'
62 colormap = 'jet'
63 plot_type = 'pcolor'
63 plot_type = 'pcolor'
64
64
65
65
66
66
67 class SnrPlot(RTIPlot):
67 class SnrPlot(RTIPlot):
68 '''
68 '''
69 Plot for SNR Data
69 Plot for SNR Data
70 '''
70 '''
71
71
72 CODE = 'snr'
72 CODE = 'snr'
73 colormap = 'jet'
73 colormap = 'jet'
74
74
75 def update(self, dataOut):
75 def update(self, dataOut):
76
76
77 data = {
77 data = {
78 'snr': 10*numpy.log10(dataOut.data_snr)
78 'snr': 10*numpy.log10(dataOut.data_snr)
79 }
79 }
80
80
81 return data, {}
81 return data, {}
82
82
83 class DopplerPlot(RTIPlot):
83 class DopplerPlot(RTIPlot):
84 '''
84 '''
85 Plot for DOPPLER Data (1st moment)
85 Plot for DOPPLER Data (1st moment)
86 '''
86 '''
87
87
88 CODE = 'dop'
88 CODE = 'dop'
89 colormap = 'jet'
89 colormap = 'jet'
90
90
91 def update(self, dataOut):
91 def update(self, dataOut):
92
92
93 data = {
93 data = {
94 'dop': 10*numpy.log10(dataOut.data_dop)
94 'dop': 10*numpy.log10(dataOut.data_dop)
95 }
95 }
96
96
97 return data, {}
97 return data, {}
98
98
99 class DopplerEEJPlot_V0(RTIPlot):
99 class DopplerEEJPlot_V0(RTIPlot):
100 '''
100 '''
101 Written by R. Flores
101 Written by R. Flores
102 '''
102 '''
103 '''
103 '''
104 Plot for EEJ
104 Plot for EEJ
105 '''
105 '''
106
106
107 CODE = 'dop'
107 CODE = 'dop'
108 colormap = 'RdBu_r'
108 colormap = 'RdBu_r'
109 colormap = 'jet'
109 colormap = 'jet'
110
110
111 def setup(self):
111 def setup(self):
112
112
113 self.xaxis = 'time'
113 self.xaxis = 'time'
114 self.ncols = 1
114 self.ncols = 1
115 self.nrows = len(self.data.channels)
115 self.nrows = len(self.data.channels)
116 self.nplots = len(self.data.channels)
116 self.nplots = len(self.data.channels)
117 self.ylabel = 'Range [km]'
117 self.ylabel = 'Range [km]'
118 self.xlabel = 'Time'
118 self.xlabel = 'Time'
119 self.cb_label = '(m/s)'
119 self.cb_label = '(m/s)'
120 self.plots_adjust.update({'hspace':0.8, 'left': 0.1, 'bottom': 0.1, 'right':0.95})
120 self.plots_adjust.update({'hspace':0.8, 'left': 0.1, 'bottom': 0.1, 'right':0.95})
121 self.titles = ['{} Channel {}'.format(
121 self.titles = ['{} Channel {}'.format(
122 self.CODE.upper(), x) for x in range(self.nrows)]
122 self.CODE.upper(), x) for x in range(self.nrows)]
123
123
124 def update(self, dataOut):
124 def update(self, dataOut):
125 #print(self.EEJtype)
125 #print(self.EEJtype)
126
126
127 if self.EEJtype == 1:
127 if self.EEJtype == 1:
128 data = {
128 data = {
129 'dop': dataOut.Oblique_params[:,-2,:]
129 'dop': dataOut.Oblique_params[:,-2,:]
130 }
130 }
131 elif self.EEJtype == 2:
131 elif self.EEJtype == 2:
132 data = {
132 data = {
133 'dop': dataOut.Oblique_params[:,-1,:]
133 'dop': dataOut.Oblique_params[:,-1,:]
134 }
134 }
135
135
136 return data, {}
136 return data, {}
137
137
138 class DopplerEEJPlot(RTIPlot):
138 class DopplerEEJPlot(RTIPlot):
139 '''
139 '''
140 Written by R. Flores
140 Written by R. Flores
141 '''
141 '''
142 '''
142 '''
143 Plot for Doppler Shift EEJ
143 Plot for Doppler Shift EEJ
144 '''
144 '''
145
145
146 CODE = 'dop'
146 CODE = 'dop'
147 colormap = 'RdBu_r'
147 colormap = 'RdBu_r'
148 #colormap = 'jet'
148 #colormap = 'jet'
149
149
150 def setup(self):
150 def setup(self):
151
151
152 self.xaxis = 'time'
152 self.xaxis = 'time'
153 self.ncols = 1
153 self.ncols = 1
154 self.nrows = 2
154 self.nrows = 2
155 self.nplots = 2
155 self.nplots = 2
156 self.ylabel = 'Range [km]'
156 self.ylabel = 'Range [km]'
157 self.xlabel = 'Time'
157 self.xlabel = 'Time'
158 self.cb_label = '(m/s)'
158 self.cb_label = '(m/s)'
159 self.plots_adjust.update({'hspace':0.8, 'left': 0.1, 'bottom': 0.1, 'right':0.95})
159 self.plots_adjust.update({'hspace':0.8, 'left': 0.1, 'bottom': 0.1, 'right':0.95})
160 self.titles = ['{} EJJ Type {} /'.format(
160 self.titles = ['{} EJJ Type {} /'.format(
161 self.CODE.upper(), x) for x in range(1,1+self.nrows)]
161 self.CODE.upper(), x) for x in range(1,1+self.nrows)]
162
162
163 def update(self, dataOut):
163 def update(self, dataOut):
164
164
165 if dataOut.mode == 11: #Double Gaussian
165 if dataOut.mode == 11: #Double Gaussian
166 doppler = numpy.append(dataOut.Oblique_params[:,1,:],dataOut.Oblique_params[:,4,:],axis=0)
166 doppler = numpy.append(dataOut.Oblique_params[:,1,:],dataOut.Oblique_params[:,4,:],axis=0)
167 elif dataOut.mode == 9: #Double Skew Gaussian
167 elif dataOut.mode == 9: #Double Skew Gaussian
168 doppler = numpy.append(dataOut.Oblique_params[:,-2,:],dataOut.Oblique_params[:,-1,:],axis=0)
168 doppler = numpy.append(dataOut.Oblique_params[:,-2,:],dataOut.Oblique_params[:,-1,:],axis=0)
169 data = {
169 data = {
170 'dop': doppler
170 'dop': doppler
171 }
171 }
172
172
173 return data, {}
173 return data, {}
174
174
175 class SpcWidthEEJPlot(RTIPlot):
175 class SpcWidthEEJPlot(RTIPlot):
176 '''
176 '''
177 Written by R. Flores
177 Written by R. Flores
178 '''
178 '''
179 '''
179 '''
180 Plot for EEJ Spectral Width
180 Plot for EEJ Spectral Width
181 '''
181 '''
182
182
183 CODE = 'width'
183 CODE = 'width'
184 colormap = 'RdBu_r'
184 colormap = 'RdBu_r'
185 colormap = 'jet'
185 colormap = 'jet'
186
186
187 def setup(self):
187 def setup(self):
188
188
189 self.xaxis = 'time'
189 self.xaxis = 'time'
190 self.ncols = 1
190 self.ncols = 1
191 self.nrows = 2
191 self.nrows = 2
192 self.nplots = 2
192 self.nplots = 2
193 self.ylabel = 'Range [km]'
193 self.ylabel = 'Range [km]'
194 self.xlabel = 'Time'
194 self.xlabel = 'Time'
195 self.cb_label = '(m/s)'
195 self.cb_label = '(m/s)'
196 self.plots_adjust.update({'hspace':0.8, 'left': 0.1, 'bottom': 0.1, 'right':0.95})
196 self.plots_adjust.update({'hspace':0.8, 'left': 0.1, 'bottom': 0.1, 'right':0.95})
197 self.titles = ['{} EJJ Type {} /'.format(
197 self.titles = ['{} EJJ Type {} /'.format(
198 self.CODE.upper(), x) for x in range(1,1+self.nrows)]
198 self.CODE.upper(), x) for x in range(1,1+self.nrows)]
199
199
200 def update(self, dataOut):
200 def update(self, dataOut):
201
201
202 if dataOut.mode == 11: #Double Gaussian
202 if dataOut.mode == 11: #Double Gaussian
203 width = numpy.append(dataOut.Oblique_params[:,2,:],dataOut.Oblique_params[:,5,:],axis=0)
203 width = numpy.append(dataOut.Oblique_params[:,2,:],dataOut.Oblique_params[:,5,:],axis=0)
204 elif dataOut.mode == 9: #Double Skew Gaussian
204 elif dataOut.mode == 9: #Double Skew Gaussian
205 width = numpy.append(dataOut.Oblique_params[:,2,:],dataOut.Oblique_params[:,6,:],axis=0)
205 width = numpy.append(dataOut.Oblique_params[:,2,:],dataOut.Oblique_params[:,6,:],axis=0)
206 data = {
206 data = {
207 'width': width
207 'width': width
208 }
208 }
209
209
210 return data, {}
210 return data, {}
211
211
212 class PowerPlot(RTIPlot):
212 class PowerPlot(RTIPlot):
213 '''
213 '''
214 Plot for Power Data (0 moment)
214 Plot for Power Data (0 moment)
215 '''
215 '''
216
216
217 CODE = 'pow'
217 CODE = 'pow'
218 colormap = 'jet'
218 colormap = 'jet'
219
219
220 def update(self, dataOut):
220 def update(self, dataOut):
221
221
222 data = {
222 data = {
223 'pow': 10*numpy.log10(dataOut.data_pow/dataOut.normFactor)
223 'pow': 10*numpy.log10(dataOut.data_pow/dataOut.normFactor)
224 }
224 }
225
225
226 return data, {}
226 return data, {}
227
227
228 class SpectralWidthPlot(RTIPlot):
228 class SpectralWidthPlot(RTIPlot):
229 '''
229 '''
230 Plot for Spectral Width Data (2nd moment)
230 Plot for Spectral Width Data (2nd moment)
231 '''
231 '''
232
232
233 CODE = 'width'
233 CODE = 'width'
234 colormap = 'jet'
234 colormap = 'jet'
235
235
236 def update(self, dataOut):
236 def update(self, dataOut):
237
237
238 data = {
238 data = {
239 'width': dataOut.data_width
239 'width': dataOut.data_width
240 }
240 }
241
241
242 return data, {}
242 return data, {}
243
243
244 class SkyMapPlot(Plot):
244 class SkyMapPlot(Plot):
245 '''
245 '''
246 Plot for meteors detection data
246 Plot for meteors detection data
247 '''
247 '''
248
248
249 CODE = 'param'
249 CODE = 'param'
250
250
251 def setup(self):
251 def setup(self):
252
252
253 self.ncols = 1
253 self.ncols = 1
254 self.nrows = 1
254 self.nrows = 1
255 self.width = 7.2
255 self.width = 7.2
256 self.height = 7.2
256 self.height = 7.2
257 self.nplots = 1
257 self.nplots = 1
258 self.xlabel = 'Zonal Zenith Angle (deg)'
258 self.xlabel = 'Zonal Zenith Angle (deg)'
259 self.ylabel = 'Meridional Zenith Angle (deg)'
259 self.ylabel = 'Meridional Zenith Angle (deg)'
260 self.polar = True
260 self.polar = True
261 self.ymin = -180
261 self.ymin = -180
262 self.ymax = 180
262 self.ymax = 180
263 self.colorbar = False
263 self.colorbar = False
264
264
265 def plot(self):
265 def plot(self):
266
266
267 arrayParameters = numpy.concatenate(self.data['param'])
267 arrayParameters = numpy.concatenate(self.data['param'])
268 error = arrayParameters[:, -1]
268 error = arrayParameters[:, -1]
269 indValid = numpy.where(error == 0)[0]
269 indValid = numpy.where(error == 0)[0]
270 finalMeteor = arrayParameters[indValid, :]
270 finalMeteor = arrayParameters[indValid, :]
271 finalAzimuth = finalMeteor[:, 3]
271 finalAzimuth = finalMeteor[:, 3]
272 finalZenith = finalMeteor[:, 4]
272 finalZenith = finalMeteor[:, 4]
273
273
274 x = finalAzimuth * numpy.pi / 180
274 x = finalAzimuth * numpy.pi / 180
275 y = finalZenith
275 y = finalZenith
276
276
277 ax = self.axes[0]
277 ax = self.axes[0]
278
278
279 if ax.firsttime:
279 if ax.firsttime:
280 ax.plot = ax.plot(x, y, 'bo', markersize=5)[0]
280 ax.plot = ax.plot(x, y, 'bo', markersize=5)[0]
281 else:
281 else:
282 ax.plot.set_data(x, y)
282 ax.plot.set_data(x, y)
283
283
284 dt1 = self.getDateTime(self.data.min_time).strftime('%y/%m/%d %H:%M:%S')
284 dt1 = self.getDateTime(self.data.min_time).strftime('%y/%m/%d %H:%M:%S')
285 dt2 = self.getDateTime(self.data.max_time).strftime('%y/%m/%d %H:%M:%S')
285 dt2 = self.getDateTime(self.data.max_time).strftime('%y/%m/%d %H:%M:%S')
286 title = 'Meteor Detection Sky Map\n %s - %s \n Number of events: %5.0f\n' % (dt1,
286 title = 'Meteor Detection Sky Map\n %s - %s \n Number of events: %5.0f\n' % (dt1,
287 dt2,
287 dt2,
288 len(x))
288 len(x))
289 self.titles[0] = title
289 self.titles[0] = title
290
290
291
291
292 class GenericRTIPlot(Plot):
292 class GenericRTIPlot(Plot):
293 '''
293 '''
294 Plot for data_xxxx object
294 Plot for data_xxxx object
295 '''
295 '''
296
296
297 CODE = 'param'
297 CODE = 'param'
298 colormap = 'viridis'
298 colormap = 'viridis'
299 plot_type = 'pcolorbuffer'
299 plot_type = 'pcolorbuffer'
300
300
301 def setup(self):
301 def setup(self):
302 self.xaxis = 'time'
302 self.xaxis = 'time'
303 self.ncols = 1
303 self.ncols = 1
304 self.nrows = self.data.shape('param')[0]
304 self.nrows = self.data.shape('param')[0]
305 self.nplots = self.nrows
305 self.nplots = self.nrows
306 self.plots_adjust.update({'hspace':0.8, 'left': 0.1, 'bottom': 0.08, 'right':0.95, 'top': 0.95})
306 self.plots_adjust.update({'hspace':0.8, 'left': 0.1, 'bottom': 0.08, 'right':0.95, 'top': 0.95})
307
307
308 if not self.xlabel:
308 if not self.xlabel:
309 self.xlabel = 'Time'
309 self.xlabel = 'Time'
310
310
311 self.ylabel = 'Range [km]'
311 self.ylabel = 'Range [km]'
312 if not self.titles:
312 if not self.titles:
313 self.titles = ['Param {}'.format(x) for x in range(self.nrows)]
313 self.titles = ['Param {}'.format(x) for x in range(self.nrows)]
314
314
315 def update(self, dataOut):
315 def update(self, dataOut):
316
316
317 data = {
317 data = {
318 'param' : numpy.concatenate([getattr(dataOut, attr) for attr in self.attr_data], axis=0)
318 'param' : numpy.concatenate([getattr(dataOut, attr) for attr in self.attr_data], axis=0)
319 }
319 }
320
320
321 meta = {}
321 meta = {}
322
322
323 return data, meta
323 return data, meta
324
324
325 def plot(self):
325 def plot(self):
326 # self.data.normalize_heights()
326 # self.data.normalize_heights()
327 self.x = self.data.times
327 self.x = self.data.times
328 self.y = self.data.yrange
328 self.y = self.data.yrange
329 self.z = self.data['param']
329 self.z = self.data['param']
330
330
331 self.z = numpy.ma.masked_invalid(self.z)
331 self.z = numpy.ma.masked_invalid(self.z)
332
332
333 if self.decimation is None:
333 if self.decimation is None:
334 x, y, z = self.fill_gaps(self.x, self.y, self.z)
334 x, y, z = self.fill_gaps(self.x, self.y, self.z)
335 else:
335 else:
336 x, y, z = self.fill_gaps(*self.decimate())
336 x, y, z = self.fill_gaps(*self.decimate())
337
337
338 for n, ax in enumerate(self.axes):
338 for n, ax in enumerate(self.axes):
339
339
340 self.zmax = self.zmax if self.zmax is not None else numpy.max(
340 self.zmax = self.zmax if self.zmax is not None else numpy.max(
341 self.z[n])
341 self.z[n])
342 self.zmin = self.zmin if self.zmin is not None else numpy.min(
342 self.zmin = self.zmin if self.zmin is not None else numpy.min(
343 self.z[n])
343 self.z[n])
344
344
345 if ax.firsttime:
345 if ax.firsttime:
346 if self.zlimits is not None:
346 if self.zlimits is not None:
347 self.zmin, self.zmax = self.zlimits[n]
347 self.zmin, self.zmax = self.zlimits[n]
348
348
349 ax.plt = ax.pcolormesh(x, y, z[n].T * self.factors[n],
349 ax.plt = ax.pcolormesh(x, y, z[n].T * self.factors[n],
350 vmin=self.zmin,
350 vmin=self.zmin,
351 vmax=self.zmax,
351 vmax=self.zmax,
352 cmap=self.cmaps[n]
352 cmap=self.cmaps[n]
353 )
353 )
354 else:
354 else:
355 if self.zlimits is not None:
355 if self.zlimits is not None:
356 self.zmin, self.zmax = self.zlimits[n]
356 self.zmin, self.zmax = self.zlimits[n]
357 ax.plt.remove()
357
358 try:
359 ax.collections.remove(ax.collections[0])
360 except:
361 pass
362 # ax.plt.remove()
358 ax.plt = ax.pcolormesh(x, y, z[n].T * self.factors[n],
363 ax.plt = ax.pcolormesh(x, y, z[n].T * self.factors[n],
359 vmin=self.zmin,
364 vmin=self.zmin,
360 vmax=self.zmax,
365 vmax=self.zmax,
361 cmap=self.cmaps[n]
366 cmap=self.cmaps[n]
362 )
367 )
363
368
364
369
365 class PolarMapPlot(Plot):
370 class PolarMapPlot(Plot):
366 '''
371 '''
367 Plot for weather radar
372 Plot for weather radar
368 '''
373 '''
369
374
370 CODE = 'param'
375 CODE = 'param'
371 colormap = 'seismic'
376 colormap = 'seismic'
372
377
373 def setup(self):
378 def setup(self):
374 self.ncols = 1
379 self.ncols = 1
375 self.nrows = 1
380 self.nrows = 1
376 self.width = 9
381 self.width = 9
377 self.height = 8
382 self.height = 8
378 self.mode = self.data.meta['mode']
383 self.mode = self.data.meta['mode']
379 if self.channels is not None:
384 if self.channels is not None:
380 self.nplots = len(self.channels)
385 self.nplots = len(self.channels)
381 self.nrows = len(self.channels)
386 self.nrows = len(self.channels)
382 else:
387 else:
383 self.nplots = self.data.shape(self.CODE)[0]
388 self.nplots = self.data.shape(self.CODE)[0]
384 self.nrows = self.nplots
389 self.nrows = self.nplots
385 self.channels = list(range(self.nplots))
390 self.channels = list(range(self.nplots))
386 if self.mode == 'E':
391 if self.mode == 'E':
387 self.xlabel = 'Longitude'
392 self.xlabel = 'Longitude'
388 self.ylabel = 'Latitude'
393 self.ylabel = 'Latitude'
389 else:
394 else:
390 self.xlabel = 'Range (km)'
395 self.xlabel = 'Range (km)'
391 self.ylabel = 'Height (km)'
396 self.ylabel = 'Height (km)'
392 self.bgcolor = 'white'
397 self.bgcolor = 'white'
393 self.cb_labels = self.data.meta['units']
398 self.cb_labels = self.data.meta['units']
394 self.lat = self.data.meta['latitude']
399 self.lat = self.data.meta['latitude']
395 self.lon = self.data.meta['longitude']
400 self.lon = self.data.meta['longitude']
396 self.xmin, self.xmax = float(
401 self.xmin, self.xmax = float(
397 km2deg(self.xmin) + self.lon), float(km2deg(self.xmax) + self.lon)
402 km2deg(self.xmin) + self.lon), float(km2deg(self.xmax) + self.lon)
398 self.ymin, self.ymax = float(
403 self.ymin, self.ymax = float(
399 km2deg(self.ymin) + self.lat), float(km2deg(self.ymax) + self.lat)
404 km2deg(self.ymin) + self.lat), float(km2deg(self.ymax) + self.lat)
400 # self.polar = True
405 # self.polar = True
401
406
402 def plot(self):
407 def plot(self):
403
408
404 for n, ax in enumerate(self.axes):
409 for n, ax in enumerate(self.axes):
405 data = self.data['param'][self.channels[n]]
410 data = self.data['param'][self.channels[n]]
406
411
407 zeniths = numpy.linspace(
412 zeniths = numpy.linspace(
408 0, self.data.meta['max_range'], data.shape[1])
413 0, self.data.meta['max_range'], data.shape[1])
409 if self.mode == 'E':
414 if self.mode == 'E':
410 azimuths = -numpy.radians(self.data.yrange)+numpy.pi/2
415 azimuths = -numpy.radians(self.data.yrange)+numpy.pi/2
411 r, theta = numpy.meshgrid(zeniths, azimuths)
416 r, theta = numpy.meshgrid(zeniths, azimuths)
412 x, y = r*numpy.cos(theta)*numpy.cos(numpy.radians(self.data.meta['elevation'])), r*numpy.sin(
417 x, y = r*numpy.cos(theta)*numpy.cos(numpy.radians(self.data.meta['elevation'])), r*numpy.sin(
413 theta)*numpy.cos(numpy.radians(self.data.meta['elevation']))
418 theta)*numpy.cos(numpy.radians(self.data.meta['elevation']))
414 x = km2deg(x) + self.lon
419 x = km2deg(x) + self.lon
415 y = km2deg(y) + self.lat
420 y = km2deg(y) + self.lat
416 else:
421 else:
417 azimuths = numpy.radians(self.data.yrange)
422 azimuths = numpy.radians(self.data.yrange)
418 r, theta = numpy.meshgrid(zeniths, azimuths)
423 r, theta = numpy.meshgrid(zeniths, azimuths)
419 x, y = r*numpy.cos(theta), r*numpy.sin(theta)
424 x, y = r*numpy.cos(theta), r*numpy.sin(theta)
420 self.y = zeniths
425 self.y = zeniths
421
426
422 if ax.firsttime:
427 if ax.firsttime:
423 if self.zlimits is not None:
428 if self.zlimits is not None:
424 self.zmin, self.zmax = self.zlimits[n]
429 self.zmin, self.zmax = self.zlimits[n]
425 ax.plt = ax.pcolormesh( # r, theta, numpy.ma.array(data, mask=numpy.isnan(data)),
430 ax.plt = ax.pcolormesh( # r, theta, numpy.ma.array(data, mask=numpy.isnan(data)),
426 x, y, numpy.ma.array(data, mask=numpy.isnan(data)),
431 x, y, numpy.ma.array(data, mask=numpy.isnan(data)),
427 vmin=self.zmin,
432 vmin=self.zmin,
428 vmax=self.zmax,
433 vmax=self.zmax,
429 cmap=self.cmaps[n])
434 cmap=self.cmaps[n])
430 else:
435 else:
431 if self.zlimits is not None:
436 if self.zlimits is not None:
432 self.zmin, self.zmax = self.zlimits[n]
437 self.zmin, self.zmax = self.zlimits[n]
433 ax.plt.remove()
438 ax.plt.remove()
434 ax.plt = ax.pcolormesh( # r, theta, numpy.ma.array(data, mask=numpy.isnan(data)),
439 ax.plt = ax.pcolormesh( # r, theta, numpy.ma.array(data, mask=numpy.isnan(data)),
435 x, y, numpy.ma.array(data, mask=numpy.isnan(data)),
440 x, y, numpy.ma.array(data, mask=numpy.isnan(data)),
436 vmin=self.zmin,
441 vmin=self.zmin,
437 vmax=self.zmax,
442 vmax=self.zmax,
438 cmap=self.cmaps[n])
443 cmap=self.cmaps[n])
439
444
440 if self.mode == 'A':
445 if self.mode == 'A':
441 continue
446 continue
442
447
443 # plot district names
448 # plot district names
444 f = open('/data/workspace/schain_scripts/distrito.csv')
449 f = open('/data/workspace/schain_scripts/distrito.csv')
445 for line in f:
450 for line in f:
446 label, lon, lat = [s.strip() for s in line.split(',') if s]
451 label, lon, lat = [s.strip() for s in line.split(',') if s]
447 lat = float(lat)
452 lat = float(lat)
448 lon = float(lon)
453 lon = float(lon)
449 # ax.plot(lon, lat, '.b', ms=2)
454 # ax.plot(lon, lat, '.b', ms=2)
450 ax.text(lon, lat, label.decode('utf8'), ha='center',
455 ax.text(lon, lat, label.decode('utf8'), ha='center',
451 va='bottom', size='8', color='black')
456 va='bottom', size='8', color='black')
452
457
453 # plot limites
458 # plot limites
454 limites = []
459 limites = []
455 tmp = []
460 tmp = []
456 for line in open('/data/workspace/schain_scripts/lima.csv'):
461 for line in open('/data/workspace/schain_scripts/lima.csv'):
457 if '#' in line:
462 if '#' in line:
458 if tmp:
463 if tmp:
459 limites.append(tmp)
464 limites.append(tmp)
460 tmp = []
465 tmp = []
461 continue
466 continue
462 values = line.strip().split(',')
467 values = line.strip().split(',')
463 tmp.append((float(values[0]), float(values[1])))
468 tmp.append((float(values[0]), float(values[1])))
464 for points in limites:
469 for points in limites:
465 ax.add_patch(
470 ax.add_patch(
466 Polygon(points, ec='k', fc='none', ls='--', lw=0.5))
471 Polygon(points, ec='k', fc='none', ls='--', lw=0.5))
467
472
468 # plot Cuencas
473 # plot Cuencas
469 for cuenca in ('rimac', 'lurin', 'mala', 'chillon', 'chilca', 'chancay-huaral'):
474 for cuenca in ('rimac', 'lurin', 'mala', 'chillon', 'chilca', 'chancay-huaral'):
470 f = open('/data/workspace/schain_scripts/{}.csv'.format(cuenca))
475 f = open('/data/workspace/schain_scripts/{}.csv'.format(cuenca))
471 values = [line.strip().split(',') for line in f]
476 values = [line.strip().split(',') for line in f]
472 points = [(float(s[0]), float(s[1])) for s in values]
477 points = [(float(s[0]), float(s[1])) for s in values]
473 ax.add_patch(Polygon(points, ec='b', fc='none'))
478 ax.add_patch(Polygon(points, ec='b', fc='none'))
474
479
475 # plot grid
480 # plot grid
476 for r in (15, 30, 45, 60):
481 for r in (15, 30, 45, 60):
477 ax.add_artist(plt.Circle((self.lon, self.lat),
482 ax.add_artist(plt.Circle((self.lon, self.lat),
478 km2deg(r), color='0.6', fill=False, lw=0.2))
483 km2deg(r), color='0.6', fill=False, lw=0.2))
479 ax.text(
484 ax.text(
480 self.lon + (km2deg(r))*numpy.cos(60*numpy.pi/180),
485 self.lon + (km2deg(r))*numpy.cos(60*numpy.pi/180),
481 self.lat + (km2deg(r))*numpy.sin(60*numpy.pi/180),
486 self.lat + (km2deg(r))*numpy.sin(60*numpy.pi/180),
482 '{}km'.format(r),
487 '{}km'.format(r),
483 ha='center', va='bottom', size='8', color='0.6', weight='heavy')
488 ha='center', va='bottom', size='8', color='0.6', weight='heavy')
484
489
485 if self.mode == 'E':
490 if self.mode == 'E':
486 title = 'El={}$^\circ$'.format(self.data.meta['elevation'])
491 title = 'El={}$^\circ$'.format(self.data.meta['elevation'])
487 label = 'E{:02d}'.format(int(self.data.meta['elevation']))
492 label = 'E{:02d}'.format(int(self.data.meta['elevation']))
488 else:
493 else:
489 title = 'Az={}$^\circ$'.format(self.data.meta['azimuth'])
494 title = 'Az={}$^\circ$'.format(self.data.meta['azimuth'])
490 label = 'A{:02d}'.format(int(self.data.meta['azimuth']))
495 label = 'A{:02d}'.format(int(self.data.meta['azimuth']))
491
496
492 self.save_labels = ['{}-{}'.format(lbl, label) for lbl in self.labels]
497 self.save_labels = ['{}-{}'.format(lbl, label) for lbl in self.labels]
493 self.titles = ['{} {}'.format(
498 self.titles = ['{} {}'.format(
494 self.data.parameters[x], title) for x in self.channels]
499 self.data.parameters[x], title) for x in self.channels]
General Comments 0
You need to be logged in to leave comments. Login now