##// END OF EJS Templates
schain-jc
RSS
formatting y raise cuando digitalrf recibe spectra
José Chávez -
r1120:4cb1c6729a0e
parent child
Show More
Side by Side Unified
Context file:
r1120:4cb1c6729a0e Loading diff...:
Show file before | Show file after | Hide comments
@@ -52,7 +52,8 RADAR_STRUCTURE = numpy.dtype([
52 52 ('sRangeTxB','<a20'),
53 53 ])
54 54
55 SAMPLING_STRUCTURE = numpy.dtype([('h0','<f4'),('dh','<f4'),('nsa','<u4')])
55 SAMPLING_STRUCTURE = numpy.dtype(
56 [('h0', '<f4'), ('dh', '<f4'), ('nsa', '<u4')])
56 57
57 58
58 59 PROCESSING_STRUCTURE = numpy.dtype([
@@ -68,6 +69,7 PROCESSING_STRUCTURE = numpy.dtype([
68 69 ('nTotalSpectra','<u4')
69 70 ])
70 71
72
71 73 class Header(object):
72 74
73 75 def __init__(self):
@@ -122,6 +124,7 class Header(object):
122 124
123 125 print message
124 126
127
125 128 class BasicHeader(Header):
126 129
127 130 size = None
@@ -180,7 +183,8 class BasicHeader(Header):
180 183
181 184 def write(self, fp):
182 185
183 headerTuple = (self.size,self.version,self.dataBlock,self.utc,self.miliSecond,self.timeZone,self.dstFlag,self.errorCount)
186 headerTuple = (self.size, self.version, self.dataBlock, self.utc,
187 self.miliSecond, self.timeZone, self.dstFlag, self.errorCount)
184 188 header = numpy.array(headerTuple, BASIC_STRUCTURE)
185 189 header.tofile(fp)
186 190
@@ -201,6 +205,7 class BasicHeader(Header):
201 205 ltc = property(get_ltc, set_ltc)
202 206 datatime = property(get_datatime)
203 207
208
204 209 class SystemHeader(Header):
205 210
206 211 size = None
@@ -244,16 +249,17 class SystemHeader(Header):
244 249 self.adcResolution = header['nADCResolution'][0]
245 250 self.pciDioBusWidth = header['nPCDIOBusWidth'][0]
246 251
247
248 252 if startFp is not None:
249 253 endFp = self.size + startFp
250 254
251 255 if fp.tell() > endFp:
252 sys.stderr.write("Warning %s: Size value read from System Header is lower than it has to be\n" %fp.name)
256 sys.stderr.write(
257 "Warning %s: Size value read from System Header is lower than it has to be\n" % fp.name)
253 258 return 0
254 259
255 260 if fp.tell() < endFp:
256 sys.stderr.write("Warning %s: Size value read from System Header size is greater than it has to be\n" %fp.name)
261 sys.stderr.write(
262 "Warning %s: Size value read from System Header size is greater than it has to be\n" % fp.name)
257 263 return 0
258 264
259 265 self.length = header.nbytes
@@ -261,12 +267,14 class SystemHeader(Header):
261 267
262 268 def write(self, fp):
263 269
264 headerTuple = (self.size,self.nSamples,self.nProfiles,self.nChannels,self.adcResolution,self.pciDioBusWidth)
270 headerTuple = (self.size, self.nSamples, self.nProfiles,
271 self.nChannels, self.adcResolution, self.pciDioBusWidth)
265 272 header = numpy.array(headerTuple,SYSTEM_STRUCTURE)
266 273 header.tofile(fp)
267 274
268 275 return 1
269 276
277
270 278 class RadarControllerHeader(Header):
271 279
272 280 expType = None
@@ -370,9 +378,11 class RadarControllerHeader(Header):
370 378
371 379 try:
372 380 if hasattr(fp, 'read'):
373 samplingWindow = numpy.fromfile(fp, SAMPLING_STRUCTURE, self.nWindows)
381 samplingWindow = numpy.fromfile(
382 fp, SAMPLING_STRUCTURE, self.nWindows)
374 383 else:
375 samplingWindow = numpy.fromstring(fp[self.length:], SAMPLING_STRUCTURE, self.nWindows)
384 samplingWindow = numpy.fromstring(
385 fp[self.length:], SAMPLING_STRUCTURE, self.nWindows)
376 386 self.length += samplingWindow.nbytes
377 387 except Exception, e:
378 388 print "RadarControllerHeader: " + str(e)
@@ -382,20 +392,17 class RadarControllerHeader(Header):
382 392 self.deltaHeight = samplingWindow['dh']
383 393 self.samplesWin = samplingWindow['nsa']
384 394
385
386
387 395 try:
388 396 if hasattr(fp, 'read'):
389 397 self.Taus = numpy.fromfile(fp, '<f4', self.numTaus)
390 398 else:
391 self.Taus = numpy.fromstring(fp[self.length:], '<f4', self.numTaus)
399 self.Taus = numpy.fromstring(
400 fp[self.length:], '<f4', self.numTaus)
392 401 self.length += self.Taus.nbytes
393 402 except Exception, e:
394 403 print "RadarControllerHeader: " + str(e)
395 404 return 0
396 405
397
398
399 406 self.code_size = 0
400 407 if self.codeType != 0:
401 408
@@ -406,9 +413,11 class RadarControllerHeader(Header):
406 413 self.nBaud = numpy.fromfile(fp, '<u4', 1)[0]
407 414 self.length += self.nBaud.nbytes
408 415 else:
409 self.nCode = numpy.fromstring(fp[self.length:], '<u4', 1)[0]
416 self.nCode = numpy.fromstring(
417 fp[self.length:], '<u4', 1)[0]
410 418 self.length += self.nCode.nbytes
411 self.nBaud = numpy.fromstring(fp[self.length:], '<u4', 1)[0]
419 self.nBaud = numpy.fromstring(
420 fp[self.length:], '<u4', 1)[0]
412 421 self.length += self.nBaud.nbytes
413 422 except Exception, e:
414 423 print "RadarControllerHeader: " + str(e)
@@ -418,9 +427,11 class RadarControllerHeader(Header):
418 427 for ic in range(self.nCode):
419 428 try:
420 429 if hasattr(fp, 'read'):
421 temp = numpy.fromfile(fp,'u4', int(numpy.ceil(self.nBaud/32.)))
430 temp = numpy.fromfile(fp, 'u4', int(
431 numpy.ceil(self.nBaud / 32.)))
422 432 else:
423 temp = numpy.fromstring(fp,'u4', int(numpy.ceil(self.nBaud/32.)))
433 temp = numpy.fromstring(
434 fp, 'u4', int(numpy.ceil(self.nBaud / 32.)))
424 435 self.length += temp.nbytes
425 436 except Exception, e:
426 437 print "RadarControllerHeader: " + str(e)
@@ -447,12 +458,13 class RadarControllerHeader(Header):
447 458 # return 0
448 459
449 460 if fp.tell() > endFp:
450 sys.stderr.write("Warning %s: Size value read from Radar Controller header is lower than it has to be\n" %fp.name)
461 sys.stderr.write(
462 "Warning %s: Size value read from Radar Controller header is lower than it has to be\n" % fp.name)
451 463 # return 0
452 464
453 465 if fp.tell() < endFp:
454 sys.stderr.write("Warning %s: Size value read from Radar Controller header is greater than it has to be\n" %fp.name)
455
466 sys.stderr.write(
467 "Warning %s: Size value read from Radar Controller header is greater than it has to be\n" % fp.name)
456 468
457 469 return 1
458 470
@@ -479,7 +491,8 class RadarControllerHeader(Header):
479 491 header = numpy.array(headerTuple,RADAR_STRUCTURE)
480 492 header.tofile(fp)
481 493
482 sampleWindowTuple = (self.firstHeight,self.deltaHeight,self.samplesWin)
494 sampleWindowTuple = (
495 self.firstHeight, self.deltaHeight, self.samplesWin)
483 496 samplingWindow = numpy.array(sampleWindowTuple,SAMPLING_STRUCTURE)
484 497 samplingWindow.tofile(fp)
485 498
@@ -501,7 +514,8 class RadarControllerHeader(Header):
501 514 code_selected = code1[ic,start:end]
502 515 for j in range(len(code_selected)-1,-1,-1):
503 516 if code_selected[j] == 1:
504 tempx[i] = tempx[i] + 2**(len(code_selected)-1-j)
517 tempx[i] = tempx[i] + \
518 2**(len(code_selected) - 1 - j)
505 519 start = start + 32
506 520 end = end + 32
507 521
@@ -536,7 +550,8 class RadarControllerHeader(Header):
536 550 self.__size = 116 + 12*self.nWindows + 4*self.numTaus
537 551
538 552 if self.codeType != 0:
539 self.__size += 4 + 4 + 4*self.nCode*numpy.ceil(self.nBaud/32.)
553 self.__size += 4 + 4 + 4 * self.nCode * \
554 numpy.ceil(self.nBaud / 32.)
540 555
541 556 return self.__size
542 557
@@ -549,6 +564,7 class RadarControllerHeader(Header):
549 564 ippSeconds = property(get_ippSeconds, set_ippSeconds)
550 565 size = property(get_size, set_size)
551 566
567
552 568 class ProcessingHeader(Header):
553 569
554 570 # size = None
@@ -628,9 +644,11 class ProcessingHeader(Header):
628 644
629 645 try:
630 646 if hasattr(fp, 'read'):
631 samplingWindow = numpy.fromfile(fp, SAMPLING_STRUCTURE, self.nWindows)
647 samplingWindow = numpy.fromfile(
648 fp, SAMPLING_STRUCTURE, self.nWindows)
632 649 else:
633 samplingWindow = numpy.fromstring(fp[self.length:], SAMPLING_STRUCTURE, self.nWindows)
650 samplingWindow = numpy.fromstring(
651 fp[self.length:], SAMPLING_STRUCTURE, self.nWindows)
634 652 self.length += samplingWindow.nbytes
635 653 except Exception, e:
636 654 print "ProcessingHeader: " + str(e)
@@ -641,12 +659,13 class ProcessingHeader(Header):
641 659 self.deltaHeight = float(samplingWindow['dh'][0])
642 660 self.samplesWin = samplingWindow['nsa'][0]
643 661
644
645 662 try:
646 663 if hasattr(fp, 'read'):
647 self.spectraComb = numpy.fromfile(fp, 'u1', 2*self.totalSpectra)
664 self.spectraComb = numpy.fromfile(
665 fp, 'u1', 2 * self.totalSpectra)
648 666 else:
649 self.spectraComb = numpy.fromstring(fp[self.length:], 'u1', 2*self.totalSpectra)
667 self.spectraComb = numpy.fromstring(
668 fp[self.length:], 'u1', 2 * self.totalSpectra)
650 669 self.length += self.spectraComb.nbytes
651 670 except Exception, e:
652 671 print "ProcessingHeader: " + str(e)
@@ -655,7 +674,8 class ProcessingHeader(Header):
655 674 if ((self.processFlags & PROCFLAG.DEFINE_PROCESS_CODE) == PROCFLAG.DEFINE_PROCESS_CODE):
656 675 self.nCode = int(numpy.fromfile(fp,'<u4',1))
657 676 self.nBaud = int(numpy.fromfile(fp,'<u4',1))
658 self.code = numpy.fromfile(fp,'<f4',self.nCode*self.nBaud).reshape(self.nCode,self.nBaud)
677 self.code = numpy.fromfile(
678 fp, '<f4', self.nCode * self.nBaud).reshape(self.nCode, self.nBaud)
659 679
660 680 if ((self.processFlags & PROCFLAG.EXP_NAME_ESP) == PROCFLAG.EXP_NAME_ESP):
661 681 exp_name_len = int(numpy.fromfile(fp,'<u4',1))
@@ -696,16 +716,16 class ProcessingHeader(Header):
696 716 if nPairs > 0:
697 717 self.flag_cspc = True
698 718
699
700
701 719 if startFp is not None:
702 720 endFp = size + startFp
703 721 if fp.tell() > endFp:
704 sys.stderr.write("Warning: Processing header size is lower than it has to be")
722 sys.stderr.write(
723 "Warning: Processing header size is lower than it has to be")
705 724 return 0
706 725
707 726 if fp.tell() < endFp:
708 sys.stderr.write("Warning: Processing header size is greater than it is considered")
727 sys.stderr.write(
728 "Warning: Processing header size is greater than it is considered")
709 729
710 730 return 1
711 731
@@ -728,7 +748,8 class ProcessingHeader(Header):
728 748 header.tofile(fp)
729 749
730 750 if self.nWindows != 0:
731 sampleWindowTuple = (self.firstHeight,self.deltaHeight,self.samplesWin)
751 sampleWindowTuple = (
752 self.firstHeight, self.deltaHeight, self.samplesWin)
732 753 samplingWindow = numpy.array(sampleWindowTuple,SAMPLING_STRUCTURE)
733 754 samplingWindow.tofile(fp)
734 755
@@ -768,6 +789,7 class ProcessingHeader(Header):
768 789
769 790 size = property(get_size, set_size)
770 791
792
771 793 class RCfunction:
772 794 NONE=0
773 795 FLIP=1
@@ -776,6 +798,7 class RCfunction:
776 798 LIN6DIV256=4
777 799 SYNCHRO=5
778 800
801
779 802 class nCodeType:
780 803 NONE=0
781 804 USERDEFINE=1
@@ -796,6 +819,7 class nCodeType:
796 819 COMPLEMENTARYCODE128=16
797 820 CODE_BINARY28=17
798 821
822
799 823 class PROCFLAG:
800 824
801 825 COHERENT_INTEGRATION = numpy.uint32(0x00000001)
@@ -834,6 +858,7 class PROCFLAG:
834 858 DATAARRANGE_MASK = numpy.uint32(0x00007000)
835 859 ACQ_SYS_MASK = numpy.uint32(0x001C0000)
836 860
861
837 862 dtype0 = numpy.dtype([('real','<i1'),('imag','<i1')])
838 863 dtype1 = numpy.dtype([('real','<i2'),('imag','<i2')])
839 864 dtype2 = numpy.dtype([('real','<i4'),('imag','<i4')])
@@ -852,6 +877,7 PROCFLAG_DTYPE_LIST = [PROCFLAG.DATATYPE_CHAR,
852 877
853 878 DTYPE_WIDTH = [1, 2, 4, 8, 4, 8]
854 879
880
855 881 def get_dtype_index(numpy_dtype):
856 882
857 883 index = None
@@ -863,14 +889,17 def get_dtype_index(numpy_dtype):
863 889
864 890 return index
865 891
892
866 893 def get_numpy_dtype(index):
867 894
868 895 return NUMPY_DTYPE_LIST[index]
869 896
897
870 898 def get_procflag_dtype(index):
871 899
872 900 return PROCFLAG_DTYPE_LIST[index]
873 901
902
874 903 def get_dtype_width(index):
875 904
876 905 return DTYPE_WIDTH[index]
Show file before | Show file after | Hide comments
@@ -657,11 +657,11 class DigitalRFWriter(Operation):
657 657 self.metadata_dict['flagDecodeData'] = self.dataOut.flagDecodeData
658 658 self.metadata_dict['flagDeflipData'] = self.dataOut.flagDeflipData
659 659 self.metadata_dict['flagShiftFFT'] = self.dataOut.flagShiftFFT
660 self.metadata_dict['flagDataAsBlock'] = self.dataOut.flagDataAsBlock
661 660 self.metadata_dict['useLocalTime'] = self.dataOut.useLocalTime
662 661 self.metadata_dict['nCohInt'] = self.dataOut.nCohInt
663
664 return
662 self.metadata_dict['type'] = self.dataOut.type
663 self.metadata_dict['flagDataAsBlock'] = getattr(
664 self.dataOut, 'flagDataAsBlock', None) # chequear
665 665
666 666 def setup(self, dataOut, path, frequency, fileCadence, dirCadence, metadataCadence, set=0, metadataFile='metadata', ext='.h5'):
667 667 '''
@@ -678,9 +678,13 class DigitalRFWriter(Operation):
678 678 self.__dtype = dataOut.dtype[0]
679 679 self.__nSamples = dataOut.systemHeaderObj.nSamples
680 680 self.__nProfiles = dataOut.nProfiles
681 self.__blocks_per_file = dataOut.processingHeaderObj.dataBlocksPerFile
682 681
683 self.arr_data = arr_data = numpy.ones((self.__nSamples, len(
682 if self.dataOut.type != 'Voltage':
683 raise 'Digital RF cannot be used with this data type'
684 self.arr_data = numpy.ones((1, dataOut.nFFTPoints * len(
685 self.dataOut.channelList)), dtype=[('r', self.__dtype), ('i', self.__dtype)])
686 else:
687 self.arr_data = numpy.ones((self.__nSamples, len(
684 688 self.dataOut.channelList)), dtype=[('r', self.__dtype), ('i', self.__dtype)])
685 689
686 690 file_cadence_millisecs = 1000
@@ -702,14 +706,11 class DigitalRFWriter(Operation):
702 706 fileCadence, start_global_index,
703 707 sample_rate_numerator, sample_rate_denominator, uuid, compression_level, checksum,
704 708 is_complex, num_subchannels, is_continuous, marching_periods)
705
706 709 metadata_dir = os.path.join(path, 'metadata')
707 710 os.system('mkdir %s' % (metadata_dir))
708
709 711 self.digitalMetadataWriteObj = digital_rf.DigitalMetadataWriter(metadata_dir, dirCadence, 1, # 236, file_cadence_millisecs / 1000
710 712 sample_rate_numerator, sample_rate_denominator,
711 713 metadataFile)
712
713 714 self.isConfig = True
714 715 self.currentSample = 0
715 716 self.oldAverage = 0
@@ -717,7 +718,6 class DigitalRFWriter(Operation):
717 718 return
718 719
719 720 def writeMetadata(self):
720 print '[Writing] - Writing metadata'
721 721 start_idx = self.__sample_rate * self.dataOut.utctime
722 722
723 723 self.metadata_dict['processingHeader'] = self.dataOut.processingHeaderObj.getAsDict(
@@ -741,6 +741,15 class DigitalRFWriter(Operation):
741 741 return
742 742
743 743 def writeData(self):
744 if self.dataOut.type != 'Voltage':
745 raise 'Digital RF cannot be used with this data type'
746 for channel in self.dataOut.channelList:
747 for i in range(self.dataOut.nFFTPoints):
748 self.arr_data[1][channel * self.dataOut.nFFTPoints +
749 i]['r'] = self.dataOut.data[channel][i].real
750 self.arr_data[1][channel * self.dataOut.nFFTPoints +
751 i]['i'] = self.dataOut.data[channel][i].imag
752 else:
744 753 for i in range(self.dataOut.systemHeaderObj.nSamples):
745 754 for channel in self.dataOut.channelList:
746 755 self.arr_data[i][channel]['r'] = self.dataOut.data[channel][i].real
Show file before | Show file after | Hide comments
@@ -6,6 +6,7 from jroproc_base import ProcessingUnit, Operation
6 6 from schainpy.model.data.jrodata import Spectra
7 7 from schainpy.model.data.jrodata import hildebrand_sekhon
8 8
9
9 10 class SpectraProc(ProcessingUnit):
10 11
11 12 def __init__(self, **kwargs):
@@ -25,7 +26,10 class SpectraProc(ProcessingUnit):
25 26 self.dataOut.dstFlag = self.dataIn.dstFlag
26 27 self.dataOut.errorCount = self.dataIn.errorCount
27 28 self.dataOut.useLocalTime = self.dataIn.useLocalTime
28
29 try:
30 self.dataOut.processingHeaderObj = self.dataIn.processingHeaderObj.copy()
31 except:
32 pass
29 33 self.dataOut.radarControllerHeaderObj = self.dataIn.radarControllerHeaderObj.copy()
30 34 self.dataOut.systemHeaderObj = self.dataIn.systemHeaderObj.copy()
31 35 self.dataOut.channelList = self.dataIn.channelList
@@ -39,8 +43,10 class SpectraProc(ProcessingUnit):
39 43
40 44 self.dataOut.flagDiscontinuousBlock = self.dataIn.flagDiscontinuousBlock
41 45 self.dataOut.utctime = self.firstdatatime
42 self.dataOut.flagDecodeData = self.dataIn.flagDecodeData #asumo q la data esta decodificada
43 self.dataOut.flagDeflipData = self.dataIn.flagDeflipData #asumo q la data esta sin flip
46 # asumo q la data esta decodificada
47 self.dataOut.flagDecodeData = self.dataIn.flagDecodeData
48 # asumo q la data esta sin flip
49 self.dataOut.flagDeflipData = self.dataIn.flagDeflipData
44 50 self.dataOut.flagShiftFFT = False
45 51
46 52 self.dataOut.nCohInt = self.dataIn.nCohInt
@@ -71,7 +77,8 class SpectraProc(ProcessingUnit):
71 77 self.buffer
72 78 self.dataOut.flagNoData
73 79 """
74 fft_volt = numpy.fft.fft(self.buffer,n=self.dataOut.nFFTPoints,axis=1)
80 fft_volt = numpy.fft.fft(
81 self.buffer, n=self.dataOut.nFFTPoints, axis=1)
75 82 fft_volt = fft_volt.astype(numpy.dtype('complex'))
76 83 dc = fft_volt[:,0,:]
77 84
@@ -88,14 +95,18 class SpectraProc(ProcessingUnit):
88 95 pairIndex = 0
89 96 if self.dataOut.pairsList != None:
90 97 #calculo de cross-spectra
91 cspc = numpy.zeros((self.dataOut.nPairs, self.dataOut.nFFTPoints, self.dataOut.nHeights), dtype='complex')
98 cspc = numpy.zeros(
99 (self.dataOut.nPairs, self.dataOut.nFFTPoints, self.dataOut.nHeights), dtype='complex')
92 100 for pair in self.dataOut.pairsList:
93 101 if pair[0] not in self.dataOut.channelList:
94 raise ValueError, "Error getting CrossSpectra: pair 0 of %s is not in channelList = %s" %(str(pair), str(self.dataOut.channelList))
102 raise ValueError, "Error getting CrossSpectra: pair 0 of %s is not in channelList = %s" % (
103 str(pair), str(self.dataOut.channelList))
95 104 if pair[1] not in self.dataOut.channelList:
96 raise ValueError, "Error getting CrossSpectra: pair 1 of %s is not in channelList = %s" %(str(pair), str(self.dataOut.channelList))
105 raise ValueError, "Error getting CrossSpectra: pair 1 of %s is not in channelList = %s" % (
106 str(pair), str(self.dataOut.channelList))
97 107
98 cspc[pairIndex,:,:] = fft_volt[pair[0],:,:] * numpy.conjugate(fft_volt[pair[1],:,:])
108 cspc[pairIndex, :, :] = fft_volt[pair[0], :, :] * \
109 numpy.conjugate(fft_volt[pair[1], :, :])
99 110 pairIndex += 1
100 111 blocksize += cspc.size
101 112
@@ -154,12 +165,14 class SpectraProc(ProcessingUnit):
154 165 self.id_min = 0
155 166 self.id_max = nVoltProfiles
156 167
157 self.buffer[:,self.id_min:self.id_max,:] = self.dataIn.data
168 self.buffer[:, self.id_min:self.id_max,
169 :] = self.dataIn.data
158 170 self.profIndex += nVoltProfiles
159 171 self.id_min += nVoltProfiles
160 172 self.id_max += nVoltProfiles
161 173 else:
162 raise ValueError, "The type object %s has %d profiles, it should just has %d profiles"%(self.dataIn.type,self.dataIn.data.shape[1],nProfiles)
174 raise ValueError, "The type object %s has %d profiles, it should just has %d profiles" % (
175 self.dataIn.type, self.dataIn.data.shape[1], nProfiles)
163 176 self.dataOut.flagNoData = True
164 177 return 0
165 178 else:
@@ -179,7 +192,8 class SpectraProc(ProcessingUnit):
179 192
180 193 return True
181 194
182 raise ValueError, "The type of input object '%s' is not valid"%(self.dataIn.type)
195 raise ValueError, "The type of input object '%s' is not valid" % (
196 self.dataIn.type)
183 197
184 198 def __selectPairs(self, pairsList):
185 199
@@ -222,7 +236,8 class SpectraProc(ProcessingUnit):
222 236 return
223 237
224 238 self.dataOut.data_cspc = self.dataOut.data_cspc[pairsIndexListSelected]
225 self.dataOut.pairsList = [self.dataOut.pairsList[i] for i in pairsIndexListSelected]
239 self.dataOut.pairsList = [self.dataOut.pairsList[i]
240 for i in pairsIndexListSelected]
226 241
227 242 return
228 243
@@ -232,7 +247,8 class SpectraProc(ProcessingUnit):
232 247
233 248 for channel in channelList:
234 249 if channel not in self.dataOut.channelList:
235 raise ValueError, "Error selecting channels, Channel %d is not valid.\nAvailable channels = %s" %(channel, str(self.dataOut.channelList))
250 raise ValueError, "Error selecting channels, Channel %d is not valid.\nAvailable channels = %s" % (
251 channel, str(self.dataOut.channelList))
236 252
237 253 index = self.dataOut.channelList.index(channel)
238 254 channelIndexList.append(index)
@@ -257,7 +273,8 class SpectraProc(ProcessingUnit):
257 273
258 274 for channelIndex in channelIndexList:
259 275 if channelIndex not in self.dataOut.channelIndexList:
260 raise ValueError, "Error selecting channels: The value %d in channelIndexList is not valid.\nAvailable channel indexes = " %(channelIndex, self.dataOut.channelIndexList)
276 raise ValueError, "Error selecting channels: The value %d in channelIndexList is not valid.\nAvailable channel indexes = " % (
277 channelIndex, self.dataOut.channelIndexList)
261 278
262 279 # nChannels = len(channelIndexList)
263 280
@@ -267,7 +284,8 class SpectraProc(ProcessingUnit):
267 284 self.dataOut.data_spc = data_spc
268 285 self.dataOut.data_dc = data_dc
269 286
270 self.dataOut.channelList = [self.dataOut.channelList[i] for i in channelIndexList]
287 self.dataOut.channelList = [
288 self.dataOut.channelList[i] for i in channelIndexList]
271 289 # self.dataOut.nChannels = nChannels
272 290
273 291 self.__selectPairsByChannel(self.dataOut.channelList)
@@ -291,7 +309,8 class SpectraProc(ProcessingUnit):
291 309 """
292 310
293 311 if (minHei > maxHei):
294 raise ValueError, "Error selecting heights: Height range (%d,%d) is not valid" % (minHei, maxHei)
312 raise ValueError, "Error selecting heights: Height range (%d,%d) is not valid" % (
313 minHei, maxHei)
295 314
296 315 if (minHei < self.dataOut.heightList[0]):
297 316 minHei = self.dataOut.heightList[0]
@@ -321,7 +340,8 class SpectraProc(ProcessingUnit):
321 340 return 1
322 341
323 342 def getBeaconSignal(self, tauindex = 0, channelindex = 0, hei_ref=None):
324 newheis = numpy.where(self.dataOut.heightList>self.dataOut.radarControllerHeaderObj.Taus[tauindex])
343 newheis = numpy.where(
344 self.dataOut.heightList > self.dataOut.radarControllerHeaderObj.Taus[tauindex])
325 345
326 346 if hei_ref != None:
327 347 newheis = numpy.where(self.dataOut.heightList>hei_ref)
@@ -332,8 +352,10 class SpectraProc(ProcessingUnit):
332 352 heightList = self.dataOut.heightList[minIndex:maxIndex+1]
333 353
334 354 # determina indices
335 nheis = int(self.dataOut.radarControllerHeaderObj.txB/(self.dataOut.heightList[1]-self.dataOut.heightList[0]))
336 avg_dB = 10*numpy.log10(numpy.sum(data_spc[channelindex,:,:],axis=0))
355 nheis = int(self.dataOut.radarControllerHeaderObj.txB /
356 (self.dataOut.heightList[1] - self.dataOut.heightList[0]))
357 avg_dB = 10 * \
358 numpy.log10(numpy.sum(data_spc[channelindex, :, :], axis=0))
337 359 beacon_dB = numpy.sort(avg_dB)[-nheis:]
338 360 beacon_heiIndexList = []
339 361 for val in avg_dB.tolist():
@@ -359,7 +381,6 class SpectraProc(ProcessingUnit):
359 381
360 382 return 1
361 383
362
363 384 def selectHeightsByIndex(self, minIndex, maxIndex):
364 385 """
365 386 Selecciona un bloque de datos en base a un grupo indices de alturas segun el rango
@@ -380,7 +401,8 class SpectraProc(ProcessingUnit):
380 401 """
381 402
382 403 if (minIndex < 0) or (minIndex > maxIndex):
383 raise ValueError, "Error selecting heights: Index range (%d,%d) is not valid" % (minIndex, maxIndex)
404 raise ValueError, "Error selecting heights: Index range (%d,%d) is not valid" % (
405 minIndex, maxIndex)
384 406
385 407 if (maxIndex >= self.dataOut.nHeights):
386 408 maxIndex = self.dataOut.nHeights-1
@@ -408,14 +430,14 class SpectraProc(ProcessingUnit):
408 430 jspectra = self.dataOut.data_spc
409 431 jcspectra = self.dataOut.data_cspc
410 432
411
412 433 num_chan = jspectra.shape[0]
413 434 num_hei = jspectra.shape[2]
414 435
415 436 if jcspectra is not None:
416 437 jcspectraExist = True
417 438 num_pairs = jcspectra.shape[0]
418 else: jcspectraExist = False
439 else:
440 jcspectraExist = False
419 441
420 442 freq_dc = jspectra.shape[1]/2
421 443 ind_vel = numpy.array([-2,-1,1,2]) + freq_dc
@@ -424,10 +446,12 class SpectraProc(ProcessingUnit):
424 446 ind_vel[range(0,1)] = ind_vel[range(0,1)] + self.num_prof
425 447
426 448 if mode == 1:
427 jspectra[:,freq_dc,:] = (jspectra[:,ind_vel[1],:] + jspectra[:,ind_vel[2],:])/2 #CORRECCION
449 jspectra[:, freq_dc, :] = (
450 jspectra[:, ind_vel[1], :] + jspectra[:, ind_vel[2], :]) / 2 # CORRECCION
428 451
429 452 if jcspectraExist:
430 jcspectra[:,freq_dc,:] = (jcspectra[:,ind_vel[1],:] + jcspectra[:,ind_vel[2],:])/2
453 jcspectra[:, freq_dc, :] = (
454 jcspectra[:, ind_vel[1], :] + jcspectra[:, ind_vel[2], :]) / 2
431 455
432 456 if mode == 2:
433 457
@@ -448,14 +472,14 class SpectraProc(ProcessingUnit):
448 472 cjunkid = sum(junkid)
449 473
450 474 if cjunkid.any():
451 jspectra[ich,freq_dc,junkid.nonzero()] = (jspectra[ich,ind_vel[1],junkid] + jspectra[ich,ind_vel[2],junkid])/2
475 jspectra[ich, freq_dc, junkid.nonzero()] = (
476 jspectra[ich, ind_vel[1], junkid] + jspectra[ich, ind_vel[2], junkid]) / 2
452 477
453 478 if jcspectraExist:
454 479 for ip in range(num_pairs):
455 480 yy = jcspectra[ip,ind_vel,:]
456 481 jcspectra[ip,freq_dc,:] = numpy.dot(xx_aux,yy)
457 482
458
459 483 self.dataOut.data_spc = jspectra
460 484 self.dataOut.data_cspc = jcspectra
461 485
@@ -494,7 +518,6 class SpectraProc(ProcessingUnit):
494 518 num_mask_prof = mask_prof.size
495 519 comp_mask_prof = [0, num_prof/2]
496 520
497
498 521 #noise_exist: Determina si la variable jnoise ha sido definida y contiene la informacion del ruido de cada canal
499 522 if (jnoise.size < num_channel or numpy.isnan(jnoise).any()):
500 523 jnoise = numpy.nan
@@ -509,7 +532,8 class SpectraProc(ProcessingUnit):
509 532 psort = power.ravel().argsort()
510 533
511 534 #Se estima la interferencia promedio en los Espectros de Potencia empleando
512 junkspc_interf = jspectra[ich,:,hei_interf[psort[range(offhei_interf, nhei_interf + offhei_interf)]]]
535 junkspc_interf = jspectra[ich, :, hei_interf[psort[range(
536 offhei_interf, nhei_interf + offhei_interf)]]]
513 537
514 538 if noise_exist:
515 539 # tmp_noise = jnoise[ich] / num_prof
@@ -520,43 +544,51 class SpectraProc(ProcessingUnit):
520 544 jspc_interf = junkspc_interf.sum(axis = 0) / nhei_interf
521 545 jspc_interf = jspc_interf.transpose()
522 546 #Calculando el espectro de interferencia promedio
523 noiseid = numpy.where(jspc_interf <= tmp_noise/ numpy.sqrt(num_incoh))
547 noiseid = numpy.where(
548 jspc_interf <= tmp_noise / numpy.sqrt(num_incoh))
524 549 noiseid = noiseid[0]
525 550 cnoiseid = noiseid.size
526 interfid = numpy.where(jspc_interf > tmp_noise/ numpy.sqrt(num_incoh))
551 interfid = numpy.where(
552 jspc_interf > tmp_noise / numpy.sqrt(num_incoh))
527 553 interfid = interfid[0]
528 554 cinterfid = interfid.size
529 555
530 if (cnoiseid > 0): jspc_interf[noiseid] = 0
556 if (cnoiseid > 0):
557 jspc_interf[noiseid] = 0
531 558
532 559 #Expandiendo los perfiles a limpiar
533 560 if (cinterfid > 0):
534 new_interfid = (numpy.r_[interfid - 1, interfid, interfid + 1] + num_prof)%num_prof
561 new_interfid = (
562 numpy.r_[interfid - 1, interfid, interfid + 1] + num_prof) % num_prof
535 563 new_interfid = numpy.asarray(new_interfid)
536 564 new_interfid = {x for x in new_interfid}
537 565 new_interfid = numpy.array(list(new_interfid))
538 566 new_cinterfid = new_interfid.size
539 else: new_cinterfid = 0
567 else:
568 new_cinterfid = 0
540 569
541 570 for ip in range(new_cinterfid):
542 571 ind = junkspc_interf[:,new_interfid[ip]].ravel().argsort()
543 jspc_interf[new_interfid[ip]] = junkspc_interf[ind[nhei_interf/2],new_interfid[ip]]
544
572 jspc_interf[new_interfid[ip]
573 ] = junkspc_interf[ind[nhei_interf / 2], new_interfid[ip]]
545 574
546 jspectra[ich,:,ind_hei] = jspectra[ich,:,ind_hei] - jspc_interf #Corregir indices
575 jspectra[ich, :, ind_hei] = jspectra[ich, :,
576 ind_hei] - jspc_interf # Corregir indices
547 577
548 578 #Removiendo la interferencia del punto de mayor interferencia
549 579 ListAux = jspc_interf[mask_prof].tolist()
550 580 maxid = ListAux.index(max(ListAux))
551 581
552
553 582 if cinterfid > 0:
554 583 for ip in range(cinterfid*(interf == 2) - 1):
555 ind = (jspectra[ich,interfid[ip],:] < tmp_noise*(1 + 1/numpy.sqrt(num_incoh))).nonzero()
584 ind = (jspectra[ich, interfid[ip], :] < tmp_noise *
585 (1 + 1 / numpy.sqrt(num_incoh))).nonzero()
556 586 cind = len(ind)
557 587
558 588 if (cind > 0):
559 jspectra[ich,interfid[ip],ind] = tmp_noise*(1 + (numpy.random.uniform(cind) - 0.5)/numpy.sqrt(num_incoh))
589 jspectra[ich, interfid[ip], ind] = tmp_noise * \
590 (1 + (numpy.random.uniform(cind) - 0.5) /
591 numpy.sqrt(num_incoh))
560 592
561 593 ind = numpy.array([-2,-1,1,2])
562 594 xx = numpy.zeros([4,4])
@@ -568,14 +600,17 class SpectraProc(ProcessingUnit):
568 600 xx = xx_inv[:,0]
569 601 ind = (ind + maxid + num_mask_prof)%num_mask_prof
570 602 yy = jspectra[ich,mask_prof[ind],:]
571 jspectra[ich,mask_prof[maxid],:] = numpy.dot(yy.transpose(),xx)
603 jspectra[ich, mask_prof[maxid], :] = numpy.dot(
604 yy.transpose(), xx)
572 605
573
574 indAux = (jspectra[ich,:,:] < tmp_noise*(1-1/numpy.sqrt(num_incoh))).nonzero()
575 jspectra[ich,indAux[0],indAux[1]] = tmp_noise * (1 - 1/numpy.sqrt(num_incoh))
606 indAux = (jspectra[ich, :, :] < tmp_noise *
607 (1 - 1 / numpy.sqrt(num_incoh))).nonzero()
608 jspectra[ich, indAux[0], indAux[1]] = tmp_noise * \
609 (1 - 1 / numpy.sqrt(num_incoh))
576 610
577 611 #Remocion de Interferencia en el Cross Spectra
578 if jcspectra is None: return jspectra, jcspectra
612 if jcspectra is None:
613 return jspectra, jcspectra
579 614 num_pairs = jcspectra.size/(num_prof*num_hei)
580 615 jcspectra = jcspectra.reshape(num_pairs, num_prof, num_hei)
581 616
@@ -588,22 +623,28 class SpectraProc(ProcessingUnit):
588 623 cspower = cspower.sum(axis = 0)
589 624
590 625 cspsort = cspower.ravel().argsort()
591 junkcspc_interf = jcspectra[ip,:,hei_interf[cspsort[range(offhei_interf, nhei_interf + offhei_interf)]]]
626 junkcspc_interf = jcspectra[ip, :, hei_interf[cspsort[range(
627 offhei_interf, nhei_interf + offhei_interf)]]]
592 628 junkcspc_interf = junkcspc_interf.transpose()
593 629 jcspc_interf = junkcspc_interf.sum(axis = 1)/nhei_interf
594 630
595 631 ind = numpy.abs(jcspc_interf[mask_prof]).ravel().argsort()
596 632
597 median_real = numpy.median(numpy.real(junkcspc_interf[mask_prof[ind[range(3*num_prof/4)]],:]))
598 median_imag = numpy.median(numpy.imag(junkcspc_interf[mask_prof[ind[range(3*num_prof/4)]],:]))
599 junkcspc_interf[comp_mask_prof,:] = numpy.complex(median_real, median_imag)
633 median_real = numpy.median(numpy.real(
634 junkcspc_interf[mask_prof[ind[range(3 * num_prof / 4)]], :]))
635 median_imag = numpy.median(numpy.imag(
636 junkcspc_interf[mask_prof[ind[range(3 * num_prof / 4)]], :]))
637 junkcspc_interf[comp_mask_prof, :] = numpy.complex(
638 median_real, median_imag)
600 639
601 640 for iprof in range(num_prof):
602 641 ind = numpy.abs(junkcspc_interf[iprof,:]).ravel().argsort()
603 jcspc_interf[iprof] = junkcspc_interf[iprof, ind[nhei_interf/2]]
642 jcspc_interf[iprof] = junkcspc_interf[iprof,
643 ind[nhei_interf / 2]]
604 644
605 645 #Removiendo la Interferencia
606 jcspectra[ip,:,ind_hei] = jcspectra[ip,:,ind_hei] - jcspc_interf
646 jcspectra[ip, :, ind_hei] = jcspectra[ip,
647 :, ind_hei] - jcspc_interf
607 648
608 649 ListAux = numpy.abs(jcspc_interf[mask_prof]).tolist()
609 650 maxid = ListAux.index(max(ListAux))
@@ -690,7 +731,8 class SpectraProc(ProcessingUnit):
690 731 maxIndex = len(heights)
691 732
692 733 if (minIndex < 0) or (minIndex > maxIndex):
693 raise ValueError, "some value in (%d,%d) is not valid" % (minIndex, maxIndex)
734 raise ValueError, "some value in (%d,%d) is not valid" % (
735 minIndex, maxIndex)
694 736
695 737 if (maxIndex >= self.dataOut.nHeights):
696 738 maxIndex = self.dataOut.nHeights-1
@@ -709,7 +751,8 class SpectraProc(ProcessingUnit):
709 751 maxIndexVel = len(velrange)
710 752
711 753 #seleccion del espectro
712 data_spc = self.dataOut.data_spc[:,minIndexVel:maxIndexVel+1,minIndex:maxIndex+1]
754 data_spc = self.dataOut.data_spc[:,
755 minIndexVel:maxIndexVel + 1, minIndex:maxIndex + 1]
713 756 #estimacion de ruido
714 757 noise = numpy.zeros(self.dataOut.nChannels)
715 758
@@ -721,8 +764,8 class SpectraProc(ProcessingUnit):
721 764
722 765 return 1
723 766
724 class IncohInt(Operation):
725 767
768 class IncohInt(Operation):
726 769
727 770 __profIndex = 0
728 771 __withOverapping = False
@@ -742,8 +785,6 class IncohInt(Operation):
742 785
743 786 n = None
744 787
745
746
747 788 def __init__(self, **kwargs):
748 789
749 790 Operation.__init__(self, **kwargs)
@@ -778,12 +819,12 class IncohInt(Operation):
778 819 if n is not None:
779 820 self.n = int(n)
780 821 else:
781 self.__integrationtime = int(timeInterval) #if (type(timeInterval)!=integer) -> change this line
822 # if (type(timeInterval)!=integer) -> change this line
823 self.__integrationtime = int(timeInterval)
782 824 self.n = None
783 825 self.__byTime = True
784 826
785 827 def putData(self, data_spc, data_cspc, data_dc):
786
787 828 """
788 829 Add a profile to the __buffer_spc and increase in one the __profileIndex
789 830
@@ -866,7 +907,8 class IncohInt(Operation):
866 907 self.__initime = datatime
867 908
868 909 if self.__byTime:
869 avgdata_spc, avgdata_cspc, avgdata_dc = self.byTime(datatime, *args)
910 avgdata_spc, avgdata_cspc, avgdata_dc = self.byTime(
911 datatime, *args)
870 912 else:
871 913 avgdata_spc, avgdata_cspc, avgdata_dc = self.byProfiles(*args)
872 914
General Comments 0
You need to be logged in to leave comments. Login now