Hydrometeorology Research Group

In [5]:
from IPython.display import HTML

    function code_toggle() {
     if (code_show){
     } else {
     code_show = !code_show
    $( document ).ready(code_toggle);
    <form action="javascript:code_toggle()"><input type="submit" value="Click here to toggle on/off the raw code."></form>''')
In [1]:
%matplotlib inline
    from rain import *
In [107]:
    path = '../input/CHARLOTTE/'

    rg = Rain(path=path, name='Charlotte_CRN_gage_{YEAR}.csv', year=range(1993, 2016), ngages=71, freq='15min')
    rad = Rain(path=path, name='Charlotte_CRN_radar_{YEAR}.csv', year=range(2001, 2016), ngages=71, freq='15min')
    rg.ll_file = 'Charlotte_CRN_lat_lon.csv'
    rad.ll = rg.ll
    p = RadarGage(gage=rg, radar=rad)
    p.save_path = '../output/Charlotte/'
In [10]:
p.year = rad.year

Rain Rate over Period of Record

This red in this plot shows mean gage rain rate (mm/hr) at each 15min timestep over the period of record and the blue shows the mean radar rain rate at each timestep (radar data is only available for half the year).

In [7]:

To get a correlation, we plot measured hourly rain rate at each gage against hourly rain rate predicted by the radar.

In [11]:

Seasonal and diurnal patterns of rain rate.

In [12]:
In [13]:

Probability of wet days occuring

Plots of the mean probability that a day will be wet using a variety of thresholds. In this instance we take the mean probability for each month.

In [14]:
for mm, c in zip([.253, 1, 5, 25], ['red','blue','green','cyan']):
        rg.thresh = mm * rg.per_hour
        rg.plot_prob_wet(interval='seasonal', time_step='24H', base=12, color=c, save=False)
    plt.legend([.253, 1, 5, 25])
    title = 'Probability of wet 24 hours 1993-2015 with varying thresholds (gage)'
    plt.savefig(rg.save_path + title + '.png')

We can choose a threshold that plot probability of wetness for radar vs gage

In [16]:
p.thresh = 5 * p.per_hour
    p.plot_prob_wet(interval='seasonal', time_step='24H', base=12)

Distribution of rain rates for wet days

Boxplots to highlight the distribution of rainfall across months we'll use a threshold of 1mm for these:

In [13]:
rg.thresh = rg.per_hour*1

Across hours:

In [14]:

Mapping Rainfall

Once we are satisfied with the temporal components of the data, we are ready to try mapping the rain in space. As a first pass, we can just plot the mean rain rate at each location.

Rain rate

In [17]:
p.plot_rate(gage=p.list_gages(), bar=False, map=True, basemap=True, sharec=True)

Probability of wet 15 min (threshold=1mm)

In [18]:
p.thresh = 1*p.per_hour
    p.plot_prob_wet(gage=p.list_gages(), bar=False, map=True,  hide_title=True, basemap=True, sharec=True)

Rain rate given that it is raining (more than 1mm in 15min)

In [19]:
p.thresh = 1*p.per_hour
    e = Event(p.ll.join((p.wet_rate.mean(axis=1))))
    e.map_rain(basemap=True, sharec=True, hide_title=True)

We can use lowess to find the time of peak rainfall within the average day.

Time of peak rainfall

In [24]:
lmg = p.get_max_lowess(df=p.rate.gage, example_plot=False)
    lmg.columns = ['Gage lowess max']
    lmr = p.get_max_lowess(df=p.rate.radar, example_plot=False)
    lmr.columns = ['Radar lowess max']
In [27]:
e = Event(p.ll.join((lmg,lmr)))
    title = 'Map of Time of Peak Rainfall using lowess (f=.25) for 15min of {y}'
    e.map_rain(latlon=True, cmap='hsv', sharec=(0,24), save_path=rg.save_path,
               title=title.format(y=p.year), hide_title=True, basemap=True)

Month of peak rainfall (just for gage data)

Since we don't have complete years of radar data, we shouldn't use lowess to find the peaks because of edge effects.

In [105]:
lm = rg.get_max_lowess(example_plot=False, interval='seasonal')
In [106]:
e = Event(p.ll.join(lm))
    title = 'Map of Time of Peak Rainfall using lowess (f=.25) for 15min of {y}'
    e.map_rain(latlon=True, cmap='hsv', sharec=(1,13), hide_title=True, basemap=True)

Particular Events

In this case we are looking for the the rainiest 4 days of the entire dataset. Using the gage data we only really get tropical cyclones.

In [22]:
fig, axes = plt.subplots(nrows=2, ncols=2, figsize=(16,10))

The 1997 storm is a nice storm so we will check that one out first.

In [54]:
    rg.get_storm('1997-07-23 3:00', storm_end='1997-07-23 15:00')
    storm = Event(rg.storm)
In [55]:
    res = storm.res
In [56]:
fig, ax = plt.subplots(1,1,figsize= (10,6))
    sc = Event(res).krige(5, animated=True, ax=ax, vmin=0, vmax=100)

    def animate(i):
        global res
        Event(res).krige(i+5, animated=True, ax=ax, vmin=0, vmax=80)
    animation.FuncAnimation(fig, animate, frames=len(res.columns[5:]), interval=300, blit=True)

Once Loop Reflect