@@ -308,6 +308,14 class ACFPlot(Figure): | |||
|
308 | 308 | y = dataOut.getHeiRange() |
|
309 | 309 | deltaHeight = dataOut.heightList[1]-dataOut.heightList[0] |
|
310 | 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 | 319 | #z = numpy.where(numpy.isfinite(z), z, numpy.NAN) |
|
312 | 320 | shape = dataOut.data_acf.shape |
|
313 | 321 | hei_index = numpy.arange(shape[2]) |
@@ -150,11 +150,50 class SpectraAFCProc(ProcessingUnit): | |||
|
150 | 150 | acf = data.imag |
|
151 | 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 | 185 | # Normalizando |
|
154 | 186 | for i in range(shape[0]): |
|
155 | 187 | for j in range(shape[2]): |
|
156 | 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 | 197 | self.dataOut.data_acf = acf |
|
159 | 198 | return True |
|
160 | 199 |
@@ -6,6 +6,70 from jroproc_base import ProcessingUnit, Operation | |||
|
6 | 6 | from schainpy.model.data.jrodata import Voltage |
|
7 | 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 | 73 | class VoltageProc(ProcessingUnit): |
|
10 | 74 | |
|
11 | 75 | |
@@ -1258,6 +1322,67 class SSheightProfiles(Operation): | |||
|
1258 | 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 | 1386 | # import collections |
|
1262 | 1387 | # from scipy.stats import mode |
|
1263 | 1388 | # |
General Comments 0
You need to be logged in to leave comments.
Login now