##// END OF EJS Templates
schain
RSS
UPDATE
avaldez -
r1246:3f5c93f33c3c
parent child
Show More
Side by Side Unified
Context file:
r1246:3f5c93f33c3c Loading diff...:
Show file before | Show file after | Hide comments
@@ -308,6 +308,14 class ACFPlot(Figure):
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])
Show file before | Show file after | Hide comments
@@ -150,11 +150,50 class SpectraAFCProc(ProcessingUnit):
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
Show file before | Show file after | Hide comments
@@ -6,6 +6,70 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
@@ -1258,6 +1322,67 class SSheightProfiles(Operation):
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 #
General Comments 0
You need to be logged in to leave comments. Login now