Python/PyArt tests - statistical methods & finding modified bins


I have been getting familiar with Python and PyArt. This is a script that adds a data field to the radar object that contains the number of Nyquist intervals added or subtracted for each bin after application of the region-based dealiasing algorithm. Output is a PPI plot and the example is based in the following case:


In [1]:
import matplotlib.pyplot as plt
import pylab as plb
import matplotlib as mpl
import pyart
import numpy as np
import numpy.ma as ma
from pylab import *

In [3]:
pyart.__version__


Out[3]:
'1.6.0.dev+84bb00c'

In [4]:
## COLORMAP #####################################################################
d_cmap = cm.get_cmap('PiYG', 11)

In [5]:
## SETTINGS #####################################################################

in_path = './data/'
out_path = './output/'
filename = 'CDV130618145623.RAWCBRF'
radar_abbr = filename[:3]
sw_sel = 2

In [10]:
## DATA ##########################################################################

in_file = in_path + filename
#radar = pyart.io.read_rsl(in_file)
radar = pyart.io.read(in_file)
radar.metadata['instrument_name'] = radar_abbr

Ny_vel = radar.instrument_parameters['nyquist_velocity']['data'][0]
sw_num = radar.nsweeps
sw_elevs = [radar.fixed_angle['data'][sw] for sw in range(0, sw_num-1)]

el_sel = sw_elevs[sw_sel]

corrV_reg = pyart.correct.dealias_region_based(radar, interval_splits=20, rays_wrap_around=True, keep_original=False)
radar.add_field('corrected_velocity_reg', corrV_reg)

In [11]:
radar.instrument_parameters


Out[11]:
{'nyquist_velocity': {'comments': 'Unambiguous velocity',
  'data': array([ 13.32499981,  13.32499981,  13.32499981, ...,  13.32499981,
          13.32499981,  13.32499981], dtype=float32),
  'long_name': 'Nyquist velocity',
  'meta_group': 'instrument_parameters',
  'units': 'meters_per_second'},
 'prt': {'comments': 'Pulse repetition time. For staggered prt, also see prt_ratio.',
  'data': array([ 0.001,  0.001,  0.001, ...,  0.001,  0.001,  0.001], dtype=float32),
  'long_name': 'Pulse repetition time',
  'meta_group': 'instrument_parameters',
  'units': 'seconds'},
 'prt_mode': {'comments': 'Pulsing mode Options are: "fixed", "staggered", "dual". Assumed "fixed" if missing.',
  'data': array(['fixed', 'fixed', 'fixed', 'fixed', 'fixed', 'fixed', 'fixed',
         'fixed', 'fixed'], 
        dtype='|S5'),
  'long_name': 'Pulsing mode',
  'meta_group': 'instrument_parameters',
  'units': 'unitless'},
 'pulse_width': {'comments': 'Pulse width',
  'data': array([  4.99999987e-06,   4.99999987e-06,   4.99999987e-06, ...,
           4.99999987e-06,   4.99999987e-06,   4.99999987e-06], dtype=float32),
  'long_name': 'Pulse width',
  'meta_group': 'instrument_parameters',
  'units': 'seconds'},
 'radar_beam_width_h': {'data': array([ 1.20000005], dtype=float32),
  'long_name': 'Antenna beam width H polarization',
  'meta_group': 'radar_parameters',
  'units': 'degrees'},
 'radar_beam_width_v': {'data': array([ 1.10000002], dtype=float32),
  'long_name': 'Antenna beam width V polarization',
  'meta_group': 'radar_parameters',
  'units': 'degrees'},
 'unambiguous_range': {'comments': 'Unambiguous range',
  'data': array([ 149896.5,  149896.5,  149896.5, ...,  149896.5,  149896.5,
          149896.5], dtype=float32),
  'long_name': 'Unambiguous range',
  'meta_group': 'instrument_parameters',
  'units': 'meters'}}

In [ ]:
## STATISTICS of MASKED ARRAYS ###################################################

dZ = radar.fields['reflectivity']['data']
dZ_mean = ma.mean(dZ)
dZ_std = ma.std(dZ)
dZ_max = ma.max(dZ)
indxs = ma.where(dZ == ma.max(dZ))

print('Statistics of reflectivity field: %.0f' % (dZ_mean) + " +/- " + '%.0f' % (dZ_std))

In [ ]:
## CHANGES AFTER DEALIASING #####################################################

diff_reg = radar.fields['corrected_velocity_reg']['data'] - radar.fields['velocity']['data']
diff_reg_sc = diff_reg/(2*Ny_vel)
diff_f = radar.fields['corrected_velocity_reg']
diff_f['data'] = diff_reg_sc
diff_f['long_name'] = 'Added Nyquist intervals'
diff_f['standard_name'] = "added_Ny_intervals"
diff_f['units']=''
radar.add_field('dealiasing_differences_sc', diff_f)

In [ ]:
## PLOTTING #####################################################################

out_file = out_path + filename.split('.', 1)[0]+ '_el%.0f' % (el_sel) + '_changes.png'

display = pyart.graph.RadarDisplay(radar)
fig = plt.figure(figsize=(8,6.5))
ax = fig.add_subplot(111)
#display.plot('dealiasing_differences_sc', sw_sel, vmin=-5, vmax=5, ax=ax, mask_outside=False, cmap=d_cmap)
display.plot('velocity', sw_sel, vmin=-Ny_vel, vmax=Ny_vel, ax=ax, mask_outside=False)
display.plot_range_rings(range(25, 125, 25))
display.plot_cross_hair(0.5)
plt.xlim((-75, 75))
plt.ylim((-75, 75))
#plt.savefig(out_file)
#plt.close()
plt.show()

In [ ]:
indxs = ma.nonzero(diff_reg_sc)
Nbin = shape(indxs)[1]