schain Commit - r1749:e3f1b2c53e4b · Jicamarca Repository

plot

joabAM -

r1749:e3f1b2c53e4b

parent child

Context file:

r1749:e3f1b2c53e4b

Collapse all files

schainpy/model/graphics/jroplot_spectra.py +1 -1

             # Copyright (c) 2012-2021 Jicamarca Radio Observatory
             # All rights reserved.
             #
             # Distributed under the terms of the BSD 3-clause license.
             """Classes to plot Spectra data
             """
             import os
             import numpy
             import datetime
             from schainpy.model.graphics.jroplot_base import Plot, plt, log
             from itertools import combinations
             from matplotlib.ticker import LinearLocator
             from schainpy.model.utils.BField import BField
             from scipy.interpolate import splrep
             from scipy.interpolate import splev
             from matplotlib import __version__ as plt_version
             if plt_version >='3.3.4':
                 EXTRA_POINTS = 0
             else:
                 EXTRA_POINTS = 1
             class SpectraPlot(Plot):
                 '''
                 Plot for Spectra data
                 '''
                 CODE = 'spc'
                 colormap = 'jet'
                 plot_type = 'pcolor'
                 buffering = False
                 channelList = []
                 elevationList = []
                 azimuthList = []
                 def setup(self):
                     self.nplots = len(self.data.channels)
                     self.ncols = int(numpy.sqrt(self.nplots) + 0.9)
                     self.nrows = int((1.0 * self.nplots / self.ncols) + 0.9)
                     self.height = 3.4 * self.nrows
                     self.cb_label = 'dB'
                     if self.showprofile:
                         self.width = 5.2 * self.ncols
                     else:
                         self.width = 4.2* self.ncols
                     self.plots_adjust.update({'wspace': 0.4, 'hspace':0.4, 'left': 0.1, 'right': 0.9, 'bottom': 0.12})
                     self.ylabel = 'Range [km]'
                 def update_list(self,dataOut):
                     if len(self.channelList) == 0:
                         self.channelList = dataOut.channelList
                     if len(self.elevationList) == 0:
                         self.elevationList = dataOut.elevationList
                     if len(self.azimuthList) == 0:
                         self.azimuthList = dataOut.azimuthList
                 def update(self, dataOut):
                     self.update_list(dataOut)
                     data = {}
                     meta = {}
                     norm = dataOut.nProfiles * dataOut.max_nIncohInt * dataOut.nCohInt  * dataOut.windowOfFilter
                     if dataOut.type == "Parameters":
                         noise = 10*numpy.log10(dataOut.getNoise()/dataOut.normFactor)
                         spc = 10*numpy.log10(dataOut.data_spc/(dataOut.nProfiles))
                     else:
                         noise = 10*numpy.log10(dataOut.getNoise()/norm)
                         z = numpy.zeros((dataOut.nChannels, dataOut.nFFTPoints, dataOut.nHeights))
                         for ch in range(dataOut.nChannels):
                             if hasattr(dataOut.normFactor,'ndim'):
                                 if dataOut.normFactor.ndim > 1:
                                     z[ch] = (numpy.divide(dataOut.data_spc[ch],dataOut.normFactor[ch]))
                                 else:
                                     z[ch] = (numpy.divide(dataOut.data_spc[ch],dataOut.normFactor))
                             else:
                                 z[ch] = (numpy.divide(dataOut.data_spc[ch],dataOut.normFactor))
                         z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
                         spc = 10*numpy.log10(z)
                     data['spc'] = spc
                     data['rti'] = spc.mean(axis=1)
                     data['noise'] = noise
                     meta['xrange'] = (dataOut.getFreqRange(EXTRA_POINTS)/1000., dataOut.getAcfRange(EXTRA_POINTS), dataOut.getVelRange(EXTRA_POINTS))
                     if self.CODE == 'spc_moments':
                         data['moments'] = dataOut.moments
                     return data, meta
                 def plot(self):
                     if self.xaxis == "frequency":
                         x = self.data.xrange[0]
                         self.xlabel = "Frequency (kHz)"
                     elif self.xaxis == "time":
                         x = self.data.xrange[1]
                         self.xlabel = "Time (ms)"
                     else:
                         x = self.data.xrange[2]
                         self.xlabel = "Velocity (m/s)"
                     if (self.CODE == 'spc_moments') | (self.CODE == 'gaussian_fit'):
                         x = self.data.xrange[2]
                         self.xlabel = "Velocity (m/s)"
                     self.titles = []
                     y = self.data.yrange
                     self.y = y
                     data = self.data[-1]
                     z = data['spc']
                     for n, ax in enumerate(self.axes):
                         noise = self.data['noise'][n][0]
                         # noise = data['noise'][n]
                         if self.CODE == 'spc_moments':
                             mean = data['moments'][n, 1]
                         if self.CODE == 'gaussian_fit':
                             gau0 = data['gaussfit'][n][2,:,0]
                             gau1 = data['gaussfit'][n][2,:,1]
                         if ax.firsttime:
                             self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
                             self.xmin = self.xmin if self.xmin else -self.xmax
                             self.zmin = self.zmin if self.zmin else numpy.nanmin(z)
                             self.zmax = self.zmax if self.zmax else numpy.nanmax(z)
                             ax.plt = ax.pcolormesh(x, y, z[n].T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
                             if self.showprofile:
                                 ax.plt_profile = self.pf_axes[n].plot(
                                     data['rti'][n], y)[0]
                                 ax.plt_noise = self.pf_axes[n].plot(numpy.repeat(noise, len(y)), y,
                                                                     color="k", linestyle="dashed", lw=1)[0]
                             if self.CODE == 'spc_moments':
                                 ax.plt_mean = ax.plot(mean, y, color='k', lw=1)[0]
                             if self.CODE == 'gaussian_fit':
                                 ax.plt_gau0 = ax.plot(gau0, y, color='r', lw=1)[0]
                                 ax.plt_gau1 = ax.plot(gau1, y, color='y', lw=1)[0]
                         else:
                             ax.plt.set_array(z[n].T.ravel())
                             if self.showprofile:
                                 ax.plt_profile.set_data(data['rti'][n], y)
                                 ax.plt_noise.set_data(numpy.repeat(noise, len(y)), y)
                             if self.CODE == 'spc_moments':
                                 ax.plt_mean.set_data(mean, y)
                             if self.CODE == 'gaussian_fit':
                                 ax.plt_gau0.set_data(gau0, y)
                                 ax.plt_gau1.set_data(gau1, y)
                         if len(self.azimuthList) > 0 and len(self.elevationList) > 0:
                             self.titles.append('CH {}: {:2.1f}elv {:2.1f}az {:3.2f}dB'.format(self.channelList[n], noise, self.elevationList[n], self.azimuthList[n]))
                         else:
                             self.titles.append('CH {}:  {:3.2f}dB'.format(self.channelList[n], noise))
             class SpectraObliquePlot(Plot):
                 '''
                 Plot for Spectra data
                 '''
                 CODE = 'spc_oblique'
                 colormap = 'jet'
                 plot_type = 'pcolor'
                 def setup(self):
                     self.xaxis = "oblique"
                     self.nplots = len(self.data.channels)
                     self.ncols = int(numpy.sqrt(self.nplots) + 0.9)
                     self.nrows = int((1.0 * self.nplots / self.ncols) + 0.9)
                     self.height = 2.6 * self.nrows
                     self.cb_label = 'dB'
                     if self.showprofile:
                         self.width = 4 * self.ncols
                     else:
                         self.width = 3.5 * self.ncols
                     self.plots_adjust.update({'wspace': 0.8, 'hspace':0.2, 'left': 0.2, 'right': 0.9, 'bottom': 0.18})
                     self.ylabel = 'Range [km]'
                 def update(self, dataOut):
                     data = {}
                     meta = {}
                     spc = 10*numpy.log10(dataOut.data_spc/dataOut.normFactor)
                     data['spc'] = spc
                     data['rti'] = dataOut.getPower()
                     data['noise'] = 10*numpy.log10(dataOut.getNoise()/dataOut.normFactor)
                     meta['xrange'] = (dataOut.getFreqRange(1)/1000., dataOut.getAcfRange(1), dataOut.getVelRange(1))
                     data['shift1'] = dataOut.Dop_EEJ_T1[0]
                     data['shift2'] = dataOut.Dop_EEJ_T2[0]
                     data['max_val_2'] = dataOut.Oblique_params[0,-1,:]
                     data['shift1_error'] = dataOut.Err_Dop_EEJ_T1[0]
                     data['shift2_error'] = dataOut.Err_Dop_EEJ_T2[0]
                     return data, meta
                 def plot(self):
                     if self.xaxis == "frequency":
                         x = self.data.xrange[0]
                         self.xlabel = "Frequency (kHz)"
                     elif self.xaxis == "time":
                         x = self.data.xrange[1]
                         self.xlabel = "Time (ms)"
                     else:
                         x = self.data.xrange[2]
                         self.xlabel = "Velocity (m/s)"
                     self.titles = []
                     y = self.data.yrange
                     self.y = y
                     data = self.data[-1]
                     z = data['spc']
                     for n, ax in enumerate(self.axes):
                         noise = self.data['noise'][n][-1]
                         shift1 = data['shift1']
                         shift2 = data['shift2']
                         max_val_2 = data['max_val_2']
                         err1 = data['shift1_error']
                         err2 = data['shift2_error']
                         if ax.firsttime:
                             self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
                             self.xmin = self.xmin if self.xmin else -self.xmax
                             self.zmin = self.zmin if self.zmin else numpy.nanmin(z)
                             self.zmax = self.zmax if self.zmax else numpy.nanmax(z)
                             ax.plt = ax.pcolormesh(x, y, z[n].T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
                             if self.showprofile:
                                 ax.plt_profile = self.pf_axes[n].plot(
                                     self.data['rti'][n][-1], y)[0]
                                 ax.plt_noise = self.pf_axes[n].plot(numpy.repeat(noise, len(y)), y,
                                                                     color="k", linestyle="dashed", lw=1)[0]
                             self.ploterr1 = ax.errorbar(shift1, y, xerr=err1, fmt='k^', elinewidth=2.2, marker='o', linestyle='None',markersize=2.5,capsize=0.3,markeredgewidth=0.2)
                             self.ploterr2 = ax.errorbar(shift2, y, xerr=err2, fmt='m^',elinewidth=2.2,marker='o',linestyle='None',markersize=2.5,capsize=0.3,markeredgewidth=0.2)
                             self.ploterr3 = ax.errorbar(max_val_2, y, xerr=0, fmt='g^',elinewidth=2.2,marker='o',linestyle='None',markersize=2.5,capsize=0.3,markeredgewidth=0.2)
                         else:
                             self.ploterr1.remove()
                             self.ploterr2.remove()
                             self.ploterr3.remove()
                             ax.plt.set_array(z[n].T.ravel())
                             if self.showprofile:
                                 ax.plt_profile.set_data(self.data['rti'][n][-1], y)
                                 ax.plt_noise.set_data(numpy.repeat(noise, len(y)), y)
                             self.ploterr1 = ax.errorbar(shift1, y, xerr=err1, fmt='k^', elinewidth=2.2, marker='o', linestyle='None',markersize=2.5,capsize=0.3,markeredgewidth=0.2)
                             self.ploterr2 = ax.errorbar(shift2, y, xerr=err2, fmt='m^',elinewidth=2.2,marker='o',linestyle='None',markersize=2.5,capsize=0.3,markeredgewidth=0.2)
                             self.ploterr3 = ax.errorbar(max_val_2, y, xerr=0, fmt='g^',elinewidth=2.2,marker='o',linestyle='None',markersize=2.5,capsize=0.3,markeredgewidth=0.2)
                         self.titles.append('CH {}: {:3.2f}dB'.format(n, noise))
             class CrossSpectraPlot(Plot):
                 CODE = 'cspc'
                 colormap = 'jet'
                 plot_type = 'pcolor'
                 zmin_coh = None
                 zmax_coh = None
                 zmin_phase = None
                 zmax_phase = None
                 realChannels = None
                 crossPairs = None
                 def setup(self):
                     self.ncols = 4
                     self.nplots = len(self.data.pairs) * 2
                     self.nrows = int((1.0 * self.nplots / self.ncols) + 0.9)
                     self.width = 3.1 * self.ncols
                     self.height = 2.6 * self.nrows
                     self.ylabel = 'Range [km]'
                     self.showprofile = False
                     self.plots_adjust.update({'left': 0.08, 'right': 0.92, 'wspace': 0.5, 'hspace':0.4, 'top':0.95, 'bottom': 0.08})
                 def update(self, dataOut):
                     data = {}
                     meta = {}
                     spc = dataOut.data_spc
                     cspc = dataOut.data_cspc
                     meta['xrange'] = (dataOut.getFreqRange(EXTRA_POINTS)/1000., dataOut.getAcfRange(EXTRA_POINTS), dataOut.getVelRange(EXTRA_POINTS))
                     rawPairs = list(combinations(list(range(dataOut.nChannels)), 2))
                     meta['pairs'] = rawPairs
                     if self.crossPairs == None:
                         self.crossPairs = dataOut.pairsList
                     tmp = []
                     for n, pair in enumerate(meta['pairs']):
                         out = cspc[n] / numpy.sqrt(spc[pair[0]] * spc[pair[1]])
                         coh = numpy.abs(out)
                         phase = numpy.arctan2(out.imag, out.real) * 180 / numpy.pi
                         tmp.append(coh)
                         tmp.append(phase)
                     data['cspc'] = numpy.array(tmp)
                     return data, meta
                 def plot(self):
                     if self.xaxis == "frequency":
                         x = self.data.xrange[0]
                         self.xlabel = "Frequency (kHz)"
                     elif self.xaxis == "time":
                         x = self.data.xrange[1]
                         self.xlabel = "Time (ms)"
                     else:
                         x = self.data.xrange[2]
                         self.xlabel = "Velocity (m/s)"
                     self.titles = []
                     y = self.data.yrange
                     self.y = y
                     data = self.data[-1]
                     cspc = data['cspc']
                     for n in range(len(self.data.pairs)):
                         pair = self.crossPairs[n]
                         coh = cspc[n*2]
                         phase = cspc[n*2+1]
                         ax = self.axes[2 * n]
                         if ax.firsttime:
                             ax.plt = ax.pcolormesh(x, y, coh.T,
                                                    vmin=self.zmin_coh,
                                                    vmax=self.zmax_coh,
                                                    cmap=plt.get_cmap(self.colormap_coh)
                                                    )
                         else:
                             ax.plt.set_array(coh.T.ravel())
                         self.titles.append(
                             'Coherence Ch{} * Ch{}'.format(pair[0], pair[1]))
                         ax = self.axes[2 * n + 1]
                         if ax.firsttime:
                             ax.plt = ax.pcolormesh(x, y, phase.T,
                                                    vmin=-180,
                                                    vmax=180,
                                                    cmap=plt.get_cmap(self.colormap_phase)
                                                    )
                         else:
                             ax.plt.set_array(phase.T.ravel())
                         self.titles.append('Phase CH{} * CH{}'.format(pair[0], pair[1]))
             class CrossSpectra4Plot(Plot):
                 CODE = 'cspc'
                 colormap = 'jet'
                 plot_type = 'pcolor'
                 zmin_coh = None
                 zmax_coh = None
                 zmin_phase = None
                 zmax_phase = None
                 def setup(self):
                     self.ncols = 4
                     self.nrows = len(self.data.pairs)
                     self.nplots = self.nrows * 4
                     self.width = 3.1 * self.ncols
                     self.height = 5 * self.nrows
                     self.ylabel = 'Range [km]'
                     self.showprofile = False
                     self.plots_adjust.update({'left': 0.08, 'right': 0.92, 'wspace': 0.5, 'hspace':0.4, 'top':0.95, 'bottom': 0.08})
                 def plot(self):
                     if self.xaxis == "frequency":
                         x = self.data.xrange[0]
                         self.xlabel = "Frequency (kHz)"
                     elif self.xaxis == "time":
                         x = self.data.xrange[1]
                         self.xlabel = "Time (ms)"
                     else:
                         x = self.data.xrange[2]
                         self.xlabel = "Velocity (m/s)"
                     self.titles = []
                     y = self.data.heights
                     self.y = y
                     nspc = self.data['spc']
                     spc = self.data['cspc'][0]
                     cspc = self.data['cspc'][1]
                     for n in range(self.nrows):
                         noise = self.data['noise'][:,-1]
                         pair = self.data.pairs[n]
                         ax = self.axes[4 * n]
                         if ax.firsttime:
                             self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
                             self.xmin = self.xmin if self.xmin else -self.xmax
                             self.zmin = self.zmin if self.zmin else numpy.nanmin(nspc)
                             self.zmax = self.zmax if self.zmax else numpy.nanmax(nspc)
                             ax.plt = ax.pcolormesh(x , y , nspc[pair[0]].T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
                         else:
                             ax.plt.set_array(nspc[pair[0]].T.ravel())
                         self.titles.append('CH {}: {:3.2f}dB'.format(pair[0], noise[pair[0]]))
                         ax = self.axes[4 * n + 1]
                         if ax.firsttime:
                             ax.plt = ax.pcolormesh(x , y, numpy.flip(nspc[pair[1]],axis=0).T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
                         else:
                             ax.plt.set_array(numpy.flip(nspc[pair[1]],axis=0).T.ravel())
                         self.titles.append('CH {}: {:3.2f}dB'.format(pair[1], noise[pair[1]]))
                         out = cspc[n] / numpy.sqrt(spc[pair[0]] * spc[pair[1]])
                         coh = numpy.abs(out)
                         phase = numpy.arctan2(out.imag, out.real) * 180 / numpy.pi
                         ax = self.axes[4 * n + 2]
                         if ax.firsttime:
                             ax.plt = ax.pcolormesh(x, y, numpy.flip(coh,axis=0).T,
                                                    vmin=0,
                                                    vmax=1,
                                                    cmap=plt.get_cmap(self.colormap_coh)
                                                    )
                         else:
                             ax.plt.set_array(numpy.flip(coh,axis=0).T.ravel())
                         self.titles.append(
                             'Coherence Ch{} * Ch{}'.format(pair[0], pair[1]))
                         ax = self.axes[4 * n + 3]
                         if ax.firsttime:
                             ax.plt = ax.pcolormesh(x, y, numpy.flip(phase,axis=0).T,
                                                    vmin=-180,
                                                    vmax=180,
                                                    cmap=plt.get_cmap(self.colormap_phase)
                                                    )
                         else:
                             ax.plt.set_array(numpy.flip(phase,axis=0).T.ravel())
                         self.titles.append('Phase CH{} * CH{}'.format(pair[0], pair[1]))
             class CrossSpectra2Plot(Plot):
                 CODE = 'cspc'
                 colormap = 'jet'
                 plot_type = 'pcolor'
                 zmin_coh = None
                 zmax_coh = None
                 zmin_phase = None
                 zmax_phase = None
                 def setup(self):
                     self.ncols = 1
                     self.nrows = len(self.data.pairs)
                     self.nplots = self.nrows * 1
                     self.width = 3.1 * self.ncols
                     self.height = 5 * self.nrows
                     self.ylabel = 'Range [km]'
                     self.showprofile = False
                     self.plots_adjust.update({'left': 0.22, 'right': .90, 'wspace': 0.5, 'hspace':0.4, 'top':0.95, 'bottom': 0.08})
                 def plot(self):
                     if self.xaxis == "frequency":
                         x = self.data.xrange[0]
                         self.xlabel = "Frequency (kHz)"
                     elif self.xaxis == "time":
                         x = self.data.xrange[1]
                         self.xlabel = "Time (ms)"
                     else:
                         x = self.data.xrange[2]
                         self.xlabel = "Velocity (m/s)"
                     self.titles = []
                     y = self.data.heights
                     self.y = y
                     cspc = self.data['cspc'][1]
                     for n in range(self.nrows):
                         noise = self.data['noise'][:,-1]
                         pair = self.data.pairs[n]
                         out = cspc[n]
                         cross = numpy.abs(out)
                         z = cross/self.data.nFactor
                         cross = 10*numpy.log10(z)
                         ax = self.axes[1 * n]
                         if ax.firsttime:
                             self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
                             self.xmin = self.xmin if self.xmin else -self.xmax
                             self.zmin = self.zmin if self.zmin else numpy.nanmin(cross)
                             self.zmax = self.zmax if self.zmax else numpy.nanmax(cross)
                             ax.plt = ax.pcolormesh(x, y, cross.T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
                         else:
                             ax.plt.set_array(cross.T.ravel())
                         self.titles.append(
                             'Cross Spectra Power Ch{} * Ch{}'.format(pair[0], pair[1]))
             class CrossSpectra3Plot(Plot):
                 CODE = 'cspc'
                 colormap = 'jet'
                 plot_type = 'pcolor'
                 zmin_coh = None
                 zmax_coh = None
                 zmin_phase = None
                 zmax_phase = None
                 def setup(self):
                     self.ncols = 3
                     self.nrows = len(self.data.pairs)
                     self.nplots = self.nrows * 3
                     self.width = 3.1 * self.ncols
                     self.height = 5 * self.nrows
                     self.ylabel = 'Range [km]'
                     self.showprofile = False
                     self.plots_adjust.update({'left': 0.22, 'right': .90, 'wspace': 0.5, 'hspace':0.4, 'top':0.95, 'bottom': 0.08})
                 def plot(self):
                     if self.xaxis == "frequency":
                         x = self.data.xrange[0]
                         self.xlabel = "Frequency (kHz)"
                     elif self.xaxis == "time":
                         x = self.data.xrange[1]
                         self.xlabel = "Time (ms)"
                     else:
                         x = self.data.xrange[2]
                         self.xlabel = "Velocity (m/s)"
                     self.titles = []
                     y = self.data.heights
                     self.y = y
                     cspc = self.data['cspc'][1]
                     for n in range(self.nrows):
                         noise = self.data['noise'][:,-1]
                         pair = self.data.pairs[n]
                         out = cspc[n]
                         cross = numpy.abs(out)
                         z = cross/self.data.nFactor
                         cross = 10*numpy.log10(z)
                         out_r= out.real/self.data.nFactor
                         out_i= out.imag/self.data.nFactor
                         ax = self.axes[3 * n]
                         if ax.firsttime:
                             self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
                             self.xmin = self.xmin if self.xmin else -self.xmax
                             self.zmin = self.zmin if self.zmin else numpy.nanmin(cross)
                             self.zmax = self.zmax if self.zmax else numpy.nanmax(cross)
                             ax.plt = ax.pcolormesh(x, y, cross.T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
                         else:
                             ax.plt.set_array(cross.T.ravel())
                         self.titles.append(
                             'Cross Spectra Power Ch{} * Ch{}'.format(pair[0], pair[1]))
                         ax = self.axes[3 * n + 1]
                         if ax.firsttime:
                             self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
                             self.xmin = self.xmin if self.xmin else -self.xmax
                             self.zmin = self.zmin if self.zmin else numpy.nanmin(cross)
                             self.zmax = self.zmax if self.zmax else numpy.nanmax(cross)
                             ax.plt = ax.pcolormesh(x, y, out_r.T,
                                                    vmin=-1.e6,
                                                    vmax=0,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
                         else:
                             ax.plt.set_array(out_r.T.ravel())
                         self.titles.append(
                             'Cross Spectra Real Ch{} * Ch{}'.format(pair[0], pair[1]))
                         ax = self.axes[3 * n + 2]
                         if ax.firsttime:
                             self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
                             self.xmin = self.xmin if self.xmin else -self.xmax
                             self.zmin = self.zmin if self.zmin else numpy.nanmin(cross)
                             self.zmax = self.zmax if self.zmax else numpy.nanmax(cross)
                             ax.plt = ax.pcolormesh(x, y, out_i.T,
                                                    vmin=-1.e6,
                                                    vmax=1.e6,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
                         else:
                             ax.plt.set_array(out_i.T.ravel())
                         self.titles.append(
                             'Cross Spectra Imag Ch{} * Ch{}'.format(pair[0], pair[1]))
             class RTIPlot(Plot):
                 '''
                 Plot for RTI data
                 '''
                 CODE = 'rti'
                 colormap = 'jet'
                 plot_type = 'pcolorbuffer'
                 titles = None
                 channelList = []
                 elevationList = []
                 azimuthList = []
                 def setup(self):
                     self.xaxis = 'time'
                     self.ncols = 1
                     self.nrows = len(self.data.channels)
                     self.nplots = len(self.data.channels)
                     self.ylabel = 'Range [km]'
                     #self.xlabel = 'Time'
                     self.cb_label = 'dB'
                     self.plots_adjust.update({'hspace':0.8, 'left': 0.1, 'bottom': 0.1, 'right':0.95})
                     self.titles = ['{} Channel {}'.format(
                         self.CODE.upper(), x) for x in range(self.nplots)]
                 def update_list(self,dataOut):
                     if len(self.channelList) == 0:
                         self.channelList = dataOut.channelList
                     if len(self.elevationList) == 0:
                         self.elevationList = dataOut.elevationList
                     if len(self.azimuthList) == 0:
                         self.azimuthList = dataOut.azimuthList
                 def update(self, dataOut):
                     if len(self.channelList) == 0:
                         self.update_list(dataOut)
                     data = {}
                     meta = {}
                     data['rti'] = dataOut.getPower()
                     norm = dataOut.nProfiles * dataOut.max_nIncohInt * dataOut.nCohInt  * dataOut.windowOfFilter
                     noise = 10*numpy.log10(dataOut.getNoise()/norm)
                     data['noise'] = noise
                     return data, meta
                 def plot(self):
                     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)
                     try:
                         if self.channelList != None:
                             if len(self.elevationList) > 0 and len(self.azimuthList) > 0:
                                 self.titles = ['{}  Channel {} ({:2.1f} Elev,  {:2.1f} Azth)'.format(
                                     self.CODE.upper(), x, self.elevationList[x], self.azimuthList[x]) for x in self.channelList]
                             else:
                                 self.titles = ['{}  Channel {}'.format(
                                     self.CODE.upper(), x) for x in self.channelList]
                     except:
                         if self.channelList.any() != None:
                             if len(self.elevationList) > 0 and len(self.azimuthList) > 0:
                                 self.titles = ['{}  Channel {} ({:2.1f} Elev,  {:2.1f} Azth)'.format(
                                     self.CODE.upper(), x, self.elevationList[x], self.azimuthList[x]) for x in self.channelList]
                             else:
                                 self.titles = ['{} Channel {}'.format(
                                     self.CODE.upper(), x) for x in self.channelList]
                     if self.decimation is None:
                         x, y, z = self.fill_gaps(self.x, self.y, self.z)
                     else:
                         x, y, z = self.fill_gaps(*self.decimate())
                     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,
                                                    vmax=self.zmax,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
                             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:
-                            # ax.collections.remove(ax.collections[0]) #  error while running
+                            ax.collections.remove(ax.collections[0]) #  error while running in 3.12
                             ax.plt = ax.pcolormesh(x, y, z[n].T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
                             if self.showprofile:
                                 ax.plot_profile.set_data(data[self.CODE][n], self.y)
                                 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]
             class SpectrogramPlot(Plot):
                 '''
                 Plot for Spectrogram data
                 '''
                 CODE = 'Spectrogram_Profile'
                 colormap = 'binary'
                 plot_type = 'pcolorbuffer'
                 def setup(self):
                     self.xaxis = 'time'
                     self.ncols = 1
                     self.nrows = len(self.data.channels)
                     self.nplots = len(self.data.channels)
                     self.xlabel = 'Time'
                     self.plots_adjust.update({'hspace':1.2, 'left': 0.1, 'bottom': 0.12, 'right':0.95})
                     self.titles = []
                     self.titles = ['{} Channel {}'.format(
                         self.CODE.upper(), x) for x in range(self.nrows)]
                 def update(self, dataOut):
                     data = {}
                     meta = {}
                     maxHei = 1620#+12000
                     indb = numpy.where(dataOut.heightList <= maxHei)
                     hei = indb[0][-1]
                     factor = dataOut.nIncohInt
                     z = dataOut.data_spc[:,:,hei] / factor
                     z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
                     meta['xrange'] = (dataOut.getFreqRange(1)/1000., dataOut.getAcfRange(1), dataOut.getVelRange(1))
                     data['Spectrogram_Profile'] = 10 * numpy.log10(z)
                     data['hei'] = hei
                     data['DH'] = (dataOut.heightList[1] - dataOut.heightList[0])/dataOut.step
                     data['nProfiles'] = dataOut.nProfiles
                     return data, meta
                 def plot(self):
                     self.x = self.data.times
                     self.z = self.data[self.CODE]
                     self.y = self.data.xrange[0]
                     hei = self.data['hei'][-1]
                     DH = self.data['DH'][-1]
                     nProfiles = self.data['nProfiles'][-1]
                     self.ylabel = "Frequency (kHz)"
                     self.z = numpy.ma.masked_invalid(self.z)
                     if self.decimation is None:
                         x, y, z = self.fill_gaps(self.x, self.y, self.z)
                     else:
                         x, y, z = self.fill_gaps(*self.decimate())
                     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,
                                                    vmax=self.zmax,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
                         else:
                             # ax.collections.remove(ax.collections[0]) #  error while running
                             ax.plt = ax.pcolormesh(x, y, z[n].T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
             class CoherencePlot(RTIPlot):
                 '''
                 Plot for Coherence data
                 '''
                 CODE = 'coh'
                 titles = None
                 def setup(self):
                     self.xaxis = 'time'
                     self.ncols = 1
                     self.nrows = len(self.data.pairs)
                     self.nplots = len(self.data.pairs)
                     self.ylabel = 'Range [km]'
                     self.xlabel = 'Time'
                     self.plots_adjust.update({'hspace':0.6, 'left': 0.1, 'bottom': 0.1,'right':0.95})
                     if self.CODE == 'coh':
                         self.cb_label = ''
                         self.titles = [
                             'Coherence Map Ch{} * Ch{}'.format(x[0], x[1]) for x in self.data.pairs]
                     else:
                         self.cb_label = 'Degrees'
                         self.titles = [
                             'Phase Map Ch{} * Ch{}'.format(x[0], x[1]) for x in self.data.pairs]
                 def update(self, dataOut):
                     data = {}
                     meta = {}
                     data['coh'] = dataOut.getCoherence()
                     meta['pairs'] = dataOut.pairsList
                     return data, meta
             class PhasePlot(CoherencePlot):
                 '''
                 Plot for Phase map data
                 '''
                 CODE = 'phase'
                 colormap = 'seismic'
                 def update(self, dataOut):
                     data = {}
                     meta = {}
                     data['phase'] = dataOut.getCoherence(phase=True)
                     meta['pairs'] = dataOut.pairsList
                     return data, meta
             class NoisePlot(Plot):
                 '''
                 Plot for noise
                 '''
                 CODE = 'noise'
                 plot_type = 'scatterbuffer'
                 def setup(self):
                     self.xaxis = 'time'
                     self.ncols = 1
                     self.nrows = 1
                     self.nplots = 1
                     self.ylabel = 'Intensity [dB]'
                     self.xlabel = 'Time'
                     self.titles = ['Noise']
                     self.colorbar = False
                     self.plots_adjust.update({'right': 0.85 })
                     self.titles = ['Noise Plot']
                 def update(self, dataOut):
                     data = {}
                     meta = {}
                     noise =  10*numpy.log10(dataOut.getNoise())
                     noise = noise.reshape(dataOut.nChannels, 1)
                     data['noise'] = noise
                     meta['yrange'] = numpy.array([])
                     return data, meta
                 def plot(self):
                     x = self.data.times
                     xmin = self.data.min_time
                     xmax = xmin + self.xrange * 60 * 60
                     Y = self.data['noise']
                     if self.axes[0].firsttime:
                         self.ymin = numpy.nanmin(Y) - 5
                         self.ymax = numpy.nanmax(Y) + 5
                         for ch in self.data.channels:
                             y = Y[ch]
                             self.axes[0].plot(x, y, lw=1, label='Ch{}'.format(ch))
                         plt.legend(bbox_to_anchor=(1.18, 1.0))
                     else:
                         for ch in self.data.channels:
                             y = Y[ch]
                             self.axes[0].lines[ch].set_data(x, y)
             class PowerProfilePlot(Plot):
                 CODE = 'pow_profile'
                 plot_type = 'scatter'
                 def setup(self):
                     self.ncols = 1
                     self.nrows = 1
                     self.nplots = 1
                     self.height = 4
                     self.width = 3
                     self.ylabel = 'Range [km]'
                     self.xlabel = 'Intensity [dB]'
                     self.titles = ['Power Profile']
                     self.colorbar = False
                 def update(self, dataOut):
                     data = {}
                     meta = {}
                     data[self.CODE] = dataOut.getPower()
                     return data, meta
                 def plot(self):
                     y = self.data.yrange
                     self.y = y
                     x = self.data[-1][self.CODE]
                     if self.xmin is None: self.xmin = numpy.nanmin(x)*0.9
                     if self.xmax is None: self.xmax = numpy.nanmax(x)*1.1
                     if self.axes[0].firsttime:
                         for ch in self.data.channels:
                             self.axes[0].plot(x[ch], y, lw=1, label='Ch{}'.format(ch))
                         plt.legend()
                     else:
                         for ch in self.data.channels:
                             self.axes[0].lines[ch].set_data(x[ch], y)
             class SpectraCutPlot(Plot):
                 CODE = 'spc_cut'
                 plot_type = 'scatter'
                 buffering = False
                 heights = []
                 channelList = []
                 maintitle = "Spectra Cuts"
                 flag_setIndex = False
                 def setup(self):
                     self.nplots = len(self.data.channels)
                     self.ncols = int(numpy.sqrt(self.nplots) + 0.9)
                     self.nrows = int((1.0 * self.nplots / self.ncols) + 0.9)
                     self.width = 4.5 * self.ncols + 2.5
                     self.height = 4.8 * self.nrows
                     self.ylabel = 'Power [dB]'
                     self.colorbar = False
                     self.plots_adjust.update({'left':0.1, 'hspace':0.3, 'right': 0.9, 'bottom':0.08})
                     if len(self.selectedHeightsList) > 0:
                         self.maintitle = "Spectra Cut"# for %d km " %(int(self.selectedHeight))
                 def update(self, dataOut):
                     if len(self.channelList) == 0:
                         self.channelList = dataOut.channelList
                     self.heights = dataOut.heightList
                     #print("sels: ",self.selectedHeightsList)
                     if len(self.selectedHeightsList)>0 and not self.flag_setIndex:
                         for sel_height in self.selectedHeightsList:
                             index_list = numpy.where(self.heights >= sel_height)
                             index_list = index_list[0]
                             self.height_index.append(index_list[0])
                         #print("sels i:"", self.height_index)
                         self.flag_setIndex = True
                         #print(self.height_index)
                     data = {}
                     meta = {}
                     norm = dataOut.nProfiles * dataOut.max_nIncohInt * dataOut.nCohInt  * dataOut.windowOfFilter#*dataOut.nFFTPoints
                     n0 = 10*numpy.log10(dataOut.getNoise()/norm)
                     noise = numpy.repeat(n0,(dataOut.nFFTPoints*dataOut.nHeights)).reshape(dataOut.nChannels,dataOut.nFFTPoints,dataOut.nHeights)
                     z = []
                     for ch in range(dataOut.nChannels):
                         if hasattr(dataOut.normFactor,'shape'):
                             z.append(numpy.divide(dataOut.data_spc[ch],dataOut.normFactor[ch]))
                         else:
                             z.append(numpy.divide(dataOut.data_spc[ch],dataOut.normFactor))
                     z = numpy.asarray(z)
                     z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
                     spc = 10*numpy.log10(z)
                     data['spc'] = spc - noise
                     meta['xrange'] = (dataOut.getFreqRange(EXTRA_POINTS)/1000., dataOut.getAcfRange(EXTRA_POINTS), dataOut.getVelRange(EXTRA_POINTS))
                     return data, meta
                 def plot(self):
                     if self.xaxis == "frequency":
                         x = self.data.xrange[0][0:]
                         self.xlabel = "Frequency (kHz)"
                     elif self.xaxis == "time":
                         x = self.data.xrange[1]
                         self.xlabel = "Time (ms)"
                     else:
                         x = self.data.xrange[2]
                         self.xlabel = "Velocity (m/s)"
                     self.titles = []
                     y = self.data.yrange
                     z = self.data[-1]['spc']
                     #print(z.shape)
                     if len(self.height_index) > 0:
                         index = self.height_index
                     else:
                         index = numpy.arange(0, len(y), int((len(y))/9))
                     #print("inde x ", index, self.axes)
                     for n, ax in enumerate(self.axes):
                         if ax.firsttime:
                             self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
                             self.xmin = self.xmin if self.xmin else -self.xmax
                             self.ymin = self.ymin if self.ymin else numpy.nanmin(z)
                             self.ymax = self.ymax if self.ymax else numpy.nanmax(z)
                             ax.plt = ax.plot(x, z[n, :, index].T)
                             labels = ['Range = {:2.1f}km'.format(y[i]) for i in index]
                             self.figures[0].legend(ax.plt, labels, loc='center right', prop={'size': 8})
                             ax.minorticks_on()
                             ax.grid(which='major', axis='both')
                             ax.grid(which='minor', axis='x')
                         else:
                             for i, line in enumerate(ax.plt):
                                 line.set_data(x, z[n, :, index[i]])
                         self.titles.append('CH {}'.format(self.channelList[n]))
                     plt.suptitle(self.maintitle,  fontsize=10)
             class BeaconPhase(Plot):
                 __isConfig = None
                 __nsubplots = None
                 PREFIX = 'beacon_phase'
                 def __init__(self):
                     Plot.__init__(self)
                     self.timerange = 24*60*60
                     self.isConfig = False
                     self.__nsubplots = 1
                     self.counter_imagwr = 0
                     self.WIDTH = 800
                     self.HEIGHT = 400
                     self.WIDTHPROF = 120
                     self.HEIGHTPROF = 0
                     self.xdata = None
                     self.ydata = None
                     self.PLOT_CODE = BEACON_CODE
                     self.FTP_WEI = None
                     self.EXP_CODE = None
                     self.SUB_EXP_CODE = None
                     self.PLOT_POS = None
                     self.filename_phase = None
                     self.figfile = None
                     self.xmin = None
                     self.xmax = None
                 def getSubplots(self):
                     ncol = 1
                     nrow = 1
                     return nrow, ncol
                 def setup(self, id, nplots, wintitle, showprofile=True, show=True):
                     self.__showprofile = showprofile
                     self.nplots = nplots
                     ncolspan = 7
                     colspan = 6
                     self.__nsubplots = 2
                     self.createFigure(id = id,
                                       wintitle = wintitle,
                                       widthplot = self.WIDTH+self.WIDTHPROF,
                                       heightplot = self.HEIGHT+self.HEIGHTPROF,
                                       show=show)
                     nrow, ncol = self.getSubplots()
                     self.addAxes(nrow, ncol*ncolspan, 0, 0, colspan, 1)
                 def save_phase(self, filename_phase):
                     f = open(filename_phase,'w+')
                     f.write('\n\n')
                     f.write('JICAMARCA RADIO OBSERVATORY - Beacon Phase \n')
                     f.write('DD MM YYYY  HH MM SS   pair(2,0) pair(2,1) pair(2,3) pair(2,4)\n\n' )
                     f.close()
                 def save_data(self, filename_phase, data, data_datetime):
                     f=open(filename_phase,'a')
                     timetuple_data = data_datetime.timetuple()
                     day = str(timetuple_data.tm_mday)
                     month = str(timetuple_data.tm_mon)
                     year = str(timetuple_data.tm_year)
                     hour = str(timetuple_data.tm_hour)
                     minute = str(timetuple_data.tm_min)
                     second = str(timetuple_data.tm_sec)
                     f.write(day+' '+month+' '+year+'  '+hour+' '+minute+' '+second+'   '+str(data[0])+'   '+str(data[1])+'   '+str(data[2])+'   '+str(data[3])+'\n')
                     f.close()
                 def plot(self):
                     log.warning('TODO: Not yet implemented...')
                 def run(self, dataOut, id, wintitle="", pairsList=None, showprofile='True',
                         xmin=None, xmax=None, ymin=None, ymax=None, hmin=None, hmax=None,
                         timerange=None,
                         save=False, figpath='./', figfile=None, show=True, ftp=False, wr_period=1,
                         server=None, folder=None, username=None, password=None,
                         ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0):
                     if dataOut.flagNoData:
                         return dataOut
                     if not isTimeInHourRange(dataOut.datatime, xmin, xmax):
                         return
                     if pairsList == None:
                         pairsIndexList = dataOut.pairsIndexList[:10]
                     else:
                         pairsIndexList = []
                         for pair in pairsList:
                             if pair not in dataOut.pairsList:
                                 raise ValueError("Pair %s is not in dataOut.pairsList" %(pair))
                             pairsIndexList.append(dataOut.pairsList.index(pair))
                     if pairsIndexList == []:
                         return
              #         if len(pairsIndexList) > 4:
              #             pairsIndexList = pairsIndexList[0:4]
                     hmin_index = None
                     hmax_index = None
                     if hmin != None and hmax != None:
                         indexes = numpy.arange(dataOut.nHeights)
                         hmin_list = indexes[dataOut.heightList >= hmin]
                         hmax_list = indexes[dataOut.heightList <= hmax]
                         if hmin_list.any():
                             hmin_index = hmin_list[0]
                         if hmax_list.any():
                             hmax_index = hmax_list[-1]+1
                     x = dataOut.getTimeRange()
                     thisDatetime = dataOut.datatime
                     title = wintitle + " Signal Phase" # : %s" %(thisDatetime.strftime("%d-%b-%Y"))
                     xlabel = "Local Time"
                     ylabel = "Phase (degrees)"
                     update_figfile = False
                     nplots = len(pairsIndexList)
                     phase_beacon = numpy.zeros(len(pairsIndexList))
                     for i in range(nplots):
                         pair = dataOut.pairsList[pairsIndexList[i]]
                         ccf = numpy.average(dataOut.data_cspc[pairsIndexList[i], :, hmin_index:hmax_index], axis=0)
                         powa = numpy.average(dataOut.data_spc[pair[0], :, hmin_index:hmax_index], axis=0)
                         powb = numpy.average(dataOut.data_spc[pair[1], :, hmin_index:hmax_index], axis=0)
                         avgcoherenceComplex = ccf/numpy.sqrt(powa*powb)
                         phase = numpy.arctan2(avgcoherenceComplex.imag, avgcoherenceComplex.real)*180/numpy.pi
                         if dataOut.beacon_heiIndexList:
                             phase_beacon[i] = numpy.average(phase[dataOut.beacon_heiIndexList])
                         else:
                             phase_beacon[i] = numpy.average(phase)
                     if not self.isConfig:
                         nplots = len(pairsIndexList)
                         self.setup(id=id,
                                    nplots=nplots,
                                    wintitle=wintitle,
                                    showprofile=showprofile,
                                    show=show)
                         if timerange != None:
                             self.timerange = timerange
                         self.xmin, self.xmax = self.getTimeLim(x, xmin, xmax, timerange)
                         if ymin == None: ymin = 0
                         if ymax == None: ymax = 360
                         self.FTP_WEI = ftp_wei
                         self.EXP_CODE = exp_code
                         self.SUB_EXP_CODE = sub_exp_code
                         self.PLOT_POS = plot_pos
                         self.name = thisDatetime.strftime("%Y%m%d_%H%M%S")
                         self.isConfig = True
                         self.figfile = figfile
                         self.xdata = numpy.array([])
                         self.ydata = numpy.array([])
                         update_figfile = True
                         #open file beacon phase
                         path = '%s%03d' %(self.PREFIX, self.id)
                         beacon_file = os.path.join(path,'%s.txt'%self.name)
                         self.filename_phase = os.path.join(figpath,beacon_file)
                     self.setWinTitle(title)
                     title = "Phase Plot %s" %(thisDatetime.strftime("%Y/%m/%d %H:%M:%S"))
                     legendlabels = ["Pair (%d,%d)"%(pair[0], pair[1]) for pair in dataOut.pairsList]
                     axes = self.axesList[0]
                     self.xdata = numpy.hstack((self.xdata, x[0:1]))
                     if len(self.ydata)==0:
                         self.ydata = phase_beacon.reshape(-1,1)
                     else:
                         self.ydata = numpy.hstack((self.ydata, phase_beacon.reshape(-1,1)))
                     axes.pmultilineyaxis(x=self.xdata, y=self.ydata,
                                 xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax,
                                 xlabel=xlabel, ylabel=ylabel, title=title, legendlabels=legendlabels, marker='x', markersize=8, linestyle="solid",
                                 XAxisAsTime=True, grid='both'
                                 )
                     self.draw()
                     if dataOut.ltctime >= self.xmax:
                         self.counter_imagwr = wr_period
                         self.isConfig = False
                         update_figfile = True
                     self.save(figpath=figpath,
                               figfile=figfile,
                               save=save,
                               ftp=ftp,
                               wr_period=wr_period,
                               thisDatetime=thisDatetime,
                               update_figfile=update_figfile)
                     return dataOut
             #####################################
             class NoiselessSpectraPlot(Plot):
                 '''
                 Plot for Spectra data, subtracting
                 the noise in all channels, using for
                 amisr-14 data
                 '''
                 CODE = 'noiseless_spc'
                 colormap = 'jet'
                 plot_type = 'pcolor'
                 buffering = False
                 channelList = []
                 last_noise = None
                 def setup(self):
                     self.nplots = len(self.data.channels)
                     self.ncols = int(numpy.sqrt(self.nplots) + 0.9)
                     self.nrows = int((1.0 * self.nplots / self.ncols) + 0.9)
                     self.height = 3.5 * self.nrows
                     self.cb_label = 'dB'
                     if self.showprofile:
                         self.width = 5.8 * self.ncols
                     else:
                         self.width = 4.8* self.ncols
                     self.plots_adjust.update({'wspace': 0.4, 'hspace':0.4, 'left': 0.1, 'right': 0.92, 'bottom': 0.12})
                     self.ylabel = 'Range [km]'
                 def update_list(self,dataOut):
                     if len(self.channelList) == 0:
                         self.channelList = dataOut.channelList
                 def update(self, dataOut):
                     self.update_list(dataOut)
                     data = {}
                     meta = {}
                     norm = dataOut.nProfiles * dataOut.max_nIncohInt * dataOut.nCohInt  * dataOut.windowOfFilter
                     n0 = (dataOut.getNoise()/norm)
                     noise = numpy.repeat(n0,(dataOut.nFFTPoints*dataOut.nHeights)).reshape(dataOut.nChannels,dataOut.nFFTPoints,dataOut.nHeights)
                     noise = 10*numpy.log10(noise)
                     z = numpy.zeros((dataOut.nChannels, dataOut.nFFTPoints, dataOut.nHeights))
                     for ch in range(dataOut.nChannels):
                         if hasattr(dataOut.normFactor,'ndim'):
                             if dataOut.normFactor.ndim > 1:
                                 z[ch] = (numpy.divide(dataOut.data_spc[ch],dataOut.normFactor[ch]))
                             else:
                                 z[ch] = (numpy.divide(dataOut.data_spc[ch],dataOut.normFactor))
                         else:
                             z[ch] = (numpy.divide(dataOut.data_spc[ch],dataOut.normFactor))
                     z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
                     spc = 10*numpy.log10(z)
                     data['spc'] = spc - noise
                     #print(spc.shape)
                     data['rti'] = spc.mean(axis=1)
                     data['noise'] = noise
                     # data['noise'] = noise
                     meta['xrange'] = (dataOut.getFreqRange(EXTRA_POINTS)/1000., dataOut.getAcfRange(EXTRA_POINTS), dataOut.getVelRange(EXTRA_POINTS))
                     return data, meta
                 def plot(self):
                     if self.xaxis == "frequency":
                         x = self.data.xrange[0]
                         self.xlabel = "Frequency (kHz)"
                     elif self.xaxis == "time":
                         x = self.data.xrange[1]
                         self.xlabel = "Time (ms)"
                     else:
                         x = self.data.xrange[2]
                         self.xlabel = "Velocity (m/s)"
                     self.titles = []
                     y = self.data.yrange
                     self.y = y
                     data = self.data[-1]
                     z = data['spc']
                     for n, ax in enumerate(self.axes):
                         #noise = data['noise'][n]
                         if ax.firsttime:
                             self.xmax = self.xmax if self.xmax else numpy.nanmax(x)
                             self.xmin = self.xmin if self.xmin else -self.xmax
                             self.zmin = self.zmin if self.zmin else numpy.nanmin(z)
                             self.zmax = self.zmax if self.zmax else numpy.nanmax(z)
                             ax.plt = ax.pcolormesh(x, y, z[n].T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
                             if self.showprofile:
                                 ax.plt_profile = self.pf_axes[n].plot(
                                     data['rti'][n], y)[0]
                         else:
                             ax.plt.set_array(z[n].T.ravel())
                             if self.showprofile:
                                 ax.plt_profile.set_data(data['rti'][n], y)
                         self.titles.append('CH {}'.format(self.channelList[n]))
             class NoiselessRTIPlot(RTIPlot):
                 '''
                 Plot for RTI data
                 '''
                 CODE = 'noiseless_rti'
                 colormap = 'jet'
                 plot_type = 'pcolorbuffer'
                 titles = None
                 channelList = []
                 elevationList = []
                 azimuthList = []
                 last_noise = None
                 def setup(self):
                     self.xaxis = 'time'
                     self.ncols = 1
                     #print("dataChannels ",self.data.channels)
                     self.nrows = len(self.data.channels)
                     self.nplots = len(self.data.channels)
                     self.ylabel = 'Range [km]'
                     #self.xlabel = 'Time'
                     self.cb_label = 'dB'
                     self.plots_adjust.update({'hspace':0.8, 'left': 0.08, 'bottom': 0.2, 'right':0.94})
                     self.titles = ['{} Channel {}'.format(
                         self.CODE.upper(), x) for x in range(self.nplots)]
                 def update_list(self,dataOut):
                     if len(self.channelList) == 0:
                         self.channelList = dataOut.channelList
                     if len(self.elevationList) == 0:
                         self.elevationList = dataOut.elevationList
                     if len(self.azimuthList) == 0:
                         self.azimuthList = dataOut.azimuthList
                 def update(self, dataOut):
                     if len(self.channelList) == 0:
                         self.update_list(dataOut)
                     data = {}
                     meta = {}
                     #print(dataOut.max_nIncohInt, dataOut.nIncohInt)
                     #print(dataOut.windowOfFilter,dataOut.nCohInt,dataOut.nProfiles,dataOut.max_nIncohInt,dataOut.nIncohInt
                     norm = dataOut.nProfiles * dataOut.max_nIncohInt * dataOut.nCohInt  * dataOut.windowOfFilter
                     n0 = 10*numpy.log10(dataOut.getNoise()/norm)
                     data['noise'] = n0
                     noise = numpy.repeat(n0,dataOut.nHeights).reshape(dataOut.nChannels,dataOut.nHeights)
                     noiseless_data = dataOut.getPower() - noise
                     #print("power, noise:", dataOut.getPower(), n0)
                     #print(noise)
                     #print(noiseless_data)
                     data['noiseless_rti'] = noiseless_data
                     return data, meta
                 def plot(self):
                     from matplotlib import pyplot as plt
                     self.x = self.data.times
                     self.y = self.data.yrange
                     self.z = self.data['noiseless_rti']
                     self.z = numpy.array(self.z, dtype=float)
                     self.z = numpy.ma.masked_invalid(self.z)
                     try:
                         if self.channelList != None:
                             if len(self.elevationList) > 0 and len(self.azimuthList) > 0:
                                 self.titles = ['{}  Channel {} ({:2.1f} Elev,  {:2.1f} Azth)'.format(
                                     self.CODE.upper(), x, self.elevationList[x], self.azimuthList[x]) for x in self.channelList]
                             else:
                                 self.titles = ['{}  Channel {}'.format(
                                     self.CODE.upper(), x) for x in self.channelList]
                     except:
                         if self.channelList.any() != None:
                             if len(self.elevationList) > 0 and len(self.azimuthList) > 0:
                                 self.titles = ['{}  Channel {} ({:2.1f} Elev,  {:2.1f} Azth)'.format(
                                     self.CODE.upper(), x, self.elevationList[x], self.azimuthList[x]) for x in self.channelList]
                             else:
                                 self.titles = ['{} Channel {}'.format(
                                     self.CODE.upper(), x) for x in self.channelList]
                     if self.decimation is None:
                         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
                     #print("plot shapes ", z.shape, x.shape, y.shape)
                     #print(self.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:
                             if (n+1) == len(self.channelList):
                                 ax.set_xlabel('Time')
                             ax.plt = ax.pcolormesh(x, y, z[n].T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
                             if self.showprofile:
                                 ax.plot_profile = self.pf_axes[n].plot(data['noiseless_rti'][n], self.y)[0]
                         else:
                             # ax.collections.remove(ax.collections[0]) #  error while running
                             ax.plt = ax.pcolormesh(x, y, z[n].T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=plt.get_cmap(self.colormap)
                                                    )
                             if self.showprofile:
                                 ax.plot_profile.set_data(data['noiseless_rti'][n], self.y)
                                 # if "noise" in self.data:
                                 #     #ax.plot_noise.set_data(numpy.repeat(data['noise'][n], len(self.y)), self.y)
                                 #     ax.plot_noise.set_data(data['noise'][n], self.y)
             class OutliersRTIPlot(Plot):
                 '''
                 Plot for data_xxxx object
                 '''
                 CODE = 'outlier_rtc' # Range Time Counts
                 colormap = 'cool'
                 plot_type = 'pcolorbuffer'
                 def setup(self):
                     self.xaxis = 'time'
                     self.ncols = 1
                     self.nrows = self.data.shape('outlier_rtc')[0]
                     self.nplots = self.nrows
                     self.plots_adjust.update({'hspace':0.8, 'left': 0.08, 'bottom': 0.2, 'right':0.94})
                     if not self.xlabel:
                         self.xlabel = 'Time'
                     self.ylabel = 'Height [km]'
                     if not self.titles:
                         self.titles = ['Outliers Ch:{}'.format(x) for x in range(self.nrows)]
                 def update(self, dataOut):
                     data = {}
                     data['outlier_rtc'] = dataOut.data_outlier
                     meta = {}
                     return data, meta
                 def plot(self):
                     # self.data.normalize_heights()
                     self.x = self.data.times
                     self.y = self.data.yrange
                     self.z = self.data['outlier_rtc']
                     #self.z = numpy.ma.masked_invalid(self.z)
                     if self.decimation is None:
                         x, y, z = self.fill_gaps(self.x, self.y, self.z)
                     else:
                         x, y, z = self.fill_gaps(*self.decimate())
                     for n, ax in enumerate(self.axes):
                         self.zmax = self.zmax if self.zmax is not None else numpy.max(
                             self.z[n])
                         self.zmin = self.zmin if self.zmin is not None else numpy.min(
                             self.z[n])
                         data = self.data[-1]
                         if ax.firsttime:
                             if self.zlimits is not None:
                                 self.zmin, self.zmax = self.zlimits[n]
                             ax.plt = ax.pcolormesh(x, y, z[n].T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=self.cmaps[n]
                                                    )
                             if self.showprofile:
                                 ax.plot_profile = self.pf_axes[n].plot(data['outlier_rtc'][n], self.y)[0]
                                 self.pf_axes[n].set_xlabel('')
                         else:
                             if self.zlimits is not None:
                                 self.zmin, self.zmax = self.zlimits[n]
                             # ax.collections.remove(ax.collections[0]) #  error while running
                             ax.plt = ax.pcolormesh(x, y, z[n].T ,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=self.cmaps[n]
                                                    )
                             if self.showprofile:
                                 ax.plot_profile.set_data(data['outlier_rtc'][n], self.y)
                                 self.pf_axes[n].set_xlabel('')
             class NIncohIntRTIPlot(Plot):
                 '''
                 Plot for data_xxxx object
                 '''
                 CODE = 'integrations_rtc' # Range Time Counts
                 colormap = 'BuGn'
                 plot_type = 'pcolorbuffer'
                 def setup(self):
                     self.xaxis = 'time'
                     self.ncols = 1
                     self.nrows = self.data.shape('integrations_rtc')[0]
                     self.nplots = self.nrows
                     self.plots_adjust.update({'hspace':0.8, 'left': 0.08, 'bottom': 0.2, 'right':0.94})
                     if not self.xlabel:
                         self.xlabel = 'Time'
                     self.ylabel = 'Height [km]'
                     if not self.titles:
                         self.titles = ['Integration Ch:{}'.format(x) for x in range(self.nrows)]
                 def update(self, dataOut):
                     data = {}
                     data['integrations_rtc'] = dataOut.nIncohInt
                     meta = {}
                     return data, meta
                 def plot(self):
                     # self.data.normalize_heights()
                     self.x = self.data.times
                     self.y = self.data.yrange
                     self.z = self.data['integrations_rtc']
                     #self.z = numpy.ma.masked_invalid(self.z)
                     if self.decimation is None:
                         x, y, z = self.fill_gaps(self.x, self.y, self.z)
                     else:
                         x, y, z = self.fill_gaps(*self.decimate())
                     for n, ax in enumerate(self.axes):
                         self.zmax = self.zmax if self.zmax is not None else numpy.max(
                             self.z[n])
                         self.zmin = self.zmin if self.zmin is not None else numpy.min(
                             self.z[n])
                         data = self.data[-1]
                         if ax.firsttime:
                             if self.zlimits is not None:
                                 self.zmin, self.zmax = self.zlimits[n]
                             ax.plt = ax.pcolormesh(x, y, z[n].T,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=self.cmaps[n]
                                                    )
                             if self.showprofile:
                                 ax.plot_profile = self.pf_axes[n].plot(data['integrations_rtc'][n], self.y)[0]
                                 self.pf_axes[n].set_xlabel('')
                         else:
                             if self.zlimits is not None:
                                 self.zmin, self.zmax = self.zlimits[n]
                             # ax.collections.remove(ax.collections[0]) #  error while running
                             ax.plt = ax.pcolormesh(x, y, z[n].T ,
                                                    vmin=self.zmin,
                                                    vmax=self.zmax,
                                                    cmap=self.cmaps[n]
                                                    )
                             if self.showprofile:
                                 ax.plot_profile.set_data(data['integrations_rtc'][n], self.y)
                                 self.pf_axes[n].set_xlabel('')
             class RTIMapPlot(Plot):
                 '''
                 Plot for RTI data
                 Example:
                 controllerObj = Project()
                 controllerObj.setup(id = '11', name='eej_proc', description=desc)
                 ##.......................................................................................
                 ##.......................................................................................
                 readUnitConfObj = controllerObj.addReadUnit(datatype='AMISRReader', path=inPath, startDate='2023/05/24',endDate='2023/05/24',
                                 startTime='12:00:00',endTime='12:45:59',walk=1,timezone='lt',margin_days=1,code = code,nCode = nCode,
                                 nBaud = nBaud,nOsamp = nosamp,nChannels=nChannels,nFFT=NFFT,
                                 syncronization=False,shiftChannels=0)
                 volts_proc = controllerObj.addProcUnit(datatype='VoltageProc', inputId=readUnitConfObj.getId())
                 opObj01 = volts_proc.addOperation(name='Decoder', optype='other')
                 opObj01.addParameter(name='code', value=code, format='floatlist')
                 opObj01.addParameter(name='nCode', value=1, format='int')
                 opObj01.addParameter(name='nBaud', value=nBaud, format='int')
                 opObj01.addParameter(name='osamp', value=nosamp, format='int')
                 opObj12 = volts_proc.addOperation(name='selectHeights', optype='self')
                 opObj12.addParameter(name='minHei', value='90', format='float')
                 opObj12.addParameter(name='maxHei', value='150', format='float')
                 proc_spc = controllerObj.addProcUnit(datatype='SpectraProc', inputId=volts_proc.getId())
                 proc_spc.addParameter(name='nFFTPoints', value='8', format='int')
                 opObj11 = proc_spc.addOperation(name='IncohInt', optype='other')
                 opObj11.addParameter(name='n', value='1', format='int')
                 beamMapFile = "/home/japaza/Documents/AMISR_sky_mapper/UMET_beamcodes.csv"
                 opObj12 = proc_spc.addOperation(name='RTIMapPlot', optype='external')
                 opObj12.addParameter(name='selectedHeightsList', value='95, 100, 105, 110 ', format='int')
                 opObj12.addParameter(name='bField', value='100', format='int')
                 opObj12.addParameter(name='filename', value=beamMapFile, format='str')
                 '''
                 CODE = 'rti_skymap'
                 plot_type = 'scatter'
                 titles = None
                 colormap = 'jet'
                 channelList = []
                 elevationList = []
                 azimuthList = []
                 last_noise = None
                 flag_setIndex = False
                 heights = []
                 dcosx = []
                 dcosy = []
                 fullDcosy = None
                 fullDcosy = None
                 hindex = []
                 mapFile = False
                 #####  BField    ####
                 flagBField = False
                 dcosxB = []
                 dcosyB = []
                 Bmarker = ['+','*','D','x','s','>','o','^']
                 def setup(self):
                     self.xaxis = 'Range (Km)'
                     if len(self.selectedHeightsList) > 0:
                         self.nplots = len(self.selectedHeightsList)
                     else:
                         self.nplots = 4
                     self.ncols = int(numpy.ceil(self.nplots/2))
                     self.nrows = int(numpy.ceil(self.nplots/self.ncols))
                     self.ylabel = 'dcosy'
                     self.xlabel = 'dcosx'
                     self.colorbar = True
                     self.width = 6 + 4.1*self.nrows
                     self.height = 3 + 3.5*self.ncols
                     if self.extFile!=None:
                         try:
                             pointings =  numpy.genfromtxt(self.extFile, delimiter=',')
                             full_azi = pointings[:,1]
                             full_elev = pointings[:,2]
                             self.fullDcosx = numpy.cos(numpy.radians(full_elev))*numpy.sin(numpy.radians(full_azi))
                             self.fullDcosy = numpy.cos(numpy.radians(full_elev))*numpy.cos(numpy.radians(full_azi))
                             mapFile = True
                         except Exception as e:
                             self.extFile = None
                             print(e)
                 def update_list(self,dataOut):
                     if len(self.channelList) == 0:
                         self.channelList = dataOut.channelList
                     if len(self.elevationList) == 0:
                         self.elevationList = dataOut.elevationList
                     if len(self.azimuthList) == 0:
                         self.azimuthList = dataOut.azimuthList
                     a = numpy.radians(numpy.asarray(self.azimuthList))
                     e = numpy.radians(numpy.asarray(self.elevationList))
                     self.heights = dataOut.heightList
                     self.dcosx = numpy.cos(e)*numpy.sin(a)
                     self.dcosy = numpy.cos(e)*numpy.cos(a)
                     if len(self.bFieldList)>0:
                         datetObj = datetime.datetime.fromtimestamp(dataOut.utctime)
                         doy = datetObj.timetuple().tm_yday
                         year = datetObj.year
                         # self.dcosxB, self.dcosyB
                         ObjB = BField(year=year,doy=doy,site=2,heights=self.bFieldList)
                         [dcos, alpha, nlon, nlat] = ObjB.getBField()
                         alpha_location = numpy.zeros((nlon,2,len(self.bFieldList)))
                         for ih in range(len(self.bFieldList)):
                             alpha_location[:,0,ih] = dcos[:,0,ih,0]
                             for ilon in numpy.arange(nlon):
                                 myx = (alpha[ilon,:,ih])[::-1]
                                 myy = (dcos[ilon,:,ih,0])[::-1]
                                 tck = splrep(myx,myy,s=0)
                                 mydcosx = splev(ObjB.alpha_i,tck,der=0)
                                 myx = (alpha[ilon,:,ih])[::-1]
                                 myy = (dcos[ilon,:,ih,1])[::-1]
                                 tck = splrep(myx,myy,s=0)
                                 mydcosy = splev(ObjB.alpha_i,tck,der=0)
                                 alpha_location[ilon,:,ih] = numpy.array([mydcosx, mydcosy])
                             self.dcosxB.append(alpha_location[:,0,ih])
                             self.dcosyB.append(alpha_location[:,1,ih])
                         self.flagBField = True
                     if len(self.celestialList)>0:
                         #getBField(self.bFieldList, date)
                         #pass = kwargs.get('celestial', [])
                         pass
                 def update(self, dataOut):
                     if len(self.channelList) == 0:
                         self.update_list(dataOut)
                     if not self.flag_setIndex:
                         if len(self.selectedHeightsList)>0:
                             for sel_height in self.selectedHeightsList:
                                 index_list = numpy.where(self.heights >= sel_height)
                                 index_list = index_list[0]
                                 self.hindex.append(index_list[0])
                         self.flag_setIndex = True
                     data = {}
                     meta = {}
                     data['rti_skymap'] = dataOut.getPower()
                     norm = dataOut.nProfiles * dataOut.max_nIncohInt * dataOut.nCohInt  * dataOut.windowOfFilter
                     noise = 10*numpy.log10(dataOut.getNoise()/norm)
                     data['noise'] = noise
                     return data, meta
                 def plot(self):
                     self.x = self.dcosx
                     self.y = self.dcosy
                     self.z = self.data[-1]['rti_skymap']
                     self.z = numpy.array(self.z, dtype=float)
                     if len(self.hindex) > 0:
                         index = self.hindex
                     else:
                         index = numpy.arange(0, len(self.heights), int((len(self.heights))/4.2))
                     self.titles = ['Height {:.2f} km '.format(self.heights[i])+" " for i in index]
                     for n, ax in enumerate(self.axes):
                         if ax.firsttime:
                             self.xmax = self.xmax if self.xmax else numpy.nanmax(self.x)
                             self.xmin = self.xmin if self.xmin else numpy.nanmin(self.x)
                             self.ymax = self.ymax if self.ymax else numpy.nanmax(self.y)
                             self.ymin = self.ymin if self.ymin else numpy.nanmin(self.y)
                             self.zmax = self.zmax if self.zmax else numpy.nanmax(self.z)
                             self.zmin = self.zmin if self.zmin else numpy.nanmin(self.z)
                             if self.extFile!=None:
                                 ax.scatter(self.fullDcosx, self.fullDcosy, marker="+", s=20)
                             ax.plt = ax.scatter(self.x, self.y, c=self.z[:,index[n]], cmap = 'jet',vmin = self.zmin,
                                 s=60, marker="s", vmax = self.zmax)
                             ax.minorticks_on()
                             ax.grid(which='major', axis='both')
                             ax.grid(which='minor', axis='x')
                             if self.flagBField :
                                 for ih in range(len(self.bFieldList)):
                                     label = str(self.bFieldList[ih]) + ' km'
                                     ax.plot(self.dcosxB[ih], self.dcosyB[ih], color='k', marker=self.Bmarker[ih % 8],
                                                      label=label,   linestyle='--',  ms=4.0,lw=0.5)
                                 handles, labels = ax.get_legend_handles_labels()
                                 a = -0.05
                                 b = 1.15 - 1.19*(self.nrows)
                                 self.axes[0].legend(handles,labels, bbox_to_anchor=(a,b), prop={'size': (5.8+ 1.1*self.nplots)}, title='B Field ⊥')
                         else:
                             ax.plt = ax.scatter(self.x, self.y, c=self.z[:,index[n]], cmap = 'jet',vmin = self.zmin,
                                  s=80, marker="s", vmax = self.zmax)
                             if self.flagBField :
                                 for ih in range(len(self.bFieldList)):
                                     ax.plot (self.dcosxB[ih], self.dcosyB[ih], color='k', marker=self.Bmarker[ih % 8],
                                             linestyle='--', ms=4.0,lw=0.5)

schainpy/model/io/jroIO_kamisr.py +1 -2

             ''''
             Created on Set 9, 2015
             @author: roj-idl71 Karim Kuyeng
             @upgrade: 2021, Joab Apaza
             '''
             import os
             import sys
             import glob
             import fnmatch
             import datetime
             import time
             import re
             import h5py
             import numpy
             try:
                 from gevent import sleep
             except:
                 from time import sleep
             from schainpy.model.data.jroheaderIO import RadarControllerHeader, SystemHeader,ProcessingHeader
             from schainpy.model.data.jrodata import Voltage
             from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation, MPDecorator
             from numpy import imag
             from schainpy.utils import log
             class AMISRReader(ProcessingUnit):
                 '''
                 classdocs
                 '''
                 def __init__(self):
                     '''
                     Constructor
                     '''
                     ProcessingUnit.__init__(self)
                     self.set = None
                     self.subset = None
                     self.extension_file = '.h5'
                     self.dtc_str = 'dtc'
                     self.dtc_id = 0
                     self.status = True
                     self.isConfig = False
                     self.dirnameList = []
                     self.filenameList = []
                     self.fileIndex = None
                     self.flagNoMoreFiles = False
                     self.flagIsNewFile = 0
                     self.filename = ''
                     self.amisrFilePointer = None
                     self.beamCodeMap = None
                     self.azimuthList = []
                     self.elevationList = []
                     self.dataShape = None
                     self.flag_old_beams = False
                     self.flagAsync = False        #Use when the experiment has no syncronization
                     self.shiftChannels = 0
                     self.profileIndex = 0
                     self.beamCodeByFrame = None
                     self.radacTimeByFrame = None
                     self.dataset = None
                     self.__firstFile = True
                     self.buffer = None
                     self.timezone = 'ut'
                     self.__waitForNewFile = 20
                     self.__filename_online = None
                     #Is really necessary create the output object in the initializer
                     self.dataOut = Voltage()
                     self.dataOut.error=False
                     self.margin_days = 1
                     self.flag_ignoreFiles = False  #to activate the ignoring Files flag
                     self.flag_standby = False # just keep waiting, use when ignoring files
                     self.ignStartDateTime=None
                     self.ignEndDateTime=None
                 def setup(self,path=None,
                                 startDate=None,
                                 endDate=None,
                                 startTime=None,
                                 endTime=None,
                                 walk=True,
                                 timezone='ut',
                                 all=0,
                                 code = 1,
                                 nCode = 1,
                                 nBaud = 0,
                                 nOsamp = 0,
                                 online=False,
                                 old_beams=False,
                                 margin_days=1,
                                 nFFT = None,
                                 nChannels = None,
                                 ignStartDate=None,
                                 ignEndDate=None,
                                 ignStartTime=None,
                                 ignEndTime=None,
                                 syncronization=True,
                                 shiftChannels=0
                                 ):
                     self.timezone = timezone
                     self.all = all
                     self.online = online
                     self.flag_old_beams = old_beams
                     self.code = code
                     self.nCode = int(nCode)
                     self.nBaud = int(nBaud)
                     self.nOsamp = int(nOsamp)
                     self.margin_days = margin_days
                     self.__sampleRate = None
                     self.flagAsync = not syncronization
                     self.shiftChannels = shiftChannels
                     self.nFFT = nFFT
                     self.nChannels = nChannels
                     if ignStartTime!=None and ignEndTime!=None:
                         if ignStartDate!=None and ignEndDate!=None:
                             self.ignStartDateTime=datetime.datetime.combine(ignStartDate,ignStartTime)
                             self.ignEndDateTime=datetime.datetime.combine(ignEndDate,ignEndTime)
                         else:
                             self.ignStartDateTime=datetime.datetime.combine(startDate,ignStartTime)
                             self.ignEndDateTime=datetime.datetime.combine(endDate,ignEndTime)
                         self.flag_ignoreFiles = True
                     #self.findFiles()
                     if not(online):
                         #Busqueda de archivos offline
                         self.searchFilesOffLine(path, startDate, endDate, startTime, endTime, walk,)
                     else:
                         self.searchFilesOnLine(path, startDate, endDate, startTime,endTime,walk)
                     if not(self.filenameList):
                         raise schainpy.admin.SchainWarning("There is no files into the folder: %s"%(path))
                         #sys.exit(0)
                         self.dataOut.error = True
                     self.fileIndex = 0
                     self.readNextFile(online)
                     '''
                     Add code
                     '''
                     self.isConfig = True
                     # print("Setup Done")
                     pass
                 def readAMISRHeader(self,fp):
                     if self.isConfig and (not self.flagNoMoreFiles):
                         newShape = fp.get('Raw11/Data/Samples/Data').shape[1:]
                         if self.dataShape != newShape and newShape != None and not self.flag_standby:
                             raise schainpy.admin.SchainError("NEW FILE HAS A DIFFERENT SHAPE: ")
                             print(self.dataShape,newShape,"\n")
                             return 0
                     else:
                         self.dataShape = fp.get('Raw11/Data/Samples/Data').shape[1:]
                     header = 'Raw11/Data/RadacHeader'
                     if self.nChannels == None:
                         expFile = fp['Setup/Experimentfile'][()].decode()
                         linesExp = expFile.split("\n")
                         a = [line for line in linesExp if "nbeamcodes" in line]
                         self.nChannels = int(a[0][11:])
                     if not self.flagAsync:  #for experiments with no syncronization
                        self.shiftChannels = 0
                     self.beamCodeByPulse = fp.get(header+'/BeamCode') # LIST OF BEAMS PER PROFILE, TO BE USED ON REARRANGE
                     if (self.startDate > datetime.date(2021, 7, 15)) or self.flag_old_beams: #Se cambió la forma de extracción de Apuntes el 17 o forzar con flag de reorganización
                         self.beamcodeFile = fp['Setup/Beamcodefile'][()].decode()
                         self.trueBeams = self.beamcodeFile.split("\n")
                         self.trueBeams.pop()#remove last
                         if self.nFFT == None:
                             log.error("FFT or number of repetitions per channels is needed",self.name)
                         beams_idx = [k*self.nFFT for k in range(self.nChannels)]
                         beams = [self.trueBeams[b] for b in beams_idx]
                         self.beamCode = [int(x, 16) for x in beams]
                         if(self.flagAsync and self.shiftChannels == 0):
                             initBeam = self.beamCodeByPulse[0, 0]
                             self.shiftChannels = numpy.argwhere(self.beamCode ==initBeam)[0,0]
                     else:
                         _beamCode= fp.get('Raw11/Data/Beamcodes') #se usa la manera previa al cambio de apuntes
                         self.beamCode = _beamCode[0,:]
                     if self.beamCodeMap == None:
                         self.beamCodeMap = fp['Setup/BeamcodeMap']
                         for beam in self.beamCode:
                             beamAziElev = numpy.where(self.beamCodeMap[:,0]==beam)
                             beamAziElev = beamAziElev[0].squeeze()
                             self.azimuthList.append(self.beamCodeMap[beamAziElev,1])
                             self.elevationList.append(self.beamCodeMap[beamAziElev,2])
                             #print("Beamssss: ",self.beamCodeMap[beamAziElev,1],self.beamCodeMap[beamAziElev,2])
                     #print(self.beamCode)
                     #self.code = fp.get(header+'/Code') # NOT USE FOR THIS
                     self.frameCount = fp.get(header+'/FrameCount')# NOT USE FOR THIS
                     self.modeGroup = fp.get(header+'/ModeGroup')# NOT USE FOR THIS
                     self.nsamplesPulse = fp.get(header+'/NSamplesPulse')# TO GET NSA OR USING DATA FOR THAT
                     self.pulseCount = fp.get(header+'/PulseCount')# NOT USE FOR THIS
                     self.radacTime = fp.get(header+'/RadacTime')# 1st TIME ON FILE ANDE CALCULATE THE REST WITH IPP*nindexprofile
                     self.timeCount = fp.get(header+'/TimeCount')# NOT USE FOR THIS
                     self.timeStatus = fp.get(header+'/TimeStatus')# NOT USE FOR THIS
                     self.rangeFromFile = fp.get('Raw11/Data/Samples/Range')
                     self.frequency =  fp.get('Rx/Frequency')
                     txAus = fp.get('Raw11/Data/Pulsewidth') #seconds
                     self.baud = fp.get('Raw11/Data/TxBaud')
                     sampleRate = fp.get('Rx/SampleRate')
                     self.__sampleRate = sampleRate[()]
                     self.nblocks = self.pulseCount.shape[0] #nblocks
                     self.profPerBlockRAW = self.pulseCount.shape[1] #profiles per block in raw data
                     self.nprofiles = self.pulseCount.shape[1] #nprofile
                     #self.nsa = self.nsamplesPulse[0,0] #ngates
                     self.nsa = len(self.rangeFromFile[0])
                     self.nchannels = len(self.beamCode)
                     self.ippSeconds = (self.radacTime[0][1] -self.radacTime[0][0]) #Ipp in seconds
                     #print("IPPS secs: ",self.ippSeconds)
                     #self.__waitForNewFile = self.nblocks  # wait depending on the number of blocks since each block is 1 sec
                     self.__waitForNewFile = self.nblocks * self.nprofiles * self.ippSeconds # wait until new file is created
                     #filling radar controller header parameters
                     self.__ippKm = self.ippSeconds *.15*1e6 # in km
                     #self.__txA = txAus[()]*.15 #(ipp[us]*.15km/1us) in km
                     self.__txA = txAus[()] #seconds
                     self.__txAKm = self.__txA*1e6*.15
                     self.__txB = 0
                     nWindows=1
                     self.__nSamples = self.nsa
                     self.__firstHeight = self.rangeFromFile[0][0]/1000 #in km
                     self.__deltaHeight = (self.rangeFromFile[0][1] - self.rangeFromFile[0][0])/1000
                     #print("amisr-ipp:",self.ippSeconds, self.__ippKm)
                     #for now until understand why the code saved is different (code included even though code not in tuf file)
                     #self.__codeType = 0
                    # self.__nCode = None
                    # self.__nBaud = None
                     self.__code = self.code
                     self.__codeType = 0
                     if self.code != None:
                         self.__codeType = 1
                     self.__nCode = self.nCode
                     self.__nBaud = self.nBaud
                     self.__baudTX = self.__txA/(self.nBaud)
                     #self.__code = 0
                     #filling system header parameters
                     self.__nSamples = self.nsa
                     self.newProfiles = self.nprofiles/self.nchannels
                     self.__channelList = [n for n in range(self.nchannels)]
                     self.__frequency = self.frequency[0][0]
                     return 1
                 def createBuffers(self):
                     pass
                 def __setParameters(self,path='', startDate='',endDate='',startTime='', endTime='', walk=''):
                     self.path = path
                     self.startDate = startDate
                     self.endDate = endDate
                     self.startTime = startTime
                     self.endTime = endTime
                     self.walk = walk
                 def __checkPath(self):
                     if os.path.exists(self.path):
                         self.status = 1
                     else:
                         self.status = 0
                         print('Path:%s does not exists'%self.path)
                     return
                 def __selDates(self, amisr_dirname_format):
                     try:
                         year = int(amisr_dirname_format[0:4])
                         month = int(amisr_dirname_format[4:6])
                         dom = int(amisr_dirname_format[6:8])
                         thisDate = datetime.date(year,month,dom)
                         #margen de un día extra, igual luego se filtra for fecha y hora
                         if (thisDate>=(self.startDate - datetime.timedelta(days=self.margin_days)) and thisDate <= (self.endDate)+ datetime.timedelta(days=1)):
                             return amisr_dirname_format
                     except:
                         return None
                 def __findDataForDates(self,online=False):
                     if not(self.status):
                         return None
                     pat = '\d+.\d+'
                     dirnameList = [re.search(pat,x) for x in os.listdir(self.path)]
                     dirnameList = [x for x in dirnameList if x!=None]
                     dirnameList = [x.string for x in dirnameList]
                     if not(online):
                         dirnameList = [self.__selDates(x) for x in dirnameList]
                         dirnameList = [x for x in dirnameList if x!=None]
                     if len(dirnameList)>0:
                         self.status = 1
                         self.dirnameList = dirnameList
                         self.dirnameList.sort()
                     else:
                         self.status = 0
                         return None
                 def __getTimeFromData(self):
                     startDateTime_Reader = datetime.datetime.combine(self.startDate,self.startTime)
                     endDateTime_Reader = datetime.datetime.combine(self.endDate,self.endTime)
                     print('Filtering Files from %s to %s'%(startDateTime_Reader, endDateTime_Reader))
                     print('........................................')
                     filter_filenameList = []
                     self.filenameList.sort()
                     total_files = len(self.filenameList)
                     #for i in range(len(self.filenameList)-1):
                     for i in range(total_files):
                         filename = self.filenameList[i]
                         #print("file-> ",filename)
                         try:
                             fp = h5py.File(filename,'r')
                             time_str = fp.get('Time/RadacTimeString')
                             startDateTimeStr_File = time_str[0][0].decode('UTF-8').split('.')[0]
                             #startDateTimeStr_File = "2019-12-16 09:21:11"
                             junk = time.strptime(startDateTimeStr_File, '%Y-%m-%d %H:%M:%S')
                             startDateTime_File = datetime.datetime(junk.tm_year,junk.tm_mon,junk.tm_mday,junk.tm_hour, junk.tm_min, junk.tm_sec)
                             #endDateTimeStr_File = "2019-12-16 11:10:11"
                             endDateTimeStr_File = time_str[-1][-1].decode('UTF-8').split('.')[0]
                             junk = time.strptime(endDateTimeStr_File, '%Y-%m-%d %H:%M:%S')
                             endDateTime_File = datetime.datetime(junk.tm_year,junk.tm_mon,junk.tm_mday,junk.tm_hour, junk.tm_min, junk.tm_sec)
                             fp.close()
                             #print("check time", startDateTime_File)
                             if self.timezone == 'lt':
                                 startDateTime_File = startDateTime_File - datetime.timedelta(minutes = 300)
                                 endDateTime_File = endDateTime_File - datetime.timedelta(minutes = 300)
                             if (startDateTime_File >=startDateTime_Reader and endDateTime_File<=endDateTime_Reader):
                                 filter_filenameList.append(filename)
                             if (startDateTime_File>endDateTime_Reader):
                                 break
                         except Exception as e:
                             log.warning("Error opening file {} -> {}".format(os.path.split(filename)[1],e))
                     filter_filenameList.sort()
                     self.filenameList = filter_filenameList
                     return 1
                 def __filterByGlob1(self, dirName):
                     filter_files = glob.glob1(dirName, '*.*%s'%self.extension_file)
                     filter_files.sort()
                     filterDict = {}
                     filterDict.setdefault(dirName)
                     filterDict[dirName] = filter_files
                     return filterDict
                 def __getFilenameList(self, fileListInKeys, dirList):
                     for value in fileListInKeys:
                         dirName = list(value.keys())[0]
                         for file in value[dirName]:
                             filename = os.path.join(dirName, file)
                             self.filenameList.append(filename)
                 def __selectDataForTimes(self, online=False):
                     #aun no esta implementado el filtro for tiempo-> implementado en readNextFile
                     if not(self.status):
                         return None
                     dirList = [os.path.join(self.path,x) for x in self.dirnameList]
                     fileListInKeys = [self.__filterByGlob1(x) for x in dirList]
                     self.__getFilenameList(fileListInKeys, dirList)
                     if not(online):
                         #filtro por tiempo
                         if not(self.all):
                             self.__getTimeFromData()
                         if len(self.filenameList)>0:
                             self.status = 1
                             self.filenameList.sort()
                         else:
                             self.status = 0
                             return None
                     else:
                         #get the last file - 1
                         self.filenameList = [self.filenameList[-2]]
                     new_dirnameList = []
                     for dirname in self.dirnameList:
                         junk = numpy.array([dirname in x for x in self.filenameList])
                         junk_sum = junk.sum()
                         if junk_sum > 0:
                             new_dirnameList.append(dirname)
                     self.dirnameList = new_dirnameList
                     return 1
                 def searchFilesOnLine(self, path, startDate, endDate, startTime=datetime.time(0,0,0),
                                         endTime=datetime.time(23,59,59),walk=True):
                     if endDate ==None:
                      startDate = datetime.datetime.utcnow().date()
                      endDate = datetime.datetime.utcnow().date()
                     self.__setParameters(path=path, startDate=startDate, endDate=endDate,startTime = startTime,endTime=endTime, walk=walk)
                     self.__checkPath()
                     self.__findDataForDates(online=True)
                     self.dirnameList = [self.dirnameList[-1]]
                     self.__selectDataForTimes(online=True)
                     return
                 def searchFilesOffLine(self,
                                         path,
                                         startDate,
                                         endDate,
                                         startTime=datetime.time(0,0,0),
                                         endTime=datetime.time(23,59,59),
                                         walk=True):
                     self.__setParameters(path, startDate, endDate, startTime, endTime, walk)
                     self.__checkPath()
                     self.__findDataForDates()
                     self.__selectDataForTimes()
                     for i in range(len(self.filenameList)):
                         print("%s" %(self.filenameList[i]))
                     return
                 def __setNextFileOffline(self):
                     try:
                         self.filename = self.filenameList[self.fileIndex]
                         self.amisrFilePointer = h5py.File(self.filename,'r')
                         self.fileIndex += 1
                     except:
                         self.flagNoMoreFiles = 1
                         raise schainpy.admin.SchainError('No more files to read')
                         return 0
                     self.flagIsNewFile = 1
                     print("Setting the file: %s"%self.filename)
                     return 1
                 def __setNextFileOnline(self):
                     filename = self.filenameList[0]
                     if self.__filename_online != None:
                         self.__selectDataForTimes(online=True)
                         filename = self.filenameList[0]
                         wait = 0
                         self.__waitForNewFile=300 ## DEBUG:
                         while self.__filename_online == filename:
                             print('waiting %d seconds to get a new file...'%(self.__waitForNewFile))
                             if wait == 5:
                                 self.flagNoMoreFiles = 1
                                 return 0
                             sleep(self.__waitForNewFile)
                             self.__selectDataForTimes(online=True)
                             filename = self.filenameList[0]
                             wait += 1
                     self.__filename_online = filename
                     self.amisrFilePointer = h5py.File(filename,'r')
                     self.flagIsNewFile = 1
                     self.filename = filename
                     print("Setting the file: %s"%self.filename)
                     return 1
                 def readData(self):
                     buffer = self.amisrFilePointer.get('Raw11/Data/Samples/Data')
                     re = buffer[:,:,:,0]
                     im = buffer[:,:,:,1]
                     dataset = re + im*1j
                     self.radacTime = self.amisrFilePointer.get('Raw11/Data/RadacHeader/RadacTime')
                     timeset = self.radacTime[:,0]
                     return dataset,timeset
                 def reshapeData(self):
                     #print(self.beamCodeByPulse, self.beamCode, self.nblocks, self.nprofiles, self.nsa)
                     channels = self.beamCodeByPulse[0,:]
                     nchan = self.nchannels
                     #self.newProfiles = self.nprofiles/nchan #must be defined on filljroheader
                     nblocks = self.nblocks
                     nsamples = self.nsa
                     #print("Channels: ",self.nChannels)
                     #Dimensions : nChannels, nProfiles, nSamples
                     new_block = numpy.empty((nblocks, nchan, numpy.int_(self.newProfiles), nsamples), dtype="complex64")
                     ############################################
                     profPerCH = int(self.profPerBlockRAW / (self.nFFT* self.nChannels))
                     #profPerCH = int(self.profPerBlockRAW /  self.nChannels)
                     for thisChannel in range(nchan):
                         ich = thisChannel
                         idx_ch = [self.nFFT*(ich + nchan*k) for k in range(profPerCH)]
                         #print(idx_ch)
                         if self.nFFT > 1:
                             aux = [numpy.arange(i, i+self.nFFT) for i in idx_ch]
                             idx_ch  = None
                             idx_ch =aux
                             idx_ch = numpy.array(idx_ch, dtype=int).flatten()
                         else:
                             idx_ch = numpy.array(idx_ch, dtype=int)
                         #print(ich,profPerCH,idx_ch)
                         #print(numpy.where(channels==self.beamCode[ich])[0])
                         #new_block[:,ich,:,:] = self.dataset[:,numpy.where(channels==self.beamCode[ich])[0],:]
                         new_block[:,ich,:,:] = self.dataset[:,idx_ch,:]
                     new_block = numpy.transpose(new_block, (1,0,2,3))
                     new_block = numpy.reshape(new_block, (nchan,-1, nsamples))
                     if  self.flagAsync:
                         new_block = numpy.roll(new_block,  self.shiftChannels, axis=0)
                     return new_block
                 def updateIndexes(self):
                     pass
                 def fillJROHeader(self):
                     #fill radar controller header
                     #fill system header
                     self.dataOut.systemHeaderObj = SystemHeader(nSamples=self.__nSamples,
                                                                 nProfiles=self.newProfiles,
                                                                 nChannels=len(self.__channelList),
                                                                 adcResolution=14,
                                                                 pciDioBusWidth=32)
                     self.dataOut.type = "Voltage"
                     self.dataOut.data = None
                     self.dataOut.dtype = numpy.dtype([('real','<i8'),('imag','<i8')])
             #        self.dataOut.nChannels = 0
             #        self.dataOut.nHeights = 0
                     self.dataOut.nProfiles = self.newProfiles*self.nblocks
                     #self.dataOut.heightList = self.__firstHeigth + numpy.arange(self.__nSamples, dtype = numpy.float)*self.__deltaHeigth
                     ranges = numpy.reshape(self.rangeFromFile[()],(-1))
                     self.dataOut.heightList =  ranges/1000.0 #km
                     self.dataOut.channelList = self.__channelList
                     self.dataOut.blocksize = self.dataOut.nChannels * self.dataOut.nHeights
             #        self.dataOut.channelIndexList = None
                     # #self.dataOut.azimuthList = numpy.roll( numpy.array(self.azimuthList) ,self.shiftChannels)
                     # #self.dataOut.elevationList = numpy.roll(numpy.array(self.elevationList) ,self.shiftChannels)
                     # #self.dataOut.codeList = numpy.roll(numpy.array(self.beamCode), self.shiftChannels)
                     self.dataOut.azimuthList = self.azimuthList
                     self.dataOut.elevationList = self.elevationList
                     self.dataOut.codeList = self.beamCode
                     #print(self.dataOut.elevationList)
                     self.dataOut.flagNoData = True
                     #Set to TRUE if the data is discontinuous
                     self.dataOut.flagDiscontinuousBlock = False
                     self.dataOut.utctime = None
                     #self.dataOut.timeZone = -5 #self.__timezone/60  #timezone like jroheader, difference in minutes between UTC and localtime
                     if self.timezone == 'lt':
                         self.dataOut.timeZone = time.timezone / 60. #get the timezone in minutes
                     else:
                         self.dataOut.timeZone = 0 #by default time is UTC
                     self.dataOut.dstFlag = 0
                     self.dataOut.errorCount = 0
                     self.dataOut.nCohInt = 1
                     self.dataOut.flagDecodeData = False #asumo que la data esta decodificada
                     self.dataOut.flagDeflipData = False #asumo que la data esta sin flip
                     self.dataOut.flagShiftFFT = False
                     self.dataOut.ippSeconds = self.ippSeconds
                     self.dataOut.ipp = self.__ippKm
                     self.dataOut.nCode = self.__nCode
                     self.dataOut.code = self.__code
                     self.dataOut.nBaud = self.__nBaud
                     self.dataOut.frequency = self.__frequency
                     self.dataOut.realtime = self.online
                     self.dataOut.radarControllerHeaderObj = RadarControllerHeader(ipp=self.__ippKm,
                                                                                   txA=self.__txAKm,
                                                                                   txB=0,
                                                                                   nWindows=1,
                                                                                   nHeights=self.__nSamples,
                                                                                   firstHeight=self.__firstHeight,
                                                                                   codeType=self.__codeType,
                                                                                   nCode=self.__nCode, nBaud=self.__nBaud,
                                                                                   code = self.__code,
                                                                                   nOsamp=self.nOsamp,
                                                                                   frequency = self.__frequency,
                                                                                   sampleRate= self.__sampleRate,
                                                                                   fClock=self.__sampleRate)
                     self.dataOut.radarControllerHeaderObj.heightList = ranges/1000.0 #km
                     self.dataOut.radarControllerHeaderObj.heightResolution = self.__deltaHeight
                     self.dataOut.radarControllerHeaderObj.rangeIpp = self.__ippKm   #km
                     self.dataOut.radarControllerHeaderObj.rangeTxA = self.__txA*1e6*.15 #km
                     self.dataOut.radarControllerHeaderObj.nChannels = self.nchannels
                     self.dataOut.radarControllerHeaderObj.channelList = self.__channelList
                     self.dataOut.radarControllerHeaderObj.azimuthList = self.azimuthList
                     self.dataOut.radarControllerHeaderObj.elevationList = self.elevationList
                     self.dataOut.radarControllerHeaderObj.dtype = "Voltage"
                     self.dataOut.ippSeconds = self.ippSeconds
                     self.dataOut.ippFactor = self.nFFT
                     pass
                 def readNextFile(self,online=False):
                     if not(online):
                         newFile = self.__setNextFileOffline()
                     else:
                         newFile = self.__setNextFileOnline()
                     if not(newFile):
                         self.dataOut.error = True
                         return 0
                     if not self.readAMISRHeader(self.amisrFilePointer):
                         self.dataOut.error = True
                         return 0
                     #self.createBuffers()
                     self.fillJROHeader()
                     #self.__firstFile = False
                     self.dataset,self.timeset = self.readData()
                     if self.endDate!=None:
                         endDateTime_Reader = datetime.datetime.combine(self.endDate,self.endTime)
                         time_str = self.amisrFilePointer.get('Time/RadacTimeString')
                         startDateTimeStr_File = time_str[0][0].decode('UTF-8').split('.')[0]
                         junk = time.strptime(startDateTimeStr_File, '%Y-%m-%d %H:%M:%S')
                         startDateTime_File = datetime.datetime(junk.tm_year,junk.tm_mon,junk.tm_mday,junk.tm_hour, junk.tm_min, junk.tm_sec)
                         if self.timezone == 'lt':
                             startDateTime_File = startDateTime_File - datetime.timedelta(minutes = 300)
                         if (startDateTime_File>endDateTime_Reader):
                             self.flag_standby = False
                             return 0
                         if self.flag_ignoreFiles and (startDateTime_File >= self.ignStartDateTime and startDateTime_File <= self.ignEndDateTime):
                             print("Ignoring...")
                             self.flag_standby = True
+                            self.profileIndex = 99999999999999999
                             return 1
                         self.flag_standby = False
                     self.jrodataset = self.reshapeData()
                     #----self.updateIndexes()
                     self.profileIndex = 0
                     return 1
                 def __hasNotDataInBuffer(self):
                     if self.profileIndex >= (self.newProfiles*self.nblocks):
                         return 1
                     return 0
                 def getData(self):
                     if self.flagNoMoreFiles:
                         self.dataOut.flagNoData = True
                         return 0
                     if self.profileIndex >= (self.newProfiles*self.nblocks): #
                     #if self.__hasNotDataInBuffer():
                         if not (self.readNextFile(self.online)):
                             print("Profile Index break...")
                             return 0
                     if self.flag_standby: #Standby mode, if files are being ignoring, just return with no error flag
                         return 0
                     if self.dataset is None: # setear esta condicion cuando no hayan datos por leer
                         self.dataOut.flagNoData = True
                         print("No more data break...")
                         return 0
                     #self.dataOut.data = numpy.reshape(self.jrodataset[self.profileIndex,:],(1,-1))
                     self.dataOut.data = self.jrodataset[:,self.profileIndex,:]
                     #print("R_t",self.timeset)
                     #self.dataOut.utctime = self.jrotimeset[self.profileIndex]
                     #verificar basic header de jro data y ver si es compatible con este valor
                     #self.dataOut.utctime = self.timeset + (self.profileIndex * self.ippSeconds * self.nchannels)
                     indexprof = numpy.mod(self.profileIndex, self.newProfiles)
                     indexblock = self.profileIndex/self.newProfiles
                     #print (indexblock, indexprof)
                     diffUTC = 0
                     t_comp = (indexprof * self.ippSeconds * self.nchannels) + diffUTC #
                     #print("utc :",indexblock," __ ",t_comp)
                     #print(numpy.shape(self.timeset))
                     self.dataOut.utctime = self.timeset[numpy.int_(indexblock)] + t_comp
                     #self.dataOut.utctime = self.timeset[self.profileIndex] + t_comp
                     self.dataOut.profileIndex = self.profileIndex
                     #print("N profile:",self.profileIndex,self.newProfiles,self.nblocks,self.dataOut.utctime)
                     self.dataOut.flagNoData = False
                     # if indexprof == 0:
                     #     print("kamisr: ",self.dataOut.utctime)
                     self.profileIndex += 1
                     return self.dataOut.data #retorno necesario??
                 def run(self, **kwargs):
                     '''
                     This method will be called many times so here you should put all your code
                     '''
                     #print("running kamisr")
                     if not self.isConfig:
                         self.setup(**kwargs)
                         self.isConfig = True
                     self.getData()

schainpy/model/proc/jroproc_spectra.py +5 -3

             # Copyright (c) 2012-2020 Jicamarca Radio Observatory
             # All rights reserved.
             #
             # Distributed under the terms of the BSD 3-clause license.
             """Spectra processing Unit and operations
             Here you will find the processing unit `SpectraProc` and several operations
             to work with Spectra data type
             """
             import time
             import itertools
             import numpy
             from schainpy.model.proc.jroproc_base import ProcessingUnit, MPDecorator, Operation
             from schainpy.model.data.jrodata import Spectra
             from schainpy.model.data.jrodata import hildebrand_sekhon
             from schainpy.model.data import _noise
             from schainpy.utils import log
             import matplotlib.pyplot as plt
             from schainpy.model.io.utilsIO import getHei_index
             import datetime
             class SpectraProc(ProcessingUnit):
                 def __init__(self):
                     ProcessingUnit.__init__(self)
                     self.buffer = None
                     self.firstdatatime = None
                     self.profIndex = 0
                     self.dataOut = Spectra()
+                    self.dataOut.error=False
                     self.id_min = None
                     self.id_max = None
                     self.setupReq = False #Agregar a todas las unidades de proc
                     self.nsamplesFFT = 0
                 def __updateSpecFromVoltage(self):
                     self.dataOut.timeZone = self.dataIn.timeZone
                     self.dataOut.dstFlag = self.dataIn.dstFlag
                     self.dataOut.errorCount = self.dataIn.errorCount
                     self.dataOut.useLocalTime = self.dataIn.useLocalTime
                     try:
                         self.dataOut.processingHeaderObj = self.dataIn.processingHeaderObj.copy()
                     except:
                         pass
                     self.dataOut.radarControllerHeaderObj = self.dataIn.radarControllerHeaderObj.copy()
                     self.dataOut.radarControllerHeaderObj = self.dataIn.radarControllerHeaderObj.copy()
                     self.dataOut.ippSeconds = self.dataIn.ippSeconds
                     self.dataOut.ipp = self.dataIn.ipp
                     self.dataOut.systemHeaderObj = self.dataIn.systemHeaderObj.copy()
                     self.dataOut.channelList = self.dataIn.channelList
                     self.dataOut.heightList = self.dataIn.heightList
                     self.dataOut.dtype = numpy.dtype([('real', '<f4'), ('imag', '<f4')])
                     self.dataOut.nProfiles = self.dataOut.nFFTPoints
                     self.dataOut.flagDiscontinuousBlock = self.dataIn.flagDiscontinuousBlock
                     self.dataOut.utctime = self.firstdatatime
                     self.dataOut.flagDecodeData = self.dataIn.flagDecodeData
                     self.dataOut.flagDeflipData = self.dataIn.flagDeflipData
                     self.dataOut.flagShiftFFT = False
                     self.dataOut.nCohInt = self.dataIn.nCohInt
                     self.dataOut.nIncohInt = 1
                     self.dataOut.deltaHeight = self.dataIn.deltaHeight
                     self.dataOut.windowOfFilter = self.dataIn.windowOfFilter
                     self.dataOut.frequency = self.dataIn.frequency
                     self.dataOut.realtime = self.dataIn.realtime
                     self.dataOut.azimuth = self.dataIn.azimuth
                     self.dataOut.zenith = self.dataIn.zenith
                     self.dataOut.codeList = self.dataIn.codeList
                     self.dataOut.azimuthList = self.dataIn.azimuthList
                     self.dataOut.elevationList = self.dataIn.elevationList
                     self.dataOut.code = self.dataIn.code
                     self.dataOut.nCode = self.dataIn.nCode
                     self.dataOut.flagProfilesByRange = self.dataIn.flagProfilesByRange
                     self.dataOut.nProfilesByRange = self.dataIn.nProfilesByRange
                     self.dataOut.runNextUnit = self.dataIn.runNextUnit
                     try:
                         self.dataOut.step = self.dataIn.step
                     except:
                         pass
                 def __getFft(self):
                     """
                     Convierte valores de Voltaje a Spectra
                     Affected:
                         self.dataOut.data_spc
                         self.dataOut.data_cspc
                         self.dataOut.data_dc
                         self.dataOut.heightList
                         self.profIndex
                         self.buffer
                         self.dataOut.flagNoData
                     """
                     fft_volt = numpy.fft.fft(
                         self.buffer, n=self.dataOut.nFFTPoints, axis=1)
                     fft_volt = fft_volt.astype(numpy.dtype('complex'))
                     dc = fft_volt[:, 0, :]
                     # calculo de self-spectra
                     fft_volt = numpy.fft.fftshift(fft_volt, axes=(1,))
                     spc = fft_volt * numpy.conjugate(fft_volt)
                     spc = spc.real
                     blocksize = 0
                     blocksize += dc.size
                     blocksize += spc.size
                     cspc = None
                     pairIndex = 0
                     if self.dataOut.pairsList != None:
                         # calculo de cross-spectra
                         cspc = numpy.zeros(
                             (self.dataOut.nPairs, self.dataOut.nFFTPoints, self.dataOut.nHeights), dtype='complex')
                         for pair in self.dataOut.pairsList:
                             if pair[0] not in self.dataOut.channelList:
                                 raise ValueError("Error getting CrossSpectra: pair 0 of %s is not in channelList = %s" % (
                                     str(pair), str(self.dataOut.channelList)))
                             if pair[1] not in self.dataOut.channelList:
                                 raise ValueError("Error getting CrossSpectra: pair 1 of %s is not in channelList = %s" % (
                                     str(pair), str(self.dataOut.channelList)))
                             cspc[pairIndex, :, :] = fft_volt[pair[0], :, :] * \
                                 numpy.conjugate(fft_volt[pair[1], :, :])
                             pairIndex += 1
                         blocksize += cspc.size
                     self.dataOut.data_spc = spc
                     self.dataOut.data_cspc = cspc
                     self.dataOut.data_dc = dc
                     self.dataOut.blockSize = blocksize
                     self.dataOut.flagShiftFFT = False
                 def run(self, nProfiles=None, nFFTPoints=None, pairsList=None, ippFactor=None, shift_fft=False,
                         zeroPad=False, zeroPoints=0, runNextUnit=0):
                     self.dataIn.runNextUnit = runNextUnit
                     try:
                         type = self.dataIn.type.decode("utf-8")
                         self.dataIn.type = type
-                    except:
+                    except Exception as e:
+                        # print("spc -> ",e)
                         pass
                     if self.dataIn.type == "Spectra":
+                        #print("AQUI")
                         try:
                             self.dataOut.copy(self.dataIn)
                             self.dataOut.radarControllerHeaderObj = self.dataIn.radarControllerHeaderObj.copy()
                             self.dataOut.processingHeaderObj = self.dataIn.processingHeaderObj.copy()
                             self.dataOut.nProfiles = self.dataOut.nFFTPoints
                             #self.dataOut.nHeights = len(self.dataOut.heightList)
                         except Exception as e:
                             print("Error dataIn ",e)
                         if shift_fft:
                             #desplaza a la derecha en el eje 2 determinadas posiciones
                             shift = int(self.dataOut.nFFTPoints/2)
                             self.dataOut.data_spc = numpy.roll(self.dataOut.data_spc, shift , axis=1)
                             if self.dataOut.data_cspc is not None:
                                 #desplaza a la derecha en el eje 2 determinadas posiciones
                                 self.dataOut.data_cspc = numpy.roll(self.dataOut.data_cspc, shift, axis=1)
                         if pairsList:
                             self.__selectPairs(pairsList)
                     elif self.dataIn.type == "Voltage":
                         self.dataOut.flagNoData = True
                         self.dataOut.radarControllerHeaderObj = self.dataIn.radarControllerHeaderObj.copy()
                         self.dataOut.processingHeaderObj = self.dataIn.processingHeaderObj.copy()
                         if nFFTPoints == None:
                             raise ValueError("This SpectraProc.run() need nFFTPoints input variable")
                         if nProfiles == None:
                             nProfiles = nFFTPoints
                         if ippFactor == None:
                             self.dataOut.ippFactor = self.dataIn.ippFactor
                         else:
                             self.dataOut.ippFactor = ippFactor
                         if self.buffer is None:
                             if not zeroPad:
                                 self.buffer = numpy.zeros((self.dataIn.nChannels,
                                                        nProfiles,
                                                        self.dataIn.nHeights),
                                                       dtype='complex')
                                 zeroPoints = 0
                             else:
                                 self.buffer = numpy.zeros((self.dataIn.nChannels,
                                                        nFFTPoints+int(zeroPoints),
                                                        self.dataIn.nHeights),
                                                       dtype='complex')
                         self.dataOut.nFFTPoints = nFFTPoints
                         if self.buffer is None:
                             self.buffer = numpy.zeros((self.dataIn.nChannels,
                                                        nProfiles,
                                                        self.dataIn.nHeights),
                                                       dtype='complex')
                         if self.dataIn.flagDataAsBlock:
                             nVoltProfiles = self.dataIn.data.shape[1]
                             zeroPoints = 0
                             if nVoltProfiles == nProfiles or zeroPad:
                                 self.buffer = self.dataIn.data.copy()
                                 self.profIndex = nVoltProfiles
                             elif nVoltProfiles < nProfiles:
                                 if self.profIndex == 0:
                                     self.id_min = 0
                                     self.id_max = nVoltProfiles
                                 self.buffer[:, self.id_min:self.id_max,
                                             :] = self.dataIn.data
                                 self.profIndex += nVoltProfiles
                                 self.id_min += nVoltProfiles
                                 self.id_max += nVoltProfiles
                             elif nVoltProfiles > nProfiles:
                                 self.reader.bypass = True
                                 if self.profIndex == 0:
                                     self.id_min = 0
                                     self.id_max = nProfiles
                                 self.buffer = self.dataIn.data[:, self.id_min:self.id_max,:]
                                 self.profIndex += nProfiles
                                 self.id_min += nProfiles
                                 self.id_max += nProfiles
                                 if self.id_max == nVoltProfiles:
                                     self.reader.bypass = False
                             else:
                                 raise ValueError("The type object %s has %d profiles, it should just has %d profiles" % (
                                     self.dataIn.type, self.dataIn.data.shape[1], nProfiles))
                                 self.dataOut.flagNoData = True
                         else:
                             self.buffer[:, self.profIndex, :] = self.dataIn.data.copy()
                             self.profIndex += 1
                         if self.firstdatatime == None:
                             self.firstdatatime = self.dataIn.utctime
                         if self.profIndex == nProfiles or (zeroPad and zeroPoints==0):
                             self.__updateSpecFromVoltage()
                             if pairsList == None:
                                 self.dataOut.pairsList = [pair for pair in itertools.combinations(self.dataOut.channelList, 2)]
                             else:
                                 self.dataOut.pairsList = pairsList
                             self.__getFft()
                             self.dataOut.flagNoData = False
                             self.firstdatatime = None
                             self.nsamplesFFT = self.profIndex
                             #if not self.reader.bypass:
                             self.profIndex = 0
                         #update Processing Header:
                         self.dataOut.processingHeaderObj.dtype = "Spectra"
                         self.dataOut.processingHeaderObj.nFFTPoints = self.dataOut.nFFTPoints
                         self.dataOut.processingHeaderObj.nSamplesFFT = self.nsamplesFFT
                         self.dataOut.processingHeaderObj.nIncohInt = 1
                     elif self.dataIn.type == "Parameters":  #when get data from h5 spc file
                         self.dataOut.data_spc = self.dataIn.data_spc
                         self.dataOut.data_cspc = self.dataIn.data_cspc
                         self.dataOut.data_outlier = self.dataIn.data_outlier
                         self.dataOut.nProfiles = self.dataIn.nProfiles
                         self.dataOut.nIncohInt = self.dataIn.nIncohInt
                         self.dataOut.nFFTPoints = self.dataIn.nFFTPoints
                         self.dataOut.ippFactor = self.dataIn.ippFactor
                         self.dataOut.max_nIncohInt = self.dataIn.max_nIncohInt
                         self.dataOut.radarControllerHeaderObj = self.dataIn.radarControllerHeaderObj.copy()
                         self.dataOut.ProcessingHeader = self.dataIn.ProcessingHeader.copy()
                         self.dataOut.ippSeconds = self.dataIn.ippSeconds
                         self.dataOut.ipp = self.dataIn.ipp
                         #self.dataOut.abscissaList = self.dataIn.getVelRange(1)
                         #self.dataOut.spc_noise = self.dataIn.getNoise()
                         #self.dataOut.spc_range = (self.dataIn.getFreqRange(1) , self.dataIn.getAcfRange(1) , self.dataIn.getVelRange(1))
                         # self.dataOut.normFactor = self.dataIn.normFactor
                         if hasattr(self.dataIn, 'channelList'):
                             self.dataOut.channelList = self.dataIn.channelList
                         if hasattr(self.dataIn, 'pairsList'):
                             self.dataOut.pairsList = self.dataIn.pairsList
                             self.dataOut.groupList = self.dataIn.pairsList
                         self.dataOut.flagNoData = False
                         if hasattr(self.dataIn, 'ChanDist'): #Distances of receiver channels
                             self.dataOut.ChanDist = self.dataIn.ChanDist
                         else: self.dataOut.ChanDist = None
                         #if hasattr(self.dataIn, 'VelRange'): #Velocities range
                         #    self.dataOut.VelRange = self.dataIn.VelRange
                         #else: self.dataOut.VelRange = None
                     else:
                         raise ValueError("The type of input object '%s' is not valid".format(
                             self.dataIn.type))
+                    # print("SPC done")
                 def __selectPairs(self, pairsList):
                     if not pairsList:
                         return
                     pairs = []
                     pairsIndex = []
                     for pair in pairsList:
                         if pair[0] not in self.dataOut.channelList or pair[1] not in self.dataOut.channelList:
                             continue
                         pairs.append(pair)
                         pairsIndex.append(pairs.index(pair))
                     self.dataOut.data_cspc = self.dataOut.data_cspc[pairsIndex]
                     self.dataOut.pairsList = pairs
                     return
                 def selectFFTs(self, minFFT, maxFFT ):
                     """
                     Selecciona un bloque de datos en base a un grupo de valores de puntos FFTs segun el rango
                     minFFT<= FFT <= maxFFT
                     """
                     if (minFFT > maxFFT):
                         raise ValueError("Error selecting heights: Height range (%d,%d) is not valid" % (minFFT, maxFFT))
                     if (minFFT < self.dataOut.getFreqRange()[0]):
                         minFFT = self.dataOut.getFreqRange()[0]
                     if (maxFFT > self.dataOut.getFreqRange()[-1]):
                         maxFFT = self.dataOut.getFreqRange()[-1]
                     minIndex = 0
                     maxIndex = 0
                     FFTs = self.dataOut.getFreqRange()
                     inda = numpy.where(FFTs >= minFFT)
                     indb = numpy.where(FFTs <= maxFFT)
                     try:
                         minIndex = inda[0][0]
                     except:
                         minIndex = 0
                     try:
                         maxIndex = indb[0][-1]
                     except:
                         maxIndex = len(FFTs)
                     self.selectFFTsByIndex(minIndex, maxIndex)
                     return 1
                 def getBeaconSignal(self, tauindex=0, channelindex=0, hei_ref=None):
                     newheis = numpy.where(
                         self.dataOut.heightList > self.dataOut.radarControllerHeaderObj.Taus[tauindex])
                     if hei_ref != None:
                         newheis = numpy.where(self.dataOut.heightList > hei_ref)
                     minIndex = min(newheis[0])
                     maxIndex = max(newheis[0])
                     data_spc = self.dataOut.data_spc[:, :, minIndex:maxIndex + 1]
                     heightList = self.dataOut.heightList[minIndex:maxIndex + 1]
                     # determina indices
                     nheis = int(self.dataOut.radarControllerHeaderObj.txB /
                                 (self.dataOut.heightList[1] - self.dataOut.heightList[0]))
                     avg_dB = 10 * \
                         numpy.log10(numpy.sum(data_spc[channelindex, :, :], axis=0))
                     beacon_dB = numpy.sort(avg_dB)[-nheis:]
                     beacon_heiIndexList = []
                     for val in avg_dB.tolist():
                         if val >= beacon_dB[0]:
                             beacon_heiIndexList.append(avg_dB.tolist().index(val))
                     data_cspc = None
                     if self.dataOut.data_cspc is not None:
                         data_cspc = self.dataOut.data_cspc[:, :, minIndex:maxIndex + 1]
                     data_dc = None
                     if self.dataOut.data_dc is not None:
                         data_dc = self.dataOut.data_dc[:, minIndex:maxIndex + 1]
                     self.dataOut.data_spc = data_spc
                     self.dataOut.data_cspc = data_cspc
                     self.dataOut.data_dc = data_dc
                     self.dataOut.heightList = heightList
                     self.dataOut.beacon_heiIndexList = beacon_heiIndexList
                     return 1
                 def selectFFTsByIndex(self, minIndex, maxIndex):
                     """
                     """
                     if (minIndex < 0) or (minIndex > maxIndex):
                         raise ValueError("Error selecting heights: Index range (%d,%d) is not valid" % (minIndex, maxIndex))
                     if (maxIndex >= self.dataOut.nProfiles):
                         maxIndex = self.dataOut.nProfiles-1
                     #Spectra
                     data_spc = self.dataOut.data_spc[:,minIndex:maxIndex+1,:]
                     data_cspc = None
                     if self.dataOut.data_cspc is not None:
                         data_cspc = self.dataOut.data_cspc[:,minIndex:maxIndex+1,:]
                     data_dc = None
                     if self.dataOut.data_dc is not None:
                         data_dc = self.dataOut.data_dc[minIndex:maxIndex+1,:]
                     self.dataOut.data_spc = data_spc
                     self.dataOut.data_cspc = data_cspc
                     self.dataOut.data_dc = data_dc
                     self.dataOut.ippSeconds = self.dataOut.ippSeconds*(self.dataOut.nFFTPoints / numpy.shape(data_cspc)[1])
                     self.dataOut.nFFTPoints = numpy.shape(data_cspc)[1]
                     self.dataOut.profilesPerBlock = numpy.shape(data_cspc)[1]
                     return 1
                 def getNoise(self, minHei=None, maxHei=None, minVel=None, maxVel=None):
                     # validacion de rango
                     if minHei == None:
                         minHei = self.dataOut.heightList[0]
                     if maxHei == None:
                         maxHei = self.dataOut.heightList[-1]
                     if (minHei < self.dataOut.heightList[0]) or (minHei > maxHei):
                         print('minHei: %.2f is out of the heights range' % (minHei))
                         print('minHei is setting to %.2f' % (self.dataOut.heightList[0]))
                         minHei = self.dataOut.heightList[0]
                     if (maxHei > self.dataOut.heightList[-1]) or (maxHei < minHei):
                         print('maxHei: %.2f is out of the heights range' % (maxHei))
                         print('maxHei is setting to %.2f' % (self.dataOut.heightList[-1]))
                         maxHei = self.dataOut.heightList[-1]
                     # validacion de velocidades
                     velrange = self.dataOut.getVelRange(1)
                     if minVel == None:
                         minVel = velrange[0]
                     if maxVel == None:
                         maxVel = velrange[-1]
                     if (minVel < velrange[0]) or (minVel > maxVel):
                         print('minVel: %.2f is out of the velocity range' % (minVel))
                         print('minVel is setting to %.2f' % (velrange[0]))
                         minVel = velrange[0]
                     if (maxVel > velrange[-1]) or (maxVel < minVel):
                         print('maxVel: %.2f is out of the velocity range' % (maxVel))
                         print('maxVel is setting to %.2f' % (velrange[-1]))
                         maxVel = velrange[-1]
                     # seleccion de indices para rango
                     minIndex = 0
                     maxIndex = 0
                     heights = self.dataOut.heightList
                     inda = numpy.where(heights >= minHei)
                     indb = numpy.where(heights <= maxHei)
                     try:
                         minIndex = inda[0][0]
                     except:
                         minIndex = 0
                     try:
                         maxIndex = indb[0][-1]
                     except:
                         maxIndex = len(heights)
                     if (minIndex < 0) or (minIndex > maxIndex):
                         raise ValueError("some value in (%d,%d) is not valid" % (
                             minIndex, maxIndex))
                     if (maxIndex >= self.dataOut.nHeights):
                         maxIndex = self.dataOut.nHeights - 1
                     # seleccion de indices para velocidades
                     indminvel = numpy.where(velrange >= minVel)
                     indmaxvel = numpy.where(velrange <= maxVel)
                     try:
                         minIndexVel = indminvel[0][0]
                     except:
                         minIndexVel = 0
                     try:
                         maxIndexVel = indmaxvel[0][-1]
                     except:
                         maxIndexVel = len(velrange)
                     # seleccion del espectro
                     data_spc = self.dataOut.data_spc[:,
                                                      minIndexVel:maxIndexVel + 1, minIndex:maxIndex + 1]
                     # estimacion de ruido
                     noise = numpy.zeros(self.dataOut.nChannels)
                     for channel in range(self.dataOut.nChannels):
                         daux = data_spc[channel, :, :]
                         sortdata = numpy.sort(daux, axis=None)
                         noise[channel] = hildebrand_sekhon(sortdata, self.dataOut.nIncohInt)
                     self.dataOut.noise_estimation = noise.copy()
                     return 1
             class GetSNR(Operation):
                 '''
                 Written by R. Flores
                 '''
                 """Operation to get SNR.
                 Parameters:
                 -----------
                 Example
                 --------
                 op = proc_unit.addOperation(name='GetSNR', optype='other')
                 """
                 def __init__(self, **kwargs):
                     Operation.__init__(self, **kwargs)
                 def run(self,dataOut):
                     noise = dataOut.getNoise(ymin_index=-10) #Región superior donde solo debería de haber ruido
                     dataOut.data_snr = (dataOut.data_spc.sum(axis=1)-noise[:,None]*dataOut.nFFTPoints)/(noise[:,None]*dataOut.nFFTPoints) #It works apparently
                     dataOut.snl = numpy.log10(dataOut.data_snr)
                     dataOut.snl = numpy.where(dataOut.data_snr<.01, numpy.nan, dataOut.snl)
                     return dataOut
             class removeDC(Operation):
                 def run(self, dataOut, mode=2):
                     self.dataOut = dataOut
                     jspectra = self.dataOut.data_spc
                     jcspectra = self.dataOut.data_cspc
                     num_chan = jspectra.shape[0]
                     num_hei = jspectra.shape[2]
                     if jcspectra is not None:
                         jcspectraExist = True
                         num_pairs = jcspectra.shape[0]
                     else:
                         jcspectraExist = False
                     freq_dc = int(jspectra.shape[1] / 2)
                     ind_vel = numpy.array([-2, -1, 1, 2]) + freq_dc
                     ind_vel = ind_vel.astype(int)
                     if ind_vel[0] < 0:
                         ind_vel[list(range(0, 1))] = ind_vel[list(range(0, 1))] + self.num_prof
                     if mode == 1:
                         jspectra[:, freq_dc, :] = (
                             jspectra[:, ind_vel[1], :] + jspectra[:, ind_vel[2], :]) / 2  # CORRECCION
                         if jcspectraExist:
                             jcspectra[:, freq_dc, :] = (
                                 jcspectra[:, ind_vel[1], :] + jcspectra[:, ind_vel[2], :]) / 2
                     if mode == 2:
                         vel = numpy.array([-2, -1, 1, 2])
                         xx = numpy.zeros([4, 4])
                         for fil in range(4):
                             xx[fil, :] = vel[fil]**numpy.asarray(list(range(4)))
                         xx_inv = numpy.linalg.inv(xx)
                         xx_aux = xx_inv[0, :]
                         for ich in range(num_chan):
                             yy = jspectra[ich, ind_vel, :]
                             jspectra[ich, freq_dc, :] = numpy.dot(xx_aux, yy)
                             junkid = jspectra[ich, freq_dc, :] <= 0
                             cjunkid = sum(junkid)
                             if cjunkid.any():
                                 jspectra[ich, freq_dc, junkid.nonzero()] = (
                                     jspectra[ich, ind_vel[1], junkid] + jspectra[ich, ind_vel[2], junkid]) / 2
                         if jcspectraExist:
                             for ip in range(num_pairs):
                                 yy = jcspectra[ip, ind_vel, :]
                                 jcspectra[ip, freq_dc, :] = numpy.dot(xx_aux, yy)
                     self.dataOut.data_spc = jspectra
                     self.dataOut.data_cspc = jcspectra
                     return self.dataOut
             class getNoiseB(Operation):
                 """
                 Get noise from custom heights and frequency ranges,
                 offset for additional manual correction
                 J. Apaza -> developed to amisr isr spectra
                 """
                 __slots__ =('offset','warnings', 'isConfig', 'minIndex','maxIndex','minIndexFFT','maxIndexFFT')
                 def __init__(self):
                     Operation.__init__(self)
                     self.isConfig = False
                 def setup(self, offset=None, minHei=None, maxHei=None,minVel=None, maxVel=None, minFreq= None, maxFreq=None, warnings=False):
                     self.warnings = warnings
                     if minHei == None:
                         minHei = self.dataOut.heightList[0]
                     if maxHei == None:
                         maxHei = self.dataOut.heightList[-1]
                     if (minHei < self.dataOut.heightList[0]) or (minHei > maxHei):
                         if self.warnings:
                             print('minHei: %.2f is out of the heights range' % (minHei))
                             print('minHei is setting to %.2f' % (self.dataOut.heightList[0]))
                         minHei = self.dataOut.heightList[0]
                     if (maxHei > self.dataOut.heightList[-1]) or (maxHei < minHei):
                         if self.warnings:
                             print('maxHei: %.2f is out of the heights range' % (maxHei))
                             print('maxHei is setting to %.2f' % (self.dataOut.heightList[-1]))
                         maxHei = self.dataOut.heightList[-1]
                     #indices relativos a los puntos de fft, puede ser de acuerdo a velocidad o frecuencia
                     minIndexFFT = 0
                     maxIndexFFT = 0
                     # validacion de velocidades
                     indminPoint = None
                     indmaxPoint = None
                     if self.dataOut.type == 'Spectra':
                         if minVel == None and maxVel == None :
                             freqrange = self.dataOut.getFreqRange(1)
                             if minFreq == None:
                                 minFreq = freqrange[0]
                             if maxFreq == None:
                                 maxFreq = freqrange[-1]
                             if (minFreq < freqrange[0]) or (minFreq > maxFreq):
                                 if self.warnings:
                                     print('minFreq: %.2f is out of the frequency range' % (minFreq))
                                     print('minFreq is setting to %.2f' % (freqrange[0]))
                                 minFreq = freqrange[0]
                             if (maxFreq > freqrange[-1]) or (maxFreq < minFreq):
                                 if self.warnings:
                                     print('maxFreq: %.2f is out of the frequency range' % (maxFreq))
                                     print('maxFreq is setting to %.2f' % (freqrange[-1]))
                                 maxFreq = freqrange[-1]
                             indminPoint = numpy.where(freqrange >= minFreq)
                             indmaxPoint = numpy.where(freqrange <= maxFreq)
                         else:
                             velrange = self.dataOut.getVelRange(1)
                             if minVel == None:
                                 minVel = velrange[0]
                             if maxVel == None:
                                 maxVel = velrange[-1]
                             if (minVel < velrange[0]) or (minVel > maxVel):
                                 if self.warnings:
                                     print('minVel: %.2f is out of the velocity range' % (minVel))
                                     print('minVel is setting to %.2f' % (velrange[0]))
                                 minVel = velrange[0]
                             if (maxVel > velrange[-1]) or (maxVel < minVel):
                                 if self.warnings:
                                     print('maxVel: %.2f is out of the velocity range' % (maxVel))
                                     print('maxVel is setting to %.2f' % (velrange[-1]))
                                 maxVel = velrange[-1]
                             indminPoint = numpy.where(velrange >= minVel)
                             indmaxPoint = numpy.where(velrange <= maxVel)
                     # seleccion de indices para rango  REEMPLAZAR FOR FUNCION EXTERNA LUEGO
                     # minIndex = 0
                     # maxIndex = 0
                     # heights = self.dataOut.heightList
                     # inda = numpy.where(heights >= minHei)
                     # indb = numpy.where(heights <= maxHei)
                     # try:
                     #     minIndex = inda[0][0]
                     # except:
                     #     minIndex = 0
                     # try:
                     #     maxIndex = indb[0][-1]
                     # except:
                     #     maxIndex = len(heights)
                     # if (minIndex < 0) or (minIndex > maxIndex):
                     #     raise ValueError("some value in (%d,%d) is not valid" % (
                     #         minIndex, maxIndex))
                     # if (maxIndex >= self.dataOut.nHeights):
                     #     maxIndex = self.dataOut.nHeights - 1
                     minIndex, maxIndex = getHei_index(minHei,maxHei,self.dataOut.heightList)
                     #############################################################3
                     # seleccion de indices para velocidades
                     if self.dataOut.type == 'Spectra':
                         try:
                             minIndexFFT = indminPoint[0][0]
                         except:
                             minIndexFFT = 0
                         try:
                             maxIndexFFT = indmaxPoint[0][-1]
                         except:
                             maxIndexFFT = len( self.dataOut.getFreqRange(1))
                     self.minIndex, self.maxIndex, self.minIndexFFT, self.maxIndexFFT = minIndex, maxIndex, minIndexFFT, maxIndexFFT
                     self.isConfig = True
                     self.offset = 1
                     if offset!=None:
                         self.offset = 10**(offset/10)
                 def run(self, dataOut, offset=None, minHei=None, maxHei=None,minVel=None, maxVel=None, minFreq= None, maxFreq=None, warnings=False):
                     self.dataOut = dataOut
                     if not self.isConfig:
                         self.setup(offset, minHei, maxHei,minVel, maxVel, minFreq, maxFreq, warnings)
                     self.dataOut.noise_estimation = None
                     noise = None
                     if self.dataOut.type == 'Voltage':
                         noise = self.dataOut.getNoise(ymin_index=self.minIndex, ymax_index=self.maxIndex)
                     elif self.dataOut.type == 'Spectra':
                         noise = numpy.zeros( self.dataOut.nChannels)
                         norm = 1
                         for channel in range( self.dataOut.nChannels):
                             if not hasattr(self.dataOut.nIncohInt,'__len__'):
                                 norm = 1
                             else:
                                 norm = self.dataOut.max_nIncohInt[channel]/self.dataOut.nIncohInt[channel, self.minIndex:self.maxIndex]
                             daux =  self.dataOut.data_spc[channel,self.minIndexFFT:self.maxIndexFFT, self.minIndex:self.maxIndex]
                             daux = numpy.multiply(daux, norm)
                             sortdata = numpy.sort(daux, axis=None)
                             noise[channel] = _noise.hildebrand_sekhon(sortdata, self.dataOut.max_nIncohInt[channel])/self.offset
                     else:
                         noise = self.dataOut.getNoise(xmin_index=self.minIndexFFT, xmax_index=self.maxIndexFFT, ymin_index=self.minIndex, ymax_index=self.maxIndex)
                     self.dataOut.noise_estimation = noise.copy() # dataOut.noise
                     return self.dataOut
                 def getNoiseByMean(self,data):
                     #data debe estar ordenado
                     data = numpy.mean(data,axis=1)
                     sortdata = numpy.sort(data, axis=None)
                     pnoise = None
                     j = 0
                     mean = numpy.mean(sortdata)
                     min = numpy.min(sortdata)
                     delta = mean - min
                     indexes = numpy.where(sortdata > (mean+delta))[0] #only array of indexes
                     #print(len(indexes))
                     if len(indexes)==0:
                         pnoise = numpy.mean(sortdata)
                     else:
                         j = indexes[0]
                         pnoise = numpy.mean(sortdata[0:j])
                     return pnoise
                 def getNoiseByHS(self,data, navg):
                     #data debe estar ordenado
                     #data = numpy.mean(data,axis=1)
                     sortdata = numpy.sort(data, axis=None)
                     lenOfData = len(sortdata)
                     nums_min = lenOfData*0.2
                     if nums_min <= 5:
                         nums_min = 5
                     sump = 0.
                     sumq = 0.
                     j = 0
                     cont = 1
                     while((cont == 1)and(j < lenOfData)):
                         sump += sortdata[j]
                         sumq += sortdata[j]**2
                         #sumq -= sump**2
                         if j > nums_min:
                             rtest = float(j)/(j-1) + 1.0/navg
                             #if ((sumq*j) > (sump**2)):
                             if ((sumq*j) > (rtest*sump**2)):
                                 j = j - 1
                                 sump = sump - sortdata[j]
                                 sumq = sumq - sortdata[j]**2
                                 cont = 0
                         j += 1
                     lnoise = sump / j
                     return lnoise
             class removeInterference(Operation):
                 def removeInterference2(self):
                     cspc = self.dataOut.data_cspc
                     spc = self.dataOut.data_spc
                     Heights = numpy.arange(cspc.shape[2])
                     realCspc = numpy.abs(cspc)
                     for i in range(cspc.shape[0]):
                         LinePower= numpy.sum(realCspc[i], axis=0)
                         Threshold = numpy.amax(LinePower)-numpy.sort(LinePower)[len(Heights)-int(len(Heights)*0.1)]
                         SelectedHeights = Heights[ numpy.where( LinePower < Threshold ) ]
                         InterferenceSum = numpy.sum( realCspc[i,:,SelectedHeights], axis=0 )
                         InterferenceThresholdMin = numpy.sort(InterferenceSum)[int(len(InterferenceSum)*0.98)]
                         InterferenceThresholdMax = numpy.sort(InterferenceSum)[int(len(InterferenceSum)*0.99)]
                         InterferenceRange = numpy.where( ([InterferenceSum > InterferenceThresholdMin]))# , InterferenceSum < InterferenceThresholdMax]) )
                         #InterferenceRange = numpy.where( ([InterferenceRange < InterferenceThresholdMax]))
                         if len(InterferenceRange)<int(cspc.shape[1]*0.3):
                             cspc[i,InterferenceRange,:] = numpy.NaN
                     self.dataOut.data_cspc = cspc
                 def removeInterference(self, interf=2, hei_interf=None, nhei_interf=None, offhei_interf=None):
                     jspectra = self.dataOut.data_spc
                     jcspectra = self.dataOut.data_cspc
                     jnoise = self.dataOut.getNoise()
                     num_incoh = self.dataOut.nIncohInt
                     num_channel = jspectra.shape[0]
                     num_prof = jspectra.shape[1]
                     num_hei = jspectra.shape[2]
                     # hei_interf
                     if hei_interf is None:
                         count_hei = int(num_hei / 2)
                         hei_interf = numpy.asmatrix(list(range(count_hei))) + num_hei - count_hei
                         hei_interf = numpy.asarray(hei_interf)[0]
                     # nhei_interf
                     if (nhei_interf == None):
                         nhei_interf = 5
                     if (nhei_interf < 1):
                         nhei_interf = 1
                     if (nhei_interf > count_hei):
                         nhei_interf = count_hei
                     if (offhei_interf == None):
                         offhei_interf = 0
                     ind_hei = list(range(num_hei))
             #         mask_prof = numpy.asarray(range(num_prof - 2)) + 1
             #         mask_prof[range(num_prof/2 - 1,len(mask_prof))] += 1
                     mask_prof = numpy.asarray(list(range(num_prof)))
                     num_mask_prof = mask_prof.size
                     comp_mask_prof = [0, num_prof / 2]
                     # noise_exist:    Determina si la variable jnoise ha sido definida y contiene la informacion del ruido de cada canal
                     if (jnoise.size < num_channel or numpy.isnan(jnoise).any()):
                         jnoise = numpy.nan
                     noise_exist = jnoise[0] < numpy.Inf
                     # Subrutina de Remocion de la Interferencia
                     for ich in range(num_channel):
                         # Se ordena los espectros segun su potencia (menor a mayor)
                         power = jspectra[ich, mask_prof, :]
                         power = power[:, hei_interf]
                         power = power.sum(axis=0)
                         psort = power.ravel().argsort()
                         # Se estima la interferencia promedio en los Espectros de Potencia empleando
                         junkspc_interf = jspectra[ich, :, hei_interf[psort[list(range(
                             offhei_interf, nhei_interf + offhei_interf))]]]
                         if noise_exist:
                             #    tmp_noise = jnoise[ich] / num_prof
                             tmp_noise = jnoise[ich]
                         junkspc_interf = junkspc_interf - tmp_noise
                         #junkspc_interf[:,comp_mask_prof] = 0
                         jspc_interf = junkspc_interf.sum(axis=0) / nhei_interf
                         jspc_interf = jspc_interf.transpose()
                         # Calculando el espectro de interferencia promedio
                         noiseid = numpy.where(
                             jspc_interf <= tmp_noise / numpy.sqrt(num_incoh))
                         noiseid = noiseid[0]
                         cnoiseid = noiseid.size
                         interfid = numpy.where(
                             jspc_interf > tmp_noise / numpy.sqrt(num_incoh))
                         interfid = interfid[0]
                         cinterfid = interfid.size
                         if (cnoiseid > 0):
                             jspc_interf[noiseid] = 0
                         # Expandiendo los perfiles a limpiar
                         if (cinterfid > 0):
                             new_interfid = (
                                 numpy.r_[interfid - 1, interfid, interfid + 1] + num_prof) % num_prof
                             new_interfid = numpy.asarray(new_interfid)
                             new_interfid = {x for x in new_interfid}
                             new_interfid = numpy.array(list(new_interfid))
                             new_cinterfid = new_interfid.size
                         else:
                             new_cinterfid = 0
                         for ip in range(new_cinterfid):
                             ind = junkspc_interf[:, new_interfid[ip]].ravel().argsort()
                             jspc_interf[new_interfid[ip]
                                         ] = junkspc_interf[ind[nhei_interf // 2], new_interfid[ip]]
                         jspectra[ich, :, ind_hei] = jspectra[ich, :,
                                                              ind_hei] - jspc_interf  # Corregir indices
                         # Removiendo la interferencia del punto de mayor interferencia
                         ListAux = jspc_interf[mask_prof].tolist()
                         maxid = ListAux.index(max(ListAux))
                         if cinterfid > 0:
                             for ip in range(cinterfid * (interf == 2) - 1):
                                 ind = (jspectra[ich, interfid[ip], :] < tmp_noise *
                                        (1 + 1 / numpy.sqrt(num_incoh))).nonzero()
                                 cind = len(ind)
                                 if (cind > 0):
                                     jspectra[ich, interfid[ip], ind] = tmp_noise * \
                                         (1 + (numpy.random.uniform(cind) - 0.5) /
                                          numpy.sqrt(num_incoh))
                             ind = numpy.array([-2, -1, 1, 2])
                             xx = numpy.zeros([4, 4])
                             for id1 in range(4):
                                 xx[:, id1] = ind[id1]**numpy.asarray(list(range(4)))
                             xx_inv = numpy.linalg.inv(xx)
                             xx = xx_inv[:, 0]
                             ind = (ind + maxid + num_mask_prof) % num_mask_prof
                             yy = jspectra[ich, mask_prof[ind], :]
                             jspectra[ich, mask_prof[maxid], :] = numpy.dot(
                                 yy.transpose(), xx)
                         indAux = (jspectra[ich, :, :] < tmp_noise *
                                   (1 - 1 / numpy.sqrt(num_incoh))).nonzero()
                         jspectra[ich, indAux[0], indAux[1]] = tmp_noise * \
                             (1 - 1 / numpy.sqrt(num_incoh))
                     # Remocion de Interferencia en el Cross Spectra
                     if jcspectra is None:
                         return jspectra, jcspectra
                     num_pairs = int(jcspectra.size / (num_prof * num_hei))
                     jcspectra = jcspectra.reshape(num_pairs, num_prof, num_hei)
                     for ip in range(num_pairs):
                         #-------------------------------------------
                         cspower = numpy.abs(jcspectra[ip, mask_prof, :])
                         cspower = cspower[:, hei_interf]
                         cspower = cspower.sum(axis=0)
                         cspsort = cspower.ravel().argsort()
                         junkcspc_interf = jcspectra[ip, :, hei_interf[cspsort[list(range(
                             offhei_interf, nhei_interf + offhei_interf))]]]
                         junkcspc_interf = junkcspc_interf.transpose()
                         jcspc_interf = junkcspc_interf.sum(axis=1) / nhei_interf
                         ind = numpy.abs(jcspc_interf[mask_prof]).ravel().argsort()
                         median_real = int(numpy.median(numpy.real(
                             junkcspc_interf[mask_prof[ind[list(range(3 * num_prof // 4))]], :])))
                         median_imag = int(numpy.median(numpy.imag(
                             junkcspc_interf[mask_prof[ind[list(range(3 * num_prof // 4))]], :])))
                         comp_mask_prof = [int(e) for e in comp_mask_prof]
                         junkcspc_interf[comp_mask_prof, :] = numpy.complex_(
                             median_real, median_imag)
                         for iprof in range(num_prof):
                             ind = numpy.abs(junkcspc_interf[iprof, :]).ravel().argsort()
                             jcspc_interf[iprof] = junkcspc_interf[iprof, ind[nhei_interf // 2]]
                         # Removiendo la Interferencia
                         jcspectra[ip, :, ind_hei] = jcspectra[ip,
                                                               :, ind_hei] - jcspc_interf
                         ListAux = numpy.abs(jcspc_interf[mask_prof]).tolist()
                         maxid = ListAux.index(max(ListAux))
                         ind = numpy.array([-2, -1, 1, 2])
                         xx = numpy.zeros([4, 4])
                         for id1 in range(4):
                             xx[:, id1] = ind[id1]**numpy.asarray(list(range(4)))
                         xx_inv = numpy.linalg.inv(xx)
                         xx = xx_inv[:, 0]
                         ind = (ind + maxid + num_mask_prof) % num_mask_prof
                         yy = jcspectra[ip, mask_prof[ind], :]
                         jcspectra[ip, mask_prof[maxid], :] = numpy.dot(yy.transpose(), xx)
                     # Guardar Resultados
                     self.dataOut.data_spc = jspectra
                     self.dataOut.data_cspc = jcspectra
                     return 1
                 def run(self, dataOut, interf=2,hei_interf=None, nhei_interf=None, offhei_interf=None, mode=1):
                     self.dataOut = dataOut
                     if mode == 1:
                         self.removeInterference(interf=2,hei_interf=None, nhei_interf=None, offhei_interf=None)
                     elif mode == 2:
                         self.removeInterference2()
                     return self.dataOut
             class deflip(Operation):
                 def run(self, dataOut):
                     # arreglo 1: (num_chan, num_profiles, num_heights)
                     self.dataOut = dataOut
                     # JULIA-oblicua, indice 2
                     # arreglo 2: (num_profiles, num_heights)
                     jspectra = self.dataOut.data_spc[2]
                     jspectra_tmp=numpy.zeros(jspectra.shape)
                     num_profiles=jspectra.shape[0]
                     freq_dc = int(num_profiles / 2)
                     # Flip con for
                     for j in range(num_profiles):
                      jspectra_tmp[num_profiles-j-1]= jspectra[j]
                     # Intercambio perfil de DC con perfil inmediato anterior
                     jspectra_tmp[freq_dc-1]= jspectra[freq_dc-1]
                     jspectra_tmp[freq_dc]= jspectra[freq_dc]
                     # canal modificado es re-escrito en el arreglo de canales
                     self.dataOut.data_spc[2] = jspectra_tmp
                     return self.dataOut
             class IncohInt(Operation):
                 __profIndex = 0
                 __withOverapping = False
                 __byTime = False
                 __initime = None
                 __lastdatatime = None
                 __integrationtime = None
                 __buffer_spc = None
                 __buffer_cspc = None
                 __buffer_dc = None
                 __dataReady = False
                 __timeInterval = None
                 incohInt = 0
                 nOutliers = 0
                 n = None
                 _flagProfilesByRange = False
                 _nProfilesByRange = 0
                 def __init__(self):
                     Operation.__init__(self)
                 def setup(self, n=None, timeInterval=None, overlapping=False):
                     """
                     Set the parameters of the integration class.
                     Inputs:
                         n        :    Number of coherent integrations
                         timeInterval   :    Time of integration. If the parameter "n" is selected this one does not work
                         overlapping    :
                     """
                     self.__initime = None
                     self.__lastdatatime = 0
                     self.__buffer_spc = 0
                     self.__buffer_cspc = 0
                     self.__buffer_dc = 0
                     self.__profIndex = 0
                     self.__dataReady = False
                     self.__byTime = False
                     self.incohInt = 0
                     self.nOutliers = 0
                     if n is None and timeInterval is None:
                         raise ValueError("n or timeInterval should be specified ...")
                     if n is not None:
                         self.n = int(n)
                     else:
                         self.__integrationtime = int(timeInterval)
                         self.n = None
                         self.__byTime = True
                 def putData(self, data_spc, data_cspc, data_dc):
                     """
                     Add a profile to the __buffer_spc and increase in one the __profileIndex
                     """
                     if data_spc.all() == numpy.nan :
                         print("nan ")
                         return
                     self.__buffer_spc += data_spc
                     if data_cspc is None:
                         self.__buffer_cspc = None
                     else:
                         self.__buffer_cspc += data_cspc
                     if data_dc is None:
                         self.__buffer_dc = None
                     else:
                         self.__buffer_dc += data_dc
                     self.__profIndex += 1
                     return
                 def pushData(self):
                     """
                     Return the sum of the last profiles and the profiles used in the sum.
                     Affected:
                     self.__profileIndex
                     """
                     data_spc = self.__buffer_spc
                     data_cspc = self.__buffer_cspc
                     data_dc = self.__buffer_dc
                     n = self.__profIndex
                     self.__buffer_spc = 0
                     self.__buffer_cspc = 0
                     self.__buffer_dc = 0
                     return data_spc, data_cspc, data_dc, n
                 def byProfiles(self, *args):
                     self.__dataReady = False
                     avgdata_spc = None
                     avgdata_cspc = None
                     avgdata_dc = None
                     self.putData(*args)
                     if self.__profIndex == self.n:
                         avgdata_spc, avgdata_cspc, avgdata_dc, n = self.pushData()
                         self.n = n
                         self.__dataReady = True
                     return avgdata_spc, avgdata_cspc, avgdata_dc
                 def byTime(self, datatime, *args):
                     self.__dataReady = False
                     avgdata_spc = None
                     avgdata_cspc = None
                     avgdata_dc = None
                     self.putData(*args)
                     if (datatime - self.__initime) >= self.__integrationtime:
                         avgdata_spc, avgdata_cspc, avgdata_dc, n = self.pushData()
                         self.n = n
                         self.__dataReady = True
                     return avgdata_spc, avgdata_cspc, avgdata_dc
                 def integrate(self, datatime, *args):
                     if self.__profIndex == 0:
                         self.__initime = datatime
                     if self.__byTime:
                         avgdata_spc, avgdata_cspc, avgdata_dc = self.byTime(
                             datatime, *args)
                     else:
                         avgdata_spc, avgdata_cspc, avgdata_dc = self.byProfiles(*args)
                     if not self.__dataReady:
                         return None, None, None, None
                     return self.__initime, avgdata_spc, avgdata_cspc, avgdata_dc
                 def run(self, dataOut, n=None, timeInterval=None, overlapping=False):
                     if n == 1:
                         return dataOut
                     if dataOut.flagNoData == True:
                         return dataOut
                     if dataOut.flagProfilesByRange == True:
                         self._flagProfilesByRange = True
                     dataOut.flagNoData = True
                     dataOut.processingHeaderObj.timeIncohInt = timeInterval
                     if not self.isConfig:
                         self._nProfilesByRange = numpy.zeros((1,len(dataOut.heightList)))
                         self.setup(n, timeInterval, overlapping)
                         self.isConfig = True
                     avgdatatime, avgdata_spc, avgdata_cspc, avgdata_dc = self.integrate(dataOut.utctime,
                                                                                         dataOut.data_spc,
                                                                                         dataOut.data_cspc,
                                                                                         dataOut.data_dc)
                     self.incohInt += dataOut.nIncohInt
                     if isinstance(dataOut.data_outlier,numpy.ndarray) or isinstance(dataOut.data_outlier,int) or isinstance(dataOut.data_outlier, float):
                         self.nOutliers += dataOut.data_outlier
                     if self._flagProfilesByRange:
                         dataOut.flagProfilesByRange = True
                         self._nProfilesByRange += dataOut.nProfilesByRange
                     if self.__dataReady:
                         #print("prof: ",dataOut.max_nIncohInt,self.__profIndex)
                         dataOut.data_spc = avgdata_spc
                         dataOut.data_cspc = avgdata_cspc
                         dataOut.data_dc = avgdata_dc
                         dataOut.nIncohInt = self.incohInt
                         dataOut.data_outlier = self.nOutliers
                         dataOut.utctime = avgdatatime
                         dataOut.flagNoData = False
                         self.incohInt = 0
                         self.nOutliers = 0
                         self.__profIndex = 0
                         dataOut.nProfilesByRange = self._nProfilesByRange
                         self._nProfilesByRange = numpy.zeros((1,len(dataOut.heightList)))
                         self._flagProfilesByRange = False
                         # print("IncohInt Done")
                     return dataOut
             class IntegrationFaradaySpectra(Operation):
                 __profIndex = 0
                 __withOverapping = False
                 __byTime = False
                 __initime = None
                 __lastdatatime = None
                 __integrationtime = None
                 __buffer_spc = None
                 __buffer_cspc = None
                 __buffer_dc = None
                 __dataReady = False
                 __timeInterval = None
                 n_ints = None #matriz de numero de integracions (CH,HEI)
                 n = None
                 minHei_ind = None
                 maxHei_ind = None
                 navg = 1.0
                 factor = 0.0
                 dataoutliers = None # (CHANNELS, HEIGHTS)
                 _flagProfilesByRange = False
                 _nProfilesByRange = 0
                 def __init__(self):
                     Operation.__init__(self)
                 def setup(self, dataOut,n=None, timeInterval=None, overlapping=False, DPL=None, minHei=None, maxHei=None, avg=1,factor=0.75):
                     """
                     Set the parameters of the integration class.
                     Inputs:
                         n        :    Number of coherent integrations
                         timeInterval   :    Time of integration. If the parameter "n" is selected this one does not work
                         overlapping    :
                     """
                     self.__initime = None
                     self.__lastdatatime = 0
                     self.__buffer_spc = []
                     self.__buffer_cspc = []
                     self.__buffer_dc = 0
                     self.__profIndex = 0
                     self.__dataReady = False
                     self.__byTime = False
                     self.factor = factor
                     self.navg = avg
                     #self.ByLags = dataOut.ByLags   ###REDEFINIR
                     self.ByLags = False
                     self.maxProfilesInt = 0
                     self.__nChannels = dataOut.nChannels
                     if DPL != None:
                         self.DPL=DPL
                     else:
                         #self.DPL=dataOut.DPL    ###REDEFINIR
                         self.DPL=0
                     if n is None and timeInterval is None:
                         raise ValueError("n or timeInterval should be specified ...")
                     if n is not None:
                         self.n = int(n)
                     else:
                         self.__integrationtime = int(timeInterval)
                         self.n = None
                         self.__byTime = True
                     if minHei == None:
                         minHei = self.dataOut.heightList[0]
                     if maxHei == None:
                         maxHei = self.dataOut.heightList[-1]
                     if (minHei < self.dataOut.heightList[0]) or (minHei > maxHei):
                         print('minHei: %.2f is out of the heights range' % (minHei))
                         print('minHei is setting to %.2f' % (self.dataOut.heightList[0]))
                         minHei = self.dataOut.heightList[0]
                     if (maxHei > self.dataOut.heightList[-1]) or (maxHei < minHei):
                         print('maxHei: %.2f is out of the heights range' % (maxHei))
                         print('maxHei is setting to %.2f' % (self.dataOut.heightList[-1]))
                         maxHei = self.dataOut.heightList[-1]
                     ind_list1 = numpy.where(self.dataOut.heightList >= minHei)
                     ind_list2 = numpy.where(self.dataOut.heightList <= maxHei)
                     self.minHei_ind = ind_list1[0][0]
                     self.maxHei_ind = ind_list2[0][-1]
                 def putData(self, data_spc, data_cspc, data_dc):
                     """
                     Add a profile to the __buffer_spc and increase in one the __profileIndex
                     """
                     self.__buffer_spc.append(data_spc)
                     if self.__nChannels < 2:
                         self.__buffer_cspc = None
                     else:
                         self.__buffer_cspc.append(data_cspc)
                     if data_dc is None:
                         self.__buffer_dc = None
                     else:
                         self.__buffer_dc += data_dc
                     self.__profIndex += 1
                     return
                 def hildebrand_sekhon_Integration(self,sortdata,navg, factor):
                     #data debe estar ordenado
                     #sortdata = numpy.sort(data, axis=None)
                     #sortID=data.argsort()
                     lenOfData = len(sortdata)
                     nums_min = lenOfData*factor
                     if nums_min <= 5:
                         nums_min = 5
                     sump = 0.
                     sumq = 0.
                     j = 0
                     cont = 1
                     while((cont == 1)and(j < lenOfData)):
                         sump += sortdata[j]
                         sumq += sortdata[j]**2
                         if j > nums_min:
                             rtest = float(j)/(j-1) + 1.0/navg
                             if ((sumq*j) > (rtest*sump**2)):
                                 j = j - 1
                                 sump = sump - sortdata[j]
                                 sumq = sumq - sortdata[j]**2
                                 cont = 0
                         j += 1
                     #lnoise = sump / j
                     #print("H S done")
                     #return j,sortID
                     return j
                 def pushData(self):
                     """
                     Return the sum of the last profiles and the profiles used in the sum.
                     Affected:
                     self.__profileIndex
                     """
                     bufferH=None
                     buffer=None
                     buffer1=None
                     buffer_cspc=None
                     #print("aes: ", self.__buffer_cspc)
                     self.__buffer_spc=numpy.array(self.__buffer_spc)
                     if self.__nChannels > 1 :
                         self.__buffer_cspc=numpy.array(self.__buffer_cspc)
                     #print("FREQ_DC",self.__buffer_spc.shape,self.__buffer_cspc.shape)
                     freq_dc = int(self.__buffer_spc.shape[2] / 2)
                     #print("FREQ_DC",freq_dc,self.__buffer_spc.shape,self.nHeights)
                     self.dataOutliers = numpy.zeros((self.nChannels,self.nHeights)) # --> almacen de outliers
                     for k in range(self.minHei_ind,self.maxHei_ind):
                         if self.__nChannels > 1:
                             buffer_cspc=numpy.copy(self.__buffer_cspc[:,:,:,k])
                         outliers_IDs_cspc=[]
                         cspc_outliers_exist=False
                         for i in range(self.nChannels):#dataOut.nChannels):
                             buffer1=numpy.copy(self.__buffer_spc[:,i,:,k])
                             indexes=[]
                             #sortIDs=[]
                             outliers_IDs=[]
                             for j in range(self.nProfiles): #frecuencias en el tiempo
                                 # if i==0 and j==freq_dc: #NOT CONSIDERING DC PROFILE AT CHANNEL 0
                                 #     continue
                                 # if i==1 and j==0: #NOT CONSIDERING DC PROFILE AT CHANNEL 1
                                 #     continue
                                 buffer=buffer1[:,j]
                                 sortdata = numpy.sort(buffer, axis=None)
                                 sortID=buffer.argsort()
                                 index = _noise.hildebrand_sekhon2(sortdata,self.navg)
                                 #index,sortID=self.hildebrand_sekhon_Integration(buffer,1,self.factor)
                                 # fig,ax = plt.subplots()
                                 # ax.set_title(str(k)+" "+str(j))
                                 # x=range(len(sortdata))
                                 # ax.scatter(x,sortdata)
                                 # ax.axvline(index)
                                 # plt.show()
                                 indexes.append(index)
                                 #sortIDs.append(sortID)
                                 outliers_IDs=numpy.append(outliers_IDs,sortID[index:])
                             #print("Outliers: ",outliers_IDs)
                             outliers_IDs=numpy.array(outliers_IDs)
                             outliers_IDs=outliers_IDs.ravel()
                             outliers_IDs=numpy.unique(outliers_IDs)
                             outliers_IDs=outliers_IDs.astype(numpy.dtype('int64'))
                             indexes=numpy.array(indexes)
                             indexmin=numpy.min(indexes)
                             #print(indexmin,buffer1.shape[0], k)
                             # fig,ax = plt.subplots()
                             # ax.plot(sortdata)
                             # ax2 = ax.twinx()
                             # x=range(len(indexes))
                             # #plt.scatter(x,indexes)
                             # ax2.scatter(x,indexes)
                             # plt.show()
                             if indexmin != buffer1.shape[0]:
                                 if self.__nChannels > 1:
                                     cspc_outliers_exist= True
                                 lt=outliers_IDs
                                 #avg=numpy.mean(buffer1[[t for t in range(buffer1.shape[0]) if t not in lt],:],axis=0)
                                 for p in list(outliers_IDs):
                                     #buffer1[p,:]=avg
                                     buffer1[p,:] = numpy.NaN
                                 self.dataOutliers[i,k] = len(outliers_IDs)
                             self.__buffer_spc[:,i,:,k]=numpy.copy(buffer1)
                             if self.__nChannels > 1:
                                 outliers_IDs_cspc=numpy.append(outliers_IDs_cspc,outliers_IDs)
                         if self.__nChannels > 1:
                             outliers_IDs_cspc=outliers_IDs_cspc.astype(numpy.dtype('int64'))
                         if cspc_outliers_exist:
                             lt=outliers_IDs_cspc
                             #avg=numpy.mean(buffer_cspc[[t for t in range(buffer_cspc.shape[0]) if t not in lt],:],axis=0)
                             for p in list(outliers_IDs_cspc):
                                 #buffer_cspc[p,:]=avg
                                 buffer_cspc[p,:] = numpy.NaN
                         if self.__nChannels > 1:
                             self.__buffer_cspc[:,:,:,k]=numpy.copy(buffer_cspc)
                     nOutliers = len(outliers_IDs)
                     #print("Outliers  n: ",self.dataOutliers,nOutliers)
                     buffer=None
                     bufferH=None
                     buffer1=None
                     buffer_cspc=None
                     buffer=None
                     #data_spc = numpy.sum(self.__buffer_spc,axis=0)
                     data_spc = numpy.nansum(self.__buffer_spc,axis=0)
                     if self.__nChannels > 1:
                         #data_cspc = numpy.sum(self.__buffer_cspc,axis=0)
                         data_cspc = numpy.nansum(self.__buffer_cspc,axis=0)
                     else:
                         data_cspc = None
                     data_dc = self.__buffer_dc
                     #(CH, HEIGH)
                     self.maxProfilesInt = self.__profIndex - 1
                     n = self.__profIndex - self.dataOutliers # n becomes a matrix
                     self.__buffer_spc = []
                     self.__buffer_cspc = []
                     self.__buffer_dc = 0
                     self.__profIndex = 0
                     #print("cleaned ",data_cspc)
                     return data_spc, data_cspc, data_dc, n
                 def byProfiles(self, *args):
                     self.__dataReady = False
                     avgdata_spc = None
                     avgdata_cspc = None
                     avgdata_dc = None
                     self.putData(*args)
                     if self.__profIndex >= self.n:
                         avgdata_spc, avgdata_cspc, avgdata_dc, n = self.pushData()
                         self.n_ints = n
                         self.__dataReady = True
                     return avgdata_spc, avgdata_cspc, avgdata_dc
                 def byTime(self, datatime, *args):
                     self.__dataReady = False
                     avgdata_spc = None
                     avgdata_cspc = None
                     avgdata_dc = None
                     self.putData(*args)
                     if (datatime - self.__initime) >= self.__integrationtime:
                         avgdata_spc, avgdata_cspc, avgdata_dc, n = self.pushData()
                         self.n_ints = n
                         self.__dataReady = True
                     return avgdata_spc, avgdata_cspc, avgdata_dc
                 def integrate(self, datatime, *args):
                     if self.__profIndex == 0:
                         self.__initime = datatime
                     if self.__byTime:
                         avgdata_spc, avgdata_cspc, avgdata_dc = self.byTime(
                             datatime, *args)
                     else:
                         avgdata_spc, avgdata_cspc, avgdata_dc = self.byProfiles(*args)
                     if not self.__dataReady:
                         return None, None, None, None
                     #print("integrate", avgdata_cspc)
                     return self.__initime, avgdata_spc, avgdata_cspc, avgdata_dc
                 def run(self, dataOut, n=None, DPL = None,timeInterval=None, overlapping=False, minHei=None, maxHei=None, avg=1, factor=0.75):
                     self.dataOut = dataOut
                     if n == 1:
                         return self.dataOut
                     self.dataOut.processingHeaderObj.timeIncohInt = timeInterval
                     if dataOut.flagProfilesByRange:
                         self._flagProfilesByRange = True
                     if self.dataOut.nChannels == 1:
                         self.dataOut.data_cspc = None #si es un solo canal no vale la pena acumular DATOS
                     #print("IN spc:", self.dataOut.data_spc.shape, self.dataOut.data_cspc)
                     if not self.isConfig:
                         self.setup(self.dataOut, n, timeInterval, overlapping,DPL ,minHei, maxHei, avg, factor)
                         self.isConfig = True
                     if not self.ByLags:
                         self.nProfiles=self.dataOut.nProfiles
                         self.nChannels=self.dataOut.nChannels
                         self.nHeights=self.dataOut.nHeights
                         avgdatatime, avgdata_spc, avgdata_cspc, avgdata_dc = self.integrate(self.dataOut.utctime,
                                                                                             self.dataOut.data_spc,
                                                                                             self.dataOut.data_cspc,
                                                                                             self.dataOut.data_dc)
                     else:
                         self.nProfiles=self.dataOut.nProfiles
                         self.nChannels=self.dataOut.nChannels
                         self.nHeights=self.dataOut.nHeights
                         avgdatatime, avgdata_spc, avgdata_cspc, avgdata_dc = self.integrate(self.dataOut.utctime,
                                                                                             self.dataOut.dataLag_spc,
                                                                                             self.dataOut.dataLag_cspc,
                                                                                             self.dataOut.dataLag_dc)
                     self.dataOut.flagNoData = True
                     if self._flagProfilesByRange:
                         dataOut.flagProfilesByRange = True
                         self._nProfilesByRange += dataOut.nProfilesByRange
                     if self.__dataReady:
                         if not self.ByLags:
                             if self.nChannels == 1:
                                 #print("f int", avgdata_spc.shape)
                                 self.dataOut.data_spc = avgdata_spc
                                 self.dataOut.data_cspc = None
                             else:
                                 self.dataOut.data_spc = numpy.squeeze(avgdata_spc)
                                 self.dataOut.data_cspc = numpy.squeeze(avgdata_cspc)
                             self.dataOut.data_dc = avgdata_dc
                             self.dataOut.data_outlier = self.dataOutliers
                         else:
                             self.dataOut.dataLag_spc = avgdata_spc
                             self.dataOut.dataLag_cspc = avgdata_cspc
                             self.dataOut.dataLag_dc = avgdata_dc
                             self.dataOut.data_spc=self.dataOut.dataLag_spc[:,:,:,self.dataOut.LagPlot]
                             self.dataOut.data_cspc=self.dataOut.dataLag_cspc[:,:,:,self.dataOut.LagPlot]
                             self.dataOut.data_dc=self.dataOut.dataLag_dc[:,:,self.dataOut.LagPlot]
                         self.dataOut.nIncohInt *= self.n_ints
                         self.dataOut.utctime = avgdatatime
                         self.dataOut.flagNoData = False
                         dataOut.nProfilesByRange = self._nProfilesByRange
                         self._nProfilesByRange = numpy.zeros((1,len(dataOut.heightList)))
                         self._flagProfilesByRange = False
                     return self.dataOut
             class dopplerFlip(Operation):
                 def run(self, dataOut, chann = None):
                     # arreglo 1: (num_chan, num_profiles, num_heights)
                     self.dataOut = dataOut
                     # JULIA-oblicua, indice 2
                     # arreglo 2: (num_profiles, num_heights)
                     jspectra = self.dataOut.data_spc[chann]
                     jspectra_tmp = numpy.zeros(jspectra.shape)
                     num_profiles = jspectra.shape[0]
                     freq_dc = int(num_profiles / 2)
                     # Flip con for
                     for j in range(num_profiles):
                         jspectra_tmp[num_profiles-j-1]= jspectra[j]
                     # Intercambio perfil de DC con perfil inmediato anterior
                     jspectra_tmp[freq_dc-1]= jspectra[freq_dc-1]
                     jspectra_tmp[freq_dc]= jspectra[freq_dc]
                     # canal modificado es re-escrito en el arreglo de canales
                     self.dataOut.data_spc[chann] = jspectra_tmp
                     return self.dataOut

General Comments 0

Write
Preview

You need to be logged in to leave comments. Login now

No TODOs yet

	Site-wide shortcuts
/	Use quick search box
g h	Goto home page
g g	Goto my private gists page
g G	Goto my public gists page
g 0-9	Goto bookmarked items from 0-9
n r	New repository page
n g	New gist page

	Repositories
g s	Goto summary page
g c	Goto changelog page
g f	Goto files page
g F	Goto files page with file search activated
g p	Goto pull requests page
g o	Goto repository settings
g O	Goto repository access permissions settings
t s	Toggle sidebar on some pages