diff --git a/schainpy/model/data/jrodata.py b/schainpy/model/data/jrodata.py
index 62752c5..3abae66 100644
--- a/schainpy/model/data/jrodata.py
+++ b/schainpy/model/data/jrodata.py
@@ -637,8 +637,6 @@ class Spectra(JROData):
def getVelRange(self, extrapoints=0):
- print 'VELMAX', self.getVmax()
- asdasdasd
deltav = self.getVmax() / (self.nFFTPoints*self.ippFactor)
velrange = deltav*(numpy.arange(self.nFFTPoints+extrapoints)-self.nFFTPoints/2.) #- deltav/2
diff --git a/schainpy/model/graphics/jroplot_parameters.py b/schainpy/model/graphics/jroplot_parameters.py
index 2691d40..7965653 100644
--- a/schainpy/model/graphics/jroplot_parameters.py
+++ b/schainpy/model/graphics/jroplot_parameters.py
@@ -6,14 +6,14 @@ from figure import Figure, isRealtime, isTimeInHourRange
from plotting_codes import *
-class FitGauPlot(Figure):
+class SpcParamPlot(Figure):
isConfig = None
__nsubplots = None
WIDTHPROF = None
HEIGHTPROF = None
- PREFIX = 'fitgau'
+ PREFIX = 'SpcParam'
def __init__(self, **kwargs):
Figure.__init__(self, **kwargs)
@@ -82,7 +82,7 @@ class FitGauPlot(Figure):
save=False, figpath='./', figfile=None, show=True, ftp=False, wr_period=1,
server=None, folder=None, username=None, password=None,
ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0, realtime=False,
- xaxis="frequency", colormap='jet', normFactor=None , GauSelector = 0):
+ xaxis="frequency", colormap='jet', normFactor=None , Selector = 0):
"""
@@ -118,23 +118,23 @@ class FitGauPlot(Figure):
# else:
# factor = normFactor
if xaxis == "frequency":
- x = dataOut.spc_range[0]
+ x = dataOut.spcparam_range[0]
xlabel = "Frequency (kHz)"
elif xaxis == "time":
- x = dataOut.spc_range[1]
+ x = dataOut.spcparam_range[1]
xlabel = "Time (ms)"
else:
- x = dataOut.spc_range[2]
+ x = dataOut.spcparam_range[2]
xlabel = "Velocity (m/s)"
ylabel = "Range (Km)"
y = dataOut.getHeiRange()
- z = dataOut.GauSPC[:,GauSelector,:,:] #GauSelector] #dataOut.data_spc/factor
- print 'GausSPC', z[0,32,10:40]
+ z = dataOut.SPCparam[Selector] #GauSelector] #dataOut.data_spc/factor
+ #print 'GausSPC', z[0,32,10:40]
z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
zdB = 10*numpy.log10(z)
diff --git a/schainpy/model/graphics/jroplot_spectra.py b/schainpy/model/graphics/jroplot_spectra.py
index a16916d..657958d 100644
--- a/schainpy/model/graphics/jroplot_spectra.py
+++ b/schainpy/model/graphics/jroplot_spectra.py
@@ -128,12 +128,6 @@ class SpectraPlot(Figure):
factor = normFactor
if xaxis == "frequency":
x = dataOut.getFreqRange(1)/1000.
- print 'FRECUENCIA MAXIMA', numpy.amax(x)
- asfasfasdfaf
- print '#######################################################'
- print 'xlen', len(x)
- print x
- print '#######################################################'
xlabel = "Frequency (kHz)"
elif xaxis == "time":
diff --git a/schainpy/model/proc/jroproc_parameters.py b/schainpy/model/proc/jroproc_parameters.py
index 76cde9a..5cbfebe 100644
--- a/schainpy/model/proc/jroproc_parameters.py
+++ b/schainpy/model/proc/jroproc_parameters.py
@@ -28,6 +28,7 @@ from scipy.optimize import curve_fit
import warnings
from numpy import NaN
from scipy.optimize.optimize import OptimizeWarning
+from IPython.parallel.controller.scheduler import numpy
warnings.filterwarnings('ignore')
@@ -122,7 +123,7 @@ class ParametersProc(ProcessingUnit):
print 'self.dataIn.data_spc', self.dataIn.data_spc.shape
self.dataOut.abscissaList = self.dataIn.getVelRange(1)
self.dataOut.spc_noise = self.dataIn.getNoise()
- self.dataOut.spc_range = (self.dataIn.getFreqRange(1) , self.dataIn.getAcfRange(1) , self.dataIn.getVelRange(1) )
+ self.dataOut.spc_range = numpy.asanyarray((self.dataIn.getFreqRange(1) , self.dataIn.getAcfRange(1) , self.dataIn.getVelRange(1) ))
self.dataOut.normFactor = self.dataIn.normFactor
#self.dataOut.outputInterval = self.dataIn.outputInterval
@@ -187,6 +188,117 @@ def target(tups):
return obj.FitGau(args)
+class SpectralFilters(Operation):
+
+ '''This class allows the Rainfall / Wind Selection for CLAIRE RADAR
+
+ LimitR : It is the limit in m/s of Rainfall
+ LimitW : It is the limit in m/s for Winds
+
+ Input:
+
+ self.dataOut.data_pre : SPC and CSPC
+ self.dataOut.spc_range : To select wind and rainfall velocities
+
+ Affected:
+
+ self.dataOut.data_pre : It is used for the new SPC and CSPC ranges of wind
+ self.dataOut.spcparam_range : Used in SpcParamPlot
+ self.dataOut.SPCparam : Used in PrecipitationProc
+
+
+ '''
+
+ def __init__(self, **kwargs):
+ Operation.__init__(self, **kwargs)
+ self.i=0
+
+ def run(self, dataOut, Rain_Velocity_Limit=1.5, Wind_Velocity_Limit=2.5):
+
+ #Limite de vientos
+ LimitR = Rain_Velocity_Limit
+ LimitW = Wind_Velocity_Limit
+
+ self.spc = dataOut.data_pre[0].copy()
+ self.cspc = dataOut.data_pre[1].copy()
+
+ self.Num_Hei = self.spc.shape[2]
+ self.Num_Bin = self.spc.shape[1]
+ self.Num_Chn = self.spc.shape[0]
+
+ VelRange = dataOut.spc_range[2]
+ TimeRange = dataOut.spc_range[1]
+ FrecRange = dataOut.spc_range[0]
+
+ Vmax= 2*numpy.max(dataOut.spc_range[2])
+ Tmax= 2*numpy.max(dataOut.spc_range[1])
+ Fmax= 2*numpy.max(dataOut.spc_range[0])
+
+ Breaker1R=VelRange[numpy.abs(VelRange-(-LimitR)).argmin()]
+ Breaker1R=numpy.where(VelRange == Breaker1R)
+
+ Breaker1W=VelRange[numpy.abs(VelRange-(-LimitW)).argmin()]
+ Breaker1W=numpy.where(VelRange == Breaker1W)
+
+ Breaker2W=VelRange[numpy.abs(VelRange-(LimitW)).argmin()]
+ Breaker2W=numpy.where(VelRange == Breaker2W)
+
+
+ '''Reacomodando SPCrange'''
+
+ VelRange=numpy.roll(VelRange,-Breaker1R[0],axis=0)
+
+ VelRange[-int(Breaker1R[0]):]+= Vmax
+
+ FrecRange=numpy.roll(FrecRange,-Breaker1R[0],axis=0)
+
+ FrecRange[-int(Breaker1R[0]):]+= Fmax
+
+ TimeRange=numpy.roll(TimeRange,-Breaker1R[0],axis=0)
+
+ TimeRange[-int(Breaker1R[0]):]+= Tmax
+
+ ''' ------------------ '''
+
+ Breaker2R=VelRange[numpy.abs(VelRange-(LimitR)).argmin()]
+ Breaker2R=numpy.where(VelRange == Breaker2R)
+
+
+
+
+ SPCroll = numpy.roll(self.spc,-Breaker1R[0],axis=1)
+
+ SPCcut = SPCroll.copy()
+ for i in range(self.Num_Chn):
+ SPCcut[i,0:int(Breaker2R[0]),:] = dataOut.noise[i]
+
+ self.spc[i, 0:int(Breaker1W[0]) ,:] = dataOut.noise[i]
+ self.spc[i, int(Breaker2W[0]):self.Num_Bin ,:] = dataOut.noise[i]
+
+ self.cspc[i, 0:int(Breaker1W[0]) ,:] = dataOut.noise[i]
+ self.cspc[i, int(Breaker2W[0]):self.Num_Bin ,:] = dataOut.noise[i]
+
+
+ SPC_ch1 = SPCroll
+
+ SPC_ch2 = SPCcut
+
+ SPCparam = (SPC_ch1, SPC_ch2, self.spc)
+ dataOut.SPCparam = numpy.asarray(SPCparam)
+
+ dataOut.data_pre= (self.spc , self.cspc)
+
+ #dataOut.data_preParam = (self.spc , self.cspc)
+
+ dataOut.spcparam_range=numpy.zeros([self.Num_Chn,self.Num_Bin+1])
+
+ dataOut.spcparam_range[2]=VelRange
+ dataOut.spcparam_range[1]=TimeRange
+ dataOut.spcparam_range[0]=FrecRange
+
+
+
+
class GaussianFit(Operation):
'''
@@ -198,7 +310,7 @@ class GaussianFit(Operation):
self.dataOut.data_pre : SelfSpectra
Output:
- self.dataOut.GauSPC : SPC_ch1, SPC_ch2
+ self.dataOut.SPCparam : SPC_ch1, SPC_ch2
'''
def __init__(self, **kwargs):
@@ -252,7 +364,7 @@ class GaussianFit(Operation):
objs = [self for __ in range(self.Num_Chn)]
attrs = zip(objs, args)
gauSPC = pool.map(target, attrs)
- dataOut.GauSPC = numpy.asarray(gauSPC)
+ dataOut.SPCparam = numpy.asarray(SPCparam)
@@ -280,7 +392,7 @@ class GaussianFit(Operation):
#print 'HEIGHTS', self.Num_Hei
- GauSPC = []
+ SPCparam = []
SPC_ch1 = numpy.empty([self.Num_Bin,self.Num_Hei])
SPC_ch2 = numpy.empty([self.Num_Bin,self.Num_Hei])
SPC_ch1[:] = 0#numpy.NaN
@@ -331,7 +443,7 @@ class GaussianFit(Operation):
snr = numpy.NaN
SPC_ch1[:,ht] = 0#numpy.NaN
SPC_ch1[:,ht] = 0#numpy.NaN
- GauSPC = (SPC_ch1,SPC_ch2)
+ SPCparam = (SPC_ch1,SPC_ch2)
continue
#print 'snr',snrdB #, sum(spcs) , tot_noise
@@ -517,7 +629,7 @@ class GaussianFit(Operation):
#print 'SPC_ch1.shape',SPC_ch1.shape
#print 'SPC_ch2.shape',SPC_ch2.shape
#dataOut.data_param = SPC_ch1
- GauSPC = (SPC_ch1,SPC_ch2)
+ SPCparam = (SPC_ch1,SPC_ch2)
#GauSPC[1] = SPC_ch2
# print 'shift0', shift0
@@ -581,13 +693,30 @@ class PrecipitationProc(Operation):
Parameters affected:
'''
-
+ def gaus(self,xSamples,Amp,Mu,Sigma):
+ return ( Amp / ((2*numpy.pi)**0.5 * Sigma) ) * numpy.exp( -( xSamples - Mu )**2 / ( 2 * (Sigma**2) ))
+
+
+
+ def Moments(self, ySamples, xSamples):
+ Pot = numpy.nansum( ySamples ) # Potencia, momento 0
+ yNorm = ySamples / Pot
+
+ Vr = numpy.nansum( yNorm * xSamples ) # Velocidad radial, mu, corrimiento doppler, primer momento
+ Sigma2 = abs(numpy.nansum( yNorm * ( xSamples - Vr )**2 )) # Segundo Momento
+ Desv = Sigma2**0.5 # Desv. Estandar, Ancho espectral
+
+ return numpy.array([Pot, Vr, Desv])
def run(self, dataOut, radar=None, Pt=5000, Gt=295.1209, Gr=70.7945, Lambda=0.6741, aL=2.5118,
tauW=4e-06, ThetaT=0.1656317, ThetaR=0.36774087, Km = 0.93, Altitude=3350):
- Velrange = dataOut.abscissaList
+ Velrange = dataOut.spc_range[2]
+ FrecRange = dataOut.spc_range[0]
+
+ dV= Velrange[1]-Velrange[0]
+ dF= FrecRange[1]-FrecRange[0]
if radar == "MIRA35C" :
@@ -599,7 +728,7 @@ class PrecipitationProc(Operation):
else:
- self.spc = dataOut.GauSPC[1] #dataOut.data_pre[0].copy()
+ self.spc = dataOut.SPCparam[1] #dataOut.data_pre[0].copy() #
self.Num_Hei = self.spc.shape[2]
self.Num_Bin = self.spc.shape[1]
self.Num_Chn = self.spc.shape[0]
@@ -619,60 +748,93 @@ class PrecipitationProc(Operation):
Numerator = ( (4*numpy.pi)**3 * aL**2 * 16 * numpy.log(2) )
Denominator = ( Pt * Gt * Gr * Lambda**2 * SPEED_OF_LIGHT * tauW * numpy.pi * ThetaT * ThetaR)
- RadarConstant = Numerator / Denominator
+ RadarConstant = 4.1396e+08# Numerator / Denominator
print '***'
print '*** RadarConstant' , RadarConstant, '****'
print '***'
''' ============================= '''
SPCmean = numpy.mean(self.spc,0)
- ETA = numpy.zeros(self.Num_Hei,self.Num_Bin)
+ ETAf = numpy.zeros([self.Num_Bin,self.Num_Hei])
+ ETAv = numpy.zeros([self.Num_Bin,self.Num_Hei])
+ ETAd = numpy.zeros([self.Num_Bin,self.Num_Hei])
+
Pr = self.spc[0,:,:]
VelMeteoro = numpy.mean(SPCmean,axis=0)
- print '==================== Vel SHAPE',VelMeteoro
- D_range = numpy.zeros(self.Num_Hei)
- SIGMA = numpy.zeros(self.Num_Hei)
- N_dist = numpy.zeros(self.Num_Hei)
- EqSec = numpy.zeros(self.Num_Hei)
+ #print '==================== Vel SHAPE',VelMeteoro
+
+ D_range = numpy.zeros([self.Num_Bin,self.Num_Hei])
+ SIGMA = numpy.zeros([self.Num_Bin,self.Num_Hei])
+ N_dist = numpy.zeros([self.Num_Bin,self.Num_Hei])
+ D_mean = numpy.zeros(self.Num_Hei)
del_V = numpy.zeros(self.Num_Hei)
+ Z = numpy.zeros(self.Num_Hei)
+ Ze = numpy.zeros(self.Num_Hei)
+ RR = numpy.zeros(self.Num_Hei)
for R in range(self.Num_Hei):
- ETA[:,R] = RadarConstant * Pr[:,R] * R**2 #Reflectivity (ETA)
- h = R + Altitude #Range from ground to radar pulse altitude
+ h = R*75 + Altitude #Range from ground to radar pulse altitude
del_V[R] = 1 + 3.68 * 10**-5 * h + 1.71 * 10**-9 * h**2 #Density change correction for velocity
- D_range[R] = numpy.log( (9.65 - (VelMeteoro[R]/del_V[R])) / 10.3 ) / -0.6 #Range of Diameter of drops related to velocity
- SIGMA[R] = numpy.pi**5 / Lambda**4 * Km * D_range[R]**6 #Equivalent Section of drops (sigma)
- #print '******* D_range ********', self.Num_Hei
- #print D_range
- N_dist[R] = ETA[R] / SIGMA[R]
+ D_range[:,R] = numpy.log( (9.65 - (Velrange[0:self.Num_Bin] / del_V[R])) / 10.3 ) / -0.6 #Range of Diameter of drops related to velocity
+
+ ETAf[:,R] = 1/RadarConstant * Pr[:,R] * (R*0.075)**2 #Reflectivity (ETA)
+
+ ETAv[:,R]=ETAf[:,R]*dF/dV
+
+ ETAd[:,R]=ETAv[:,R]*6.18*exp(-0.6*D_range[:,R])
+
+ SIGMA[:,R] = numpy.pi**5 / Lambda**4 * Km * D_range[:,R]**6 #Equivalent Section of drops (sigma)
+
+ N_dist[:,R] = ETAd[:,R] / SIGMA[:,R]
+
+ DMoments = self.Moments(Pr[:,R], D_range[:,R])
+ try:
+ popt01,pcov = curve_fit(self.gaus, D_range[:,R] , Pr[:,R] , p0=DMoments)
+ except:
+ popt01=numpy.zeros(3)
+ popt01[1]= DMoments[1]
+ D_mean[R]=popt01[1]
+
+ Z[R] = numpy.nansum( N_dist[:,R] * D_range[:,R]**6 )
+ RR[R] = 6*10**-4.*numpy.pi * numpy.nansum( D_range[:,R]**3 * N_dist[:,R] * Velrange[0:self.Num_Bin] ) #Rainfall rate
- Ze = (ETA * Lambda**4) / (numpy.pi * Km)
- Z = numpy.sum( N_dist * D_range**6 )
- #print 'Velrange',len(Velrange), 'D_range', len(D_range)
- RR = 6*10**-4*numpy.pi * numpy.sum( D_range[R]**3 * N_dist[R] * VelMeteoro[R] ) #Rainfall rate
+ Ze[R] = (numpy.nansum(ETAd[:,R]) * Lambda**4) / (numpy.pi * Km)
- RR2 = (Ze/200)**(1/1.6)
+ RR2 = (Z/200)**(1/1.6)
dBRR = 10*numpy.log10(RR)
+ dBRR2 = 10*numpy.log10(RR2)
dBZe = 10*numpy.log10(Ze)
- dataOut.data_output = Ze
+ dBZ = 10*numpy.log10(Z)
+
+ dataOut.data_output = Z
dataOut.data_param = numpy.ones([3,self.Num_Hei])
dataOut.channelList = [0,1,2]
- print 'channelList', dataOut.channelList
- dataOut.data_param[0]=dBZe
- dataOut.data_param[1]=dBRR
+
+ dataOut.data_param[0]=dBZ
+ dataOut.data_param[1]=dBZe
+ dataOut.data_param[2]=dBRR2
+
print 'RR SHAPE', dBRR.shape
- print 'Ze SHAPE', dBZe.shape
- print 'dataOut.data_param SHAPE', dataOut.data_param.shape
+ #print 'Ze ', dBZe
+ print 'Z ', dBZ
+ print 'RR ', dBRR
+ #print 'RR2 ', dBRR2
+ #print 'D_mean', D_mean
+ #print 'del_V', del_V
+ #print 'D_range',D_range.shape, D_range[:,30]
+ #print 'Velrange', Velrange
+ #print 'numpy.nansum( N_dist[:,R]', numpy.nansum( N_dist, axis=0)
+ #print 'dataOut.data_param SHAPE', dataOut.data_param.shape
def dBZeMODE2(self, dataOut): # Processing for MIRA35C
@@ -747,7 +909,7 @@ class FullSpectralAnalysis(Operation):
self.indice=int(numpy.random.rand()*1000)
- spc = dataOut.data_pre[0]
+ spc = dataOut.data_pre[0].copy()
cspc = dataOut.data_pre[1]
nChannel = spc.shape[0]
@@ -760,12 +922,7 @@ class FullSpectralAnalysis(Operation):
else:
ChanDist = numpy.array([[E01, N01],[E02,N02],[E12,N12]])
- #print 'ChanDist', ChanDist
-
- #if dataOut.VelRange is not None:
- AbbsisaRange= dataOut.spc_range#dataOut.VelRange
- #else:
- # AbbsisaRange= dataOut.spc_range#dataOut.abscissaList
+ FrecRange = dataOut.spc_range[0]
ySamples=numpy.ones([nChannel,nProfiles])
phase=numpy.ones([nChannel,nProfiles])
@@ -773,38 +930,34 @@ class FullSpectralAnalysis(Operation):
coherence=numpy.ones([nChannel,nProfiles])
PhaseSlope=numpy.ones(nChannel)
PhaseInter=numpy.ones(nChannel)
- dataSNR = dataOut.data_SNR
-
-
+ data_SNR=numpy.zeros([nProfiles])
data = dataOut.data_pre
noise = dataOut.noise
- #print 'noise',noise
- #SNRdB = 10*numpy.log10(dataOut.data_SNR)
- FirstMoment = dataOut.data_param[0,1,:]#numpy.average(dataOut.data_param[:,1,:],0)
- SecondMoment = numpy.average(dataOut.data_param[:,2,:],0)
+ dataOut.data_SNR = (numpy.mean(spc,axis=1)- noise[0]) / noise[0]
- #SNRdBMean = []
-
- #for j in range(nHeights):
- # FirstMoment = numpy.append(FirstMoment,numpy.mean([dataOut.data_param[0,1,j],dataOut.data_param[1,1,j],dataOut.data_param[2,1,j]]))
- # SNRdBMean = numpy.append(SNRdBMean,numpy.mean([SNRdB[0,j],SNRdB[1,j],SNRdB[2,j]]))
- data_output=numpy.ones([3,spc.shape[2]])*numpy.NaN
+ print dataOut.data_SNR.shape
+ #FirstMoment = dataOut.moments[0,1,:]#numpy.average(dataOut.data_param[:,1,:],0)
+ #SecondMoment = numpy.average(dataOut.moments[:,2,:],0)
+
+ #SNRdBMean = []
+
+ data_output=numpy.ones([spc.shape[0],spc.shape[2]])*numpy.NaN
velocityX=[]
velocityY=[]
velocityV=[]
PhaseLine=[]
- dbSNR = 10*numpy.log10(dataSNR)
+ dbSNR = 10*numpy.log10(dataOut.data_SNR)
dbSNR = numpy.average(dbSNR,0)
for Height in range(nHeights):
- [Vzon,Vmer,Vver, GaussCenter, PhaseSlope, FitGaussCSPC]= self.WindEstimation(dataOut.data_param, spc, cspc, pairsList, ChanDist, Height, noise, AbbsisaRange, dbSNR[Height], SNRlimit)
+ [Vzon,Vmer,Vver, GaussCenter, PhaseSlope, FitGaussCSPC]= self.WindEstimation(spc, cspc, pairsList, ChanDist, Height, noise, dataOut.spc_range.copy(), dbSNR[Height], SNRlimit)
PhaseLine = numpy.append(PhaseLine, PhaseSlope)
if abs(Vzon)<100. and abs(Vzon)> 0.:
@@ -822,7 +975,7 @@ class FullSpectralAnalysis(Operation):
velocityY=numpy.append(velocityY, numpy.NaN)
if dbSNR[Height] > SNRlimit:
- velocityV=numpy.append(velocityV, FirstMoment[Height])
+ velocityV=numpy.append(velocityV, -Vver)#FirstMoment[Height])
else:
velocityV=numpy.append(velocityV, numpy.NaN)
#FirstMoment[Height]= numpy.NaN
@@ -836,24 +989,23 @@ class FullSpectralAnalysis(Operation):
data_output[0] = numpy.array(velocityX) #self.moving_average(numpy.array(velocityX) , N=1)
data_output[1] = numpy.array(velocityY) #self.moving_average(numpy.array(velocityY) , N=1)
data_output[2] = -velocityV#FirstMoment
-
- print ' '
- #print 'FirstMoment'
+
+ print 'FirstMoment', data_output[2]
#print FirstMoment
- print 'velocityX',numpy.shape(data_output[0])
- print 'velocityX',data_output[0]
- print ' '
- print 'velocityY',numpy.shape(data_output[1])
- print 'velocityY',data_output[1]
- print 'velocityV',data_output[2]
- print 'PhaseLine',PhaseLine
+# print 'velocityX',numpy.shape(data_output[0])
+# print 'velocityX',data_output[0]
+# print ' '
+# print 'velocityY',numpy.shape(data_output[1])
+# print 'velocityY',data_output[1]
+# print 'velocityV',data_output[2]
+# print 'PhaseLine',PhaseLine
#print numpy.array(velocityY)
#print 'SNR'
#print 10*numpy.log10(dataOut.data_SNR)
#print numpy.shape(10*numpy.log10(dataOut.data_SNR))
print ' '
- xFrec=AbbsisaRange[0][0:spc.shape[1]]
+ xFrec=FrecRange[0:spc.shape[1]]
dataOut.data_output=data_output
@@ -878,11 +1030,11 @@ class FullSpectralAnalysis(Operation):
return numpy.array([Pot, Vr, Desv])
- def WindEstimation(self, data_param, spc, cspc, pairsList, ChanDist, Height, noise, AbbsisaRange, dbSNR, SNRlimit):
+ def WindEstimation(self, spc, cspc, pairsList, ChanDist, Height, noise, AbbsisaRange, dbSNR, SNRlimit):
- print ' '
- print '######################## Height',Height, (1000 + 75*Height), '##############################'
- print ' '
+# print ' '
+# print '######################## Height',Height, (1000 + 75*Height), '##############################'
+# print ' '
ySamples=numpy.ones([spc.shape[0],spc.shape[1]])
phase=numpy.ones([spc.shape[0],spc.shape[1]])
@@ -956,10 +1108,10 @@ class FullSpectralAnalysis(Operation):
#print '##### SUMA de CSPC #####', len(coherence)
#print numpy.sum(numpy.abs(CSPCNorm))
#print numpy.sum(coherence[0])
- print 'len',len(xSamples)
- print 'CSPCmoments', numpy.shape(CSPCmoments)
- print CSPCmoments
- print '#######################'
+# print 'len',len(xSamples)
+# print 'CSPCmoments', numpy.shape(CSPCmoments)
+# print CSPCmoments
+# print '#######################'
popt=[1e-10,1e-10,1e-10]
popt01, popt02, popt12 = [1e-10,1e-10,1e-10], [1e-10,1e-10,1e-10] ,[1e-10,1e-10,1e-10]
@@ -981,8 +1133,6 @@ class FullSpectralAnalysis(Operation):
-
-
'''***Fit Gauss CSPC01***'''
if dbSNR > SNRlimit:
try:
@@ -1044,56 +1194,25 @@ class FullSpectralAnalysis(Operation):
else:
Fij = xSamples[PointGauCenter] - xSamples[PointFij]
- print 'CSPCopt'
- print CSPCopt
- print 'popt'
- print popt
- print '#######################################'
+# print 'CSPCopt'
+# print CSPCopt
+# print 'popt'
+# print popt
+# print '#######################################'
#print 'dataOut.data_param', numpy.shape(data_param)
#print 'dataOut.data_param0', data_param[0,0,Height]
#print 'dataOut.data_param1', data_param[0,1,Height]
#print 'dataOut.data_param2', data_param[0,2,Height]
- print 'yMoments', yMoments
- print 'Moments', SPCmoments
- print 'Fij2 Moment', Fij
- #print 'Fij', Fij, 'popt[2]/2',popt[2]/2
- print 'Fijcspc',Fijcspc
- print '#######################################'
-
-
-# Range = numpy.empty([3,2])
-# GaussCenter = numpy.empty(3)
-# FrecRange, VelRange = [[],[],[]],[[],[],[]]
-# FrecRange01, FrecRange02, FrecRange12 = numpy.empty(3), numpy.empty(3), numpy.empty(3)
-# VelRange01, VelRange02, VelRange12 = numpy.empty(3), numpy.empty(3), numpy.empty(3)
-#
-# GauWidth = numpy.empty(3)
-# ClosRangeMin, ClosRangeMax = numpy.empty(3), numpy.empty(3)
-# PointRangeMin, PointRangeMax = numpy.empty(3), numpy.empty(3)
-# '''****** Taking frequency ranges from CSPCs ******'''
-# for i in range(spc.shape[0]):
-#
-# GaussCenter[i] = CSPCopt[i][1] #Primer momento 01
-# GauWidth[i] = CSPCopt[i][2]*2/2 #Ancho de banda de Gau01
-#
-# Range[i][0] = GaussCenter[i] - GauWidth[i]
-# Range[i][1] = GaussCenter[i] + GauWidth[i]
-# #Punto en Eje X de la Gaussiana donde se encuentra ancho de banda (min:max)
-# ClosRangeMin[i] = xSamples[numpy.abs(xSamples-Range[i][0]).argmin()]
-# ClosRangeMax[i] = xSamples[numpy.abs(xSamples-Range[i][1]).argmin()]
-#
-# PointRangeMin[i] = numpy.where(xSamples==ClosRangeMin[i])[0][0]
-# PointRangeMax[i] = numpy.where(xSamples==ClosRangeMax[i])[0][0]
-#
-# Range[i]=numpy.array([ PointRangeMin[i], PointRangeMax[i] ])
-#
-# FrecRange[i] = xFrec[ Range[i][0] : Range[i][1] ]
-# VelRange[i] = xVel[ Range[i][0] : Range[i][1] ]
-
-
+# print 'yMoments', yMoments
+# print 'Moments', SPCmoments
+# print 'Fij2 Moment', Fij
+# #print 'Fij', Fij, 'popt[2]/2',popt[2]/2
+# print 'Fijcspc',Fijcspc
+# print '#######################################'
+
'''****** Taking frequency ranges from SPCs ******'''
@@ -1178,9 +1297,9 @@ class FullSpectralAnalysis(Operation):
(cA,cB,cH) = numpy.linalg.solve(WijEijNij, WijResults)
VxVy=numpy.array([[cA,cH],[cH,cB]])
-
VxVyResults=numpy.array([-cF,-cG])
(Vx,Vy) = numpy.linalg.solve(VxVy, VxVyResults)
+
#print 'MijResults, cC, PhaseSlope', MijResults, cC, PhaseSlope
#print 'W01,02,12', W01, W02, W12
#print 'WijResult0,1,2',WijResult0, WijResult1, WijResult2, 'Results', WijResults
@@ -1230,7 +1349,7 @@ class FullSpectralAnalysis(Operation):
- print 'vzon y vmer', Vzon, Vmer
+# print 'vzon y vmer', Vzon, Vmer
return Vzon, Vmer, Vver, GaussCenter, PhaseSlope, FitGaussCSPC
class SpectralMoments(Operation):
@@ -1257,7 +1376,7 @@ class SpectralMoments(Operation):
self.dataOut.noise : Noise level per channel
Affected:
- self.dataOut.data_param : Parameters per channel
+ self.dataOut.moments : Parameters per channel
self.dataOut.data_SNR : SNR per channel
'''
@@ -1274,7 +1393,7 @@ class SpectralMoments(Operation):
for ind in range(nChannel):
data_param[ind,:,:] = self.__calculateMoments( data[ind,:,:] , absc , noise[ind] )
- dataOut.data_param = data_param[:,1:,:]
+ dataOut.moments = data_param[:,1:,:]
dataOut.data_SNR = data_param[:,0]
return
diff --git a/schainpy/model/proc/jroproc_spectra.py b/schainpy/model/proc/jroproc_spectra.py
index 2154351..fba373f 100644
--- a/schainpy/model/proc/jroproc_spectra.py
+++ b/schainpy/model/proc/jroproc_spectra.py
@@ -162,8 +162,7 @@ class SpectraProc(ProcessingUnit):
self.dataOut.flagNoData = True
return 0
else:
- print 'DATA shape', self.dataIn.data.shape
- sadsdf
+
self.buffer[:,self.profIndex,:] = self.dataIn.data.copy()
self.profIndex += 1
diff --git a/schainpy/scripts/schain.xml b/schainpy/scripts/schain.xml
index 8979c10..0f9977c 100644
--- a/schainpy/scripts/schain.xml
+++ b/schainpy/scripts/schain.xml
@@ -1 +1 @@
-
\ No newline at end of file
+
\ No newline at end of file
diff --git a/schainpy/scripts/testRawData.py b/schainpy/scripts/testRawData.py
index c8ee5f2..e2286c0 100644
--- a/schainpy/scripts/testRawData.py
+++ b/schainpy/scripts/testRawData.py
@@ -7,8 +7,8 @@ if __name__ == '__main__':
desc = "Segundo Test"
filename = "schain.xml"
- pathW='/media/erick/6F60F7113095A154/CLAIRE/CLAIRE_WINDS_2MHZ/DATA/pdatatest/test1024'
- figpath = '/media/erick/6F60F7113095A154/CLAIRE/CLAIRE_WINDS_2MHZ/Images/test1024'
+ pathW='/media/erick/6F60F7113095A154/CLAIRE/CLAIRE_WINDS_2MHZ/DATA/pdataCLAIRE/Extra'
+ figpath = '/media/erick/6F60F7113095A154/CLAIRE/CLAIRE_WINDS_2MHZ/DATA/pdataCLAIRE/Extra'
controllerObj = Project()
@@ -17,10 +17,10 @@ if __name__ == '__main__':
readUnitConfObj = controllerObj.addReadUnit(datatype='VoltageReader',
path='/media/erick/6F60F7113095A154/CLAIRE/CLAIRE_WINDS_2MHZ/DATA/',
#path='/home/erick/Documents/Data/Claire_Data/raw',
- startDate='2017/08/22',
- endDate='2017/08/22',
- startTime='01:00:00',
- endTime='06:00:00',
+ startDate='2018/02/01',
+ endDate='2018/02/01',
+ startTime='12:00:00',
+ endTime='20:00:00',
online=0,
walk=1)
@@ -54,8 +54,8 @@ if __name__ == '__main__':
opObj10 = procUnitConfObj0.addOperation(name='setRadarFrequency')
opObj10.addParameter(name='frequency', value='445.09e6', format='float')
- #opObj10 = procUnitConfObj0.addOperation(name='CohInt', optype='external')
- #opObj10.addParameter(name='n', value='1', format='float')
+ opObj10 = procUnitConfObj0.addOperation(name='CohInt', optype='external')
+ opObj10.addParameter(name='n', value='2', format='float')
#opObj10 = procUnitConfObj0.addOperation(name='selectHeights')
#opObj10.addParameter(name='minHei', value='1', format='float')
@@ -71,13 +71,13 @@ if __name__ == '__main__':
# Creating a processing object with its parameters
# schainpy.model.proc.jroproc_spectra.SpectraProc.run()
# If you need to add more parameters can use the "addParameter method"
- procUnitConfObj1.addParameter(name='nFFTPoints', value='1024', format='int')
+ procUnitConfObj1.addParameter(name='nFFTPoints', value='256', format='int')
opObj10 = procUnitConfObj1.addOperation(name='removeDC')
#opObj10 = procUnitConfObj1.addOperation(name='removeInterference')
- #opObj10 = procUnitConfObj1.addOperation(name='IncohInt', optype='external')
- #opObj10.addParameter(name='n', value='30', format='float')
+ opObj10 = procUnitConfObj1.addOperation(name='IncohInt', optype='external')
+ opObj10.addParameter(name='n', value='10', format='float')
@@ -151,9 +151,10 @@ if __name__ == '__main__':
opObj11.addParameter(name='ymax', value='7', format='int')
opObj11.addParameter(name='showprofile', value='1', format='int')
# opObj11.addParameter(name='timerange', value=str(5*60*60*60), format='int')
- opObj11.addParameter(name='xmin', value='1', format='float')
- opObj11.addParameter(name='xmax', value='6', format='float')
+ #opObj11.addParameter(name='xmin', value='1', format='float')
+ #opObj11.addParameter(name='xmax', value='6', format='float')
opObj11.addParameter(name='save', value='1', format='int')
+ opObj11.addParameter(name='figpath', value=figpath, format='str')
# '''#########################################################################################'''