##// END OF EJS Templates
schain
RSS
update for weather radar options
avaldez -
r1296:162740cfc930
parent child
Show More
Side by Side Unified
Context file:
r1296:162740cfc930 Loading diff...:
Show file before | Show file after | Hide comments
@@ -0,0 +1,473
1 import numpy,math,random,time
2 #---------------1 Heredamos JRODatareader
3 from schainpy.model.io.jroIO_base import *
4 #---------------2 Heredamos las propiedades de ProcessingUnit
5 from schainpy.model.proc.jroproc_base import ProcessingUnit,Operation,MPDecorator
6 #---------------3 Importaremos las clases BascicHeader, SystemHeader, RadarControlHeader, ProcessingHeader
7 from schainpy.model.data.jroheaderIO import PROCFLAG, BasicHeader,SystemHeader,RadarControllerHeader, ProcessingHeader
8 #---------------4 Importaremos el objeto Voltge
9 from schainpy.model.data.jrodata import Voltage
10
11 class SimulatorReader(JRODataReader, ProcessingUnit):
12 incIntFactor = 1
13 nFFTPoints = 0
14 FixPP_IncInt = 1
15 FixRCP_IPP = 1000
16 FixPP_CohInt = 1
17 Tau_0 = 250
18 AcqH0_0 = 70
19 H0 = AcqH0_0
20 AcqDH_0 = 1.25
21 DH0 = AcqDH_0
22 Bauds = 32
23 BaudWidth = None
24 FixRCP_TXA = 40
25 FixRCP_TXB = 70
26 fAngle = 2.0*math.pi*(1/16)
27 DC_level = 500
28 stdev = 8
29 Num_Codes = 2
30 #code0 = numpy.array([1,1,1,0,1,1,0,1,1,1,1,0,0,0,1,0,1,1,1,0,1,1,0,1,0,0,0,1,1,1,0,1])
31 #code1 = numpy.array([1,1,1,0,1,1,0,1,1,1,1,0,0,0,1,0,0,0,0,1,0,0,1,0,1,1,1,0,0,0,1,0])
32 #Dyn_snCode = numpy.array([Num_Codes,Bauds])
33 Dyn_snCode = None
34 Samples = 200
35 channels = 5
36 pulses = None
37 Reference = None
38 pulse_size = None
39 prof_gen = None
40 Fdoppler = 100
41 Hdoppler = 36
42 Adoppler = 300
43 frequency = 9345
44 nTotalReadFiles = 1000
45
46 def __init__(self):
47 """
48 Inicializador de la clases SimulatorReader para
49 generar datos de voltage simulados.
50 Input:
51 dataOut: Objeto de la clase Voltage.
52 Este Objeto sera utilizado apra almacenar
53 un perfil de datos cada vez qe se haga
54 un requerimiento (getData)
55 """
56 ProcessingUnit.__init__(self)
57 print(" [ START ] init - Metodo Simulator Reader")
58
59 self.isConfig = False
60 self.basicHeaderObj = BasicHeader(LOCALTIME)
61 self.systemHeaderObj = SystemHeader()
62 self.radarControllerHeaderObj = RadarControllerHeader()
63 self.processingHeaderObj = ProcessingHeader()
64 self.profileIndex = 2**32-1
65 self.dataOut = Voltage()
66 #code0 = numpy.array([1,1,1,0,1,1,0,1,1,1,1,0,0,0,1,0,1,1,1,0,1,1,0,1,0,0,0,1,1,1,0,1])
67 code0 = numpy.array([1,1,1,-1,1,1,-1,1,1,1,1,-1,-1,-1,1,-1,1,1,1,-1,1,1,-1,1,-1,-1,-1,1,1,1,-1,1])
68 #code1 = numpy.array([1,1,1,0,1,1,0,1,1,1,1,0,0,0,1,0,0,0,0,1,0,0,1,0,1,1,1,0,0,0,1,0])
69 code1 = numpy.array([1,1,1,-1,1,1,-1,1,1,1,1,-1,-1,-1,1,-1,-1,-1,-1,1,-1,-1,1,-1,1,1,1,-1,-1,-1,1,-1])
70 #self.Dyn_snCode = numpy.array([code0,code1])
71 self.Dyn_snCode = None
72
73 def set_kwargs(self, **kwargs):
74 for key, value in kwargs.items():
75 setattr(self, key, value)
76
77 def __hasNotDataInBuffer(self):
78
79 if self.profileIndex >= self.processingHeaderObj.profilesPerBlock* self.nTxs:
80 if self.nReadBlocks>0:
81 tmp = self.dataOut.utctime
82 tmp_utc = int(self.dataOut.utctime)
83 tmp_milisecond = int((tmp-tmp_utc)*1000)
84 self.basicHeaderObj.utc = tmp_utc
85 self.basicHeaderObj.miliSecond= tmp_milisecond
86 return 1
87 return 0
88
89 def setNextFile(self):
90 """Set the next file to be readed open it and parse de file header"""
91
92 if (self.nReadBlocks >= self.processingHeaderObj.dataBlocksPerFile):
93 self.nReadFiles=self.nReadFiles+1
94 if self.nReadFiles ==self.nTotalReadFiles:
95 self.flagNoMoreFiles=1
96 raise schainpy.admin.SchainWarning('No more files to read')
97
98 print('------------------- [Opening file] ------------------------------',self.nReadFiles)
99 self.nReadBlocks = 0
100
101 def __setNewBlock(self):
102 self.setNextFile()
103 if self.flagIsNewFile:
104 return 1
105
106 def readNextBlock(self):
107 while True:
108 self.__setNewBlock()
109 if not(self.readBlock()):
110 return 0
111 self.getBasicHeader()
112 break
113 if self.verbose:
114 print("[Reading] Block No. %d/%d -> %s" %(self.nReadBlocks,
115 self.processingHeaderObj.dataBlocksPerFile,
116 self.dataOut.datatime.ctime()) )
117 return 1
118
119 def getFirstHeader(self):
120 self.getBasicHeader()
121 self.dataOut.processingHeaderObj = self.processingHeaderObj.copy()
122 self.dataOut.systemHeaderObj = self.systemHeaderObj.copy()
123 self.dataOut.radarControllerHeaderObj = self.radarControllerHeaderObj.copy()
124
125 self.dataOut.nProfiles = self.processingHeaderObj.profilesPerBlock
126 self.dataOut.heightList = numpy.arange(self.processingHeaderObj.nHeights) * self.processingHeaderObj.deltaHeight + self.processingHeaderObj.firstHeight
127 self.dataOut.channelList = list(range(self.systemHeaderObj.nChannels))
128 self.dataOut.nCohInt = self.processingHeaderObj.nCohInt
129 # asumo q la data no esta decodificada
130 self.dataOut.flagDecodeData = self.processingHeaderObj.flag_decode
131 # asumo q la data no esta sin flip
132 self.dataOut.flagDeflipData = self.processingHeaderObj.flag_deflip
133 self.dataOut.flagShiftFFT = self.processingHeaderObj.shif_fft
134 self.dataOut.frequency = self.frequency
135
136 def getBasicHeader(self):
137 self.dataOut.utctime = self.basicHeaderObj.utc + self.basicHeaderObj.miliSecond / \
138 1000. + self.profileIndex * self.radarControllerHeaderObj.ippSeconds
139
140 self.dataOut.flagDiscontinuousBlock = self.flagDiscontinuousBlock
141 self.dataOut.timeZone = self.basicHeaderObj.timeZone
142 self.dataOut.dstFlag = self.basicHeaderObj.dstFlag
143 self.dataOut.errorCount = self.basicHeaderObj.errorCount
144 self.dataOut.useLocalTime = self.basicHeaderObj.useLocalTime
145 self.dataOut.ippSeconds = self.radarControllerHeaderObj.ippSeconds / self.nTxs
146
147 def set_RCH(self, expType=2, nTx=1,ipp=None, txA=0, txB=0,
148 nWindows=None, nHeights=None, firstHeight=None, deltaHeight=None,
149 numTaus=0, line6Function=0, line5Function=0, fClock=None,
150 prePulseBefore=0, prePulseAfter=0,
151 codeType=0, nCode=0, nBaud=0, code=None,
152 flip1=0, flip2=0):
153 self.radarControllerHeaderObj.expType = expType
154 self.radarControllerHeaderObj.nTx = nTx
155 self.radarControllerHeaderObj.ipp = float(ipp)
156 self.radarControllerHeaderObj.txA = float(txA)
157 self.radarControllerHeaderObj.txB = float(txB)
158 self.radarControllerHeaderObj.rangeIPP = ipp
159 self.radarControllerHeaderObj.rangeTxA = txA
160 self.radarControllerHeaderObj.rangeTxB = txB
161
162 self.radarControllerHeaderObj.nHeights = int(nHeights)
163 self.radarControllerHeaderObj.firstHeight = numpy.array([firstHeight])
164 self.radarControllerHeaderObj.deltaHeight = numpy.array([deltaHeight])
165 self.radarControllerHeaderObj.samplesWin = numpy.array([nHeights])
166
167
168 self.radarControllerHeaderObj.nWindows = nWindows
169 self.radarControllerHeaderObj.numTaus = numTaus
170 self.radarControllerHeaderObj.codeType = codeType
171 self.radarControllerHeaderObj.line6Function = line6Function
172 self.radarControllerHeaderObj.line5Function = line5Function
173 self.radarControllerHeaderObj.fclock = fClock
174 self.radarControllerHeaderObj.prePulseBefore= prePulseBefore
175 self.radarControllerHeaderObj.prePulseAfter = prePulseAfter
176
177 self.radarControllerHeaderObj.nCode = nCode
178 self.radarControllerHeaderObj.nBaud = nBaud
179 self.radarControllerHeaderObj.code = code
180 self.radarControllerHeaderObj.flip1 = flip1
181 self.radarControllerHeaderObj.flip2 = flip2
182
183 self.radarControllerHeaderObj.code_size = int(numpy.ceil(nBaud / 32.)) * nCode * 4
184
185 if fClock is None and deltaHeight is not None:
186 self.fClock = 0.15 / (deltaHeight * 1e-6)
187
188 def set_PH(self, dtype=0, blockSize=0, profilesPerBlock=0,
189 dataBlocksPerFile=0, nWindows=0, processFlags=0, nCohInt=0,
190 nIncohInt=0, totalSpectra=0, nHeights=0, firstHeight=0,
191 deltaHeight=0, samplesWin=0, spectraComb=0, nCode=0,
192 code=0, nBaud=None, shif_fft=False, flag_dc=False,
193 flag_cspc=False, flag_decode=False, flag_deflip=False):
194
195 self.processingHeaderObj.profilesPerBlock = profilesPerBlock
196 self.processingHeaderObj.dataBlocksPerFile = dataBlocksPerFile
197 self.processingHeaderObj.nWindows = nWindows
198 self.processingHeaderObj.nCohInt = nCohInt
199 self.processingHeaderObj.nIncohInt = nIncohInt
200 self.processingHeaderObj.totalSpectra = totalSpectra
201 self.processingHeaderObj.nHeights = int(nHeights)
202 self.processingHeaderObj.firstHeight = firstHeight
203 self.processingHeaderObj.deltaHeight = deltaHeight
204 self.processingHeaderObj.samplesWin = nHeights
205
206 def set_BH(self, utc = 0, miliSecond = 0, timeZone = 0):
207 self.basicHeaderObj.utc = utc
208 self.basicHeaderObj.miliSecond = miliSecond
209 self.basicHeaderObj.timeZone = timeZone
210
211 def set_SH(self, nSamples=0, nProfiles=0, nChannels=0, adcResolution=14, pciDioBusWidth=0):
212 self.systemHeaderObj.nSamples = nSamples
213 self.systemHeaderObj.nProfiles = nProfiles
214 self.systemHeaderObj.nChannels = nChannels
215 self.systemHeaderObj.adcResolution = adcResolution
216 self.systemHeaderObj.pciDioBusWidth = pciDioBusWidth
217
218 def init_acquisition(self):
219
220 if self.nFFTPoints != 0:
221 self.incIntFactor = m_nProfilesperBlock/self.nFFTPoints
222 if (self.FixPP_IncInt > self.incIntFactor):
223 self.incIntFactor = self.FixPP_IncInt/ self.incIntFactor
224 elif(self.FixPP_IncInt< self.incIntFactor):
225 print("False alert...")
226
227 ProfilesperBlock = self.processingHeaderObj.profilesPerBlock
228
229 self.timeperblock =int(((self.FixRCP_IPP
230 *ProfilesperBlock
231 *self.FixPP_CohInt
232 *self.incIntFactor)
233 /150.0)
234 *0.9
235 +0.5)
236 # para cada canal
237 self.profiles = ProfilesperBlock*self.FixPP_CohInt
238 self.profiles = ProfilesperBlock
239 self.Reference = int((self.Tau_0-self.AcqH0_0)/(self.AcqDH_0)+0.5)
240 self.BaudWidth = int((self.FixRCP_TXA/self.AcqDH_0)/self.Bauds + 0.5 )
241
242 if (self.BaudWidth==0):
243 self.BaudWidth=1
244
245 def init_pulse(self,Num_Codes=Num_Codes,Bauds=Bauds,BaudWidth=BaudWidth,Dyn_snCode=Dyn_snCode):
246
247 Num_Codes = Num_Codes
248 Bauds = Bauds
249 BaudWidth = BaudWidth
250 Dyn_snCode = Dyn_snCode
251
252 if Dyn_snCode:
253 print("EXISTE")
254 else:
255 print("No existe")
256
257 if Dyn_snCode: # if Bauds:
258 pulses = list(range(0,Num_Codes))
259 num_codes = Num_Codes
260 for i in range(num_codes):
261 pulse_size = Bauds*BaudWidth
262 pulses[i] = numpy.zeros(pulse_size)
263 for j in range(Bauds):
264 for k in range(BaudWidth):
265 pulses[i][j*BaudWidth+k] = int(Dyn_snCode[i][j]*600)
266 else:
267 print("sin code")
268 pulses = list(range(1))
269 if self.AcqDH_0>0.149:
270 pulse_size = int(self.FixRCP_TXB/0.15+0.5)
271 else:
272 pulse_size = int((self.FixRCP_TXB/self.AcqDH_0)+0.5) #0.0375
273 pulses[0] = numpy.ones(pulse_size)
274 pulses = 600*pulses[0]
275
276 return pulses,pulse_size
277
278 def jro_GenerateBlockOfData(self,Samples=Samples,DC_level= DC_level,stdev=stdev,
279 Reference= Reference,pulses= pulses,
280 Num_Codes= Num_Codes,pulse_size=pulse_size,
281 prof_gen= prof_gen,H0 = H0,DH0=DH0,
282 Adoppler=Adoppler,Fdoppler= Fdoppler,Hdoppler=Hdoppler):
283 Samples = Samples
284 DC_level = DC_level
285 stdev = stdev
286 m_nR = Reference
287 pulses = pulses
288 num_codes = Num_Codes
289 ps = pulse_size
290 prof_gen = prof_gen
291 channels = self.channels
292 H0 = H0
293 DH0 = DH0
294 ippSec = self.radarControllerHeaderObj.ippSeconds
295 Fdoppler = self.Fdoppler
296 Hdoppler = self.Hdoppler
297 Adoppler = self.Adoppler
298
299 self.datablock = numpy.zeros([channels,prof_gen,Samples],dtype= numpy.complex64)
300 for i in range(channels):
301 for k in range(prof_gen):
302 #Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·NOISEΒ·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·
303 Noise_r = numpy.random.normal(DC_level,stdev,Samples)
304 Noise_i = numpy.random.normal(DC_level,stdev,Samples)
305 Noise = numpy.zeros(Samples,dtype=complex)
306 Noise.real = Noise_r
307 Noise.imag = Noise_i
308 #Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·PULSOSΒ·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·
309 Pulso = numpy.zeros(pulse_size,dtype=complex)
310 Pulso.real = pulses[k%num_codes]
311 Pulso.imag = pulses[k%num_codes]
312 #Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β· PULSES+NOISEΒ·Β·Β·Β·Β·Β·Β·Β·Β·Β·
313 InBuffer = numpy.zeros(Samples,dtype=complex)
314 InBuffer[m_nR:m_nR+ps] = Pulso
315 InBuffer = InBuffer+Noise
316 #Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β· ANGLE Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·
317 InBuffer.real[m_nR:m_nR+ps] = InBuffer.real[m_nR:m_nR+ps]*(math.cos( self.fAngle)*5)
318 InBuffer.imag[m_nR:m_nR+ps] = InBuffer.imag[m_nR:m_nR+ps]*(math.sin( self.fAngle)*5)
319 InBuffer=InBuffer
320 self.datablock[i][k]= InBuffer
321 #plot_cts(InBuffer,H0=H0,DH0=DH0)
322 #wave_fft(x=InBuffer,plot_show=True)
323 #time.sleep(1)
324 #Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·DOPPLER SIGNAL...............................................
325 time_vec = numpy.linspace(0,(prof_gen-1)*ippSec,int(prof_gen))+self.nReadBlocks*ippSec*prof_gen+(self.nReadFiles-1)*ippSec*prof_gen
326 fd = Fdoppler #+(600.0/120)*self.nReadBlocks
327 d_signal = Adoppler*numpy.array(numpy.exp(1.0j*2.0*math.pi*fd*time_vec),dtype=numpy.complex64)
328 #Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β· DATABLOCK + DOPPLERΒ·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·Β·...........................
329 HD=int(Hdoppler/self.AcqDH_0)
330 for i in range(12):
331 self.datablock[:,:,HD+i]=self.datablock[:,:,HD+i]+ d_signal # RESULT
332 '''
333 a= numpy.zeros(10)
334 for i in range(10):
335 a[i]=i+self.nReadBlocks+20
336 for i in a:
337 self.datablock[0,:,int(i)]=self.datablock[0,:,int(i)]+ d_signal # RESULT
338 '''
339
340 def readBlock(self):
341
342 self.jro_GenerateBlockOfData(Samples= self.samples,DC_level=self.DC_level,
343 stdev=self.stdev,Reference= self.Reference,
344 pulses = self.pulses,Num_Codes=self.Num_Codes,
345 pulse_size=self.pulse_size,prof_gen=self.profiles,
346 H0=self.H0,DH0=self.DH0)
347
348 self.profileIndex = 0
349 self.flagIsNewFile = 0
350 self.flagIsNewBlock = 1
351 self.nTotalBlocks += 1
352 self.nReadBlocks += 1
353
354 return 1
355
356
357 def getData(self):
358 if self.flagNoMoreFiles:
359 self.dataOut.flagNodata = True
360 return 0
361 self.flagDiscontinuousBlock = 0
362 self.flagIsNewBlock = 0
363 if self.__hasNotDataInBuffer(): # aqui es verdad
364 if not(self.readNextBlock()): # return 1 y por eso el if not salta a getBasic Header
365 return 0
366 self.getFirstHeader() # atributo
367
368 if not self.getByBlock:
369 self.dataOut.flagDataAsBlock = False
370 self.dataOut.data = self.datablock[:, self.profileIndex, :]
371 self.dataOut.profileIndex = self.profileIndex
372 self.profileIndex += 1
373 else:
374 pass
375 self.dataOut.flagNoData = False
376 self.getBasicHeader()
377 self.dataOut.realtime = self.online
378 return self.dataOut.data
379
380
381 def setup(self,frequency=49.92e6,incIntFactor= 1, nFFTPoints = 0, FixPP_IncInt=1,FixRCP_IPP=1000,
382 FixPP_CohInt= 1,Tau_0= 250,AcqH0_0 = 70 ,AcqDH_0=1.25, Bauds= 32,
383 FixRCP_TXA = 40, FixRCP_TXB = 50, fAngle = 2.0*math.pi*(1/16),DC_level= 50,
384 stdev= 8,Num_Codes = 1 , Dyn_snCode = None, samples=200,
385 channels=2,Fdoppler=20,Hdoppler=36,Adoppler=500,nTotalReadFiles=10000,
386 **kwargs):
387
388 self.set_kwargs(**kwargs)
389 self.nReadBlocks = 0
390 self.nReadFiles = 1
391 print('------------------- [Opening file: ] ------------------------------',self.nReadFiles)
392
393 tmp = time.time()
394 tmp_utc = int(tmp)
395 tmp_milisecond = int((tmp-tmp_utc)*1000)
396 print(" SETUP -basicHeaderObj.utc",datetime.datetime.utcfromtimestamp(tmp))
397 if Dyn_snCode is None:
398 Num_Codes=1
399 Bauds =1
400
401
402
403 self.set_BH(utc= tmp_utc,miliSecond= tmp_milisecond,timeZone=300 )
404
405 self.set_RCH( expType=0, nTx=150,ipp=FixRCP_IPP, txA=FixRCP_TXA, txB= FixRCP_TXB,
406 nWindows=1 , nHeights=samples, firstHeight=AcqH0_0, deltaHeight=AcqDH_0,
407 numTaus=1, line6Function=0, line5Function=0, fClock=None,
408 prePulseBefore=0, prePulseAfter=0,
409 codeType=14, nCode=Num_Codes, nBaud=32, code=Dyn_snCode,
410 flip1=0, flip2=0)
411
412 self.set_PH(dtype=0, blockSize=0, profilesPerBlock=300,
413 dataBlocksPerFile=120, nWindows=1, processFlags=0, nCohInt=1,
414 nIncohInt=1, totalSpectra=0, nHeights=samples, firstHeight=AcqH0_0,
415 deltaHeight=AcqDH_0, samplesWin=samples, spectraComb=0, nCode=0,
416 code=0, nBaud=None, shif_fft=False, flag_dc=False,
417 flag_cspc=False, flag_decode=False, flag_deflip=False)
418
419 self.set_SH(nSamples=samples, nProfiles=300, nChannels=channels)
420
421
422 self.frequency = frequency
423 self.incIntFactor = incIntFactor
424 self.nFFTPoints = nFFTPoints
425 self.FixPP_IncInt = FixPP_IncInt
426 self.FixRCP_IPP = FixRCP_IPP
427 self.FixPP_CohInt = FixPP_CohInt
428 self.Tau_0 = Tau_0
429 self.AcqH0_0 = AcqH0_0
430 self.H0 = AcqH0_0
431 self.AcqDH_0 = AcqDH_0
432 self.DH0 = AcqDH_0
433 self.Bauds = Bauds
434 self.FixRCP_TXA = FixRCP_TXA
435 self.FixRCP_TXB = FixRCP_TXB
436 self.fAngle = fAngle
437 self.DC_level = DC_level
438 self.stdev = stdev
439 self.Num_Codes = Num_Codes
440 self.Dyn_snCode = Dyn_snCode
441 self.samples = samples
442 self.channels = channels
443 self.profiles = None
444 self.m_nReference = None
445 self.Baudwidth = None
446 self.Fdoppler = Fdoppler
447 self.Hdoppler = Hdoppler
448 self.Adoppler = Adoppler
449 self.nTotalReadFiles = int(nTotalReadFiles)
450
451 print("IPP ", self.FixRCP_IPP)
452 print("Tau_0 ",self.Tau_0)
453 print("AcqH0_0",self.AcqH0_0)
454 print("samples,window ",self.samples)
455 print("AcqDH_0",AcqDH_0)
456 print("FixRCP_TXA",self.FixRCP_TXA)
457 print("FixRCP_TXB",self.FixRCP_TXB)
458 print("Dyn_snCode",Dyn_snCode)
459 print("Fdoppler", Fdoppler)
460 print("Hdoppler",Hdoppler)
461 print("Vdopplermax",Fdoppler*(3.0e8/self.frequency)/2.0)
462 print("nTotalReadFiles", nTotalReadFiles)
463
464 self.init_acquisition()
465 self.pulses,self.pulse_size=self.init_pulse(Num_Codes=self.Num_Codes,Bauds=self.Bauds,BaudWidth=self.BaudWidth,Dyn_snCode=Dyn_snCode)
466 print(" [ END ] - SETUP metodo")
467 return
468
469 def run(self,**kwargs): # metodo propio
470 if not(self.isConfig):
471 self.setup(**kwargs)
472 self.isConfig = True
473 self.getData()
Show file before | Show file after | Hide comments
@@ -360,7 +360,8 class Voltage(JROData):
360
360
361 # data es un numpy array de 2 dmensiones (canales, alturas)
361 # data es un numpy array de 2 dmensiones (canales, alturas)
362 data = None
362 data = None
363
363 data_intensity = None
364 data_velocity = None
364 def __init__(self):
365 def __init__(self):
365 '''
366 '''
366 Constructor
367 Constructor
@@ -1139,6 +1140,7 class PlotterData(object):
1139 for plot in self.plottypes:
1140 for plot in self.plottypes:
1140 self.data[plot] = {}
1141 self.data[plot] = {}
1141
1142
1143
1142 def __str__(self):
1144 def __str__(self):
1143 dum = ['{}{}'.format(key, self.shape(key)) for key in self.data]
1145 dum = ['{}{}'.format(key, self.shape(key)) for key in self.data]
1144 return 'Data[{}][{}]'.format(';'.join(dum), len(self.__times))
1146 return 'Data[{}][{}]'.format(';'.join(dum), len(self.__times))
@@ -1167,7 +1169,6 class PlotterData(object):
1167 '''
1169 '''
1168 Configure object
1170 Configure object
1169 '''
1171 '''
1170
1171 self.type = ''
1172 self.type = ''
1172 self.ready = False
1173 self.ready = False
1173 self.data = {}
1174 self.data = {}
@@ -1204,7 +1205,6 class PlotterData(object):
1204 '''
1205 '''
1205 Update data object with new dataOut
1206 Update data object with new dataOut
1206 '''
1207 '''
1207
1208 if tm in self.__times:
1208 if tm in self.__times:
1209 return
1209 return
1210 self.profileIndex = dataOut.profileIndex
1210 self.profileIndex = dataOut.profileIndex
@@ -1227,7 +1227,6 class PlotterData(object):
1227 self.__heights.append(dataOut.heightList)
1227 self.__heights.append(dataOut.heightList)
1228 self.__all_heights.update(dataOut.heightList)
1228 self.__all_heights.update(dataOut.heightList)
1229 self.__times.append(tm)
1229 self.__times.append(tm)
1230
1231 for plot in self.plottypes:
1230 for plot in self.plottypes:
1232 if plot in ('spc', 'spc_moments', 'spc_cut'):
1231 if plot in ('spc', 'spc_moments', 'spc_cut'):
1233 z = dataOut.data_spc/dataOut.normFactor
1232 z = dataOut.data_spc/dataOut.normFactor
@@ -1261,6 +1260,14 class PlotterData(object):
1261 buffer = dataOut.data
1260 buffer = dataOut.data
1262 self.flagDataAsBlock = dataOut.flagDataAsBlock
1261 self.flagDataAsBlock = dataOut.flagDataAsBlock
1263 self.nProfiles = dataOut.nProfiles
1262 self.nProfiles = dataOut.nProfiles
1263 if plot == 'pp_power':
1264 buffer = dataOut.data_intensity
1265 self.flagDataAsBlock = dataOut.flagDataAsBlock
1266 self.nProfiles = dataOut.nProfiles
1267 if plot == 'pp_velocity':
1268 buffer = dataOut.data_velocity
1269 self.flagDataAsBlock = dataOut.flagDataAsBlock
1270 self.nProfiles = dataOut.nProfiles
1264
1271
1265 if plot == 'spc':
1272 if plot == 'spc':
1266 self.data['spc'] = buffer
1273 self.data['spc'] = buffer
Show file before | Show file after | Hide comments
@@ -60,6 +60,7 class ScopePlot(Plot):
60 def plot_power(self, x, y, channelIndexList, thisDatetime, wintitle):
60 def plot_power(self, x, y, channelIndexList, thisDatetime, wintitle):
61 y = y[channelIndexList,:] * numpy.conjugate(y[channelIndexList,:])
61 y = y[channelIndexList,:] * numpy.conjugate(y[channelIndexList,:])
62 yreal = y.real
62 yreal = y.real
63 yreal = 10*numpy.log10(yreal)
63 self.y = yreal
64 self.y = yreal
64 title = wintitle + " Scope: %s" %(thisDatetime.strftime("%d-%b-%Y"))
65 title = wintitle + " Scope: %s" %(thisDatetime.strftime("%d-%b-%Y"))
65 self.xlabel = "Range (Km)"
66 self.xlabel = "Range (Km)"
@@ -81,24 +82,71 class ScopePlot(Plot):
81 #pass
82 #pass
82 ax.plt_r.set_data(x, ychannel)
83 ax.plt_r.set_data(x, ychannel)
83
84
85 def plot_weatherpower(self, x, y, channelIndexList, thisDatetime, wintitle):
84
86
85 def plot(self):
86
87
88 y = y[channelIndexList,:]
89 yreal = y.real
90 yreal = 10*numpy.log10(yreal)
91 self.y = yreal
92 title = wintitle + " Scope: %s" %(thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
93 self.xlabel = "Range (Km)"
94 self.ylabel = "Intensity"
95 self.xmin = min(x)
96 self.xmax = max(x)
97
98 self.titles[0] =title
99 for i,ax in enumerate(self.axes):
100 title = "Channel %d" %(i)
101
102 ychannel = yreal[i,:]
103
104 if ax.firsttime:
105 ax.plt_r = ax.plot(x, ychannel)[0]
106 else:
107 #pass
108 ax.plt_r.set_data(x, ychannel)
109
110 def plot_weathervelocity(self, x, y, channelIndexList, thisDatetime, wintitle):
111
112 x = x[channelIndexList,:]
113 yreal = y
114 self.y = yreal
115 title = wintitle + " Scope: %s" %(thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
116 self.xlabel = "Velocity (m/s)"
117 self.ylabel = "Range (Km)"
118 self.xmin = numpy.min(x)
119 self.xmax = numpy.max(x)
120 self.titles[0] =title
121 for i,ax in enumerate(self.axes):
122 title = "Channel %d" %(i)
123 xchannel = x[i,:]
124 if ax.firsttime:
125 ax.plt_r = ax.plot(xchannel, yreal)[0]
126 else:
127 #pass
128 ax.plt_r.set_data(xchannel, yreal)
129
130 def plot(self):
87 if self.channels:
131 if self.channels:
88 channels = self.channels
132 channels = self.channels
89 else:
133 else:
90 channels = self.data.channels
134 channels = self.data.channels
91
135
92 thisDatetime = datetime.datetime.utcfromtimestamp(self.data.times[-1])
136 thisDatetime = datetime.datetime.utcfromtimestamp(self.data.times[-1])
93
137 if self.CODE == "pp_power":
94 scope = self.data['scope']
138 scope = self.data['pp_power']
139 elif self.CODE == "pp_velocity":
140 scope = self.data["pp_velocity"]
141 else:
142 scope =self.data["scope"]
95
143
96 if self.data.flagDataAsBlock:
144 if self.data.flagDataAsBlock:
97
145
98 for i in range(self.data.nProfiles):
146 for i in range(self.data.nProfiles):
99
147
100 wintitle1 = " [Profile = %d] " %i
148 wintitle1 = " [Profile = %d] " %i
101
149 if self.CODE =="scope":
102 if self.type == "power":
150 if self.type == "power":
103 self.plot_power(self.data.heights,
151 self.plot_power(self.data.heights,
104 scope[:,i,:],
152 scope[:,i,:],
@@ -114,9 +162,23 class ScopePlot(Plot):
114 thisDatetime,
162 thisDatetime,
115 wintitle1
163 wintitle1
116 )
164 )
165 if self.CODE=="pp_power":
166 self.plot_weatherpower(self.data.heights,
167 scope[:,i,:],
168 channels,
169 thisDatetime,
170 wintitle
171 )
172 if self.CODE=="pp_velocity":
173 self.plot_weathervelocity(scope[:,i,:],
174 self.data.heights,
175 channels,
176 thisDatetime,
177 wintitle
178 )
117 else:
179 else:
118 wintitle = " [Profile = %d] " %self.data.profileIndex
180 wintitle = " [Profile = %d] " %self.data.profileIndex
119
181 if self.CODE== "scope":
120 if self.type == "power":
182 if self.type == "power":
121 self.plot_power(self.data.heights,
183 self.plot_power(self.data.heights,
122 scope,
184 scope,
@@ -132,3 +194,40 class ScopePlot(Plot):
132 thisDatetime,
194 thisDatetime,
133 wintitle
195 wintitle
134 )
196 )
197 if self.CODE=="pp_power":
198 self.plot_weatherpower(self.data.heights,
199 scope,
200 channels,
201 thisDatetime,
202 wintitle
203 )
204 if self.CODE=="pp_velocity":
205 self.plot_weathervelocity(scope,
206 self.data.heights,
207 channels,
208 thisDatetime,
209 wintitle
210 )
211
212
213
214 class PulsepairPowerPlot(ScopePlot):
215 '''
216 Plot for
217 '''
218
219 CODE = 'pp_power'
220 plot_name = 'PulsepairPower'
221 plot_type = 'scatter'
222 buffering = False
223
224
225
226 class PulsepairVelocityPlot(ScopePlot):
227 '''
228 Plot for
229 '''
230 CODE = 'pp_velocity'
231 plot_name = 'PulsepairVelocity'
232 plot_type = 'scatter'
233 buffering = False
Show file before | Show file after | Hide comments
@@ -1,5 +1,5
1 import sys
1 import sys
2 import numpy
2 import numpy,math
3 from scipy import interpolate
3 from scipy import interpolate
4 from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation, MPDecorator
4 from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation, MPDecorator
5 from schainpy.model.data.jrodata import Voltage
5 from schainpy.model.data.jrodata import Voltage
@@ -59,14 +59,12 class selectChannels(Operation):
59
59
60 channelIndexList = []
60 channelIndexList = []
61 self.dataOut = dataOut
61 self.dataOut = dataOut
62
63 for channel in channelList:
62 for channel in channelList:
64 if channel not in self.dataOut.channelList:
63 if channel not in self.dataOut.channelList:
65 raise ValueError("Channel %d is not in %s" %(channel, str(self.dataOut.channelList)))
64 raise ValueError("Channel %d is not in %s" %(channel, str(self.dataOut.channelList)))
66
65
67 index = self.dataOut.channelList.index(channel)
66 index = self.dataOut.channelList.index(channel)
68 channelIndexList.append(index)
67 channelIndexList.append(index)
69
70 self.selectChannelsByIndex(channelIndexList)
68 self.selectChannelsByIndex(channelIndexList)
71 return self.dataOut
69 return self.dataOut
72
70
@@ -105,6 +103,7 class selectChannels(Operation):
105 self.dataOut.data = data
103 self.dataOut.data = data
106 # self.dataOut.channelList = [self.dataOut.channelList[i] for i in channelIndexList]
104 # self.dataOut.channelList = [self.dataOut.channelList[i] for i in channelIndexList]
107 self.dataOut.channelList = range(len(channelIndexList))
105 self.dataOut.channelList = range(len(channelIndexList))
106
108 elif self.dataOut.type == 'Spectra':
107 elif self.dataOut.type == 'Spectra':
109 data_spc = self.dataOut.data_spc[channelIndexList, :]
108 data_spc = self.dataOut.data_spc[channelIndexList, :]
110 data_dc = self.dataOut.data_dc[channelIndexList, :]
109 data_dc = self.dataOut.data_dc[channelIndexList, :]
@@ -1272,6 +1271,155 class CombineProfiles(Operation):
1272 dataOut.ippSeconds *= n
1271 dataOut.ippSeconds *= n
1273
1272
1274 return dataOut
1273 return dataOut
1274
1275 class PulsePairVoltage(Operation):
1276 '''
1277 Function PulsePair(Signal Power, Velocity)
1278 The real component of Lag[0] provides Intensity Information
1279 The imag component of Lag[1] Phase provides Velocity Information
1280
1281 Configuration Parameters:
1282 nPRF = Number of Several PRF
1283 theta = Degree Azimuth angel Boundaries
1284
1285 Input:
1286 self.dataOut
1287 lag[N]
1288 Affected:
1289 self.dataOut.spc
1290 '''
1291 isConfig = False
1292 __profIndex = 0
1293 __initime = None
1294 __lastdatatime = None
1295 __buffer = None
1296 __buffer2 = []
1297 __buffer3 = None
1298 __dataReady = False
1299 n = None
1300 __nch = 0
1301 __nHeis = 0
1302 removeDC = False
1303 ipp = None
1304 lambda_ = 0
1305
1306 def __init__(self,**kwargs):
1307 Operation.__init__(self,**kwargs)
1308
1309 def setup(self, dataOut, n = None, removeDC=False):
1310 '''
1311 n= Numero de PRF's de entrada
1312 '''
1313 self.__initime = None
1314 self.__lastdatatime = 0
1315 self.__dataReady = False
1316 self.__buffer = 0
1317 self.__buffer2 = []
1318 self.__buffer3 = 0
1319 self.__profIndex = 0
1320
1321 self.__nch = dataOut.nChannels
1322 self.__nHeis = dataOut.nHeights
1323 self.removeDC = removeDC
1324 self.lambda_ = 3.0e8/(9345.0e6)
1325 self.ippSec = dataOut.ippSeconds
1326 self.nCohInt = dataOut.nCohInt
1327 print("IPPseconds",dataOut.ippSeconds)
1328
1329 print("ELVALOR DE n es:", n)
1330 if n == None:
1331 raise ValueError("n should be specified.")
1332
1333 if n != None:
1334 if n<2:
1335 raise ValueError("n should be greater than 2")
1336
1337 self.n = n
1338 self.__nProf = n
1339
1340 self.__buffer = numpy.zeros((dataOut.nChannels,
1341 n,
1342 dataOut.nHeights),
1343 dtype='complex')
1344
1345 def putData(self,data):
1346 '''
1347 Add a profile to he __buffer and increase in one the __profiel Index
1348 '''
1349 self.__buffer[:,self.__profIndex,:]= data
1350 self.__profIndex += 1
1351 return
1352
1353 def pushData(self):
1354 '''
1355 Return the PULSEPAIR and the profiles used in the operation
1356 Affected : self.__profileIndex
1357 '''
1358
1359 if self.removeDC==True:
1360 mean = numpy.mean(self.__buffer,1)
1361 tmp = mean.reshape(self.__nch,1,self.__nHeis)
1362 dc= numpy.tile(tmp,[1,self.__nProf,1])
1363 self.__buffer = self.__buffer - dc
1364
1365 data_intensity = numpy.sum(self.__buffer*numpy.conj(self.__buffer),1)/(self.n*self.nCohInt)#*self.nCohInt)
1366 pair1 = self.__buffer[:,1:,:]*numpy.conjugate(self.__buffer[:,:-1,:])
1367 angle = numpy.angle(numpy.sum(pair1,1))*180/(math.pi)
1368 #print(angle.shape)#print("__ANGLE__") #print("angle",angle[:,:10])
1369 data_velocity = (self.lambda_/(4*math.pi*self.ippSec))*numpy.angle(numpy.sum(pair1,1))#self.ippSec*self.nCohInt
1370 n = self.__profIndex
1371
1372 self.__buffer = numpy.zeros((self.__nch, self.__nProf,self.__nHeis), dtype='complex')
1373 self.__profIndex = 0
1374 return data_intensity,data_velocity,n
1375
1376 def pulsePairbyProfiles(self,data):
1377
1378 self.__dataReady = False
1379 data_intensity = None
1380 data_velocity = None
1381 self.putData(data)
1382 if self.__profIndex == self.n:
1383 data_intensity, data_velocity, n = self.pushData()
1384 self.__dataReady = True
1385
1386 return data_intensity, data_velocity
1387
1388 def pulsePairOp(self, data, datatime= None):
1389
1390 if self.__initime == None:
1391 self.__initime = datatime
1392
1393 data_intensity, data_velocity = self.pulsePairbyProfiles(data)
1394 self.__lastdatatime = datatime
1395
1396 if data_intensity is None:
1397 return None, None, None
1398
1399 avgdatatime = self.__initime
1400 deltatime = datatime - self.__lastdatatime
1401 self.__initime = datatime
1402
1403 return data_intensity, data_velocity, avgdatatime
1404
1405 def run(self, dataOut,n = None,removeDC= False, overlapping= False,**kwargs):
1406
1407 if not self.isConfig:
1408 self.setup(dataOut = dataOut, n = n , removeDC=removeDC , **kwargs)
1409 self.isConfig = True
1410 data_intensity, data_velocity, avgdatatime = self.pulsePairOp(dataOut.data, dataOut.utctime)
1411 dataOut.flagNoData = True
1412
1413 if self.__dataReady:
1414 dataOut.nCohInt *= self.n
1415 dataOut.data_intensity = data_intensity #valor para intensidad
1416 dataOut.data_velocity = data_velocity #valor para velocidad
1417 dataOut.PRFbyAngle = self.n #numero de PRF*cada angulo rotado que equivale a un tiempo.
1418 dataOut.utctime = avgdatatime
1419 dataOut.flagNoData = False
1420 return dataOut
1421
1422
1275 # import collections
1423 # import collections
1276 # from scipy.stats import mode
1424 # from scipy.stats import mode
1277 #
1425 #
General Comments 0
You need to be logged in to leave comments. Login now