sophy_proc_rev006.py
229 lines
| 8.8 KiB
| text/x-python
|
PythonLexer
|
r1516 | #!/usr/bin/env python | ||
|
r1520 | # Ing. AVP | ||
# 22/06/2022 | ||||
# ARCHIVO DE LECTURA y PLOT | ||||
|
r1516 | import matplotlib.pyplot as pl | ||
import matplotlib | ||||
import wradlib | ||||
import numpy | ||||
import warnings | ||||
import argparse | ||||
from wradlib.io import read_generic_hdf5 | ||||
from wradlib.util import get_wradlib_data_file | ||||
from plotting_codes import sophy_cb_tables | ||||
import os,time | ||||
for name, cb_table in sophy_cb_tables: | ||||
ncmap = matplotlib.colors.ListedColormap(cb_table, name=name) | ||||
matplotlib.pyplot.register_cmap(cmap=ncmap) | ||||
#LINUX bash: export WRADLIB_DATA=/path/to/wradlib-data | ||||
|
r1525 | #example | ||
#export WRADLIB_DATA="/home/soporte/Documents/EVENTO/HYO_PM@2022-06-09T15-05-12/paramC0N36.0/2022-06-09T18-00-00/" | ||||
|
r1516 | warnings.filterwarnings('ignore') | ||
PARAM = { | ||||
|
r1519 | 'S': {'var': 'power','vmin': -45, 'vmax': -15, 'cmap': 'jet', 'label': 'Power','unit': 'dBm'}, | ||
'V': {'var': 'velocity', 'vmin': -10, 'vmax': 10 , 'cmap': 'sophy_v', 'label': 'Velocity','unit': 'm/s'}, | ||||
'Z': {'var': 'reflectivity','vmin': -30, 'vmax': 80 , 'cmap': 'sophy_r','label': 'Reflectivity','unit': 'dBZ'}, | ||||
|
r1529 | 'W': {'var': 'spectral_width', 'vmin': 0 , 'vmax': 12 , 'cmap': 'sophy_w','label': 'Spectral Width','unit': 'm/s'} | ||
|
r1516 | } | ||
class Readsophy(): | ||||
def __init__(self): | ||||
self.list_file = None | ||||
self.grado = None | ||||
self.variable = None | ||||
self.save = None | ||||
|
r1519 | self.range = None | ||
|
r1516 | |||
def read_files(self,path_file,grado=None, variable=None): | ||||
filter= "_E"+str(grado)+".0_"+variable | ||||
validFilelist = [] | ||||
fileList= os.listdir(path_file) | ||||
for thisFile in fileList: | ||||
if (os.path.splitext(thisFile)[0][-7:] != filter): | ||||
#print("s_:",os.path.splitext(thisFile)[0][-7:]) | ||||
continue | ||||
validFilelist.append(thisFile) | ||||
validFilelist.sort() | ||||
return validFilelist | ||||
|
r1519 | def setup(self, path_file,grado,range,variable,save): | ||
|
r1516 | self.path_file = path_file | ||
|
r1519 | self.range = range | ||
|
r1516 | self.grado = grado | ||
self.variable = variable | ||||
self.save = save | ||||
self.list_file = self.read_files(path_file=self.path_file,grado=self.grado, variable=self.variable) | ||||
|
r1519 | def selectHeights(self,heightList,minHei,maxHei): | ||
if minHei and maxHei: | ||||
if (minHei < heightList[0]): | ||||
minHei = heightList[0] | ||||
if (maxHei > heightList[-1]): | ||||
maxHei = heightList[-1] | ||||
minIndex = 0 | ||||
maxIndex = 0 | ||||
heights = heightList | ||||
inda = numpy.where(heights >= minHei) | ||||
indb = numpy.where(heights <= maxHei) | ||||
try: | ||||
minIndex = inda[0][0] | ||||
except: | ||||
minIndex = 0 | ||||
try: | ||||
maxIndex = indb[0][-1] | ||||
except: | ||||
maxIndex = len(heights) | ||||
new_heightList= self.selectHeightsByIndex(heightList=heightList,minIndex=minIndex, maxIndex=maxIndex) | ||||
return new_heightList, minIndex,maxIndex | ||||
def selectHeightsByIndex(self,heightList,minIndex, maxIndex): | ||||
if (minIndex < 0) or (minIndex > maxIndex): | ||||
raise ValueError("Height index range (%d,%d) is not valid" % (minIndex, maxIndex)) | ||||
if (maxIndex >= len(heightList)): | ||||
maxIndex = len(heightList) | ||||
new_h = heightList[minIndex:maxIndex] | ||||
return new_h | ||||
|
r1525 | def readAttributes(self,obj,variable): | ||
var = PARAM[variable]['var'] | ||||
unit = PARAM[variable]['unit'] | ||||
cmap = PARAM[variable]['cmap'] | ||||
vmin = PARAM[variable]['vmin'] | ||||
vmax = PARAM[variable]['vmax'] | ||||
label = PARAM[variable]['label'] | ||||
var_ = 'Data/'+var+'/H' | ||||
data_arr = numpy.array(obj[var_]['data']) # data | ||||
utc_time = numpy.array(obj['Data/time']['data']) | ||||
data_azi = numpy.array(obj['Metadata/azimuth']['data']) # th | ||||
data_ele = numpy.array(obj["Metadata/elevation"]['data']) | ||||
heightList = numpy.array(obj["Metadata/range"]['data']) # r | ||||
return data_arr, utc_time, data_azi,data_ele, heightList,unit,cmap,vmin,vmax,label | ||||
|
r1516 | def run(self): | ||
count= 0 | ||||
len_files = len(self.list_file) | ||||
for thisFile in self.list_file: | ||||
count = count +1 | ||||
fullpathfile = self.path_file + thisFile | ||||
filename = get_wradlib_data_file(fullpathfile) | ||||
test_hdf5 = read_generic_hdf5(filename) | ||||
|
r1525 | # LECTURA | ||
data_arr, utc_time, data_azi,data_ele, heightList,unit,cmap,vmin,vmax,label = self.readAttributes(obj= test_hdf5,variable=self.variable) | ||||
|
r1519 | |||
if self.range==0: | ||||
self.range == heightList[-1] | ||||
new_heightList,minIndex,maxIndex = self.selectHeights(heightList,0.06,self.range) | ||||
|
r1525 | # TIEMPO | ||
|
r1516 | my_time = time.strftime('%Y-%m-%d %H:%M:%S', time.localtime(utc_time[0])) | ||
time_save = time.strftime('%Y%m%d_%H%M%S',time.localtime(utc_time[0])) | ||||
# PLOT DATA WITH ANNOTATION | ||||
if count ==1: | ||||
fig = pl.figure(figsize=(10,8)) | ||||
|
r1519 | cgax, pm = wradlib.vis.plot_ppi(data_arr[:,minIndex:maxIndex],r=new_heightList,az=data_azi,rf=1,fig=fig, ax=111,proj='cg',cmap=cmap,vmin=vmin, vmax=vmax) | ||
|
r1516 | caax = cgax.parasites[0] | ||
|
r1521 | title = 'Simple PPI'+"-"+ my_time+" E."+self.grado | ||
|
r1516 | t = pl.title(title, fontsize=12,y=1.05) | ||
cbar = pl.gcf().colorbar(pm, pad=0.075) | ||||
pl.text(1.0, 1.05, 'azimuth', transform=caax.transAxes, va='bottom',ha='right') | ||||
cbar.set_label(label+'[' + unit + ']') | ||||
gh = cgax.get_grid_helper() | ||||
else: | ||||
|
r1519 | cgax, pm = wradlib.vis.plot_ppi(data_arr[:,minIndex:maxIndex],r=new_heightList,az=data_azi,rf=1,fig=fig, ax=111,proj='cg',cmap=cmap,vmin=vmin, vmax=vmax) | ||
|
r1516 | caax = cgax.parasites[0] | ||
|
r1521 | title = 'Simple PPI'+"-"+my_time+" E."+self.grado | ||
|
r1516 | t = pl.title(title, fontsize=12,y=1.05) | ||
cbar = pl.gcf().colorbar(pm, pad=0.075) | ||||
pl.text(1.0, 1.05, 'azimuth', transform=caax.transAxes, va='bottom',ha='right') | ||||
cbar.set_label(label+'[' + unit + ']') | ||||
gh = cgax.get_grid_helper() | ||||
if self.save == 1: | ||||
if count ==1: | ||||
filename = "SOPHY"+"_"+time_save+"_"+"E."+self.grado+"_"+self.variable+".png" | ||||
dir =self.variable+"_"+"E."+self.grado+"CH0/" | ||||
filesavepath = os.path.join(self.path_file,dir) | ||||
|
r1519 | try: | ||
os.mkdir(filesavepath) | ||||
except: | ||||
pass | ||||
|
r1516 | else: | ||
filename = "SOPHY"+"_"+time_save+"_"+"E."+self.grado+"_"+self.variable+".png" | ||||
pl.savefig(filesavepath+filename) | ||||
|
r1519 | pl.pause(1) | ||
|
r1516 | pl.clf() | ||
if count==len_files: | ||||
pl.close() | ||||
pl.show() | ||||
|
r1519 | PATH = "/home/soporte/Documents/EVENTO/HYO_PM@2022-06-09T15-05-12/paramC0N36.0/2022-06-09T18-00-00/" | ||
|
r1521 | #PATH = "/home/soporte/Documents/EVENTO/HYO_PM@2022-06-09T15-05-12/paramC0N36.0/2022-06-09T19-00-00/" | ||
|
r1519 | |||
#PATH = "/home/soporte/Documents/EVENTO/HYO_PM@2022-05-31T12-00-17/paramC0N36.0/2022-05-31T16-00-00/" | ||||
|
r1516 | |||
def main(args): | ||||
grado = args.grado | ||||
parameters = args.parameters | ||||
save = args.save | ||||
|
r1519 | range = args.range | ||
|
r1516 | obj = Readsophy() | ||
for param in parameters: | ||||
|
r1519 | obj.setup(path_file = PATH,grado = grado,range=range, variable=param,save=int(save)) | ||
|
r1516 | obj.run() | ||
if __name__ == '__main__': | ||||
parser = argparse.ArgumentParser(description='Script to process SOPHy data.') | ||||
parser.add_argument('--parameters', nargs='*', default=['S'], | ||||
help='Variables to process: P, Z, V ,W') | ||||
parser.add_argument('--grado', default=2, | ||||
help='Angle in Elev to plot') | ||||
parser.add_argument('--save', default=0, | ||||
help='Save plot') | ||||
|
r1519 | parser.add_argument('--range', default=0, type=float, | ||
help='Max range to plot') | ||||
|
r1516 | args = parser.parse_args() | ||
main(args) | ||||
|
r1519 | #python sophy_proc_rev006.py --parameters Z --grado 8 --save 1 --range 28 | ||
|
r1516 | ''' | ||
def read_and_overview(filename): | ||||
"""Read HDF5 using read_generic_hdf5 and print upper level dictionary keys | ||||
""" | ||||
test = read_generic_hdf5(filename) | ||||
print("\nPrint keys for file %s" % os.path.basename(filename)) | ||||
for key in test.keys(): | ||||
print("\t%s" % key) | ||||
return test | ||||
file__ = "/home/soporte/Documents/EVENTO/HYO_PM@2022-06-09T15-05-12/paramC0_FD_PL_R15.0km/2022-06-09T18-00-00/SOPHY_20220609_180229_E2.0_Z.hdf5" | ||||
filename = get_wradlib_data_file(file__) | ||||
print("filename:\n",filename ) | ||||
test= read_and_overview(filename) | ||||
# ANADIR INFORMACION | ||||
# informacion de los pulsos de TX | ||||
# informacion de los ruidos | ||||
# informacion de los SNR ¿? | ||||
# Aumentar la amplitud de la USRP | ||||
|
r1522 | LAST_UPDATE | ||
---- Noise | ||||
---- Mapas | ||||
|
r1516 | ''' | ||