In [1]:
#first we do some imports and check the version of Py-ART for consistency
import pyart
from matplotlib import pyplot as plt
import netCDF4
import numpy as np
from copy import deepcopy
import os
import copy
%matplotlib inline
print pyart.__version__


1.0.0.dev-e4273ce

In [2]:
data_file = '/data/sgpcsaprppi_20110520095101.nc'
radar = pyart.io.read(data_file)

In [26]:
display = pyart.graph.RadarDisplay(radar)
fig = plt.figure(figsize = [15,10])
display.plot_ppi('reflectivity', sweep = 1, vmin = -8, vmax =64)
bbox_props = dict(boxstyle="rarrow,pad=0.3", fc="white", ec="b", lw=2)

display.set_limits(xlim = [-120, 120], ylim = [-120, 120])

t = plt.gca().text(-70, -70, "Severe Attenuation", ha="center", va="center", rotation=45,
            size=25,
            bbox=bbox_props)

t2 = plt.gca().text(20, -50, "Multi \n Trip", ha="center", va="center", rotation=0,
            size=25,
            bbox=bbox_props)

plt.savefig('figures/ze.png', bbox_inches='tight')



In [31]:
fig = plt.figure(figsize = [15,10])
display.plot_ppi('differential_phase', vmin = -180, vmax =180)


display.set_limits(xlim = [-120, 120], ylim = [-120, 120])

t = plt.gca().text(-80, -70, "Phase \n folding", ha="center", va="center", rotation=0,
            size=25,
            bbox=bbox_props)

t2 = plt.gca().text(20, -50, "Multi \n Trip", ha="center", va="center", rotation=0,
            size=25,
            bbox=bbox_props)


plt.savefig('figures/phidp.png', bbox_inches='tight')



In [5]:
# Perform a Linear programming phase processing to get the conditionally
# fitted phase profile and sobel filtered KDP field

c_offset = - 3.42

gates = radar.range['data'][1] - radar.range['data'][0]
rge = 30.0 * 1000.0
ng = rge / gates


reproc_phase, sob_kdp = pyart.correct.phase_proc.phase_proc_lp(radar, c_offset, debug=True, 
                                                               nowrap=ng, min_ncp = 0.6, min_rhv = 0.8)

#in some sectors where there is not enough data the LP method fails.. find these and mask


Unfolding
Exec time: 
/Users/scollis/anaconda/lib/python2.7/site-packages/pyart/correct/phase_proc.py:186: RuntimeWarning: invalid value encountered in sqrt
  noise = smooth_and_trim(np.sqrt((line - signal) ** 2), window_len=wl)
/Users/scollis/anaconda/lib/python2.7/site-packages/cylp-0.7.1-py2.7-macosx-10.6-x86_64.egg/cylp/py/modeling/CyLPModel.py:372: FutureWarning: comparison to `None` will result in an elementwise object comparison in the future.
  if left == None:
/Users/scollis/anaconda/lib/python2.7/site-packages/cylp-0.7.1-py2.7-macosx-10.6-x86_64.egg/cylp/py/utils/sparseUtil.py:388: FutureWarning: comparison to `None` will result in an elementwise object comparison in the future.
  if b == None:
/Users/scollis/anaconda/lib/python2.7/site-packages/pyart/correct/phase_proc.py:923: FutureWarning: comparison to `None` will result in an elementwise object comparison in the future.
  s = CyClpSimplex(model)
 2.69928312302
Doing  0
Doing  1

In [12]:
reproc_phase['data'][np.where(reproc_phase['data'] < 0.0)] = 0.0

radar.fields.update({'recalculated_diff_phase': sob_kdp,
                         'proc_dp_phase_shift': reproc_phase})

In [33]:
fig = plt.figure(figsize = [15,10])
display.plot_ppi('recalculated_diff_phase', vmin = -1, vmax =10)
display.set_limits(xlim = [-120, 120], ylim = [-120, 120])
t = plt.gca().text(-60, -60, "High \n KDP", ha="center", va="center", rotation=0,
            size=25,
            bbox=bbox_props)

plt.savefig('figures/kdp.png', bbox_inches='tight')



In [35]:
display = pyart.graph.RadarDisplay(radar)
fig = plt.figure(figsize = [15,10])
display.plot_ppi('proc_dp_phase_shift', vmin = 0, vmax =400)
display.set_limits(xlim = [-120, 120], ylim = [-120, 120])
t = plt.gca().text(-80, -60, "High \n Gradient", ha="center", va="center", rotation=0,
            size=25,
            bbox=bbox_props)
plt.savefig('figures/phidp_f.png', bbox_inches='tight')



In [8]:
spec_at, cor_z = pyart.correct.attenuation.calculate_attenuation(radar, c_offset, 
                                                                 debug=True, ncp_min=0.4)

#again filter where the algorithm has failed due to lack of data

spec_at['data'][np.where(spec_at['data'] < 0.0)] = 0.0

radar.fields.update({'specific_attenuation': spec_at})
radar.fields.update({'corrected_reflectivity': cor_z})


Doing  0
Doing  1

In [36]:
display = pyart.graph.RadarDisplay(radar)
fig = plt.figure(figsize = [15,10])
display.plot_ppi('specific_attenuation', vmin = 0, vmax =1)
display.set_limits(xlim = [-120, 120], ylim = [-120, 120])
t = plt.gca().text(-80, -60, "High \n Attenuation", ha="center", va="center", rotation=0,
            size=25,
            bbox=bbox_props)

plt.savefig('figures/speca.png', bbox_inches='tight')



In [ ]: