##// END OF EJS Templates
New Updates
rflores -
r1561:3446c343bac9
parent child
Show More
@@ -557,7 +557,7 class Plot(Operation):
557
557
558 self.sender_time = last_time
558 self.sender_time = last_time
559
559
560 attrs = ['titles', 'zmin', 'zmax', 'tag', 'ymin', 'ymax']
560 attrs = ['titles', 'zmin', 'zmax', 'tag', 'ymin', 'ymax', 'zlimits']
561 for attr in attrs:
561 for attr in attrs:
562 value = getattr(self, attr)
562 value = getattr(self, attr)
563 if value:
563 if value:
@@ -660,11 +660,12 class Plot(Operation):
660 self.poll.register(self.socket, zmq.POLLIN)
660 self.poll.register(self.socket, zmq.POLLIN)
661
661
662 tm = getattr(dataOut, self.attr_time)
662 tm = getattr(dataOut, self.attr_time)
663
664 if self.data and 'time' in self.xaxis and (tm - self.tmin) >= self.xrange*60*60:
663 if self.data and 'time' in self.xaxis and (tm - self.tmin) >= self.xrange*60*60:
665 self.save_time = tm
664 self.save_time = tm
666 self.__plot()
665 self.__plot()
667 self.tmin += self.xrange*60*60
666 #self.tmin += self.xrange*60*60 #Modified by R. Flores
667 self.tmin += 24*60*60 #Modified by R. Flores
668
668 self.data.setup()
669 self.data.setup()
669 self.clear_figures()
670 self.clear_figures()
670
671
@@ -675,6 +676,7 class Plot(Operation):
675 self.isPlotConfig = True
676 self.isPlotConfig = True
676 if self.xaxis == 'time':
677 if self.xaxis == 'time':
677 dt = self.getDateTime(tm)
678 dt = self.getDateTime(tm)
679
678 if self.xmin is None:
680 if self.xmin is None:
679 self.tmin = tm
681 self.tmin = tm
680 self.xmin = dt.hour
682 self.xmin = dt.hour
@@ -8,7 +8,7
8
8
9 import os
9 import os
10 import numpy
10 import numpy
11 import collections.abc
11 #import collections.abc
12
12
13 from schainpy.model.graphics.jroplot_base import Plot, plt, log
13 from schainpy.model.graphics.jroplot_base import Plot, plt, log
14
14
@@ -186,6 +186,7 class SpectraObliquePlot(Plot):
186 '''
186 '''
187 data['shift1'] = dataOut.Dop_EEJ_T1[0]
187 data['shift1'] = dataOut.Dop_EEJ_T1[0]
188 data['shift2'] = dataOut.Dop_EEJ_T2[0]
188 data['shift2'] = dataOut.Dop_EEJ_T2[0]
189 data['max_val_2'] = dataOut.Oblique_params[0,-1,:]
189 data['shift1_error'] = dataOut.Err_Dop_EEJ_T1[0]
190 data['shift1_error'] = dataOut.Err_Dop_EEJ_T1[0]
190 data['shift2_error'] = dataOut.Err_Dop_EEJ_T2[0]
191 data['shift2_error'] = dataOut.Err_Dop_EEJ_T2[0]
191
192
@@ -216,6 +217,7 class SpectraObliquePlot(Plot):
216 shift1 = data['shift1']
217 shift1 = data['shift1']
217 #print(shift1)
218 #print(shift1)
218 shift2 = data['shift2']
219 shift2 = data['shift2']
220 max_val_2 = data['max_val_2']
219 err1 = data['shift1_error']
221 err1 = data['shift1_error']
220 err2 = data['shift2_error']
222 err2 = data['shift2_error']
221 if ax.firsttime:
223 if ax.firsttime:
@@ -238,18 +240,22 class SpectraObliquePlot(Plot):
238
240
239 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)
241 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)
240 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)
242 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)
243 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)
244
241 #print("plotter1: ", self.ploterr1,shift1)
245 #print("plotter1: ", self.ploterr1,shift1)
242
246
243 else:
247 else:
244 #print("else plotter1: ", self.ploterr1,shift1)
248 #print("else plotter1: ", self.ploterr1,shift1)
245 self.ploterr1.remove()
249 self.ploterr1.remove()
246 self.ploterr2.remove()
250 self.ploterr2.remove()
251 self.ploterr3.remove()
247 ax.plt.set_array(z[n].T.ravel())
252 ax.plt.set_array(z[n].T.ravel())
248 if self.showprofile:
253 if self.showprofile:
249 ax.plt_profile.set_data(self.data['rti'][n][-1], y)
254 ax.plt_profile.set_data(self.data['rti'][n][-1], y)
250 ax.plt_noise.set_data(numpy.repeat(noise, len(y)), y)
255 ax.plt_noise.set_data(numpy.repeat(noise, len(y)), y)
251 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)
256 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)
252 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)
257 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)
258 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)
253
259
254 self.titles.append('CH {}: {:3.2f}dB'.format(n, noise))
260 self.titles.append('CH {}: {:3.2f}dB'.format(n, noise))
255
261
@@ -4,7 +4,7 import time
4 import math
4 import math
5 import datetime
5 import datetime
6 import numpy
6 import numpy
7 import collections.abc
7
8 from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation, MPDecorator #YONG
8 from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation, MPDecorator #YONG
9
9
10 from .jroplot_spectra import RTIPlot, NoisePlot
10 from .jroplot_spectra import RTIPlot, NoisePlot
@@ -749,9 +749,9 class FaradayAnglePlot(Plot):
749
749
750 data['angle'] = numpy.degrees(dataOut.phi)
750 data['angle'] = numpy.degrees(dataOut.phi)
751 #'''
751 #'''
752 print(dataOut.phi_uwrp)
752 #print(dataOut.phi_uwrp)
753 print(data['angle'])
753 #print(data['angle'])
754 exit(1)
754 #exit(1)
755 #'''
755 #'''
756 data['dphi'] = dataOut.dphi_uc*10
756 data['dphi'] = dataOut.dphi_uc*10
757 #print(dataOut.dphi)
757 #print(dataOut.dphi)
@@ -522,9 +522,8 class Reader(object):
522
522
523 def find_files(self, folders, ext, filefmt, startDate=None, endDate=None,
523 def find_files(self, folders, ext, filefmt, startDate=None, endDate=None,
524 expLabel='', last=False):
524 expLabel='', last=False):
525
526 for path in folders:
525 for path in folders:
527 files = glob.glob1(path, '*{}'.format(ext))
526 files = glob.glob1(path+'/'+expLabel, '*{}'.format(ext))
528 files.sort()
527 files.sort()
529 if last:
528 if last:
530 if files:
529 if files:
@@ -567,6 +566,7 class Reader(object):
567 if walk:
566 if walk:
568 folders = self.find_folders(
567 folders = self.find_folders(
569 path, startDate, endDate, folderfmt)
568 path, startDate, endDate, folderfmt)
569 #print("folders: ", folders)
570 else:
570 else:
571 folders = path.split(',')
571 folders = path.split(',')
572
572
@@ -928,7 +928,6 class JRODataReader(Reader):
928 self.lastUTTime = self.basicHeaderObj.utc
928 self.lastUTTime = self.basicHeaderObj.utc
929
929
930 self.flagDiscontinuousBlock = 0
930 self.flagDiscontinuousBlock = 0
931
932 if deltaTime > self.maxTimeStep:
931 if deltaTime > self.maxTimeStep:
933 self.flagDiscontinuousBlock = 1
932 self.flagDiscontinuousBlock = 1
934
933
@@ -384,6 +384,7 Inputs:
384
384
385 __attrs__ = ['path', 'oneDDict', 'ind2DList', 'twoDDict','metadata', 'format', 'blocks']
385 __attrs__ = ['path', 'oneDDict', 'ind2DList', 'twoDDict','metadata', 'format', 'blocks']
386 missing = -32767
386 missing = -32767
387 currentDay = None
387
388
388 def __init__(self):
389 def __init__(self):
389
390
@@ -608,12 +609,29 Inputs:
608 header.createHeader(**self.header)
609 header.createHeader(**self.header)
609 header.write()
610 header.write()
610
611
612 def timeFlag(self):
613 currentTime = self.dataOut.utctime
614 timeTuple = time.localtime(currentTime)
615 dataDay = timeTuple.tm_yday
616
617 if self.currentDay is None:
618 self.currentDay = dataDay
619 return False
620
621 #Si el dia es diferente
622 if dataDay != self.currentDay:
623 self.currentDay = dataDay
624 return True
625
626 else:
627 return False
628
611 def putData(self):
629 def putData(self):
612
630
613 if self.dataOut.flagNoData:
631 if self.dataOut.flagNoData:
614 return 0
632 return 0
615
633
616 if self.dataOut.flagDiscontinuousBlock or self.counter == self.blocks:
634 if self.dataOut.flagDiscontinuousBlock or self.counter == self.blocks or self.timeFlag():
617 if self.counter > 0:
635 if self.counter > 0:
618 self.setHeader()
636 self.setHeader()
619 self.counter = 0
637 self.counter = 0
@@ -75,7 +75,7 class SpectraReader(JRODataReader, ProcessingUnit):
75
75
76 self.pts2read_SelfSpectra = 0
76 self.pts2read_SelfSpectra = 0
77 self.pts2read_CrossSpectra = 0
77 self.pts2read_CrossSpectra = 0
78 self.pts2read_DCchannels = 0
78 self.pts2read_DCchannels = 0
79 self.ext = ".pdata"
79 self.ext = ".pdata"
80 self.optchar = "P"
80 self.optchar = "P"
81 self.basicHeaderObj = BasicHeader(LOCALTIME)
81 self.basicHeaderObj = BasicHeader(LOCALTIME)
@@ -162,7 +162,7 class SpectraReader(JRODataReader, ProcessingUnit):
162 Exceptions:
162 Exceptions:
163 Si un bloque leido no es un bloque valido
163 Si un bloque leido no es un bloque valido
164 """
164 """
165
165
166 fpointer = self.fp.tell()
166 fpointer = self.fp.tell()
167
167
168 spc = numpy.fromfile( self.fp, self.dtype[0], self.pts2read_SelfSpectra )
168 spc = numpy.fromfile( self.fp, self.dtype[0], self.pts2read_SelfSpectra )
@@ -364,7 +364,7 class SpectraWriter(JRODataWriter, Operation):
364 data.tofile(self.fp)
364 data.tofile(self.fp)
365
365
366 if self.data_cspc is not None:
366 if self.data_cspc is not None:
367
367
368 cspc = numpy.transpose( self.data_cspc, (0,2,1) )
368 cspc = numpy.transpose( self.data_cspc, (0,2,1) )
369 data = numpy.zeros( numpy.shape(cspc), self.dtype )
369 data = numpy.zeros( numpy.shape(cspc), self.dtype )
370 #print 'data.shape', self.shape_cspc_Buffer
370 #print 'data.shape', self.shape_cspc_Buffer
@@ -376,7 +376,7 class SpectraWriter(JRODataWriter, Operation):
376 data.tofile(self.fp)
376 data.tofile(self.fp)
377
377
378 if self.data_dc is not None:
378 if self.data_dc is not None:
379
379
380 dc = self.data_dc
380 dc = self.data_dc
381 data = numpy.zeros( numpy.shape(dc), self.dtype )
381 data = numpy.zeros( numpy.shape(dc), self.dtype )
382 data['real'] = dc.real
382 data['real'] = dc.real
@@ -524,4 +524,4 class SpectraWriter(JRODataWriter, Operation):
524
524
525 self.processingHeaderObj.processFlags = self.getProcessFlags()
525 self.processingHeaderObj.processFlags = self.getProcessFlags()
526
526
527 self.setBasicHeader() No newline at end of file
527 self.setBasicHeader()
@@ -1159,11 +1159,12 class Oblique_Gauss_Fit(Operation):
1159 freq_max = numpy.max(numpy.abs(freq))
1159 freq_max = numpy.max(numpy.abs(freq))
1160 spc_max = numpy.max(spc)
1160 spc_max = numpy.max(spc)
1161
1161
1162 from scipy.signal import medfilt
1162 #from scipy.signal import medfilt
1163 Nincoh = 20
1163 #Nincoh = 20
1164 Nincoh = 80
1164 #Nincoh = 80
1165 Nincoh = Nincoh
1165 Nincoh = Nincoh
1166 spcm = medfilt(spc,11)/numpy.sqrt(Nincoh)
1166 #spcm = medfilt(spc,11)/numpy.sqrt(Nincoh)
1167 spcm = spc/numpy.sqrt(Nincoh)
1167
1168
1168 # define a least squares function to optimize
1169 # define a least squares function to optimize
1169 def lsq_func(params):
1170 def lsq_func(params):
@@ -1174,7 +1175,9 class Oblique_Gauss_Fit(Operation):
1174 # bounds=([0,-numpy.inf,0,0,-numpy.inf,0,-numpy.inf,0],[numpy.inf,-200,numpy.inf,numpy.inf,0,numpy.inf,0,numpy.inf])
1175 # bounds=([0,-numpy.inf,0,0,-numpy.inf,0,-numpy.inf,0],[numpy.inf,-200,numpy.inf,numpy.inf,0,numpy.inf,0,numpy.inf])
1175 #print(a1,b1,c1,a2,b2,c2,k2,d)
1176 #print(a1,b1,c1,a2,b2,c2,k2,d)
1176 #bounds=([0,-numpy.inf,0,-numpy.inf,0,-400,0,0,0],[numpy.inf,-340,numpy.inf,0,numpy.inf,0,numpy.inf,numpy.inf,numpy.inf])
1177 #bounds=([0,-numpy.inf,0,-numpy.inf,0,-400,0,0,0],[numpy.inf,-340,numpy.inf,0,numpy.inf,0,numpy.inf,numpy.inf,numpy.inf])
1177 bounds=([0,-numpy.inf,0,-numpy.inf,0,-400,0,0,0],[numpy.inf,-140,numpy.inf,0,numpy.inf,0,numpy.inf,numpy.inf,numpy.inf])
1178 #bounds=([0,-numpy.inf,0,-numpy.inf,0,-400,0,0,0],[numpy.inf,-140,numpy.inf,0,numpy.inf,0,numpy.inf,numpy.inf,numpy.inf])
1179 bounds=([0,-numpy.inf,0,-5,0,-400,0,0,0],[numpy.inf,-300,numpy.inf,5,numpy.inf,0,numpy.inf,numpy.inf,numpy.inf])
1180
1178 #print(bounds)
1181 #print(bounds)
1179 #bounds=([0,-numpy.inf,0,0,-numpy.inf,0,0,0],[numpy.inf,-200,numpy.inf,numpy.inf,0,numpy.inf,numpy.inf,numpy.inf])
1182 #bounds=([0,-numpy.inf,0,0,-numpy.inf,0,0,0],[numpy.inf,-200,numpy.inf,numpy.inf,0,numpy.inf,numpy.inf,numpy.inf])
1180 params_scale = [spc_max,freq_max,freq_max,1,spc_max,freq_max,freq_max,1,spc_max]
1183 params_scale = [spc_max,freq_max,freq_max,1,spc_max,freq_max,freq_max,1,spc_max]
@@ -1278,9 +1281,6 class Oblique_Gauss_Fit(Operation):
1278 #print("before return")
1281 #print("before return")
1279 return A1f, B1f, C1f, A2f, B2f, C2f, Df, error
1282 return A1f, B1f, C1f, A2f, B2f, C2f, Df, error
1280
1283
1281
1282
1283
1284 def Double_Gauss_Double_Skew_fit_weight_bound_with_inputs(self, spc, freq, a1, b1, c1, a2, b2, c2, k2, d):
1284 def Double_Gauss_Double_Skew_fit_weight_bound_with_inputs(self, spc, freq, a1, b1, c1, a2, b2, c2, k2, d):
1285
1285
1286 from scipy.optimize import least_squares
1286 from scipy.optimize import least_squares
@@ -1469,7 +1469,7 class Oblique_Gauss_Fit(Operation):
1469 return A1f, B1f, C1f, Df, error
1469 return A1f, B1f, C1f, Df, error
1470
1470
1471
1471
1472 def run(self, dataOut, mode = 0, Hmin1 = None, Hmax1 = None, Hmin2 = None, Hmax2 = None):
1472 def run(self, dataOut, mode = 0, Hmin1 = None, Hmax1 = None, Hmin2 = None, Hmax2 = None, Dop = 'Shift'):
1473
1473
1474 pwcode = 1
1474 pwcode = 1
1475
1475
@@ -1589,16 +1589,18 class Oblique_Gauss_Fit(Operation):
1589
1589
1590 elif mode == 9: #Double Skewed Weighted Bounded no inputs
1590 elif mode == 9: #Double Skewed Weighted Bounded no inputs
1591 #if numpy.max(spc) <= 0:
1591 #if numpy.max(spc) <= 0:
1592 if x[numpy.argmax(spc)] <= 0:
1592 from scipy.signal import medfilt
1593 spcm = medfilt(spc,11)
1594 if x[numpy.argmax(spcm)] <= 0:
1593 #print("EEJ")
1595 #print("EEJ")
1594 dataOut.Oblique_params[0,0,hei],dataOut.Oblique_params[0,1,hei],dataOut.Oblique_params[0,2,hei],dataOut.Oblique_params[0,3,hei],dataOut.Oblique_params[0,4,hei],dataOut.Oblique_params[0,5,hei],dataOut.Oblique_params[0,6,hei],dataOut.Oblique_params[0,7,hei],dataOut.Oblique_params[0,8,hei],dataOut.Oblique_params[0,9,hei],dataOut.Oblique_params[0,10,hei],dataOut.Oblique_param_errors[0,:,hei] = self.Double_Gauss_Double_Skew_fit_weight_bound_no_inputs(spc,x,dataOut.nIncohInt)
1596 dataOut.Oblique_params[0,0,hei],dataOut.Oblique_params[0,1,hei],dataOut.Oblique_params[0,2,hei],dataOut.Oblique_params[0,3,hei],dataOut.Oblique_params[0,4,hei],dataOut.Oblique_params[0,5,hei],dataOut.Oblique_params[0,6,hei],dataOut.Oblique_params[0,7,hei],dataOut.Oblique_params[0,8,hei],dataOut.Oblique_params[0,9,hei],dataOut.Oblique_params[0,10,hei],dataOut.Oblique_param_errors[0,:,hei] = self.Double_Gauss_Double_Skew_fit_weight_bound_no_inputs(spcm,x,dataOut.nIncohInt)
1595 #if dataOut.Oblique_params[0,-2,hei] < -500 or dataOut.Oblique_params[0,-2,hei] > 500 or dataOut.Oblique_params[0,-1,hei] < -500 or dataOut.Oblique_params[0,-1,hei] > 500:
1597 #if dataOut.Oblique_params[0,-2,hei] < -500 or dataOut.Oblique_params[0,-2,hei] > 500 or dataOut.Oblique_params[0,-1,hei] < -500 or dataOut.Oblique_params[0,-1,hei] > 500:
1596 # dataOut.Oblique_params[0,:,hei] *= numpy.NAN
1598 # dataOut.Oblique_params[0,:,hei] *= numpy.NAN
1597 dataOut.dplr_2_u[0,0,hei] = dataOut.Oblique_params[0,10,hei]/numpy.sin(numpy.arccos(100./dataOut.heightList[hei]))
1599 dataOut.dplr_2_u[0,0,hei] = dataOut.Oblique_params[0,10,hei]/numpy.sin(numpy.arccos(100./dataOut.heightList[hei]))
1598
1600
1599 else:
1601 else:
1600 #print("CEEJ")
1602 #print("CEEJ")
1601 dataOut.Oblique_params[0,0,hei],dataOut.Oblique_params[0,1,hei],dataOut.Oblique_params[0,2,hei],dataOut.Oblique_params[0,3,hei],dataOut.Oblique_params[0,4,hei],dataOut.Oblique_params[0,5,hei],dataOut.Oblique_params[0,6,hei],dataOut.Oblique_params[0,7,hei],dataOut.Oblique_params[0,8,hei],dataOut.Oblique_params[0,9,hei],dataOut.Oblique_params[0,10,hei],dataOut.Oblique_param_errors[0,:,hei] = self.CEEJ_Skew_fit_weight_bound_no_inputs(spc,x,dataOut.nIncohInt)
1603 dataOut.Oblique_params[0,0,hei],dataOut.Oblique_params[0,1,hei],dataOut.Oblique_params[0,2,hei],dataOut.Oblique_params[0,3,hei],dataOut.Oblique_params[0,4,hei],dataOut.Oblique_params[0,5,hei],dataOut.Oblique_params[0,6,hei],dataOut.Oblique_params[0,7,hei],dataOut.Oblique_params[0,8,hei],dataOut.Oblique_params[0,9,hei],dataOut.Oblique_params[0,10,hei],dataOut.Oblique_param_errors[0,:,hei] = self.CEEJ_Skew_fit_weight_bound_no_inputs(spcm,x,dataOut.nIncohInt)
1602 #if dataOut.Oblique_params[0,-2,hei] < -500 or dataOut.Oblique_params[0,-2,hei] > 500 or dataOut.Oblique_params[0,-1,hei] < -500 or dataOut.Oblique_params[0,-1,hei] > 500:
1604 #if dataOut.Oblique_params[0,-2,hei] < -500 or dataOut.Oblique_params[0,-2,hei] > 500 or dataOut.Oblique_params[0,-1,hei] < -500 or dataOut.Oblique_params[0,-1,hei] > 500:
1603 # dataOut.Oblique_params[0,:,hei] *= numpy.NAN
1605 # dataOut.Oblique_params[0,:,hei] *= numpy.NAN
1604 dataOut.dplr_2_u[0,0,hei] = dataOut.Oblique_params[0,10,hei]/numpy.sin(numpy.arccos(100./dataOut.heightList[hei]))
1606 dataOut.dplr_2_u[0,0,hei] = dataOut.Oblique_params[0,10,hei]/numpy.sin(numpy.arccos(100./dataOut.heightList[hei]))
@@ -1672,10 +1674,18 class Oblique_Gauss_Fit(Operation):
1672 dataOut.lon=-76.87
1674 dataOut.lon=-76.87
1673
1675
1674 if mode == 9: #Double Skew Gaussian
1676 if mode == 9: #Double Skew Gaussian
1675 dataOut.Dop_EEJ_T1 = dataOut.Oblique_params[:,-2,:]
1677 #dataOut.Dop_EEJ_T1 = dataOut.Oblique_params[:,-2,:] #Pos[Max_value]
1678 #dataOut.Dop_EEJ_T1 = dataOut.Oblique_params[:,1,:] #Shift
1676 dataOut.Spec_W_T1 = dataOut.Oblique_params[:,2,:]
1679 dataOut.Spec_W_T1 = dataOut.Oblique_params[:,2,:]
1677 dataOut.Dop_EEJ_T2 = dataOut.Oblique_params[:,-1,:]
1680 #dataOut.Dop_EEJ_T2 = dataOut.Oblique_params[:,-1,:] #Pos[Max_value]
1681 #dataOut.Dop_EEJ_T2 = dataOut.Oblique_params[:,5,:] #Shift
1678 dataOut.Spec_W_T2 = dataOut.Oblique_params[:,6,:]
1682 dataOut.Spec_W_T2 = dataOut.Oblique_params[:,6,:]
1683 if Dop == 'Shift':
1684 dataOut.Dop_EEJ_T1 = dataOut.Oblique_params[:,1,:] #Shift
1685 dataOut.Dop_EEJ_T2 = dataOut.Oblique_params[:,5,:] #Shift
1686 elif Dop == 'Max':
1687 dataOut.Dop_EEJ_T1 = dataOut.Oblique_params[:,-2,:] #Pos[Max_value]
1688 dataOut.Dop_EEJ_T2 = dataOut.Oblique_params[:,-1,:] #Pos[Max_value]
1679
1689
1680 dataOut.Err_Dop_EEJ_T1 = dataOut.Oblique_param_errors[:,1,:] #En realidad este es el error?
1690 dataOut.Err_Dop_EEJ_T1 = dataOut.Oblique_param_errors[:,1,:] #En realidad este es el error?
1681 dataOut.Err_Spec_W_T1 = dataOut.Oblique_param_errors[:,2,:]
1691 dataOut.Err_Spec_W_T1 = dataOut.Oblique_param_errors[:,2,:]
@@ -1694,6 +1704,8 class Oblique_Gauss_Fit(Operation):
1694 dataOut.Err_Spec_W_T2 = dataOut.Oblique_param_errors[:,5,:]
1704 dataOut.Err_Spec_W_T2 = dataOut.Oblique_param_errors[:,5,:]
1695
1705
1696 dataOut.mode = mode
1706 dataOut.mode = mode
1707 dataOut.flagNoData = numpy.all(numpy.isnan(dataOut.Dop_EEJ_T1)) #Si todos los valores son NaN no se prosigue
1708 #dataOut.flagNoData = False #Descomentar solo para ploteo sino mantener comentado
1697
1709
1698 return dataOut
1710 return dataOut
1699
1711
@@ -6601,8 +6613,9 class IGRFModel(Operation):
6601 dataOut.ut=dataOut.bd_time.tm_hour+dataOut.bd_time.tm_min/60.0+dataOut.bd_time.tm_sec/3600.0
6613 dataOut.ut=dataOut.bd_time.tm_hour+dataOut.bd_time.tm_min/60.0+dataOut.bd_time.tm_sec/3600.0
6602
6614
6603 self.aux=0
6615 self.aux=0
6604
6616 dh = dataOut.heightList[1]-dataOut.heightList[0]
6605 dataOut.h=numpy.arange(0.0,15.0*dataOut.MAXNRANGENDT,15.0,dtype='float32')
6617 #dataOut.h=numpy.arange(0.0,15.0*dataOut.MAXNRANGENDT,15.0,dtype='float32')
6618 dataOut.h=numpy.arange(0.0,dh*dataOut.MAXNRANGENDT,dh,dtype='float32')
6606 dataOut.bfm=numpy.zeros(dataOut.MAXNRANGENDT,dtype='float32')
6619 dataOut.bfm=numpy.zeros(dataOut.MAXNRANGENDT,dtype='float32')
6607 dataOut.bfm=numpy.array(dataOut.bfm,order='F')
6620 dataOut.bfm=numpy.array(dataOut.bfm,order='F')
6608 dataOut.thb=numpy.zeros(dataOut.MAXNRANGENDT,dtype='float32')
6621 dataOut.thb=numpy.zeros(dataOut.MAXNRANGENDT,dtype='float32')
@@ -6627,7 +6640,7 class MergeProc(ProcessingUnit):
6627 #print(data_inputs)
6640 #print(data_inputs)
6628 #print("Run: ",self.dataOut.runNextUnit)
6641 #print("Run: ",self.dataOut.runNextUnit)
6629 #exit(1)
6642 #exit(1)
6630 #print(numpy.shape([getattr(data, attr_data) for data in data_inputs][1]))
6643 #print("a:", [getattr(data, attr_data) for data in data_inputs][1])
6631 #exit(1)
6644 #exit(1)
6632 if mode==0:
6645 if mode==0:
6633 data = numpy.concatenate([getattr(data, attr_data) for data in data_inputs])
6646 data = numpy.concatenate([getattr(data, attr_data) for data in data_inputs])
@@ -6747,3 +6760,13 class MergeProc(ProcessingUnit):
6747 #print(numpy.shape(self.dataOut.data_spc))
6760 #print(numpy.shape(self.dataOut.data_spc))
6748 #print("*************************GOOD*************************")
6761 #print("*************************GOOD*************************")
6749 #exit(1)
6762 #exit(1)
6763
6764 if mode==11: #MST ISR
6765 #data = numpy.concatenate([getattr(data, attr_data) for data in data_inputs],axis=1)
6766 #setattr(self.dataOut, attr_data, data)
6767 setattr(self.dataOut, 'ph2', [getattr(data, attr_data) for data in data_inputs][1])
6768 print("MST Density", numpy.shape(self.dataOut.ph2))
6769 print("cf MST: ", self.dataOut.cf)
6770 exit(1)
6771 self.dataOut.ph2 *= self.dataOut.cf
6772 self.dataOut.sdp2 *= self.dataOut.cf
@@ -945,6 +945,7 class SpectraAFCProc(ProcessingUnit):
945 #data = numpy.fft.ifft(data, axis=1, n = 32)
945 #data = numpy.fft.ifft(data, axis=1, n = 32)
946 #data = numpy.fft.fftshift( data, axes=(1,))
946 #data = numpy.fft.fftshift( data, axes=(1,))
947 #acf = numpy.abs(data)
947 #acf = numpy.abs(data)
948 #print("data", data[0,0,0])
948 acf = data[:,:16,:]
949 acf = data[:,:16,:]
949 #acf = data[:,16:,:]
950 #acf = data[:,16:,:]
950 #print("SUM: ",numpy.sum(acf))
951 #print("SUM: ",numpy.sum(acf))
@@ -216,6 +216,7 class SpectraLagProc(ProcessingUnit):
216 self.dataOut.LagPlot=LagPlot
216 self.dataOut.LagPlot=LagPlot
217
217
218 #print(self.dataIn.data.shape)
218 #print(self.dataIn.data.shape)
219 #exit(1)
219 '''
220 '''
220 try:
221 try:
221 print(self.dataIn.data.shape)
222 print(self.dataIn.data.shape)
@@ -243,6 +244,10 class SpectraLagProc(ProcessingUnit):
243
244
244 if not self.dataOut.ByLags:
245 if not self.dataOut.ByLags:
245 #self.dataOut.data = self.dataIn.data
246 #self.dataOut.data = self.dataIn.data
247 try:
248 self.dataOut.FlipChannels=self.dataIn.FlipChannels
249 except: pass
250 self.dataOut.TimeBlockSeconds=self.dataIn.TimeBlockSeconds
246 self.VoltageType(nFFTPoints,nProfiles,ippFactor,pairsList)
251 self.VoltageType(nFFTPoints,nProfiles,ippFactor,pairsList)
247 else:
252 else:
248 self.dataOut.nLags = nLags
253 self.dataOut.nLags = nLags
@@ -1810,7 +1815,7 class IntegrationFaradaySpectra2(Operation):
1810 '''
1815 '''
1811 #print(self.nHeights)
1816 #print(self.nHeights)
1812 #exit(1)
1817 #exit(1)
1813 for l in range(self.nLags):#dataOut.DPL):
1818 for l in range(self.nLags):#dataOut.DPL): #if DP --> nLags=11, elif HP --> nLags=16
1814 #breakFlag=False
1819 #breakFlag=False
1815 for k in range(7,self.nHeights):
1820 for k in range(7,self.nHeights):
1816 if self.__buffer_cspc is not None:
1821 if self.__buffer_cspc is not None:
@@ -1818,7 +1823,7 class IntegrationFaradaySpectra2(Operation):
1818 outliers_IDs_cspc=[]
1823 outliers_IDs_cspc=[]
1819 cspc_outliers_exist=False
1824 cspc_outliers_exist=False
1820 #indexmin_cspc=0
1825 #indexmin_cspc=0
1821 for i in range(2):
1826 for i in range(2): #Solo nos interesa los 2 primeros canales que son los canales con seΓ±al
1822 #for i in range(self.nChannels):#dataOut.nChannels):
1827 #for i in range(self.nChannels):#dataOut.nChannels):
1823 #if self.TrueLags:
1828 #if self.TrueLags:
1824 #print("HERE")
1829 #print("HERE")
@@ -2871,7 +2876,8 class IntegrationFaradaySpectraNoLags(Operation):
2871 self.__lastdatatime = 0
2876 self.__lastdatatime = 0
2872
2877
2873 self.__buffer_spc = []
2878 self.__buffer_spc = []
2874 self.__buffer_cspc = []
2879 #self.__buffer_cspc = []
2880 self.__buffer_cspc = None
2875 self.__buffer_dc = 0
2881 self.__buffer_dc = 0
2876
2882
2877 self.__profIndex = 0
2883 self.__profIndex = 0
@@ -2949,7 +2955,7 class IntegrationFaradaySpectraNoLags(Operation):
2949
2955
2950 return j,sortID
2956 return j,sortID
2951 #'''
2957 #'''
2952 def pushData(self):
2958 def pushData_original_09_11_22(self):
2953 """
2959 """
2954 Return the sum of the last profiles and the profiles used in the sum.
2960 Return the sum of the last profiles and the profiles used in the sum.
2955
2961
@@ -2963,16 +2969,20 class IntegrationFaradaySpectraNoLags(Operation):
2963 buffer1=None
2969 buffer1=None
2964 buffer_cspc=None
2970 buffer_cspc=None
2965 self.__buffer_spc=numpy.array(self.__buffer_spc)
2971 self.__buffer_spc=numpy.array(self.__buffer_spc)
2966 self.__buffer_cspc=numpy.array(self.__buffer_cspc)
2972 #self.__buffer_cspc=numpy.array(self.__buffer_cspc)
2973 if self.__buffer_cspc is not None:
2974 self.__buffer_cspc=numpy.array(self.__buffer_cspc)
2967 freq_dc = int(self.__buffer_spc.shape[2] / 2)
2975 freq_dc = int(self.__buffer_spc.shape[2] / 2)
2968 #print("FREQ_DC",freq_dc,self.__buffer_spc.shape,self.nHeights)
2976 #print("FREQ_DC",freq_dc,self.__buffer_spc.shape,self.nHeights)
2969 for k in range(7,self.nHeights):
2977 for k in range(7,self.nHeights):
2970 buffer_cspc=numpy.copy(self.__buffer_cspc[:,:,:,k])
2978 if self.__buffer_cspc is not None:
2971 outliers_IDs_cspc=[]
2979 buffer_cspc=numpy.copy(self.__buffer_cspc[:,:,:,k])
2972 cspc_outliers_exist=False
2980 outliers_IDs_cspc=[]
2973 for i in range(self.nChannels):#dataOut.nChannels):
2981 cspc_outliers_exist=False
2974
2982 #for i in range(self.nChannels):#dataOut.nChannels):
2983 for i in range(2):#dataOut.nChannels):
2975 buffer1=numpy.copy(self.__buffer_spc[:,i,:,k])
2984 buffer1=numpy.copy(self.__buffer_spc[:,i,:,k])
2985
2976 indexes=[]
2986 indexes=[]
2977 #sortIDs=[]
2987 #sortIDs=[]
2978 outliers_IDs=[]
2988 outliers_IDs=[]
@@ -2984,8 +2994,9 class IntegrationFaradaySpectraNoLags(Operation):
2984 continue
2994 continue
2985 buffer=buffer1[:,j]
2995 buffer=buffer1[:,j]
2986 #if k != 100:
2996 #if k != 100:
2987 index=int(_HS_algorithm.HS_algorithm(numpy.sort(buffer, axis=None),1))
2997 #index=int(_HS_algorithm.HS_algorithm(numpy.sort(buffer, axis=None),1))
2988 sortID = buffer.argsort()
2998 index,sortID=self.hildebrand_sekhon_Integration(buffer,1)
2999 #sortID = buffer.argsort()
2989 #else:
3000 #else:
2990 #index,sortID=self.hildebrand_sekhon_Integration(buffer,1)
3001 #index,sortID=self.hildebrand_sekhon_Integration(buffer,1)
2991 #if k == 100:
3002 #if k == 100:
@@ -3008,29 +3019,35 class IntegrationFaradaySpectraNoLags(Operation):
3008 cspc_outliers_exist=True
3019 cspc_outliers_exist=True
3009 ###sortdata=numpy.sort(buffer1,axis=0)
3020 ###sortdata=numpy.sort(buffer1,axis=0)
3010 ###avg2=numpy.mean(sortdata[:indexmin,:],axis=0)
3021 ###avg2=numpy.mean(sortdata[:indexmin,:],axis=0)
3011 lt=outliers_IDs
3012 avg=numpy.mean(buffer1[[t for t in range(buffer1.shape[0]) if t not in lt],:],axis=0)
3013
3022
3023 lt=outliers_IDs
3024 #print("buffer1: ", numpy.sum(buffer1))
3025 #print("outliers: ", buffer1[lt])
3026 #print("outliers_IDs: ", outliers_IDs)
3027 avg=numpy.nanmean(buffer1[[t for t in range(buffer1.shape[0]) if t not in lt],:],axis=0)
3028 #print("avg: ", avg)
3014 for p in list(outliers_IDs):
3029 for p in list(outliers_IDs):
3015 buffer1[p,:]=avg
3030 buffer1[p,:]=avg
3016
3031
3017 self.__buffer_spc[:,i,:,k]=numpy.copy(buffer1)
3032 self.__buffer_spc[:,i,:,k]=numpy.copy(buffer1)
3018 ###cspc IDs
3033 ###cspc IDs
3019 #indexmin_cspc+=indexmin_cspc
3034 #indexmin_cspc+=indexmin_cspc
3020 outliers_IDs_cspc=numpy.append(outliers_IDs_cspc,outliers_IDs)
3035 if self.__buffer_cspc is not None:
3036 outliers_IDs_cspc=numpy.append(outliers_IDs_cspc,outliers_IDs)
3021
3037
3022 #if not breakFlag:
3038 #if not breakFlag:
3023 outliers_IDs_cspc=outliers_IDs_cspc.astype(numpy.dtype('int64'))
3039 if self.__buffer_cspc is not None:
3024 if cspc_outliers_exist:
3040 outliers_IDs_cspc=outliers_IDs_cspc.astype(numpy.dtype('int64'))
3025 #sortdata=numpy.sort(buffer_cspc,axis=0)
3041 if cspc_outliers_exist:
3026 #avg=numpy.mean(sortdata[:indexmin_cpsc,:],axis=0)
3042 #sortdata=numpy.sort(buffer_cspc,axis=0)
3027 lt=outliers_IDs_cspc
3043 #avg=numpy.mean(sortdata[:indexmin_cpsc,:],axis=0)
3044 lt=outliers_IDs_cspc
3028
3045
3029 avg=numpy.mean(buffer_cspc[[t for t in range(buffer_cspc.shape[0]) if t not in lt],:],axis=0)
3046 avg=numpy.mean(buffer_cspc[[t for t in range(buffer_cspc.shape[0]) if t not in lt],:],axis=0)
3030 for p in list(outliers_IDs_cspc):
3047 for p in list(outliers_IDs_cspc):
3031 buffer_cspc[p,:]=avg
3048 buffer_cspc[p,:]=avg
3032
3049
3033 self.__buffer_cspc[:,:,:,k]=numpy.copy(buffer_cspc)
3050 self.__buffer_cspc[:,:,:,k]=numpy.copy(buffer_cspc)
3034 #else:
3051 #else:
3035 #break
3052 #break
3036
3053
@@ -3050,7 +3067,15 class IntegrationFaradaySpectraNoLags(Operation):
3050 #print(self.__buffer_spc[:,1,3,20,0])
3067 #print(self.__buffer_spc[:,1,3,20,0])
3051 #print(self.__buffer_spc[:,1,5,37,0])
3068 #print(self.__buffer_spc[:,1,5,37,0])
3052 data_spc = numpy.sum(self.__buffer_spc,axis=0)
3069 data_spc = numpy.sum(self.__buffer_spc,axis=0)
3053 data_cspc = numpy.sum(self.__buffer_cspc,axis=0)
3070 print("data_spc: ", data_spc[0,:,0])
3071 print("data_spc: ", data_spc[0,:,7])
3072 print("shape: ", numpy.shape(data_spc))
3073 #exit(1)
3074 #data_cspc = numpy.sum(self.__buffer_cspc,axis=0)
3075 if self.__buffer_cspc is not None:
3076 data_cspc = numpy.sum(self.__buffer_cspc,axis=0)
3077 else:
3078 data_cspc = None
3054
3079
3055 #print(numpy.shape(data_spc))
3080 #print(numpy.shape(data_spc))
3056 #data_spc[1,4,20,0]=numpy.nan
3081 #data_spc[1,4,20,0]=numpy.nan
@@ -3060,7 +3085,150 class IntegrationFaradaySpectraNoLags(Operation):
3060 n = self.__profIndex
3085 n = self.__profIndex
3061
3086
3062 self.__buffer_spc = []
3087 self.__buffer_spc = []
3063 self.__buffer_cspc = []
3088 #self.__buffer_cspc = []
3089 self.__buffer_cspc = None
3090 self.__buffer_dc = 0
3091 self.__profIndex = 0
3092
3093 return data_spc, data_cspc, data_dc, n
3094
3095 def pushData(self):
3096 """
3097 Return the sum of the last profiles and the profiles used in the sum.
3098
3099 Affected:
3100
3101 self.__profileIndex
3102
3103 """
3104 bufferH=None
3105 buffer=None
3106 buffer1=None
3107 buffer_cspc=None
3108 self.__buffer_spc=numpy.array(self.__buffer_spc)
3109 #self.__buffer_cspc=numpy.array(self.__buffer_cspc)
3110 if self.__buffer_cspc is not None:
3111 self.__buffer_cspc=numpy.array(self.__buffer_cspc)
3112 freq_dc = int(self.__buffer_spc.shape[2] / 2)
3113 #print("FREQ_DC",freq_dc,self.__buffer_spc.shape,self.nHeights)
3114 for k in range(7,self.nHeights):
3115 if self.__buffer_cspc is not None:
3116 buffer_cspc=numpy.copy(self.__buffer_cspc[:,:,:,k])
3117 outliers_IDs_cspc=[]
3118 cspc_outliers_exist=False
3119 #for i in range(self.nChannels):#dataOut.nChannels):
3120 for i in range(2):#dataOut.nChannels):
3121 buffer1=numpy.copy(self.__buffer_spc[:,i,:,k])
3122
3123 indexes=[]
3124 #sortIDs=[]
3125 outliers_IDs=[]
3126
3127 for j in range(self.nProfiles):
3128 if i==0 and j==freq_dc: #NOT CONSIDERING DC PROFILE AT CHANNEL 0
3129 continue
3130 if i==1 and j==0: #NOT CONSIDERING DC PROFILE AT CHANNEL 1
3131 continue
3132 buffer=buffer1[:,j]
3133 #if k != 100:
3134 index=int(_HS_algorithm.HS_algorithm(numpy.sort(buffer, axis=None),1))
3135 #index,sortID=self.hildebrand_sekhon_Integration(buffer,1)
3136 sortID = buffer.argsort()
3137 #else:
3138 #index,sortID=self.hildebrand_sekhon_Integration(buffer,1)
3139 #if k == 100:
3140 # print(k,index,sortID)
3141 # exit(1)
3142 #print("index: ", index)
3143 indexes.append(index)
3144 out_IDs = sortID[index:]
3145 avg=numpy.nanmean(buffer1[[t for t in range(buffer1.shape[0]) if t not in out_IDs],:],axis=0)
3146 #print("avg: ", avg)
3147 for p in list(out_IDs):
3148 buffer1[p]=avg
3149 #sortIDs.append(sortID)
3150 outliers_IDs=numpy.append(outliers_IDs,sortID[index:])
3151 #if k == 100:
3152 # exit(1)
3153 outliers_IDs=numpy.array(outliers_IDs)
3154 outliers_IDs=outliers_IDs.ravel()
3155 outliers_IDs=numpy.unique(outliers_IDs)
3156 outliers_IDs=outliers_IDs.astype(numpy.dtype('int64'))
3157 indexes=numpy.array(indexes)
3158 indexmin=numpy.min(indexes)
3159 '''
3160 if indexmin != buffer1.shape[0]:
3161 cspc_outliers_exist=True
3162 ###sortdata=numpy.sort(buffer1,axis=0)
3163 ###avg2=numpy.mean(sortdata[:indexmin,:],axis=0)
3164
3165 lt=outliers_IDs
3166 #print("buffer1: ", numpy.sum(buffer1))
3167 #print("outliers: ", buffer1[lt])
3168 #print("outliers_IDs: ", outliers_IDs)
3169 avg=numpy.nanmean(buffer1[[t for t in range(buffer1.shape[0]) if t not in lt],:],axis=0)
3170 #print("avg: ", avg)
3171 for p in list(outliers_IDs):
3172 buffer1[p,:]=avg
3173 '''
3174 self.__buffer_spc[:,i,:,k]=numpy.copy(buffer1)
3175 ###cspc IDs
3176 #indexmin_cspc+=indexmin_cspc
3177 if self.__buffer_cspc is not None:
3178 outliers_IDs_cspc=numpy.append(outliers_IDs_cspc,outliers_IDs)
3179
3180 #if not breakFlag:
3181 if self.__buffer_cspc is not None:
3182 outliers_IDs_cspc=outliers_IDs_cspc.astype(numpy.dtype('int64'))
3183 if cspc_outliers_exist:
3184 #sortdata=numpy.sort(buffer_cspc,axis=0)
3185 #avg=numpy.mean(sortdata[:indexmin_cpsc,:],axis=0)
3186 lt=outliers_IDs_cspc
3187
3188 avg=numpy.mean(buffer_cspc[[t for t in range(buffer_cspc.shape[0]) if t not in lt],:],axis=0)
3189 for p in list(outliers_IDs_cspc):
3190 buffer_cspc[p,:]=avg
3191
3192 self.__buffer_cspc[:,:,:,k]=numpy.copy(buffer_cspc)
3193 #else:
3194 #break
3195
3196
3197
3198
3199 buffer=None
3200 bufferH=None
3201 buffer1=None
3202 buffer_cspc=None
3203
3204 #print("cpsc",self.__buffer_cspc[:,0,0,0,0])
3205 #print(self.__profIndex)
3206 #exit()
3207
3208 buffer=None
3209 #print(self.__buffer_spc[:,1,3,20,0])
3210 #print(self.__buffer_spc[:,1,5,37,0])
3211 data_spc = numpy.sum(self.__buffer_spc,axis=0)
3212 print("data_spc: ", data_spc[0,:,0])
3213 print("data_spc: ", data_spc[0,:,7])
3214 print("shape: ", numpy.shape(data_spc))
3215 #exit(1)
3216 #data_cspc = numpy.sum(self.__buffer_cspc,axis=0)
3217 if self.__buffer_cspc is not None:
3218 data_cspc = numpy.sum(self.__buffer_cspc,axis=0)
3219 else:
3220 data_cspc = None
3221
3222 #print(numpy.shape(data_spc))
3223 #data_spc[1,4,20,0]=numpy.nan
3224
3225 #data_cspc = self.__buffer_cspc
3226 data_dc = self.__buffer_dc
3227 n = self.__profIndex
3228
3229 self.__buffer_spc = []
3230 #self.__buffer_cspc = []
3231 self.__buffer_cspc = None
3064 self.__buffer_dc = 0
3232 self.__buffer_dc = 0
3065 self.__profIndex = 0
3233 self.__profIndex = 0
3066
3234
@@ -3120,7 +3288,10 class IntegrationFaradaySpectraNoLags(Operation):
3120 return dataOut
3288 return dataOut
3121
3289
3122 dataOut.flagNoData = True
3290 dataOut.flagNoData = True
3123
3291 #print(numpy.shape(dataOut.data_spc))
3292 #print(numpy.shape(dataOut.data_cspc))
3293 #exit(1)
3294 dataOut.data_cspc = None
3124 if not self.isConfig:
3295 if not self.isConfig:
3125 self.setup(dataOut, n, timeInterval, overlapping,DPL )
3296 self.setup(dataOut, n, timeInterval, overlapping,DPL )
3126 self.isConfig = True
3297 self.isConfig = True
@@ -3751,7 +3922,7 class SnrFaraday(Operation):
3751
3922
3752 return dataOut
3923 return dataOut
3753
3924
3754 class SpectraDataToFaraday(Operation):
3925 class SpectraDataToFaraday_07_11_2022(Operation):
3755 '''
3926 '''
3756 Written by R. Flores
3927 Written by R. Flores
3757 '''
3928 '''
@@ -3990,23 +4161,34 class SpectraDataToFaraday(Operation):
3990 data_to_remov_eej = dataOut.dataLag_spc[:,:,:,0]
4161 data_to_remov_eej = dataOut.dataLag_spc[:,:,:,0]
3991
4162
3992 self.normFactor(dataOut)
4163 self.normFactor(dataOut)
3993
4164 #print(dataOut.NDP)
3994 dataOut.NDP=dataOut.nHeights
4165 dataOut.NDP=dataOut.nHeights
3995 dataOut.NR=len(dataOut.channelList)
4166 dataOut.NR=len(dataOut.channelList)
3996 dataOut.DH=dataOut.heightList[1]-dataOut.heightList[0]
4167 dataOut.DH=dataOut.heightList[1]-dataOut.heightList[0]
3997 dataOut.H0=int(dataOut.heightList[0])
4168 dataOut.H0=int(dataOut.heightList[0])
3998
4169
3999 self.ConvertData(dataOut)
4170 self.ConvertData(dataOut)
4000
4171 #print(dataOut.NDP)
4172 #exit(1)
4001 dataOut.NAVG=16#dataOut.rnint2[0] #CHECK THIS!
4173 dataOut.NAVG=16#dataOut.rnint2[0] #CHECK THIS!
4002 if hasattr(dataOut, 'NRANGE'):
4174 if hasattr(dataOut, 'NRANGE'):
4003 dataOut.MAXNRANGENDT = max(dataOut.NRANGE,dataOut.NDT)
4175 dataOut.MAXNRANGENDT = max(dataOut.NRANGE,dataOut.NDT)
4004 else:
4176 else:
4005 dataOut.MAXNRANGENDT = dataOut.NDP
4177 dataOut.MAXNRANGENDT = dataOut.NDP
4006
4178
4179
4180 #if hasattr(dataOut, 'HP'):
4181 #pass
4182 #else:
4183 self.noise(dataOut)
4184
4185 '''
4007 if not hasattr(dataOut, 'tnoise'):
4186 if not hasattr(dataOut, 'tnoise'):
4187 print("noise")
4008 self.noise(dataOut)
4188 self.noise(dataOut)
4009
4189 else:
4190 delattr(dataOut, 'tnoise')
4191 '''
4010 #dataOut.pan = numpy.mean(dataOut.pan)
4192 #dataOut.pan = numpy.mean(dataOut.pan)
4011 #dataOut.pbn = numpy.mean(dataOut.pbn)
4193 #dataOut.pbn = numpy.mean(dataOut.pbn)
4012 #print(dataOut.pan)
4194 #print(dataOut.pan)
@@ -4023,7 +4205,197 class SpectraDataToFaraday(Operation):
4023 #exit(1)
4205 #exit(1)
4024 return dataOut
4206 return dataOut
4025
4207
4208 class SpectraDataToFaraday(Operation):
4209 """Operation to use spectra data in Faraday processing.
4210
4211 Parameters:
4212 -----------
4213 nint : int
4214 Number of integrations.
4215
4216 Example
4217 --------
4218
4219 op = proc_unit.addOperation(name='SpectraDataToFaraday', optype='other')
4220
4221 """
4222
4223 def __init__(self, **kwargs):
4224
4225 Operation.__init__(self, **kwargs)
4226
4227 self.dataLag_spc=None
4228 self.dataLag_cspc=None
4229 self.dataLag_dc=None
4230
4231 def noise_original(self,dataOut):
4232
4233 dataOut.noise_lag = numpy.zeros((dataOut.nChannels,dataOut.DPL),'float32')
4234 #print("Lags")
4235 '''
4236 for lag in range(dataOut.DPL):
4237 #print(lag)
4238 dataOut.data_spc = dataOut.dataLag_spc[:,:,:,lag]
4239 dataOut.noise_lag[:,lag] = dataOut.getNoise(ymin_index=46)
4240 #dataOut.noise_lag[:,lag] = dataOut.getNoise(ymin_index=33,ymax_index=46)
4241 '''
4242 #print(dataOut.NDP)
4243 #exit(1)
4244 #Channel B
4245 for lag in range(dataOut.DPL):
4246 #print(lag)
4247 dataOut.data_spc = dataOut.dataLag_spc[:,:,:,lag]
4248 max_hei_id = dataOut.NDP - 2*lag
4249 #if lag < 6:
4250 dataOut.noise_lag[1,lag] = dataOut.getNoise(ymin_index=53,ymax_index=max_hei_id)[1]
4251 #else:
4252 #dataOut.noise_lag[1,lag] = numpy.mean(dataOut.noise_lag[1,:6])
4253 #dataOut.noise_lag[:,lag] = dataOut.getNoise(ymin_index=33,ymax_index=46)
4254 #Channel A
4255 for lag in range(dataOut.DPL):
4256 #print(lag)
4257 dataOut.data_spc = dataOut.dataLag_spc[:,:,:,lag]
4258 dataOut.noise_lag[0,lag] = dataOut.getNoise(ymin_index=53)[0]
4026
4259
4260 nanindex = numpy.argwhere(numpy.isnan(numpy.log10(dataOut.noise_lag[1,:])))
4261 i1 = nanindex[0][0]
4262 dataOut.noise_lag[1,i1:] = numpy.mean(dataOut.noise_lag[1,:i1]) #El ruido de lags contaminados se
4263 #determina a partir del promedio del
4264 #ruido de los lags limpios
4265 '''
4266 dataOut.noise_lag[1,:] = dataOut.noise_lag[1,0] #El ruido de los lags diferentes de cero para
4267 #el canal B es contaminado por el Tx y EEJ
4268 #del siguiente perfil, por ello se asigna el ruido
4269 #del lag 0 a todos los lags
4270 '''
4271 #print("Noise lag: ", 10*numpy.log10(dataOut.noise_lag/dataOut.normFactor))
4272 #exit(1)
4273 '''
4274 dataOut.tnoise = dataOut.getNoise(ymin_index=46)
4275 dataOut.tnoise /= float(dataOut.nProfiles*dataOut.nIncohInt)
4276 dataOut.pan = dataOut.tnoise[0]
4277 dataOut.pbn = dataOut.tnoise[1]
4278 '''
4279
4280 dataOut.tnoise = dataOut.noise_lag/float(dataOut.nProfiles*dataOut.nIncohInt)
4281 #dataOut.tnoise /= float(dataOut.nProfiles*dataOut.nIncohInt)
4282 dataOut.pan = dataOut.tnoise[0]
4283 dataOut.pbn = dataOut.tnoise[1]
4284
4285 def noise(self,dataOut):
4286
4287 dataOut.noise_lag = numpy.zeros((dataOut.nChannels),'float32')
4288 #print("Lags")
4289 '''
4290 for lag in range(dataOut.DPL):
4291 #print(lag)
4292 dataOut.data_spc = dataOut.dataLag_spc[:,:,:,lag]
4293 dataOut.noise_lag[:,lag] = dataOut.getNoise(ymin_index=46)
4294 #dataOut.noise_lag[:,lag] = dataOut.getNoise(ymin_index=33,ymax_index=46)
4295 '''
4296 #print(dataOut.NDP)
4297 #exit(1)
4298 #Channel B
4299
4300 #print(lag)
4301 dataOut.data_spc = dataOut.dataLag_spc[:,:,:,0]
4302 max_hei_id = dataOut.NDP - 2*0
4303 #if lag < 6:
4304 #dataOut.noise_lag[1] = dataOut.getNoise(ymin_index=80,ymax_index=106)[1]
4305 dataOut.noise_lag[1] = dataOut.getNoise()[1]
4306 #else:
4307 #dataOut.noise_lag[1,lag] = numpy.mean(dataOut.noise_lag[1,:6])
4308 #dataOut.noise_lag[:,lag] = dataOut.getNoise(ymin_index=33,ymax_index=46)
4309 #Channel A
4310
4311 #print(lag)
4312 dataOut.data_spc = dataOut.dataLag_spc[:,:,:,0]
4313 dataOut.noise_lag[0] = dataOut.getNoise()[0]
4314
4315 dataOut.tnoise = dataOut.noise_lag/float(dataOut.nProfiles*dataOut.nIncohInt)
4316 #dataOut.tnoise /= float(dataOut.nProfiles*dataOut.nIncohInt)
4317 dataOut.pan = dataOut.tnoise[0]#*.95
4318 dataOut.pbn = dataOut.tnoise[1]#*.95
4319
4320 def ConvertData(self,dataOut):
4321
4322 dataOut.TimeBlockSeconds_for_dp_power=dataOut.utctime
4323 dataOut.bd_time=gmtime(dataOut.TimeBlockSeconds_for_dp_power)
4324 dataOut.year=dataOut.bd_time.tm_year+(dataOut.bd_time.tm_yday-1)/364.0
4325 dataOut.ut_Faraday=dataOut.bd_time.tm_hour+dataOut.bd_time.tm_min/60.0+dataOut.bd_time.tm_sec/3600.0
4326
4327
4328 tmpx=numpy.zeros((dataOut.nHeights,dataOut.DPL,2),'float32')
4329 tmpx_a2=numpy.zeros((dataOut.nHeights,dataOut.DPL,2),'float32')
4330 tmpx_b2=numpy.zeros((dataOut.nHeights,dataOut.DPL,2),'float32')
4331 tmpx_abr=numpy.zeros((dataOut.nHeights,dataOut.DPL,2),'float32')
4332 tmpx_abi=numpy.zeros((dataOut.nHeights,dataOut.DPL,2),'float32')
4333 dataOut.kabxys_integrated=[tmpx,tmpx,tmpx,tmpx,tmpx_a2,tmpx,tmpx_b2,tmpx,tmpx_abr,tmpx,tmpx_abi,tmpx,tmpx,tmpx]
4334
4335 dataOut.rnint2=numpy.zeros(dataOut.DPL,'float32')
4336 for l in range(dataOut.DPL):
4337 dataOut.rnint2[l]=1.0/(dataOut.nIncohInt*dataOut.nProfiles)#*dataOut.nProfiles
4338
4339
4340
4341 self.dataLag_spc=(dataOut.dataLag_spc.sum(axis=1))*(dataOut.rnint2[0]/dataOut.nProfiles)
4342 self.dataLag_cspc=(dataOut.dataLag_cspc.sum(axis=1))*(dataOut.rnint2[0]/dataOut.nProfiles)
4343 #self.dataLag_dc=dataOut.dataLag_dc.sum(axis=1)/dataOut.rnint2[0]
4344
4345
4346 dataOut.kabxys_integrated[4][:,:,0]=self.dataLag_spc[0,:,:].real
4347 #dataOut.kabxys_integrated[5][:,:,0]+=self.dataLag_spc[0,:,:].imag
4348 dataOut.kabxys_integrated[6][:,:,0]=self.dataLag_spc[1,:,:].real
4349 #dataOut.kabxys_integrated[7][:,:,0]+=self.dataLag_spc[1,:,:].imag
4350
4351 dataOut.kabxys_integrated[8][:,:,0]=self.dataLag_cspc[0,:,:].real
4352 dataOut.kabxys_integrated[10][:,:,0]=self.dataLag_cspc[0,:,:].imag
4353
4354 '''
4355 print(dataOut.kabxys_integrated[4][:,0,0])
4356 print(dataOut.kabxys_integrated[6][:,0,0])
4357 print("times 12")
4358 print(dataOut.kabxys_integrated[4][:,0,0]*dataOut.nProfiles)
4359 print(dataOut.kabxys_integrated[6][:,0,0]*dataOut.nProfiles)
4360 print(dataOut.rnint2[0])
4361 input()
4362 '''
4363
4364
4365
4366
4367
4368
4369 def run(self,dataOut):
4370
4371 dataOut.paramInterval=0#int(dataOut.nint*dataOut.header[7][0]*2 )
4372 dataOut.lat=-11.95
4373 dataOut.lon=-76.87
4374
4375 dataOut.NDP=dataOut.nHeights
4376 dataOut.NR=len(dataOut.channelList)
4377 dataOut.DH=dataOut.heightList[1]-dataOut.heightList[0]
4378 dataOut.H0=int(dataOut.heightList[0])
4379 #print(dataOut.data_spc.shape)
4380 dataOut.dataLag_spc = numpy.stack((dataOut.data_spc, dataOut.data_spc), axis=-1)
4381 dataOut.dataLag_cspc = numpy.stack((dataOut.data_cspc, dataOut.data_cspc), axis=-1)
4382 #print(dataOut.dataLag_spc.shape)
4383 dataOut.DPL = numpy.shape(dataOut.dataLag_spc)[-1]
4384 #exit(1)
4385 self.ConvertData(dataOut)
4386 self.noise(dataOut)
4387 dataOut.NAVG=16#dataOut.rnint2[0] #CHECK THIS!
4388 dataOut.MAXNRANGENDT=dataOut.NDP
4389 '''
4390 print(dataOut.kabxys_integrated[4][:,0,0])
4391 import matplotlib.pyplot as plt
4392 plt.plot(dataOut.kabxys_integrated[4][:,0,0],dataOut.heightList)
4393 plt.axvline(dataOut.pan)
4394 plt.xlim(0,1.e3)
4395 plt.show()
4396 '''
4397 dataOut.DPL = 1
4398 return dataOut
4027
4399
4028 class SpectraDataToHybrid(SpectraDataToFaraday):
4400 class SpectraDataToHybrid(SpectraDataToFaraday):
4029 '''
4401 '''
@@ -4352,7 +4724,7 class SpectraDataToHybrid_V2(SpectraDataToFaraday):
4352 dataOut.pbn_LP=dataOut.tnoise[1]/float(dataOut.nProfiles_LP*dataOut.nIncohInt_LP)
4724 dataOut.pbn_LP=dataOut.tnoise[1]/float(dataOut.nProfiles_LP*dataOut.nIncohInt_LP)
4353
4725
4354 def ConvertDataLP(self,dataOut):
4726 def ConvertDataLP(self,dataOut):
4355
4727 #print(numpy.shape(dataOut.data_acf))
4356 #print(dataOut.dataLag_spc[:,:,:,1]/dataOut.data_spc)
4728 #print(dataOut.dataLag_spc[:,:,:,1]/dataOut.data_spc)
4357 #exit(1)
4729 #exit(1)
4358 self.normfactor_LP=1.0/(dataOut.nIncohInt_LP*dataOut.nProfiles_LP)#*dataOut.nProfiles
4730 self.normfactor_LP=1.0/(dataOut.nIncohInt_LP*dataOut.nProfiles_LP)#*dataOut.nProfiles
@@ -4360,10 +4732,10 class SpectraDataToHybrid_V2(SpectraDataToFaraday):
4360 #print("Power: ",numpy.mean(dataOut.dataLag_spc_LP[0,:,100]))
4732 #print("Power: ",numpy.mean(dataOut.dataLag_spc_LP[0,:,100]))
4361 #buffer = dataOut.data_acf*(1./(normfactor*dataOut.nProfiles_LP))
4733 #buffer = dataOut.data_acf*(1./(normfactor*dataOut.nProfiles_LP))
4362 #buffer = dataOut.data_acf*(1./(normfactor))
4734 #buffer = dataOut.data_acf*(1./(normfactor))
4363 buffer = dataOut.data_acf#*(self.normfactor_LP)
4735 buffer = dataOut.data_acf#*(self.normfactor_LP) #nChannels x nProfiles (nLags) x nHeights
4364 #print("acf: ",numpy.sum(buffer))
4736 #print("acf: ",numpy.sum(buffer))
4365
4737
4366 dataOut.output_LP_integrated = numpy.transpose(buffer,(1,2,0))
4738 dataOut.output_LP_integrated = numpy.transpose(buffer,(1,2,0)) #nProfiles (nLags) x nHeights x nChannels
4367
4739
4368 def normFactor(self,dataOut):
4740 def normFactor(self,dataOut):
4369 dataOut.rnint2=numpy.zeros(dataOut.DPL,'float32')
4741 dataOut.rnint2=numpy.zeros(dataOut.DPL,'float32')
@@ -4394,7 +4766,9 class SpectraDataToHybrid_V2(SpectraDataToFaraday):
4394 #dataOut.output_LP_integrated[:,:,3] *= float(dataOut.NSCAN/22)#(dataOut.nNoiseProfiles) #Corrects the zero padding
4766 #dataOut.output_LP_integrated[:,:,3] *= float(dataOut.NSCAN/22)#(dataOut.nNoiseProfiles) #Corrects the zero padding
4395
4767
4396 dataOut.nis=dataOut.NSCAN*dataOut.nIncohInt_LP*10
4768 dataOut.nis=dataOut.NSCAN*dataOut.nIncohInt_LP*10
4769 #print("nis/10: ", dataOut.NSCAN,dataOut.nIncohInt_LP,dataOut.nProfiles_LP)
4397 dataOut.nis=dataOut.NSCAN*dataOut.nIncohInt_LP*dataOut.nProfiles_LP*10
4770 dataOut.nis=dataOut.NSCAN*dataOut.nIncohInt_LP*dataOut.nProfiles_LP*10
4771 dataOut.nis=dataOut.nIncohInt_LP*dataOut.nProfiles_LP*10 #Removemos NSCAN debido a que estΓ‘ incluido en nProfiles_LP
4398
4772
4399 self.ConvertData(dataOut)
4773 self.ConvertData(dataOut)
4400
4774
This diff has been collapsed as it changes many lines, (988 lines changed) Show them Hide them
@@ -2543,12 +2543,11 class CleanCohEchoes(Operation):
2543 dataOut.flagSpreadF = True
2543 dataOut.flagSpreadF = True
2544
2544
2545 #Removing echoes greater than 35 dB
2545 #Removing echoes greater than 35 dB
2546 if isinstance(dataOut.pbn, collections.abc.Sequence):
2546 if hasattr(dataOut.pbn, "__len__"):
2547 maxdB = 10*numpy.log10(dataOut.pbn[0]) + 10 #Lag 0 NOise
2547 maxdB = 10*numpy.log10(dataOut.pbn[0]) + 10 #Lag 0 Noise
2548 else:
2548 else:
2549 maxdB = 10*numpy.log10(dataOut.pbn) + 10 #Lag 0 NOise
2549 maxdB = 10*numpy.log10(dataOut.pbn) + 10
2550 #maxdB = 35 #DEBERÍA SER NOISE+ALGO!!!!!!!!!!!!!!!!!!!!!!
2550
2551 #print("noise: ",maxdB - 10)
2552 #print(dataOut.kabxys_integrated[6][:,0,0])
2551 #print(dataOut.kabxys_integrated[6][:,0,0])
2553 data = numpy.copy(10*numpy.log10(dataOut.kabxys_integrated[6][:,0,0])) #Lag0 ChB
2552 data = numpy.copy(10*numpy.log10(dataOut.kabxys_integrated[6][:,0,0])) #Lag0 ChB
2554 #print("data: ",data)
2553 #print("data: ",data)
@@ -3052,8 +3051,8 class NoisePower(Operation):
3052 dataOut.pnoise[i]=dataOut.pnoise[i]/dataOut.DPL
3051 dataOut.pnoise[i]=dataOut.pnoise[i]/dataOut.DPL
3053
3052
3054
3053
3055 dataOut.pan=1.0*dataOut.pnoise[0] # weights could change
3054 dataOut.pan=.8*dataOut.pnoise[0] # weights could change
3056 dataOut.pbn=1.0*dataOut.pnoise[1] # weights could change
3055 dataOut.pbn=.8*dataOut.pnoise[1] # weights could change
3057 '''
3056 '''
3058 print("pan: ",dataOut.pan)
3057 print("pan: ",dataOut.pan)
3059 print("pbn: ",dataOut.pbn)
3058 print("pbn: ",dataOut.pbn)
@@ -3107,11 +3106,12 class DoublePulseACFs(Operation):
3107 panrm=numpy.zeros((dataOut.NDP,dataOut.DPL), dtype=float)
3106 panrm=numpy.zeros((dataOut.NDP,dataOut.DPL), dtype=float)
3108
3107
3109 id = numpy.where(dataOut.heightList>700)[0]
3108 id = numpy.where(dataOut.heightList>700)[0]
3110
3109 #print("kabxys: ", numpy.shape(dataOut.kabxys_integrated))
3111 for i in range(dataOut.NDP):
3110 for i in range(dataOut.NDP):
3112 for j in range(dataOut.DPL):
3111 for j in range(dataOut.DPL):
3113 ################# Total power
3112 ################# Total power
3114 pa=numpy.abs(dataOut.kabxys_integrated[4][i,j,0]+dataOut.kabxys_integrated[5][i,j,0])
3113 pa=numpy.abs(dataOut.kabxys_integrated[4][i,j,0]+dataOut.kabxys_integrated[5][i,j,0])
3114 #print("pa::",pa)
3115 pb=numpy.abs(dataOut.kabxys_integrated[6][i,j,0]+dataOut.kabxys_integrated[7][i,j,0])
3115 pb=numpy.abs(dataOut.kabxys_integrated[6][i,j,0]+dataOut.kabxys_integrated[7][i,j,0])
3116 st4=pa*pb
3116 st4=pa*pb
3117
3117
@@ -3200,17 +3200,19 class DoublePulseACFs(Operation):
3200 exit(1)
3200 exit(1)
3201 '''
3201 '''
3202 #print(pa)
3202 #print(pa)
3203 #print("pa: ", numpy.shape(pa))
3204 #print(numpy.shape(dataOut.heightList))
3203 '''
3205 '''
3204 import matplotlib.pyplot as plt
3206 import matplotlib.pyplot as plt
3205 #plt.plot(dataOut.p[:,-1],dataOut.heightList)
3207 plt.plot(dataOut.p[:,-1],dataOut.heightList)
3206 plt.plot(pa/dataOut.pan-1.,dataOut.heightList)
3208 #plt.plot(pa/dataOut.pan-1.,dataOut.heightList)
3207 plt.plot(pb/dataOut.pbn-1.,dataOut.heightList)
3209 #plt.plot(pb/dataOut.pbn-1.,dataOut.heightList)
3208 plt.grid()
3210 plt.grid()
3209 #plt.xlim(0,1e5)
3211 plt.xlim(0,1e5)
3210 plt.show()
3212 plt.show()
3211 #print("p: ",dataOut.p[33,:])
3213 #print("p: ",dataOut.p[33,:])
3212 #exit(1)
3214 #exit(1)
3213 '''
3215 #'''
3214 #print(numpy.sum(dataOut.rhor))
3216 #print(numpy.sum(dataOut.rhor))
3215 #exit(1)
3217 #exit(1)
3216 return dataOut
3218 return dataOut
@@ -3346,7 +3348,7 class DoublePulseACFs_PerLag(Operation):
3346 exit(1)
3348 exit(1)
3347 '''
3349 '''
3348 #print(pa)
3350 #print(pa)
3349 '''
3351 #'''
3350 import matplotlib.pyplot as plt
3352 import matplotlib.pyplot as plt
3351 #plt.plot(dataOut.p[:,-1],dataOut.heightList)
3353 #plt.plot(dataOut.p[:,-1],dataOut.heightList)
3352 plt.plot(pa/dataOut.pan-1.,dataOut.heightList)
3354 plt.plot(pa/dataOut.pan-1.,dataOut.heightList)
@@ -3356,7 +3358,7 class DoublePulseACFs_PerLag(Operation):
3356 plt.show()
3358 plt.show()
3357 #print("p: ",dataOut.p[33,:])
3359 #print("p: ",dataOut.p[33,:])
3358 #exit(1)
3360 #exit(1)
3359 '''
3361 #'''
3360 return dataOut
3362 return dataOut
3361
3363
3362 class FaradayAngleAndDPPower(Operation):
3364 class FaradayAngleAndDPPower(Operation):
@@ -3446,6 +3448,7 class FaradayAngleAndDPPower(Operation):
3446
3448
3447 dataOut.flagTeTiCorrection = False
3449 dataOut.flagTeTiCorrection = False
3448 #print(dataOut.ph2)
3450 #print(dataOut.ph2)
3451
3449 #exit(1)
3452 #exit(1)
3450
3453
3451 return dataOut
3454 return dataOut
@@ -3519,10 +3522,11 class ElectronDensityFaraday(Operation):
3519 #exit(1)
3522 #exit(1)
3520 '''
3523 '''
3521 import matplotlib.pyplot as plt
3524 import matplotlib.pyplot as plt
3522 plt.plot(dataOut.bki)
3525 plt.plot(dataOut.phi,dataOut.heightList)
3523 plt.show()
3526 plt.show()
3524 '''
3527 '''
3525
3528 #print(dataOut.bki)
3529 print(dataOut.NDP)
3526 for i in range(2,dataOut.NSHTS-2):
3530 for i in range(2,dataOut.NSHTS-2):
3527 fact=(-0.5/(dataOut.RATE*dataOut.DH))*dataOut.bki[i]
3531 fact=(-0.5/(dataOut.RATE*dataOut.DH))*dataOut.bki[i]
3528 #four-point derivative, no phase unwrapping necessary
3532 #four-point derivative, no phase unwrapping necessary
@@ -3535,8 +3539,15 class ElectronDensityFaraday(Operation):
3535 dataOut.dphi[i]=abs(dataOut.dphi[i]*fact)
3539 dataOut.dphi[i]=abs(dataOut.dphi[i]*fact)
3536 dataOut.sdn1[i]=(4.*(dataOut.sdn2[i-2]+dataOut.sdn2[i+2])+dataOut.sdn2[i-1]+dataOut.sdn2[i+1])
3540 dataOut.sdn1[i]=(4.*(dataOut.sdn2[i-2]+dataOut.sdn2[i+2])+dataOut.sdn2[i-1]+dataOut.sdn2[i+1])
3537 dataOut.sdn1[i]=numpy.sqrt(dataOut.sdn1[i])*fact
3541 dataOut.sdn1[i]=numpy.sqrt(dataOut.sdn1[i])*fact
3542 '''
3538 #print(dataOut.dphi)
3543 #print(dataOut.dphi)
3539 #exit(1)
3544 #exit(1)
3545 import matplotlib.pyplot as plt
3546 plt.plot(dataOut.dphi,dataOut.heightList)
3547 plt.grid()
3548 plt.xlim(0,1e7)
3549 plt.show()
3550 '''
3540 return dataOut
3551 return dataOut
3541
3552
3542
3553
@@ -3962,6 +3973,7 class NormalizeDPPowerRoberto_V2(Operation):
3962 #print(dataOut.ph2)
3973 #print(dataOut.ph2)
3963 #input()
3974 #input()
3964 # in case of spread F, normalize much higher
3975 # in case of spread F, normalize much higher
3976 #print("dens: ", dataOut.dphi,dataOut.ph2)
3965 if(dataOut.cf<dataOut.cflast[0]/10.0):
3977 if(dataOut.cf<dataOut.cflast[0]/10.0):
3966 i1=(night_first1+100.-dataOut.range1[0])/dataOut.DH
3978 i1=(night_first1+100.-dataOut.range1[0])/dataOut.DH
3967 i2=(night_end+100.0-dataOut.range1[0])/dataOut.DH
3979 i2=(night_end+100.0-dataOut.range1[0])/dataOut.DH
@@ -3977,21 +3989,22 class NormalizeDPPowerRoberto_V2(Operation):
3977 #'''
3989 #'''
3978 #if (time_text.hour == 5 and time_text.minute == 32): #Year: 2022, DOY:104
3990 #if (time_text.hour == 5 and time_text.minute == 32): #Year: 2022, DOY:104
3979 #if (time_text.hour == 0 and time_text.minute == 12): #Year: 2022, DOY:93
3991 #if (time_text.hour == 0 and time_text.minute == 12): #Year: 2022, DOY:93
3980 #if (time_text.hour == 0 and time_text.minute == 22) or (time_text.hour == 0 and time_text.minute == 54) or (time_text.hour == 1 and time_text.minute == 48): #Year: 2022, DOY:242
3992 if (time_text.hour == 0 and time_text.minute == 22) or (time_text.hour == 0 and time_text.minute == 54) or (time_text.hour == 1 and time_text.minute == 48): #Year: 2022, DOY:242
3981 #if (time_text.hour == 1 and time_text.minute == 23) or (time_text.hour == 1 and time_text.minute == 44): #Year: 2022, DOY:243
3993 #if (time_text.hour == 1 and time_text.minute == 23) or (time_text.hour == 1 and time_text.minute == 44): #Year: 2022, DOY:243
3982 if (time_text.hour == 0 and time_text.minute == 4): #Year: 2022, DOY:244
3983 dataOut.cf = dataOut.cflast[0]
3994 dataOut.cf = dataOut.cflast[0]
3995 #if (time_text.hour == 0 and time_text.minute == 4): #Year: 2022, DOY:244
3984 #dataOut.cf = 0.08
3996 #dataOut.cf = 0.08
3985 #print("here")
3997 #print("here")
3986 if (time_text.hour == 2 and time_text.minute == 23): #Year: 2022, DOY:244
3998 #if (time_text.hour == 2 and time_text.minute == 23): #Year: 2022, DOY:244
3987 dataOut.cf = 0.08
3999 #dataOut.cf = 0.08
3988 if (time_text.hour == 2 and time_text.minute == 33): #Year: 2022, DOY:244
4000 #if (time_text.hour == 2 and time_text.minute == 33): #Year: 2022, DOY:244
3989 dataOut.cf = 0.09
4001 #dataOut.cf = 0.09
3990 if (time_text.hour == 3 and time_text.minute == 59) or (time_text.hour == 4 and time_text.minute == 20): #Year: 2022, DOY:244
4002 #if (time_text.hour == 3 and time_text.minute == 59) or (time_text.hour == 4 and time_text.minute == 20): #Year: 2022, DOY:244
3991 dataOut.cf = 0.09
4003 #dataOut.cf = 0.09
3992 #'''
4004 #'''
4005
3993 dataOut.cflast[0]=dataOut.cf
4006 dataOut.cflast[0]=dataOut.cf
3994 #print(dataOut.cf)
4007 print("cf: ", dataOut.cf)
3995
4008
3996 #print(dataOut.ph2)
4009 #print(dataOut.ph2)
3997 #print(dataOut.sdp2)
4010 #print(dataOut.sdp2)
@@ -4337,30 +4350,41 class DenCorrection(NormalizeDPPowerRoberto_V2):
4337 plt.title("{}".format(datetime.datetime.fromtimestamp(dataOut.utctime)))
4350 plt.title("{}".format(datetime.datetime.fromtimestamp(dataOut.utctime)))
4338 plt.xlim(.99,3)
4351 plt.xlim(.99,3)
4339 plt.grid()
4352 plt.grid()
4340 plt.savefig("/home/roberto/Pictures/Density_Comparison/TeTi_from_temps/{}.png".format(dataOut.utctime))
4353 plt.savefig("/home/roberto/Pictures/Density_Comparison/V2/TeTi_from_temps/{}.png".format(dataOut.utctime))
4341 '''
4354 '''
4342
4355 #dataOut.ti2 *= 5
4343 my_te2 = dataOut.ti2*ratio2
4356 my_te2 = dataOut.ti2*ratio2
4344 #'''
4357 #'''
4358
4359 te2_aux = dataOut.te2.copy()
4360 te2_aux[26:] = 1000
4361 te2_aux[:12] = 1000
4362 te2_aux -= 1000
4345 def func(params):
4363 def func(params):
4346 return (dataOut.te2-self.gaussian(dataOut.heightList[:dataOut.NSHTS],params[0],params[1],params[2]))
4364 #return (dataOut.te2-self.gaussian(dataOut.heightList[:dataOut.NSHTS],params[0],params[1],params[2]))
4347 x0_value = numpy.array([2000,250,20])
4365 return (te2_aux-self.gaussian(dataOut.heightList[:dataOut.NSHTS],params[0],params[1],params[2]))
4366 x0_value = numpy.array([1000,250,20])
4348 popt = least_squares(func,x0=x0_value,verbose=0)
4367 popt = least_squares(func,x0=x0_value,verbose=0)
4349 A = popt.x[0]; B = popt.x[1]; C = popt.x[2]
4368 A = popt.x[0]; B = popt.x[1]; C = popt.x[2]
4350 te2_smooth = self.gaussian(dataOut.heightList[:dataOut.NSHTS], A, B, C)
4369 te2_smooth = self.gaussian(dataOut.heightList[:dataOut.NSHTS], A, B, C)+1000
4370 te2_aux += 1000
4351 #'''
4371 #'''
4352
4372 ti2_smooth = te2_smooth/ratio2
4353 '''
4373 '''
4354 import matplotlib.pyplot as plt
4374 import matplotlib.pyplot as plt
4355 plt.clf()
4375 plt.clf()
4356 plt.plot(te2_smooth,dataOut.heightList[:dataOut.NSHTS],'*-',label = 'My Te')
4376 plt.plot(te2_smooth,dataOut.heightList[:dataOut.NSHTS],'*-',label = 'My Te')
4377 plt.plot(ti2_smooth,dataOut.heightList[:dataOut.NSHTS],'*-',label = 'My Ti')
4357 plt.plot(dataOut.te2,dataOut.heightList[:dataOut.NSHTS],'*-',label = 'Te')
4378 plt.plot(dataOut.te2,dataOut.heightList[:dataOut.NSHTS],'*-',label = 'Te')
4379 plt.plot(te2_aux,dataOut.heightList[:dataOut.NSHTS],'*-',label = 'Te_aux')
4380 plt.plot(ratio2*1000,dataOut.heightList[:dataOut.NSHTS],'*-',label = 'ratio*1000')
4358 #plt.plot(signal.medfilt(dataOut.te2),dataOut.heightList[:dataOut.NSHTS],'*-',label = 'Te')
4381 #plt.plot(signal.medfilt(dataOut.te2),dataOut.heightList[:dataOut.NSHTS],'*-',label = 'Te')
4359 plt.title("{}".format(datetime.datetime.fromtimestamp(dataOut.utctime)))
4382 plt.title("{}".format(datetime.datetime.fromtimestamp(dataOut.utctime)))
4360 plt.xlim(-50,3000)
4383 plt.xlim(-50,3000)
4361 plt.grid()
4384 plt.grid()
4362 plt.legend()
4385 plt.legend()
4363 plt.savefig("/home/roberto/Pictures/Density_Comparison/Te/{}.png".format(dataOut.utctime))
4386 plt.savefig("/home/roberto/Pictures/Density_Comparison/V3/Temps+{}.png".format(dataOut.utctime))
4387 #plt.show()
4364 '''
4388 '''
4365 #print("**** ACF2 WRAPPER ***** ",fitacf_acf2.acf2.__doc__ )
4389 #print("**** ACF2 WRAPPER ***** ",fitacf_acf2.acf2.__doc__ )
4366
4390
@@ -4370,18 +4394,21 class DenCorrection(NormalizeDPPowerRoberto_V2):
4370 nue=nui[0]=nui[1]=nui[2]=0.0#nui[3]=0.0
4394 nue=nui[0]=nui[1]=nui[2]=0.0#nui[3]=0.0
4371 wion[0]=16 #O
4395 wion[0]=16 #O
4372 wion[1]=1 #H
4396 wion[1]=1 #H
4373 wion[2]=4
4397 wion[2]=4 #He
4374 tion[0]=tion[1]=tion[2]=dataOut.ti2[i]
4398 #tion[0]=tion[1]=tion[2]=dataOut.ti2[i]
4399 tion[0]=tion[1]=tion[2]=ti2_smooth[i]
4375 fion[0]=1.0-dataOut.phy2[i] #1
4400 fion[0]=1.0-dataOut.phy2[i] #1
4376 fion[1]=dataOut.phy2[i] #0
4401 fion[1]=dataOut.phy2[i] #0
4377 fion[2]=0.0 #0
4402 fion[2]=0.0 #0
4378 for j in range(dataOut.DPL):
4403 for j in range(dataOut.DPL):
4379 tau=dataOut.alag[j]*1.0e-3
4404 tau=dataOut.alag[j]*1.0e-3
4380
4405
4381 with suppress_stdout_stderr():
4406 with suppress_stdout_stderr():#The smoothness in range of "y" depends on the smoothness of the input parameters
4382 y[j]=fitacf_acf2.acf2(wl,tau,dataOut.te2[i],tion,fion,nue,nui,wion,angle,dataOut.ph2[i],dataOut.bfm[i],y[j],three)
4407 y[j]=fitacf_acf2.acf2(wl,tau,dataOut.te2[i],tion,fion,nue,nui,wion,angle,dataOut.ph2[i],dataOut.bfm[i],y[j],three)
4383 #y[j]=fitacf_acf2.acf2(wl,tau,my_te2[i],tion,fion,nue,nui,wion,angle,dataOut.ph2[i],dataOut.bfm[i],y[j],three)
4408 #y[j]=fitacf_acf2.acf2(wl,tau,te2_smooth[i],tion,fion,nue,nui,wion,angle,dataOut.ph2[i],dataOut.bfm[i],y[j],three)
4384
4409
4410 #y[j]=fitacf_acf2.acf2(wl,tau,my_te2[i],tion,fion,nue,nui,wion,angle,dataOut.ph2[i],dataOut.bfm[i],y[j],three)
4411 #exit(1)
4385 #if dataOut.ut_Faraday>11.0 and dataOut.range1[i]>150.0 and dataOut.range1[i]<400.0:
4412 #if dataOut.ut_Faraday>11.0 and dataOut.range1[i]>150.0 and dataOut.range1[i]<400.0:
4386
4413
4387 if dataOut.ut_Faraday>11.0 and dataOut.range1[i]>150.0 and dataOut.range1[i]<300.0:
4414 if dataOut.ut_Faraday>11.0 and dataOut.range1[i]>150.0 and dataOut.range1[i]<300.0:
@@ -4419,7 +4446,7 class DenCorrection(NormalizeDPPowerRoberto_V2):
4419 '''
4446 '''
4420 import matplotlib.pyplot as plt
4447 import matplotlib.pyplot as plt
4421 plt.clf()
4448 plt.clf()
4422 plt.plot(aux,dataOut.heightList[:dataOut.NSHTS],'*:',label='Fitting')
4449 #plt.plot(aux,dataOut.heightList[:dataOut.NSHTS],'*:',label='Fitting')
4423 plt.plot(my_aux,dataOut.heightList[:dataOut.NSHTS],'*:',label='Ratio')
4450 plt.plot(my_aux,dataOut.heightList[:dataOut.NSHTS],'*:',label='Ratio')
4424 #plt.plot(acf_Temps,dataOut.heightList[:dataOut.NSHTS],'b*:',label='Temps')
4451 #plt.plot(acf_Temps,dataOut.heightList[:dataOut.NSHTS],'b*:',label='Temps')
4425 #plt.plot(acf_no_Temps,dataOut.heightList[:dataOut.NSHTS],'k*:',label='No Temps')
4452 #plt.plot(acf_no_Temps,dataOut.heightList[:dataOut.NSHTS],'k*:',label='No Temps')
@@ -4431,7 +4458,9 class DenCorrection(NormalizeDPPowerRoberto_V2):
4431 #plt.xlim(.99,1.25)
4458 #plt.xlim(.99,1.25)
4432 #plt.show()
4459 #plt.show()
4433 #plt.savefig("/home/roberto/Pictures/Density_Comparison/FactorEf_NoLimits/{}.png".format(dataOut.utctime))
4460 #plt.savefig("/home/roberto/Pictures/Density_Comparison/FactorEf_NoLimits/{}.png".format(dataOut.utctime))
4434 plt.savefig("/home/roberto/Pictures/Faraday/2022/08/Density_Comparison/FactorEf/{}.png".format(dataOut.utctime))
4461 #plt.savefig("/home/roberto/Pictures/Density_Comparison/V2/TeTi_cte/Te2/{}.png".format(dataOut.utctime))
4462 #plt.savefig("/home/roberto/Pictures/Density_Comparison/V2/bline/Te2/{}.png".format(dataOut.utctime))
4463 plt.savefig("/home/roberto/Pictures/Density_Comparison/V3/{}.png".format(dataOut.utctime))
4435 #plt.savefig("/home/roberto/Pictures/Faraday_TeTi_Test/Ratio/{}.png".format(dataOut.utctime))
4464 #plt.savefig("/home/roberto/Pictures/Faraday_TeTi_Test/Ratio/{}.png".format(dataOut.utctime))
4436 '''
4465 '''
4437 #print("inside correction",dataOut.ph2)
4466 #print("inside correction",dataOut.ph2)
@@ -4775,10 +4804,10 class DataSaveCleaner(Operation):
4775
4804
4776 #print(time_text.hour,time_text.minute)
4805 #print(time_text.hour,time_text.minute)
4777 #if (time_text.hour == 16 and time_text.minute==48) or (time_text.hour == 19 and time_text.minute ==49 ) or (time_text.hour >= 0 and time_text.hour < 5): #Year: 2022, DOY:241
4806 #if (time_text.hour == 16 and time_text.minute==48) or (time_text.hour == 19 and time_text.minute ==49 ) or (time_text.hour >= 0 and time_text.hour < 5): #Year: 2022, DOY:241
4778 #if (time_text.hour == 5 and time_text.minute==21) or (time_text.hour == 19 and time_text.minute ==49 ) or (time_text.hour == 7 and time_text.minute==40) or (time_text.hour == 7 and time_text.minute==50) or (time_text.hour >= 8 and time_text.hour < 11) or (time_text.hour == 11 and time_text.minute==2) or (time_text.hour == 11 and time_text.minute==13): #Year: 2022, DOY:242
4807 if (time_text.hour == 5 and time_text.minute==21) or (time_text.hour == 19 and time_text.minute ==49 ) or (time_text.hour == 7 and time_text.minute==40) or (time_text.hour == 7 and time_text.minute==50) or (time_text.hour >= 8 and time_text.hour < 11) or (time_text.hour == 11 and time_text.minute==2) or (time_text.hour == 11 and time_text.minute==13): #Year: 2022, DOY:242
4779 #if (time_text.hour >= 8 and time_text.hour < 11) or (time_text.hour == 11 and time_text.minute==2) or (time_text.hour == 11 and time_text.minute==13) or (time_text.hour == 11 and time_text.minute==24): #Year: 2022, DOY:243
4808 #if (time_text.hour >= 8 and time_text.hour < 11) or (time_text.hour == 11 and time_text.minute==2) or (time_text.hour == 11 and time_text.minute==13) or (time_text.hour == 11 and time_text.minute==24): #Year: 2022, DOY:243
4780 #if (time_text.hour >= 9 and time_text.hour < 11) or (time_text.hour == 8 and time_text.minute==12) or (time_text.hour == 8 and time_text.minute==22) or (time_text.hour == 8 and time_text.minute==33) or (time_text.hour == 8 and time_text.minute==44) or (time_text.hour == 8 and time_text.minute==54) or (time_text.hour == 11 and time_text.minute==2) or (time_text.hour == 11 and time_text.minute==13): #Year: 2022, DOY:245
4809 #if (time_text.hour >= 9 and time_text.hour < 11) or (time_text.hour == 8 and time_text.minute==12) or (time_text.hour == 8 and time_text.minute==22) or (time_text.hour == 8 and time_text.minute==33) or (time_text.hour == 8 and time_text.minute==44) or (time_text.hour == 8 and time_text.minute==54) or (time_text.hour == 11 and time_text.minute==2) or (time_text.hour == 11 and time_text.minute==13): #Year: 2022, DOY:245
4781 if (time_text.hour >= 8 and time_text.hour < 11) or (time_text.hour == 1) or (time_text.hour == 0 and time_text.minute==25) or (time_text.hour == 0 and time_text.minute==36) or (time_text.hour == 0 and time_text.minute==47) or (time_text.hour == 0 and time_text.minute==57) or (time_text.hour == 2 and time_text.minute==1) or (time_text.hour == 11 and time_text.minute==2) or (time_text.hour == 11 and time_text.minute==13) or (time_text.hour == 11 and time_text.minute==24) or (time_text.hour == 7 and time_text.minute==40) or (time_text.hour == 7 and time_text.minute==50) or (time_text.hour == 3 and time_text.minute==5): #Year: 2022, DOY:244
4810 #if (time_text.hour >= 8 and time_text.hour < 11) or (time_text.hour == 1) or (time_text.hour == 0 and time_text.minute==15) or (time_text.hour == 0 and time_text.minute==25) or (time_text.hour == 0 and time_text.minute==36) or (time_text.hour == 0 and time_text.minute==47) or (time_text.hour == 0 and time_text.minute==57) or (time_text.hour == 2 and time_text.minute==1) or (time_text.hour == 11 and time_text.minute==2) or (time_text.hour == 11 and time_text.minute==13) or (time_text.hour == 11 and time_text.minute==24) or (time_text.hour == 7 and time_text.minute==40) or (time_text.hour == 7 and time_text.minute==50) or (time_text.hour == 3 and time_text.minute==5) or (time_text.hour == 3 and time_text.minute==16) or (time_text.hour == 3 and time_text.minute==27): #Year: 2022, DOY:244
4782
4811
4783 dataOut.DensityFinal[0,:]=missing
4812 dataOut.DensityFinal[0,:]=missing
4784 dataOut.EDensityFinal[0,:]=missing
4813 dataOut.EDensityFinal[0,:]=missing
@@ -4813,7 +4842,17 class DataSaveCleaner(Operation):
4813 dataOut.EIonTempFinal[0,id_aux:]=missing
4842 dataOut.EIonTempFinal[0,id_aux:]=missing
4814 dataOut.PhyFinal[0,id_aux:]=missing
4843 dataOut.PhyFinal[0,id_aux:]=missing
4815 dataOut.EPhyFinal[0,id_aux:]=missing
4844 dataOut.EPhyFinal[0,id_aux:]=missing
4816 if (time_text.hour == 20 and time_text.minute == 29) or (time_text.hour == 20 and time_text.minute == 44): #Year: 2022, DOY:243
4845 if (time_text.hour == 20 and time_text.minute == 29): #Year: 2022, DOY:243
4846 id_aux = 30
4847 dataOut.DensityFinal[0,id_aux:]=missing
4848 dataOut.EDensityFinal[0,id_aux:]=missing
4849 dataOut.ElecTempFinal[0,id_aux:]=missing
4850 dataOut.EElecTempFinal[0,id_aux:]=missing
4851 dataOut.IonTempFinal[0,id_aux:]=missing
4852 dataOut.EIonTempFinal[0,id_aux:]=missing
4853 dataOut.PhyFinal[0,id_aux:]=missing
4854 dataOut.EPhyFinal[0,id_aux:]=missing
4855 if (time_text.hour == 20 and time_text.minute == 44): #Year: 2022, DOY:243
4817 id_aux = 31
4856 id_aux = 31
4818 dataOut.DensityFinal[0,id_aux:]=missing
4857 dataOut.DensityFinal[0,id_aux:]=missing
4819 dataOut.EDensityFinal[0,id_aux:]=missing
4858 dataOut.EDensityFinal[0,id_aux:]=missing
@@ -4823,8 +4862,187 class DataSaveCleaner(Operation):
4823 dataOut.EIonTempFinal[0,id_aux:]=missing
4862 dataOut.EIonTempFinal[0,id_aux:]=missing
4824 dataOut.PhyFinal[0,id_aux:]=missing
4863 dataOut.PhyFinal[0,id_aux:]=missing
4825 dataOut.EPhyFinal[0,id_aux:]=missing
4864 dataOut.EPhyFinal[0,id_aux:]=missing
4865 if (time_text.hour <= 8): #Year: 2022, DOY:243
4866 id_aux = 11
4867 dataOut.DensityFinal[0,:id_aux]=missing
4868 dataOut.EDensityFinal[0,:id_aux]=missing
4869 dataOut.ElecTempFinal[0,:id_aux]=missing
4870 dataOut.EElecTempFinal[0,:id_aux]=missing
4871 dataOut.IonTempFinal[0,:id_aux]=missing
4872 dataOut.EIonTempFinal[0,:id_aux]=missing
4873 dataOut.PhyFinal[0,:id_aux]=missing
4874 dataOut.EPhyFinal[0,:id_aux]=missing
4875 if (time_text.hour == 23): #Year: 2022, DOY:243
4876 id_aux = 12
4877 dataOut.DensityFinal[0,:id_aux]=missing
4878 dataOut.EDensityFinal[0,:id_aux]=missing
4879 dataOut.ElecTempFinal[0,:id_aux]=missing
4880 dataOut.EElecTempFinal[0,:id_aux]=missing
4881 dataOut.IonTempFinal[0,:id_aux]=missing
4882 dataOut.EIonTempFinal[0,:id_aux]=missing
4883 dataOut.PhyFinal[0,:id_aux]=missing
4884 dataOut.EPhyFinal[0,:id_aux]=missing
4885 if (time_text.hour == 5 and time_text.minute == 21): #Year: 2022, DOY:243
4886 id_aux = (36,37)
4887 dataOut.DensityFinal[0,id_aux]=missing
4888 dataOut.EDensityFinal[0,id_aux]=missing
4889 dataOut.ElecTempFinal[0,id_aux]=missing
4890 dataOut.EElecTempFinal[0,id_aux]=missing
4891 dataOut.IonTempFinal[0,id_aux]=missing
4892 dataOut.EIonTempFinal[0,id_aux]=missing
4893 dataOut.PhyFinal[0,id_aux]=missing
4894 dataOut.EPhyFinal[0,id_aux]=missing
4895 if (time_text.hour == 5 and time_text.minute == 53): #Year: 2022, DOY:243
4896 id_aux = (37,38)
4897 dataOut.DensityFinal[0,id_aux]=missing
4898 dataOut.EDensityFinal[0,id_aux]=missing
4899 dataOut.ElecTempFinal[0,id_aux]=missing
4900 dataOut.EElecTempFinal[0,id_aux]=missing
4901 dataOut.IonTempFinal[0,id_aux]=missing
4902 dataOut.EIonTempFinal[0,id_aux]=missing
4903 dataOut.PhyFinal[0,id_aux]=missing
4904 dataOut.EPhyFinal[0,id_aux]=missing
4905 if (time_text.hour == 6 and time_text.minute == 4): #Year: 2022, DOY:243
4906 id_aux = (38,39)
4907 dataOut.DensityFinal[0,id_aux]=missing
4908 dataOut.EDensityFinal[0,id_aux]=missing
4909 dataOut.ElecTempFinal[0,id_aux]=missing
4910 dataOut.EElecTempFinal[0,id_aux]=missing
4911 dataOut.IonTempFinal[0,id_aux]=missing
4912 dataOut.EIonTempFinal[0,id_aux]=missing
4913 dataOut.PhyFinal[0,id_aux]=missing
4914 dataOut.EPhyFinal[0,id_aux]=missing
4915 if (time_text.hour == 12 and time_text.minute == 6): #Year: 2022, DOY:243
4916 id_aux = (29,30)
4917 dataOut.DensityFinal[0,id_aux]=missing
4918 dataOut.EDensityFinal[0,id_aux]=missing
4919 dataOut.ElecTempFinal[0,id_aux]=missing
4920 dataOut.EElecTempFinal[0,id_aux]=missing
4921 dataOut.IonTempFinal[0,id_aux]=missing
4922 dataOut.EIonTempFinal[0,id_aux]=missing
4923 dataOut.PhyFinal[0,id_aux]=missing
4924 dataOut.EPhyFinal[0,id_aux]=missing
4925 if (time_text.hour == 14 and time_text.minute == 14): #Year: 2022, DOY:243
4926 id_aux = (35,36)
4927 dataOut.DensityFinal[0,id_aux]=missing
4928 dataOut.EDensityFinal[0,id_aux]=missing
4929 dataOut.ElecTempFinal[0,id_aux]=missing
4930 dataOut.EElecTempFinal[0,id_aux]=missing
4931 dataOut.IonTempFinal[0,id_aux]=missing
4932 dataOut.EIonTempFinal[0,id_aux]=missing
4933 dataOut.PhyFinal[0,id_aux]=missing
4934 dataOut.EPhyFinal[0,id_aux]=missing
4935 if (time_text.hour == 23 and time_text.minute == 2): #Year: 2022, DOY:243
4936 id_aux = (41,42)
4937 dataOut.DensityFinal[0,id_aux]=missing
4938 dataOut.EDensityFinal[0,id_aux]=missing
4939 dataOut.ElecTempFinal[0,id_aux]=missing
4940 dataOut.EElecTempFinal[0,id_aux]=missing
4941 dataOut.IonTempFinal[0,id_aux]=missing
4942 dataOut.EIonTempFinal[0,id_aux]=missing
4943 dataOut.PhyFinal[0,id_aux]=missing
4944 dataOut.EPhyFinal[0,id_aux]=missing
4945 if (time_text.hour == 0 and time_text.minute == 8): #Year: 2022, DOY:243
4946 id_aux = 33
4947 dataOut.DensityFinal[0,id_aux:]=missing
4948 dataOut.EDensityFinal[0,id_aux:]=missing
4949 dataOut.ElecTempFinal[0,id_aux:]=missing
4950 dataOut.EElecTempFinal[0,id_aux:]=missing
4951 dataOut.IonTempFinal[0,id_aux:]=missing
4952 dataOut.EIonTempFinal[0,id_aux:]=missing
4953 dataOut.PhyFinal[0,id_aux:]=missing
4954 dataOut.EPhyFinal[0,id_aux:]=missing
4955 id_aux = 18
4956 dataOut.DensityFinal[0,id_aux]=missing
4957 dataOut.EDensityFinal[0,id_aux]=missing
4958 dataOut.ElecTempFinal[0,id_aux]=missing
4959 dataOut.EElecTempFinal[0,id_aux]=missing
4960 dataOut.IonTempFinal[0,id_aux]=missing
4961 dataOut.EIonTempFinal[0,id_aux]=missing
4962 dataOut.PhyFinal[0,id_aux]=missing
4963 dataOut.EPhyFinal[0,id_aux]=missing
4964 if (time_text.hour == 23 and time_text.minute == 26): #Year: 2022, DOY:243
4965 id_aux = (12,13,14)
4966 dataOut.DensityFinal[0,id_aux]=missing
4967 dataOut.EDensityFinal[0,id_aux]=missing
4968 dataOut.ElecTempFinal[0,id_aux]=missing
4969 dataOut.EElecTempFinal[0,id_aux]=missing
4970 dataOut.IonTempFinal[0,id_aux]=missing
4971 dataOut.EIonTempFinal[0,id_aux]=missing
4972 dataOut.PhyFinal[0,id_aux]=missing
4973 dataOut.EPhyFinal[0,id_aux]=missing
4974 if (time_text.hour == 23 and time_text.minute == 36): #Year: 2022, DOY:243
4975 id_aux = (14,15,16)
4976 dataOut.DensityFinal[0,id_aux]=missing
4977 dataOut.EDensityFinal[0,id_aux]=missing
4978 dataOut.ElecTempFinal[0,id_aux]=missing
4979 dataOut.EElecTempFinal[0,id_aux]=missing
4980 dataOut.IonTempFinal[0,id_aux]=missing
4981 dataOut.EIonTempFinal[0,id_aux]=missing
4982 dataOut.PhyFinal[0,id_aux]=missing
4983 dataOut.EPhyFinal[0,id_aux]=missing
4984 if (time_text.hour == 2 and time_text.minute == 6): #Year: 2022, DOY:243
4985 id_aux = (36,37,38)
4986 dataOut.DensityFinal[0,id_aux]=missing
4987 dataOut.EDensityFinal[0,id_aux]=missing
4988 dataOut.ElecTempFinal[0,id_aux]=missing
4989 dataOut.EElecTempFinal[0,id_aux]=missing
4990 dataOut.IonTempFinal[0,id_aux]=missing
4991 dataOut.EIonTempFinal[0,id_aux]=missing
4992 dataOut.PhyFinal[0,id_aux]=missing
4993 dataOut.EPhyFinal[0,id_aux]=missing
4994 if (time_text.hour == 2 and time_text.minute == 16): #Year: 2022, DOY:243
4995 id_aux = (34,35)
4996 dataOut.DensityFinal[0,id_aux]=missing
4997 dataOut.EDensityFinal[0,id_aux]=missing
4998 dataOut.ElecTempFinal[0,id_aux]=missing
4999 dataOut.EElecTempFinal[0,id_aux]=missing
5000 dataOut.IonTempFinal[0,id_aux]=missing
5001 dataOut.EIonTempFinal[0,id_aux]=missing
5002 dataOut.PhyFinal[0,id_aux]=missing
5003 dataOut.EPhyFinal[0,id_aux]=missing
5004 if (time_text.hour == 2 and time_text.minute == 38): #Year: 2022, DOY:243
5005 id_aux = (35,36)
5006 dataOut.DensityFinal[0,id_aux]=missing
5007 dataOut.EDensityFinal[0,id_aux]=missing
5008 dataOut.ElecTempFinal[0,id_aux]=missing
5009 dataOut.EElecTempFinal[0,id_aux]=missing
5010 dataOut.IonTempFinal[0,id_aux]=missing
5011 dataOut.EIonTempFinal[0,id_aux]=missing
5012 dataOut.PhyFinal[0,id_aux]=missing
5013 dataOut.EPhyFinal[0,id_aux]=missing
5014 if (time_text.hour == 3 and time_text.minute == 20): #Year: 2022, DOY:243
5015 id_aux = (33,34)
5016 dataOut.DensityFinal[0,id_aux]=missing
5017 dataOut.EDensityFinal[0,id_aux]=missing
5018 dataOut.ElecTempFinal[0,id_aux]=missing
5019 dataOut.EElecTempFinal[0,id_aux]=missing
5020 dataOut.IonTempFinal[0,id_aux]=missing
5021 dataOut.EIonTempFinal[0,id_aux]=missing
5022 dataOut.PhyFinal[0,id_aux]=missing
5023 dataOut.EPhyFinal[0,id_aux]=missing
5024 if (time_text.hour == 3 and time_text.minute == 42): #Year: 2022, DOY:243
5025 id_aux = 34
5026 dataOut.DensityFinal[0,id_aux:]=missing
5027 dataOut.EDensityFinal[0,id_aux:]=missing
5028 dataOut.ElecTempFinal[0,id_aux:]=missing
5029 dataOut.EElecTempFinal[0,id_aux:]=missing
5030 dataOut.IonTempFinal[0,id_aux:]=missing
5031 dataOut.EIonTempFinal[0,id_aux:]=missing
5032 dataOut.PhyFinal[0,id_aux:]=missing
5033 dataOut.EPhyFinal[0,id_aux:]=missing
5034 if (time_text.hour == 4 and time_text.minute == 35): #Year: 2022, DOY:243
5035 id_aux = (36,37)
5036 dataOut.DensityFinal[0,id_aux]=missing
5037 dataOut.EDensityFinal[0,id_aux]=missing
5038 dataOut.ElecTempFinal[0,id_aux]=missing
5039 dataOut.EElecTempFinal[0,id_aux]=missing
5040 dataOut.IonTempFinal[0,id_aux]=missing
5041 dataOut.EIonTempFinal[0,id_aux]=missing
5042 dataOut.PhyFinal[0,id_aux]=missing
5043 dataOut.EPhyFinal[0,id_aux]=missing
4826 '''
5044 '''
4827 #'''
5045 '''
4828 if (time_text.hour == 2 and time_text.minute == 23): #Year: 2022, DOY:244
5046 if (time_text.hour == 2 and time_text.minute == 23): #Year: 2022, DOY:244
4829 id_aux = 12
5047 id_aux = 12
4830 dataOut.DensityFinal[0,:id_aux]=missing
5048 dataOut.DensityFinal[0,:id_aux]=missing
@@ -4905,16 +5123,367 class DataSaveCleaner(Operation):
4905 dataOut.EIonTempFinal[0,id_aux:]=missing
5123 dataOut.EIonTempFinal[0,id_aux:]=missing
4906 dataOut.PhyFinal[0,id_aux:]=missing
5124 dataOut.PhyFinal[0,id_aux:]=missing
4907 dataOut.EPhyFinal[0,id_aux:]=missing
5125 dataOut.EPhyFinal[0,id_aux:]=missing
5126 if (time_text.hour <= 8): #Year: 2022, DOY:244
5127 id_aux = 12
5128 dataOut.DensityFinal[0,:id_aux]=missing
5129 dataOut.EDensityFinal[0,:id_aux]=missing
5130 dataOut.ElecTempFinal[0,:id_aux]=missing
5131 dataOut.EElecTempFinal[0,:id_aux]=missing
5132 dataOut.IonTempFinal[0,:id_aux]=missing
5133 dataOut.EIonTempFinal[0,:id_aux]=missing
5134 dataOut.PhyFinal[0,:id_aux]=missing
5135 dataOut.EPhyFinal[0,:id_aux]=missing
5136 if (time_text.hour == 23): #Year: 2022, DOY:244
5137 id_aux = 12
5138 dataOut.DensityFinal[0,:id_aux]=missing
5139 dataOut.EDensityFinal[0,:id_aux]=missing
5140 dataOut.ElecTempFinal[0,:id_aux]=missing
5141 dataOut.EElecTempFinal[0,:id_aux]=missing
5142 dataOut.IonTempFinal[0,:id_aux]=missing
5143 dataOut.EIonTempFinal[0,:id_aux]=missing
5144 dataOut.PhyFinal[0,:id_aux]=missing
5145 dataOut.EPhyFinal[0,:id_aux]=missing
5146 if (time_text.hour == 5 and time_text.minute == 42): #Year: 2022, DOY:244
5147 id_aux = (32,33)
5148 dataOut.DensityFinal[0,id_aux]=missing
5149 dataOut.EDensityFinal[0,id_aux]=missing
5150 dataOut.ElecTempFinal[0,id_aux]=missing
5151 dataOut.EElecTempFinal[0,id_aux]=missing
5152 dataOut.IonTempFinal[0,id_aux]=missing
5153 dataOut.EIonTempFinal[0,id_aux]=missing
5154 dataOut.PhyFinal[0,id_aux]=missing
5155 dataOut.EPhyFinal[0,id_aux]=missing
5156 if (time_text.hour == 11 and time_text.minute == 56): #Year: 2022, DOY:244
5157 id_aux = (39,40)
5158 dataOut.DensityFinal[0,id_aux]=missing
5159 dataOut.EDensityFinal[0,id_aux]=missing
5160 dataOut.ElecTempFinal[0,id_aux]=missing
5161 dataOut.EElecTempFinal[0,id_aux]=missing
5162 dataOut.IonTempFinal[0,id_aux]=missing
5163 dataOut.EIonTempFinal[0,id_aux]=missing
5164 dataOut.PhyFinal[0,id_aux]=missing
5165 dataOut.EPhyFinal[0,id_aux]=missing
5166 if (time_text.hour == 12 and time_text.minute == 52): #Year: 2022, DOY:244
5167 id_aux = (36,37)
5168 dataOut.DensityFinal[0,id_aux]=missing
5169 dataOut.EDensityFinal[0,id_aux]=missing
5170 dataOut.ElecTempFinal[0,id_aux]=missing
5171 dataOut.EElecTempFinal[0,id_aux]=missing
5172 dataOut.IonTempFinal[0,id_aux]=missing
5173 dataOut.EIonTempFinal[0,id_aux]=missing
5174 dataOut.PhyFinal[0,id_aux]=missing
5175 dataOut.EPhyFinal[0,id_aux]=missing
5176 if (time_text.hour == 13 and time_text.minute == 3): #Year: 2022, DOY:244
5177 id_aux = (37,38)
5178 dataOut.DensityFinal[0,id_aux]=missing
5179 dataOut.EDensityFinal[0,id_aux]=missing
5180 dataOut.ElecTempFinal[0,id_aux]=missing
5181 dataOut.EElecTempFinal[0,id_aux]=missing
5182 dataOut.IonTempFinal[0,id_aux]=missing
5183 dataOut.EIonTempFinal[0,id_aux]=missing
5184 dataOut.PhyFinal[0,id_aux]=missing
5185 dataOut.EPhyFinal[0,id_aux]=missing
5186 if (time_text.hour == 23 and time_text.minute == 11): #Year: 2022, DOY:244
5187 id_aux = (40,41)
5188 dataOut.DensityFinal[0,id_aux]=missing
5189 dataOut.EDensityFinal[0,id_aux]=missing
5190 dataOut.ElecTempFinal[0,id_aux]=missing
5191 dataOut.EElecTempFinal[0,id_aux]=missing
5192 dataOut.IonTempFinal[0,id_aux]=missing
5193 dataOut.EIonTempFinal[0,id_aux]=missing
5194 dataOut.PhyFinal[0,id_aux]=missing
5195 dataOut.EPhyFinal[0,id_aux]=missing
5196 if (time_text.hour == 23 and time_text.minute == 21): #Year: 2022, DOY:244
5197 id_aux = (12,13,39,40,41)
5198 dataOut.DensityFinal[0,id_aux]=missing
5199 dataOut.EDensityFinal[0,id_aux]=missing
5200 dataOut.ElecTempFinal[0,id_aux]=missing
5201 dataOut.EElecTempFinal[0,id_aux]=missing
5202 dataOut.IonTempFinal[0,id_aux]=missing
5203 dataOut.EIonTempFinal[0,id_aux]=missing
5204 dataOut.PhyFinal[0,id_aux]=missing
5205 dataOut.EPhyFinal[0,id_aux]=missing
5206 if (time_text.hour == 23 and time_text.minute == 53): #Year: 2022, DOY:244
5207 id_aux = (15,16,17,18)
5208 dataOut.DensityFinal[0,id_aux]=missing
5209 dataOut.EDensityFinal[0,id_aux]=missing
5210 dataOut.ElecTempFinal[0,id_aux]=missing
5211 dataOut.EElecTempFinal[0,id_aux]=missing
5212 dataOut.IonTempFinal[0,id_aux]=missing
5213 dataOut.EIonTempFinal[0,id_aux]=missing
5214 dataOut.PhyFinal[0,id_aux]=missing
5215 dataOut.EPhyFinal[0,id_aux]=missing
5216 if (time_text.hour == 2 and time_text.minute == 44): #Year: 2022, DOY:244
5217 id_aux = (40,41,42)
5218 dataOut.DensityFinal[0,id_aux]=missing
5219 dataOut.EDensityFinal[0,id_aux]=missing
5220 dataOut.ElecTempFinal[0,id_aux]=missing
5221 dataOut.EElecTempFinal[0,id_aux]=missing
5222 dataOut.IonTempFinal[0,id_aux]=missing
5223 dataOut.EIonTempFinal[0,id_aux]=missing
5224 dataOut.PhyFinal[0,id_aux]=missing
5225 dataOut.EPhyFinal[0,id_aux]=missing
5226 if (time_text.hour == 3 and time_text.minute == 37): #Year: 2022, DOY:244
5227 id_aux = (36,37)
5228 dataOut.DensityFinal[0,id_aux]=missing
5229 dataOut.EDensityFinal[0,id_aux]=missing
5230 dataOut.ElecTempFinal[0,id_aux]=missing
5231 dataOut.EElecTempFinal[0,id_aux]=missing
5232 dataOut.IonTempFinal[0,id_aux]=missing
5233 dataOut.EIonTempFinal[0,id_aux]=missing
5234 dataOut.PhyFinal[0,id_aux]=missing
5235 dataOut.EPhyFinal[0,id_aux]=missing
5236 if (time_text.hour == 4 and time_text.minute == 9): #Year: 2022, DOY:244
5237 id_aux = (32,33,34)
5238 dataOut.DensityFinal[0,id_aux]=missing
5239 dataOut.EDensityFinal[0,id_aux]=missing
5240 dataOut.ElecTempFinal[0,id_aux]=missing
5241 dataOut.EElecTempFinal[0,id_aux]=missing
5242 dataOut.IonTempFinal[0,id_aux]=missing
5243 dataOut.EIonTempFinal[0,id_aux]=missing
5244 dataOut.PhyFinal[0,id_aux]=missing
5245 dataOut.EPhyFinal[0,id_aux]=missing
5246 if (time_text.hour == 4 and time_text.minute == 20): #Year: 2022, DOY:244
5247 id_aux = 37
5248 dataOut.DensityFinal[0,id_aux:]=missing
5249 dataOut.EDensityFinal[0,id_aux:]=missing
5250 dataOut.ElecTempFinal[0,id_aux:]=missing
5251 dataOut.EElecTempFinal[0,id_aux:]=missing
5252 dataOut.IonTempFinal[0,id_aux:]=missing
5253 dataOut.EIonTempFinal[0,id_aux:]=missing
5254 dataOut.PhyFinal[0,id_aux:]=missing
5255 dataOut.EPhyFinal[0,id_aux:]=missing
5256 if (time_text.hour == 4 and time_text.minute == 31): #Year: 2022, DOY:244
5257 id_aux = 33
5258 dataOut.DensityFinal[0,id_aux:]=missing
5259 dataOut.EDensityFinal[0,id_aux:]=missing
5260 dataOut.ElecTempFinal[0,id_aux:]=missing
5261 dataOut.EElecTempFinal[0,id_aux:]=missing
5262 dataOut.IonTempFinal[0,id_aux:]=missing
5263 dataOut.EIonTempFinal[0,id_aux:]=missing
5264 dataOut.PhyFinal[0,id_aux:]=missing
5265 dataOut.EPhyFinal[0,id_aux:]=missing
5266 '''
5267 '''
5268 if (time_text.hour <= 10): #Year: 2022, DOY:245
5269 id_aux = 11
5270 dataOut.DensityFinal[0,:id_aux]=missing
5271 dataOut.EDensityFinal[0,:id_aux]=missing
5272 dataOut.ElecTempFinal[0,:id_aux]=missing
5273 dataOut.EElecTempFinal[0,:id_aux]=missing
5274 dataOut.IonTempFinal[0,:id_aux]=missing
5275 dataOut.EIonTempFinal[0,:id_aux]=missing
5276 dataOut.PhyFinal[0,:id_aux]=missing
5277 dataOut.EPhyFinal[0,:id_aux]=missing
5278 if (time_text.hour == 5 and time_text.minute == 10): #Year: 2022, DOY:245
5279 id_aux = 35
5280 dataOut.DensityFinal[0,id_aux:]=missing
5281 dataOut.EDensityFinal[0,id_aux:]=missing
5282 dataOut.ElecTempFinal[0,id_aux:]=missing
5283 dataOut.EElecTempFinal[0,id_aux:]=missing
5284 dataOut.IonTempFinal[0,id_aux:]=missing
5285 dataOut.EIonTempFinal[0,id_aux:]=missing
5286 dataOut.PhyFinal[0,id_aux:]=missing
5287 dataOut.EPhyFinal[0,id_aux:]=missing
5288 if (time_text.hour == 5 and time_text.minute == 21): #Year: 2022, DOY:245
5289 id_aux = 36
5290 dataOut.DensityFinal[0,id_aux:]=missing
5291 dataOut.EDensityFinal[0,id_aux:]=missing
5292 dataOut.ElecTempFinal[0,id_aux:]=missing
5293 dataOut.EElecTempFinal[0,id_aux:]=missing
5294 dataOut.IonTempFinal[0,id_aux:]=missing
5295 dataOut.EIonTempFinal[0,id_aux:]=missing
5296 dataOut.PhyFinal[0,id_aux:]=missing
5297 dataOut.EPhyFinal[0,id_aux:]=missing
5298 if (time_text.hour == 11 and time_text.minute == 45): #Year: 2022, DOY:245
5299 id_aux = 7
5300 dataOut.DensityFinal[0,id_aux]=missing
5301 dataOut.EDensityFinal[0,id_aux]=missing
5302 dataOut.ElecTempFinal[0,id_aux]=missing
5303 dataOut.EElecTempFinal[0,id_aux]=missing
5304 dataOut.IonTempFinal[0,id_aux]=missing
5305 dataOut.EIonTempFinal[0,id_aux]=missing
5306 dataOut.PhyFinal[0,id_aux]=missing
5307 dataOut.EPhyFinal[0,id_aux]=missing
5308 '''
5309 '''
5310 if (time_text.hour == 23 and time_text.minute > 30): #Year: 2022, DOY:241
5311 id_aux = 17
5312 dataOut.DensityFinal[0,:id_aux]=missing
5313 dataOut.EDensityFinal[0,:id_aux]=missing
5314 dataOut.ElecTempFinal[0,:id_aux]=missing
5315 dataOut.EElecTempFinal[0,:id_aux]=missing
5316 dataOut.IonTempFinal[0,:id_aux]=missing
5317 dataOut.EIonTempFinal[0,:id_aux]=missing
5318 dataOut.PhyFinal[0,:id_aux]=missing
5319 dataOut.EPhyFinal[0,:id_aux]=missing
5320
5321 if (time_text.hour == 13 and time_text.minute == 36): #Year: 2022, DOY:241
5322 id_aux = 33
5323 dataOut.DensityFinal[0,id_aux:]=missing
5324 dataOut.EDensityFinal[0,id_aux:]=missing
5325 dataOut.ElecTempFinal[0,id_aux:]=missing
5326 dataOut.EElecTempFinal[0,id_aux:]=missing
5327 dataOut.IonTempFinal[0,id_aux:]=missing
5328 dataOut.EIonTempFinal[0,id_aux:]=missing
5329 dataOut.PhyFinal[0,id_aux:]=missing
5330 dataOut.EPhyFinal[0,id_aux:]=missing
5331 if (time_text.hour == 13 and time_text.minute == 47): #Year: 2022, DOY:241
5332 id_aux = 36
5333 dataOut.DensityFinal[0,id_aux:]=missing
5334 dataOut.EDensityFinal[0,id_aux:]=missing
5335 dataOut.ElecTempFinal[0,id_aux:]=missing
5336 dataOut.EElecTempFinal[0,id_aux:]=missing
5337 dataOut.IonTempFinal[0,id_aux:]=missing
5338 dataOut.EIonTempFinal[0,id_aux:]=missing
5339 dataOut.PhyFinal[0,id_aux:]=missing
5340 dataOut.EPhyFinal[0,id_aux:]=missing
5341
5342 if (time_text.hour == 13 and time_text.minute == 57): #Year: 2022, DOY:241
5343 id_aux = 36
5344 dataOut.DensityFinal[0,id_aux:]=missing
5345 dataOut.EDensityFinal[0,id_aux:]=missing
5346 dataOut.ElecTempFinal[0,id_aux:]=missing
5347 dataOut.EElecTempFinal[0,id_aux:]=missing
5348 dataOut.IonTempFinal[0,id_aux:]=missing
5349 dataOut.EIonTempFinal[0,id_aux:]=missing
5350 dataOut.PhyFinal[0,id_aux:]=missing
5351 dataOut.EPhyFinal[0,id_aux:]=missing
5352 '''
5353 #'''
5354 #print("den: ", dataOut.DensityFinal[0,27])
5355 if (time_text.hour == 5 and time_text.minute == 42): #Year: 2022, DOY:242
5356 id_aux = 16
5357 dataOut.DensityFinal[0,:id_aux]=missing
5358 dataOut.EDensityFinal[0,:id_aux]=missing
5359 dataOut.ElecTempFinal[0,:id_aux]=missing
5360 dataOut.EElecTempFinal[0,:id_aux]=missing
5361 dataOut.IonTempFinal[0,:id_aux]=missing
5362 dataOut.EIonTempFinal[0,:id_aux]=missing
5363 dataOut.PhyFinal[0,:id_aux]=missing
5364 dataOut.EPhyFinal[0,:id_aux]=missing
5365 if (time_text.hour == 5 and time_text.minute == 53): #Year: 2022, DOY:242
5366 id_aux = 9
5367 dataOut.DensityFinal[0,id_aux]=missing
5368 dataOut.EDensityFinal[0,id_aux]=missing
5369 dataOut.ElecTempFinal[0,id_aux]=missing
5370 dataOut.EElecTempFinal[0,id_aux]=missing
5371 dataOut.IonTempFinal[0,id_aux]=missing
5372 dataOut.EIonTempFinal[0,id_aux]=missing
5373 dataOut.PhyFinal[0,id_aux]=missing
5374 dataOut.EPhyFinal[0,id_aux]=missing
5375 if (time_text.hour == 6): #Year: 2022, DOY:242
5376 id_aux = 9
5377 dataOut.DensityFinal[0,id_aux]=missing
5378 dataOut.EDensityFinal[0,id_aux]=missing
5379 dataOut.ElecTempFinal[0,id_aux]=missing
5380 dataOut.EElecTempFinal[0,id_aux]=missing
5381 dataOut.IonTempFinal[0,id_aux]=missing
5382 dataOut.EIonTempFinal[0,id_aux]=missing
5383 dataOut.PhyFinal[0,id_aux]=missing
5384 dataOut.EPhyFinal[0,id_aux]=missing
5385 if (time_text.hour == 6 and time_text.minute == 36): #Year: 2022, DOY:242
5386 id_aux = (10,36,37)
5387 dataOut.DensityFinal[0,id_aux]=missing
5388 dataOut.EDensityFinal[0,id_aux]=missing
5389 dataOut.ElecTempFinal[0,id_aux]=missing
5390 dataOut.EElecTempFinal[0,id_aux]=missing
5391 dataOut.IonTempFinal[0,id_aux]=missing
5392 dataOut.EIonTempFinal[0,id_aux]=missing
5393 dataOut.PhyFinal[0,id_aux]=missing
5394 dataOut.EPhyFinal[0,id_aux]=missing
5395 if (time_text.hour == 7): #Year: 2022, DOY:242
5396 id_aux = 9
5397 dataOut.DensityFinal[0,id_aux]=missing
5398 dataOut.EDensityFinal[0,id_aux]=missing
5399 dataOut.ElecTempFinal[0,id_aux]=missing
5400 dataOut.EElecTempFinal[0,id_aux]=missing
5401 dataOut.IonTempFinal[0,id_aux]=missing
5402 dataOut.EIonTempFinal[0,id_aux]=missing
5403 dataOut.PhyFinal[0,id_aux]=missing
5404 dataOut.EPhyFinal[0,id_aux]=missing
5405 if (time_text.hour == 13 and time_text.minute == 32): #Year: 2022, DOY:242
5406 id_aux = (36,37)
5407 dataOut.DensityFinal[0,id_aux]=missing
5408 dataOut.EDensityFinal[0,id_aux]=missing
5409 dataOut.ElecTempFinal[0,id_aux]=missing
5410 dataOut.EElecTempFinal[0,id_aux]=missing
5411 dataOut.IonTempFinal[0,id_aux]=missing
5412 dataOut.EIonTempFinal[0,id_aux]=missing
5413 dataOut.PhyFinal[0,id_aux]=missing
5414 dataOut.EPhyFinal[0,id_aux]=missing
5415 if (time_text.hour == 23 or time_text.hour <= 4): #Year: 2022, DOY:242
5416 id_aux = 15
5417 dataOut.DensityFinal[0,:id_aux]=missing
5418 dataOut.EDensityFinal[0,:id_aux]=missing
5419 dataOut.ElecTempFinal[0,:id_aux]=missing
5420 dataOut.EElecTempFinal[0,:id_aux]=missing
5421 dataOut.IonTempFinal[0,:id_aux]=missing
5422 dataOut.EIonTempFinal[0,:id_aux]=missing
5423 dataOut.PhyFinal[0,:id_aux]=missing
5424 dataOut.EPhyFinal[0,:id_aux]=missing
5425 if (time_text.hour == 3 and time_text.minute == 13): #Year: 2022, DOY:242
5426 id_aux = (37,38)
5427 dataOut.DensityFinal[0,id_aux]=missing
5428 dataOut.EDensityFinal[0,id_aux]=missing
5429 dataOut.ElecTempFinal[0,id_aux]=missing
5430 dataOut.EElecTempFinal[0,id_aux]=missing
5431 dataOut.IonTempFinal[0,id_aux]=missing
5432 dataOut.EIonTempFinal[0,id_aux]=missing
5433 dataOut.PhyFinal[0,id_aux]=missing
5434 dataOut.EPhyFinal[0,id_aux]=missing
5435 if (time_text.hour == 3 and time_text.minute == 34): #Year: 2022, DOY:242
5436 id_aux = (35,36)
5437 dataOut.DensityFinal[0,id_aux]=missing
5438 dataOut.EDensityFinal[0,id_aux]=missing
5439 dataOut.ElecTempFinal[0,id_aux]=missing
5440 dataOut.EElecTempFinal[0,id_aux]=missing
5441 dataOut.IonTempFinal[0,id_aux]=missing
5442 dataOut.EIonTempFinal[0,id_aux]=missing
5443 dataOut.PhyFinal[0,id_aux]=missing
5444 dataOut.EPhyFinal[0,id_aux]=missing
5445 if (time_text.hour == 4 and time_text.minute == 17): #Year: 2022, DOY:242
5446 id_aux = (34,35)
5447 dataOut.DensityFinal[0,id_aux]=missing
5448 dataOut.EDensityFinal[0,id_aux]=missing
5449 dataOut.ElecTempFinal[0,id_aux]=missing
5450 dataOut.EElecTempFinal[0,id_aux]=missing
5451 dataOut.IonTempFinal[0,id_aux]=missing
5452 dataOut.EIonTempFinal[0,id_aux]=missing
5453 dataOut.PhyFinal[0,id_aux]=missing
5454 dataOut.EPhyFinal[0,id_aux]=missing
5455 if (time_text.hour == 18 and time_text.minute == 30): #Year: 2022, DOY:242
5456 id_aux = (26,27)
5457 dataOut.DensityFinal[0,id_aux]=missing
5458 dataOut.EDensityFinal[0,id_aux]=missing
5459 dataOut.ElecTempFinal[0,id_aux]=missing
5460 dataOut.EElecTempFinal[0,id_aux]=missing
5461 dataOut.IonTempFinal[0,id_aux]=missing
5462 dataOut.EIonTempFinal[0,id_aux]=missing
5463 dataOut.PhyFinal[0,id_aux]=missing
5464 dataOut.EPhyFinal[0,id_aux]=missing
5465 if (time_text.hour == 14 and time_text.minute == 14): #Year: 2022, DOY:242
5466 id_aux = (35,36)
5467 dataOut.DensityFinal[0,id_aux]=missing
5468 dataOut.EDensityFinal[0,id_aux]=missing
5469 dataOut.ElecTempFinal[0,id_aux]=missing
5470 dataOut.EElecTempFinal[0,id_aux]=missing
5471 dataOut.IonTempFinal[0,id_aux]=missing
5472 dataOut.EIonTempFinal[0,id_aux]=missing
5473 dataOut.PhyFinal[0,id_aux]=missing
5474 dataOut.EPhyFinal[0,id_aux]=missing
4908 #'''
5475 #'''
4909 #print("den_final",dataOut.DensityFinal)
5476 #print("den_final",dataOut.DensityFinal)
4910
5477
5478
4911 dataOut.flagNoData = numpy.all(numpy.isnan(dataOut.DensityFinal)) #Si todos los valores son NaN no se prosigue
5479 dataOut.flagNoData = numpy.all(numpy.isnan(dataOut.DensityFinal)) #Si todos los valores son NaN no se prosigue
4912
5480
4913 ####dataOut.flagNoData = False #Solo para ploteo
5481 ####dataOut.flagNoData = False #Solo para ploteo
4914
5482
4915 dataOut.DensityFinal *= 1.e6 #Convert units to m^⁻3
5483 dataOut.DensityFinal *= 1.e6 #Convert units to m^⁻3
4916 dataOut.EDensityFinal *= 1.e6 #Convert units to m^⁻3
5484 dataOut.EDensityFinal *= 1.e6 #Convert units to m^⁻3
4917 #print(dataOut.flagNoData)
5485 print("Save Cleaner: ", dataOut.flagNoData)
5486 #print("den: ", dataOut.DensityFinal[0,27])
4918 return dataOut
5487 return dataOut
4919
5488
4920
5489
@@ -5704,6 +6273,8 class SSheightProfiles(Operation):
5704 #exit(1)
6273 #exit(1)
5705 #print(dataOut.data[0,:,150])
6274 #print(dataOut.data[0,:,150])
5706 #exit(1)
6275 #exit(1)
6276 #print(dataOut.data[0,:,0]*numpy.conjugate(dataOut.data[0,0,0]))
6277 #exit(1)
5707
6278
5708 return dataOut
6279 return dataOut
5709
6280
@@ -6905,6 +7476,7 class CrossProdHybrid(CrossProdDP):
6905 #self.dataOut.nptsfft2=150
7476 #self.dataOut.nptsfft2=150
6906 self.cnorm=float((dataOut.nProfiles-dataOut.NSCAN)/dataOut.NSCAN)
7477 self.cnorm=float((dataOut.nProfiles-dataOut.NSCAN)/dataOut.NSCAN)
6907 self.lagp0=numpy.zeros((dataOut.NLAG,dataOut.NRANGE,dataOut.NAVG),'complex128')
7478 self.lagp0=numpy.zeros((dataOut.NLAG,dataOut.NRANGE,dataOut.NAVG),'complex128')
7479 ww=numpy.zeros((dataOut.NLAG,dataOut.NRANGE,dataOut.NSCAN,dataOut.NAVG),'complex128')
6908 self.lagp1=numpy.zeros((dataOut.NLAG,dataOut.NRANGE,dataOut.NAVG),'complex128')
7480 self.lagp1=numpy.zeros((dataOut.NLAG,dataOut.NRANGE,dataOut.NAVG),'complex128')
6909 self.lagp2=numpy.zeros((dataOut.NLAG,dataOut.NRANGE,dataOut.NAVG),'complex128')
7481 self.lagp2=numpy.zeros((dataOut.NLAG,dataOut.NRANGE,dataOut.NAVG),'complex128')
6910 self.lagp3=numpy.zeros((dataOut.NLAG,dataOut.NRANGE,dataOut.NAVG),'complex128')
7482 self.lagp3=numpy.zeros((dataOut.NLAG,dataOut.NRANGE,dataOut.NAVG),'complex128')
@@ -6927,6 +7499,7 class CrossProdHybrid(CrossProdDP):
6927 #exit(1)
7499 #exit(1)
6928 if i==0:
7500 if i==0:
6929 self.lagp0[n][j][self.bcounter-1]=numpy.sum(c[:dataOut.NSCAN])
7501 self.lagp0[n][j][self.bcounter-1]=numpy.sum(c[:dataOut.NSCAN])
7502 ww[n,j,:,self.bcounter-1]=c[:dataOut.NSCAN]
6930 self.lagp3[n][j][self.bcounter-1]=numpy.sum(c[dataOut.NSCAN:]/self.cnorm)
7503 self.lagp3[n][j][self.bcounter-1]=numpy.sum(c[dataOut.NSCAN:]/self.cnorm)
6931 elif i==1:
7504 elif i==1:
6932 self.lagp1[n][j][self.bcounter-1]=numpy.sum(c[:dataOut.NSCAN])
7505 self.lagp1[n][j][self.bcounter-1]=numpy.sum(c[:dataOut.NSCAN])
@@ -6945,6 +7518,9 class CrossProdHybrid(CrossProdDP):
6945 #print(sum(self.buffer[3,:,199,2]))
7518 #print(sum(self.buffer[3,:,199,2]))
6946 #print(self.cnorm)
7519 #print(self.cnorm)
6947 #exit(1)
7520 #exit(1)
7521 #print("self,lagp0: ", self.lagp0[0,0,self.bcounter-1])
7522 print(ww[:,0,0,self.bcounter-1])
7523 exit(1)
6948
7524
6949
7525
6950 def LP_median_estimates_original(self,dataOut):
7526 def LP_median_estimates_original(self,dataOut):
@@ -6955,7 +7531,7 class CrossProdHybrid(CrossProdDP):
6955 self.output=numpy.zeros((dataOut.NLAG,dataOut.NRANGE,dataOut.NR),'complex128')
7531 self.output=numpy.zeros((dataOut.NLAG,dataOut.NRANGE,dataOut.NR),'complex128')
6956 self.lag_products_LP_median_estimates_aux=0
7532 self.lag_products_LP_median_estimates_aux=0
6957
7533
6958
7534 #print("self,lagp0: ", numpy.sum(self.lagp0[0,0,:]))
6959 for i in range(dataOut.NLAG):
7535 for i in range(dataOut.NLAG):
6960 for j in range(dataOut.NRANGE):
7536 for j in range(dataOut.NRANGE):
6961 for l in range(4): #four outputs
7537 for l in range(4): #four outputs
@@ -6995,7 +7571,7 class CrossProdHybrid(CrossProdDP):
6995 if self.lag_products_LP_median_estimates_aux==1:
7571 if self.lag_products_LP_median_estimates_aux==1:
6996 self.output=numpy.zeros((dataOut.NLAG,dataOut.NRANGE,dataOut.NR),'complex128')
7572 self.output=numpy.zeros((dataOut.NLAG,dataOut.NRANGE,dataOut.NR),'complex128')
6997 self.lag_products_LP_median_estimates_aux=0
7573 self.lag_products_LP_median_estimates_aux=0
6998
7574 #print("self,lagp0: ", numpy.sum(self.lagp0[0,0,:]))
6999
7575
7000 for i in range(dataOut.NLAG):
7576 for i in range(dataOut.NLAG):
7001 #my_list = ([0,1,2,3,4,5,6,7]) #hasta 7 funciona, en 6 ya no
7577 #my_list = ([0,1,2,3,4,5,6,7]) #hasta 7 funciona, en 6 ya no
@@ -7388,12 +7964,12 class LongPulseAnalysis(Operation):
7388 self.aux=0
7964 self.aux=0
7389
7965
7390 dataOut.cut=30
7966 dataOut.cut=30
7391 for i in range(30,15,-1):
7967 for i in range(30,15,-1): #AquΓ­ se calcula en donde se unirΓ‘ DP y LP en la parte final
7392 if numpy.nanmax(dataOut.acfs_error_to_plot[i,:])>=10 or dataOut.info2[i]==0:
7968 if numpy.nanmax(dataOut.acfs_error_to_plot[i,:])>=10 or dataOut.info2[i]==0:
7393 dataOut.cut=i-1
7969 dataOut.cut=i-1
7394
7970
7395 for i in range(dataOut.NLAG):
7971 for i in range(dataOut.NLAG):
7396 self.cal[i]=sum(dataOut.output_LP_integrated[i,:,3].real)
7972 self.cal[i]=sum(dataOut.output_LP_integrated[i,:,3].real) #Lag x Height x Channel
7397
7973
7398 #print(numpy.sum(self.cal)) #Coinciden
7974 #print(numpy.sum(self.cal)) #Coinciden
7399 #exit(1)
7975 #exit(1)
@@ -7401,7 +7977,18 class LongPulseAnalysis(Operation):
7401 #print(anoise0)
7977 #print(anoise0)
7402 #print(anoise1)
7978 #print(anoise1)
7403 #exit(1)
7979 #exit(1)
7404
7980 print("nis: ", dataOut.nis)
7981 print("pan: ", dataOut.pan)
7982 print("pbn: ", dataOut.pbn)
7983 #print(numpy.sum(dataOut.output_LP_integrated[0,:,0]))
7984 '''
7985 import matplotlib.pyplot as plt
7986 plt.plot(dataOut.output_LP_integrated[:,40,0])
7987 plt.show()
7988 '''
7989 #print(dataOut.output_LP_integrated[0,40,0])
7990 print(numpy.sum(dataOut.output_LP_integrated[:,0,0]))
7991 exit(1)
7405
7992
7406 #################### PROBAR MÁS INTEGRACIΓ“N, SINO MODIFICAR VALOR DE "NIS" ####################
7993 #################### PROBAR MÁS INTEGRACIΓ“N, SINO MODIFICAR VALOR DE "NIS" ####################
7407 # VER dataOut.nProfiles_LP #
7994 # VER dataOut.nProfiles_LP #
@@ -7427,7 +8014,7 class LongPulseAnalysis(Operation):
7427 for i in range(1,dataOut.NLAG): #remove cal data from certain lags
8014 for i in range(1,dataOut.NLAG): #remove cal data from certain lags
7428 dataOut.output_LP_integrated.real[i,j,0]-=self.cal[i]
8015 dataOut.output_LP_integrated.real[i,j,0]-=self.cal[i]
7429 k=max(j,26) #constant power below range 26
8016 k=max(j,26) #constant power below range 26
7430 self.powera[j]=dataOut.output_LP_integrated.real[0,k,0]
8017 self.powera[j]=dataOut.output_LP_integrated.real[0,k,0] #Lag0 and Channel 0
7431
8018
7432 ## examine drifts here - based on 60 'indep.' estimates
8019 ## examine drifts here - based on 60 'indep.' estimates
7433 #print(numpy.sum(self.powera))
8020 #print(numpy.sum(self.powera))
@@ -7520,7 +8107,7 class LongPulseAnalysis(Operation):
7520 self.perror[:range2_nnls]=1.00/self.perror[:range2_nnls]
8107 self.perror[:range2_nnls]=1.00/self.perror[:range2_nnls]
7521
8108
7522 b_nnlswrap=numpy.zeros(range2_nnls,'float64')
8109 b_nnlswrap=numpy.zeros(range2_nnls,'float64')
7523 b_nnlswrap[:]=numpy.matmul(self.powera[dataOut.IBITS+ioff:range1_nnls+dataOut.IBITS+ioff],g_nnlswrap)
8110 b_nnlswrap[:]=numpy.matmul(self.powera[dataOut.IBITS+ioff:range1_nnls+dataOut.IBITS+ioff],g_nnlswrap) #match filter alturas
7524
8111
7525 x_nnlswrap=numpy.zeros(range2_nnls,'float64')
8112 x_nnlswrap=numpy.zeros(range2_nnls,'float64')
7526 x_nnlswrap[:]=nnls(a_nnlswrap,b_nnlswrap)[0]
8113 x_nnlswrap[:]=nnls(a_nnlswrap,b_nnlswrap)[0]
@@ -7630,6 +8217,8 class LongPulseAnalysis(Operation):
7630 exit(1)
8217 exit(1)
7631 '''
8218 '''
7632 print("Success")
8219 print("Success")
8220 ###################Correlation pulse and itself
8221
7633 #print(dataOut.NRANGE)
8222 #print(dataOut.NRANGE)
7634 with suppress_stdout_stderr():
8223 with suppress_stdout_stderr():
7635 #pass
8224 #pass
@@ -7649,6 +8238,303 class LongPulseAnalysis(Operation):
7649
8238
7650 return dataOut
8239 return dataOut
7651
8240
8241 class LongPulseAnalysisSpectra(Operation):
8242 """Operation to estimate ACFs, temperatures, total electron density and Hydrogen/Helium fractions from the Long Pulse data.
8243
8244 Parameters:
8245 -----------
8246 NACF : int
8247 .*
8248
8249 Example
8250 --------
8251
8252 op = proc_unit.addOperation(name='LongPulseAnalysis', optype='other')
8253 op.addParameter(name='NACF', value='16', format='int')
8254
8255 """
8256
8257 def __init__(self, **kwargs):
8258
8259 Operation.__init__(self, **kwargs)
8260 self.aux=1
8261
8262 def run(self,dataOut,NACF):
8263
8264 dataOut.NACF=NACF
8265 dataOut.heightList=dataOut.DH*(numpy.arange(dataOut.NACF))
8266 anoise0=dataOut.tnoise[0]
8267 anoise1=anoise0*0.0 #seems to be noise in 1st lag 0.015 before '14
8268 #print(anoise0)
8269 #exit(1)
8270 if self.aux:
8271 #dataOut.cut=31#26#height=31*15=465
8272 self.cal=numpy.zeros((dataOut.NLAG),'float32')
8273 self.drift=numpy.zeros((200),'float32')
8274 self.rdrift=numpy.zeros((200),'float32')
8275 self.ddrift=numpy.zeros((200),'float32')
8276 self.sigma=numpy.zeros((dataOut.NRANGE),order='F',dtype='float32')
8277 self.powera=numpy.zeros((dataOut.NRANGE),order='F',dtype='float32')
8278 self.powerb=numpy.zeros((dataOut.NRANGE),order='F',dtype='float32')
8279 self.perror=numpy.zeros((dataOut.NRANGE),order='F',dtype='float32')
8280 dataOut.ene=numpy.zeros((dataOut.NRANGE),'float32')
8281 self.dpulse=numpy.zeros((dataOut.NACF),'float32')
8282 self.lpulse=numpy.zeros((dataOut.NACF),'float32')
8283 dataOut.lags_LP=numpy.zeros((dataOut.IBITS),order='F',dtype='float32')
8284 self.lagp=numpy.zeros((dataOut.NACF),'float32')
8285 self.u=numpy.zeros((2*dataOut.NACF,2*dataOut.NACF),'float32')
8286 dataOut.ne=numpy.zeros((dataOut.NRANGE),order='F',dtype='float32')
8287 dataOut.te=numpy.zeros((dataOut.NACF),order='F',dtype='float32')
8288 dataOut.ete=numpy.zeros((dataOut.NACF),order='F',dtype='float32')
8289 dataOut.ti=numpy.zeros((dataOut.NACF),order='F',dtype='float32')
8290 dataOut.eti=numpy.zeros((dataOut.NACF),order='F',dtype='float32')
8291 dataOut.ph=numpy.zeros((dataOut.NACF),order='F',dtype='float32')
8292 dataOut.eph=numpy.zeros((dataOut.NACF),order='F',dtype='float32')
8293 dataOut.phe=numpy.zeros((dataOut.NACF),order='F',dtype='float32')
8294 dataOut.ephe=numpy.zeros((dataOut.NACF),order='F',dtype='float32')
8295 dataOut.errors=numpy.zeros((dataOut.IBITS,max(dataOut.NRANGE,dataOut.NSHTS)),order='F',dtype='float32')
8296 dataOut.fit_array_real=numpy.zeros((max(dataOut.NRANGE,dataOut.NSHTS),dataOut.NLAG),order='F',dtype='float32')
8297 dataOut.status=numpy.zeros(1,'float32')
8298 dataOut.tx=240.0 #deberΓ­a provenir del header #hybrid
8299
8300 for i in range(dataOut.IBITS):
8301 dataOut.lags_LP[i]=float(i)*(dataOut.tx/150.0)/float(dataOut.IBITS) # (float)i*(header.tx/150.0)/(float)IBITS;
8302
8303 self.aux=0
8304
8305 dataOut.cut=30
8306 for i in range(30,15,-1): #AquΓ­ se calcula en donde se unirΓ‘ DP y LP en la parte final
8307 if numpy.nanmax(dataOut.acfs_error_to_plot[i,:])>=10 or dataOut.info2[i]==0:
8308 dataOut.cut=i-1
8309
8310 for i in range(dataOut.NLAG):
8311 self.cal[i]=sum(dataOut.output_LP_integrated[i,:,3].real) #Lag x Height x Channel
8312
8313 #print(numpy.sum(self.cal)) #Coinciden
8314 #exit(1)
8315 self.cal/=float(dataOut.NRANGE)
8316
8317
8318 #################### PROBAR MÁS INTEGRACIΓ“N, SINO MODIFICAR VALOR DE "NIS" ####################
8319 # VER dataOut.nProfiles_LP #
8320
8321 '''
8322 #PLOTEAR POTENCIA VS RUIDO, QUIZA SE ESTA REMOVIENDO MUCHA SEΓ‘AL
8323 #print(dataOut.heightList)
8324 import matplotlib.pyplot as plt
8325 plt.plot(10*numpy.log10(dataOut.output_LP_integrated.real[0,:,0]),dataOut.range1)
8326 #plt.plot(10*numpy.log10(dataOut.output_LP_integrated.real[0,:,0]/dataOut.nProfiles_LP),dataOut.range1)
8327 plt.axvline(10*numpy.log10(anoise0),color='k',linestyle='dashed')
8328 plt.grid()
8329 plt.xlim(20,100)
8330 plt.show()
8331 '''
8332
8333
8334 for j in range(dataOut.NACF+2*dataOut.IBITS+2):
8335
8336 dataOut.output_LP_integrated.real[0,j,0]-=anoise0 #lag0 ch0
8337 dataOut.output_LP_integrated.real[1,j,0]-=anoise1 #lag1 ch0
8338
8339 for i in range(1,dataOut.NLAG): #remove cal data from certain lags
8340 dataOut.output_LP_integrated.real[i,j,0]-=self.cal[i]
8341 k=max(j,26) #constant power below range 26
8342 self.powera[j]=dataOut.output_LP_integrated.real[0,k,0] #Lag0 and Channel 0
8343
8344 ## examine drifts here - based on 60 'indep.' estimates
8345 #print(numpy.sum(self.powera))
8346 #exit(1)
8347 #nis=dataOut.NSCAN*dataOut.NAVG*dataOut.nint*10
8348 nis = dataOut.nis
8349 #print("nis",nis)
8350 alpha=beta=delta=0.0
8351 nest=0
8352 gamma=3.0/(2.0*numpy.pi*dataOut.lags_LP[1]*1.0e-3)
8353 beta=gamma*(math.atan2(dataOut.output_LP_integrated.imag[14,0,2],dataOut.output_LP_integrated.real[14,0,2])-math.atan2(dataOut.output_LP_integrated.imag[1,0,2],dataOut.output_LP_integrated.real[1,0,2]))/13.0
8354 #print(gamma,beta)
8355 #exit(1)
8356 for i in range(1,3):
8357 gamma=3.0/(2.0*numpy.pi*dataOut.lags_LP[i]*1.0e-3)
8358 #print("gamma",gamma)
8359 for j in range(34,44):
8360 rho2=numpy.abs(dataOut.output_LP_integrated[i,j,0])/numpy.abs(dataOut.output_LP_integrated[0,j,0])
8361 dataOut.dphi2=(1.0/rho2-1.0)/(float(2*nis))
8362 dataOut.dphi2*=gamma**2
8363 pest=gamma*math.atan(dataOut.output_LP_integrated.imag[i,j,0]/dataOut.output_LP_integrated.real[i,j,0])
8364 #print("1",dataOut.output_LP_integrated.imag[i,j,0])
8365 #print("2",dataOut.output_LP_integrated.real[i,j,0])
8366 self.drift[nest]=pest
8367 self.ddrift[nest]=dataOut.dphi2
8368 self.rdrift[nest]=float(nest)
8369 nest+=1
8370
8371 sorted(self.drift[:nest])
8372
8373 #print(dataOut.dphi2)
8374 #exit(1)
8375
8376 for j in range(int(nest/4),int(3*nest/4)):
8377 #i=int(self.rdrift[j])
8378 alpha+=self.drift[j]/self.ddrift[j]
8379 delta+=1.0/self.ddrift[j]
8380
8381 alpha/=delta
8382 delta=1./numpy.sqrt(delta)
8383 vdrift=alpha-beta
8384 dvdrift=delta
8385
8386 #need to develop estimate of complete density profile using all
8387 #available data
8388
8389 #estimate sample variances for long-pulse power profile
8390
8391 #nis=dataOut.NSCAN*dataOut.NAVG*dataOut.nint
8392 nis = dataOut.nis/10
8393 #print("nis",nis)
8394
8395 self.sigma[:dataOut.NACF+2*dataOut.IBITS+2]=((anoise0+self.powera[:dataOut.NACF+2*dataOut.IBITS+2])**2)/float(nis)
8396 #print(self.sigma)
8397 #exit(1)
8398 ioff=1
8399
8400 #deconvolve rectangular pulse shape from profile ==> powerb, perror
8401
8402 '''
8403 ############# START nnlswrap#############
8404
8405 if dataOut.ut_Faraday>14.0:
8406 alpha_nnlswrap=20.0
8407 else:
8408 alpha_nnlswrap=30.0
8409
8410 range1_nnls=dataOut.NACF
8411 range2_nnls=dataOut.NACF+dataOut.IBITS-1
8412
8413 g_nnlswrap=numpy.zeros((range1_nnls,range2_nnls),'float32')
8414 a_nnlswrap=numpy.zeros((range2_nnls,range2_nnls),'float64')
8415
8416 for i in range(range1_nnls):
8417 for j in range(range2_nnls):
8418 if j>=i and j<i+dataOut.IBITS:
8419 g_nnlswrap[i,j]=1.0
8420 else:
8421 g_nnlswrap[i,j]=0.0
8422
8423 a_nnlswrap[:]=numpy.matmul(numpy.transpose(g_nnlswrap),g_nnlswrap)
8424
8425 numpy.fill_diagonal(a_nnlswrap,a_nnlswrap.diagonal()+alpha_nnlswrap**2)
8426
8427 #ERROR ANALYSIS#
8428
8429 self.perror[:range2_nnls]=0.0
8430 self.perror[:range2_nnls]=numpy.matmul(1./(self.sigma[dataOut.IBITS+ioff:range1_nnls+dataOut.IBITS+ioff]),g_nnlswrap**2)
8431 self.perror[:range1_nnls]+=(alpha_nnlswrap**2)/(self.sigma[dataOut.IBITS+ioff:range1_nnls+dataOut.IBITS+ioff])
8432 self.perror[:range2_nnls]=1.00/self.perror[:range2_nnls]
8433
8434 b_nnlswrap=numpy.zeros(range2_nnls,'float64')
8435 b_nnlswrap[:]=numpy.matmul(self.powera[dataOut.IBITS+ioff:range1_nnls+dataOut.IBITS+ioff],g_nnlswrap)
8436
8437 x_nnlswrap=numpy.zeros(range2_nnls,'float64')
8438 x_nnlswrap[:]=nnls(a_nnlswrap,b_nnlswrap)[0]
8439
8440 self.powerb[:range2_nnls]=x_nnlswrap
8441 #print(self.powerb[40])
8442 #print(self.powerb[66])
8443 #exit(1)
8444 #############END nnlswrap#############
8445 '''
8446 self.powerb[:] = self.powera
8447 self.perror[:] = 0.
8448 #print(numpy.sum(numpy.sqrt(self.perror[0:dataOut.NACF])))
8449 #print(self.powerb[0:dataOut.NACF])
8450 #exit(1)
8451 #estimate relative error for deconvolved profile (scaling irrelevant)
8452 #print(dataOut.NACF)
8453 dataOut.ene[0:dataOut.NACF]=numpy.sqrt(self.perror[0:dataOut.NACF])/self.powerb[0:dataOut.NACF]
8454 #print(numpy.sum(dataOut.ene))
8455 #exit(1)
8456 aux=0
8457
8458 for i in range(dataOut.IBITS,dataOut.NACF):
8459 self.dpulse[i]=self.lpulse[i]=0.0
8460 for j in range(dataOut.IBITS):
8461 k=int(i-j)
8462 if k<36-aux and k>16:
8463 self.dpulse[i]+=dataOut.ph2[k]/dataOut.h2[k]
8464 elif k>=36-aux:
8465 self.lpulse[i]+=self.powerb[k]
8466 self.lagp[i]=self.powera[i]
8467
8468 #find scale factor that best merges profiles
8469
8470 qi=sum(self.dpulse[32:dataOut.NACF]**2/(self.lagp[32:dataOut.NACF]+anoise0)**2)
8471 ri=sum((self.dpulse[32:dataOut.NACF]*self.lpulse[32:dataOut.NACF])/(self.lagp[32:dataOut.NACF]+anoise0)**2)
8472 si=sum((self.dpulse[32:dataOut.NACF]*self.lagp[32:dataOut.NACF])/(self.lagp[32:dataOut.NACF]+anoise0)**2)
8473 ui=sum(self.lpulse[32:dataOut.NACF]**2/(self.lagp[32:dataOut.NACF]+anoise0)**2)
8474 vi=sum((self.lpulse[32:dataOut.NACF]*self.lagp[32:dataOut.NACF])/(self.lagp[32:dataOut.NACF]+anoise0)**2)
8475
8476 alpha=(si*ui-vi*ri)/(qi*ui-ri*ri)
8477 beta=(qi*vi-ri*si)/(qi*ui-ri*ri)
8478
8479 #form density profile estimate, merging rescaled power profiles
8480 #print(dataOut.h2)
8481 #print(numpy.sum(alpha))
8482 #print(numpy.sum(dataOut.ph2))
8483 self.powerb[16:36-aux]=alpha*dataOut.ph2[16:36-aux]/dataOut.h2[16:36-aux]
8484 self.powerb[36-aux:dataOut.NACF]*=beta
8485
8486 #form Ne estimate, fill in error estimate at low altitudes
8487
8488 dataOut.ene[0:36-aux]=dataOut.sdp2[0:36-aux]/dataOut.ph2[0:36-aux]
8489 dataOut.ne[:dataOut.NACF]=self.powerb[:dataOut.NACF]*dataOut.h2[:dataOut.NACF]/alpha
8490 #print(numpy.sum(self.powerb))
8491 #print(numpy.sum(dataOut.ene))
8492 #print(numpy.sum(dataOut.ne))
8493 #exit(1)
8494 #now do error propagation: store zero lag error covariance in u
8495
8496 nis=dataOut.NSCAN*dataOut.NAVG*dataOut.nint/1 # DLH serious debris removal
8497
8498 for i in range(dataOut.NACF):
8499 for j in range(i,dataOut.NACF):
8500 if j-i>=dataOut.IBITS:
8501 self.u[i,j]=0.0
8502 else:
8503 self.u[i,j]=dataOut.output_LP_integrated.real[j-i,i,0]**2/float(nis)
8504 self.u[i,j]*=(anoise0+dataOut.output_LP_integrated.real[0,i,0])/dataOut.output_LP_integrated.real[0,i,0]
8505 self.u[i,j]*=(anoise0+dataOut.output_LP_integrated.real[0,j,0])/dataOut.output_LP_integrated.real[0,j,0]
8506
8507 self.u[j,i]=self.u[i,j]
8508
8509 #now error analyis for lag product matrix (diag), place in acf_err
8510
8511 for i in range(dataOut.NACF):
8512 for j in range(dataOut.IBITS):
8513 if j==0:
8514 dataOut.errors[0,i]=numpy.sqrt(self.u[i,i])
8515 else:
8516 dataOut.errors[j,i]=numpy.sqrt(((dataOut.output_LP_integrated.real[0,i,0]+anoise0)*(dataOut.output_LP_integrated.real[0,i+j,0]+anoise0)+dataOut.output_LP_integrated.real[j,i,0]**2)/float(2*nis))
8517
8518 print("Success")
8519 #print(dataOut.NRANGE)
8520 with suppress_stdout_stderr():
8521 pass
8522 #full_profile_profile.profile(numpy.transpose(dataOut.output_LP_integrated,(2,1,0)),numpy.transpose(dataOut.errors),self.powerb,dataOut.ne,dataOut.lags_LP,dataOut.thb,dataOut.bfm,dataOut.te,dataOut.ete,dataOut.ti,dataOut.eti,dataOut.ph,dataOut.eph,dataOut.phe,dataOut.ephe,dataOut.range1,dataOut.ut,dataOut.NACF,dataOut.fit_array_real,dataOut.status,dataOut.NRANGE,dataOut.IBITS)
8523
8524 print("status: ",dataOut.status)
8525
8526 if dataOut.status>=3.5:
8527 dataOut.te[:]=numpy.nan
8528 dataOut.ete[:]=numpy.nan
8529 dataOut.ti[:]=numpy.nan
8530 dataOut.eti[:]=numpy.nan
8531 dataOut.ph[:]=numpy.nan
8532 dataOut.eph[:]=numpy.nan
8533 dataOut.phe[:]=numpy.nan
8534 dataOut.ephe[:]=numpy.nan
8535
8536 return dataOut
8537
7652 class LongPulseAnalysis_V2(Operation):
8538 class LongPulseAnalysis_V2(Operation):
7653 """Operation to estimate ACFs, temperatures, total electron density and Hydrogen/Helium fractions from the Long Pulse data.
8539 """Operation to estimate ACFs, temperatures, total electron density and Hydrogen/Helium fractions from the Long Pulse data.
7654
8540
General Comments 0
You need to be logged in to leave comments. Login now