From c5c71a9421125fb0885610abc8bab64e5009e4a3 2022-05-25 22:26:06 From: joabAM Date: 2022-05-25 22:26:06 Subject: [PATCH] nuevo --- diff --git a/schainc/_noise.c b/schainc/_noise.c index 32b4e30..97e07a0 100644 --- a/schainc/_noise.c +++ b/schainc/_noise.c @@ -49,7 +49,55 @@ static PyObject *hildebrand_sekhon(PyObject *self, PyObject *args) { // return PyLong_FromLong(lnoise); return PyFloat_FromDouble(lnoise); } +/* +static PyObject *hildebrand_sekhon2(PyObject *self, PyObject *args) { + double navg; + double th; + PyObject *data_obj, *data_array; + + if (!PyArg_ParseTuple(args, "Od", &data_obj, &navg, &th)) { + return NULL; + } + + data_array = PyArray_FROM_OTF(data_obj, NPY_FLOAT64, NPY_IN_ARRAY); + + if (data_array == NULL) { + Py_XDECREF(data_array); + Py_XDECREF(data_obj); + return NULL; + } + double *sortdata = (double*)PyArray_DATA(data_array); + int lenOfData = (int)PyArray_SIZE(data_array) ; + double nums_min = lenOfData*th; + if (nums_min <= 5) nums_min = 5; + double sump = 0; + double sumq = 0; + long j = 0; + int cont = 1; + double rtest = 0; + while ((cont == 1) && (j < lenOfData)) { + sump = sump + sortdata[j]; + sumq = sumq + pow(sortdata[j], 2); + if (j > nums_min) { + rtest = (double)j/(j-1) + 1/navg; + if ((sumq*j) > (rtest*pow(sump, 2))) { + j = j - 1; + sump = sump - sortdata[j]; + sumq = sumq - pow(sortdata[j],2); + cont = 0; + } + } + j = j + 1; + } + //double lnoise = sump / j; + + Py_DECREF(data_array); + + // return PyLong_FromLong(lnoise); + return PyFloat_FromDouble(j,sortID); +} +*/ static PyMethodDef noiseMethods[] = { { "hildebrand_sekhon", hildebrand_sekhon, METH_VARARGS, "Get noise with hildebrand_sekhon algorithm" }, diff --git a/schainpy/model/graphics/jroplot_parameters.py b/schainpy/model/graphics/jroplot_parameters.py index b4c4216..93e9230 100644 --- a/schainpy/model/graphics/jroplot_parameters.py +++ b/schainpy/model/graphics/jroplot_parameters.py @@ -50,10 +50,12 @@ class SnrPlot(RTIPlot): def update(self, dataOut): if len(self.channelList) == 0: self.update_list(dataOut) - data = {} + meta = {} - data['snr'] = 10*numpy.log10(dataOut.data_snr) - + data = { + 'snr': 10 * numpy.log10(dataOut.data_snr) + } + #print(data['snr']) return data, meta class DopplerPlot(RTIPlot): diff --git a/schainpy/model/graphics/jroplot_spectra.py b/schainpy/model/graphics/jroplot_spectra.py index 0da983d..16d851e 100644 --- a/schainpy/model/graphics/jroplot_spectra.py +++ b/schainpy/model/graphics/jroplot_spectra.py @@ -258,6 +258,7 @@ class RTIPlot(Plot): self.x = self.data.times self.y = self.data.yrange + self.z = self.data[self.CODE] self.z = numpy.array(self.z, dtype=float) self.z = numpy.ma.masked_invalid(self.z) @@ -274,11 +275,13 @@ class RTIPlot(Plot): x, y, z = self.fill_gaps(self.x, self.y, self.z) else: x, y, z = self.fill_gaps(*self.decimate()) - dummy_var = self.axes #Extrañamente esto actualiza el valor axes + + #dummy_var = self.axes #Extrañamente esto actualiza el valor axes for n, ax in enumerate(self.axes): self.zmin = self.zmin if self.zmin else numpy.min(self.z) self.zmax = self.zmax if self.zmax else numpy.max(self.z) data = self.data[-1] + if ax.firsttime: ax.plt = ax.pcolormesh(x, y, z[n].T, vmin=self.zmin, @@ -287,8 +290,8 @@ class RTIPlot(Plot): ) if self.showprofile: ax.plot_profile = self.pf_axes[n].plot(data[self.CODE][n], self.y)[0] - if "noise" in self.data: + ax.plot_noise = self.pf_axes[n].plot(numpy.repeat(data['noise'][n], len(self.y)), self.y, color="k", linestyle="dashed", lw=1)[0] else: @@ -301,6 +304,7 @@ class RTIPlot(Plot): if self.showprofile: ax.plot_profile.set_data(data[self.CODE][n], self.y) if "noise" in self.data: + ax.plot_noise.set_data(numpy.repeat( data['noise'][n], len(self.y)), self.y) diff --git a/schainpy/model/io/jroIO_param.py b/schainpy/model/io/jroIO_param.py index 0d2d818..d3db918 100644 --- a/schainpy/model/io/jroIO_param.py +++ b/schainpy/model/io/jroIO_param.py @@ -308,8 +308,8 @@ class HDFReader(Reader, ProcessingUnit): self.getData() -#@MPDecorator -class HDFWrite(Operation): +@MPDecorator +class HDFWriter(Operation): """Operation to write HDF5 files. The HDF5 file contains by default two groups Data and Metadata where @@ -465,7 +465,7 @@ class HDFWrite(Operation): def run(self, dataOut,**kwargs): - self.dataOut = dataOut + self.dataOut = dataOut.copy() if not(self.isConfig): self.setup(**kwargs) @@ -474,7 +474,7 @@ class HDFWrite(Operation): self.putData() - return self.dataOut + #return self.dataOut def setNextFile(self): diff --git a/schainpy/model/proc/jroproc_base.py b/schainpy/model/proc/jroproc_base.py index a4aa194..3a7ebba 100644 --- a/schainpy/model/proc/jroproc_base.py +++ b/schainpy/model/proc/jroproc_base.py @@ -89,7 +89,7 @@ class ProcessingUnit(object): elif optype == 'external' and self.dataOut.error: op.queue.put(copy.deepcopy(self.dataOut)) - return 'Error' if self.dataOut.error else True#self.dataOut.isReady() + return 'Error' if self.dataOut.error else self.dataOut.isReady() def setup(self): diff --git a/schainpy/model/proc/jroproc_spectra.py b/schainpy/model/proc/jroproc_spectra.py index 57dc8a9..03c38e0 100644 --- a/schainpy/model/proc/jroproc_spectra.py +++ b/schainpy/model/proc/jroproc_spectra.py @@ -202,7 +202,7 @@ class SpectraProc(ProcessingUnit): self.dataOut.flagNoData = False self.firstdatatime = None self.profIndex = 0 - self.dataOut.noise_estimation = None + else: raise ValueError("The type of input object '%s' is not valid".format( self.dataIn.type)) @@ -490,13 +490,13 @@ class removeDC(Operation): return self.dataOut class getNoise(Operation): + def __init__(self): Operation.__init__(self) - def run(self, dataOut, minHei=None, maxHei=None, minVel=None, maxVel=None, minFreq= None, maxFreq=None,): - self.dataOut = dataOut.copy() - print("1: ",dataOut.noise_estimation, dataOut.normFactor) + def run(self, dataOut, minHei=None, maxHei=None,minVel=None, maxVel=None, minFreq= None, maxFreq=None): + self.dataOut = dataOut if minHei == None: minHei = self.dataOut.heightList[0] @@ -606,9 +606,10 @@ class getNoise(Operation): maxIndexFFT = len( self.dataOut.getFreqRange(1)) #print(minIndex, maxIndex,minIndexVel, maxIndexVel) + self.dataOut.noise_estimation = None noise = self.dataOut.getNoise(xmin_index=minIndexFFT, xmax_index=maxIndexFFT, ymin_index=minIndex, ymax_index=maxIndex) - self.dataOut.noise_estimation = noise.copy() + self.dataOut.noise_estimation = noise.copy() # dataOut.noise #print("2: ",10*numpy.log10(self.dataOut.noise_estimation/64)) return self.dataOut @@ -1201,10 +1202,14 @@ class IntegrationFaradaySpectra(Operation): buffer_cspc=None self.__buffer_spc=numpy.array(self.__buffer_spc) self.__buffer_cspc=numpy.array(self.__buffer_cspc) + freq_dc = int(self.__buffer_spc.shape[2] / 2) #print("FREQ_DC",freq_dc,self.__buffer_spc.shape,self.nHeights) for k in range(7,self.nHeights): - buffer_cspc=numpy.copy(self.__buffer_cspc[:,:,:,k]) + try: + buffer_cspc=numpy.copy(self.__buffer_cspc[:,:,:,k]) + except Exception as e: + print(e) outliers_IDs_cspc=[] cspc_outliers_exist=False for i in range(self.nChannels):#dataOut.nChannels): @@ -1292,7 +1297,7 @@ class IntegrationFaradaySpectra(Operation): self.__buffer_cspc = [] self.__buffer_dc = 0 self.__profIndex = 0 - + #print("pushData Done") return data_spc, data_cspc, data_dc, n def byProfiles(self, *args): @@ -1374,10 +1379,14 @@ class IntegrationFaradaySpectra(Operation): if self.__dataReady: if not self.ByLags: - - dataOut.data_spc = numpy.squeeze(avgdata_spc) + if self.nChannels == 1: + dataOut.data_spc = avgdata_spc + else: + dataOut.data_spc = numpy.squeeze(avgdata_spc) dataOut.data_cspc = numpy.squeeze(avgdata_cspc) dataOut.data_dc = avgdata_dc + + else: dataOut.dataLag_spc = avgdata_spc dataOut.dataLag_cspc = avgdata_cspc