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