##// END OF EJS Templates
nuevo
joabAM -
r1472:c5c71a942112
parent child
Show More
@@ -49,7 +49,55 static PyObject *hildebrand_sekhon(PyObject *self, PyObject *args) {
49 49 // return PyLong_FromLong(lnoise);
50 50 return PyFloat_FromDouble(lnoise);
51 51 }
52 /*
53 static PyObject *hildebrand_sekhon2(PyObject *self, PyObject *args) {
54 double navg;
55 double th;
56 PyObject *data_obj, *data_array;
57
58 if (!PyArg_ParseTuple(args, "Od", &data_obj, &navg, &th)) {
59 return NULL;
60 }
61
62 data_array = PyArray_FROM_OTF(data_obj, NPY_FLOAT64, NPY_IN_ARRAY);
63
64 if (data_array == NULL) {
65 Py_XDECREF(data_array);
66 Py_XDECREF(data_obj);
67 return NULL;
68 }
69 double *sortdata = (double*)PyArray_DATA(data_array);
70 int lenOfData = (int)PyArray_SIZE(data_array) ;
71 double nums_min = lenOfData*th;
72 if (nums_min <= 5) nums_min = 5;
73 double sump = 0;
74 double sumq = 0;
75 long j = 0;
76 int cont = 1;
77 double rtest = 0;
78 while ((cont == 1) && (j < lenOfData)) {
79 sump = sump + sortdata[j];
80 sumq = sumq + pow(sortdata[j], 2);
81 if (j > nums_min) {
82 rtest = (double)j/(j-1) + 1/navg;
83 if ((sumq*j) > (rtest*pow(sump, 2))) {
84 j = j - 1;
85 sump = sump - sortdata[j];
86 sumq = sumq - pow(sortdata[j],2);
87 cont = 0;
88 }
89 }
90 j = j + 1;
91 }
52 92
93 //double lnoise = sump / j;
94
95 Py_DECREF(data_array);
96
97 // return PyLong_FromLong(lnoise);
98 return PyFloat_FromDouble(j,sortID);
99 }
100 */
53 101
54 102 static PyMethodDef noiseMethods[] = {
55 103 { "hildebrand_sekhon", hildebrand_sekhon, METH_VARARGS, "Get noise with hildebrand_sekhon algorithm" },
@@ -50,10 +50,12 class SnrPlot(RTIPlot):
50 50 def update(self, dataOut):
51 51 if len(self.channelList) == 0:
52 52 self.update_list(dataOut)
53 data = {}
53
54 54 meta = {}
55 data['snr'] = 10*numpy.log10(dataOut.data_snr)
56
55 data = {
56 'snr': 10 * numpy.log10(dataOut.data_snr)
57 }
58 #print(data['snr'])
57 59 return data, meta
58 60
59 61 class DopplerPlot(RTIPlot):
@@ -258,6 +258,7 class RTIPlot(Plot):
258 258
259 259 self.x = self.data.times
260 260 self.y = self.data.yrange
261
261 262 self.z = self.data[self.CODE]
262 263 self.z = numpy.array(self.z, dtype=float)
263 264 self.z = numpy.ma.masked_invalid(self.z)
@@ -274,11 +275,13 class RTIPlot(Plot):
274 275 x, y, z = self.fill_gaps(self.x, self.y, self.z)
275 276 else:
276 277 x, y, z = self.fill_gaps(*self.decimate())
277 dummy_var = self.axes #ExtraΓ±amente esto actualiza el valor axes
278
279 #dummy_var = self.axes #ExtraΓ±amente esto actualiza el valor axes
278 280 for n, ax in enumerate(self.axes):
279 281 self.zmin = self.zmin if self.zmin else numpy.min(self.z)
280 282 self.zmax = self.zmax if self.zmax else numpy.max(self.z)
281 283 data = self.data[-1]
284
282 285 if ax.firsttime:
283 286 ax.plt = ax.pcolormesh(x, y, z[n].T,
284 287 vmin=self.zmin,
@@ -287,8 +290,8 class RTIPlot(Plot):
287 290 )
288 291 if self.showprofile:
289 292 ax.plot_profile = self.pf_axes[n].plot(data[self.CODE][n], self.y)[0]
290
291 293 if "noise" in self.data:
294
292 295 ax.plot_noise = self.pf_axes[n].plot(numpy.repeat(data['noise'][n], len(self.y)), self.y,
293 296 color="k", linestyle="dashed", lw=1)[0]
294 297 else:
@@ -301,6 +304,7 class RTIPlot(Plot):
301 304 if self.showprofile:
302 305 ax.plot_profile.set_data(data[self.CODE][n], self.y)
303 306 if "noise" in self.data:
307
304 308 ax.plot_noise.set_data(numpy.repeat(
305 309 data['noise'][n], len(self.y)), self.y)
306 310
@@ -308,8 +308,8 class HDFReader(Reader, ProcessingUnit):
308 308
309 309 self.getData()
310 310
311 #@MPDecorator
312 class HDFWrite(Operation):
311 @MPDecorator
312 class HDFWriter(Operation):
313 313 """Operation to write HDF5 files.
314 314
315 315 The HDF5 file contains by default two groups Data and Metadata where
@@ -465,7 +465,7 class HDFWrite(Operation):
465 465
466 466 def run(self, dataOut,**kwargs):
467 467
468 self.dataOut = dataOut
468 self.dataOut = dataOut.copy()
469 469 if not(self.isConfig):
470 470 self.setup(**kwargs)
471 471
@@ -474,7 +474,7 class HDFWrite(Operation):
474 474
475 475 self.putData()
476 476
477 return self.dataOut
477 #return self.dataOut
478 478
479 479 def setNextFile(self):
480 480
@@ -89,7 +89,7 class ProcessingUnit(object):
89 89 elif optype == 'external' and self.dataOut.error:
90 90 op.queue.put(copy.deepcopy(self.dataOut))
91 91
92 return 'Error' if self.dataOut.error else True#self.dataOut.isReady()
92 return 'Error' if self.dataOut.error else self.dataOut.isReady()
93 93
94 94 def setup(self):
95 95
@@ -202,7 +202,7 class SpectraProc(ProcessingUnit):
202 202 self.dataOut.flagNoData = False
203 203 self.firstdatatime = None
204 204 self.profIndex = 0
205 self.dataOut.noise_estimation = None
205
206 206 else:
207 207 raise ValueError("The type of input object '%s' is not valid".format(
208 208 self.dataIn.type))
@@ -490,13 +490,13 class removeDC(Operation):
490 490 return self.dataOut
491 491
492 492 class getNoise(Operation):
493
493 494 def __init__(self):
494 495
495 496 Operation.__init__(self)
496 497
497 def run(self, dataOut, minHei=None, maxHei=None, minVel=None, maxVel=None, minFreq= None, maxFreq=None,):
498 self.dataOut = dataOut.copy()
499 print("1: ",dataOut.noise_estimation, dataOut.normFactor)
498 def run(self, dataOut, minHei=None, maxHei=None,minVel=None, maxVel=None, minFreq= None, maxFreq=None):
499 self.dataOut = dataOut
500 500
501 501 if minHei == None:
502 502 minHei = self.dataOut.heightList[0]
@@ -606,9 +606,10 class getNoise(Operation):
606 606 maxIndexFFT = len( self.dataOut.getFreqRange(1))
607 607
608 608 #print(minIndex, maxIndex,minIndexVel, maxIndexVel)
609 self.dataOut.noise_estimation = None
609 610 noise = self.dataOut.getNoise(xmin_index=minIndexFFT, xmax_index=maxIndexFFT, ymin_index=minIndex, ymax_index=maxIndex)
610 611
611 self.dataOut.noise_estimation = noise.copy()
612 self.dataOut.noise_estimation = noise.copy() # dataOut.noise
612 613 #print("2: ",10*numpy.log10(self.dataOut.noise_estimation/64))
613 614 return self.dataOut
614 615
@@ -1201,10 +1202,14 class IntegrationFaradaySpectra(Operation):
1201 1202 buffer_cspc=None
1202 1203 self.__buffer_spc=numpy.array(self.__buffer_spc)
1203 1204 self.__buffer_cspc=numpy.array(self.__buffer_cspc)
1205
1204 1206 freq_dc = int(self.__buffer_spc.shape[2] / 2)
1205 1207 #print("FREQ_DC",freq_dc,self.__buffer_spc.shape,self.nHeights)
1206 1208 for k in range(7,self.nHeights):
1207 buffer_cspc=numpy.copy(self.__buffer_cspc[:,:,:,k])
1209 try:
1210 buffer_cspc=numpy.copy(self.__buffer_cspc[:,:,:,k])
1211 except Exception as e:
1212 print(e)
1208 1213 outliers_IDs_cspc=[]
1209 1214 cspc_outliers_exist=False
1210 1215 for i in range(self.nChannels):#dataOut.nChannels):
@@ -1292,7 +1297,7 class IntegrationFaradaySpectra(Operation):
1292 1297 self.__buffer_cspc = []
1293 1298 self.__buffer_dc = 0
1294 1299 self.__profIndex = 0
1295
1300 #print("pushData Done")
1296 1301 return data_spc, data_cspc, data_dc, n
1297 1302
1298 1303 def byProfiles(self, *args):
@@ -1374,10 +1379,14 class IntegrationFaradaySpectra(Operation):
1374 1379 if self.__dataReady:
1375 1380
1376 1381 if not self.ByLags:
1377
1378 dataOut.data_spc = numpy.squeeze(avgdata_spc)
1382 if self.nChannels == 1:
1383 dataOut.data_spc = avgdata_spc
1384 else:
1385 dataOut.data_spc = numpy.squeeze(avgdata_spc)
1379 1386 dataOut.data_cspc = numpy.squeeze(avgdata_cspc)
1380 1387 dataOut.data_dc = avgdata_dc
1388
1389
1381 1390 else:
1382 1391 dataOut.dataLag_spc = avgdata_spc
1383 1392 dataOut.dataLag_cspc = avgdata_cspc
General Comments 0
You need to be logged in to leave comments. Login now