##// END OF EJS Templates
schain
RSS
Correcion de plot eje, graficos centrados, seleccion de ploteo real,imag, abs por defecto, pruebas integrando5 min
Correcion de plot eje, graficos centrados, seleccion de ploteo real,imag, abs por defecto, pruebas integrando5 min
Juan C. Espinoza - - Load All Authors

File last commit:

r1085:0d7a50b6a171
r1244:464f48b8a1f9
Show More
bltrproc_parameters.py
402 lines | 15.6 KiB | text/x-python | PythonLexer
/ schainpy / model / proc / bltrproc_parameters.py
History | Source | Raw |Copy content |Copy permalink
ebocanegra
first commit
 r965 '''
Created on Oct 24, 2016
@author: roj- LouVD
'''
import numpy
import copy
import datetime
import time
from time import gmtime
from numpy import transpose
Juan C. Espinoza
BLTR ok
 r1018 from jroproc_base import ProcessingUnit, Operation
from schainpy.model.data.jrodata import Parameters
ebocanegra
first commit
 r965
Juan C. Espinoza
BLTRParamreader ready
 r1010 class BLTRParametersProc(ProcessingUnit):
ebocanegra
first commit
 r965 '''
Juan C. Espinoza
BLTRParamreader ready
 r1010 Processing unit for BLTR parameters data (winds)
ebocanegra
first commit
 r965 Inputs:
self.dataOut.nmodes - Number of operation modes
self.dataOut.nchannels - Number of channels
self.dataOut.nranges - Number of ranges
Juan C. Espinoza
BLTRParamreader ready
 r1010
ebocanegra
first commit
 r965 self.dataOut.data_SNR - SNR array
self.dataOut.data_output - Zonal, Vertical and Meridional velocity array
self.dataOut.height - Height array (km)
self.dataOut.time - Time array (seconds)
Juan C. Espinoza
BLTRParamreader ready
 r1010
ebocanegra
first commit
 r965 self.dataOut.fileIndex -Index of the file currently read
self.dataOut.lat - Latitude coordinate of BLTR location
Juan C. Espinoza
BLTRParamreader ready
 r1010
ebocanegra
first commit
 r965 self.dataOut.doy - Experiment doy (number of the day in the current year)
self.dataOut.month - Experiment month
self.dataOut.day - Experiment day
self.dataOut.year - Experiment year
'''
Juan C. Espinoza
BLTRParamreader ready
 r1010
def __init__(self, **kwargs):
ebocanegra
first commit
 r965 '''
Juan C. Espinoza
BLTRParamreader ready
 r1010 Inputs: None
ebocanegra
first commit
 r965 '''
ebocanegra
22/08/2017
 r1006 ProcessingUnit.__init__(self, **kwargs)
ebocanegra
first commit
 r965 self.dataOut = Parameters()
Juan C. Espinoza
Operation MAD2Writer done Task #343
 r1021 self.isConfig = False
ebocanegra
first commit
 r965
Juan C. Espinoza
Operation MAD2Writer done Task #343
 r1021 def setup(self, mode):
'''
Juan C. Espinoza
BLTRParamreader ready
 r1010 '''
Juan C. Espinoza
Operation MAD2Writer done Task #343
 r1021 self.dataOut.mode = mode
Juan C. Espinoza
BLTR ok
 r1018
Juan C. Espinoza
Operation MAD2Writer done Task #343
 r1021 def run(self, mode, snr_threshold=None):
'''
Juan C. Espinoza
BLTR ok
 r1018 Inputs:
mode = High resolution (0) or Low resolution (1) data
snr_threshold = snr filter value
Juan C. Espinoza
BLTRParamreader ready
 r1010 '''
Juan C. Espinoza
Operation MAD2Writer done Task #343
 r1021
if not self.isConfig:
self.setup(mode)
self.isConfig = True
Juan C. Espinoza
BLTR ok
 r1018 if self.dataIn.type == 'Parameters':
ebocanegra
first commit
 r965 self.dataOut.copy(self.dataIn)
Juan C. Espinoza
Fix utc times and bugs in BLTR modules
 r1085
self.dataOut.data_param = self.dataOut.data[mode]
Juan C. Espinoza
BLTR ok
 r1018 self.dataOut.heightList = self.dataOut.height[0]
self.dataOut.data_SNR = self.dataOut.data_SNR[mode]
ebocanegra
first commit
 r965
Juan C. Espinoza
BLTR ok
 r1018 if snr_threshold is not None:
SNRavg = numpy.average(self.dataOut.data_SNR, axis=0)
SNRavgdB = 10*numpy.log10(SNRavg)
for i in range(3):
Juan C. Espinoza
Fix utc times and bugs in BLTR modules
 r1085 self.dataOut.data_param[i][SNRavgdB <= snr_threshold] = numpy.nan
Juan C. Espinoza
BLTR ok
 r1018
# TODO
class OutliersFilter(Operation):
def __init__(self, **kwargs):
ebocanegra
first commit
 r965 '''
'''
Juan C. Espinoza
BLTR ok
 r1018 Operation.__init__(self, **kwargs)
ebocanegra
first commit
 r965
Juan C. Espinoza
BLTR ok
 r1018 def run(self, svalue2, method, factor, filter, npoints=9):
ebocanegra
first commit
 r965 '''
Inputs:
svalue - string to select array velocity
svalue2 - string to choose axis filtering
method - 0 for SMOOTH or 1 for MEDIAN
Juan C. Espinoza
BLTR ok
 r1018 factor - number used to set threshold
ebocanegra
first commit
 r965 filter - 1 for data filtering using the standard deviation criteria else 0
npoints - number of points for mask filter
Juan C. Espinoza
BLTR ok
 r1018 '''
print ' Outliers Filter {} {} / threshold = {}'.format(svalue, svalue, factor)
ebocanegra
first commit
 r965
Juan C. Espinoza
BLTR ok
 r1018 yaxis = self.dataOut.heightList
xaxis = numpy.array([[self.dataOut.utctime]])
ebocanegra
first commit
 r965
Juan C. Espinoza
BLTR ok
 r1018 # Zonal
value_temp = self.dataOut.data_output[0]
ebocanegra
first commit
 r965
Juan C. Espinoza
BLTR ok
 r1018 # Zonal
value_temp = self.dataOut.data_output[1]
ebocanegra
first commit
 r965
Juan C. Espinoza
BLTR ok
 r1018 # Vertical
value_temp = numpy.transpose(self.dataOut.data_output[2])
htemp = yaxis
ebocanegra
first commit
 r965 std = value_temp
for h in range(len(htemp)):
nvalues_valid = len(numpy.where(numpy.isfinite(value_temp[h]))[0])
Juan C. Espinoza
BLTR ok
 r1018 minvalid = npoints
ebocanegra
first commit
 r965
#only if valid values greater than the minimum required (10%)
if nvalues_valid > minvalid:
if method == 0:
#SMOOTH
w = value_temp[h] - self.Smooth(input=value_temp[h], width=npoints, edge_truncate=1)
if method == 1:
#MEDIAN
w = value_temp[h] - self.Median(input=value_temp[h], width = npoints)
dw = numpy.std(w[numpy.where(numpy.isfinite(w))],ddof = 1)
threshold = dw*factor
value_temp[numpy.where(w > threshold),h] = numpy.nan
value_temp[numpy.where(w < -1*threshold),h] = numpy.nan
#At the end
if svalue2 == 'inHeight':
value_temp = numpy.transpose(value_temp)
output_array[:,m] = value_temp
if svalue == 'zonal':
self.dataOut.data_output[0] = output_array
elif svalue == 'meridional':
self.dataOut.data_output[1] = output_array
elif svalue == 'vertical':
self.dataOut.data_output[2] = output_array
return self.dataOut.data_output
def Median(self,input,width):
'''
Inputs:
input - Velocity array
width - Number of points for mask filter
'''
if numpy.mod(width,2) == 1:
pc = int((width - 1) / 2)
cont = 0
output = []
for i in range(len(input)):
if i >= pc and i < len(input) - pc:
new2 = input[i-pc:i+pc+1]
temp = numpy.where(numpy.isfinite(new2))
new = new2[temp]
value = numpy.median(new)
output.append(value)
output = numpy.array(output)
output = numpy.hstack((input[0:pc],output))
output = numpy.hstack((output,input[-pc:len(input)]))
return output
def Smooth(self,input,width,edge_truncate = None):
'''
Inputs:
input - Velocity array
width - Number of points for mask filter
edge_truncate - 1 for truncate the convolution product else
'''
if numpy.mod(width,2) == 0:
real_width = width + 1
nzeros = width / 2
else:
real_width = width
nzeros = (width - 1) / 2
half_width = int(real_width)/2
length = len(input)
gate = numpy.ones(real_width,dtype='float')
norm_of_gate = numpy.sum(gate)
nan_process = 0
nan_id = numpy.where(numpy.isnan(input))
if len(nan_id[0]) > 0:
nan_process = 1
pb = numpy.zeros(len(input))
pb[nan_id] = 1.
input[nan_id] = 0.
if edge_truncate == True:
output = numpy.convolve(input/norm_of_gate,gate,mode='same')
elif edge_truncate == False or edge_truncate == None:
output = numpy.convolve(input/norm_of_gate,gate,mode='valid')
output = numpy.hstack((input[0:half_width],output))
output = numpy.hstack((output,input[len(input)-half_width:len(input)]))
if nan_process:
pb = numpy.convolve(pb/norm_of_gate,gate,mode='valid')
pb = numpy.hstack((numpy.zeros(half_width),pb))
pb = numpy.hstack((pb,numpy.zeros(half_width)))
output[numpy.where(pb > 0.9999)] = numpy.nan
input[nan_id] = numpy.nan
return output
def Average(self,aver=0,nhaver=1):
'''
Inputs:
aver - Indicates the time period over which is averaged or consensus data
nhaver - Indicates the decimation factor in heights
'''
nhpoints = 48
lat_piura = -5.17
lat_huancayo = -12.04
lat_porcuya = -5.8
if '%2.2f'%self.dataOut.lat == '%2.2f'%lat_piura:
hcm = 3.
if self.dataOut.year == 2003 :
if self.dataOut.doy >= 25 and self.dataOut.doy < 64:
nhpoints = 12
elif '%2.2f'%self.dataOut.lat == '%2.2f'%lat_huancayo:
hcm = 3.
if self.dataOut.year == 2003 :
if self.dataOut.doy >= 25 and self.dataOut.doy < 64:
nhpoints = 12
elif '%2.2f'%self.dataOut.lat == '%2.2f'%lat_porcuya:
hcm = 5.#2
pdata = 0.2
taver = [1,2,3,4,6,8,12,24]
t0 = 0
tf = 24
ntime =(tf-t0)/taver[aver]
ti = numpy.arange(ntime)
tf = numpy.arange(ntime) + taver[aver]
old_height = self.dataOut.heightList
if nhaver > 1:
num_hei = len(self.dataOut.heightList)/nhaver/self.dataOut.nmodes
deltha = 0.05*nhaver
minhvalid = pdata*nhaver
for im in range(self.dataOut.nmodes):
new_height = numpy.arange(num_hei)*deltha + self.dataOut.height[im,0] + deltha/2.
data_fHeigths_List = []
data_fZonal_List = []
data_fMeridional_List = []
data_fVertical_List = []
startDTList = []
for i in range(ntime):
height = old_height
start = datetime.datetime(self.dataOut.year,self.dataOut.month,self.dataOut.day) + datetime.timedelta(hours = int(ti[i])) - datetime.timedelta(hours = 5)
stop = datetime.datetime(self.dataOut.year,self.dataOut.month,self.dataOut.day) + datetime.timedelta(hours = int(tf[i])) - datetime.timedelta(hours = 5)
limit_sec1 = time.mktime(start.timetuple())
limit_sec2 = time.mktime(stop.timetuple())
t1 = numpy.where(self.f_timesec >= limit_sec1)
t2 = numpy.where(self.f_timesec < limit_sec2)
time_select = []
for val_sec in t1[0]:
if val_sec in t2[0]:
time_select.append(val_sec)
time_select = numpy.array(time_select,dtype = 'int')
minvalid = numpy.ceil(pdata*nhpoints)
zon_aver = numpy.zeros([self.dataOut.nranges,self.dataOut.nmodes],dtype='f4') + numpy.nan
mer_aver = numpy.zeros([self.dataOut.nranges,self.dataOut.nmodes],dtype='f4') + numpy.nan
ver_aver = numpy.zeros([self.dataOut.nranges,self.dataOut.nmodes],dtype='f4') + numpy.nan
if nhaver > 1:
new_zon_aver = numpy.zeros([num_hei,self.dataOut.nmodes],dtype='f4') + numpy.nan
new_mer_aver = numpy.zeros([num_hei,self.dataOut.nmodes],dtype='f4') + numpy.nan
new_ver_aver = numpy.zeros([num_hei,self.dataOut.nmodes],dtype='f4') + numpy.nan
if len(time_select) > minvalid:
time_average = self.f_timesec[time_select]
for im in range(self.dataOut.nmodes):
for ih in range(self.dataOut.nranges):
if numpy.sum(numpy.isfinite(self.f_zon[time_select,ih,im])) >= minvalid:
zon_aver[ih,im] = numpy.nansum(self.f_zon[time_select,ih,im]) / numpy.sum(numpy.isfinite(self.f_zon[time_select,ih,im]))
if numpy.sum(numpy.isfinite(self.f_mer[time_select,ih,im])) >= minvalid:
mer_aver[ih,im] = numpy.nansum(self.f_mer[time_select,ih,im]) / numpy.sum(numpy.isfinite(self.f_mer[time_select,ih,im]))
if numpy.sum(numpy.isfinite(self.f_ver[time_select,ih,im])) >= minvalid:
ver_aver[ih,im] = numpy.nansum(self.f_ver[time_select,ih,im]) / numpy.sum(numpy.isfinite(self.f_ver[time_select,ih,im]))
if nhaver > 1:
for ih in range(num_hei):
hvalid = numpy.arange(nhaver) + nhaver*ih
if numpy.sum(numpy.isfinite(zon_aver[hvalid,im])) >= minvalid:
new_zon_aver[ih,im] = numpy.nansum(zon_aver[hvalid,im]) / numpy.sum(numpy.isfinite(zon_aver[hvalid,im]))
if numpy.sum(numpy.isfinite(mer_aver[hvalid,im])) >= minvalid:
new_mer_aver[ih,im] = numpy.nansum(mer_aver[hvalid,im]) / numpy.sum(numpy.isfinite(mer_aver[hvalid,im]))
if numpy.sum(numpy.isfinite(ver_aver[hvalid,im])) >= minvalid:
new_ver_aver[ih,im] = numpy.nansum(ver_aver[hvalid,im]) / numpy.sum(numpy.isfinite(ver_aver[hvalid,im]))
if nhaver > 1:
zon_aver = new_zon_aver
mer_aver = new_mer_aver
ver_aver = new_ver_aver
height = new_height
tstart = time_average[0]
tend = time_average[-1]
startTime = time.gmtime(tstart)
year = startTime.tm_year
month = startTime.tm_mon
day = startTime.tm_mday
hour = startTime.tm_hour
minute = startTime.tm_min
second = startTime.tm_sec
startDTList.append(datetime.datetime(year,month,day,hour,minute,second))
o_height = numpy.array([])
o_zon_aver = numpy.array([])
o_mer_aver = numpy.array([])
o_ver_aver = numpy.array([])
if self.dataOut.nmodes > 1:
for im in range(self.dataOut.nmodes):
if im == 0:
h_select = numpy.where(numpy.bitwise_and(height[0,:] >=0,height[0,:] <= hcm,numpy.isfinite(height[0,:])))
else:
h_select = numpy.where(numpy.bitwise_and(height[1,:] > hcm,height[1,:] < 20,numpy.isfinite(height[1,:])))
ht = h_select[0]
o_height = numpy.hstack((o_height,height[im,ht]))
o_zon_aver = numpy.hstack((o_zon_aver,zon_aver[ht,im]))
o_mer_aver = numpy.hstack((o_mer_aver,mer_aver[ht,im]))
o_ver_aver = numpy.hstack((o_ver_aver,ver_aver[ht,im]))
data_fHeigths_List.append(o_height)
data_fZonal_List.append(o_zon_aver)
data_fMeridional_List.append(o_mer_aver)
data_fVertical_List.append(o_ver_aver)
else:
h_select = numpy.where(numpy.bitwise_and(height[0,:] <= hcm,numpy.isfinite(height[0,:])))
ht = h_select[0]
o_height = numpy.hstack((o_height,height[im,ht]))
o_zon_aver = numpy.hstack((o_zon_aver,zon_aver[ht,im]))
o_mer_aver = numpy.hstack((o_mer_aver,mer_aver[ht,im]))
o_ver_aver = numpy.hstack((o_ver_aver,ver_aver[ht,im]))
data_fHeigths_List.append(o_height)
data_fZonal_List.append(o_zon_aver)
data_fMeridional_List.append(o_mer_aver)
data_fVertical_List.append(o_ver_aver)
return startDTList, data_fHeigths_List, data_fZonal_List, data_fMeridional_List, data_fVertical_List