Radars that transmit and receive signals in polarizations aligned both horizontal and vertical to the horizon collect a number of measurements. The relation both between these measurements and between measurements and desired microphysical quantities (such as rainfall rate) is complicated due to a number of scattering mechanisms. The result is that there ends up being an intractable number of often incompatible techniques for extracting geophysical insight. This presentation will discuss methods developed by the Atmospheric Measurement Climate (ARM) Research Facility to streamline the creation of application chains for retrieving rainfall properties for the purposes of fine scale model evaluation. By using a Common Data Model (CDM) approach and working in the popular open source Python scientific environment analysis techniques such as Linear Programming (LP) can be bought to bear on the task of retrieving insight from radar signals. This presentation will outline how we have used these techniques to detangle polarimetric phase signals, estimate a three-dimensional precipitation field and then objectively compare to cloud resolving model derived rainfall fields from the NASA/DoE Mid-Latitude Continental Convective Clouds Experiment (MC3E). All techniques show will be available, open source, in the Python-ARM Radar Toolkit (Py-ART).


In [2]:
#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
import copy
import scipy
from matplotlib import rc
from skimage.transform import warp, AffineTransform
import os
from __future__ import print_function


%matplotlib inline
print(pyart._debug_info())


Py-ART version: 1.0.0.dev-f9f55d4

---- Dependencies ----
Python version: 2.7.8
NumPy version: 1.9.1
SciPy version: 0.14.0
matplotlib version: 1.4.0
netCDF4 version: 1.1.1

---- Optional dependencies ----
TRMM RSL version: v1.46
CyLP: MISSING
PyGLPK version: MISSING
CVXOPT version: MISSING
basemap version: 1.0.7
nose version: 1.3.4
None

In [3]:
data_dir = '/data/qpe/20110520/'
files = os.listdir(data_dir)
files.sort()

In [4]:
print(files[80:90])


['sgpcsaprqpeI7.c0.20110520.091602.nc', 'sgpcsaprqpeI7.c0.20110520.092302.nc', 'sgpcsaprqpeI7.c0.20110520.093001.nc', 'sgpcsaprqpeI7.c0.20110520.093701.nc', 'sgpcsaprqpeI7.c0.20110520.094401.nc', 'sgpcsaprqpeI7.c0.20110520.095101.nc', 'sgpcsaprqpeI7.c0.20110520.095801.nc', 'sgpcsaprqpeI7.c0.20110520.100501.nc', 'sgpcsaprqpeI7.c0.20110520.101200.nc', 'sgpcsaprqpeI7.c0.20110520.101901.nc']

In [5]:
list_of_grids = [pyart.io.read_grid(data_dir + files[i]) for i in range(80,90)]

In [6]:
print(list_of_grids[0].fields.keys())


[u'rain_rate_A']

In [43]:
nw_loc = [36.7675300687551 , -97.5476898998022]
sw_loc = [36.4910300262272 , -97.5943199358881]
se_loc = [36.5787000395358 , -97.3637798801064]

cf_loc = [bad_way.axes['lat']['data'], bad_way.axes['lon']['data']]
csapr_loc = [36.7961578369141, -97.4505462646484 ]

In [49]:
max_lat = 39
min_lat =35
min_lon = -100
max_lon = -94.5
plt.figure(figsize = [15,15])
display = pyart.graph.GridMapDisplay(list_of_grids[5])
display.plot_basemap(lat_lines=np.arange(min_lat,max_lat,.1),
                         lon_lines=np.arange(min_lon, max_lon, .3),
                         resolution='h')
display.plot_grid('rain_rate_A', vmin=0, vmax=150)
#xcf,ycf = display.basemap(-87.9706,41.6815)
#display.basemap.plot(xcf,ycf,'ro')
#plt.text(xcf+2000.,ycf+2000., 'Argonne Lab')
display.basemap.drawcounties()
display.basemap.drawrivers()


xc,yc = display.basemap(csapr_loc[1], csapr_loc[0])
display.basemap.plot(xc,yc,'go')
plt.text(xc-1500.,yc-6000., 'C', fontdict={'color':'white', 'weight':'bold'})


xsw,ysw = display.basemap(sw_loc[1], sw_loc[0])
display.basemap.plot(xsw,ysw,'go')
plt.text(xsw+2000.,ysw+2000., 'X', fontdict={'color':'white', 'weight':'bold'})

xse,yse = display.basemap(se_loc[1], se_loc[0])
display.basemap.plot(xse,yse,'go')
plt.text(xse+2000.,yse+2000., 'X', fontdict={'color':'white', 'weight':'bold'})

xnw,ynw = display.basemap(nw_loc[1], nw_loc[0])
display.basemap.plot(xnw,ynw,'go')
plt.text(xnw+2000.,ynw+2000., 'X', fontdict={'color':'white', 'weight':'bold'})

display.basemap.drawmapscale(-97.0, 35.9, -98, 
                             36, 100.0, barstyle='fancy')

display.plot_colorbar()
plt.savefig('figures/qpe.png', bbox_inches='tight')



In [19]:
Delta_t = (netCDF4.num2date( list_of_grids[-1].axes['time']['data'], list_of_grids[-1].axes['time']['units']) -\
           netCDF4.num2date( list_of_grids[0].axes['time']['data'],  list_of_grids[0].axes['time']['units']))[0].seconds/(60.*60.)

n = len(list_of_grids)
wts = np.ones(n)*2.0
wts[0] = 1.
wts[-1] = 1.
wts = wts * (Delta_t/(2.0*(len(list_of_grids)-1.)))

In [50]:
bad_way = pyart.retrieve.add_grids(list_of_grids, weights = wts)
bad_way_trimmed = pyart.retrieve.grid_shift(bad_way, [0.0, 0.0], trim_edges = 30)
bad_way_trimmed.fields['rain_rate_A']['units'] = 'mm'
bad_way_trimmed.fields['rain_rate_A']['standard_name'] = 'Rainfall'
bad_way_trimmed.fields['rain_rate_A']['long_name'] = 'Rainfall'

In [51]:
max_lat = 39
min_lat =35
min_lon = -100
max_lon = -94.5
plt.figure(figsize = [15,15])
display = pyart.graph.GridMapDisplay(bad_way_trimmed)
display.plot_basemap(lat_lines=np.arange(min_lat,max_lat,.1),
                         lon_lines=np.arange(min_lon, max_lon, .3),
                         resolution='h')
display.plot_grid('rain_rate_A', vmin=0, vmax=100)


xcf,ycf = display.basemap(cf_loc[1], cf_loc[0])
display.basemap.plot(xcf,ycf,'ro')
plt.text(xcf+2000.,ycf+2000., 'CF', fontdict={'color':'white', 'weight':'bold'})


xc,yc = display.basemap(csapr_loc[1], csapr_loc[0])
display.basemap.plot(xc,yc,'go')
plt.text(xc-1500.,yc-6000., 'C', fontdict={'color':'white', 'weight':'bold'})


xsw,ysw = display.basemap(sw_loc[1], sw_loc[0])
display.basemap.plot(xsw,ysw,'go')
plt.text(xsw+2000.,ysw+2000., 'X', fontdict={'color':'white', 'weight':'bold'})

xse,yse = display.basemap(se_loc[1], se_loc[0])
display.basemap.plot(xse,yse,'go')
plt.text(xse+2000.,yse+2000., 'X', fontdict={'color':'white', 'weight':'bold'})

xnw,ynw = display.basemap(nw_loc[1], nw_loc[0])
display.basemap.plot(xnw,ynw,'go')
plt.text(xnw+2000.,ynw+2000., 'X', fontdict={'color':'white', 'weight':'bold'})

display.basemap.drawcounties()
display.basemap.drawrivers()
display.basemap.drawmapscale(-97.0, 37.2, -98, 
                             37, 100.0, barstyle='fancy')

display.plot_colorbar()
plt.savefig('figures/basic_accumulation.png', bbox_inches='tight')



In [31]:
agg_list = []
for i in range(len(list_of_grids)-1):
    grid1 = list_of_grids[i]
    grid2 = list_of_grids[i+1]
    Delta_t = (netCDF4.num2date( grid2.axes['time']['data'], 
                                grid2.axes['time']['units']) -\
               netCDF4.num2date( grid1.axes['time']['data'],  
                                grid1.axes['time']['units']))[0].seconds/(60.*60.)
    
    displacement = pyart.retrieve.grid_displacement_pc(grid1, grid2, 'rain_rate_A', 0)
    print(displacement)
    n = 5
    grid_list = pyart.retrieve.create_substep_grids(grid1, grid2, displacement, n, 30 )
    swts = np.ones(n)*2.0
    swts[0] = 1.
    swts[-1] = 1.
    swts = swts * (Delta_t/(2.0*(len(grid_list)-1.)))
    good_way = pyart.retrieve.add_grids(grid_list, weights=swts)
    agg_list.append(good_way)


(1, 8)
(2, 8)
(2, 9)
(2, 7)
(2, 8)
(3, 6)
(3, 8)
(3, 8)
(2, 7)

In [32]:
n = len(agg_list)
awts = np.ones(n)*2.0
awts[0] = 1.
awts[-1] = 1.
awts = awts * (Delta_t/(2.0*(len(agg_list)-1.)))

displaced_agg =  pyart.retrieve.add_grids(agg_list, weights=1.0*np.ones(len(agg_list)))

In [52]:
displaced_agg.fields['rain_rate_A']['units'] = 'mm'
displaced_agg.fields['rain_rate_A']['standard_name'] = 'Rainfall'
displaced_agg.fields['rain_rate_A']['long_name'] = 'Rainfall'

In [53]:
max_lat = 39
min_lat =35
min_lon = -100
max_lon = -94.5
plt.figure(figsize = [15,15])
display = pyart.graph.GridMapDisplay(displaced_agg)
display.plot_basemap(lat_lines=np.arange(min_lat,max_lat,.1),
                         lon_lines=np.arange(min_lon, max_lon, .3),
                         resolution='h')
display.plot_grid('rain_rate_A', vmin=0, vmax=100)
#xcf,ycf = display.basemap(-87.9706,41.6815)
#display.basemap.plot(xcf,ycf,'ro')
#plt.text(xcf+2000.,ycf+2000., 'Argonne Lab')
display.basemap.drawcounties()
display.basemap.drawrivers()

xcf,ycf = display.basemap(cf_loc[1], cf_loc[0])
display.basemap.plot(xcf,ycf,'ro')
plt.text(xcf+2000.,ycf+2000., 'CF', fontdict={'color':'white', 'weight':'bold'})


xc,yc = display.basemap(csapr_loc[1], csapr_loc[0])
display.basemap.plot(xc,yc,'go')
plt.text(xc-1500.,yc-6000., 'C', fontdict={'color':'white', 'weight':'bold'})


xsw,ysw = display.basemap(sw_loc[1], sw_loc[0])
display.basemap.plot(xsw,ysw,'go')
plt.text(xsw+2000.,ysw+2000., 'X', fontdict={'color':'white', 'weight':'bold'})

xse,yse = display.basemap(se_loc[1], se_loc[0])
display.basemap.plot(xse,yse,'go')
plt.text(xse+2000.,yse+2000., 'X', fontdict={'color':'white', 'weight':'bold'})

xnw,ynw = display.basemap(nw_loc[1], nw_loc[0])
display.basemap.plot(xnw,ynw,'go')
plt.text(xnw+2000.,ynw+2000., 'X', fontdict={'color':'white', 'weight':'bold'})

display.basemap.drawmapscale(-97.0, 37.2, -98, 
                             37, 100.0, barstyle='fancy')

display.plot_colorbar()
plt.savefig('figures/advective_accum.png', bbox_inches='tight')



In [ ]: