##// END OF EJS Templates
update y fix
avaldez -
r1257:f1a6137b3999
parent child
Show More
@@ -1,1736 +1,1737
1 '''
1 '''
2 Created on Jul 9, 2014
2 Created on Jul 9, 2014
3
3
4 @author: roj-idl71
4 @author: roj-idl71
5 '''
5 '''
6 import os
6 import os
7 import datetime
7 import datetime
8 import numpy
8 import numpy
9
9
10 from figure import Figure, isRealtime, isTimeInHourRange
10 from figure import Figure, isRealtime, isTimeInHourRange
11 from plotting_codes import *
11 from plotting_codes import *
12
12
13
13
14 class SpectraPlot(Figure):
14 class SpectraPlot(Figure):
15
15
16 isConfig = None
16 isConfig = None
17 __nsubplots = None
17 __nsubplots = None
18
18
19 WIDTHPROF = None
19 WIDTHPROF = None
20 HEIGHTPROF = None
20 HEIGHTPROF = None
21 PREFIX = 'spc'
21 PREFIX = 'spc'
22
22
23 def __init__(self, **kwargs):
23 def __init__(self, **kwargs):
24 Figure.__init__(self, **kwargs)
24 Figure.__init__(self, **kwargs)
25 self.isConfig = False
25 self.isConfig = False
26 self.__nsubplots = 1
26 self.__nsubplots = 1
27
27
28 self.WIDTH = 250
28 self.WIDTH = 250
29 self.HEIGHT = 250
29 self.HEIGHT = 250
30 self.WIDTHPROF = 120
30 self.WIDTHPROF = 120
31 self.HEIGHTPROF = 0
31 self.HEIGHTPROF = 0
32 self.counter_imagwr = 0
32 self.counter_imagwr = 0
33
33
34 self.PLOT_CODE = SPEC_CODE
34 self.PLOT_CODE = SPEC_CODE
35
35
36 self.FTP_WEI = None
36 self.FTP_WEI = None
37 self.EXP_CODE = None
37 self.EXP_CODE = None
38 self.SUB_EXP_CODE = None
38 self.SUB_EXP_CODE = None
39 self.PLOT_POS = None
39 self.PLOT_POS = None
40
40
41 self.__xfilter_ena = False
41 self.__xfilter_ena = False
42 self.__yfilter_ena = False
42 self.__yfilter_ena = False
43
43
44 def getSubplots(self):
44 def getSubplots(self):
45
45
46 ncol = int(numpy.sqrt(self.nplots)+0.9)
46 ncol = int(numpy.sqrt(self.nplots)+0.9)
47 nrow = int(self.nplots*1./ncol + 0.9)
47 nrow = int(self.nplots*1./ncol + 0.9)
48
48
49 return nrow, ncol
49 return nrow, ncol
50
50
51 def setup(self, id, nplots, wintitle, showprofile=True, show=True):
51 def setup(self, id, nplots, wintitle, showprofile=True, show=True):
52
52
53 self.__showprofile = showprofile
53 self.__showprofile = showprofile
54 self.nplots = nplots
54 self.nplots = nplots
55
55
56 ncolspan = 1
56 ncolspan = 1
57 colspan = 1
57 colspan = 1
58 if showprofile:
58 if showprofile:
59 ncolspan = 3
59 ncolspan = 3
60 colspan = 2
60 colspan = 2
61 self.__nsubplots = 2
61 self.__nsubplots = 2
62
62
63 self.createFigure(id = id,
63 self.createFigure(id = id,
64 wintitle = wintitle,
64 wintitle = wintitle,
65 widthplot = self.WIDTH + self.WIDTHPROF,
65 widthplot = self.WIDTH + self.WIDTHPROF,
66 heightplot = self.HEIGHT + self.HEIGHTPROF,
66 heightplot = self.HEIGHT + self.HEIGHTPROF,
67 show=show)
67 show=show)
68
68
69 nrow, ncol = self.getSubplots()
69 nrow, ncol = self.getSubplots()
70
70
71 counter = 0
71 counter = 0
72 for y in range(nrow):
72 for y in range(nrow):
73 for x in range(ncol):
73 for x in range(ncol):
74
74
75 if counter >= self.nplots:
75 if counter >= self.nplots:
76 break
76 break
77
77
78 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan, colspan, 1)
78 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan, colspan, 1)
79
79
80 if showprofile:
80 if showprofile:
81 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan+colspan, 1, 1)
81 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan+colspan, 1, 1)
82
82
83 counter += 1
83 counter += 1
84
84
85 def run(self, dataOut, id, wintitle="", channelList=None, showprofile=True,
85 def run(self, dataOut, id, wintitle="", channelList=None, showprofile=True,
86 xmin=None, xmax=None, ymin=None, ymax=None, zmin=None, zmax=None,
86 xmin=None, xmax=None, ymin=None, ymax=None, zmin=None, zmax=None,
87 save=False, figpath='./', figfile=None, show=True, ftp=False, wr_period=1,
87 save=False, figpath='./', figfile=None, show=True, ftp=False, wr_period=1,
88 server=None, folder=None, username=None, password=None,
88 server=None, folder=None, username=None, password=None,
89 ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0, realtime=False,
89 ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0, realtime=False,
90 xaxis="frequency", colormap='jet', normFactor=None):
90 xaxis="frequency", colormap='jet', normFactor=None):
91
91
92 """
92 """
93
93
94 Input:
94 Input:
95 dataOut :
95 dataOut :
96 id :
96 id :
97 wintitle :
97 wintitle :
98 channelList :
98 channelList :
99 showProfile :
99 showProfile :
100 xmin : None,
100 xmin : None,
101 xmax : None,
101 xmax : None,
102 ymin : None,
102 ymin : None,
103 ymax : None,
103 ymax : None,
104 zmin : None,
104 zmin : None,
105 zmax : None
105 zmax : None
106 """
106 """
107 if realtime:
107 if realtime:
108 if not(isRealtime(utcdatatime = dataOut.utctime)):
108 if not(isRealtime(utcdatatime = dataOut.utctime)):
109 print 'Skipping this plot function'
109 print 'Skipping this plot function'
110 return
110 return
111
111
112 if channelList == None:
112 if channelList == None:
113 channelIndexList = dataOut.channelIndexList
113 channelIndexList = dataOut.channelIndexList
114 else:
114 else:
115 channelIndexList = []
115 channelIndexList = []
116 for channel in channelList:
116 for channel in channelList:
117 if channel not in dataOut.channelList:
117 if channel not in dataOut.channelList:
118 raise ValueError, "Channel %d is not in dataOut.channelList" %channel
118 raise ValueError, "Channel %d is not in dataOut.channelList" %channel
119 channelIndexList.append(dataOut.channelList.index(channel))
119 channelIndexList.append(dataOut.channelList.index(channel))
120
120
121 if normFactor is None:
121 if normFactor is None:
122 factor = dataOut.normFactor
122 factor = dataOut.normFactor
123 else:
123 else:
124 factor = normFactor
124 factor = normFactor
125 if xaxis == "frequency":
125 if xaxis == "frequency":
126 x = dataOut.getFreqRange(1)/1000.
126 x = dataOut.getFreqRange(1)/1000.
127 xlabel = "Frequency (kHz)"
127 xlabel = "Frequency (kHz)"
128
128
129 elif xaxis == "time":
129 elif xaxis == "time":
130 x = dataOut.getAcfRange(1)
130 x = dataOut.getAcfRange(1)
131 xlabel = "Time (ms)"
131 xlabel = "Time (ms)"
132
132
133 else:
133 else:
134 x = dataOut.getVelRange(1)
134 x = dataOut.getVelRange(1)
135 xlabel = "Velocity (m/s)"
135 xlabel = "Velocity (m/s)"
136
136
137 ylabel = "Range (Km)"
137 ylabel = "Range (Km)"
138
138
139 y = dataOut.getHeiRange()
139 y = dataOut.getHeiRange()
140
140
141 z = dataOut.data_spc/factor
141 z = dataOut.data_spc/factor
142 z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
142 z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
143 zdB = 10*numpy.log10(z)
143 zdB = 10*numpy.log10(z)
144
144
145 #print "a000",dataOut.data_spc.dtype
145 #print "a000",dataOut.data_spc.dtype
146 avg = numpy.average(z, axis=1)
146 avg = numpy.average(z, axis=1)
147 avgdB = 10*numpy.log10(avg)
147 avgdB = 10*numpy.log10(avg)
148 #print "before plot"
148 #print "before plot"
149 noise = dataOut.getNoise()/factor
149 noise = dataOut.getNoise()/factor
150 noisedB = 10*numpy.log10(noise)
150 noisedB = 10*numpy.log10(noise)
151
151
152 thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[0])
152 thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[0])
153 title = wintitle + " Spectra"
153 title = wintitle + " Spectra"
154 if ((dataOut.azimuth!=None) and (dataOut.zenith!=None)):
154 if ((dataOut.azimuth!=None) and (dataOut.zenith!=None)):
155 title = title + '_' + 'azimuth,zenith=%2.2f,%2.2f'%(dataOut.azimuth, dataOut.zenith)
155 title = title + '_' + 'azimuth,zenith=%2.2f,%2.2f'%(dataOut.azimuth, dataOut.zenith)
156
156
157 if not self.isConfig:
157 if not self.isConfig:
158
158
159 nplots = len(channelIndexList)
159 nplots = len(channelIndexList)
160
160
161 self.setup(id=id,
161 self.setup(id=id,
162 nplots=nplots,
162 nplots=nplots,
163 wintitle=wintitle,
163 wintitle=wintitle,
164 showprofile=showprofile,
164 showprofile=showprofile,
165 show=show)
165 show=show)
166
166
167 if xmin == None: xmin = numpy.nanmin(x)
167 if xmin == None: xmin = numpy.nanmin(x)
168 if xmax == None: xmax = numpy.nanmax(x)
168 if xmax == None: xmax = numpy.nanmax(x)
169 if ymin == None: ymin = numpy.nanmin(y)
169 if ymin == None: ymin = numpy.nanmin(y)
170 if ymax == None: ymax = numpy.nanmax(y)
170 if ymax == None: ymax = numpy.nanmax(y)
171 if zmin == None: zmin = numpy.floor(numpy.nanmin(noisedB)) - 3
171 if zmin == None: zmin = numpy.floor(numpy.nanmin(noisedB)) - 3
172 if zmax == None: zmax = numpy.ceil(numpy.nanmax(avgdB)) + 3
172 if zmax == None: zmax = numpy.ceil(numpy.nanmax(avgdB)) + 3
173
173
174 self.FTP_WEI = ftp_wei
174 self.FTP_WEI = ftp_wei
175 self.EXP_CODE = exp_code
175 self.EXP_CODE = exp_code
176 self.SUB_EXP_CODE = sub_exp_code
176 self.SUB_EXP_CODE = sub_exp_code
177 self.PLOT_POS = plot_pos
177 self.PLOT_POS = plot_pos
178
178
179 self.isConfig = True
179 self.isConfig = True
180
180
181 self.setWinTitle(title)
181 self.setWinTitle(title)
182
182
183 for i in range(self.nplots):
183 for i in range(self.nplots):
184 index = channelIndexList[i]
184 index = channelIndexList[i]
185 str_datetime = '%s %s'%(thisDatetime.strftime("%Y/%m/%d"),thisDatetime.strftime("%H:%M:%S"))
185 str_datetime = '%s %s'%(thisDatetime.strftime("%Y/%m/%d"),thisDatetime.strftime("%H:%M:%S"))
186 title = "Channel %d: %4.2fdB: %s" %(dataOut.channelList[index], noisedB[index], str_datetime)
186 title = "Channel %d: %4.2fdB: %s" %(dataOut.channelList[index], noisedB[index], str_datetime)
187 if len(dataOut.beam.codeList) != 0:
187 if len(dataOut.beam.codeList) != 0:
188 title = "Ch%d:%4.2fdB,%2.2f,%2.2f:%s" %(dataOut.channelList[index], noisedB[index], dataOut.beam.azimuthList[index], dataOut.beam.zenithList[index], str_datetime)
188 title = "Ch%d:%4.2fdB,%2.2f,%2.2f:%s" %(dataOut.channelList[index], noisedB[index], dataOut.beam.azimuthList[index], dataOut.beam.zenithList[index], str_datetime)
189
189
190 axes = self.axesList[i*self.__nsubplots]
190 axes = self.axesList[i*self.__nsubplots]
191 axes.pcolor(x, y, zdB[index,:,:],
191 axes.pcolor(x, y, zdB[index,:,:],
192 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=zmin, zmax=zmax,
192 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=zmin, zmax=zmax,
193 xlabel=xlabel, ylabel=ylabel, title=title, colormap=colormap,
193 xlabel=xlabel, ylabel=ylabel, title=title, colormap=colormap,
194 ticksize=9, cblabel='')
194 ticksize=9, cblabel='')
195
195
196 if self.__showprofile:
196 if self.__showprofile:
197 axes = self.axesList[i*self.__nsubplots +1]
197 axes = self.axesList[i*self.__nsubplots +1]
198 axes.pline(avgdB[index,:], y,
198 axes.pline(avgdB[index,:], y,
199 xmin=zmin, xmax=zmax, ymin=ymin, ymax=ymax,
199 xmin=zmin, xmax=zmax, ymin=ymin, ymax=ymax,
200 xlabel='dB', ylabel='', title='',
200 xlabel='dB', ylabel='', title='',
201 ytick_visible=False,
201 ytick_visible=False,
202 grid='x')
202 grid='x')
203
203
204 noiseline = numpy.repeat(noisedB[index], len(y))
204 noiseline = numpy.repeat(noisedB[index], len(y))
205 axes.addpline(noiseline, y, idline=1, color="black", linestyle="dashed", lw=2)
205 axes.addpline(noiseline, y, idline=1, color="black", linestyle="dashed", lw=2)
206
206
207 self.draw()
207 self.draw()
208
208
209 if figfile == None:
209 if figfile == None:
210 str_datetime = thisDatetime.strftime("%Y%m%d_%H%M%S")
210 str_datetime = thisDatetime.strftime("%Y%m%d_%H%M%S")
211 name = str_datetime
211 name = str_datetime
212 if ((dataOut.azimuth!=None) and (dataOut.zenith!=None)):
212 if ((dataOut.azimuth!=None) and (dataOut.zenith!=None)):
213 name = name + '_az' + '_%2.2f'%(dataOut.azimuth) + '_zn' + '_%2.2f'%(dataOut.zenith)
213 name = name + '_az' + '_%2.2f'%(dataOut.azimuth) + '_zn' + '_%2.2f'%(dataOut.zenith)
214 figfile = self.getFilename(name)
214 figfile = self.getFilename(name)
215
215
216 self.save(figpath=figpath,
216 self.save(figpath=figpath,
217 figfile=figfile,
217 figfile=figfile,
218 save=save,
218 save=save,
219 ftp=ftp,
219 ftp=ftp,
220 wr_period=wr_period,
220 wr_period=wr_period,
221 thisDatetime=thisDatetime)
221 thisDatetime=thisDatetime)
222
222
223 class ACFPlot(Figure):
223 class ACFPlot(Figure):
224
224
225 isConfig = None
225 isConfig = None
226 __nsubplots = None
226 __nsubplots = None
227
227
228 WIDTHPROF = None
228 WIDTHPROF = None
229 HEIGHTPROF = None
229 HEIGHTPROF = None
230 PREFIX = 'acf'
230 PREFIX = 'acf'
231
231
232 def __init__(self, **kwargs):
232 def __init__(self, **kwargs):
233 Figure.__init__(self, **kwargs)
233 Figure.__init__(self, **kwargs)
234 self.isConfig = False
234 self.isConfig = False
235 self.__nsubplots = 1
235 self.__nsubplots = 1
236
236
237 self.PLOT_CODE = ACF_CODE
237 self.PLOT_CODE = ACF_CODE
238
238
239 self.WIDTH = 900
239 self.WIDTH = 900
240 self.HEIGHT = 700
240 self.HEIGHT = 700
241 self.counter_imagwr = 0
241 self.counter_imagwr = 0
242
242
243 def getSubplots(self):
243 def getSubplots(self):
244 ncol = 1
244 ncol = 1
245 nrow = 1
245 nrow = 1
246
246
247 return nrow, ncol
247 return nrow, ncol
248
248
249 def setup(self, id, nplots, wintitle, show):
249 def setup(self, id, nplots, wintitle, show):
250
250
251 self.nplots = nplots
251 self.nplots = nplots
252
252
253 ncolspan = 1
253 ncolspan = 1
254 colspan = 1
254 colspan = 1
255
255
256 self.createFigure(id = id,
256 self.createFigure(id = id,
257 wintitle = wintitle,
257 wintitle = wintitle,
258 widthplot = self.WIDTH,
258 widthplot = self.WIDTH,
259 heightplot = self.HEIGHT,
259 heightplot = self.HEIGHT,
260 show=show)
260 show=show)
261
261
262 nrow, ncol = self.getSubplots()
262 nrow, ncol = self.getSubplots()
263
263
264 counter = 0
264 counter = 0
265 for y in range(nrow):
265 for y in range(nrow):
266 for x in range(ncol):
266 for x in range(ncol):
267 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan, colspan, 1)
267 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan, colspan, 1)
268
268
269 def run(self, dataOut, id, wintitle="", channelList=None,channel=None,nSamples=None,
269 def run(self, dataOut, id, wintitle="", channelList=None,channel=None,nSamples=None,
270 nSampleList= None,resolutionFactor=None,
270 nSampleList= None,resolutionFactor=None,
271 xmin=None, xmax=None, ymin=None, ymax=None,
271 xmin=None, xmax=None, ymin=None, ymax=None,
272 save=False, figpath='./', figfile=None, show=True,
272 save=False, figpath='./', figfile=None, show=True,
273 ftp=False, wr_period=1, server=None,
273 ftp=False, wr_period=1, server=None,
274 folder=None, username=None, password=None,
274 folder=None, username=None, password=None,
275 xaxis="frequency"):
275 xaxis="frequency"):
276
276
277 channel0 = channel
277 channel0 = channel
278 nSamples = nSamples
278 nSamples = nSamples
279 resFactor = resolutionFactor
279 resFactor = resolutionFactor
280
280
281 if nSamples == None:
281 if nSamples == None:
282 nSamples = 20
282 nSamples = 20
283
283
284 if resFactor == None:
284 if resFactor == None:
285 resFactor = 5
285 resFactor = 5
286 #else:
286 #else:
287 # if nSamples not in dataOut.channelList:
287 # if nSamples not in dataOut.channelList:
288 # raise ValueError, "Channel %d is not in %s dataOut.channelList"%(channel0, dataOut.channelList)
288 # raise ValueError, "Channel %d is not in %s dataOut.channelList"%(channel0, dataOut.channelList)
289
289
290 if channel0 == None:
290 if channel0 == None:
291 channel0 = 0
291 channel0 = 0
292 else:
292 else:
293 if channel0 not in dataOut.channelList:
293 if channel0 not in dataOut.channelList:
294 raise ValueError, "Channel %d is not in %s dataOut.channelList"%(channel0, dataOut.channelList)
294 raise ValueError, "Channel %d is not in %s dataOut.channelList"%(channel0, dataOut.channelList)
295
295
296 if channelList == None:
296 if channelList == None:
297 channelIndexList = dataOut.channelIndexList
297 channelIndexList = dataOut.channelIndexList
298 channelList = dataOut.channelList
298 channelList = dataOut.channelList
299 else:
299 else:
300 channelIndexList = []
300 channelIndexList = []
301 for channel in channelList:
301 for channel in channelList:
302 if channel not in dataOut.channelList:
302 if channel not in dataOut.channelList:
303 raise ValueError, "Channel %d is not in dataOut.channelList"
303 raise ValueError, "Channel %d is not in dataOut.channelList"
304 channelIndexList.append(dataOut.channelList.index(channel))
304 channelIndexList.append(dataOut.channelList.index(channel))
305
305
306 #z = dataOut.data_spc/factor
306 #z = dataOut.data_spc/factor
307 factor = dataOut.normFactor
307 factor = dataOut.normFactor
308 y = dataOut.getHeiRange()
308 y = dataOut.getHeiRange()
309 deltaHeight = dataOut.heightList[1]-dataOut.heightList[0]
309 deltaHeight = dataOut.heightList[1]-dataOut.heightList[0]
310 z = dataOut.data_acf
310 z = dataOut.data_acf
311
311
312 #import matplotlib.pyplot as plt
312 #import matplotlib.pyplot as plt
313 #plt.plot(z[0,:,85])
313 #plt.plot(z[0,:,85])
314 #plt.show()
314 #plt.show()
315 #import time
315 #import time
316 #time.sleep(20)
316 #time.sleep(20)
317
317
318
318
319 #z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
319 #z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
320 shape = dataOut.data_acf.shape
320 shape = dataOut.data_acf.shape
321 print "shape_plot",shape
321 hei_index = numpy.arange(shape[2])
322 hei_index = numpy.arange(shape[2])
322 hei_plot = numpy.arange(nSamples)*resFactor
323 hei_plot = numpy.arange(nSamples)*resFactor
323
324
324 if nSampleList is not None:
325 if nSampleList is not None:
325 for nsample in nSampleList:
326 for nsample in nSampleList:
326 if nsample not in dataOut.heightList/deltaHeight:
327 if nsample not in dataOut.heightList/deltaHeight:
327 raise ValueError, "nsample %d is not in %s dataOut.heightList"%(nsample,dataOut.heightList)
328 raise ValueError, "nsample %d is not in %s dataOut.heightList"%(nsample,dataOut.heightList)
328
329
329 if nSampleList is not None:
330 if nSampleList is not None:
330 hei_plot = numpy.array(nSampleList)*resFactor
331 hei_plot = numpy.array(nSampleList)*resFactor
331
332
332 if hei_plot[-1] >= hei_index[-1]:
333 if hei_plot[-1] >= hei_index[-1]:
333 print ("La cantidad de puntos en altura es %d y la resolucion es %d Km"%(hei_plot.shape[0],deltaHeight*resFactor ))
334 print ("La cantidad de puntos en altura es %d y la resolucion es %f Km"%(hei_plot.shape[0],deltaHeight*resFactor ))
334 raise ValueError, "resFactor %d multiplicado por el valor de %d nSamples es mayor a %d cantidad total de puntos"%(resFactor,nSamples,hei_index[-1])
335 raise ValueError, "resFactor %d multiplicado por el valor de %d nSamples es mayor a %d cantidad total de puntos"%(resFactor,nSamples,hei_index[-1])
335
336
336 #escalamiento -1 a 1 a resolucion (factor de resolucion en altura)* deltaHeight
337 #escalamiento -1 a 1 a resolucion (factor de resolucion en altura)* deltaHeight
337 min = numpy.min(z[0,:,0])
338 min = numpy.min(z[0,:,0])
338 max =numpy.max(z[0,:,0])
339 max =numpy.max(z[0,:,0])
339 for i in range(shape[0]):
340 for i in range(shape[0]):
340 for j in range(shape[2]):
341 for j in range(shape[2]):
341 z[i,:,j]= (((z[i,:,j]-min)/(max-min))*deltaHeight*resFactor + j*deltaHeight)
342 z[i,:,j]= (((z[i,:,j]-min)/(max-min))*deltaHeight*resFactor + j*deltaHeight)
342
343
343 if xaxis == "frequency":
344 if xaxis == "frequency":
344 x = dataOut.getFreqRange()/1000.
345 x = dataOut.getFreqRange()/1000.
345 zdB = 10*numpy.log10(z[channel0,:,hei_plot])
346 zdB = 10*numpy.log10(z[channel0,:,hei_plot])
346 xlabel = "Frequency (kHz)"
347 xlabel = "Frequency (kHz)"
347 ylabel = "Power (dB)"
348 ylabel = "Power (dB)"
348
349
349 elif xaxis == "time":
350 elif xaxis == "time":
350 delta= dataOut.getAcfRange()[1]-dataOut.getAcfRange()[0]
351 delta= dataOut.getAcfRange()[1]-dataOut.getAcfRange()[0]
351 x = dataOut.getAcfRange()+delta/2.0
352 x = dataOut.getAcfRange()+delta/2.0
352 zdB = z[channel0,:,hei_plot]
353 zdB = z[channel0,:,hei_plot]
353 xlabel = "Time (ms)"
354 xlabel = "Time (ms)"
354 ylabel = "ACF"
355 ylabel = "ACF"
355
356
356 else:
357 else:
357 x = dataOut.getVelRange()
358 x = dataOut.getVelRange()
358 zdB = 10*numpy.log10(z[channel0,:,hei_plot])
359 zdB = 10*numpy.log10(z[channel0,:,hei_plot])
359 xlabel = "Velocity (m/s)"
360 xlabel = "Velocity (m/s)"
360 ylabel = "Power (dB)"
361 ylabel = "Power (dB)"
361
362
362 thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[0])
363 thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[0])
363 title = wintitle + " ACF Plot Ch %s %s" %(channel0,thisDatetime.strftime("%d-%b-%Y"))
364 title = wintitle + " ACF Plot Ch %s %s" %(channel0,thisDatetime.strftime("%d-%b-%Y"))
364
365
365 if not self.isConfig:
366 if not self.isConfig:
366
367
367 nplots = 1
368 nplots = 1
368
369
369 self.setup(id=id,
370 self.setup(id=id,
370 nplots=nplots,
371 nplots=nplots,
371 wintitle=wintitle,
372 wintitle=wintitle,
372 show=show)
373 show=show)
373
374
374 if xmin == None: xmin = numpy.nanmin(x)*0.9
375 if xmin == None: xmin = numpy.nanmin(x)*0.9
375 if xmax == None: xmax = numpy.nanmax(x)*1.1
376 if xmax == None: xmax = numpy.nanmax(x)*1.1
376 if ymin == None: ymin = numpy.nanmin(zdB)
377 if ymin == None: ymin = numpy.nanmin(zdB)
377 if ymax == None: ymax = numpy.nanmax(zdB)
378 if ymax == None: ymax = numpy.nanmax(zdB)
378
379
379 print ("El parametro resFactor es %d y la resolucion en altura es %d"%(resFactor,deltaHeight ))
380 print ("El parametro resFactor es %d y la resolucion en altura es %f"%(resFactor,deltaHeight ))
380 print ("La cantidad de puntos en altura es %d y la nueva resolucion es %d Km"%(hei_plot.shape[0],deltaHeight*resFactor ))
381 print ("La cantidad de puntos en altura es %d y la nueva resolucion es %f Km"%(hei_plot.shape[0],deltaHeight*resFactor ))
381 print ("La altura maxima es %d Km"%(hei_plot[-1]*deltaHeight ))
382 print ("La altura maxima es %d Km"%(hei_plot[-1]*deltaHeight ))
382
383
383 self.isConfig = True
384 self.isConfig = True
384
385
385 self.setWinTitle(title)
386 self.setWinTitle(title)
386
387
387 title = "ACF Plot: %s" %(thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
388 title = "ACF Plot: %s" %(thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
388 axes = self.axesList[0]
389 axes = self.axesList[0]
389
390
390 legendlabels = ["Range = %dKm" %y[i] for i in hei_plot]
391 legendlabels = ["Range = %dKm" %y[i] for i in hei_plot]
391
392
392 axes.pmultilineyaxis( x, zdB,
393 axes.pmultilineyaxis( x, zdB,
393 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax,
394 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax,
394 xlabel=xlabel, ylabel=ylabel, title=title, legendlabels=legendlabels,
395 xlabel=xlabel, ylabel=ylabel, title=title, legendlabels=legendlabels,
395 ytick_visible=True, nxticks=5,
396 ytick_visible=True, nxticks=5,
396 grid='x')
397 grid='x')
397
398
398 self.draw()
399 self.draw()
399
400
400 if figfile == None:
401 if figfile == None:
401 str_datetime = thisDatetime.strftime("%Y%m%d_%H%M%S")
402 str_datetime = thisDatetime.strftime("%Y%m%d_%H%M%S")
402 name = str_datetime
403 name = str_datetime
403 if ((dataOut.azimuth!=None) and (dataOut.zenith!=None)):
404 if ((dataOut.azimuth!=None) and (dataOut.zenith!=None)):
404 name = name + '_az' + '_%2.2f'%(dataOut.azimuth) + '_zn' + '_%2.2f'%(dataOut.zenith)
405 name = name + '_az' + '_%2.2f'%(dataOut.azimuth) + '_zn' + '_%2.2f'%(dataOut.zenith)
405 figfile = self.getFilename(name)
406 figfile = self.getFilename(name)
406
407
407 self.save(figpath=figpath,
408 self.save(figpath=figpath,
408 figfile=figfile,
409 figfile=figfile,
409 save=save,
410 save=save,
410 ftp=ftp,
411 ftp=ftp,
411 wr_period=wr_period,
412 wr_period=wr_period,
412 thisDatetime=thisDatetime)
413 thisDatetime=thisDatetime)
413
414
414
415
415
416
416 class CrossSpectraPlot(Figure):
417 class CrossSpectraPlot(Figure):
417
418
418 isConfig = None
419 isConfig = None
419 __nsubplots = None
420 __nsubplots = None
420
421
421 WIDTH = None
422 WIDTH = None
422 HEIGHT = None
423 HEIGHT = None
423 WIDTHPROF = None
424 WIDTHPROF = None
424 HEIGHTPROF = None
425 HEIGHTPROF = None
425 PREFIX = 'cspc'
426 PREFIX = 'cspc'
426
427
427 def __init__(self, **kwargs):
428 def __init__(self, **kwargs):
428 Figure.__init__(self, **kwargs)
429 Figure.__init__(self, **kwargs)
429 self.isConfig = False
430 self.isConfig = False
430 self.__nsubplots = 4
431 self.__nsubplots = 4
431 self.counter_imagwr = 0
432 self.counter_imagwr = 0
432 self.WIDTH = 250
433 self.WIDTH = 250
433 self.HEIGHT = 250
434 self.HEIGHT = 250
434 self.WIDTHPROF = 0
435 self.WIDTHPROF = 0
435 self.HEIGHTPROF = 0
436 self.HEIGHTPROF = 0
436
437
437 self.PLOT_CODE = CROSS_CODE
438 self.PLOT_CODE = CROSS_CODE
438 self.FTP_WEI = None
439 self.FTP_WEI = None
439 self.EXP_CODE = None
440 self.EXP_CODE = None
440 self.SUB_EXP_CODE = None
441 self.SUB_EXP_CODE = None
441 self.PLOT_POS = None
442 self.PLOT_POS = None
442
443
443 def getSubplots(self):
444 def getSubplots(self):
444
445
445 ncol = 4
446 ncol = 4
446 nrow = self.nplots
447 nrow = self.nplots
447
448
448 return nrow, ncol
449 return nrow, ncol
449
450
450 def setup(self, id, nplots, wintitle, showprofile=True, show=True):
451 def setup(self, id, nplots, wintitle, showprofile=True, show=True):
451
452
452 self.__showprofile = showprofile
453 self.__showprofile = showprofile
453 self.nplots = nplots
454 self.nplots = nplots
454
455
455 ncolspan = 1
456 ncolspan = 1
456 colspan = 1
457 colspan = 1
457
458
458 self.createFigure(id = id,
459 self.createFigure(id = id,
459 wintitle = wintitle,
460 wintitle = wintitle,
460 widthplot = self.WIDTH + self.WIDTHPROF,
461 widthplot = self.WIDTH + self.WIDTHPROF,
461 heightplot = self.HEIGHT + self.HEIGHTPROF,
462 heightplot = self.HEIGHT + self.HEIGHTPROF,
462 show=True)
463 show=True)
463
464
464 nrow, ncol = self.getSubplots()
465 nrow, ncol = self.getSubplots()
465
466
466 counter = 0
467 counter = 0
467 for y in range(nrow):
468 for y in range(nrow):
468 for x in range(ncol):
469 for x in range(ncol):
469 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan, colspan, 1)
470 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan, colspan, 1)
470
471
471 counter += 1
472 counter += 1
472
473
473 def run(self, dataOut, id, wintitle="", pairsList=None,
474 def run(self, dataOut, id, wintitle="", pairsList=None,
474 xmin=None, xmax=None, ymin=None, ymax=None, zmin=None, zmax=None,
475 xmin=None, xmax=None, ymin=None, ymax=None, zmin=None, zmax=None,
475 coh_min=None, coh_max=None, phase_min=None, phase_max=None,
476 coh_min=None, coh_max=None, phase_min=None, phase_max=None,
476 save=False, figpath='./', figfile=None, ftp=False, wr_period=1,
477 save=False, figpath='./', figfile=None, ftp=False, wr_period=1,
477 power_cmap='jet', coherence_cmap='jet', phase_cmap='RdBu_r', show=True,
478 power_cmap='jet', coherence_cmap='jet', phase_cmap='RdBu_r', show=True,
478 server=None, folder=None, username=None, password=None,
479 server=None, folder=None, username=None, password=None,
479 ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0, normFactor=None,
480 ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0, normFactor=None,
480 xaxis='frequency'):
481 xaxis='frequency'):
481
482
482 """
483 """
483
484
484 Input:
485 Input:
485 dataOut :
486 dataOut :
486 id :
487 id :
487 wintitle :
488 wintitle :
488 channelList :
489 channelList :
489 showProfile :
490 showProfile :
490 xmin : None,
491 xmin : None,
491 xmax : None,
492 xmax : None,
492 ymin : None,
493 ymin : None,
493 ymax : None,
494 ymax : None,
494 zmin : None,
495 zmin : None,
495 zmax : None
496 zmax : None
496 """
497 """
497
498
498 if pairsList == None:
499 if pairsList == None:
499 pairsIndexList = dataOut.pairsIndexList
500 pairsIndexList = dataOut.pairsIndexList
500 else:
501 else:
501 pairsIndexList = []
502 pairsIndexList = []
502 for pair in pairsList:
503 for pair in pairsList:
503 if pair not in dataOut.pairsList:
504 if pair not in dataOut.pairsList:
504 raise ValueError, "Pair %s is not in dataOut.pairsList" %str(pair)
505 raise ValueError, "Pair %s is not in dataOut.pairsList" %str(pair)
505 pairsIndexList.append(dataOut.pairsList.index(pair))
506 pairsIndexList.append(dataOut.pairsList.index(pair))
506
507
507 if not pairsIndexList:
508 if not pairsIndexList:
508 return
509 return
509
510
510 if len(pairsIndexList) > 4:
511 if len(pairsIndexList) > 4:
511 pairsIndexList = pairsIndexList[0:4]
512 pairsIndexList = pairsIndexList[0:4]
512
513
513 if normFactor is None:
514 if normFactor is None:
514 factor = dataOut.normFactor
515 factor = dataOut.normFactor
515 else:
516 else:
516 factor = normFactor
517 factor = normFactor
517 x = dataOut.getVelRange(1)
518 x = dataOut.getVelRange(1)
518 y = dataOut.getHeiRange()
519 y = dataOut.getHeiRange()
519 z = dataOut.data_spc[:,:,:]/factor
520 z = dataOut.data_spc[:,:,:]/factor
520 z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
521 z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
521
522
522 noise = dataOut.noise/factor
523 noise = dataOut.noise/factor
523
524
524 zdB = 10*numpy.log10(z)
525 zdB = 10*numpy.log10(z)
525 noisedB = 10*numpy.log10(noise)
526 noisedB = 10*numpy.log10(noise)
526
527
527 if coh_min == None:
528 if coh_min == None:
528 coh_min = 0.0
529 coh_min = 0.0
529 if coh_max == None:
530 if coh_max == None:
530 coh_max = 1.0
531 coh_max = 1.0
531
532
532 if phase_min == None:
533 if phase_min == None:
533 phase_min = -180
534 phase_min = -180
534 if phase_max == None:
535 if phase_max == None:
535 phase_max = 180
536 phase_max = 180
536
537
537 #thisDatetime = dataOut.datatime
538 #thisDatetime = dataOut.datatime
538 thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[0])
539 thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[0])
539 title = wintitle + " Cross-Spectra: %s" %(thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
540 title = wintitle + " Cross-Spectra: %s" %(thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
540 # xlabel = "Velocity (m/s)"
541 # xlabel = "Velocity (m/s)"
541 ylabel = "Range (Km)"
542 ylabel = "Range (Km)"
542
543
543 if xaxis == "frequency":
544 if xaxis == "frequency":
544 x = dataOut.getFreqRange(1)/1000.
545 x = dataOut.getFreqRange(1)/1000.
545 xlabel = "Frequency (kHz)"
546 xlabel = "Frequency (kHz)"
546
547
547 elif xaxis == "time":
548 elif xaxis == "time":
548 x = dataOut.getAcfRange(1)
549 x = dataOut.getAcfRange(1)
549 xlabel = "Time (ms)"
550 xlabel = "Time (ms)"
550
551
551 else:
552 else:
552 x = dataOut.getVelRange(1)
553 x = dataOut.getVelRange(1)
553 xlabel = "Velocity (m/s)"
554 xlabel = "Velocity (m/s)"
554
555
555 if not self.isConfig:
556 if not self.isConfig:
556
557
557 nplots = len(pairsIndexList)
558 nplots = len(pairsIndexList)
558
559
559 self.setup(id=id,
560 self.setup(id=id,
560 nplots=nplots,
561 nplots=nplots,
561 wintitle=wintitle,
562 wintitle=wintitle,
562 showprofile=False,
563 showprofile=False,
563 show=show)
564 show=show)
564
565
565 avg = numpy.abs(numpy.average(z, axis=1))
566 avg = numpy.abs(numpy.average(z, axis=1))
566 avgdB = 10*numpy.log10(avg)
567 avgdB = 10*numpy.log10(avg)
567
568
568 if xmin == None: xmin = numpy.nanmin(x)
569 if xmin == None: xmin = numpy.nanmin(x)
569 if xmax == None: xmax = numpy.nanmax(x)
570 if xmax == None: xmax = numpy.nanmax(x)
570 if ymin == None: ymin = numpy.nanmin(y)
571 if ymin == None: ymin = numpy.nanmin(y)
571 if ymax == None: ymax = numpy.nanmax(y)
572 if ymax == None: ymax = numpy.nanmax(y)
572 if zmin == None: zmin = numpy.floor(numpy.nanmin(noisedB)) - 3
573 if zmin == None: zmin = numpy.floor(numpy.nanmin(noisedB)) - 3
573 if zmax == None: zmax = numpy.ceil(numpy.nanmax(avgdB)) + 3
574 if zmax == None: zmax = numpy.ceil(numpy.nanmax(avgdB)) + 3
574
575
575 self.FTP_WEI = ftp_wei
576 self.FTP_WEI = ftp_wei
576 self.EXP_CODE = exp_code
577 self.EXP_CODE = exp_code
577 self.SUB_EXP_CODE = sub_exp_code
578 self.SUB_EXP_CODE = sub_exp_code
578 self.PLOT_POS = plot_pos
579 self.PLOT_POS = plot_pos
579
580
580 self.isConfig = True
581 self.isConfig = True
581
582
582 self.setWinTitle(title)
583 self.setWinTitle(title)
583
584
584 for i in range(self.nplots):
585 for i in range(self.nplots):
585 pair = dataOut.pairsList[pairsIndexList[i]]
586 pair = dataOut.pairsList[pairsIndexList[i]]
586
587
587 chan_index0 = dataOut.channelList.index(pair[0])
588 chan_index0 = dataOut.channelList.index(pair[0])
588 chan_index1 = dataOut.channelList.index(pair[1])
589 chan_index1 = dataOut.channelList.index(pair[1])
589
590
590 str_datetime = '%s %s'%(thisDatetime.strftime("%Y/%m/%d"),thisDatetime.strftime("%H:%M:%S"))
591 str_datetime = '%s %s'%(thisDatetime.strftime("%Y/%m/%d"),thisDatetime.strftime("%H:%M:%S"))
591 title = "Ch%d: %4.2fdB: %s" %(pair[0], noisedB[chan_index0], str_datetime)
592 title = "Ch%d: %4.2fdB: %s" %(pair[0], noisedB[chan_index0], str_datetime)
592 zdB = 10.*numpy.log10(dataOut.data_spc[chan_index0,:,:]/factor)
593 zdB = 10.*numpy.log10(dataOut.data_spc[chan_index0,:,:]/factor)
593 axes0 = self.axesList[i*self.__nsubplots]
594 axes0 = self.axesList[i*self.__nsubplots]
594 axes0.pcolor(x, y, zdB,
595 axes0.pcolor(x, y, zdB,
595 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=zmin, zmax=zmax,
596 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=zmin, zmax=zmax,
596 xlabel=xlabel, ylabel=ylabel, title=title,
597 xlabel=xlabel, ylabel=ylabel, title=title,
597 ticksize=9, colormap=power_cmap, cblabel='')
598 ticksize=9, colormap=power_cmap, cblabel='')
598
599
599 title = "Ch%d: %4.2fdB: %s" %(pair[1], noisedB[chan_index1], str_datetime)
600 title = "Ch%d: %4.2fdB: %s" %(pair[1], noisedB[chan_index1], str_datetime)
600 zdB = 10.*numpy.log10(dataOut.data_spc[chan_index1,:,:]/factor)
601 zdB = 10.*numpy.log10(dataOut.data_spc[chan_index1,:,:]/factor)
601 axes0 = self.axesList[i*self.__nsubplots+1]
602 axes0 = self.axesList[i*self.__nsubplots+1]
602 axes0.pcolor(x, y, zdB,
603 axes0.pcolor(x, y, zdB,
603 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=zmin, zmax=zmax,
604 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=zmin, zmax=zmax,
604 xlabel=xlabel, ylabel=ylabel, title=title,
605 xlabel=xlabel, ylabel=ylabel, title=title,
605 ticksize=9, colormap=power_cmap, cblabel='')
606 ticksize=9, colormap=power_cmap, cblabel='')
606
607
607 coherenceComplex = dataOut.data_cspc[pairsIndexList[i],:,:]/numpy.sqrt(dataOut.data_spc[chan_index0,:,:]*dataOut.data_spc[chan_index1,:,:])
608 coherenceComplex = dataOut.data_cspc[pairsIndexList[i],:,:]/numpy.sqrt(dataOut.data_spc[chan_index0,:,:]*dataOut.data_spc[chan_index1,:,:])
608 coherence = numpy.abs(coherenceComplex)
609 coherence = numpy.abs(coherenceComplex)
609 # phase = numpy.arctan(-1*coherenceComplex.imag/coherenceComplex.real)*180/numpy.pi
610 # phase = numpy.arctan(-1*coherenceComplex.imag/coherenceComplex.real)*180/numpy.pi
610 phase = numpy.arctan2(coherenceComplex.imag, coherenceComplex.real)*180/numpy.pi
611 phase = numpy.arctan2(coherenceComplex.imag, coherenceComplex.real)*180/numpy.pi
611
612
612 title = "Coherence Ch%d * Ch%d" %(pair[0], pair[1])
613 title = "Coherence Ch%d * Ch%d" %(pair[0], pair[1])
613 axes0 = self.axesList[i*self.__nsubplots+2]
614 axes0 = self.axesList[i*self.__nsubplots+2]
614 axes0.pcolor(x, y, coherence,
615 axes0.pcolor(x, y, coherence,
615 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=coh_min, zmax=coh_max,
616 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=coh_min, zmax=coh_max,
616 xlabel=xlabel, ylabel=ylabel, title=title,
617 xlabel=xlabel, ylabel=ylabel, title=title,
617 ticksize=9, colormap=coherence_cmap, cblabel='')
618 ticksize=9, colormap=coherence_cmap, cblabel='')
618
619
619 title = "Phase Ch%d * Ch%d" %(pair[0], pair[1])
620 title = "Phase Ch%d * Ch%d" %(pair[0], pair[1])
620 axes0 = self.axesList[i*self.__nsubplots+3]
621 axes0 = self.axesList[i*self.__nsubplots+3]
621 axes0.pcolor(x, y, phase,
622 axes0.pcolor(x, y, phase,
622 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=phase_min, zmax=phase_max,
623 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=phase_min, zmax=phase_max,
623 xlabel=xlabel, ylabel=ylabel, title=title,
624 xlabel=xlabel, ylabel=ylabel, title=title,
624 ticksize=9, colormap=phase_cmap, cblabel='')
625 ticksize=9, colormap=phase_cmap, cblabel='')
625
626
626
627
627
628
628 self.draw()
629 self.draw()
629
630
630 self.save(figpath=figpath,
631 self.save(figpath=figpath,
631 figfile=figfile,
632 figfile=figfile,
632 save=save,
633 save=save,
633 ftp=ftp,
634 ftp=ftp,
634 wr_period=wr_period,
635 wr_period=wr_period,
635 thisDatetime=thisDatetime)
636 thisDatetime=thisDatetime)
636
637
637
638
638 class RTIPlot(Figure):
639 class RTIPlot(Figure):
639
640
640 __isConfig = None
641 __isConfig = None
641 __nsubplots = None
642 __nsubplots = None
642
643
643 WIDTHPROF = None
644 WIDTHPROF = None
644 HEIGHTPROF = None
645 HEIGHTPROF = None
645 PREFIX = 'rti'
646 PREFIX = 'rti'
646
647
647 def __init__(self, **kwargs):
648 def __init__(self, **kwargs):
648
649
649 Figure.__init__(self, **kwargs)
650 Figure.__init__(self, **kwargs)
650 self.timerange = None
651 self.timerange = None
651 self.isConfig = False
652 self.isConfig = False
652 self.__nsubplots = 1
653 self.__nsubplots = 1
653
654
654 self.WIDTH = 800
655 self.WIDTH = 800
655 self.HEIGHT = 180
656 self.HEIGHT = 180
656 self.WIDTHPROF = 120
657 self.WIDTHPROF = 120
657 self.HEIGHTPROF = 0
658 self.HEIGHTPROF = 0
658 self.counter_imagwr = 0
659 self.counter_imagwr = 0
659
660
660 self.PLOT_CODE = RTI_CODE
661 self.PLOT_CODE = RTI_CODE
661
662
662 self.FTP_WEI = None
663 self.FTP_WEI = None
663 self.EXP_CODE = None
664 self.EXP_CODE = None
664 self.SUB_EXP_CODE = None
665 self.SUB_EXP_CODE = None
665 self.PLOT_POS = None
666 self.PLOT_POS = None
666 self.tmin = None
667 self.tmin = None
667 self.tmax = None
668 self.tmax = None
668
669
669 self.xmin = None
670 self.xmin = None
670 self.xmax = None
671 self.xmax = None
671
672
672 self.figfile = None
673 self.figfile = None
673
674
674 def getSubplots(self):
675 def getSubplots(self):
675
676
676 ncol = 1
677 ncol = 1
677 nrow = self.nplots
678 nrow = self.nplots
678
679
679 return nrow, ncol
680 return nrow, ncol
680
681
681 def setup(self, id, nplots, wintitle, showprofile=True, show=True):
682 def setup(self, id, nplots, wintitle, showprofile=True, show=True):
682
683
683 self.__showprofile = showprofile
684 self.__showprofile = showprofile
684 self.nplots = nplots
685 self.nplots = nplots
685
686
686 ncolspan = 1
687 ncolspan = 1
687 colspan = 1
688 colspan = 1
688 if showprofile:
689 if showprofile:
689 ncolspan = 7
690 ncolspan = 7
690 colspan = 6
691 colspan = 6
691 self.__nsubplots = 2
692 self.__nsubplots = 2
692
693
693 self.createFigure(id = id,
694 self.createFigure(id = id,
694 wintitle = wintitle,
695 wintitle = wintitle,
695 widthplot = self.WIDTH + self.WIDTHPROF,
696 widthplot = self.WIDTH + self.WIDTHPROF,
696 heightplot = self.HEIGHT + self.HEIGHTPROF,
697 heightplot = self.HEIGHT + self.HEIGHTPROF,
697 show=show)
698 show=show)
698
699
699 nrow, ncol = self.getSubplots()
700 nrow, ncol = self.getSubplots()
700
701
701 counter = 0
702 counter = 0
702 for y in range(nrow):
703 for y in range(nrow):
703 for x in range(ncol):
704 for x in range(ncol):
704
705
705 if counter >= self.nplots:
706 if counter >= self.nplots:
706 break
707 break
707
708
708 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan, colspan, 1)
709 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan, colspan, 1)
709
710
710 if showprofile:
711 if showprofile:
711 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan+colspan, 1, 1)
712 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan+colspan, 1, 1)
712
713
713 counter += 1
714 counter += 1
714
715
715 def run(self, dataOut, id, wintitle="", channelList=None, showprofile='True',
716 def run(self, dataOut, id, wintitle="", channelList=None, showprofile='True',
716 xmin=None, xmax=None, ymin=None, ymax=None, zmin=None, zmax=None,
717 xmin=None, xmax=None, ymin=None, ymax=None, zmin=None, zmax=None,
717 timerange=None, colormap='jet',
718 timerange=None, colormap='jet',
718 save=False, figpath='./', lastone=0,figfile=None, ftp=False, wr_period=1, show=True,
719 save=False, figpath='./', lastone=0,figfile=None, ftp=False, wr_period=1, show=True,
719 server=None, folder=None, username=None, password=None,
720 server=None, folder=None, username=None, password=None,
720 ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0, normFactor=None, HEIGHT=None):
721 ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0, normFactor=None, HEIGHT=None):
721
722
722 """
723 """
723
724
724 Input:
725 Input:
725 dataOut :
726 dataOut :
726 id :
727 id :
727 wintitle :
728 wintitle :
728 channelList :
729 channelList :
729 showProfile :
730 showProfile :
730 xmin : None,
731 xmin : None,
731 xmax : None,
732 xmax : None,
732 ymin : None,
733 ymin : None,
733 ymax : None,
734 ymax : None,
734 zmin : None,
735 zmin : None,
735 zmax : None
736 zmax : None
736 """
737 """
737
738
738 #colormap = kwargs.get('colormap', 'jet')
739 #colormap = kwargs.get('colormap', 'jet')
739 if HEIGHT is not None:
740 if HEIGHT is not None:
740 self.HEIGHT = HEIGHT
741 self.HEIGHT = HEIGHT
741
742
742 if not isTimeInHourRange(dataOut.datatime, xmin, xmax):
743 if not isTimeInHourRange(dataOut.datatime, xmin, xmax):
743 return
744 return
744
745
745 if channelList == None:
746 if channelList == None:
746 channelIndexList = dataOut.channelIndexList
747 channelIndexList = dataOut.channelIndexList
747 else:
748 else:
748 channelIndexList = []
749 channelIndexList = []
749 for channel in channelList:
750 for channel in channelList:
750 if channel not in dataOut.channelList:
751 if channel not in dataOut.channelList:
751 raise ValueError, "Channel %d is not in dataOut.channelList"
752 raise ValueError, "Channel %d is not in dataOut.channelList"
752 channelIndexList.append(dataOut.channelList.index(channel))
753 channelIndexList.append(dataOut.channelList.index(channel))
753
754
754 if normFactor is None:
755 if normFactor is None:
755 factor = dataOut.normFactor
756 factor = dataOut.normFactor
756 else:
757 else:
757 factor = normFactor
758 factor = normFactor
758
759
759 # factor = dataOut.normFactor
760 # factor = dataOut.normFactor
760 x = dataOut.getTimeRange()
761 x = dataOut.getTimeRange()
761 y = dataOut.getHeiRange()
762 y = dataOut.getHeiRange()
762
763
763 z = dataOut.data_spc/factor
764 z = dataOut.data_spc/factor
764 z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
765 z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
765 avg = numpy.average(z, axis=1)
766 avg = numpy.average(z, axis=1)
766 avgdB = 10.*numpy.log10(avg)
767 avgdB = 10.*numpy.log10(avg)
767 # avgdB = dataOut.getPower()
768 # avgdB = dataOut.getPower()
768
769
769
770
770 thisDatetime = dataOut.datatime
771 thisDatetime = dataOut.datatime
771 # thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[0])
772 # thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[0])
772 title = wintitle + " RTI" #: %s" %(thisDatetime.strftime("%d-%b-%Y"))
773 title = wintitle + " RTI" #: %s" %(thisDatetime.strftime("%d-%b-%Y"))
773 xlabel = ""
774 xlabel = ""
774 ylabel = "Range (Km)"
775 ylabel = "Range (Km)"
775
776
776 update_figfile = False
777 update_figfile = False
777
778
778 if dataOut.ltctime >= self.xmax:
779 if dataOut.ltctime >= self.xmax:
779 self.counter_imagwr = wr_period
780 self.counter_imagwr = wr_period
780 self.isConfig = False
781 self.isConfig = False
781 update_figfile = True
782 update_figfile = True
782
783
783 if not self.isConfig:
784 if not self.isConfig:
784
785
785 nplots = len(channelIndexList)
786 nplots = len(channelIndexList)
786
787
787 self.setup(id=id,
788 self.setup(id=id,
788 nplots=nplots,
789 nplots=nplots,
789 wintitle=wintitle,
790 wintitle=wintitle,
790 showprofile=showprofile,
791 showprofile=showprofile,
791 show=show)
792 show=show)
792
793
793 if timerange != None:
794 if timerange != None:
794 self.timerange = timerange
795 self.timerange = timerange
795
796
796 self.xmin, self.xmax = self.getTimeLim(x, xmin, xmax, timerange)
797 self.xmin, self.xmax = self.getTimeLim(x, xmin, xmax, timerange)
797
798
798 noise = dataOut.noise/factor
799 noise = dataOut.noise/factor
799 noisedB = 10*numpy.log10(noise)
800 noisedB = 10*numpy.log10(noise)
800
801
801 if ymin == None: ymin = numpy.nanmin(y)
802 if ymin == None: ymin = numpy.nanmin(y)
802 if ymax == None: ymax = numpy.nanmax(y)
803 if ymax == None: ymax = numpy.nanmax(y)
803 if zmin == None: zmin = numpy.floor(numpy.nanmin(noisedB)) - 3
804 if zmin == None: zmin = numpy.floor(numpy.nanmin(noisedB)) - 3
804 if zmax == None: zmax = numpy.ceil(numpy.nanmax(avgdB)) + 3
805 if zmax == None: zmax = numpy.ceil(numpy.nanmax(avgdB)) + 3
805
806
806 self.FTP_WEI = ftp_wei
807 self.FTP_WEI = ftp_wei
807 self.EXP_CODE = exp_code
808 self.EXP_CODE = exp_code
808 self.SUB_EXP_CODE = sub_exp_code
809 self.SUB_EXP_CODE = sub_exp_code
809 self.PLOT_POS = plot_pos
810 self.PLOT_POS = plot_pos
810
811
811 self.name = thisDatetime.strftime("%Y%m%d_%H%M%S")
812 self.name = thisDatetime.strftime("%Y%m%d_%H%M%S")
812 self.isConfig = True
813 self.isConfig = True
813 self.figfile = figfile
814 self.figfile = figfile
814 update_figfile = True
815 update_figfile = True
815
816
816 self.setWinTitle(title)
817 self.setWinTitle(title)
817
818
818 for i in range(self.nplots):
819 for i in range(self.nplots):
819 index = channelIndexList[i]
820 index = channelIndexList[i]
820 title = "Channel %d: %s" %(dataOut.channelList[index], thisDatetime.strftime("%Y/%m/%d %H:%M:%S"))
821 title = "Channel %d: %s" %(dataOut.channelList[index], thisDatetime.strftime("%Y/%m/%d %H:%M:%S"))
821 if ((dataOut.azimuth!=None) and (dataOut.zenith!=None)):
822 if ((dataOut.azimuth!=None) and (dataOut.zenith!=None)):
822 title = title + '_' + 'azimuth,zenith=%2.2f,%2.2f'%(dataOut.azimuth, dataOut.zenith)
823 title = title + '_' + 'azimuth,zenith=%2.2f,%2.2f'%(dataOut.azimuth, dataOut.zenith)
823 axes = self.axesList[i*self.__nsubplots]
824 axes = self.axesList[i*self.__nsubplots]
824 zdB = avgdB[index].reshape((1,-1))
825 zdB = avgdB[index].reshape((1,-1))
825 axes.pcolorbuffer(x, y, zdB,
826 axes.pcolorbuffer(x, y, zdB,
826 xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax, zmin=zmin, zmax=zmax,
827 xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax, zmin=zmin, zmax=zmax,
827 xlabel=xlabel, ylabel=ylabel, title=title, rti=True, XAxisAsTime=True,
828 xlabel=xlabel, ylabel=ylabel, title=title, rti=True, XAxisAsTime=True,
828 ticksize=9, cblabel='', cbsize="1%", colormap=colormap)
829 ticksize=9, cblabel='', cbsize="1%", colormap=colormap)
829
830
830 if self.__showprofile:
831 if self.__showprofile:
831 axes = self.axesList[i*self.__nsubplots +1]
832 axes = self.axesList[i*self.__nsubplots +1]
832 axes.pline(avgdB[index], y,
833 axes.pline(avgdB[index], y,
833 xmin=zmin, xmax=zmax, ymin=ymin, ymax=ymax,
834 xmin=zmin, xmax=zmax, ymin=ymin, ymax=ymax,
834 xlabel='dB', ylabel='', title='',
835 xlabel='dB', ylabel='', title='',
835 ytick_visible=False,
836 ytick_visible=False,
836 grid='x')
837 grid='x')
837
838
838 self.draw()
839 self.draw()
839
840
840 self.save(figpath=figpath,
841 self.save(figpath=figpath,
841 figfile=figfile,
842 figfile=figfile,
842 save=save,
843 save=save,
843 ftp=ftp,
844 ftp=ftp,
844 wr_period=wr_period,
845 wr_period=wr_period,
845 thisDatetime=thisDatetime,
846 thisDatetime=thisDatetime,
846 update_figfile=update_figfile)
847 update_figfile=update_figfile)
847
848
848 class CoherenceMap(Figure):
849 class CoherenceMap(Figure):
849 isConfig = None
850 isConfig = None
850 __nsubplots = None
851 __nsubplots = None
851
852
852 WIDTHPROF = None
853 WIDTHPROF = None
853 HEIGHTPROF = None
854 HEIGHTPROF = None
854 PREFIX = 'cmap'
855 PREFIX = 'cmap'
855
856
856 def __init__(self, **kwargs):
857 def __init__(self, **kwargs):
857 Figure.__init__(self, **kwargs)
858 Figure.__init__(self, **kwargs)
858 self.timerange = 2*60*60
859 self.timerange = 2*60*60
859 self.isConfig = False
860 self.isConfig = False
860 self.__nsubplots = 1
861 self.__nsubplots = 1
861
862
862 self.WIDTH = 800
863 self.WIDTH = 800
863 self.HEIGHT = 180
864 self.HEIGHT = 180
864 self.WIDTHPROF = 120
865 self.WIDTHPROF = 120
865 self.HEIGHTPROF = 0
866 self.HEIGHTPROF = 0
866 self.counter_imagwr = 0
867 self.counter_imagwr = 0
867
868
868 self.PLOT_CODE = COH_CODE
869 self.PLOT_CODE = COH_CODE
869
870
870 self.FTP_WEI = None
871 self.FTP_WEI = None
871 self.EXP_CODE = None
872 self.EXP_CODE = None
872 self.SUB_EXP_CODE = None
873 self.SUB_EXP_CODE = None
873 self.PLOT_POS = None
874 self.PLOT_POS = None
874 self.counter_imagwr = 0
875 self.counter_imagwr = 0
875
876
876 self.xmin = None
877 self.xmin = None
877 self.xmax = None
878 self.xmax = None
878
879
879 def getSubplots(self):
880 def getSubplots(self):
880 ncol = 1
881 ncol = 1
881 nrow = self.nplots*2
882 nrow = self.nplots*2
882
883
883 return nrow, ncol
884 return nrow, ncol
884
885
885 def setup(self, id, nplots, wintitle, showprofile=True, show=True):
886 def setup(self, id, nplots, wintitle, showprofile=True, show=True):
886 self.__showprofile = showprofile
887 self.__showprofile = showprofile
887 self.nplots = nplots
888 self.nplots = nplots
888
889
889 ncolspan = 1
890 ncolspan = 1
890 colspan = 1
891 colspan = 1
891 if showprofile:
892 if showprofile:
892 ncolspan = 7
893 ncolspan = 7
893 colspan = 6
894 colspan = 6
894 self.__nsubplots = 2
895 self.__nsubplots = 2
895
896
896 self.createFigure(id = id,
897 self.createFigure(id = id,
897 wintitle = wintitle,
898 wintitle = wintitle,
898 widthplot = self.WIDTH + self.WIDTHPROF,
899 widthplot = self.WIDTH + self.WIDTHPROF,
899 heightplot = self.HEIGHT + self.HEIGHTPROF,
900 heightplot = self.HEIGHT + self.HEIGHTPROF,
900 show=True)
901 show=True)
901
902
902 nrow, ncol = self.getSubplots()
903 nrow, ncol = self.getSubplots()
903
904
904 for y in range(nrow):
905 for y in range(nrow):
905 for x in range(ncol):
906 for x in range(ncol):
906
907
907 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan, colspan, 1)
908 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan, colspan, 1)
908
909
909 if showprofile:
910 if showprofile:
910 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan+colspan, 1, 1)
911 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan+colspan, 1, 1)
911
912
912 def run(self, dataOut, id, wintitle="", pairsList=None, showprofile='True',
913 def run(self, dataOut, id, wintitle="", pairsList=None, showprofile='True',
913 xmin=None, xmax=None, ymin=None, ymax=None, zmin=None, zmax=None,
914 xmin=None, xmax=None, ymin=None, ymax=None, zmin=None, zmax=None,
914 timerange=None, phase_min=None, phase_max=None,
915 timerange=None, phase_min=None, phase_max=None,
915 save=False, figpath='./', figfile=None, ftp=False, wr_period=1,
916 save=False, figpath='./', figfile=None, ftp=False, wr_period=1,
916 coherence_cmap='jet', phase_cmap='RdBu_r', show=True,
917 coherence_cmap='jet', phase_cmap='RdBu_r', show=True,
917 server=None, folder=None, username=None, password=None,
918 server=None, folder=None, username=None, password=None,
918 ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0):
919 ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0):
919
920
920 if not isTimeInHourRange(dataOut.datatime, xmin, xmax):
921 if not isTimeInHourRange(dataOut.datatime, xmin, xmax):
921 return
922 return
922
923
923 if pairsList == None:
924 if pairsList == None:
924 pairsIndexList = dataOut.pairsIndexList
925 pairsIndexList = dataOut.pairsIndexList
925 else:
926 else:
926 pairsIndexList = []
927 pairsIndexList = []
927 for pair in pairsList:
928 for pair in pairsList:
928 if pair not in dataOut.pairsList:
929 if pair not in dataOut.pairsList:
929 raise ValueError, "Pair %s is not in dataOut.pairsList" %(pair)
930 raise ValueError, "Pair %s is not in dataOut.pairsList" %(pair)
930 pairsIndexList.append(dataOut.pairsList.index(pair))
931 pairsIndexList.append(dataOut.pairsList.index(pair))
931
932
932 if pairsIndexList == []:
933 if pairsIndexList == []:
933 return
934 return
934
935
935 if len(pairsIndexList) > 4:
936 if len(pairsIndexList) > 4:
936 pairsIndexList = pairsIndexList[0:4]
937 pairsIndexList = pairsIndexList[0:4]
937
938
938 if phase_min == None:
939 if phase_min == None:
939 phase_min = -180
940 phase_min = -180
940 if phase_max == None:
941 if phase_max == None:
941 phase_max = 180
942 phase_max = 180
942
943
943 x = dataOut.getTimeRange()
944 x = dataOut.getTimeRange()
944 y = dataOut.getHeiRange()
945 y = dataOut.getHeiRange()
945
946
946 thisDatetime = dataOut.datatime
947 thisDatetime = dataOut.datatime
947
948
948 title = wintitle + " CoherenceMap" #: %s" %(thisDatetime.strftime("%d-%b-%Y"))
949 title = wintitle + " CoherenceMap" #: %s" %(thisDatetime.strftime("%d-%b-%Y"))
949 xlabel = ""
950 xlabel = ""
950 ylabel = "Range (Km)"
951 ylabel = "Range (Km)"
951 update_figfile = False
952 update_figfile = False
952
953
953 if not self.isConfig:
954 if not self.isConfig:
954 nplots = len(pairsIndexList)
955 nplots = len(pairsIndexList)
955 self.setup(id=id,
956 self.setup(id=id,
956 nplots=nplots,
957 nplots=nplots,
957 wintitle=wintitle,
958 wintitle=wintitle,
958 showprofile=showprofile,
959 showprofile=showprofile,
959 show=show)
960 show=show)
960
961
961 if timerange != None:
962 if timerange != None:
962 self.timerange = timerange
963 self.timerange = timerange
963
964
964 self.xmin, self.xmax = self.getTimeLim(x, xmin, xmax, timerange)
965 self.xmin, self.xmax = self.getTimeLim(x, xmin, xmax, timerange)
965
966
966 if ymin == None: ymin = numpy.nanmin(y)
967 if ymin == None: ymin = numpy.nanmin(y)
967 if ymax == None: ymax = numpy.nanmax(y)
968 if ymax == None: ymax = numpy.nanmax(y)
968 if zmin == None: zmin = 0.
969 if zmin == None: zmin = 0.
969 if zmax == None: zmax = 1.
970 if zmax == None: zmax = 1.
970
971
971 self.FTP_WEI = ftp_wei
972 self.FTP_WEI = ftp_wei
972 self.EXP_CODE = exp_code
973 self.EXP_CODE = exp_code
973 self.SUB_EXP_CODE = sub_exp_code
974 self.SUB_EXP_CODE = sub_exp_code
974 self.PLOT_POS = plot_pos
975 self.PLOT_POS = plot_pos
975
976
976 self.name = thisDatetime.strftime("%Y%m%d_%H%M%S")
977 self.name = thisDatetime.strftime("%Y%m%d_%H%M%S")
977
978
978 self.isConfig = True
979 self.isConfig = True
979 update_figfile = True
980 update_figfile = True
980
981
981 self.setWinTitle(title)
982 self.setWinTitle(title)
982
983
983 for i in range(self.nplots):
984 for i in range(self.nplots):
984
985
985 pair = dataOut.pairsList[pairsIndexList[i]]
986 pair = dataOut.pairsList[pairsIndexList[i]]
986
987
987 ccf = numpy.average(dataOut.data_cspc[pairsIndexList[i],:,:],axis=0)
988 ccf = numpy.average(dataOut.data_cspc[pairsIndexList[i],:,:],axis=0)
988 powa = numpy.average(dataOut.data_spc[pair[0],:,:],axis=0)
989 powa = numpy.average(dataOut.data_spc[pair[0],:,:],axis=0)
989 powb = numpy.average(dataOut.data_spc[pair[1],:,:],axis=0)
990 powb = numpy.average(dataOut.data_spc[pair[1],:,:],axis=0)
990
991
991
992
992 avgcoherenceComplex = ccf/numpy.sqrt(powa*powb)
993 avgcoherenceComplex = ccf/numpy.sqrt(powa*powb)
993 coherence = numpy.abs(avgcoherenceComplex)
994 coherence = numpy.abs(avgcoherenceComplex)
994
995
995 z = coherence.reshape((1,-1))
996 z = coherence.reshape((1,-1))
996
997
997 counter = 0
998 counter = 0
998
999
999 title = "Coherence Ch%d * Ch%d: %s" %(pair[0], pair[1], thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
1000 title = "Coherence Ch%d * Ch%d: %s" %(pair[0], pair[1], thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
1000 axes = self.axesList[i*self.__nsubplots*2]
1001 axes = self.axesList[i*self.__nsubplots*2]
1001 axes.pcolorbuffer(x, y, z,
1002 axes.pcolorbuffer(x, y, z,
1002 xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax, zmin=zmin, zmax=zmax,
1003 xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax, zmin=zmin, zmax=zmax,
1003 xlabel=xlabel, ylabel=ylabel, title=title, rti=True, XAxisAsTime=True,
1004 xlabel=xlabel, ylabel=ylabel, title=title, rti=True, XAxisAsTime=True,
1004 ticksize=9, cblabel='', colormap=coherence_cmap, cbsize="1%")
1005 ticksize=9, cblabel='', colormap=coherence_cmap, cbsize="1%")
1005
1006
1006 if self.__showprofile:
1007 if self.__showprofile:
1007 counter += 1
1008 counter += 1
1008 axes = self.axesList[i*self.__nsubplots*2 + counter]
1009 axes = self.axesList[i*self.__nsubplots*2 + counter]
1009 axes.pline(coherence, y,
1010 axes.pline(coherence, y,
1010 xmin=zmin, xmax=zmax, ymin=ymin, ymax=ymax,
1011 xmin=zmin, xmax=zmax, ymin=ymin, ymax=ymax,
1011 xlabel='', ylabel='', title='', ticksize=7,
1012 xlabel='', ylabel='', title='', ticksize=7,
1012 ytick_visible=False, nxticks=5,
1013 ytick_visible=False, nxticks=5,
1013 grid='x')
1014 grid='x')
1014
1015
1015 counter += 1
1016 counter += 1
1016
1017
1017 phase = numpy.arctan2(avgcoherenceComplex.imag, avgcoherenceComplex.real)*180/numpy.pi
1018 phase = numpy.arctan2(avgcoherenceComplex.imag, avgcoherenceComplex.real)*180/numpy.pi
1018
1019
1019 z = phase.reshape((1,-1))
1020 z = phase.reshape((1,-1))
1020
1021
1021 title = "Phase Ch%d * Ch%d: %s" %(pair[0], pair[1], thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
1022 title = "Phase Ch%d * Ch%d: %s" %(pair[0], pair[1], thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
1022 axes = self.axesList[i*self.__nsubplots*2 + counter]
1023 axes = self.axesList[i*self.__nsubplots*2 + counter]
1023 axes.pcolorbuffer(x, y, z,
1024 axes.pcolorbuffer(x, y, z,
1024 xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax, zmin=phase_min, zmax=phase_max,
1025 xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax, zmin=phase_min, zmax=phase_max,
1025 xlabel=xlabel, ylabel=ylabel, title=title, rti=True, XAxisAsTime=True,
1026 xlabel=xlabel, ylabel=ylabel, title=title, rti=True, XAxisAsTime=True,
1026 ticksize=9, cblabel='', colormap=phase_cmap, cbsize="1%")
1027 ticksize=9, cblabel='', colormap=phase_cmap, cbsize="1%")
1027
1028
1028 if self.__showprofile:
1029 if self.__showprofile:
1029 counter += 1
1030 counter += 1
1030 axes = self.axesList[i*self.__nsubplots*2 + counter]
1031 axes = self.axesList[i*self.__nsubplots*2 + counter]
1031 axes.pline(phase, y,
1032 axes.pline(phase, y,
1032 xmin=phase_min, xmax=phase_max, ymin=ymin, ymax=ymax,
1033 xmin=phase_min, xmax=phase_max, ymin=ymin, ymax=ymax,
1033 xlabel='', ylabel='', title='', ticksize=7,
1034 xlabel='', ylabel='', title='', ticksize=7,
1034 ytick_visible=False, nxticks=4,
1035 ytick_visible=False, nxticks=4,
1035 grid='x')
1036 grid='x')
1036
1037
1037 self.draw()
1038 self.draw()
1038
1039
1039 if dataOut.ltctime >= self.xmax:
1040 if dataOut.ltctime >= self.xmax:
1040 self.counter_imagwr = wr_period
1041 self.counter_imagwr = wr_period
1041 self.isConfig = False
1042 self.isConfig = False
1042 update_figfile = True
1043 update_figfile = True
1043
1044
1044 self.save(figpath=figpath,
1045 self.save(figpath=figpath,
1045 figfile=figfile,
1046 figfile=figfile,
1046 save=save,
1047 save=save,
1047 ftp=ftp,
1048 ftp=ftp,
1048 wr_period=wr_period,
1049 wr_period=wr_period,
1049 thisDatetime=thisDatetime,
1050 thisDatetime=thisDatetime,
1050 update_figfile=update_figfile)
1051 update_figfile=update_figfile)
1051
1052
1052 class PowerProfilePlot(Figure):
1053 class PowerProfilePlot(Figure):
1053
1054
1054 isConfig = None
1055 isConfig = None
1055 __nsubplots = None
1056 __nsubplots = None
1056
1057
1057 WIDTHPROF = None
1058 WIDTHPROF = None
1058 HEIGHTPROF = None
1059 HEIGHTPROF = None
1059 PREFIX = 'spcprofile'
1060 PREFIX = 'spcprofile'
1060
1061
1061 def __init__(self, **kwargs):
1062 def __init__(self, **kwargs):
1062 Figure.__init__(self, **kwargs)
1063 Figure.__init__(self, **kwargs)
1063 self.isConfig = False
1064 self.isConfig = False
1064 self.__nsubplots = 1
1065 self.__nsubplots = 1
1065
1066
1066 self.PLOT_CODE = POWER_CODE
1067 self.PLOT_CODE = POWER_CODE
1067
1068
1068 self.WIDTH = 300
1069 self.WIDTH = 300
1069 self.HEIGHT = 500
1070 self.HEIGHT = 500
1070 self.counter_imagwr = 0
1071 self.counter_imagwr = 0
1071
1072
1072 def getSubplots(self):
1073 def getSubplots(self):
1073 ncol = 1
1074 ncol = 1
1074 nrow = 1
1075 nrow = 1
1075
1076
1076 return nrow, ncol
1077 return nrow, ncol
1077
1078
1078 def setup(self, id, nplots, wintitle, show):
1079 def setup(self, id, nplots, wintitle, show):
1079
1080
1080 self.nplots = nplots
1081 self.nplots = nplots
1081
1082
1082 ncolspan = 1
1083 ncolspan = 1
1083 colspan = 1
1084 colspan = 1
1084
1085
1085 self.createFigure(id = id,
1086 self.createFigure(id = id,
1086 wintitle = wintitle,
1087 wintitle = wintitle,
1087 widthplot = self.WIDTH,
1088 widthplot = self.WIDTH,
1088 heightplot = self.HEIGHT,
1089 heightplot = self.HEIGHT,
1089 show=show)
1090 show=show)
1090
1091
1091 nrow, ncol = self.getSubplots()
1092 nrow, ncol = self.getSubplots()
1092
1093
1093 counter = 0
1094 counter = 0
1094 for y in range(nrow):
1095 for y in range(nrow):
1095 for x in range(ncol):
1096 for x in range(ncol):
1096 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan, colspan, 1)
1097 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan, colspan, 1)
1097
1098
1098 def run(self, dataOut, id, wintitle="", channelList=None,
1099 def run(self, dataOut, id, wintitle="", channelList=None,
1099 xmin=None, xmax=None, ymin=None, ymax=None,
1100 xmin=None, xmax=None, ymin=None, ymax=None,
1100 save=False, figpath='./', figfile=None, show=True,
1101 save=False, figpath='./', figfile=None, show=True,
1101 ftp=False, wr_period=1, server=None,
1102 ftp=False, wr_period=1, server=None,
1102 folder=None, username=None, password=None):
1103 folder=None, username=None, password=None):
1103
1104
1104
1105
1105 if channelList == None:
1106 if channelList == None:
1106 channelIndexList = dataOut.channelIndexList
1107 channelIndexList = dataOut.channelIndexList
1107 channelList = dataOut.channelList
1108 channelList = dataOut.channelList
1108 else:
1109 else:
1109 channelIndexList = []
1110 channelIndexList = []
1110 for channel in channelList:
1111 for channel in channelList:
1111 if channel not in dataOut.channelList:
1112 if channel not in dataOut.channelList:
1112 raise ValueError, "Channel %d is not in dataOut.channelList"
1113 raise ValueError, "Channel %d is not in dataOut.channelList"
1113 channelIndexList.append(dataOut.channelList.index(channel))
1114 channelIndexList.append(dataOut.channelList.index(channel))
1114
1115
1115 factor = dataOut.normFactor
1116 factor = dataOut.normFactor
1116
1117
1117 y = dataOut.getHeiRange()
1118 y = dataOut.getHeiRange()
1118
1119
1119 #for voltage
1120 #for voltage
1120 if dataOut.type == 'Voltage':
1121 if dataOut.type == 'Voltage':
1121 x = dataOut.data[channelIndexList,:] * numpy.conjugate(dataOut.data[channelIndexList,:])
1122 x = dataOut.data[channelIndexList,:] * numpy.conjugate(dataOut.data[channelIndexList,:])
1122 x = x.real
1123 x = x.real
1123 x = numpy.where(numpy.isfinite(x), x, numpy.NAN)
1124 x = numpy.where(numpy.isfinite(x), x, numpy.NAN)
1124
1125
1125 #for spectra
1126 #for spectra
1126 if dataOut.type == 'Spectra':
1127 if dataOut.type == 'Spectra':
1127 x = dataOut.data_spc[channelIndexList,:,:]/factor
1128 x = dataOut.data_spc[channelIndexList,:,:]/factor
1128 x = numpy.where(numpy.isfinite(x), x, numpy.NAN)
1129 x = numpy.where(numpy.isfinite(x), x, numpy.NAN)
1129 x = numpy.average(x, axis=1)
1130 x = numpy.average(x, axis=1)
1130
1131
1131
1132
1132 xdB = 10*numpy.log10(x)
1133 xdB = 10*numpy.log10(x)
1133
1134
1134 thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[0])
1135 thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[0])
1135 title = wintitle + " Power Profile %s" %(thisDatetime.strftime("%d-%b-%Y"))
1136 title = wintitle + " Power Profile %s" %(thisDatetime.strftime("%d-%b-%Y"))
1136 xlabel = "dB"
1137 xlabel = "dB"
1137 ylabel = "Range (Km)"
1138 ylabel = "Range (Km)"
1138
1139
1139 if not self.isConfig:
1140 if not self.isConfig:
1140
1141
1141 nplots = 1
1142 nplots = 1
1142
1143
1143 self.setup(id=id,
1144 self.setup(id=id,
1144 nplots=nplots,
1145 nplots=nplots,
1145 wintitle=wintitle,
1146 wintitle=wintitle,
1146 show=show)
1147 show=show)
1147
1148
1148 if ymin == None: ymin = numpy.nanmin(y)
1149 if ymin == None: ymin = numpy.nanmin(y)
1149 if ymax == None: ymax = numpy.nanmax(y)
1150 if ymax == None: ymax = numpy.nanmax(y)
1150 if xmin == None: xmin = numpy.nanmin(xdB)*0.9
1151 if xmin == None: xmin = numpy.nanmin(xdB)*0.9
1151 if xmax == None: xmax = numpy.nanmax(xdB)*1.1
1152 if xmax == None: xmax = numpy.nanmax(xdB)*1.1
1152
1153
1153 self.isConfig = True
1154 self.isConfig = True
1154
1155
1155 self.setWinTitle(title)
1156 self.setWinTitle(title)
1156
1157
1157 title = "Power Profile: %s" %(thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
1158 title = "Power Profile: %s" %(thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
1158 axes = self.axesList[0]
1159 axes = self.axesList[0]
1159
1160
1160 legendlabels = ["channel %d"%x for x in channelList]
1161 legendlabels = ["channel %d"%x for x in channelList]
1161 axes.pmultiline(xdB, y,
1162 axes.pmultiline(xdB, y,
1162 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax,
1163 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax,
1163 xlabel=xlabel, ylabel=ylabel, title=title, legendlabels=legendlabels,
1164 xlabel=xlabel, ylabel=ylabel, title=title, legendlabels=legendlabels,
1164 ytick_visible=True, nxticks=5,
1165 ytick_visible=True, nxticks=5,
1165 grid='x')
1166 grid='x')
1166
1167
1167 self.draw()
1168 self.draw()
1168
1169
1169 self.save(figpath=figpath,
1170 self.save(figpath=figpath,
1170 figfile=figfile,
1171 figfile=figfile,
1171 save=save,
1172 save=save,
1172 ftp=ftp,
1173 ftp=ftp,
1173 wr_period=wr_period,
1174 wr_period=wr_period,
1174 thisDatetime=thisDatetime)
1175 thisDatetime=thisDatetime)
1175
1176
1176 class SpectraCutPlot(Figure):
1177 class SpectraCutPlot(Figure):
1177
1178
1178 isConfig = None
1179 isConfig = None
1179 __nsubplots = None
1180 __nsubplots = None
1180
1181
1181 WIDTHPROF = None
1182 WIDTHPROF = None
1182 HEIGHTPROF = None
1183 HEIGHTPROF = None
1183 PREFIX = 'spc_cut'
1184 PREFIX = 'spc_cut'
1184
1185
1185 def __init__(self, **kwargs):
1186 def __init__(self, **kwargs):
1186 Figure.__init__(self, **kwargs)
1187 Figure.__init__(self, **kwargs)
1187 self.isConfig = False
1188 self.isConfig = False
1188 self.__nsubplots = 1
1189 self.__nsubplots = 1
1189
1190
1190 self.PLOT_CODE = POWER_CODE
1191 self.PLOT_CODE = POWER_CODE
1191
1192
1192 self.WIDTH = 700
1193 self.WIDTH = 700
1193 self.HEIGHT = 500
1194 self.HEIGHT = 500
1194 self.counter_imagwr = 0
1195 self.counter_imagwr = 0
1195
1196
1196 def getSubplots(self):
1197 def getSubplots(self):
1197 ncol = 1
1198 ncol = 1
1198 nrow = 1
1199 nrow = 1
1199
1200
1200 return nrow, ncol
1201 return nrow, ncol
1201
1202
1202 def setup(self, id, nplots, wintitle, show):
1203 def setup(self, id, nplots, wintitle, show):
1203
1204
1204 self.nplots = nplots
1205 self.nplots = nplots
1205
1206
1206 ncolspan = 1
1207 ncolspan = 1
1207 colspan = 1
1208 colspan = 1
1208
1209
1209 self.createFigure(id = id,
1210 self.createFigure(id = id,
1210 wintitle = wintitle,
1211 wintitle = wintitle,
1211 widthplot = self.WIDTH,
1212 widthplot = self.WIDTH,
1212 heightplot = self.HEIGHT,
1213 heightplot = self.HEIGHT,
1213 show=show)
1214 show=show)
1214
1215
1215 nrow, ncol = self.getSubplots()
1216 nrow, ncol = self.getSubplots()
1216
1217
1217 counter = 0
1218 counter = 0
1218 for y in range(nrow):
1219 for y in range(nrow):
1219 for x in range(ncol):
1220 for x in range(ncol):
1220 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan, colspan, 1)
1221 self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan, colspan, 1)
1221
1222
1222 def run(self, dataOut, id, wintitle="", channelList=None,
1223 def run(self, dataOut, id, wintitle="", channelList=None,
1223 xmin=None, xmax=None, ymin=None, ymax=None,
1224 xmin=None, xmax=None, ymin=None, ymax=None,
1224 save=False, figpath='./', figfile=None, show=True,
1225 save=False, figpath='./', figfile=None, show=True,
1225 ftp=False, wr_period=1, server=None,
1226 ftp=False, wr_period=1, server=None,
1226 folder=None, username=None, password=None,
1227 folder=None, username=None, password=None,
1227 xaxis="frequency"):
1228 xaxis="frequency"):
1228
1229
1229
1230
1230 if channelList == None:
1231 if channelList == None:
1231 channelIndexList = dataOut.channelIndexList
1232 channelIndexList = dataOut.channelIndexList
1232 channelList = dataOut.channelList
1233 channelList = dataOut.channelList
1233 else:
1234 else:
1234 channelIndexList = []
1235 channelIndexList = []
1235 for channel in channelList:
1236 for channel in channelList:
1236 if channel not in dataOut.channelList:
1237 if channel not in dataOut.channelList:
1237 raise ValueError, "Channel %d is not in dataOut.channelList"
1238 raise ValueError, "Channel %d is not in dataOut.channelList"
1238 channelIndexList.append(dataOut.channelList.index(channel))
1239 channelIndexList.append(dataOut.channelList.index(channel))
1239
1240
1240 factor = dataOut.normFactor
1241 factor = dataOut.normFactor
1241
1242
1242 y = dataOut.getHeiRange()
1243 y = dataOut.getHeiRange()
1243
1244
1244 z = dataOut.data_spc/factor
1245 z = dataOut.data_spc/factor
1245 z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
1246 z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
1246
1247
1247 hei_index = numpy.arange(25)*3 + 20
1248 hei_index = numpy.arange(25)*3 + 20
1248
1249
1249 if xaxis == "frequency":
1250 if xaxis == "frequency":
1250 x = dataOut.getFreqRange()/1000.
1251 x = dataOut.getFreqRange()/1000.
1251 zdB = 10*numpy.log10(z[0,:,hei_index])
1252 zdB = 10*numpy.log10(z[0,:,hei_index])
1252 xlabel = "Frequency (kHz)"
1253 xlabel = "Frequency (kHz)"
1253 ylabel = "Power (dB)"
1254 ylabel = "Power (dB)"
1254
1255
1255 elif xaxis == "time":
1256 elif xaxis == "time":
1256 x = dataOut.getAcfRange()
1257 x = dataOut.getAcfRange()
1257 zdB = z[0,:,hei_index]
1258 zdB = z[0,:,hei_index]
1258 xlabel = "Time (ms)"
1259 xlabel = "Time (ms)"
1259 ylabel = "ACF"
1260 ylabel = "ACF"
1260
1261
1261 else:
1262 else:
1262 x = dataOut.getVelRange()
1263 x = dataOut.getVelRange()
1263 zdB = 10*numpy.log10(z[0,:,hei_index])
1264 zdB = 10*numpy.log10(z[0,:,hei_index])
1264 xlabel = "Velocity (m/s)"
1265 xlabel = "Velocity (m/s)"
1265 ylabel = "Power (dB)"
1266 ylabel = "Power (dB)"
1266
1267
1267 thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[0])
1268 thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[0])
1268 title = wintitle + " Range Cuts %s" %(thisDatetime.strftime("%d-%b-%Y"))
1269 title = wintitle + " Range Cuts %s" %(thisDatetime.strftime("%d-%b-%Y"))
1269
1270
1270 if not self.isConfig:
1271 if not self.isConfig:
1271
1272
1272 nplots = 1
1273 nplots = 1
1273
1274
1274 self.setup(id=id,
1275 self.setup(id=id,
1275 nplots=nplots,
1276 nplots=nplots,
1276 wintitle=wintitle,
1277 wintitle=wintitle,
1277 show=show)
1278 show=show)
1278
1279
1279 if xmin == None: xmin = numpy.nanmin(x)*0.9
1280 if xmin == None: xmin = numpy.nanmin(x)*0.9
1280 if xmax == None: xmax = numpy.nanmax(x)*1.1
1281 if xmax == None: xmax = numpy.nanmax(x)*1.1
1281 if ymin == None: ymin = numpy.nanmin(zdB)
1282 if ymin == None: ymin = numpy.nanmin(zdB)
1282 if ymax == None: ymax = numpy.nanmax(zdB)
1283 if ymax == None: ymax = numpy.nanmax(zdB)
1283
1284
1284 self.isConfig = True
1285 self.isConfig = True
1285
1286
1286 self.setWinTitle(title)
1287 self.setWinTitle(title)
1287
1288
1288 title = "Spectra Cuts: %s" %(thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
1289 title = "Spectra Cuts: %s" %(thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
1289 axes = self.axesList[0]
1290 axes = self.axesList[0]
1290
1291
1291 legendlabels = ["Range = %dKm" %y[i] for i in hei_index]
1292 legendlabels = ["Range = %dKm" %y[i] for i in hei_index]
1292
1293
1293 axes.pmultilineyaxis( x, zdB,
1294 axes.pmultilineyaxis( x, zdB,
1294 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax,
1295 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax,
1295 xlabel=xlabel, ylabel=ylabel, title=title, legendlabels=legendlabels,
1296 xlabel=xlabel, ylabel=ylabel, title=title, legendlabels=legendlabels,
1296 ytick_visible=True, nxticks=5,
1297 ytick_visible=True, nxticks=5,
1297 grid='x')
1298 grid='x')
1298
1299
1299 self.draw()
1300 self.draw()
1300
1301
1301 self.save(figpath=figpath,
1302 self.save(figpath=figpath,
1302 figfile=figfile,
1303 figfile=figfile,
1303 save=save,
1304 save=save,
1304 ftp=ftp,
1305 ftp=ftp,
1305 wr_period=wr_period,
1306 wr_period=wr_period,
1306 thisDatetime=thisDatetime)
1307 thisDatetime=thisDatetime)
1307
1308
1308 class Noise(Figure):
1309 class Noise(Figure):
1309
1310
1310 isConfig = None
1311 isConfig = None
1311 __nsubplots = None
1312 __nsubplots = None
1312
1313
1313 PREFIX = 'noise'
1314 PREFIX = 'noise'
1314
1315
1315
1316
1316 def __init__(self, **kwargs):
1317 def __init__(self, **kwargs):
1317 Figure.__init__(self, **kwargs)
1318 Figure.__init__(self, **kwargs)
1318 self.timerange = 24*60*60
1319 self.timerange = 24*60*60
1319 self.isConfig = False
1320 self.isConfig = False
1320 self.__nsubplots = 1
1321 self.__nsubplots = 1
1321 self.counter_imagwr = 0
1322 self.counter_imagwr = 0
1322 self.WIDTH = 800
1323 self.WIDTH = 800
1323 self.HEIGHT = 400
1324 self.HEIGHT = 400
1324 self.WIDTHPROF = 120
1325 self.WIDTHPROF = 120
1325 self.HEIGHTPROF = 0
1326 self.HEIGHTPROF = 0
1326 self.xdata = None
1327 self.xdata = None
1327 self.ydata = None
1328 self.ydata = None
1328
1329
1329 self.PLOT_CODE = NOISE_CODE
1330 self.PLOT_CODE = NOISE_CODE
1330
1331
1331 self.FTP_WEI = None
1332 self.FTP_WEI = None
1332 self.EXP_CODE = None
1333 self.EXP_CODE = None
1333 self.SUB_EXP_CODE = None
1334 self.SUB_EXP_CODE = None
1334 self.PLOT_POS = None
1335 self.PLOT_POS = None
1335 self.figfile = None
1336 self.figfile = None
1336
1337
1337 self.xmin = None
1338 self.xmin = None
1338 self.xmax = None
1339 self.xmax = None
1339
1340
1340 def getSubplots(self):
1341 def getSubplots(self):
1341
1342
1342 ncol = 1
1343 ncol = 1
1343 nrow = 1
1344 nrow = 1
1344
1345
1345 return nrow, ncol
1346 return nrow, ncol
1346
1347
1347 def openfile(self, filename):
1348 def openfile(self, filename):
1348 dirname = os.path.dirname(filename)
1349 dirname = os.path.dirname(filename)
1349
1350
1350 if not os.path.exists(dirname):
1351 if not os.path.exists(dirname):
1351 os.mkdir(dirname)
1352 os.mkdir(dirname)
1352
1353
1353 f = open(filename,'w+')
1354 f = open(filename,'w+')
1354 f.write('\n\n')
1355 f.write('\n\n')
1355 f.write('JICAMARCA RADIO OBSERVATORY - Noise \n')
1356 f.write('JICAMARCA RADIO OBSERVATORY - Noise \n')
1356 f.write('DD MM YYYY HH MM SS Channel0 Channel1 Channel2 Channel3\n\n' )
1357 f.write('DD MM YYYY HH MM SS Channel0 Channel1 Channel2 Channel3\n\n' )
1357 f.close()
1358 f.close()
1358
1359
1359 def save_data(self, filename_phase, data, data_datetime):
1360 def save_data(self, filename_phase, data, data_datetime):
1360
1361
1361 f=open(filename_phase,'a')
1362 f=open(filename_phase,'a')
1362
1363
1363 timetuple_data = data_datetime.timetuple()
1364 timetuple_data = data_datetime.timetuple()
1364 day = str(timetuple_data.tm_mday)
1365 day = str(timetuple_data.tm_mday)
1365 month = str(timetuple_data.tm_mon)
1366 month = str(timetuple_data.tm_mon)
1366 year = str(timetuple_data.tm_year)
1367 year = str(timetuple_data.tm_year)
1367 hour = str(timetuple_data.tm_hour)
1368 hour = str(timetuple_data.tm_hour)
1368 minute = str(timetuple_data.tm_min)
1369 minute = str(timetuple_data.tm_min)
1369 second = str(timetuple_data.tm_sec)
1370 second = str(timetuple_data.tm_sec)
1370
1371
1371 data_msg = ''
1372 data_msg = ''
1372 for i in range(len(data)):
1373 for i in range(len(data)):
1373 data_msg += str(data[i]) + ' '
1374 data_msg += str(data[i]) + ' '
1374
1375
1375 f.write(day+' '+month+' '+year+' '+hour+' '+minute+' '+second+' ' + data_msg + '\n')
1376 f.write(day+' '+month+' '+year+' '+hour+' '+minute+' '+second+' ' + data_msg + '\n')
1376 f.close()
1377 f.close()
1377
1378
1378
1379
1379 def setup(self, id, nplots, wintitle, showprofile=True, show=True):
1380 def setup(self, id, nplots, wintitle, showprofile=True, show=True):
1380
1381
1381 self.__showprofile = showprofile
1382 self.__showprofile = showprofile
1382 self.nplots = nplots
1383 self.nplots = nplots
1383
1384
1384 ncolspan = 7
1385 ncolspan = 7
1385 colspan = 6
1386 colspan = 6
1386 self.__nsubplots = 2
1387 self.__nsubplots = 2
1387
1388
1388 self.createFigure(id = id,
1389 self.createFigure(id = id,
1389 wintitle = wintitle,
1390 wintitle = wintitle,
1390 widthplot = self.WIDTH+self.WIDTHPROF,
1391 widthplot = self.WIDTH+self.WIDTHPROF,
1391 heightplot = self.HEIGHT+self.HEIGHTPROF,
1392 heightplot = self.HEIGHT+self.HEIGHTPROF,
1392 show=show)
1393 show=show)
1393
1394
1394 nrow, ncol = self.getSubplots()
1395 nrow, ncol = self.getSubplots()
1395
1396
1396 self.addAxes(nrow, ncol*ncolspan, 0, 0, colspan, 1)
1397 self.addAxes(nrow, ncol*ncolspan, 0, 0, colspan, 1)
1397
1398
1398
1399
1399 def run(self, dataOut, id, wintitle="", channelList=None, showprofile='True',
1400 def run(self, dataOut, id, wintitle="", channelList=None, showprofile='True',
1400 xmin=None, xmax=None, ymin=None, ymax=None,
1401 xmin=None, xmax=None, ymin=None, ymax=None,
1401 timerange=None,
1402 timerange=None,
1402 save=False, figpath='./', figfile=None, show=True, ftp=False, wr_period=1,
1403 save=False, figpath='./', figfile=None, show=True, ftp=False, wr_period=1,
1403 server=None, folder=None, username=None, password=None,
1404 server=None, folder=None, username=None, password=None,
1404 ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0):
1405 ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0):
1405
1406
1406 if not isTimeInHourRange(dataOut.datatime, xmin, xmax):
1407 if not isTimeInHourRange(dataOut.datatime, xmin, xmax):
1407 return
1408 return
1408
1409
1409 if channelList == None:
1410 if channelList == None:
1410 channelIndexList = dataOut.channelIndexList
1411 channelIndexList = dataOut.channelIndexList
1411 channelList = dataOut.channelList
1412 channelList = dataOut.channelList
1412 else:
1413 else:
1413 channelIndexList = []
1414 channelIndexList = []
1414 for channel in channelList:
1415 for channel in channelList:
1415 if channel not in dataOut.channelList:
1416 if channel not in dataOut.channelList:
1416 raise ValueError, "Channel %d is not in dataOut.channelList"
1417 raise ValueError, "Channel %d is not in dataOut.channelList"
1417 channelIndexList.append(dataOut.channelList.index(channel))
1418 channelIndexList.append(dataOut.channelList.index(channel))
1418
1419
1419 x = dataOut.getTimeRange()
1420 x = dataOut.getTimeRange()
1420 #y = dataOut.getHeiRange()
1421 #y = dataOut.getHeiRange()
1421 factor = dataOut.normFactor
1422 factor = dataOut.normFactor
1422 noise = dataOut.noise[channelIndexList]/factor
1423 noise = dataOut.noise[channelIndexList]/factor
1423 noisedB = 10*numpy.log10(noise)
1424 noisedB = 10*numpy.log10(noise)
1424
1425
1425 thisDatetime = dataOut.datatime
1426 thisDatetime = dataOut.datatime
1426
1427
1427 title = wintitle + " Noise" # : %s" %(thisDatetime.strftime("%d-%b-%Y"))
1428 title = wintitle + " Noise" # : %s" %(thisDatetime.strftime("%d-%b-%Y"))
1428 xlabel = ""
1429 xlabel = ""
1429 ylabel = "Intensity (dB)"
1430 ylabel = "Intensity (dB)"
1430 update_figfile = False
1431 update_figfile = False
1431
1432
1432 if not self.isConfig:
1433 if not self.isConfig:
1433
1434
1434 nplots = 1
1435 nplots = 1
1435
1436
1436 self.setup(id=id,
1437 self.setup(id=id,
1437 nplots=nplots,
1438 nplots=nplots,
1438 wintitle=wintitle,
1439 wintitle=wintitle,
1439 showprofile=showprofile,
1440 showprofile=showprofile,
1440 show=show)
1441 show=show)
1441
1442
1442 if timerange != None:
1443 if timerange != None:
1443 self.timerange = timerange
1444 self.timerange = timerange
1444
1445
1445 self.xmin, self.xmax = self.getTimeLim(x, xmin, xmax, timerange)
1446 self.xmin, self.xmax = self.getTimeLim(x, xmin, xmax, timerange)
1446
1447
1447 if ymin == None: ymin = numpy.floor(numpy.nanmin(noisedB)) - 10.0
1448 if ymin == None: ymin = numpy.floor(numpy.nanmin(noisedB)) - 10.0
1448 if ymax == None: ymax = numpy.nanmax(noisedB) + 10.0
1449 if ymax == None: ymax = numpy.nanmax(noisedB) + 10.0
1449
1450
1450 self.FTP_WEI = ftp_wei
1451 self.FTP_WEI = ftp_wei
1451 self.EXP_CODE = exp_code
1452 self.EXP_CODE = exp_code
1452 self.SUB_EXP_CODE = sub_exp_code
1453 self.SUB_EXP_CODE = sub_exp_code
1453 self.PLOT_POS = plot_pos
1454 self.PLOT_POS = plot_pos
1454
1455
1455
1456
1456 self.name = thisDatetime.strftime("%Y%m%d_%H%M%S")
1457 self.name = thisDatetime.strftime("%Y%m%d_%H%M%S")
1457 self.isConfig = True
1458 self.isConfig = True
1458 self.figfile = figfile
1459 self.figfile = figfile
1459 self.xdata = numpy.array([])
1460 self.xdata = numpy.array([])
1460 self.ydata = numpy.array([])
1461 self.ydata = numpy.array([])
1461
1462
1462 update_figfile = True
1463 update_figfile = True
1463
1464
1464 #open file beacon phase
1465 #open file beacon phase
1465 path = '%s%03d' %(self.PREFIX, self.id)
1466 path = '%s%03d' %(self.PREFIX, self.id)
1466 noise_file = os.path.join(path,'%s.txt'%self.name)
1467 noise_file = os.path.join(path,'%s.txt'%self.name)
1467 self.filename_noise = os.path.join(figpath,noise_file)
1468 self.filename_noise = os.path.join(figpath,noise_file)
1468
1469
1469 self.setWinTitle(title)
1470 self.setWinTitle(title)
1470
1471
1471 title = "Noise %s" %(thisDatetime.strftime("%Y/%m/%d %H:%M:%S"))
1472 title = "Noise %s" %(thisDatetime.strftime("%Y/%m/%d %H:%M:%S"))
1472
1473
1473 legendlabels = ["channel %d"%(idchannel) for idchannel in channelList]
1474 legendlabels = ["channel %d"%(idchannel) for idchannel in channelList]
1474 axes = self.axesList[0]
1475 axes = self.axesList[0]
1475
1476
1476 self.xdata = numpy.hstack((self.xdata, x[0:1]))
1477 self.xdata = numpy.hstack((self.xdata, x[0:1]))
1477
1478
1478 if len(self.ydata)==0:
1479 if len(self.ydata)==0:
1479 self.ydata = noisedB.reshape(-1,1)
1480 self.ydata = noisedB.reshape(-1,1)
1480 else:
1481 else:
1481 self.ydata = numpy.hstack((self.ydata, noisedB.reshape(-1,1)))
1482 self.ydata = numpy.hstack((self.ydata, noisedB.reshape(-1,1)))
1482
1483
1483
1484
1484 axes.pmultilineyaxis(x=self.xdata, y=self.ydata,
1485 axes.pmultilineyaxis(x=self.xdata, y=self.ydata,
1485 xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax,
1486 xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax,
1486 xlabel=xlabel, ylabel=ylabel, title=title, legendlabels=legendlabels, marker='x', markersize=8, linestyle="solid",
1487 xlabel=xlabel, ylabel=ylabel, title=title, legendlabels=legendlabels, marker='x', markersize=8, linestyle="solid",
1487 XAxisAsTime=True, grid='both'
1488 XAxisAsTime=True, grid='both'
1488 )
1489 )
1489
1490
1490 self.draw()
1491 self.draw()
1491
1492
1492 if dataOut.ltctime >= self.xmax:
1493 if dataOut.ltctime >= self.xmax:
1493 self.counter_imagwr = wr_period
1494 self.counter_imagwr = wr_period
1494 self.isConfig = False
1495 self.isConfig = False
1495 update_figfile = True
1496 update_figfile = True
1496
1497
1497 self.save(figpath=figpath,
1498 self.save(figpath=figpath,
1498 figfile=figfile,
1499 figfile=figfile,
1499 save=save,
1500 save=save,
1500 ftp=ftp,
1501 ftp=ftp,
1501 wr_period=wr_period,
1502 wr_period=wr_period,
1502 thisDatetime=thisDatetime,
1503 thisDatetime=thisDatetime,
1503 update_figfile=update_figfile)
1504 update_figfile=update_figfile)
1504
1505
1505 #store data beacon phase
1506 #store data beacon phase
1506 if save:
1507 if save:
1507 self.save_data(self.filename_noise, noisedB, thisDatetime)
1508 self.save_data(self.filename_noise, noisedB, thisDatetime)
1508
1509
1509 class BeaconPhase(Figure):
1510 class BeaconPhase(Figure):
1510
1511
1511 __isConfig = None
1512 __isConfig = None
1512 __nsubplots = None
1513 __nsubplots = None
1513
1514
1514 PREFIX = 'beacon_phase'
1515 PREFIX = 'beacon_phase'
1515
1516
1516 def __init__(self, **kwargs):
1517 def __init__(self, **kwargs):
1517 Figure.__init__(self, **kwargs)
1518 Figure.__init__(self, **kwargs)
1518 self.timerange = 24*60*60
1519 self.timerange = 24*60*60
1519 self.isConfig = False
1520 self.isConfig = False
1520 self.__nsubplots = 1
1521 self.__nsubplots = 1
1521 self.counter_imagwr = 0
1522 self.counter_imagwr = 0
1522 self.WIDTH = 800
1523 self.WIDTH = 800
1523 self.HEIGHT = 400
1524 self.HEIGHT = 400
1524 self.WIDTHPROF = 120
1525 self.WIDTHPROF = 120
1525 self.HEIGHTPROF = 0
1526 self.HEIGHTPROF = 0
1526 self.xdata = None
1527 self.xdata = None
1527 self.ydata = None
1528 self.ydata = None
1528
1529
1529 self.PLOT_CODE = BEACON_CODE
1530 self.PLOT_CODE = BEACON_CODE
1530
1531
1531 self.FTP_WEI = None
1532 self.FTP_WEI = None
1532 self.EXP_CODE = None
1533 self.EXP_CODE = None
1533 self.SUB_EXP_CODE = None
1534 self.SUB_EXP_CODE = None
1534 self.PLOT_POS = None
1535 self.PLOT_POS = None
1535
1536
1536 self.filename_phase = None
1537 self.filename_phase = None
1537
1538
1538 self.figfile = None
1539 self.figfile = None
1539
1540
1540 self.xmin = None
1541 self.xmin = None
1541 self.xmax = None
1542 self.xmax = None
1542
1543
1543 def getSubplots(self):
1544 def getSubplots(self):
1544
1545
1545 ncol = 1
1546 ncol = 1
1546 nrow = 1
1547 nrow = 1
1547
1548
1548 return nrow, ncol
1549 return nrow, ncol
1549
1550
1550 def setup(self, id, nplots, wintitle, showprofile=True, show=True):
1551 def setup(self, id, nplots, wintitle, showprofile=True, show=True):
1551
1552
1552 self.__showprofile = showprofile
1553 self.__showprofile = showprofile
1553 self.nplots = nplots
1554 self.nplots = nplots
1554
1555
1555 ncolspan = 7
1556 ncolspan = 7
1556 colspan = 6
1557 colspan = 6
1557 self.__nsubplots = 2
1558 self.__nsubplots = 2
1558
1559
1559 self.createFigure(id = id,
1560 self.createFigure(id = id,
1560 wintitle = wintitle,
1561 wintitle = wintitle,
1561 widthplot = self.WIDTH+self.WIDTHPROF,
1562 widthplot = self.WIDTH+self.WIDTHPROF,
1562 heightplot = self.HEIGHT+self.HEIGHTPROF,
1563 heightplot = self.HEIGHT+self.HEIGHTPROF,
1563 show=show)
1564 show=show)
1564
1565
1565 nrow, ncol = self.getSubplots()
1566 nrow, ncol = self.getSubplots()
1566
1567
1567 self.addAxes(nrow, ncol*ncolspan, 0, 0, colspan, 1)
1568 self.addAxes(nrow, ncol*ncolspan, 0, 0, colspan, 1)
1568
1569
1569 def save_phase(self, filename_phase):
1570 def save_phase(self, filename_phase):
1570 f = open(filename_phase,'w+')
1571 f = open(filename_phase,'w+')
1571 f.write('\n\n')
1572 f.write('\n\n')
1572 f.write('JICAMARCA RADIO OBSERVATORY - Beacon Phase \n')
1573 f.write('JICAMARCA RADIO OBSERVATORY - Beacon Phase \n')
1573 f.write('DD MM YYYY HH MM SS pair(2,0) pair(2,1) pair(2,3) pair(2,4)\n\n' )
1574 f.write('DD MM YYYY HH MM SS pair(2,0) pair(2,1) pair(2,3) pair(2,4)\n\n' )
1574 f.close()
1575 f.close()
1575
1576
1576 def save_data(self, filename_phase, data, data_datetime):
1577 def save_data(self, filename_phase, data, data_datetime):
1577 f=open(filename_phase,'a')
1578 f=open(filename_phase,'a')
1578 timetuple_data = data_datetime.timetuple()
1579 timetuple_data = data_datetime.timetuple()
1579 day = str(timetuple_data.tm_mday)
1580 day = str(timetuple_data.tm_mday)
1580 month = str(timetuple_data.tm_mon)
1581 month = str(timetuple_data.tm_mon)
1581 year = str(timetuple_data.tm_year)
1582 year = str(timetuple_data.tm_year)
1582 hour = str(timetuple_data.tm_hour)
1583 hour = str(timetuple_data.tm_hour)
1583 minute = str(timetuple_data.tm_min)
1584 minute = str(timetuple_data.tm_min)
1584 second = str(timetuple_data.tm_sec)
1585 second = str(timetuple_data.tm_sec)
1585 f.write(day+' '+month+' '+year+' '+hour+' '+minute+' '+second+' '+str(data[0])+' '+str(data[1])+' '+str(data[2])+' '+str(data[3])+'\n')
1586 f.write(day+' '+month+' '+year+' '+hour+' '+minute+' '+second+' '+str(data[0])+' '+str(data[1])+' '+str(data[2])+' '+str(data[3])+'\n')
1586 f.close()
1587 f.close()
1587
1588
1588
1589
1589 def run(self, dataOut, id, wintitle="", pairsList=None, showprofile='True',
1590 def run(self, dataOut, id, wintitle="", pairsList=None, showprofile='True',
1590 xmin=None, xmax=None, ymin=None, ymax=None, hmin=None, hmax=None,
1591 xmin=None, xmax=None, ymin=None, ymax=None, hmin=None, hmax=None,
1591 timerange=None,
1592 timerange=None,
1592 save=False, figpath='./', figfile=None, show=True, ftp=False, wr_period=1,
1593 save=False, figpath='./', figfile=None, show=True, ftp=False, wr_period=1,
1593 server=None, folder=None, username=None, password=None,
1594 server=None, folder=None, username=None, password=None,
1594 ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0):
1595 ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0):
1595
1596
1596 if not isTimeInHourRange(dataOut.datatime, xmin, xmax):
1597 if not isTimeInHourRange(dataOut.datatime, xmin, xmax):
1597 return
1598 return
1598
1599
1599 if pairsList == None:
1600 if pairsList == None:
1600 pairsIndexList = dataOut.pairsIndexList[:10]
1601 pairsIndexList = dataOut.pairsIndexList[:10]
1601 else:
1602 else:
1602 pairsIndexList = []
1603 pairsIndexList = []
1603 for pair in pairsList:
1604 for pair in pairsList:
1604 if pair not in dataOut.pairsList:
1605 if pair not in dataOut.pairsList:
1605 raise ValueError, "Pair %s is not in dataOut.pairsList" %(pair)
1606 raise ValueError, "Pair %s is not in dataOut.pairsList" %(pair)
1606 pairsIndexList.append(dataOut.pairsList.index(pair))
1607 pairsIndexList.append(dataOut.pairsList.index(pair))
1607
1608
1608 if pairsIndexList == []:
1609 if pairsIndexList == []:
1609 return
1610 return
1610
1611
1611 # if len(pairsIndexList) > 4:
1612 # if len(pairsIndexList) > 4:
1612 # pairsIndexList = pairsIndexList[0:4]
1613 # pairsIndexList = pairsIndexList[0:4]
1613
1614
1614 hmin_index = None
1615 hmin_index = None
1615 hmax_index = None
1616 hmax_index = None
1616
1617
1617 if hmin != None and hmax != None:
1618 if hmin != None and hmax != None:
1618 indexes = numpy.arange(dataOut.nHeights)
1619 indexes = numpy.arange(dataOut.nHeights)
1619 hmin_list = indexes[dataOut.heightList >= hmin]
1620 hmin_list = indexes[dataOut.heightList >= hmin]
1620 hmax_list = indexes[dataOut.heightList <= hmax]
1621 hmax_list = indexes[dataOut.heightList <= hmax]
1621
1622
1622 if hmin_list.any():
1623 if hmin_list.any():
1623 hmin_index = hmin_list[0]
1624 hmin_index = hmin_list[0]
1624
1625
1625 if hmax_list.any():
1626 if hmax_list.any():
1626 hmax_index = hmax_list[-1]+1
1627 hmax_index = hmax_list[-1]+1
1627
1628
1628 x = dataOut.getTimeRange()
1629 x = dataOut.getTimeRange()
1629 #y = dataOut.getHeiRange()
1630 #y = dataOut.getHeiRange()
1630
1631
1631
1632
1632 thisDatetime = dataOut.datatime
1633 thisDatetime = dataOut.datatime
1633
1634
1634 title = wintitle + " Signal Phase" # : %s" %(thisDatetime.strftime("%d-%b-%Y"))
1635 title = wintitle + " Signal Phase" # : %s" %(thisDatetime.strftime("%d-%b-%Y"))
1635 xlabel = "Local Time"
1636 xlabel = "Local Time"
1636 ylabel = "Phase (degrees)"
1637 ylabel = "Phase (degrees)"
1637
1638
1638 update_figfile = False
1639 update_figfile = False
1639
1640
1640 nplots = len(pairsIndexList)
1641 nplots = len(pairsIndexList)
1641 #phase = numpy.zeros((len(pairsIndexList),len(dataOut.beacon_heiIndexList)))
1642 #phase = numpy.zeros((len(pairsIndexList),len(dataOut.beacon_heiIndexList)))
1642 phase_beacon = numpy.zeros(len(pairsIndexList))
1643 phase_beacon = numpy.zeros(len(pairsIndexList))
1643 for i in range(nplots):
1644 for i in range(nplots):
1644 pair = dataOut.pairsList[pairsIndexList[i]]
1645 pair = dataOut.pairsList[pairsIndexList[i]]
1645 ccf = numpy.average(dataOut.data_cspc[pairsIndexList[i], :, hmin_index:hmax_index], axis=0)
1646 ccf = numpy.average(dataOut.data_cspc[pairsIndexList[i], :, hmin_index:hmax_index], axis=0)
1646 powa = numpy.average(dataOut.data_spc[pair[0], :, hmin_index:hmax_index], axis=0)
1647 powa = numpy.average(dataOut.data_spc[pair[0], :, hmin_index:hmax_index], axis=0)
1647 powb = numpy.average(dataOut.data_spc[pair[1], :, hmin_index:hmax_index], axis=0)
1648 powb = numpy.average(dataOut.data_spc[pair[1], :, hmin_index:hmax_index], axis=0)
1648 avgcoherenceComplex = ccf/numpy.sqrt(powa*powb)
1649 avgcoherenceComplex = ccf/numpy.sqrt(powa*powb)
1649 phase = numpy.arctan2(avgcoherenceComplex.imag, avgcoherenceComplex.real)*180/numpy.pi
1650 phase = numpy.arctan2(avgcoherenceComplex.imag, avgcoherenceComplex.real)*180/numpy.pi
1650
1651
1651 #print "Phase %d%d" %(pair[0], pair[1])
1652 #print "Phase %d%d" %(pair[0], pair[1])
1652 #print phase[dataOut.beacon_heiIndexList]
1653 #print phase[dataOut.beacon_heiIndexList]
1653
1654
1654 if dataOut.beacon_heiIndexList:
1655 if dataOut.beacon_heiIndexList:
1655 phase_beacon[i] = numpy.average(phase[dataOut.beacon_heiIndexList])
1656 phase_beacon[i] = numpy.average(phase[dataOut.beacon_heiIndexList])
1656 else:
1657 else:
1657 phase_beacon[i] = numpy.average(phase)
1658 phase_beacon[i] = numpy.average(phase)
1658
1659
1659 if not self.isConfig:
1660 if not self.isConfig:
1660
1661
1661 nplots = len(pairsIndexList)
1662 nplots = len(pairsIndexList)
1662
1663
1663 self.setup(id=id,
1664 self.setup(id=id,
1664 nplots=nplots,
1665 nplots=nplots,
1665 wintitle=wintitle,
1666 wintitle=wintitle,
1666 showprofile=showprofile,
1667 showprofile=showprofile,
1667 show=show)
1668 show=show)
1668
1669
1669 if timerange != None:
1670 if timerange != None:
1670 self.timerange = timerange
1671 self.timerange = timerange
1671
1672
1672 self.xmin, self.xmax = self.getTimeLim(x, xmin, xmax, timerange)
1673 self.xmin, self.xmax = self.getTimeLim(x, xmin, xmax, timerange)
1673
1674
1674 if ymin == None: ymin = 0
1675 if ymin == None: ymin = 0
1675 if ymax == None: ymax = 360
1676 if ymax == None: ymax = 360
1676
1677
1677 self.FTP_WEI = ftp_wei
1678 self.FTP_WEI = ftp_wei
1678 self.EXP_CODE = exp_code
1679 self.EXP_CODE = exp_code
1679 self.SUB_EXP_CODE = sub_exp_code
1680 self.SUB_EXP_CODE = sub_exp_code
1680 self.PLOT_POS = plot_pos
1681 self.PLOT_POS = plot_pos
1681
1682
1682 self.name = thisDatetime.strftime("%Y%m%d_%H%M%S")
1683 self.name = thisDatetime.strftime("%Y%m%d_%H%M%S")
1683 self.isConfig = True
1684 self.isConfig = True
1684 self.figfile = figfile
1685 self.figfile = figfile
1685 self.xdata = numpy.array([])
1686 self.xdata = numpy.array([])
1686 self.ydata = numpy.array([])
1687 self.ydata = numpy.array([])
1687
1688
1688 update_figfile = True
1689 update_figfile = True
1689
1690
1690 #open file beacon phase
1691 #open file beacon phase
1691 path = '%s%03d' %(self.PREFIX, self.id)
1692 path = '%s%03d' %(self.PREFIX, self.id)
1692 beacon_file = os.path.join(path,'%s.txt'%self.name)
1693 beacon_file = os.path.join(path,'%s.txt'%self.name)
1693 self.filename_phase = os.path.join(figpath,beacon_file)
1694 self.filename_phase = os.path.join(figpath,beacon_file)
1694 #self.save_phase(self.filename_phase)
1695 #self.save_phase(self.filename_phase)
1695
1696
1696
1697
1697 #store data beacon phase
1698 #store data beacon phase
1698 #self.save_data(self.filename_phase, phase_beacon, thisDatetime)
1699 #self.save_data(self.filename_phase, phase_beacon, thisDatetime)
1699
1700
1700 self.setWinTitle(title)
1701 self.setWinTitle(title)
1701
1702
1702
1703
1703 title = "Phase Plot %s" %(thisDatetime.strftime("%Y/%m/%d %H:%M:%S"))
1704 title = "Phase Plot %s" %(thisDatetime.strftime("%Y/%m/%d %H:%M:%S"))
1704
1705
1705 legendlabels = ["Pair (%d,%d)"%(pair[0], pair[1]) for pair in dataOut.pairsList]
1706 legendlabels = ["Pair (%d,%d)"%(pair[0], pair[1]) for pair in dataOut.pairsList]
1706
1707
1707 axes = self.axesList[0]
1708 axes = self.axesList[0]
1708
1709
1709 self.xdata = numpy.hstack((self.xdata, x[0:1]))
1710 self.xdata = numpy.hstack((self.xdata, x[0:1]))
1710
1711
1711 if len(self.ydata)==0:
1712 if len(self.ydata)==0:
1712 self.ydata = phase_beacon.reshape(-1,1)
1713 self.ydata = phase_beacon.reshape(-1,1)
1713 else:
1714 else:
1714 self.ydata = numpy.hstack((self.ydata, phase_beacon.reshape(-1,1)))
1715 self.ydata = numpy.hstack((self.ydata, phase_beacon.reshape(-1,1)))
1715
1716
1716
1717
1717 axes.pmultilineyaxis(x=self.xdata, y=self.ydata,
1718 axes.pmultilineyaxis(x=self.xdata, y=self.ydata,
1718 xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax,
1719 xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax,
1719 xlabel=xlabel, ylabel=ylabel, title=title, legendlabels=legendlabels, marker='x', markersize=8, linestyle="solid",
1720 xlabel=xlabel, ylabel=ylabel, title=title, legendlabels=legendlabels, marker='x', markersize=8, linestyle="solid",
1720 XAxisAsTime=True, grid='both'
1721 XAxisAsTime=True, grid='both'
1721 )
1722 )
1722
1723
1723 self.draw()
1724 self.draw()
1724
1725
1725 if dataOut.ltctime >= self.xmax:
1726 if dataOut.ltctime >= self.xmax:
1726 self.counter_imagwr = wr_period
1727 self.counter_imagwr = wr_period
1727 self.isConfig = False
1728 self.isConfig = False
1728 update_figfile = True
1729 update_figfile = True
1729
1730
1730 self.save(figpath=figpath,
1731 self.save(figpath=figpath,
1731 figfile=figfile,
1732 figfile=figfile,
1732 save=save,
1733 save=save,
1733 ftp=ftp,
1734 ftp=ftp,
1734 wr_period=wr_period,
1735 wr_period=wr_period,
1735 thisDatetime=thisDatetime,
1736 thisDatetime=thisDatetime,
1736 update_figfile=update_figfile)
1737 update_figfile=update_figfile)
@@ -1,793 +1,797
1 import numpy
1 import numpy
2
2
3 from jroproc_base import ProcessingUnit, Operation
3 from jroproc_base import ProcessingUnit, Operation
4 from schainpy.model.data.jrodata import Spectra
4 from schainpy.model.data.jrodata import Spectra
5 from schainpy.model.data.jrodata import hildebrand_sekhon
5 from schainpy.model.data.jrodata import hildebrand_sekhon
6
6
7 class SpectraAFCProc(ProcessingUnit):
7 class SpectraAFCProc(ProcessingUnit):
8
8
9 def __init__(self, **kwargs):
9 def __init__(self, **kwargs):
10
10
11 ProcessingUnit.__init__(self, **kwargs)
11 ProcessingUnit.__init__(self, **kwargs)
12
12
13 self.buffer = None
13 self.buffer = None
14 self.firstdatatime = None
14 self.firstdatatime = None
15 self.profIndex = 0
15 self.profIndex = 0
16 self.dataOut = Spectra()
16 self.dataOut = Spectra()
17 self.id_min = None
17 self.id_min = None
18 self.id_max = None
18 self.id_max = None
19
19
20 def __updateSpecFromVoltage(self):
20 def __updateSpecFromVoltage(self):
21
21
22 self.dataOut.plotting = "spectra_acf"
22 self.dataOut.plotting = "spectra_acf"
23
23
24 self.dataOut.timeZone = self.dataIn.timeZone
24 self.dataOut.timeZone = self.dataIn.timeZone
25 self.dataOut.dstFlag = self.dataIn.dstFlag
25 self.dataOut.dstFlag = self.dataIn.dstFlag
26 self.dataOut.errorCount = self.dataIn.errorCount
26 self.dataOut.errorCount = self.dataIn.errorCount
27 self.dataOut.useLocalTime = self.dataIn.useLocalTime
27 self.dataOut.useLocalTime = self.dataIn.useLocalTime
28
28
29 self.dataOut.radarControllerHeaderObj = self.dataIn.radarControllerHeaderObj.copy()
29 self.dataOut.radarControllerHeaderObj = self.dataIn.radarControllerHeaderObj.copy()
30 self.dataOut.systemHeaderObj = self.dataIn.systemHeaderObj.copy()
30 self.dataOut.systemHeaderObj = self.dataIn.systemHeaderObj.copy()
31 self.dataOut.ippSeconds = self.dataIn.getDeltaH()*(10**-6)/0.15
31 self.dataOut.ippSeconds = self.dataIn.getDeltaH()*(10**-6)/0.15
32
32
33 self.dataOut.channelList = self.dataIn.channelList
33 self.dataOut.channelList = self.dataIn.channelList
34 self.dataOut.heightList = self.dataIn.heightList
34 self.dataOut.heightList = self.dataIn.heightList
35 self.dataOut.dtype = numpy.dtype([('real','<f4'),('imag','<f4')])
35 self.dataOut.dtype = numpy.dtype([('real','<f4'),('imag','<f4')])
36
36
37 self.dataOut.nBaud = self.dataIn.nBaud
37 self.dataOut.nBaud = self.dataIn.nBaud
38 self.dataOut.nCode = self.dataIn.nCode
38 self.dataOut.nCode = self.dataIn.nCode
39 self.dataOut.code = self.dataIn.code
39 self.dataOut.code = self.dataIn.code
40 # self.dataOut.nProfiles = self.dataOut.nFFTPoints
40 # self.dataOut.nProfiles = self.dataOut.nFFTPoints
41
41
42 self.dataOut.flagDiscontinuousBlock = self.dataIn.flagDiscontinuousBlock
42 self.dataOut.flagDiscontinuousBlock = self.dataIn.flagDiscontinuousBlock
43 self.dataOut.utctime = self.firstdatatime
43 self.dataOut.utctime = self.firstdatatime
44 self.dataOut.flagDecodeData = self.dataIn.flagDecodeData #asumo q la data esta decodificada
44 self.dataOut.flagDecodeData = self.dataIn.flagDecodeData #asumo q la data esta decodificada
45 self.dataOut.flagDeflipData = self.dataIn.flagDeflipData #asumo q la data esta sin flip
45 self.dataOut.flagDeflipData = self.dataIn.flagDeflipData #asumo q la data esta sin flip
46 self.dataOut.flagShiftFFT = False
46 self.dataOut.flagShiftFFT = False
47
47
48 self.dataOut.nCohInt = self.dataIn.nCohInt
48 self.dataOut.nCohInt = self.dataIn.nCohInt
49 self.dataOut.nIncohInt = 1
49 self.dataOut.nIncohInt = 1
50
50
51 self.dataOut.windowOfFilter = self.dataIn.windowOfFilter
51 self.dataOut.windowOfFilter = self.dataIn.windowOfFilter
52
52
53 self.dataOut.frequency = self.dataIn.frequency
53 self.dataOut.frequency = self.dataIn.frequency
54 self.dataOut.realtime = self.dataIn.realtime
54 self.dataOut.realtime = self.dataIn.realtime
55
55
56 self.dataOut.azimuth = self.dataIn.azimuth
56 self.dataOut.azimuth = self.dataIn.azimuth
57 self.dataOut.zenith = self.dataIn.zenith
57 self.dataOut.zenith = self.dataIn.zenith
58
58
59 self.dataOut.beam.codeList = self.dataIn.beam.codeList
59 self.dataOut.beam.codeList = self.dataIn.beam.codeList
60 self.dataOut.beam.azimuthList = self.dataIn.beam.azimuthList
60 self.dataOut.beam.azimuthList = self.dataIn.beam.azimuthList
61 self.dataOut.beam.zenithList = self.dataIn.beam.zenithList
61 self.dataOut.beam.zenithList = self.dataIn.beam.zenithList
62
62
63 def __decodeData(self, nProfiles, code):
63 def __decodeData(self, nProfiles, code):
64
64
65 if code is None:
65 if code is None:
66 return
66 return
67
67
68 for i in range(nProfiles):
68 for i in range(nProfiles):
69 self.buffer[:,i,:] = self.buffer[:,i,:]*code[0][i]
69 self.buffer[:,i,:] = self.buffer[:,i,:]*code[0][i]
70
70
71 def __getFft(self):
71 def __getFft(self):
72 """
72 """
73 Convierte valores de Voltaje a Spectra
73 Convierte valores de Voltaje a Spectra
74
74
75 Affected:
75 Affected:
76 self.dataOut.data_spc
76 self.dataOut.data_spc
77 self.dataOut.data_cspc
77 self.dataOut.data_cspc
78 self.dataOut.data_dc
78 self.dataOut.data_dc
79 self.dataOut.heightList
79 self.dataOut.heightList
80 self.profIndex
80 self.profIndex
81 self.buffer
81 self.buffer
82 self.dataOut.flagNoData
82 self.dataOut.flagNoData
83 """
83 """
84 nsegments = self.dataOut.nHeights
84 nsegments = self.dataOut.nHeights
85
85
86 _fft_buffer = numpy.zeros((self.dataOut.nChannels, self.dataOut.nProfiles, nsegments), dtype='complex')
86 _fft_buffer = numpy.zeros((self.dataOut.nChannels, self.dataOut.nProfiles, nsegments), dtype='complex')
87
87
88 for i in range(nsegments):
88 for i in range(nsegments):
89 try:
89 try:
90 _fft_buffer[:,:,i] = self.buffer[:,i:i+self.dataOut.nProfiles]
90 _fft_buffer[:,:,i] = self.buffer[:,i:i+self.dataOut.nProfiles]
91
91
92 if self.code is not None:
92 if self.code is not None:
93 _fft_buffer[:,:,i] = _fft_buffer[:,:,i]*self.code[0]
93 _fft_buffer[:,:,i] = _fft_buffer[:,:,i]*self.code[0]
94 except:
94 except:
95 pass
95 pass
96
96
97 fft_volt = numpy.fft.fft(_fft_buffer, n=self.dataOut.nFFTPoints, axis=1)
97 fft_volt = numpy.fft.fft(_fft_buffer, n=self.dataOut.nFFTPoints, axis=1)
98 fft_volt = fft_volt.astype(numpy.dtype('complex'))
98 fft_volt = fft_volt.astype(numpy.dtype('complex'))
99 dc = fft_volt[:,0,:]
99 dc = fft_volt[:,0,:]
100
100
101 #calculo de self-spectra
101 #calculo de self-spectra
102 spc = fft_volt * numpy.conjugate(fft_volt)
102 spc = fft_volt * numpy.conjugate(fft_volt)
103 data = numpy.fft.ifft(spc, axis=1)
103 data = numpy.fft.ifft(spc, axis=1)
104 data = numpy.fft.fftshift(data, axes=(1,))
104 data = numpy.fft.fftshift(data, axes=(1,))
105
105
106 spc = data
106 spc = data
107
107
108 blocksize = 0
108 blocksize = 0
109 blocksize += dc.size
109 blocksize += dc.size
110 blocksize += spc.size
110 blocksize += spc.size
111
111
112 cspc = None
112 cspc = None
113 pairIndex = 0
113 pairIndex = 0
114
114
115 if self.dataOut.pairsList != None:
115 if self.dataOut.pairsList != None:
116 #calculo de cross-spectra
116 #calculo de cross-spectra
117 cspc = numpy.zeros((self.dataOut.nPairs, self.dataOut.nFFTPoints, self.dataOut.nHeights), dtype='complex')
117 cspc = numpy.zeros((self.dataOut.nPairs, self.dataOut.nFFTPoints, self.dataOut.nHeights), dtype='complex')
118 for pair in self.dataOut.pairsList:
118 for pair in self.dataOut.pairsList:
119 if pair[0] not in self.dataOut.channelList:
119 if pair[0] not in self.dataOut.channelList:
120 raise ValueError, "Error getting CrossSpectra: pair 0 of %s is not in channelList = %s" %(str(pair), str(self.dataOut.channelList))
120 raise ValueError, "Error getting CrossSpectra: pair 0 of %s is not in channelList = %s" %(str(pair), str(self.dataOut.channelList))
121 if pair[1] not in self.dataOut.channelList:
121 if pair[1] not in self.dataOut.channelList:
122 raise ValueError, "Error getting CrossSpectra: pair 1 of %s is not in channelList = %s" %(str(pair), str(self.dataOut.channelList))
122 raise ValueError, "Error getting CrossSpectra: pair 1 of %s is not in channelList = %s" %(str(pair), str(self.dataOut.channelList))
123
123
124 chan_index0 = self.dataOut.channelList.index(pair[0])
124 chan_index0 = self.dataOut.channelList.index(pair[0])
125 chan_index1 = self.dataOut.channelList.index(pair[1])
125 chan_index1 = self.dataOut.channelList.index(pair[1])
126
126
127 cspc[pairIndex,:,:] = fft_volt[chan_index0,:,:] * numpy.conjugate(fft_volt[chan_index1,:,:])
127 cspc[pairIndex,:,:] = fft_volt[chan_index0,:,:] * numpy.conjugate(fft_volt[chan_index1,:,:])
128 pairIndex += 1
128 pairIndex += 1
129 blocksize += cspc.size
129 blocksize += cspc.size
130
130
131 self.dataOut.data_spc = spc
131 self.dataOut.data_spc = spc
132 self.dataOut.data_cspc = cspc
132 self.dataOut.data_cspc = cspc
133 self.dataOut.data_dc = dc
133 self.dataOut.data_dc = dc
134 self.dataOut.blockSize = blocksize
134 self.dataOut.blockSize = blocksize
135 self.dataOut.flagShiftFFT = True
135 self.dataOut.flagShiftFFT = True
136
136
137 def run(self, nProfiles=None, nFFTPoints=None, pairsList=[], code=None, nCode=1, nBaud=1,real= None, imag=None):
137 def run(self, nProfiles=None, nFFTPoints=None, pairsList=[], code=None, nCode=1, nBaud=1,real= None, imag=None):
138
138
139 self.dataOut.flagNoData = True
139 self.dataOut.flagNoData = True
140
140
141 if self.dataIn.type == "Spectra":
141 if self.dataIn.type == "Spectra":
142 self.dataOut.copy(self.dataIn)
142 self.dataOut.copy(self.dataIn)
143 spc = self.dataOut.data_spc
143 spc = self.dataOut.data_spc
144 data = numpy.fft.ifft(spc, axis=1)
144 data = numpy.fft.ifft(spc, axis=1)
145 data = numpy.fft.fftshift(data, axes=(1,))
145 data = numpy.fft.fftshift(data, axes=(1,))
146 acf = numpy.abs(data) # Autocorrelacion LLAMAR A ESTE VALOR ACF
146 acf = numpy.abs(data) # Autocorrelacion LLAMAR A ESTE VALOR ACF
147 if real:
147 if real:
148 acf = data.real
148 acf = data.real
149 if imag:
149 if imag:
150 acf = data.imag
150 acf = data.imag
151 shape = acf.shape # nchannels, nprofiles, nsamples
151 shape = acf.shape # nchannels, nprofiles, nsamples
152
152
153 #import matplotlib.pyplot as plt
153 #import matplotlib.pyplot as plt
154 #acf_tmp=acf[0,:,85]
154 #print "test",acf.shape
155 #acf_tmp=acf[0,:,]
155 #plt.plot(acf_tmp)
156 #plt.plot(acf_tmp)
156 #plt.show()
157 #plt.show()
158 #import time
159 #time.sleep(10)
157
160
161 for j in range(shape[0]):
158 for i in range(shape[1]):
162 for i in range(shape[1]):
159 tmp = numpy.argmax(acf[0,:,i])
163 tmp = numpy.argmax(acf[j,:,i])
160 if i>30:
164 if i>30:
161 value = (acf[0,:,i][tmp+3]+acf[0,:,i][tmp+4])/2.0
165 value = (acf[j,:,i][tmp+3]+acf[j,:,i][tmp+4])/2.0
162 acf[0,:,i][tmp] = value
166 acf[j,:,i][tmp] = value
163 acf[0,:,i][tmp-1] = value
167 acf[j,:,i][tmp-1] = value
164 acf[0,:,i][tmp+1] = value
168 acf[j,:,i][tmp+1] = value
165 acf[0,:,i][tmp-2] = value
169 acf[j,:,i][tmp-2] = value
166 acf[0,:,i][tmp+2] = value
170 acf[j,:,i][tmp+2] = value
167
171
168 import scipy as sp
172 import scipy as sp
169 from scipy import signal
173 from scipy import signal
170 acf[0,:,i] = sp.signal.medfilt(acf[0,:,i],21)
174 #acf[3,:,i] = sp.signal.medfilt(acf[3,:,i],21)
171
175
172
176
173
177
174 #print numpy.argmax(acf[0,:,85])
178 #print numpy.argmax(acf[0,:,85])
175 #import matplotlib.pyplot as plt
179 #import matplotlib.pyplot as plt
176 #plt.plot(acf[0,:,85])
180 #plt.plot(acf[0,:,85])
177 #a= acf[0,:,85]
181 #a= acf[0,:,85]
178 #b= acf[0,:,0]
182 #b= acf[0,:,0]
179 #print a[200],a[198],a[199], a[201],a[202],a[203]
183 #print a[200],a[198],a[199], a[201],a[202],a[203]
180 #plt.show()
184 #plt.show()
181 #import time
185 #import time
182 #time.sleep(10)
186 #time.sleep(10)
183
187
184
188
185 # Normalizando
189 # Normalizando
186 for i in range(shape[0]):
190 for i in range(shape[0]):
187 for j in range(shape[2]):
191 for j in range(shape[2]):
188 acf[i,:,j]= acf[i,:,j] / numpy.max(numpy.abs(acf[i,:,j]))
192 acf[i,:,j]= acf[i,:,j] / numpy.max(numpy.abs(acf[i,:,j]))
189
193
190 #import matplotlib.pyplot as plt
194 #import matplotlib.pyplot as plt
191 #plt.plot(acf[0,:,85])
195 #plt.plot(acf[0,:,85])
192 #plt.show()
196 #plt.show()
193 #import time
197 #import time
194 #time.sleep(20)
198 #time.sleep(20)
195
199
196
200
197 self.dataOut.data_acf = acf
201 self.dataOut.data_acf = acf
198 return True
202 return True
199
203
200 if code is not None:
204 if code is not None:
201 self.code = numpy.array(code).reshape(nCode,nBaud)
205 self.code = numpy.array(code).reshape(nCode,nBaud)
202 else:
206 else:
203 self.code = None
207 self.code = None
204
208
205 if self.dataIn.type == "Voltage":
209 if self.dataIn.type == "Voltage":
206
210
207 if nFFTPoints == None:
211 if nFFTPoints == None:
208 raise ValueError, "This SpectraProc.run() need nFFTPoints input variable"
212 raise ValueError, "This SpectraProc.run() need nFFTPoints input variable"
209
213
210 if nProfiles == None:
214 if nProfiles == None:
211 nProfiles = nFFTPoints
215 nProfiles = nFFTPoints
212
216
213 self.dataOut.ippFactor = 1
217 self.dataOut.ippFactor = 1
214
218
215 self.dataOut.nFFTPoints = nFFTPoints
219 self.dataOut.nFFTPoints = nFFTPoints
216 self.dataOut.nProfiles = nProfiles
220 self.dataOut.nProfiles = nProfiles
217 self.dataOut.pairsList = pairsList
221 self.dataOut.pairsList = pairsList
218
222
219 # if self.buffer is None:
223 # if self.buffer is None:
220 # self.buffer = numpy.zeros( (self.dataIn.nChannels, nProfiles, self.dataIn.nHeights),
224 # self.buffer = numpy.zeros( (self.dataIn.nChannels, nProfiles, self.dataIn.nHeights),
221 # dtype='complex')
225 # dtype='complex')
222
226
223 if not self.dataIn.flagDataAsBlock:
227 if not self.dataIn.flagDataAsBlock:
224 self.buffer = self.dataIn.data.copy()
228 self.buffer = self.dataIn.data.copy()
225
229
226 # for i in range(self.dataIn.nHeights):
230 # for i in range(self.dataIn.nHeights):
227 # self.buffer[:, self.profIndex, self.profIndex:] = voltage_data[:,:self.dataIn.nHeights - self.profIndex]
231 # self.buffer[:, self.profIndex, self.profIndex:] = voltage_data[:,:self.dataIn.nHeights - self.profIndex]
228 #
232 #
229 # self.profIndex += 1
233 # self.profIndex += 1
230
234
231 else:
235 else:
232 raise ValueError, ""
236 raise ValueError, ""
233
237
234 self.firstdatatime = self.dataIn.utctime
238 self.firstdatatime = self.dataIn.utctime
235
239
236 self.profIndex == nProfiles
240 self.profIndex == nProfiles
237
241
238 self.__updateSpecFromVoltage()
242 self.__updateSpecFromVoltage()
239
243
240 self.__getFft()
244 self.__getFft()
241
245
242 self.dataOut.flagNoData = False
246 self.dataOut.flagNoData = False
243
247
244 return True
248 return True
245
249
246 raise ValueError, "The type of input object '%s' is not valid"%(self.dataIn.type)
250 raise ValueError, "The type of input object '%s' is not valid"%(self.dataIn.type)
247
251
248 def __selectPairs(self, pairsList):
252 def __selectPairs(self, pairsList):
249
253
250 if channelList == None:
254 if channelList == None:
251 return
255 return
252
256
253 pairsIndexListSelected = []
257 pairsIndexListSelected = []
254
258
255 for thisPair in pairsList:
259 for thisPair in pairsList:
256
260
257 if thisPair not in self.dataOut.pairsList:
261 if thisPair not in self.dataOut.pairsList:
258 continue
262 continue
259
263
260 pairIndex = self.dataOut.pairsList.index(thisPair)
264 pairIndex = self.dataOut.pairsList.index(thisPair)
261
265
262 pairsIndexListSelected.append(pairIndex)
266 pairsIndexListSelected.append(pairIndex)
263
267
264 if not pairsIndexListSelected:
268 if not pairsIndexListSelected:
265 self.dataOut.data_cspc = None
269 self.dataOut.data_cspc = None
266 self.dataOut.pairsList = []
270 self.dataOut.pairsList = []
267 return
271 return
268
272
269 self.dataOut.data_cspc = self.dataOut.data_cspc[pairsIndexListSelected]
273 self.dataOut.data_cspc = self.dataOut.data_cspc[pairsIndexListSelected]
270 self.dataOut.pairsList = [self.dataOut.pairsList[i] for i in pairsIndexListSelected]
274 self.dataOut.pairsList = [self.dataOut.pairsList[i] for i in pairsIndexListSelected]
271
275
272 return
276 return
273
277
274 def __selectPairsByChannel(self, channelList=None):
278 def __selectPairsByChannel(self, channelList=None):
275
279
276 if channelList == None:
280 if channelList == None:
277 return
281 return
278
282
279 pairsIndexListSelected = []
283 pairsIndexListSelected = []
280 for pairIndex in self.dataOut.pairsIndexList:
284 for pairIndex in self.dataOut.pairsIndexList:
281 #First pair
285 #First pair
282 if self.dataOut.pairsList[pairIndex][0] not in channelList:
286 if self.dataOut.pairsList[pairIndex][0] not in channelList:
283 continue
287 continue
284 #Second pair
288 #Second pair
285 if self.dataOut.pairsList[pairIndex][1] not in channelList:
289 if self.dataOut.pairsList[pairIndex][1] not in channelList:
286 continue
290 continue
287
291
288 pairsIndexListSelected.append(pairIndex)
292 pairsIndexListSelected.append(pairIndex)
289
293
290 if not pairsIndexListSelected:
294 if not pairsIndexListSelected:
291 self.dataOut.data_cspc = None
295 self.dataOut.data_cspc = None
292 self.dataOut.pairsList = []
296 self.dataOut.pairsList = []
293 return
297 return
294
298
295 self.dataOut.data_cspc = self.dataOut.data_cspc[pairsIndexListSelected]
299 self.dataOut.data_cspc = self.dataOut.data_cspc[pairsIndexListSelected]
296 self.dataOut.pairsList = [self.dataOut.pairsList[i] for i in pairsIndexListSelected]
300 self.dataOut.pairsList = [self.dataOut.pairsList[i] for i in pairsIndexListSelected]
297
301
298 return
302 return
299
303
300 def selectChannels(self, channelList):
304 def selectChannels(self, channelList):
301
305
302 channelIndexList = []
306 channelIndexList = []
303
307
304 for channel in channelList:
308 for channel in channelList:
305 if channel not in self.dataOut.channelList:
309 if channel not in self.dataOut.channelList:
306 raise ValueError, "Error selecting channels, Channel %d is not valid.\nAvailable channels = %s" %(channel, str(self.dataOut.channelList))
310 raise ValueError, "Error selecting channels, Channel %d is not valid.\nAvailable channels = %s" %(channel, str(self.dataOut.channelList))
307
311
308 index = self.dataOut.channelList.index(channel)
312 index = self.dataOut.channelList.index(channel)
309 channelIndexList.append(index)
313 channelIndexList.append(index)
310
314
311 self.selectChannelsByIndex(channelIndexList)
315 self.selectChannelsByIndex(channelIndexList)
312
316
313 def selectChannelsByIndex(self, channelIndexList):
317 def selectChannelsByIndex(self, channelIndexList):
314 """
318 """
315 Selecciona un bloque de datos en base a canales segun el channelIndexList
319 Selecciona un bloque de datos en base a canales segun el channelIndexList
316
320
317 Input:
321 Input:
318 channelIndexList : lista sencilla de canales a seleccionar por ej. [2,3,7]
322 channelIndexList : lista sencilla de canales a seleccionar por ej. [2,3,7]
319
323
320 Affected:
324 Affected:
321 self.dataOut.data_spc
325 self.dataOut.data_spc
322 self.dataOut.channelIndexList
326 self.dataOut.channelIndexList
323 self.dataOut.nChannels
327 self.dataOut.nChannels
324
328
325 Return:
329 Return:
326 None
330 None
327 """
331 """
328
332
329 for channelIndex in channelIndexList:
333 for channelIndex in channelIndexList:
330 if channelIndex not in self.dataOut.channelIndexList:
334 if channelIndex not in self.dataOut.channelIndexList:
331 raise ValueError, "Error selecting channels: The value %d in channelIndexList is not valid.\nAvailable channel indexes = " %(channelIndex, self.dataOut.channelIndexList)
335 raise ValueError, "Error selecting channels: The value %d in channelIndexList is not valid.\nAvailable channel indexes = " %(channelIndex, self.dataOut.channelIndexList)
332
336
333 # nChannels = len(channelIndexList)
337 # nChannels = len(channelIndexList)
334
338
335 data_spc = self.dataOut.data_spc[channelIndexList,:]
339 data_spc = self.dataOut.data_spc[channelIndexList,:]
336 data_dc = self.dataOut.data_dc[channelIndexList,:]
340 data_dc = self.dataOut.data_dc[channelIndexList,:]
337
341
338 self.dataOut.data_spc = data_spc
342 self.dataOut.data_spc = data_spc
339 self.dataOut.data_dc = data_dc
343 self.dataOut.data_dc = data_dc
340
344
341 self.dataOut.channelList = [self.dataOut.channelList[i] for i in channelIndexList]
345 self.dataOut.channelList = [self.dataOut.channelList[i] for i in channelIndexList]
342 # self.dataOut.nChannels = nChannels
346 # self.dataOut.nChannels = nChannels
343
347
344 self.__selectPairsByChannel(self.dataOut.channelList)
348 self.__selectPairsByChannel(self.dataOut.channelList)
345
349
346 return 1
350 return 1
347
351
348 def selectHeights(self, minHei, maxHei):
352 def selectHeights(self, minHei, maxHei):
349 """
353 """
350 Selecciona un bloque de datos en base a un grupo de valores de alturas segun el rango
354 Selecciona un bloque de datos en base a un grupo de valores de alturas segun el rango
351 minHei <= height <= maxHei
355 minHei <= height <= maxHei
352
356
353 Input:
357 Input:
354 minHei : valor minimo de altura a considerar
358 minHei : valor minimo de altura a considerar
355 maxHei : valor maximo de altura a considerar
359 maxHei : valor maximo de altura a considerar
356
360
357 Affected:
361 Affected:
358 Indirectamente son cambiados varios valores a travez del metodo selectHeightsByIndex
362 Indirectamente son cambiados varios valores a travez del metodo selectHeightsByIndex
359
363
360 Return:
364 Return:
361 1 si el metodo se ejecuto con exito caso contrario devuelve 0
365 1 si el metodo se ejecuto con exito caso contrario devuelve 0
362 """
366 """
363
367
364 if (minHei > maxHei):
368 if (minHei > maxHei):
365 raise ValueError, "Error selecting heights: Height range (%d,%d) is not valid" % (minHei, maxHei)
369 raise ValueError, "Error selecting heights: Height range (%d,%d) is not valid" % (minHei, maxHei)
366
370
367 if (minHei < self.dataOut.heightList[0]):
371 if (minHei < self.dataOut.heightList[0]):
368 minHei = self.dataOut.heightList[0]
372 minHei = self.dataOut.heightList[0]
369
373
370 if (maxHei > self.dataOut.heightList[-1]):
374 if (maxHei > self.dataOut.heightList[-1]):
371 maxHei = self.dataOut.heightList[-1]
375 maxHei = self.dataOut.heightList[-1]
372
376
373 minIndex = 0
377 minIndex = 0
374 maxIndex = 0
378 maxIndex = 0
375 heights = self.dataOut.heightList
379 heights = self.dataOut.heightList
376
380
377 inda = numpy.where(heights >= minHei)
381 inda = numpy.where(heights >= minHei)
378 indb = numpy.where(heights <= maxHei)
382 indb = numpy.where(heights <= maxHei)
379
383
380 try:
384 try:
381 minIndex = inda[0][0]
385 minIndex = inda[0][0]
382 except:
386 except:
383 minIndex = 0
387 minIndex = 0
384
388
385 try:
389 try:
386 maxIndex = indb[0][-1]
390 maxIndex = indb[0][-1]
387 except:
391 except:
388 maxIndex = len(heights)
392 maxIndex = len(heights)
389
393
390 self.selectHeightsByIndex(minIndex, maxIndex)
394 self.selectHeightsByIndex(minIndex, maxIndex)
391
395
392 return 1
396 return 1
393
397
394 def getBeaconSignal(self, tauindex = 0, channelindex = 0, hei_ref=None):
398 def getBeaconSignal(self, tauindex = 0, channelindex = 0, hei_ref=None):
395 newheis = numpy.where(self.dataOut.heightList>self.dataOut.radarControllerHeaderObj.Taus[tauindex])
399 newheis = numpy.where(self.dataOut.heightList>self.dataOut.radarControllerHeaderObj.Taus[tauindex])
396
400
397 if hei_ref != None:
401 if hei_ref != None:
398 newheis = numpy.where(self.dataOut.heightList>hei_ref)
402 newheis = numpy.where(self.dataOut.heightList>hei_ref)
399
403
400 minIndex = min(newheis[0])
404 minIndex = min(newheis[0])
401 maxIndex = max(newheis[0])
405 maxIndex = max(newheis[0])
402 data_spc = self.dataOut.data_spc[:,:,minIndex:maxIndex+1]
406 data_spc = self.dataOut.data_spc[:,:,minIndex:maxIndex+1]
403 heightList = self.dataOut.heightList[minIndex:maxIndex+1]
407 heightList = self.dataOut.heightList[minIndex:maxIndex+1]
404
408
405 # determina indices
409 # determina indices
406 nheis = int(self.dataOut.radarControllerHeaderObj.txB/(self.dataOut.heightList[1]-self.dataOut.heightList[0]))
410 nheis = int(self.dataOut.radarControllerHeaderObj.txB/(self.dataOut.heightList[1]-self.dataOut.heightList[0]))
407 avg_dB = 10*numpy.log10(numpy.sum(data_spc[channelindex,:,:],axis=0))
411 avg_dB = 10*numpy.log10(numpy.sum(data_spc[channelindex,:,:],axis=0))
408 beacon_dB = numpy.sort(avg_dB)[-nheis:]
412 beacon_dB = numpy.sort(avg_dB)[-nheis:]
409 beacon_heiIndexList = []
413 beacon_heiIndexList = []
410 for val in avg_dB.tolist():
414 for val in avg_dB.tolist():
411 if val >= beacon_dB[0]:
415 if val >= beacon_dB[0]:
412 beacon_heiIndexList.append(avg_dB.tolist().index(val))
416 beacon_heiIndexList.append(avg_dB.tolist().index(val))
413
417
414 #data_spc = data_spc[:,:,beacon_heiIndexList]
418 #data_spc = data_spc[:,:,beacon_heiIndexList]
415 data_cspc = None
419 data_cspc = None
416 if self.dataOut.data_cspc is not None:
420 if self.dataOut.data_cspc is not None:
417 data_cspc = self.dataOut.data_cspc[:,:,minIndex:maxIndex+1]
421 data_cspc = self.dataOut.data_cspc[:,:,minIndex:maxIndex+1]
418 #data_cspc = data_cspc[:,:,beacon_heiIndexList]
422 #data_cspc = data_cspc[:,:,beacon_heiIndexList]
419
423
420 data_dc = None
424 data_dc = None
421 if self.dataOut.data_dc is not None:
425 if self.dataOut.data_dc is not None:
422 data_dc = self.dataOut.data_dc[:,minIndex:maxIndex+1]
426 data_dc = self.dataOut.data_dc[:,minIndex:maxIndex+1]
423 #data_dc = data_dc[:,beacon_heiIndexList]
427 #data_dc = data_dc[:,beacon_heiIndexList]
424
428
425 self.dataOut.data_spc = data_spc
429 self.dataOut.data_spc = data_spc
426 self.dataOut.data_cspc = data_cspc
430 self.dataOut.data_cspc = data_cspc
427 self.dataOut.data_dc = data_dc
431 self.dataOut.data_dc = data_dc
428 self.dataOut.heightList = heightList
432 self.dataOut.heightList = heightList
429 self.dataOut.beacon_heiIndexList = beacon_heiIndexList
433 self.dataOut.beacon_heiIndexList = beacon_heiIndexList
430
434
431 return 1
435 return 1
432
436
433
437
434 def selectHeightsByIndex(self, minIndex, maxIndex):
438 def selectHeightsByIndex(self, minIndex, maxIndex):
435 """
439 """
436 Selecciona un bloque de datos en base a un grupo indices de alturas segun el rango
440 Selecciona un bloque de datos en base a un grupo indices de alturas segun el rango
437 minIndex <= index <= maxIndex
441 minIndex <= index <= maxIndex
438
442
439 Input:
443 Input:
440 minIndex : valor de indice minimo de altura a considerar
444 minIndex : valor de indice minimo de altura a considerar
441 maxIndex : valor de indice maximo de altura a considerar
445 maxIndex : valor de indice maximo de altura a considerar
442
446
443 Affected:
447 Affected:
444 self.dataOut.data_spc
448 self.dataOut.data_spc
445 self.dataOut.data_cspc
449 self.dataOut.data_cspc
446 self.dataOut.data_dc
450 self.dataOut.data_dc
447 self.dataOut.heightList
451 self.dataOut.heightList
448
452
449 Return:
453 Return:
450 1 si el metodo se ejecuto con exito caso contrario devuelve 0
454 1 si el metodo se ejecuto con exito caso contrario devuelve 0
451 """
455 """
452
456
453 if (minIndex < 0) or (minIndex > maxIndex):
457 if (minIndex < 0) or (minIndex > maxIndex):
454 raise ValueError, "Error selecting heights: Index range (%d,%d) is not valid" % (minIndex, maxIndex)
458 raise ValueError, "Error selecting heights: Index range (%d,%d) is not valid" % (minIndex, maxIndex)
455
459
456 if (maxIndex >= self.dataOut.nHeights):
460 if (maxIndex >= self.dataOut.nHeights):
457 maxIndex = self.dataOut.nHeights-1
461 maxIndex = self.dataOut.nHeights-1
458
462
459 #Spectra
463 #Spectra
460 data_spc = self.dataOut.data_spc[:,:,minIndex:maxIndex+1]
464 data_spc = self.dataOut.data_spc[:,:,minIndex:maxIndex+1]
461
465
462 data_cspc = None
466 data_cspc = None
463 if self.dataOut.data_cspc is not None:
467 if self.dataOut.data_cspc is not None:
464 data_cspc = self.dataOut.data_cspc[:,:,minIndex:maxIndex+1]
468 data_cspc = self.dataOut.data_cspc[:,:,minIndex:maxIndex+1]
465
469
466 data_dc = None
470 data_dc = None
467 if self.dataOut.data_dc is not None:
471 if self.dataOut.data_dc is not None:
468 data_dc = self.dataOut.data_dc[:,minIndex:maxIndex+1]
472 data_dc = self.dataOut.data_dc[:,minIndex:maxIndex+1]
469
473
470 self.dataOut.data_spc = data_spc
474 self.dataOut.data_spc = data_spc
471 self.dataOut.data_cspc = data_cspc
475 self.dataOut.data_cspc = data_cspc
472 self.dataOut.data_dc = data_dc
476 self.dataOut.data_dc = data_dc
473
477
474 self.dataOut.heightList = self.dataOut.heightList[minIndex:maxIndex+1]
478 self.dataOut.heightList = self.dataOut.heightList[minIndex:maxIndex+1]
475
479
476 return 1
480 return 1
477
481
478 def removeDC(self, mode = 2):
482 def removeDC(self, mode = 2):
479 jspectra = self.dataOut.data_spc
483 jspectra = self.dataOut.data_spc
480 jcspectra = self.dataOut.data_cspc
484 jcspectra = self.dataOut.data_cspc
481
485
482
486
483 num_chan = jspectra.shape[0]
487 num_chan = jspectra.shape[0]
484 num_hei = jspectra.shape[2]
488 num_hei = jspectra.shape[2]
485
489
486 if jcspectra is not None:
490 if jcspectra is not None:
487 jcspectraExist = True
491 jcspectraExist = True
488 num_pairs = jcspectra.shape[0]
492 num_pairs = jcspectra.shape[0]
489 else: jcspectraExist = False
493 else: jcspectraExist = False
490
494
491 freq_dc = jspectra.shape[1]/2
495 freq_dc = jspectra.shape[1]/2
492 ind_vel = numpy.array([-2,-1,1,2]) + freq_dc
496 ind_vel = numpy.array([-2,-1,1,2]) + freq_dc
493
497
494 if ind_vel[0]<0:
498 if ind_vel[0]<0:
495 ind_vel[range(0,1)] = ind_vel[range(0,1)] + self.num_prof
499 ind_vel[range(0,1)] = ind_vel[range(0,1)] + self.num_prof
496
500
497 if mode == 1:
501 if mode == 1:
498 jspectra[:,freq_dc,:] = (jspectra[:,ind_vel[1],:] + jspectra[:,ind_vel[2],:])/2 #CORRECCION
502 jspectra[:,freq_dc,:] = (jspectra[:,ind_vel[1],:] + jspectra[:,ind_vel[2],:])/2 #CORRECCION
499
503
500 if jcspectraExist:
504 if jcspectraExist:
501 jcspectra[:,freq_dc,:] = (jcspectra[:,ind_vel[1],:] + jcspectra[:,ind_vel[2],:])/2
505 jcspectra[:,freq_dc,:] = (jcspectra[:,ind_vel[1],:] + jcspectra[:,ind_vel[2],:])/2
502
506
503 if mode == 2:
507 if mode == 2:
504
508
505 vel = numpy.array([-2,-1,1,2])
509 vel = numpy.array([-2,-1,1,2])
506 xx = numpy.zeros([4,4])
510 xx = numpy.zeros([4,4])
507
511
508 for fil in range(4):
512 for fil in range(4):
509 xx[fil,:] = vel[fil]**numpy.asarray(range(4))
513 xx[fil,:] = vel[fil]**numpy.asarray(range(4))
510
514
511 xx_inv = numpy.linalg.inv(xx)
515 xx_inv = numpy.linalg.inv(xx)
512 xx_aux = xx_inv[0,:]
516 xx_aux = xx_inv[0,:]
513
517
514 for ich in range(num_chan):
518 for ich in range(num_chan):
515 yy = jspectra[ich,ind_vel,:]
519 yy = jspectra[ich,ind_vel,:]
516 jspectra[ich,freq_dc,:] = numpy.dot(xx_aux,yy)
520 jspectra[ich,freq_dc,:] = numpy.dot(xx_aux,yy)
517
521
518 junkid = jspectra[ich,freq_dc,:]<=0
522 junkid = jspectra[ich,freq_dc,:]<=0
519 cjunkid = sum(junkid)
523 cjunkid = sum(junkid)
520
524
521 if cjunkid.any():
525 if cjunkid.any():
522 jspectra[ich,freq_dc,junkid.nonzero()] = (jspectra[ich,ind_vel[1],junkid] + jspectra[ich,ind_vel[2],junkid])/2
526 jspectra[ich,freq_dc,junkid.nonzero()] = (jspectra[ich,ind_vel[1],junkid] + jspectra[ich,ind_vel[2],junkid])/2
523
527
524 if jcspectraExist:
528 if jcspectraExist:
525 for ip in range(num_pairs):
529 for ip in range(num_pairs):
526 yy = jcspectra[ip,ind_vel,:]
530 yy = jcspectra[ip,ind_vel,:]
527 jcspectra[ip,freq_dc,:] = numpy.dot(xx_aux,yy)
531 jcspectra[ip,freq_dc,:] = numpy.dot(xx_aux,yy)
528
532
529
533
530 self.dataOut.data_spc = jspectra
534 self.dataOut.data_spc = jspectra
531 self.dataOut.data_cspc = jcspectra
535 self.dataOut.data_cspc = jcspectra
532
536
533 return 1
537 return 1
534
538
535 def removeInterference(self, interf = 2,hei_interf = None, nhei_interf = None, offhei_interf = None):
539 def removeInterference(self, interf = 2,hei_interf = None, nhei_interf = None, offhei_interf = None):
536
540
537 jspectra = self.dataOut.data_spc
541 jspectra = self.dataOut.data_spc
538 jcspectra = self.dataOut.data_cspc
542 jcspectra = self.dataOut.data_cspc
539 jnoise = self.dataOut.getNoise()
543 jnoise = self.dataOut.getNoise()
540 num_incoh = self.dataOut.nIncohInt
544 num_incoh = self.dataOut.nIncohInt
541
545
542 num_channel = jspectra.shape[0]
546 num_channel = jspectra.shape[0]
543 num_prof = jspectra.shape[1]
547 num_prof = jspectra.shape[1]
544 num_hei = jspectra.shape[2]
548 num_hei = jspectra.shape[2]
545
549
546 #hei_interf
550 #hei_interf
547 if hei_interf is None:
551 if hei_interf is None:
548 count_hei = num_hei/2 #Como es entero no importa
552 count_hei = num_hei/2 #Como es entero no importa
549 hei_interf = numpy.asmatrix(range(count_hei)) + num_hei - count_hei
553 hei_interf = numpy.asmatrix(range(count_hei)) + num_hei - count_hei
550 hei_interf = numpy.asarray(hei_interf)[0]
554 hei_interf = numpy.asarray(hei_interf)[0]
551 #nhei_interf
555 #nhei_interf
552 if (nhei_interf == None):
556 if (nhei_interf == None):
553 nhei_interf = 5
557 nhei_interf = 5
554 if (nhei_interf < 1):
558 if (nhei_interf < 1):
555 nhei_interf = 1
559 nhei_interf = 1
556 if (nhei_interf > count_hei):
560 if (nhei_interf > count_hei):
557 nhei_interf = count_hei
561 nhei_interf = count_hei
558 if (offhei_interf == None):
562 if (offhei_interf == None):
559 offhei_interf = 0
563 offhei_interf = 0
560
564
561 ind_hei = range(num_hei)
565 ind_hei = range(num_hei)
562 # mask_prof = numpy.asarray(range(num_prof - 2)) + 1
566 # mask_prof = numpy.asarray(range(num_prof - 2)) + 1
563 # mask_prof[range(num_prof/2 - 1,len(mask_prof))] += 1
567 # mask_prof[range(num_prof/2 - 1,len(mask_prof))] += 1
564 mask_prof = numpy.asarray(range(num_prof))
568 mask_prof = numpy.asarray(range(num_prof))
565 num_mask_prof = mask_prof.size
569 num_mask_prof = mask_prof.size
566 comp_mask_prof = [0, num_prof/2]
570 comp_mask_prof = [0, num_prof/2]
567
571
568
572
569 #noise_exist: Determina si la variable jnoise ha sido definida y contiene la informacion del ruido de cada canal
573 #noise_exist: Determina si la variable jnoise ha sido definida y contiene la informacion del ruido de cada canal
570 if (jnoise.size < num_channel or numpy.isnan(jnoise).any()):
574 if (jnoise.size < num_channel or numpy.isnan(jnoise).any()):
571 jnoise = numpy.nan
575 jnoise = numpy.nan
572 noise_exist = jnoise[0] < numpy.Inf
576 noise_exist = jnoise[0] < numpy.Inf
573
577
574 #Subrutina de Remocion de la Interferencia
578 #Subrutina de Remocion de la Interferencia
575 for ich in range(num_channel):
579 for ich in range(num_channel):
576 #Se ordena los espectros segun su potencia (menor a mayor)
580 #Se ordena los espectros segun su potencia (menor a mayor)
577 power = jspectra[ich,mask_prof,:]
581 power = jspectra[ich,mask_prof,:]
578 power = power[:,hei_interf]
582 power = power[:,hei_interf]
579 power = power.sum(axis = 0)
583 power = power.sum(axis = 0)
580 psort = power.ravel().argsort()
584 psort = power.ravel().argsort()
581
585
582 #Se estima la interferencia promedio en los Espectros de Potencia empleando
586 #Se estima la interferencia promedio en los Espectros de Potencia empleando
583 junkspc_interf = jspectra[ich,:,hei_interf[psort[range(offhei_interf, nhei_interf + offhei_interf)]]]
587 junkspc_interf = jspectra[ich,:,hei_interf[psort[range(offhei_interf, nhei_interf + offhei_interf)]]]
584
588
585 if noise_exist:
589 if noise_exist:
586 # tmp_noise = jnoise[ich] / num_prof
590 # tmp_noise = jnoise[ich] / num_prof
587 tmp_noise = jnoise[ich]
591 tmp_noise = jnoise[ich]
588 junkspc_interf = junkspc_interf - tmp_noise
592 junkspc_interf = junkspc_interf - tmp_noise
589 #junkspc_interf[:,comp_mask_prof] = 0
593 #junkspc_interf[:,comp_mask_prof] = 0
590
594
591 jspc_interf = junkspc_interf.sum(axis = 0) / nhei_interf
595 jspc_interf = junkspc_interf.sum(axis = 0) / nhei_interf
592 jspc_interf = jspc_interf.transpose()
596 jspc_interf = jspc_interf.transpose()
593 #Calculando el espectro de interferencia promedio
597 #Calculando el espectro de interferencia promedio
594 noiseid = numpy.where(jspc_interf <= tmp_noise/ numpy.sqrt(num_incoh))
598 noiseid = numpy.where(jspc_interf <= tmp_noise/ numpy.sqrt(num_incoh))
595 noiseid = noiseid[0]
599 noiseid = noiseid[0]
596 cnoiseid = noiseid.size
600 cnoiseid = noiseid.size
597 interfid = numpy.where(jspc_interf > tmp_noise/ numpy.sqrt(num_incoh))
601 interfid = numpy.where(jspc_interf > tmp_noise/ numpy.sqrt(num_incoh))
598 interfid = interfid[0]
602 interfid = interfid[0]
599 cinterfid = interfid.size
603 cinterfid = interfid.size
600
604
601 if (cnoiseid > 0): jspc_interf[noiseid] = 0
605 if (cnoiseid > 0): jspc_interf[noiseid] = 0
602
606
603 #Expandiendo los perfiles a limpiar
607 #Expandiendo los perfiles a limpiar
604 if (cinterfid > 0):
608 if (cinterfid > 0):
605 new_interfid = (numpy.r_[interfid - 1, interfid, interfid + 1] + num_prof)%num_prof
609 new_interfid = (numpy.r_[interfid - 1, interfid, interfid + 1] + num_prof)%num_prof
606 new_interfid = numpy.asarray(new_interfid)
610 new_interfid = numpy.asarray(new_interfid)
607 new_interfid = {x for x in new_interfid}
611 new_interfid = {x for x in new_interfid}
608 new_interfid = numpy.array(list(new_interfid))
612 new_interfid = numpy.array(list(new_interfid))
609 new_cinterfid = new_interfid.size
613 new_cinterfid = new_interfid.size
610 else: new_cinterfid = 0
614 else: new_cinterfid = 0
611
615
612 for ip in range(new_cinterfid):
616 for ip in range(new_cinterfid):
613 ind = junkspc_interf[:,new_interfid[ip]].ravel().argsort()
617 ind = junkspc_interf[:,new_interfid[ip]].ravel().argsort()
614 jspc_interf[new_interfid[ip]] = junkspc_interf[ind[nhei_interf/2],new_interfid[ip]]
618 jspc_interf[new_interfid[ip]] = junkspc_interf[ind[nhei_interf/2],new_interfid[ip]]
615
619
616
620
617 jspectra[ich,:,ind_hei] = jspectra[ich,:,ind_hei] - jspc_interf #Corregir indices
621 jspectra[ich,:,ind_hei] = jspectra[ich,:,ind_hei] - jspc_interf #Corregir indices
618
622
619 #Removiendo la interferencia del punto de mayor interferencia
623 #Removiendo la interferencia del punto de mayor interferencia
620 ListAux = jspc_interf[mask_prof].tolist()
624 ListAux = jspc_interf[mask_prof].tolist()
621 maxid = ListAux.index(max(ListAux))
625 maxid = ListAux.index(max(ListAux))
622
626
623
627
624 if cinterfid > 0:
628 if cinterfid > 0:
625 for ip in range(cinterfid*(interf == 2) - 1):
629 for ip in range(cinterfid*(interf == 2) - 1):
626 ind = (jspectra[ich,interfid[ip],:] < tmp_noise*(1 + 1/numpy.sqrt(num_incoh))).nonzero()
630 ind = (jspectra[ich,interfid[ip],:] < tmp_noise*(1 + 1/numpy.sqrt(num_incoh))).nonzero()
627 cind = len(ind)
631 cind = len(ind)
628
632
629 if (cind > 0):
633 if (cind > 0):
630 jspectra[ich,interfid[ip],ind] = tmp_noise*(1 + (numpy.random.uniform(cind) - 0.5)/numpy.sqrt(num_incoh))
634 jspectra[ich,interfid[ip],ind] = tmp_noise*(1 + (numpy.random.uniform(cind) - 0.5)/numpy.sqrt(num_incoh))
631
635
632 ind = numpy.array([-2,-1,1,2])
636 ind = numpy.array([-2,-1,1,2])
633 xx = numpy.zeros([4,4])
637 xx = numpy.zeros([4,4])
634
638
635 for id1 in range(4):
639 for id1 in range(4):
636 xx[:,id1] = ind[id1]**numpy.asarray(range(4))
640 xx[:,id1] = ind[id1]**numpy.asarray(range(4))
637
641
638 xx_inv = numpy.linalg.inv(xx)
642 xx_inv = numpy.linalg.inv(xx)
639 xx = xx_inv[:,0]
643 xx = xx_inv[:,0]
640 ind = (ind + maxid + num_mask_prof)%num_mask_prof
644 ind = (ind + maxid + num_mask_prof)%num_mask_prof
641 yy = jspectra[ich,mask_prof[ind],:]
645 yy = jspectra[ich,mask_prof[ind],:]
642 jspectra[ich,mask_prof[maxid],:] = numpy.dot(yy.transpose(),xx)
646 jspectra[ich,mask_prof[maxid],:] = numpy.dot(yy.transpose(),xx)
643
647
644
648
645 indAux = (jspectra[ich,:,:] < tmp_noise*(1-1/numpy.sqrt(num_incoh))).nonzero()
649 indAux = (jspectra[ich,:,:] < tmp_noise*(1-1/numpy.sqrt(num_incoh))).nonzero()
646 jspectra[ich,indAux[0],indAux[1]] = tmp_noise * (1 - 1/numpy.sqrt(num_incoh))
650 jspectra[ich,indAux[0],indAux[1]] = tmp_noise * (1 - 1/numpy.sqrt(num_incoh))
647
651
648 #Remocion de Interferencia en el Cross Spectra
652 #Remocion de Interferencia en el Cross Spectra
649 if jcspectra is None: return jspectra, jcspectra
653 if jcspectra is None: return jspectra, jcspectra
650 num_pairs = jcspectra.size/(num_prof*num_hei)
654 num_pairs = jcspectra.size/(num_prof*num_hei)
651 jcspectra = jcspectra.reshape(num_pairs, num_prof, num_hei)
655 jcspectra = jcspectra.reshape(num_pairs, num_prof, num_hei)
652
656
653 for ip in range(num_pairs):
657 for ip in range(num_pairs):
654
658
655 #-------------------------------------------
659 #-------------------------------------------
656
660
657 cspower = numpy.abs(jcspectra[ip,mask_prof,:])
661 cspower = numpy.abs(jcspectra[ip,mask_prof,:])
658 cspower = cspower[:,hei_interf]
662 cspower = cspower[:,hei_interf]
659 cspower = cspower.sum(axis = 0)
663 cspower = cspower.sum(axis = 0)
660
664
661 cspsort = cspower.ravel().argsort()
665 cspsort = cspower.ravel().argsort()
662 junkcspc_interf = jcspectra[ip,:,hei_interf[cspsort[range(offhei_interf, nhei_interf + offhei_interf)]]]
666 junkcspc_interf = jcspectra[ip,:,hei_interf[cspsort[range(offhei_interf, nhei_interf + offhei_interf)]]]
663 junkcspc_interf = junkcspc_interf.transpose()
667 junkcspc_interf = junkcspc_interf.transpose()
664 jcspc_interf = junkcspc_interf.sum(axis = 1)/nhei_interf
668 jcspc_interf = junkcspc_interf.sum(axis = 1)/nhei_interf
665
669
666 ind = numpy.abs(jcspc_interf[mask_prof]).ravel().argsort()
670 ind = numpy.abs(jcspc_interf[mask_prof]).ravel().argsort()
667
671
668 median_real = numpy.median(numpy.real(junkcspc_interf[mask_prof[ind[range(3*num_prof/4)]],:]))
672 median_real = numpy.median(numpy.real(junkcspc_interf[mask_prof[ind[range(3*num_prof/4)]],:]))
669 median_imag = numpy.median(numpy.imag(junkcspc_interf[mask_prof[ind[range(3*num_prof/4)]],:]))
673 median_imag = numpy.median(numpy.imag(junkcspc_interf[mask_prof[ind[range(3*num_prof/4)]],:]))
670 junkcspc_interf[comp_mask_prof,:] = numpy.complex(median_real, median_imag)
674 junkcspc_interf[comp_mask_prof,:] = numpy.complex(median_real, median_imag)
671
675
672 for iprof in range(num_prof):
676 for iprof in range(num_prof):
673 ind = numpy.abs(junkcspc_interf[iprof,:]).ravel().argsort()
677 ind = numpy.abs(junkcspc_interf[iprof,:]).ravel().argsort()
674 jcspc_interf[iprof] = junkcspc_interf[iprof, ind[nhei_interf/2]]
678 jcspc_interf[iprof] = junkcspc_interf[iprof, ind[nhei_interf/2]]
675
679
676 #Removiendo la Interferencia
680 #Removiendo la Interferencia
677 jcspectra[ip,:,ind_hei] = jcspectra[ip,:,ind_hei] - jcspc_interf
681 jcspectra[ip,:,ind_hei] = jcspectra[ip,:,ind_hei] - jcspc_interf
678
682
679 ListAux = numpy.abs(jcspc_interf[mask_prof]).tolist()
683 ListAux = numpy.abs(jcspc_interf[mask_prof]).tolist()
680 maxid = ListAux.index(max(ListAux))
684 maxid = ListAux.index(max(ListAux))
681
685
682 ind = numpy.array([-2,-1,1,2])
686 ind = numpy.array([-2,-1,1,2])
683 xx = numpy.zeros([4,4])
687 xx = numpy.zeros([4,4])
684
688
685 for id1 in range(4):
689 for id1 in range(4):
686 xx[:,id1] = ind[id1]**numpy.asarray(range(4))
690 xx[:,id1] = ind[id1]**numpy.asarray(range(4))
687
691
688 xx_inv = numpy.linalg.inv(xx)
692 xx_inv = numpy.linalg.inv(xx)
689 xx = xx_inv[:,0]
693 xx = xx_inv[:,0]
690
694
691 ind = (ind + maxid + num_mask_prof)%num_mask_prof
695 ind = (ind + maxid + num_mask_prof)%num_mask_prof
692 yy = jcspectra[ip,mask_prof[ind],:]
696 yy = jcspectra[ip,mask_prof[ind],:]
693 jcspectra[ip,mask_prof[maxid],:] = numpy.dot(yy.transpose(),xx)
697 jcspectra[ip,mask_prof[maxid],:] = numpy.dot(yy.transpose(),xx)
694
698
695 #Guardar Resultados
699 #Guardar Resultados
696 self.dataOut.data_spc = jspectra
700 self.dataOut.data_spc = jspectra
697 self.dataOut.data_cspc = jcspectra
701 self.dataOut.data_cspc = jcspectra
698
702
699 return 1
703 return 1
700
704
701 def setRadarFrequency(self, frequency=None):
705 def setRadarFrequency(self, frequency=None):
702
706
703 if frequency != None:
707 if frequency != None:
704 self.dataOut.frequency = frequency
708 self.dataOut.frequency = frequency
705
709
706 return 1
710 return 1
707
711
708 def getNoise(self, minHei=None, maxHei=None, minVel=None, maxVel=None):
712 def getNoise(self, minHei=None, maxHei=None, minVel=None, maxVel=None):
709 #validacion de rango
713 #validacion de rango
710 if minHei == None:
714 if minHei == None:
711 minHei = self.dataOut.heightList[0]
715 minHei = self.dataOut.heightList[0]
712
716
713 if maxHei == None:
717 if maxHei == None:
714 maxHei = self.dataOut.heightList[-1]
718 maxHei = self.dataOut.heightList[-1]
715
719
716 if (minHei < self.dataOut.heightList[0]) or (minHei > maxHei):
720 if (minHei < self.dataOut.heightList[0]) or (minHei > maxHei):
717 print 'minHei: %.2f is out of the heights range'%(minHei)
721 print 'minHei: %.2f is out of the heights range'%(minHei)
718 print 'minHei is setting to %.2f'%(self.dataOut.heightList[0])
722 print 'minHei is setting to %.2f'%(self.dataOut.heightList[0])
719 minHei = self.dataOut.heightList[0]
723 minHei = self.dataOut.heightList[0]
720
724
721 if (maxHei > self.dataOut.heightList[-1]) or (maxHei < minHei):
725 if (maxHei > self.dataOut.heightList[-1]) or (maxHei < minHei):
722 print 'maxHei: %.2f is out of the heights range'%(maxHei)
726 print 'maxHei: %.2f is out of the heights range'%(maxHei)
723 print 'maxHei is setting to %.2f'%(self.dataOut.heightList[-1])
727 print 'maxHei is setting to %.2f'%(self.dataOut.heightList[-1])
724 maxHei = self.dataOut.heightList[-1]
728 maxHei = self.dataOut.heightList[-1]
725
729
726 # validacion de velocidades
730 # validacion de velocidades
727 velrange = self.dataOut.getVelRange(1)
731 velrange = self.dataOut.getVelRange(1)
728
732
729 if minVel == None:
733 if minVel == None:
730 minVel = velrange[0]
734 minVel = velrange[0]
731
735
732 if maxVel == None:
736 if maxVel == None:
733 maxVel = velrange[-1]
737 maxVel = velrange[-1]
734
738
735 if (minVel < velrange[0]) or (minVel > maxVel):
739 if (minVel < velrange[0]) or (minVel > maxVel):
736 print 'minVel: %.2f is out of the velocity range'%(minVel)
740 print 'minVel: %.2f is out of the velocity range'%(minVel)
737 print 'minVel is setting to %.2f'%(velrange[0])
741 print 'minVel is setting to %.2f'%(velrange[0])
738 minVel = velrange[0]
742 minVel = velrange[0]
739
743
740 if (maxVel > velrange[-1]) or (maxVel < minVel):
744 if (maxVel > velrange[-1]) or (maxVel < minVel):
741 print 'maxVel: %.2f is out of the velocity range'%(maxVel)
745 print 'maxVel: %.2f is out of the velocity range'%(maxVel)
742 print 'maxVel is setting to %.2f'%(velrange[-1])
746 print 'maxVel is setting to %.2f'%(velrange[-1])
743 maxVel = velrange[-1]
747 maxVel = velrange[-1]
744
748
745 # seleccion de indices para rango
749 # seleccion de indices para rango
746 minIndex = 0
750 minIndex = 0
747 maxIndex = 0
751 maxIndex = 0
748 heights = self.dataOut.heightList
752 heights = self.dataOut.heightList
749
753
750 inda = numpy.where(heights >= minHei)
754 inda = numpy.where(heights >= minHei)
751 indb = numpy.where(heights <= maxHei)
755 indb = numpy.where(heights <= maxHei)
752
756
753 try:
757 try:
754 minIndex = inda[0][0]
758 minIndex = inda[0][0]
755 except:
759 except:
756 minIndex = 0
760 minIndex = 0
757
761
758 try:
762 try:
759 maxIndex = indb[0][-1]
763 maxIndex = indb[0][-1]
760 except:
764 except:
761 maxIndex = len(heights)
765 maxIndex = len(heights)
762
766
763 if (minIndex < 0) or (minIndex > maxIndex):
767 if (minIndex < 0) or (minIndex > maxIndex):
764 raise ValueError, "some value in (%d,%d) is not valid" % (minIndex, maxIndex)
768 raise ValueError, "some value in (%d,%d) is not valid" % (minIndex, maxIndex)
765
769
766 if (maxIndex >= self.dataOut.nHeights):
770 if (maxIndex >= self.dataOut.nHeights):
767 maxIndex = self.dataOut.nHeights-1
771 maxIndex = self.dataOut.nHeights-1
768
772
769 # seleccion de indices para velocidades
773 # seleccion de indices para velocidades
770 indminvel = numpy.where(velrange >= minVel)
774 indminvel = numpy.where(velrange >= minVel)
771 indmaxvel = numpy.where(velrange <= maxVel)
775 indmaxvel = numpy.where(velrange <= maxVel)
772 try:
776 try:
773 minIndexVel = indminvel[0][0]
777 minIndexVel = indminvel[0][0]
774 except:
778 except:
775 minIndexVel = 0
779 minIndexVel = 0
776
780
777 try:
781 try:
778 maxIndexVel = indmaxvel[0][-1]
782 maxIndexVel = indmaxvel[0][-1]
779 except:
783 except:
780 maxIndexVel = len(velrange)
784 maxIndexVel = len(velrange)
781
785
782 #seleccion del espectro
786 #seleccion del espectro
783 data_spc = self.dataOut.data_spc[:,minIndexVel:maxIndexVel+1,minIndex:maxIndex+1]
787 data_spc = self.dataOut.data_spc[:,minIndexVel:maxIndexVel+1,minIndex:maxIndex+1]
784 #estimacion de ruido
788 #estimacion de ruido
785 noise = numpy.zeros(self.dataOut.nChannels)
789 noise = numpy.zeros(self.dataOut.nChannels)
786
790
787 for channel in range(self.dataOut.nChannels):
791 for channel in range(self.dataOut.nChannels):
788 daux = data_spc[channel,:,:]
792 daux = data_spc[channel,:,:]
789 noise[channel] = hildebrand_sekhon(daux, self.dataOut.nIncohInt)
793 noise[channel] = hildebrand_sekhon(daux, self.dataOut.nIncohInt)
790
794
791 self.dataOut.noise_estimation = noise.copy()
795 self.dataOut.noise_estimation = noise.copy()
792
796
793 return 1
797 return 1
@@ -1,1530 +1,1540
1 import sys
1 import sys
2 import numpy
2 import numpy
3 from scipy import interpolate
3 from scipy import interpolate
4 from schainpy import cSchain
4 from schainpy import cSchain
5 from jroproc_base import ProcessingUnit, Operation
5 from jroproc_base import ProcessingUnit, Operation
6 from schainpy.model.data.jrodata import Voltage
6 from schainpy.model.data.jrodata import Voltage
7 from time import time
7 from time import time
8
8
9 import math
9 import math
10
10
11 def rep_seq(x, rep=10):
11 def rep_seq(x, rep=10):
12 L = len(x) * rep
12 L = len(x) * rep
13 res = numpy.zeros(L, dtype=x.dtype)
13 res = numpy.zeros(L, dtype=x.dtype)
14 idx = numpy.arange(len(x)) * rep
14 idx = numpy.arange(len(x)) * rep
15 for i in numpy.arange(rep):
15 for i in numpy.arange(rep):
16 res[idx + i] = x
16 res[idx + i] = x
17 return(res)
17 return(res)
18
18
19
19
20 def create_pseudo_random_code(clen=10000, seed=0):
20 def create_pseudo_random_code(clen=10000, seed=0):
21 """
21 """
22 seed is a way of reproducing the random code without
22 seed is a way of reproducing the random code without
23 having to store all actual codes. the seed can then
23 having to store all actual codes. the seed can then
24 act as a sort of station_id.
24 act as a sort of station_id.
25
25
26 """
26 """
27 numpy.random.seed(seed)
27 numpy.random.seed(seed)
28 phases = numpy.array(
28 phases = numpy.array(
29 numpy.exp(1.0j * 2.0 * math.pi * numpy.random.random(clen)),
29 numpy.exp(1.0j * 2.0 * math.pi * numpy.random.random(clen)),
30 dtype=numpy.complex64,
30 dtype=numpy.complex64,
31 )
31 )
32 return(phases)
32 return(phases)
33
33
34
34
35 def periodic_convolution_matrix(envelope, rmin=0, rmax=100):
35 def periodic_convolution_matrix(envelope, rmin=0, rmax=100):
36 """
36 """
37 we imply that the number of measurements is equal to the number of elements
37 we imply that the number of measurements is equal to the number of elements
38 in code
38 in code
39
39
40 """
40 """
41 L = len(envelope)
41 L = len(envelope)
42 ridx = numpy.arange(rmin, rmax)
42 ridx = numpy.arange(rmin, rmax)
43 A = numpy.zeros([L, rmax-rmin], dtype=numpy.complex64)
43 A = numpy.zeros([L, rmax-rmin], dtype=numpy.complex64)
44 for i in numpy.arange(L):
44 for i in numpy.arange(L):
45 A[i, :] = envelope[(i-ridx) % L]
45 A[i, :] = envelope[(i-ridx) % L]
46 result = {}
46 result = {}
47 result['A'] = A
47 result['A'] = A
48 result['ridx'] = ridx
48 result['ridx'] = ridx
49 return(result)
49 return(result)
50
50
51
51
52 B_cache = 0
52 B_cache = 0
53 r_cache = 0
53 r_cache = 0
54 B_cached = False
54 B_cached = False
55 def create_estimation_matrix(code, rmin=0, rmax=1000, cache=True):
55 def create_estimation_matrix(code, rmin=0, rmax=1000, cache=True):
56 global B_cache
56 global B_cache
57 global r_cache
57 global r_cache
58 global B_cached
58 global B_cached
59
59
60 if not cache or not B_cached:
60 if not cache or not B_cached:
61 r_cache = periodic_convolution_matrix(
61 r_cache = periodic_convolution_matrix(
62 envelope=code, rmin=rmin, rmax=rmax,
62 envelope=code, rmin=rmin, rmax=rmax,
63 )
63 )
64 A = r_cache['A']
64 A = r_cache['A']
65 Ah = numpy.transpose(numpy.conjugate(A))
65 Ah = numpy.transpose(numpy.conjugate(A))
66 B_cache = numpy.dot(numpy.linalg.inv(numpy.dot(Ah, A)), Ah)
66 B_cache = numpy.dot(numpy.linalg.inv(numpy.dot(Ah, A)), Ah)
67 r_cache['B'] = B_cache
67 r_cache['B'] = B_cache
68 B_cached = True
68 B_cached = True
69 return(r_cache)
69 return(r_cache)
70 else:
70 else:
71 return(r_cache)
71 return(r_cache)
72
72
73 class VoltageProc(ProcessingUnit):
73 class VoltageProc(ProcessingUnit):
74
74
75
75
76 def __init__(self, **kwargs):
76 def __init__(self, **kwargs):
77
77
78 ProcessingUnit.__init__(self, **kwargs)
78 ProcessingUnit.__init__(self, **kwargs)
79
79
80 # self.objectDict = {}
80 # self.objectDict = {}
81 self.dataOut = Voltage()
81 self.dataOut = Voltage()
82 self.flip = 1
82 self.flip = 1
83
83
84 def run(self):
84 def run(self):
85 if self.dataIn.type == 'AMISR':
85 if self.dataIn.type == 'AMISR':
86 self.__updateObjFromAmisrInput()
86 self.__updateObjFromAmisrInput()
87
87
88 if self.dataIn.type == 'Voltage':
88 if self.dataIn.type == 'Voltage':
89 self.dataOut.copy(self.dataIn)
89 self.dataOut.copy(self.dataIn)
90
90
91 # self.dataOut.copy(self.dataIn)
91 # self.dataOut.copy(self.dataIn)
92
92
93 def __updateObjFromAmisrInput(self):
93 def __updateObjFromAmisrInput(self):
94
94
95 self.dataOut.timeZone = self.dataIn.timeZone
95 self.dataOut.timeZone = self.dataIn.timeZone
96 self.dataOut.dstFlag = self.dataIn.dstFlag
96 self.dataOut.dstFlag = self.dataIn.dstFlag
97 self.dataOut.errorCount = self.dataIn.errorCount
97 self.dataOut.errorCount = self.dataIn.errorCount
98 self.dataOut.useLocalTime = self.dataIn.useLocalTime
98 self.dataOut.useLocalTime = self.dataIn.useLocalTime
99
99
100 self.dataOut.flagNoData = self.dataIn.flagNoData
100 self.dataOut.flagNoData = self.dataIn.flagNoData
101 self.dataOut.data = self.dataIn.data
101 self.dataOut.data = self.dataIn.data
102 self.dataOut.utctime = self.dataIn.utctime
102 self.dataOut.utctime = self.dataIn.utctime
103 self.dataOut.channelList = self.dataIn.channelList
103 self.dataOut.channelList = self.dataIn.channelList
104 # self.dataOut.timeInterval = self.dataIn.timeInterval
104 # self.dataOut.timeInterval = self.dataIn.timeInterval
105 self.dataOut.heightList = self.dataIn.heightList
105 self.dataOut.heightList = self.dataIn.heightList
106 self.dataOut.nProfiles = self.dataIn.nProfiles
106 self.dataOut.nProfiles = self.dataIn.nProfiles
107
107
108 self.dataOut.nCohInt = self.dataIn.nCohInt
108 self.dataOut.nCohInt = self.dataIn.nCohInt
109 self.dataOut.ippSeconds = self.dataIn.ippSeconds
109 self.dataOut.ippSeconds = self.dataIn.ippSeconds
110 self.dataOut.frequency = self.dataIn.frequency
110 self.dataOut.frequency = self.dataIn.frequency
111
111
112 self.dataOut.azimuth = self.dataIn.azimuth
112 self.dataOut.azimuth = self.dataIn.azimuth
113 self.dataOut.zenith = self.dataIn.zenith
113 self.dataOut.zenith = self.dataIn.zenith
114
114
115 self.dataOut.beam.codeList = self.dataIn.beam.codeList
115 self.dataOut.beam.codeList = self.dataIn.beam.codeList
116 self.dataOut.beam.azimuthList = self.dataIn.beam.azimuthList
116 self.dataOut.beam.azimuthList = self.dataIn.beam.azimuthList
117 self.dataOut.beam.zenithList = self.dataIn.beam.zenithList
117 self.dataOut.beam.zenithList = self.dataIn.beam.zenithList
118 #
118 #
119 # pass#
119 # pass#
120 #
120 #
121 # def init(self):
121 # def init(self):
122 #
122 #
123 #
123 #
124 # if self.dataIn.type == 'AMISR':
124 # if self.dataIn.type == 'AMISR':
125 # self.__updateObjFromAmisrInput()
125 # self.__updateObjFromAmisrInput()
126 #
126 #
127 # if self.dataIn.type == 'Voltage':
127 # if self.dataIn.type == 'Voltage':
128 # self.dataOut.copy(self.dataIn)
128 # self.dataOut.copy(self.dataIn)
129 # # No necesita copiar en cada init() los atributos de dataIn
129 # # No necesita copiar en cada init() los atributos de dataIn
130 # # la copia deberia hacerse por cada nuevo bloque de datos
130 # # la copia deberia hacerse por cada nuevo bloque de datos
131
131
132 def selectChannels(self, channelList):
132 def selectChannels(self, channelList):
133
133
134 channelIndexList = []
134 channelIndexList = []
135
135
136 for channel in channelList:
136 for channel in channelList:
137 if channel not in self.dataOut.channelList:
137 if channel not in self.dataOut.channelList:
138 raise ValueError, "Channel %d is not in %s" %(channel, str(self.dataOut.channelList))
138 raise ValueError, "Channel %d is not in %s" %(channel, str(self.dataOut.channelList))
139
139
140 index = self.dataOut.channelList.index(channel)
140 index = self.dataOut.channelList.index(channel)
141 channelIndexList.append(index)
141 channelIndexList.append(index)
142
142
143 self.selectChannelsByIndex(channelIndexList)
143 self.selectChannelsByIndex(channelIndexList)
144
144
145 def selectChannelsByIndex(self, channelIndexList):
145 def selectChannelsByIndex(self, channelIndexList):
146 """
146 """
147 Selecciona un bloque de datos en base a canales segun el channelIndexList
147 Selecciona un bloque de datos en base a canales segun el channelIndexList
148
148
149 Input:
149 Input:
150 channelIndexList : lista sencilla de canales a seleccionar por ej. [2,3,7]
150 channelIndexList : lista sencilla de canales a seleccionar por ej. [2,3,7]
151
151
152 Affected:
152 Affected:
153 self.dataOut.data
153 self.dataOut.data
154 self.dataOut.channelIndexList
154 self.dataOut.channelIndexList
155 self.dataOut.nChannels
155 self.dataOut.nChannels
156 self.dataOut.m_ProcessingHeader.totalSpectra
156 self.dataOut.m_ProcessingHeader.totalSpectra
157 self.dataOut.systemHeaderObj.numChannels
157 self.dataOut.systemHeaderObj.numChannels
158 self.dataOut.m_ProcessingHeader.blockSize
158 self.dataOut.m_ProcessingHeader.blockSize
159
159
160 Return:
160 Return:
161 None
161 None
162 """
162 """
163
163
164 for channelIndex in channelIndexList:
164 for channelIndex in channelIndexList:
165 if channelIndex not in self.dataOut.channelIndexList:
165 if channelIndex not in self.dataOut.channelIndexList:
166 print channelIndexList
166 print channelIndexList
167 raise ValueError, "The value %d in channelIndexList is not valid" %channelIndex
167 raise ValueError, "The value %d in channelIndexList is not valid" %channelIndex
168
168
169 if self.dataOut.flagDataAsBlock:
169 if self.dataOut.flagDataAsBlock:
170 """
170 """
171 Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
171 Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
172 """
172 """
173 data = self.dataOut.data[channelIndexList,:,:]
173 data = self.dataOut.data[channelIndexList,:,:]
174 else:
174 else:
175 data = self.dataOut.data[channelIndexList,:]
175 data = self.dataOut.data[channelIndexList,:]
176
176
177 self.dataOut.data = data
177 self.dataOut.data = data
178 self.dataOut.channelList = [self.dataOut.channelList[i] for i in channelIndexList]
178 self.dataOut.channelList = [self.dataOut.channelList[i] for i in channelIndexList]
179 # self.dataOut.nChannels = nChannels
179 # self.dataOut.nChannels = nChannels
180
180
181 return 1
181 return 1
182
182
183 def selectHeights(self, minHei=None, maxHei=None):
183 def selectHeights(self, minHei=None, maxHei=None):
184 """
184 """
185 Selecciona un bloque de datos en base a un grupo de valores de alturas segun el rango
185 Selecciona un bloque de datos en base a un grupo de valores de alturas segun el rango
186 minHei <= height <= maxHei
186 minHei <= height <= maxHei
187
187
188 Input:
188 Input:
189 minHei : valor minimo de altura a considerar
189 minHei : valor minimo de altura a considerar
190 maxHei : valor maximo de altura a considerar
190 maxHei : valor maximo de altura a considerar
191
191
192 Affected:
192 Affected:
193 Indirectamente son cambiados varios valores a travez del metodo selectHeightsByIndex
193 Indirectamente son cambiados varios valores a travez del metodo selectHeightsByIndex
194
194
195 Return:
195 Return:
196 1 si el metodo se ejecuto con exito caso contrario devuelve 0
196 1 si el metodo se ejecuto con exito caso contrario devuelve 0
197 """
197 """
198
198
199 if minHei == None:
199 if minHei == None:
200 minHei = self.dataOut.heightList[0]
200 minHei = self.dataOut.heightList[0]
201
201
202 if maxHei == None:
202 if maxHei == None:
203 maxHei = self.dataOut.heightList[-1]
203 maxHei = self.dataOut.heightList[-1]
204
204
205 if (minHei < self.dataOut.heightList[0]):
205 if (minHei < self.dataOut.heightList[0]):
206 minHei = self.dataOut.heightList[0]
206 minHei = self.dataOut.heightList[0]
207
207
208 if (maxHei > self.dataOut.heightList[-1]):
208 if (maxHei > self.dataOut.heightList[-1]):
209 maxHei = self.dataOut.heightList[-1]
209 maxHei = self.dataOut.heightList[-1]
210
210
211 minIndex = 0
211 minIndex = 0
212 maxIndex = 0
212 maxIndex = 0
213 heights = self.dataOut.heightList
213 heights = self.dataOut.heightList
214
214
215 inda = numpy.where(heights >= minHei)
215 inda = numpy.where(heights >= minHei)
216 indb = numpy.where(heights <= maxHei)
216 indb = numpy.where(heights <= maxHei)
217
217
218 try:
218 try:
219 minIndex = inda[0][0]
219 minIndex = inda[0][0]
220 except:
220 except:
221 minIndex = 0
221 minIndex = 0
222
222
223 try:
223 try:
224 maxIndex = indb[0][-1]
224 maxIndex = indb[0][-1]
225 except:
225 except:
226 maxIndex = len(heights)
226 maxIndex = len(heights)
227
227
228 self.selectHeightsByIndex(minIndex, maxIndex)
228 self.selectHeightsByIndex(minIndex, maxIndex)
229
229
230 return 1
230 return 1
231
231
232
232
233 def selectHeightsByIndex(self, minIndex, maxIndex):
233 def selectHeightsByIndex(self, minIndex, maxIndex):
234 """
234 """
235 Selecciona un bloque de datos en base a un grupo indices de alturas segun el rango
235 Selecciona un bloque de datos en base a un grupo indices de alturas segun el rango
236 minIndex <= index <= maxIndex
236 minIndex <= index <= maxIndex
237
237
238 Input:
238 Input:
239 minIndex : valor de indice minimo de altura a considerar
239 minIndex : valor de indice minimo de altura a considerar
240 maxIndex : valor de indice maximo de altura a considerar
240 maxIndex : valor de indice maximo de altura a considerar
241
241
242 Affected:
242 Affected:
243 self.dataOut.data
243 self.dataOut.data
244 self.dataOut.heightList
244 self.dataOut.heightList
245
245
246 Return:
246 Return:
247 1 si el metodo se ejecuto con exito caso contrario devuelve 0
247 1 si el metodo se ejecuto con exito caso contrario devuelve 0
248 """
248 """
249
249
250 if (minIndex < 0) or (minIndex > maxIndex):
250 if (minIndex < 0) or (minIndex > maxIndex):
251 raise ValueError, "Height index range (%d,%d) is not valid" % (minIndex, maxIndex)
251 raise ValueError, "Height index range (%d,%d) is not valid" % (minIndex, maxIndex)
252
252
253 if (maxIndex >= self.dataOut.nHeights):
253 if (maxIndex >= self.dataOut.nHeights):
254 maxIndex = self.dataOut.nHeights
254 maxIndex = self.dataOut.nHeights
255
255
256 #voltage
256 #voltage
257 if self.dataOut.flagDataAsBlock:
257 if self.dataOut.flagDataAsBlock:
258 """
258 """
259 Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
259 Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
260 """
260 """
261 data = self.dataOut.data[:,:, minIndex:maxIndex]
261 data = self.dataOut.data[:,:, minIndex:maxIndex]
262 else:
262 else:
263 data = self.dataOut.data[:, minIndex:maxIndex]
263 data = self.dataOut.data[:, minIndex:maxIndex]
264
264
265 # firstHeight = self.dataOut.heightList[minIndex]
265 # firstHeight = self.dataOut.heightList[minIndex]
266
266
267 self.dataOut.data = data
267 self.dataOut.data = data
268 self.dataOut.heightList = self.dataOut.heightList[minIndex:maxIndex]
268 self.dataOut.heightList = self.dataOut.heightList[minIndex:maxIndex]
269
269
270 if self.dataOut.nHeights <= 1:
270 if self.dataOut.nHeights <= 1:
271 raise ValueError, "selectHeights: Too few heights. Current number of heights is %d" %(self.dataOut.nHeights)
271 raise ValueError, "selectHeights: Too few heights. Current number of heights is %d" %(self.dataOut.nHeights)
272
272
273 return 1
273 return 1
274
274
275
275
276 def filterByHeights(self, window):
276 def filterByHeights(self, window):
277
277
278 deltaHeight = self.dataOut.heightList[1] - self.dataOut.heightList[0]
278 deltaHeight = self.dataOut.heightList[1] - self.dataOut.heightList[0]
279
279
280 if window == None:
280 if window == None:
281 window = (self.dataOut.radarControllerHeaderObj.txA/self.dataOut.radarControllerHeaderObj.nBaud) / deltaHeight
281 window = (self.dataOut.radarControllerHeaderObj.txA/self.dataOut.radarControllerHeaderObj.nBaud) / deltaHeight
282
282
283 newdelta = deltaHeight * window
283 newdelta = deltaHeight * window
284 r = self.dataOut.nHeights % window
284 r = self.dataOut.nHeights % window
285 newheights = (self.dataOut.nHeights-r)/window
285 newheights = (self.dataOut.nHeights-r)/window
286
286
287 if newheights <= 1:
287 if newheights <= 1:
288 raise ValueError, "filterByHeights: Too few heights. Current number of heights is %d and window is %d" %(self.dataOut.nHeights, window)
288 raise ValueError, "filterByHeights: Too few heights. Current number of heights is %d and window is %d" %(self.dataOut.nHeights, window)
289
289
290 if self.dataOut.flagDataAsBlock:
290 if self.dataOut.flagDataAsBlock:
291 """
291 """
292 Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
292 Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis]
293 """
293 """
294 buffer = self.dataOut.data[:, :, 0:self.dataOut.nHeights-r]
294 buffer = self.dataOut.data[:, :, 0:self.dataOut.nHeights-r]
295 buffer = buffer.reshape(self.dataOut.nChannels,self.dataOut.nProfiles,self.dataOut.nHeights/window,window)
295 buffer = buffer.reshape(self.dataOut.nChannels,self.dataOut.nProfiles,self.dataOut.nHeights/window,window)
296 buffer = numpy.sum(buffer,3)
296 buffer = numpy.sum(buffer,3)
297
297
298 else:
298 else:
299 buffer = self.dataOut.data[:,0:self.dataOut.nHeights-r]
299 buffer = self.dataOut.data[:,0:self.dataOut.nHeights-r]
300 buffer = buffer.reshape(self.dataOut.nChannels,self.dataOut.nHeights/window,window)
300 buffer = buffer.reshape(self.dataOut.nChannels,self.dataOut.nHeights/window,window)
301 buffer = numpy.sum(buffer,2)
301 buffer = numpy.sum(buffer,2)
302
302
303 self.dataOut.data = buffer
303 self.dataOut.data = buffer
304 self.dataOut.heightList = self.dataOut.heightList[0] + numpy.arange( newheights )*newdelta
304 self.dataOut.heightList = self.dataOut.heightList[0] + numpy.arange( newheights )*newdelta
305 self.dataOut.windowOfFilter = window
305 self.dataOut.windowOfFilter = window
306
306
307 def setH0(self, h0, deltaHeight = None):
307 def setH0(self, h0, deltaHeight = None):
308
308
309 if not deltaHeight:
309 if not deltaHeight:
310 deltaHeight = self.dataOut.heightList[1] - self.dataOut.heightList[0]
310 deltaHeight = self.dataOut.heightList[1] - self.dataOut.heightList[0]
311
311
312 nHeights = self.dataOut.nHeights
312 nHeights = self.dataOut.nHeights
313
313
314 newHeiRange = h0 + numpy.arange(nHeights)*deltaHeight
314 newHeiRange = h0 + numpy.arange(nHeights)*deltaHeight
315
315
316 self.dataOut.heightList = newHeiRange
316 self.dataOut.heightList = newHeiRange
317
317
318 def deFlip(self, channelList = []):
318 def deFlip(self, channelList = []):
319
319
320 data = self.dataOut.data.copy()
320 data = self.dataOut.data.copy()
321
321
322 if self.dataOut.flagDataAsBlock:
322 if self.dataOut.flagDataAsBlock:
323 flip = self.flip
323 flip = self.flip
324 profileList = range(self.dataOut.nProfiles)
324 profileList = range(self.dataOut.nProfiles)
325
325
326 if not channelList:
326 if not channelList:
327 for thisProfile in profileList:
327 for thisProfile in profileList:
328 data[:,thisProfile,:] = data[:,thisProfile,:]*flip
328 data[:,thisProfile,:] = data[:,thisProfile,:]*flip
329 flip *= -1.0
329 flip *= -1.0
330 else:
330 else:
331 for thisChannel in channelList:
331 for thisChannel in channelList:
332 if thisChannel not in self.dataOut.channelList:
332 if thisChannel not in self.dataOut.channelList:
333 continue
333 continue
334
334
335 for thisProfile in profileList:
335 for thisProfile in profileList:
336 data[thisChannel,thisProfile,:] = data[thisChannel,thisProfile,:]*flip
336 data[thisChannel,thisProfile,:] = data[thisChannel,thisProfile,:]*flip
337 flip *= -1.0
337 flip *= -1.0
338
338
339 self.flip = flip
339 self.flip = flip
340
340
341 else:
341 else:
342 if not channelList:
342 if not channelList:
343 data[:,:] = data[:,:]*self.flip
343 data[:,:] = data[:,:]*self.flip
344 else:
344 else:
345 for thisChannel in channelList:
345 for thisChannel in channelList:
346 if thisChannel not in self.dataOut.channelList:
346 if thisChannel not in self.dataOut.channelList:
347 continue
347 continue
348
348
349 data[thisChannel,:] = data[thisChannel,:]*self.flip
349 data[thisChannel,:] = data[thisChannel,:]*self.flip
350
350
351 self.flip *= -1.
351 self.flip *= -1.
352
352
353 self.dataOut.data = data
353 self.dataOut.data = data
354
354
355 def setRadarFrequency(self, frequency=None):
355 def setRadarFrequency(self, frequency=None):
356
356
357 if frequency != None:
357 if frequency != None:
358 self.dataOut.frequency = frequency
358 self.dataOut.frequency = frequency
359
359
360 return 1
360 return 1
361
361
362 def interpolateHeights(self, topLim, botLim):
362 def interpolateHeights(self, topLim, botLim):
363 #69 al 72 para julia
363 #69 al 72 para julia
364 #82-84 para meteoros
364 #82-84 para meteoros
365 if len(numpy.shape(self.dataOut.data))==2:
365 if len(numpy.shape(self.dataOut.data))==2:
366 sampInterp = (self.dataOut.data[:,botLim-1] + self.dataOut.data[:,topLim+1])/2
366 sampInterp = (self.dataOut.data[:,botLim-1] + self.dataOut.data[:,topLim+1])/2
367 sampInterp = numpy.transpose(numpy.tile(sampInterp,(topLim-botLim + 1,1)))
367 sampInterp = numpy.transpose(numpy.tile(sampInterp,(topLim-botLim + 1,1)))
368 #self.dataOut.data[:,botLim:limSup+1] = sampInterp
368 #self.dataOut.data[:,botLim:limSup+1] = sampInterp
369 self.dataOut.data[:,botLim:topLim+1] = sampInterp
369 self.dataOut.data[:,botLim:topLim+1] = sampInterp
370 else:
370 else:
371 nHeights = self.dataOut.data.shape[2]
371 nHeights = self.dataOut.data.shape[2]
372 x = numpy.hstack((numpy.arange(botLim),numpy.arange(topLim+1,nHeights)))
372 x = numpy.hstack((numpy.arange(botLim),numpy.arange(topLim+1,nHeights)))
373 y = self.dataOut.data[:,:,range(botLim)+range(topLim+1,nHeights)]
373 y = self.dataOut.data[:,:,range(botLim)+range(topLim+1,nHeights)]
374 f = interpolate.interp1d(x, y, axis = 2)
374 f = interpolate.interp1d(x, y, axis = 2)
375 xnew = numpy.arange(botLim,topLim+1)
375 xnew = numpy.arange(botLim,topLim+1)
376 ynew = f(xnew)
376 ynew = f(xnew)
377
377
378 self.dataOut.data[:,:,botLim:topLim+1] = ynew
378 self.dataOut.data[:,:,botLim:topLim+1] = ynew
379
379
380 # import collections
380 # import collections
381
381
382 class CohInt(Operation):
382 class CohInt(Operation):
383
383
384 isConfig = False
384 isConfig = False
385 __profIndex = 0
385 __profIndex = 0
386 __byTime = False
386 __byTime = False
387 __initime = None
387 __initime = None
388 __lastdatatime = None
388 __lastdatatime = None
389 __integrationtime = None
389 __integrationtime = None
390 __buffer = None
390 __buffer = None
391 __bufferStride = []
391 __bufferStride = []
392 __dataReady = False
392 __dataReady = False
393 __profIndexStride = 0
393 __profIndexStride = 0
394 __dataToPutStride = False
394 __dataToPutStride = False
395 n = None
395 n = None
396
396
397 def __init__(self, **kwargs):
397 def __init__(self, **kwargs):
398
398
399 Operation.__init__(self, **kwargs)
399 Operation.__init__(self, **kwargs)
400
400
401 # self.isConfig = False
401 # self.isConfig = False
402
402
403 def setup(self, n=None, timeInterval=None, stride=None, overlapping=False, byblock=False):
403 def setup(self, n=None, timeInterval=None, stride=None, overlapping=False, byblock=False):
404 """
404 """
405 Set the parameters of the integration class.
405 Set the parameters of the integration class.
406
406
407 Inputs:
407 Inputs:
408
408
409 n : Number of coherent integrations
409 n : Number of coherent integrations
410 timeInterval : Time of integration. If the parameter "n" is selected this one does not work
410 timeInterval : Time of integration. If the parameter "n" is selected this one does not work
411 overlapping :
411 overlapping :
412 """
412 """
413
413
414 self.__initime = None
414 self.__initime = None
415 self.__lastdatatime = 0
415 self.__lastdatatime = 0
416 self.__buffer = None
416 self.__buffer = None
417 self.__dataReady = False
417 self.__dataReady = False
418 self.byblock = byblock
418 self.byblock = byblock
419 self.stride = stride
419 self.stride = stride
420
420
421 if n == None and timeInterval == None:
421 if n == None and timeInterval == None:
422 raise ValueError, "n or timeInterval should be specified ..."
422 raise ValueError, "n or timeInterval should be specified ..."
423
423
424 if n != None:
424 if n != None:
425 self.n = n
425 self.n = n
426 self.__byTime = False
426 self.__byTime = False
427 else:
427 else:
428 self.__integrationtime = timeInterval #* 60. #if (type(timeInterval)!=integer) -> change this line
428 self.__integrationtime = timeInterval #* 60. #if (type(timeInterval)!=integer) -> change this line
429 self.n = 9999
429 self.n = 9999
430 self.__byTime = True
430 self.__byTime = True
431
431
432 if overlapping:
432 if overlapping:
433 self.__withOverlapping = True
433 self.__withOverlapping = True
434 self.__buffer = None
434 self.__buffer = None
435 else:
435 else:
436 self.__withOverlapping = False
436 self.__withOverlapping = False
437 self.__buffer = 0
437 self.__buffer = 0
438
438
439 self.__profIndex = 0
439 self.__profIndex = 0
440
440
441 def putData(self, data):
441 def putData(self, data):
442
442
443 """
443 """
444 Add a profile to the __buffer and increase in one the __profileIndex
444 Add a profile to the __buffer and increase in one the __profileIndex
445
445
446 """
446 """
447
447
448 if not self.__withOverlapping:
448 if not self.__withOverlapping:
449 self.__buffer += data.copy()
449 self.__buffer += data.copy()
450 self.__profIndex += 1
450 self.__profIndex += 1
451 return
451 return
452
452
453 #Overlapping data
453 #Overlapping data
454 nChannels, nHeis = data.shape
454 nChannels, nHeis = data.shape
455 data = numpy.reshape(data, (1, nChannels, nHeis))
455 data = numpy.reshape(data, (1, nChannels, nHeis))
456
456
457 #If the buffer is empty then it takes the data value
457 #If the buffer is empty then it takes the data value
458 if self.__buffer is None:
458 if self.__buffer is None:
459 self.__buffer = data
459 self.__buffer = data
460 self.__profIndex += 1
460 self.__profIndex += 1
461 return
461 return
462
462
463 #If the buffer length is lower than n then stakcing the data value
463 #If the buffer length is lower than n then stakcing the data value
464 if self.__profIndex < self.n:
464 if self.__profIndex < self.n:
465 self.__buffer = numpy.vstack((self.__buffer, data))
465 self.__buffer = numpy.vstack((self.__buffer, data))
466 self.__profIndex += 1
466 self.__profIndex += 1
467 return
467 return
468
468
469 #If the buffer length is equal to n then replacing the last buffer value with the data value
469 #If the buffer length is equal to n then replacing the last buffer value with the data value
470 self.__buffer = numpy.roll(self.__buffer, -1, axis=0)
470 self.__buffer = numpy.roll(self.__buffer, -1, axis=0)
471 self.__buffer[self.n-1] = data
471 self.__buffer[self.n-1] = data
472 self.__profIndex = self.n
472 self.__profIndex = self.n
473 return
473 return
474
474
475
475
476 def pushData(self):
476 def pushData(self):
477 """
477 """
478 Return the sum of the last profiles and the profiles used in the sum.
478 Return the sum of the last profiles and the profiles used in the sum.
479
479
480 Affected:
480 Affected:
481
481
482 self.__profileIndex
482 self.__profileIndex
483
483
484 """
484 """
485
485
486 if not self.__withOverlapping:
486 if not self.__withOverlapping:
487 data = self.__buffer
487 data = self.__buffer
488 n = self.__profIndex
488 n = self.__profIndex
489
489
490 self.__buffer = 0
490 self.__buffer = 0
491 self.__profIndex = 0
491 self.__profIndex = 0
492
492
493 return data, n
493 return data, n
494
494
495 #Integration with Overlapping
495 #Integration with Overlapping
496 data = numpy.sum(self.__buffer, axis=0)
496 data = numpy.sum(self.__buffer, axis=0)
497 # print data
497 # print data
498 # raise
498 # raise
499 n = self.__profIndex
499 n = self.__profIndex
500
500
501 return data, n
501 return data, n
502
502
503 def byProfiles(self, data):
503 def byProfiles(self, data):
504
504
505 self.__dataReady = False
505 self.__dataReady = False
506 avgdata = None
506 avgdata = None
507 # n = None
507 # n = None
508 # print data
508 # print data
509 # raise
509 # raise
510 self.putData(data)
510 self.putData(data)
511
511
512 if self.__profIndex == self.n:
512 if self.__profIndex == self.n:
513 avgdata, n = self.pushData()
513 avgdata, n = self.pushData()
514 self.__dataReady = True
514 self.__dataReady = True
515
515
516 return avgdata
516 return avgdata
517
517
518 def byTime(self, data, datatime):
518 def byTime(self, data, datatime):
519
519
520 self.__dataReady = False
520 self.__dataReady = False
521 avgdata = None
521 avgdata = None
522 n = None
522 n = None
523
523
524 self.putData(data)
524 self.putData(data)
525
525
526 if (datatime - self.__initime) >= self.__integrationtime:
526 if (datatime - self.__initime) >= self.__integrationtime:
527 avgdata, n = self.pushData()
527 avgdata, n = self.pushData()
528 self.n = n
528 self.n = n
529 self.__dataReady = True
529 self.__dataReady = True
530
530
531 return avgdata
531 return avgdata
532
532
533 def integrateByStride(self, data, datatime):
533 def integrateByStride(self, data, datatime):
534 # print data
534 # print data
535 if self.__profIndex == 0:
535 if self.__profIndex == 0:
536 self.__buffer = [[data.copy(), datatime]]
536 self.__buffer = [[data.copy(), datatime]]
537 else:
537 else:
538 self.__buffer.append([data.copy(),datatime])
538 self.__buffer.append([data.copy(),datatime])
539 self.__profIndex += 1
539 self.__profIndex += 1
540 self.__dataReady = False
540 self.__dataReady = False
541
541
542 if self.__profIndex == self.n * self.stride :
542 if self.__profIndex == self.n * self.stride :
543 self.__dataToPutStride = True
543 self.__dataToPutStride = True
544 self.__profIndexStride = 0
544 self.__profIndexStride = 0
545 self.__profIndex = 0
545 self.__profIndex = 0
546 self.__bufferStride = []
546 self.__bufferStride = []
547 for i in range(self.stride):
547 for i in range(self.stride):
548 current = self.__buffer[i::self.stride]
548 current = self.__buffer[i::self.stride]
549 data = numpy.sum([t[0] for t in current], axis=0)
549 data = numpy.sum([t[0] for t in current], axis=0)
550 avgdatatime = numpy.average([t[1] for t in current])
550 avgdatatime = numpy.average([t[1] for t in current])
551 # print data
551 # print data
552 self.__bufferStride.append((data, avgdatatime))
552 self.__bufferStride.append((data, avgdatatime))
553
553
554 if self.__dataToPutStride:
554 if self.__dataToPutStride:
555 self.__dataReady = True
555 self.__dataReady = True
556 self.__profIndexStride += 1
556 self.__profIndexStride += 1
557 if self.__profIndexStride == self.stride:
557 if self.__profIndexStride == self.stride:
558 self.__dataToPutStride = False
558 self.__dataToPutStride = False
559 # print self.__bufferStride[self.__profIndexStride - 1]
559 # print self.__bufferStride[self.__profIndexStride - 1]
560 # raise
560 # raise
561 return self.__bufferStride[self.__profIndexStride - 1]
561 return self.__bufferStride[self.__profIndexStride - 1]
562
562
563
563
564 return None, None
564 return None, None
565
565
566 def integrate(self, data, datatime=None):
566 def integrate(self, data, datatime=None):
567
567
568 if self.__initime == None:
568 if self.__initime == None:
569 self.__initime = datatime
569 self.__initime = datatime
570
570
571 if self.__byTime:
571 if self.__byTime:
572 avgdata = self.byTime(data, datatime)
572 avgdata = self.byTime(data, datatime)
573 else:
573 else:
574 avgdata = self.byProfiles(data)
574 avgdata = self.byProfiles(data)
575
575
576
576
577 self.__lastdatatime = datatime
577 self.__lastdatatime = datatime
578
578
579 if avgdata is None:
579 if avgdata is None:
580 return None, None
580 return None, None
581
581
582 avgdatatime = self.__initime
582 avgdatatime = self.__initime
583
583
584 deltatime = datatime - self.__lastdatatime
584 deltatime = datatime - self.__lastdatatime
585
585
586 if not self.__withOverlapping:
586 if not self.__withOverlapping:
587 self.__initime = datatime
587 self.__initime = datatime
588 else:
588 else:
589 self.__initime += deltatime
589 self.__initime += deltatime
590
590
591 return avgdata, avgdatatime
591 return avgdata, avgdatatime
592
592
593 def integrateByBlock(self, dataOut):
593 def integrateByBlock(self, dataOut):
594
594
595 times = int(dataOut.data.shape[1]/self.n)
595 times = int(dataOut.data.shape[1]/self.n)
596 avgdata = numpy.zeros((dataOut.nChannels, times, dataOut.nHeights), dtype=numpy.complex)
596 avgdata = numpy.zeros((dataOut.nChannels, times, dataOut.nHeights), dtype=numpy.complex)
597
597
598 id_min = 0
598 id_min = 0
599 id_max = self.n
599 id_max = self.n
600
600
601 for i in range(times):
601 for i in range(times):
602 junk = dataOut.data[:,id_min:id_max,:]
602 junk = dataOut.data[:,id_min:id_max,:]
603 avgdata[:,i,:] = junk.sum(axis=1)
603 avgdata[:,i,:] = junk.sum(axis=1)
604 id_min += self.n
604 id_min += self.n
605 id_max += self.n
605 id_max += self.n
606
606
607 timeInterval = dataOut.ippSeconds*self.n
607 timeInterval = dataOut.ippSeconds*self.n
608 avgdatatime = (times - 1) * timeInterval + dataOut.utctime
608 avgdatatime = (times - 1) * timeInterval + dataOut.utctime
609 self.__dataReady = True
609 self.__dataReady = True
610 return avgdata, avgdatatime
610 return avgdata, avgdatatime
611
611
612 def run(self, dataOut, n=None, timeInterval=None, stride=None, overlapping=False, byblock=False, **kwargs):
612 def run(self, dataOut, n=None, timeInterval=None, stride=None, overlapping=False, byblock=False, **kwargs):
613 if not self.isConfig:
613 if not self.isConfig:
614 self.setup(n=n, stride=stride, timeInterval=timeInterval, overlapping=overlapping, byblock=byblock, **kwargs)
614 self.setup(n=n, stride=stride, timeInterval=timeInterval, overlapping=overlapping, byblock=byblock, **kwargs)
615 self.isConfig = True
615 self.isConfig = True
616
616
617 if dataOut.flagDataAsBlock:
617 if dataOut.flagDataAsBlock:
618 """
618 """
619 Si la data es leida por bloques, dimension = [nChannels, nProfiles, nHeis]
619 Si la data es leida por bloques, dimension = [nChannels, nProfiles, nHeis]
620 """
620 """
621 avgdata, avgdatatime = self.integrateByBlock(dataOut)
621 avgdata, avgdatatime = self.integrateByBlock(dataOut)
622 dataOut.nProfiles /= self.n
622 dataOut.nProfiles /= self.n
623 else:
623 else:
624 if stride is None:
624 if stride is None:
625 avgdata, avgdatatime = self.integrate(dataOut.data, dataOut.utctime)
625 avgdata, avgdatatime = self.integrate(dataOut.data, dataOut.utctime)
626 else:
626 else:
627 avgdata, avgdatatime = self.integrateByStride(dataOut.data, dataOut.utctime)
627 avgdata, avgdatatime = self.integrateByStride(dataOut.data, dataOut.utctime)
628
628
629
629
630 # dataOut.timeInterval *= n
630 # dataOut.timeInterval *= n
631 dataOut.flagNoData = True
631 dataOut.flagNoData = True
632
632
633 if self.__dataReady:
633 if self.__dataReady:
634 dataOut.data = avgdata
634 dataOut.data = avgdata
635 dataOut.nCohInt *= self.n
635 dataOut.nCohInt *= self.n
636 dataOut.utctime = avgdatatime
636 dataOut.utctime = avgdatatime
637 # print avgdata, avgdatatime
637 # print avgdata, avgdatatime
638 # raise
638 # raise
639 # dataOut.timeInterval = dataOut.ippSeconds * dataOut.nCohInt
639 # dataOut.timeInterval = dataOut.ippSeconds * dataOut.nCohInt
640 dataOut.flagNoData = False
640 dataOut.flagNoData = False
641
641
642 class Decoder(Operation):
642 class Decoder(Operation):
643
643
644 isConfig = False
644 isConfig = False
645 __profIndex = 0
645 __profIndex = 0
646
646
647 code = None
647 code = None
648
648
649 nCode = None
649 nCode = None
650 nBaud = None
650 nBaud = None
651
651
652 def __init__(self, **kwargs):
652 def __init__(self, **kwargs):
653
653
654 Operation.__init__(self, **kwargs)
654 Operation.__init__(self, **kwargs)
655
655
656 self.times = None
656 self.times = None
657 self.osamp = None
657 self.osamp = None
658 # self.__setValues = False
658 # self.__setValues = False
659 self.isConfig = False
659 self.isConfig = False
660
660
661 def setup(self, code, osamp, dataOut):
661 def setup(self, code, osamp, dataOut):
662
662
663 self.__profIndex = 0
663 self.__profIndex = 0
664
664
665 self.code = code
665 self.code = code
666
666
667 self.nCode = len(code)
667 self.nCode = len(code)
668 self.nBaud = len(code[0])
668 self.nBaud = len(code[0])
669
669
670 if (osamp != None) and (osamp >1):
670 if (osamp != None) and (osamp >1):
671 self.osamp = osamp
671 self.osamp = osamp
672 self.code = numpy.repeat(code, repeats=self.osamp, axis=1)
672 self.code = numpy.repeat(code, repeats=self.osamp, axis=1)
673 self.nBaud = self.nBaud*self.osamp
673 self.nBaud = self.nBaud*self.osamp
674
674
675 self.__nChannels = dataOut.nChannels
675 self.__nChannels = dataOut.nChannels
676 self.__nProfiles = dataOut.nProfiles
676 self.__nProfiles = dataOut.nProfiles
677 self.__nHeis = dataOut.nHeights
677 self.__nHeis = dataOut.nHeights
678
678
679 if self.__nHeis < self.nBaud:
679 if self.__nHeis < self.nBaud:
680 raise ValueError, 'Number of heights (%d) should be greater than number of bauds (%d)' %(self.__nHeis, self.nBaud)
680 raise ValueError, 'Number of heights (%d) should be greater than number of bauds (%d)' %(self.__nHeis, self.nBaud)
681
681
682 #Frequency
682 #Frequency
683 __codeBuffer = numpy.zeros((self.nCode, self.__nHeis), dtype=numpy.complex)
683 __codeBuffer = numpy.zeros((self.nCode, self.__nHeis), dtype=numpy.complex)
684
684
685 __codeBuffer[:,0:self.nBaud] = self.code
685 __codeBuffer[:,0:self.nBaud] = self.code
686
686
687 self.fft_code = numpy.conj(numpy.fft.fft(__codeBuffer, axis=1))
687 self.fft_code = numpy.conj(numpy.fft.fft(__codeBuffer, axis=1))
688
688
689 if dataOut.flagDataAsBlock:
689 if dataOut.flagDataAsBlock:
690
690
691 self.ndatadec = self.__nHeis #- self.nBaud + 1
691 self.ndatadec = self.__nHeis #- self.nBaud + 1
692
692
693 self.datadecTime = numpy.zeros((self.__nChannels, self.__nProfiles, self.ndatadec), dtype=numpy.complex)
693 self.datadecTime = numpy.zeros((self.__nChannels, self.__nProfiles, self.ndatadec), dtype=numpy.complex)
694
694
695 else:
695 else:
696
696
697 #Time
697 #Time
698 self.ndatadec = self.__nHeis #- self.nBaud + 1
698 self.ndatadec = self.__nHeis #- self.nBaud + 1
699
699
700 self.datadecTime = numpy.zeros((self.__nChannels, self.ndatadec), dtype=numpy.complex)
700 self.datadecTime = numpy.zeros((self.__nChannels, self.ndatadec), dtype=numpy.complex)
701
701
702 def __convolutionInFreq(self, data):
702 def __convolutionInFreq(self, data):
703
703
704 fft_code = self.fft_code[self.__profIndex].reshape(1,-1)
704 fft_code = self.fft_code[self.__profIndex].reshape(1,-1)
705
705
706 fft_data = numpy.fft.fft(data, axis=1)
706 fft_data = numpy.fft.fft(data, axis=1)
707
707
708 conv = fft_data*fft_code
708 conv = fft_data*fft_code
709
709
710 data = numpy.fft.ifft(conv,axis=1)
710 data = numpy.fft.ifft(conv,axis=1)
711
711
712 return data
712 return data
713
713
714 def __convolutionInFreqOpt(self, data):
714 def __convolutionInFreqOpt(self, data):
715
715
716 raise NotImplementedError
716 raise NotImplementedError
717
717
718 def __convolutionInTime(self, data):
718 def __convolutionInTime(self, data):
719
719
720 code = self.code[self.__profIndex]
720 code = self.code[self.__profIndex]
721 for i in range(self.__nChannels):
721 for i in range(self.__nChannels):
722 self.datadecTime[i,:] = numpy.correlate(data[i,:], code, mode='full')[self.nBaud-1:]
722 self.datadecTime[i,:] = numpy.correlate(data[i,:], code, mode='full')[self.nBaud-1:]
723
723
724 return self.datadecTime
724 return self.datadecTime
725
725
726 def __convolutionByBlockInTime(self, data):
726 def __convolutionByBlockInTime(self, data):
727
727
728 repetitions = self.__nProfiles / self.nCode
728 repetitions = self.__nProfiles / self.nCode
729
729
730 junk = numpy.lib.stride_tricks.as_strided(self.code, (repetitions, self.code.size), (0, self.code.itemsize))
730 junk = numpy.lib.stride_tricks.as_strided(self.code, (repetitions, self.code.size), (0, self.code.itemsize))
731 junk = junk.flatten()
731 junk = junk.flatten()
732 code_block = numpy.reshape(junk, (self.nCode*repetitions, self.nBaud))
732 code_block = numpy.reshape(junk, (self.nCode*repetitions, self.nBaud))
733 profilesList = xrange(self.__nProfiles)
733 profilesList = xrange(self.__nProfiles)
734
734
735 for i in range(self.__nChannels):
735 for i in range(self.__nChannels):
736 for j in profilesList:
736 for j in profilesList:
737 self.datadecTime[i,j,:] = numpy.correlate(data[i,j,:], code_block[j,:], mode='full')[self.nBaud-1:]
737 self.datadecTime[i,j,:] = numpy.correlate(data[i,j,:], code_block[j,:], mode='full')[self.nBaud-1:]
738 return self.datadecTime
738 return self.datadecTime
739
739
740 def __convolutionByBlockInFreq(self, data):
740 def __convolutionByBlockInFreq(self, data):
741
741
742 raise NotImplementedError, "Decoder by frequency fro Blocks not implemented"
742 raise NotImplementedError, "Decoder by frequency fro Blocks not implemented"
743
743
744
744
745 fft_code = self.fft_code[self.__profIndex].reshape(1,-1)
745 fft_code = self.fft_code[self.__profIndex].reshape(1,-1)
746
746
747 fft_data = numpy.fft.fft(data, axis=2)
747 fft_data = numpy.fft.fft(data, axis=2)
748
748
749 conv = fft_data*fft_code
749 conv = fft_data*fft_code
750
750
751 data = numpy.fft.ifft(conv,axis=2)
751 data = numpy.fft.ifft(conv,axis=2)
752
752
753 return data
753 return data
754
754
755
755
756 def run(self, dataOut, code=None, nCode=None, nBaud=None, mode = 0, osamp=None, times=None):
756 def run(self, dataOut, code=None, nCode=None, nBaud=None, mode = 0, osamp=None, times=None):
757
757
758 if dataOut.flagDecodeData:
758 if dataOut.flagDecodeData:
759 print "This data is already decoded, recoding again ..."
759 print "This data is already decoded, recoding again ..."
760
760
761 if not self.isConfig:
761 if not self.isConfig:
762
762
763 if code is None:
763 if code is None:
764 if dataOut.code is None:
764 if dataOut.code is None:
765 raise ValueError, "Code could not be read from %s instance. Enter a value in Code parameter" %dataOut.type
765 raise ValueError, "Code could not be read from %s instance. Enter a value in Code parameter" %dataOut.type
766
766
767 code = dataOut.code
767 code = dataOut.code
768 else:
768 else:
769 code = numpy.array(code).reshape(nCode,nBaud)
769 code = numpy.array(code).reshape(nCode,nBaud)
770 self.setup(code, osamp, dataOut)
770 self.setup(code, osamp, dataOut)
771
771
772 self.isConfig = True
772 self.isConfig = True
773
773
774 if mode == 3:
774 if mode == 3:
775 sys.stderr.write("Decoder Warning: mode=%d is not valid, using mode=0\n" %mode)
775 sys.stderr.write("Decoder Warning: mode=%d is not valid, using mode=0\n" %mode)
776
776
777 if times != None:
777 if times != None:
778 sys.stderr.write("Decoder Warning: Argument 'times' in not used anymore\n")
778 sys.stderr.write("Decoder Warning: Argument 'times' in not used anymore\n")
779
779
780 if self.code is None:
780 if self.code is None:
781 print "Fail decoding: Code is not defined."
781 print "Fail decoding: Code is not defined."
782 return
782 return
783
783
784 self.__nProfiles = dataOut.nProfiles
784 self.__nProfiles = dataOut.nProfiles
785 datadec = None
785 datadec = None
786
786
787 if mode == 3:
787 if mode == 3:
788 mode = 0
788 mode = 0
789
789
790 if dataOut.flagDataAsBlock:
790 if dataOut.flagDataAsBlock:
791 """
791 """
792 Decoding when data have been read as block,
792 Decoding when data have been read as block,
793 """
793 """
794
794
795 if mode == 0:
795 if mode == 0:
796 datadec = self.__convolutionByBlockInTime(dataOut.data)
796 datadec = self.__convolutionByBlockInTime(dataOut.data)
797 if mode == 1:
797 if mode == 1:
798 datadec = self.__convolutionByBlockInFreq(dataOut.data)
798 datadec = self.__convolutionByBlockInFreq(dataOut.data)
799 else:
799 else:
800 """
800 """
801 Decoding when data have been read profile by profile
801 Decoding when data have been read profile by profile
802 """
802 """
803 if mode == 0:
803 if mode == 0:
804 datadec = self.__convolutionInTime(dataOut.data)
804 datadec = self.__convolutionInTime(dataOut.data)
805
805
806 if mode == 1:
806 if mode == 1:
807 datadec = self.__convolutionInFreq(dataOut.data)
807 datadec = self.__convolutionInFreq(dataOut.data)
808
808
809 if mode == 2:
809 if mode == 2:
810 datadec = self.__convolutionInFreqOpt(dataOut.data)
810 datadec = self.__convolutionInFreqOpt(dataOut.data)
811
811
812 if datadec is None:
812 if datadec is None:
813 raise ValueError, "Codification mode selected is not valid: mode=%d. Try selecting 0 or 1" %mode
813 raise ValueError, "Codification mode selected is not valid: mode=%d. Try selecting 0 or 1" %mode
814
814
815 dataOut.code = self.code
815 dataOut.code = self.code
816 dataOut.nCode = self.nCode
816 dataOut.nCode = self.nCode
817 dataOut.nBaud = self.nBaud
817 dataOut.nBaud = self.nBaud
818
818
819 dataOut.data = datadec
819 dataOut.data = datadec
820
820
821 dataOut.heightList = dataOut.heightList[0:datadec.shape[-1]]
821 dataOut.heightList = dataOut.heightList[0:datadec.shape[-1]]
822
822
823 dataOut.flagDecodeData = True #asumo q la data esta decodificada
823 dataOut.flagDecodeData = True #asumo q la data esta decodificada
824
824
825 if self.__profIndex == self.nCode-1:
825 if self.__profIndex == self.nCode-1:
826 self.__profIndex = 0
826 self.__profIndex = 0
827 return 1
827 return 1
828
828
829 self.__profIndex += 1
829 self.__profIndex += 1
830
830
831 return 1
831 return 1
832 # dataOut.flagDeflipData = True #asumo q la data no esta sin flip
832 # dataOut.flagDeflipData = True #asumo q la data no esta sin flip
833
833
834
834
835 class ProfileConcat(Operation):
835 class ProfileConcat(Operation):
836
836
837 isConfig = False
837 isConfig = False
838 buffer = None
838 buffer = None
839 concat_m =None
839 concat_m =None
840
840
841 def __init__(self, **kwargs):
841 def __init__(self, **kwargs):
842
842
843 Operation.__init__(self, **kwargs)
843 Operation.__init__(self, **kwargs)
844 self.profileIndex = 0
844 self.profileIndex = 0
845
845
846 def reset(self):
846 def reset(self):
847 self.buffer = numpy.zeros_like(self.buffer)
847 self.buffer = numpy.zeros_like(self.buffer)
848 self.start_index = 0
848 self.start_index = 0
849 self.times = 1
849 self.times = 1
850
850
851 def setup(self, data, m, n=1):
851 def setup(self, data, m, n=1):
852 self.buffer = numpy.zeros((data.shape[0],data.shape[1]*m),dtype=type(data[0,0]))
852 self.buffer = numpy.zeros((data.shape[0],data.shape[1]*m),dtype=type(data[0,0]))
853 self.nHeights = data.shape[1]#.nHeights
853 self.nHeights = data.shape[1]#.nHeights
854 self.start_index = 0
854 self.start_index = 0
855 self.times = 1
855 self.times = 1
856
856
857 def concat(self, data):
857 def concat(self, data):
858
858
859 self.buffer[:,self.start_index:self.nHeights*self.times] = data.copy()
859 self.buffer[:,self.start_index:self.nHeights*self.times] = data.copy()
860 self.start_index = self.start_index + self.nHeights
860 self.start_index = self.start_index + self.nHeights
861
861
862 def run(self, dataOut, m):
862 def run(self, dataOut, m):
863
863
864 self.concat_m= m
864 self.concat_m= m
865 dataOut.flagNoData = True
865 dataOut.flagNoData = True
866
866
867 if not self.isConfig:
867 if not self.isConfig:
868 self.setup(dataOut.data, m, 1)
868 self.setup(dataOut.data, m, 1)
869 self.isConfig = True
869 self.isConfig = True
870
870
871 if dataOut.flagDataAsBlock:
871 if dataOut.flagDataAsBlock:
872 raise ValueError, "ProfileConcat can only be used when voltage have been read profile by profile, getBlock = False"
872 raise ValueError, "ProfileConcat can only be used when voltage have been read profile by profile, getBlock = False"
873
873
874 else:
874 else:
875 self.concat(dataOut.data)
875 self.concat(dataOut.data)
876 self.times += 1
876 self.times += 1
877 if self.times > m:
877 if self.times > m:
878 dataOut.data = self.buffer
878 dataOut.data = self.buffer
879 self.reset()
879 self.reset()
880 dataOut.flagNoData = False
880 dataOut.flagNoData = False
881 # se deben actualizar mas propiedades del header y del objeto dataOut, por ejemplo, las alturas
881 # se deben actualizar mas propiedades del header y del objeto dataOut, por ejemplo, las alturas
882 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
882 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
883 xf = dataOut.heightList[0] + dataOut.nHeights * deltaHeight * m
883 xf = dataOut.heightList[0] + dataOut.nHeights * deltaHeight * m
884 dataOut.heightList = numpy.arange(dataOut.heightList[0], xf, deltaHeight)
884 dataOut.heightList = numpy.arange(dataOut.heightList[0], xf, deltaHeight)
885 dataOut.ippSeconds *= m
885 dataOut.ippSeconds *= m
886 dataOut.concat_m = int(m)
886 dataOut.concat_m = int(m)
887
887
888 class ProfileSelector(Operation):
888 class ProfileSelector(Operation):
889
889
890 profileIndex = None
890 profileIndex = None
891 # Tamanho total de los perfiles
891 # Tamanho total de los perfiles
892 nProfiles = None
892 nProfiles = None
893
893
894 def __init__(self, **kwargs):
894 def __init__(self, **kwargs):
895
895
896 Operation.__init__(self, **kwargs)
896 Operation.__init__(self, **kwargs)
897 self.profileIndex = 0
897 self.profileIndex = 0
898
898
899 def incProfileIndex(self):
899 def incProfileIndex(self):
900
900
901 self.profileIndex += 1
901 self.profileIndex += 1
902
902
903 if self.profileIndex >= self.nProfiles:
903 if self.profileIndex >= self.nProfiles:
904 self.profileIndex = 0
904 self.profileIndex = 0
905
905
906 def isThisProfileInRange(self, profileIndex, minIndex, maxIndex):
906 def isThisProfileInRange(self, profileIndex, minIndex, maxIndex):
907
907
908 if profileIndex < minIndex:
908 if profileIndex < minIndex:
909 return False
909 return False
910
910
911 if profileIndex > maxIndex:
911 if profileIndex > maxIndex:
912 return False
912 return False
913
913
914 return True
914 return True
915
915
916 def isThisProfileInList(self, profileIndex, profileList):
916 def isThisProfileInList(self, profileIndex, profileList):
917
917
918 if profileIndex not in profileList:
918 if profileIndex not in profileList:
919 return False
919 return False
920
920
921 return True
921 return True
922
922
923 def run(self, dataOut, profileList=None, profileRangeList=None, beam=None, byblock=False, rangeList = None, nProfiles=None):
923 def run(self, dataOut, profileList=None, profileRangeList=None, beam=None, byblock=False, rangeList = None, nProfiles=None):
924
924
925 """
925 """
926 ProfileSelector:
926 ProfileSelector:
927
927
928 Inputs:
928 Inputs:
929 profileList : Index of profiles selected. Example: profileList = (0,1,2,7,8)
929 profileList : Index of profiles selected. Example: profileList = (0,1,2,7,8)
930
930
931 profileRangeList : Minimum and maximum profile indexes. Example: profileRangeList = (4, 30)
931 profileRangeList : Minimum and maximum profile indexes. Example: profileRangeList = (4, 30)
932
932
933 rangeList : List of profile ranges. Example: rangeList = ((4, 30), (32, 64), (128, 256))
933 rangeList : List of profile ranges. Example: rangeList = ((4, 30), (32, 64), (128, 256))
934
934
935 """
935 """
936
936
937 if rangeList is not None:
937 if rangeList is not None:
938 if type(rangeList[0]) not in (tuple, list):
938 if type(rangeList[0]) not in (tuple, list):
939 rangeList = [rangeList]
939 rangeList = [rangeList]
940
940
941 dataOut.flagNoData = True
941 dataOut.flagNoData = True
942
942
943 if dataOut.flagDataAsBlock:
943 if dataOut.flagDataAsBlock:
944 """
944 """
945 data dimension = [nChannels, nProfiles, nHeis]
945 data dimension = [nChannels, nProfiles, nHeis]
946 """
946 """
947 if profileList != None:
947 if profileList != None:
948 dataOut.data = dataOut.data[:,profileList,:]
948 dataOut.data = dataOut.data[:,profileList,:]
949
949
950 if profileRangeList != None:
950 if profileRangeList != None:
951 minIndex = profileRangeList[0]
951 minIndex = profileRangeList[0]
952 maxIndex = profileRangeList[1]
952 maxIndex = profileRangeList[1]
953 profileList = range(minIndex, maxIndex+1)
953 profileList = range(minIndex, maxIndex+1)
954
954
955 dataOut.data = dataOut.data[:,minIndex:maxIndex+1,:]
955 dataOut.data = dataOut.data[:,minIndex:maxIndex+1,:]
956
956
957 if rangeList != None:
957 if rangeList != None:
958
958
959 profileList = []
959 profileList = []
960
960
961 for thisRange in rangeList:
961 for thisRange in rangeList:
962 minIndex = thisRange[0]
962 minIndex = thisRange[0]
963 maxIndex = thisRange[1]
963 maxIndex = thisRange[1]
964
964
965 profileList.extend(range(minIndex, maxIndex+1))
965 profileList.extend(range(minIndex, maxIndex+1))
966
966
967 dataOut.data = dataOut.data[:,profileList,:]
967 dataOut.data = dataOut.data[:,profileList,:]
968
968
969 dataOut.nProfiles = len(profileList)
969 dataOut.nProfiles = len(profileList)
970 dataOut.profileIndex = dataOut.nProfiles - 1
970 dataOut.profileIndex = dataOut.nProfiles - 1
971 dataOut.flagNoData = False
971 dataOut.flagNoData = False
972
972
973 return True
973 return True
974
974
975 """
975 """
976 data dimension = [nChannels, nHeis]
976 data dimension = [nChannels, nHeis]
977 """
977 """
978
978
979 if profileList != None:
979 if profileList != None:
980
980
981 if self.isThisProfileInList(dataOut.profileIndex, profileList):
981 if self.isThisProfileInList(dataOut.profileIndex, profileList):
982
982
983 self.nProfiles = len(profileList)
983 self.nProfiles = len(profileList)
984 dataOut.nProfiles = self.nProfiles
984 dataOut.nProfiles = self.nProfiles
985 dataOut.profileIndex = self.profileIndex
985 dataOut.profileIndex = self.profileIndex
986 dataOut.flagNoData = False
986 dataOut.flagNoData = False
987
987
988 self.incProfileIndex()
988 self.incProfileIndex()
989 return True
989 return True
990
990
991 if profileRangeList != None:
991 if profileRangeList != None:
992
992
993 minIndex = profileRangeList[0]
993 minIndex = profileRangeList[0]
994 maxIndex = profileRangeList[1]
994 maxIndex = profileRangeList[1]
995
995
996 if self.isThisProfileInRange(dataOut.profileIndex, minIndex, maxIndex):
996 if self.isThisProfileInRange(dataOut.profileIndex, minIndex, maxIndex):
997
997
998 self.nProfiles = maxIndex - minIndex + 1
998 self.nProfiles = maxIndex - minIndex + 1
999 dataOut.nProfiles = self.nProfiles
999 dataOut.nProfiles = self.nProfiles
1000 dataOut.profileIndex = self.profileIndex
1000 dataOut.profileIndex = self.profileIndex
1001 dataOut.flagNoData = False
1001 dataOut.flagNoData = False
1002
1002
1003 self.incProfileIndex()
1003 self.incProfileIndex()
1004 return True
1004 return True
1005
1005
1006 if rangeList != None:
1006 if rangeList != None:
1007
1007
1008 nProfiles = 0
1008 nProfiles = 0
1009
1009
1010 for thisRange in rangeList:
1010 for thisRange in rangeList:
1011 minIndex = thisRange[0]
1011 minIndex = thisRange[0]
1012 maxIndex = thisRange[1]
1012 maxIndex = thisRange[1]
1013
1013
1014 nProfiles += maxIndex - minIndex + 1
1014 nProfiles += maxIndex - minIndex + 1
1015
1015
1016 for thisRange in rangeList:
1016 for thisRange in rangeList:
1017
1017
1018 minIndex = thisRange[0]
1018 minIndex = thisRange[0]
1019 maxIndex = thisRange[1]
1019 maxIndex = thisRange[1]
1020
1020
1021 if self.isThisProfileInRange(dataOut.profileIndex, minIndex, maxIndex):
1021 if self.isThisProfileInRange(dataOut.profileIndex, minIndex, maxIndex):
1022
1022
1023 self.nProfiles = nProfiles
1023 self.nProfiles = nProfiles
1024 dataOut.nProfiles = self.nProfiles
1024 dataOut.nProfiles = self.nProfiles
1025 dataOut.profileIndex = self.profileIndex
1025 dataOut.profileIndex = self.profileIndex
1026 dataOut.flagNoData = False
1026 dataOut.flagNoData = False
1027
1027
1028 self.incProfileIndex()
1028 self.incProfileIndex()
1029
1029
1030 break
1030 break
1031
1031
1032 return True
1032 return True
1033
1033
1034
1034
1035 if beam != None: #beam is only for AMISR data
1035 if beam != None: #beam is only for AMISR data
1036 if self.isThisProfileInList(dataOut.profileIndex, dataOut.beamRangeDict[beam]):
1036 if self.isThisProfileInList(dataOut.profileIndex, dataOut.beamRangeDict[beam]):
1037 dataOut.flagNoData = False
1037 dataOut.flagNoData = False
1038 dataOut.profileIndex = self.profileIndex
1038 dataOut.profileIndex = self.profileIndex
1039
1039
1040 self.incProfileIndex()
1040 self.incProfileIndex()
1041
1041
1042 return True
1042 return True
1043
1043
1044 raise ValueError, "ProfileSelector needs profileList, profileRangeList or rangeList parameter"
1044 raise ValueError, "ProfileSelector needs profileList, profileRangeList or rangeList parameter"
1045
1045
1046 return False
1046 return False
1047
1047
1048 class Reshaper(Operation):
1048 class Reshaper(Operation):
1049
1049
1050 def __init__(self, **kwargs):
1050 def __init__(self, **kwargs):
1051
1051
1052 Operation.__init__(self, **kwargs)
1052 Operation.__init__(self, **kwargs)
1053
1053
1054 self.__buffer = None
1054 self.__buffer = None
1055 self.__nitems = 0
1055 self.__nitems = 0
1056
1056
1057 def __appendProfile(self, dataOut, nTxs):
1057 def __appendProfile(self, dataOut, nTxs):
1058
1058
1059 if self.__buffer is None:
1059 if self.__buffer is None:
1060 shape = (dataOut.nChannels, int(dataOut.nHeights/nTxs) )
1060 shape = (dataOut.nChannels, int(dataOut.nHeights/nTxs) )
1061 self.__buffer = numpy.empty(shape, dtype = dataOut.data.dtype)
1061 self.__buffer = numpy.empty(shape, dtype = dataOut.data.dtype)
1062
1062
1063 ini = dataOut.nHeights * self.__nitems
1063 ini = dataOut.nHeights * self.__nitems
1064 end = ini + dataOut.nHeights
1064 end = ini + dataOut.nHeights
1065
1065
1066 self.__buffer[:, ini:end] = dataOut.data
1066 self.__buffer[:, ini:end] = dataOut.data
1067
1067
1068 self.__nitems += 1
1068 self.__nitems += 1
1069
1069
1070 return int(self.__nitems*nTxs)
1070 return int(self.__nitems*nTxs)
1071
1071
1072 def __getBuffer(self):
1072 def __getBuffer(self):
1073
1073
1074 if self.__nitems == int(1./self.__nTxs):
1074 if self.__nitems == int(1./self.__nTxs):
1075
1075
1076 self.__nitems = 0
1076 self.__nitems = 0
1077
1077
1078 return self.__buffer.copy()
1078 return self.__buffer.copy()
1079
1079
1080 return None
1080 return None
1081
1081
1082 def __checkInputs(self, dataOut, shape, nTxs):
1082 def __checkInputs(self, dataOut, shape, nTxs):
1083
1083
1084 if shape is None and nTxs is None:
1084 if shape is None and nTxs is None:
1085 raise ValueError, "Reshaper: shape of factor should be defined"
1085 raise ValueError, "Reshaper: shape of factor should be defined"
1086
1086
1087 if nTxs:
1087 if nTxs:
1088 if nTxs < 0:
1088 if nTxs < 0:
1089 raise ValueError, "nTxs should be greater than 0"
1089 raise ValueError, "nTxs should be greater than 0"
1090
1090
1091 if nTxs < 1 and dataOut.nProfiles % (1./nTxs) != 0:
1091 if nTxs < 1 and dataOut.nProfiles % (1./nTxs) != 0:
1092 raise ValueError, "nProfiles= %d is not divisibled by (1./nTxs) = %f" %(dataOut.nProfiles, (1./nTxs))
1092 raise ValueError, "nProfiles= %d is not divisibled by (1./nTxs) = %f" %(dataOut.nProfiles, (1./nTxs))
1093
1093
1094 shape = [dataOut.nChannels, dataOut.nProfiles*nTxs, dataOut.nHeights/nTxs]
1094 shape = [dataOut.nChannels, dataOut.nProfiles*nTxs, dataOut.nHeights/nTxs]
1095
1095
1096 return shape, nTxs
1096 return shape, nTxs
1097
1097
1098 if len(shape) != 2 and len(shape) != 3:
1098 if len(shape) != 2 and len(shape) != 3:
1099 raise ValueError, "shape dimension should be equal to 2 or 3. shape = (nProfiles, nHeis) or (nChannels, nProfiles, nHeis). Actually shape = (%d, %d, %d)" %(dataOut.nChannels, dataOut.nProfiles, dataOut.nHeights)
1099 raise ValueError, "shape dimension should be equal to 2 or 3. shape = (nProfiles, nHeis) or (nChannels, nProfiles, nHeis). Actually shape = (%d, %d, %d)" %(dataOut.nChannels, dataOut.nProfiles, dataOut.nHeights)
1100
1100
1101 if len(shape) == 2:
1101 if len(shape) == 2:
1102 shape_tuple = [dataOut.nChannels]
1102 shape_tuple = [dataOut.nChannels]
1103 shape_tuple.extend(shape)
1103 shape_tuple.extend(shape)
1104 else:
1104 else:
1105 shape_tuple = list(shape)
1105 shape_tuple = list(shape)
1106
1106
1107 nTxs = 1.0*shape_tuple[1]/dataOut.nProfiles
1107 nTxs = 1.0*shape_tuple[1]/dataOut.nProfiles
1108
1108
1109 return shape_tuple, nTxs
1109 return shape_tuple, nTxs
1110
1110
1111 def run(self, dataOut, shape=None, nTxs=None):
1111 def run(self, dataOut, shape=None, nTxs=None):
1112
1112
1113 shape_tuple, self.__nTxs = self.__checkInputs(dataOut, shape, nTxs)
1113 shape_tuple, self.__nTxs = self.__checkInputs(dataOut, shape, nTxs)
1114
1114
1115 dataOut.flagNoData = True
1115 dataOut.flagNoData = True
1116 profileIndex = None
1116 profileIndex = None
1117
1117
1118 if dataOut.flagDataAsBlock:
1118 if dataOut.flagDataAsBlock:
1119
1119
1120 dataOut.data = numpy.reshape(dataOut.data, shape_tuple)
1120 dataOut.data = numpy.reshape(dataOut.data, shape_tuple)
1121 dataOut.flagNoData = False
1121 dataOut.flagNoData = False
1122
1122
1123 profileIndex = int(dataOut.nProfiles*self.__nTxs) - 1
1123 profileIndex = int(dataOut.nProfiles*self.__nTxs) - 1
1124
1124
1125 else:
1125 else:
1126
1126
1127 if self.__nTxs < 1:
1127 if self.__nTxs < 1:
1128
1128
1129 self.__appendProfile(dataOut, self.__nTxs)
1129 self.__appendProfile(dataOut, self.__nTxs)
1130 new_data = self.__getBuffer()
1130 new_data = self.__getBuffer()
1131
1131
1132 if new_data is not None:
1132 if new_data is not None:
1133 dataOut.data = new_data
1133 dataOut.data = new_data
1134 dataOut.flagNoData = False
1134 dataOut.flagNoData = False
1135
1135
1136 profileIndex = dataOut.profileIndex*nTxs
1136 profileIndex = dataOut.profileIndex*nTxs
1137
1137
1138 else:
1138 else:
1139 raise ValueError, "nTxs should be greater than 0 and lower than 1, or use VoltageReader(..., getblock=True)"
1139 raise ValueError, "nTxs should be greater than 0 and lower than 1, or use VoltageReader(..., getblock=True)"
1140
1140
1141 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1141 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1142
1142
1143 dataOut.heightList = numpy.arange(dataOut.nHeights/self.__nTxs) * deltaHeight + dataOut.heightList[0]
1143 dataOut.heightList = numpy.arange(dataOut.nHeights/self.__nTxs) * deltaHeight + dataOut.heightList[0]
1144
1144
1145 dataOut.nProfiles = int(dataOut.nProfiles*self.__nTxs)
1145 dataOut.nProfiles = int(dataOut.nProfiles*self.__nTxs)
1146
1146
1147 dataOut.profileIndex = profileIndex
1147 dataOut.profileIndex = profileIndex
1148
1148
1149 dataOut.ippSeconds /= self.__nTxs
1149 dataOut.ippSeconds /= self.__nTxs
1150
1150
1151 class SplitProfiles(Operation):
1151 class SplitProfiles(Operation):
1152
1152
1153 def __init__(self, **kwargs):
1153 def __init__(self, **kwargs):
1154
1154
1155 Operation.__init__(self, **kwargs)
1155 Operation.__init__(self, **kwargs)
1156
1156
1157 def run(self, dataOut, n):
1157 def run(self, dataOut, n):
1158
1158
1159 dataOut.flagNoData = True
1159 dataOut.flagNoData = True
1160 profileIndex = None
1160 profileIndex = None
1161
1161
1162 if dataOut.flagDataAsBlock:
1162 if dataOut.flagDataAsBlock:
1163
1163
1164 #nchannels, nprofiles, nsamples
1164 #nchannels, nprofiles, nsamples
1165 shape = dataOut.data.shape
1165 shape = dataOut.data.shape
1166
1166
1167 if shape[2] % n != 0:
1167 if shape[2] % n != 0:
1168 raise ValueError, "Could not split the data, n=%d has to be multiple of %d" %(n, shape[2])
1168 raise ValueError, "Could not split the data, n=%d has to be multiple of %d" %(n, shape[2])
1169
1169
1170 new_shape = shape[0], shape[1]*n, shape[2]/n
1170 new_shape = shape[0], shape[1]*n, shape[2]/n
1171
1171
1172 dataOut.data = numpy.reshape(dataOut.data, new_shape)
1172 dataOut.data = numpy.reshape(dataOut.data, new_shape)
1173 dataOut.flagNoData = False
1173 dataOut.flagNoData = False
1174
1174
1175 profileIndex = int(dataOut.nProfiles/n) - 1
1175 profileIndex = int(dataOut.nProfiles/n) - 1
1176
1176
1177 else:
1177 else:
1178
1178
1179 raise ValueError, "Could not split the data when is read Profile by Profile. Use VoltageReader(..., getblock=True)"
1179 raise ValueError, "Could not split the data when is read Profile by Profile. Use VoltageReader(..., getblock=True)"
1180
1180
1181 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1181 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1182
1182
1183 dataOut.heightList = numpy.arange(dataOut.nHeights/n) * deltaHeight + dataOut.heightList[0]
1183 dataOut.heightList = numpy.arange(dataOut.nHeights/n) * deltaHeight + dataOut.heightList[0]
1184
1184
1185 dataOut.nProfiles = int(dataOut.nProfiles*n)
1185 dataOut.nProfiles = int(dataOut.nProfiles*n)
1186
1186
1187 dataOut.profileIndex = profileIndex
1187 dataOut.profileIndex = profileIndex
1188
1188
1189 dataOut.ippSeconds /= n
1189 dataOut.ippSeconds /= n
1190
1190
1191 class CombineProfiles(Operation):
1191 class CombineProfiles(Operation):
1192
1192
1193 def __init__(self, **kwargs):
1193 def __init__(self, **kwargs):
1194
1194
1195 Operation.__init__(self, **kwargs)
1195 Operation.__init__(self, **kwargs)
1196
1196
1197 self.__remData = None
1197 self.__remData = None
1198 self.__profileIndex = 0
1198 self.__profileIndex = 0
1199
1199
1200 def run(self, dataOut, n):
1200 def run(self, dataOut, n):
1201
1201
1202 dataOut.flagNoData = True
1202 dataOut.flagNoData = True
1203 profileIndex = None
1203 profileIndex = None
1204
1204
1205 if dataOut.flagDataAsBlock:
1205 if dataOut.flagDataAsBlock:
1206
1206
1207 #nchannels, nprofiles, nsamples
1207 #nchannels, nprofiles, nsamples
1208 shape = dataOut.data.shape
1208 shape = dataOut.data.shape
1209 new_shape = shape[0], shape[1]/n, shape[2]*n
1209 new_shape = shape[0], shape[1]/n, shape[2]*n
1210
1210
1211 if shape[1] % n != 0:
1211 if shape[1] % n != 0:
1212 raise ValueError, "Could not split the data, n=%d has to be multiple of %d" %(n, shape[1])
1212 raise ValueError, "Could not split the data, n=%d has to be multiple of %d" %(n, shape[1])
1213
1213
1214 dataOut.data = numpy.reshape(dataOut.data, new_shape)
1214 dataOut.data = numpy.reshape(dataOut.data, new_shape)
1215 dataOut.flagNoData = False
1215 dataOut.flagNoData = False
1216
1216
1217 profileIndex = int(dataOut.nProfiles*n) - 1
1217 profileIndex = int(dataOut.nProfiles*n) - 1
1218
1218
1219 else:
1219 else:
1220
1220
1221 #nchannels, nsamples
1221 #nchannels, nsamples
1222 if self.__remData is None:
1222 if self.__remData is None:
1223 newData = dataOut.data
1223 newData = dataOut.data
1224 else:
1224 else:
1225 newData = numpy.concatenate((self.__remData, dataOut.data), axis=1)
1225 newData = numpy.concatenate((self.__remData, dataOut.data), axis=1)
1226
1226
1227 self.__profileIndex += 1
1227 self.__profileIndex += 1
1228
1228
1229 if self.__profileIndex < n:
1229 if self.__profileIndex < n:
1230 self.__remData = newData
1230 self.__remData = newData
1231 #continue
1231 #continue
1232 return
1232 return
1233
1233
1234 self.__profileIndex = 0
1234 self.__profileIndex = 0
1235 self.__remData = None
1235 self.__remData = None
1236
1236
1237 dataOut.data = newData
1237 dataOut.data = newData
1238 dataOut.flagNoData = False
1238 dataOut.flagNoData = False
1239
1239
1240 profileIndex = dataOut.profileIndex/n
1240 profileIndex = dataOut.profileIndex/n
1241
1241
1242
1242
1243 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1243 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1244
1244
1245 dataOut.heightList = numpy.arange(dataOut.nHeights*n) * deltaHeight + dataOut.heightList[0]
1245 dataOut.heightList = numpy.arange(dataOut.nHeights*n) * deltaHeight + dataOut.heightList[0]
1246
1246
1247 dataOut.nProfiles = int(dataOut.nProfiles/n)
1247 dataOut.nProfiles = int(dataOut.nProfiles/n)
1248
1248
1249 dataOut.profileIndex = profileIndex
1249 dataOut.profileIndex = profileIndex
1250
1250
1251 dataOut.ippSeconds *= n
1251 dataOut.ippSeconds *= n
1252
1252
1253
1253
1254 class SSheightProfiles(Operation):
1254 class SSheightProfiles(Operation):
1255
1255
1256 step = None
1256 step = None
1257 nsamples = None
1257 nsamples = None
1258 bufferShape = None
1258 bufferShape = None
1259 profileShape = None
1259 profileShape = None
1260 sshProfiles = None
1260 sshProfiles = None
1261 profileIndex = None
1261 profileIndex = None
1262
1262
1263 def __init__(self, **kwargs):
1263 def __init__(self, **kwargs):
1264
1264
1265 Operation.__init__(self, **kwargs)
1265 Operation.__init__(self, **kwargs)
1266 self.isConfig = False
1266 self.isConfig = False
1267
1267
1268 def setup(self,dataOut ,step = None , nsamples = None):
1268 def setup(self,dataOut ,step = None , nsamples = None):
1269
1269
1270 if step == None and nsamples == None:
1270 if step == None and nsamples == None:
1271 raise ValueError, "step or nheights should be specified ..."
1271 raise ValueError, "step or nheights should be specified ..."
1272
1272
1273 self.step = step
1273 self.step = step
1274 self.nsamples = nsamples
1274 self.nsamples = nsamples
1275 self.__nChannels = dataOut.nChannels
1275 self.__nChannels = dataOut.nChannels
1276 self.__nProfiles = dataOut.nProfiles
1276 self.__nProfiles = dataOut.nProfiles
1277 self.__nHeis = dataOut.nHeights
1277 self.__nHeis = dataOut.nHeights
1278 shape = dataOut.data.shape #nchannels, nprofiles, nsamples
1278 shape = dataOut.data.shape #nchannels, nprofiles, nsamples
1279 print "input nChannels",self.__nChannels
1280 print "input nProfiles",self.__nProfiles
1281 print "input nHeis",self.__nHeis
1282 print "input Shape",shape
1283
1279
1284
1280
1285
1281 residue = (shape[1] - self.nsamples) % self.step
1286 residue = (shape[1] - self.nsamples) % self.step
1282 if residue != 0:
1287 if residue != 0:
1283 print "The residue is %d, step=%d should be multiple of %d to avoid loss of %d samples"%(residue,step,shape[1] - self.nsamples,residue)
1288 print "The residue is %d, step=%d should be multiple of %d to avoid loss of %d samples"%(residue,step,shape[1] - self.nsamples,residue)
1284
1289
1285 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1290 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1286 numberProfile = self.nsamples
1291 numberProfile = self.nsamples
1287 numberSamples = (shape[1] - self.nsamples)/self.step
1292 numberSamples = (shape[1] - self.nsamples)/self.step
1288
1293
1289 print "New number of profile: %d, number of height: %d, Resolution %d Km"%(numberProfile,numberSamples,deltaHeight*self.step)
1294 print "new numberProfile",numberProfile
1295 print "new numberSamples",numberSamples
1296
1297 print "New number of profile: %d, number of height: %d, Resolution %f Km"%(numberProfile,numberSamples,deltaHeight*self.step)
1290
1298
1291 self.bufferShape = shape[0], numberSamples, numberProfile # nchannels, nsamples , nprofiles
1299 self.bufferShape = shape[0], numberSamples, numberProfile # nchannels, nsamples , nprofiles
1292 self.profileShape = shape[0], numberProfile, numberSamples # nchannels, nprofiles, nsamples
1300 self.profileShape = shape[0], numberProfile, numberSamples # nchannels, nprofiles, nsamples
1293
1301
1294 self.buffer = numpy.zeros(self.bufferShape , dtype=numpy.complex)
1302 self.buffer = numpy.zeros(self.bufferShape , dtype=numpy.complex)
1295 self.sshProfiles = numpy.zeros(self.profileShape, dtype=numpy.complex)
1303 self.sshProfiles = numpy.zeros(self.profileShape, dtype=numpy.complex)
1296
1304
1297 def run(self, dataOut, step, nsamples):
1305 def run(self, dataOut, step, nsamples):
1298
1306
1299 dataOut.flagNoData = True
1307 dataOut.flagNoData = True
1300 dataOut.flagDataAsBlock = False
1308 dataOut.flagDataAsBlock = False
1301 profileIndex = None
1309 profileIndex = None
1302
1310
1303
1311
1304 if not self.isConfig:
1312 if not self.isConfig:
1305 self.setup(dataOut, step=step , nsamples=nsamples)
1313 self.setup(dataOut, step=step , nsamples=nsamples)
1306 self.isConfig = True
1314 self.isConfig = True
1307
1315
1308 for i in range(self.buffer.shape[1]):
1316 for i in range(self.buffer.shape[1]):
1309 self.buffer[:,i] = numpy.flip(dataOut.data[:,i*self.step:i*self.step + self.nsamples])
1317 #self.buffer[:,i] = numpy.flip(dataOut.data[:,i*self.step:i*self.step + self.nsamples])
1318 self.buffer[:,i] = dataOut.data[:,i*self.step:i*self.step + self.nsamples]
1319
1310 #self.buffer[:,j,self.__nHeis-j*self.step - self.nheights:self.__nHeis-j*self.step] = numpy.flip(dataOut.data[:,j*self.step:j*self.step + self.nheights])
1320 #self.buffer[:,j,self.__nHeis-j*self.step - self.nheights:self.__nHeis-j*self.step] = numpy.flip(dataOut.data[:,j*self.step:j*self.step + self.nheights])
1311
1321
1312 for j in range(self.buffer.shape[0]):
1322 for j in range(self.buffer.shape[0]):
1313 self.sshProfiles[j] = numpy.transpose(self.buffer[j])
1323 self.sshProfiles[j] = numpy.transpose(self.buffer[j])
1314
1324
1315 profileIndex = self.nsamples
1325 profileIndex = self.nsamples
1316 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1326 deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1317 ippSeconds = (deltaHeight*1.0e-6)/(0.15)
1327 ippSeconds = (deltaHeight*1.0e-6)/(0.15)
1318 try:
1328 try:
1319 if dataOut.concat_m is not None:
1329 if dataOut.concat_m is not None:
1320 ippSeconds= ippSeconds/float(dataOut.concat_m)
1330 ippSeconds= ippSeconds/float(dataOut.concat_m)
1321 #print "Profile concat %d"%dataOut.concat_m
1331 #print "Profile concat %d"%dataOut.concat_m
1322 except:
1332 except:
1323 pass
1333 pass
1324
1334
1325 dataOut.data = self.sshProfiles
1335 dataOut.data = self.sshProfiles
1326 dataOut.flagNoData = False
1336 dataOut.flagNoData = False
1327 dataOut.heightList = numpy.arange(self.buffer.shape[1]) *self.step*deltaHeight + dataOut.heightList[0]
1337 dataOut.heightList = numpy.arange(self.buffer.shape[1]) *self.step*deltaHeight + dataOut.heightList[0]
1328 dataOut.nProfiles = int(dataOut.nProfiles*self.nsamples)
1338 dataOut.nProfiles = int(dataOut.nProfiles*self.nsamples)
1329 dataOut.profileIndex = profileIndex
1339 dataOut.profileIndex = profileIndex
1330 dataOut.flagDataAsBlock = True
1340 dataOut.flagDataAsBlock = True
1331 dataOut.ippSeconds = ippSeconds
1341 dataOut.ippSeconds = ippSeconds
1332 dataOut.step = self.step
1342 dataOut.step = self.step
1333
1343
1334
1344
1335 import time
1345 import time
1336 #################################################
1346 #################################################
1337
1347
1338 class decoPseudorandom(Operation):
1348 class decoPseudorandom(Operation):
1339
1349
1340 nProfiles= 0
1350 nProfiles= 0
1341 buffer= None
1351 buffer= None
1342 isConfig = False
1352 isConfig = False
1343
1353
1344 def setup(self, clen= 10000,seed= 0,Nranges= 1000,oversample=1):
1354 def setup(self, clen= 10000,seed= 0,Nranges= 1000,oversample=1):
1345 #code = create_pseudo_random_code(clen=clen, seed=seed)
1355 #code = create_pseudo_random_code(clen=clen, seed=seed)
1346 code= rep_seq(create_pseudo_random_code(clen=clen, seed=seed),rep=oversample)
1356 code= rep_seq(create_pseudo_random_code(clen=clen, seed=seed),rep=oversample)
1347 #print ("code_rx", code.shape)
1357 #print ("code_rx", code.shape)
1348 #N = int(an_len/clen) # 100
1358 #N = int(an_len/clen) # 100
1349 B_cache = 0
1359 B_cache = 0
1350 r_cache = 0
1360 r_cache = 0
1351 B_cached = False
1361 B_cached = False
1352 r = create_estimation_matrix(code=code, cache=True, rmax=Nranges)
1362 r = create_estimation_matrix(code=code, cache=True, rmax=Nranges)
1353 #print ("code shape", code.shape)
1363 #print ("code shape", code.shape)
1354 #print ("seed",seed)
1364 #print ("seed",seed)
1355 #print ("Code", code[0:10])
1365 #print ("Code", code[0:10])
1356 self.B = r['B']
1366 self.B = r['B']
1357
1367
1358
1368
1359 def run (self,dataOut,length_code= 10000,seed= 0,Nranges= 1000,oversample=1):
1369 def run (self,dataOut,length_code= 10000,seed= 0,Nranges= 1000,oversample=1):
1360 #print((dataOut.data.shape))
1370 #print((dataOut.data.shape))
1361 if not self.isConfig:
1371 if not self.isConfig:
1362 self.setup(clen= length_code,seed= seed,Nranges= Nranges,oversample=oversample)
1372 self.setup(clen= length_code,seed= seed,Nranges= Nranges,oversample=oversample)
1363 self.isConfig = True
1373 self.isConfig = True
1364
1374
1365 dataOut.flagNoData = True
1375 dataOut.flagNoData = True
1366 data =dataOut.data
1376 data =dataOut.data
1367 #print "length_CODE",length_code
1377 #print "length_CODE",length_code
1368 data_shape = (data.shape[1])
1378 data_shape = (data.shape[1])
1369 #print "data_shape",data_shape
1379 #print "data_shape",data_shape
1370 n = (length_code /data_shape)
1380 n = (length_code /data_shape)
1371 #print "we need this number of sample",n
1381 #print "we need this number of sample",n
1372
1382
1373 if n>0 and self.buffer is None:
1383 if n>0 and self.buffer is None:
1374 self.buffer = numpy.zeros([1, length_code], dtype=numpy.complex64)
1384 self.buffer = numpy.zeros([1, length_code], dtype=numpy.complex64)
1375 self.buffer[0][0:data_shape] = data[0]
1385 self.buffer[0][0:data_shape] = data[0]
1376 #print "FIRST CREATION",self.buffer.shape
1386 #print "FIRST CREATION",self.buffer.shape
1377
1387
1378 else:
1388 else:
1379 self.buffer[0][self.nProfiles*data_shape:(self.nProfiles+1)*data_shape]=data[0]
1389 self.buffer[0][self.nProfiles*data_shape:(self.nProfiles+1)*data_shape]=data[0]
1380
1390
1381 #print "buffer_shape",(self.buffer.shape)
1391 #print "buffer_shape",(self.buffer.shape)
1382 self.nProfiles += 1
1392 self.nProfiles += 1
1383 #print "count",self.nProfiles
1393 #print "count",self.nProfiles
1384
1394
1385 if self.nProfiles== n:
1395 if self.nProfiles== n:
1386 temporal = numpy.dot(self.B, numpy.transpose(self.buffer))
1396 temporal = numpy.dot(self.B, numpy.transpose(self.buffer))
1387 #print temporal.shape
1397 #print temporal.shape
1388 #import time
1398 #import time
1389 #time.sleep(40)
1399 #time.sleep(40)
1390 dataOut.data=numpy.transpose(temporal)
1400 dataOut.data=numpy.transpose(temporal)
1391
1401
1392 dataOut.flagNoData = False
1402 dataOut.flagNoData = False
1393 self.buffer= None
1403 self.buffer= None
1394 self.nProfiles = 0
1404 self.nProfiles = 0
1395
1405
1396 # import collections
1406 # import collections
1397 # from scipy.stats import mode
1407 # from scipy.stats import mode
1398 #
1408 #
1399 # class Synchronize(Operation):
1409 # class Synchronize(Operation):
1400 #
1410 #
1401 # isConfig = False
1411 # isConfig = False
1402 # __profIndex = 0
1412 # __profIndex = 0
1403 #
1413 #
1404 # def __init__(self, **kwargs):
1414 # def __init__(self, **kwargs):
1405 #
1415 #
1406 # Operation.__init__(self, **kwargs)
1416 # Operation.__init__(self, **kwargs)
1407 # # self.isConfig = False
1417 # # self.isConfig = False
1408 # self.__powBuffer = None
1418 # self.__powBuffer = None
1409 # self.__startIndex = 0
1419 # self.__startIndex = 0
1410 # self.__pulseFound = False
1420 # self.__pulseFound = False
1411 #
1421 #
1412 # def __findTxPulse(self, dataOut, channel=0, pulse_with = None):
1422 # def __findTxPulse(self, dataOut, channel=0, pulse_with = None):
1413 #
1423 #
1414 # #Read data
1424 # #Read data
1415 #
1425 #
1416 # powerdB = dataOut.getPower(channel = channel)
1426 # powerdB = dataOut.getPower(channel = channel)
1417 # noisedB = dataOut.getNoise(channel = channel)[0]
1427 # noisedB = dataOut.getNoise(channel = channel)[0]
1418 #
1428 #
1419 # self.__powBuffer.extend(powerdB.flatten())
1429 # self.__powBuffer.extend(powerdB.flatten())
1420 #
1430 #
1421 # dataArray = numpy.array(self.__powBuffer)
1431 # dataArray = numpy.array(self.__powBuffer)
1422 #
1432 #
1423 # filteredPower = numpy.correlate(dataArray, dataArray[0:self.__nSamples], "same")
1433 # filteredPower = numpy.correlate(dataArray, dataArray[0:self.__nSamples], "same")
1424 #
1434 #
1425 # maxValue = numpy.nanmax(filteredPower)
1435 # maxValue = numpy.nanmax(filteredPower)
1426 #
1436 #
1427 # if maxValue < noisedB + 10:
1437 # if maxValue < noisedB + 10:
1428 # #No se encuentra ningun pulso de transmision
1438 # #No se encuentra ningun pulso de transmision
1429 # return None
1439 # return None
1430 #
1440 #
1431 # maxValuesIndex = numpy.where(filteredPower > maxValue - 0.1*abs(maxValue))[0]
1441 # maxValuesIndex = numpy.where(filteredPower > maxValue - 0.1*abs(maxValue))[0]
1432 #
1442 #
1433 # if len(maxValuesIndex) < 2:
1443 # if len(maxValuesIndex) < 2:
1434 # #Solo se encontro un solo pulso de transmision de un baudio, esperando por el siguiente TX
1444 # #Solo se encontro un solo pulso de transmision de un baudio, esperando por el siguiente TX
1435 # return None
1445 # return None
1436 #
1446 #
1437 # phasedMaxValuesIndex = maxValuesIndex - self.__nSamples
1447 # phasedMaxValuesIndex = maxValuesIndex - self.__nSamples
1438 #
1448 #
1439 # #Seleccionar solo valores con un espaciamiento de nSamples
1449 # #Seleccionar solo valores con un espaciamiento de nSamples
1440 # pulseIndex = numpy.intersect1d(maxValuesIndex, phasedMaxValuesIndex)
1450 # pulseIndex = numpy.intersect1d(maxValuesIndex, phasedMaxValuesIndex)
1441 #
1451 #
1442 # if len(pulseIndex) < 2:
1452 # if len(pulseIndex) < 2:
1443 # #Solo se encontro un pulso de transmision con ancho mayor a 1
1453 # #Solo se encontro un pulso de transmision con ancho mayor a 1
1444 # return None
1454 # return None
1445 #
1455 #
1446 # spacing = pulseIndex[1:] - pulseIndex[:-1]
1456 # spacing = pulseIndex[1:] - pulseIndex[:-1]
1447 #
1457 #
1448 # #remover senales que se distancien menos de 10 unidades o muestras
1458 # #remover senales que se distancien menos de 10 unidades o muestras
1449 # #(No deberian existir IPP menor a 10 unidades)
1459 # #(No deberian existir IPP menor a 10 unidades)
1450 #
1460 #
1451 # realIndex = numpy.where(spacing > 10 )[0]
1461 # realIndex = numpy.where(spacing > 10 )[0]
1452 #
1462 #
1453 # if len(realIndex) < 2:
1463 # if len(realIndex) < 2:
1454 # #Solo se encontro un pulso de transmision con ancho mayor a 1
1464 # #Solo se encontro un pulso de transmision con ancho mayor a 1
1455 # return None
1465 # return None
1456 #
1466 #
1457 # #Eliminar pulsos anchos (deja solo la diferencia entre IPPs)
1467 # #Eliminar pulsos anchos (deja solo la diferencia entre IPPs)
1458 # realPulseIndex = pulseIndex[realIndex]
1468 # realPulseIndex = pulseIndex[realIndex]
1459 #
1469 #
1460 # period = mode(realPulseIndex[1:] - realPulseIndex[:-1])[0][0]
1470 # period = mode(realPulseIndex[1:] - realPulseIndex[:-1])[0][0]
1461 #
1471 #
1462 # print "IPP = %d samples" %period
1472 # print "IPP = %d samples" %period
1463 #
1473 #
1464 # self.__newNSamples = dataOut.nHeights #int(period)
1474 # self.__newNSamples = dataOut.nHeights #int(period)
1465 # self.__startIndex = int(realPulseIndex[0])
1475 # self.__startIndex = int(realPulseIndex[0])
1466 #
1476 #
1467 # return 1
1477 # return 1
1468 #
1478 #
1469 #
1479 #
1470 # def setup(self, nSamples, nChannels, buffer_size = 4):
1480 # def setup(self, nSamples, nChannels, buffer_size = 4):
1471 #
1481 #
1472 # self.__powBuffer = collections.deque(numpy.zeros( buffer_size*nSamples,dtype=numpy.float),
1482 # self.__powBuffer = collections.deque(numpy.zeros( buffer_size*nSamples,dtype=numpy.float),
1473 # maxlen = buffer_size*nSamples)
1483 # maxlen = buffer_size*nSamples)
1474 #
1484 #
1475 # bufferList = []
1485 # bufferList = []
1476 #
1486 #
1477 # for i in range(nChannels):
1487 # for i in range(nChannels):
1478 # bufferByChannel = collections.deque(numpy.zeros( buffer_size*nSamples, dtype=numpy.complex) + numpy.NAN,
1488 # bufferByChannel = collections.deque(numpy.zeros( buffer_size*nSamples, dtype=numpy.complex) + numpy.NAN,
1479 # maxlen = buffer_size*nSamples)
1489 # maxlen = buffer_size*nSamples)
1480 #
1490 #
1481 # bufferList.append(bufferByChannel)
1491 # bufferList.append(bufferByChannel)
1482 #
1492 #
1483 # self.__nSamples = nSamples
1493 # self.__nSamples = nSamples
1484 # self.__nChannels = nChannels
1494 # self.__nChannels = nChannels
1485 # self.__bufferList = bufferList
1495 # self.__bufferList = bufferList
1486 #
1496 #
1487 # def run(self, dataOut, channel = 0):
1497 # def run(self, dataOut, channel = 0):
1488 #
1498 #
1489 # if not self.isConfig:
1499 # if not self.isConfig:
1490 # nSamples = dataOut.nHeights
1500 # nSamples = dataOut.nHeights
1491 # nChannels = dataOut.nChannels
1501 # nChannels = dataOut.nChannels
1492 # self.setup(nSamples, nChannels)
1502 # self.setup(nSamples, nChannels)
1493 # self.isConfig = True
1503 # self.isConfig = True
1494 #
1504 #
1495 # #Append new data to internal buffer
1505 # #Append new data to internal buffer
1496 # for thisChannel in range(self.__nChannels):
1506 # for thisChannel in range(self.__nChannels):
1497 # bufferByChannel = self.__bufferList[thisChannel]
1507 # bufferByChannel = self.__bufferList[thisChannel]
1498 # bufferByChannel.extend(dataOut.data[thisChannel])
1508 # bufferByChannel.extend(dataOut.data[thisChannel])
1499 #
1509 #
1500 # if self.__pulseFound:
1510 # if self.__pulseFound:
1501 # self.__startIndex -= self.__nSamples
1511 # self.__startIndex -= self.__nSamples
1502 #
1512 #
1503 # #Finding Tx Pulse
1513 # #Finding Tx Pulse
1504 # if not self.__pulseFound:
1514 # if not self.__pulseFound:
1505 # indexFound = self.__findTxPulse(dataOut, channel)
1515 # indexFound = self.__findTxPulse(dataOut, channel)
1506 #
1516 #
1507 # if indexFound == None:
1517 # if indexFound == None:
1508 # dataOut.flagNoData = True
1518 # dataOut.flagNoData = True
1509 # return
1519 # return
1510 #
1520 #
1511 # self.__arrayBuffer = numpy.zeros((self.__nChannels, self.__newNSamples), dtype = numpy.complex)
1521 # self.__arrayBuffer = numpy.zeros((self.__nChannels, self.__newNSamples), dtype = numpy.complex)
1512 # self.__pulseFound = True
1522 # self.__pulseFound = True
1513 # self.__startIndex = indexFound
1523 # self.__startIndex = indexFound
1514 #
1524 #
1515 # #If pulse was found ...
1525 # #If pulse was found ...
1516 # for thisChannel in range(self.__nChannels):
1526 # for thisChannel in range(self.__nChannels):
1517 # bufferByChannel = self.__bufferList[thisChannel]
1527 # bufferByChannel = self.__bufferList[thisChannel]
1518 # #print self.__startIndex
1528 # #print self.__startIndex
1519 # x = numpy.array(bufferByChannel)
1529 # x = numpy.array(bufferByChannel)
1520 # self.__arrayBuffer[thisChannel] = x[self.__startIndex:self.__startIndex+self.__newNSamples]
1530 # self.__arrayBuffer[thisChannel] = x[self.__startIndex:self.__startIndex+self.__newNSamples]
1521 #
1531 #
1522 # deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1532 # deltaHeight = dataOut.heightList[1] - dataOut.heightList[0]
1523 # dataOut.heightList = numpy.arange(self.__newNSamples)*deltaHeight
1533 # dataOut.heightList = numpy.arange(self.__newNSamples)*deltaHeight
1524 # # dataOut.ippSeconds = (self.__newNSamples / deltaHeight)/1e6
1534 # # dataOut.ippSeconds = (self.__newNSamples / deltaHeight)/1e6
1525 #
1535 #
1526 # dataOut.data = self.__arrayBuffer
1536 # dataOut.data = self.__arrayBuffer
1527 #
1537 #
1528 # self.__startIndex += self.__newNSamples
1538 # self.__startIndex += self.__newNSamples
1529 #
1539 #
1530 # return
1540 # return
General Comments 0
You need to be logged in to leave comments. Login now