Hydrometeorology Research Group



In [5]:
from IPython.display import HTML

    HTML('''<script>
    code_show=true; 
    function code_toggle() {
     if (code_show){
     $('div.input').hide();
     } else {
     $('div.input').show();
     }
     code_show = !code_show
    } 
    $( document ).ready(code_toggle);
    </script>
    <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]:
%%capture
    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'
    rg.get_ll(cols=['lat','lon'])
    rad.ll = rg.ll
    p = RadarGage(gage=rg, radar=rad)
    p.get_nonan()
    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]:
rg.rate.mean(axis=1).plot(figsize=(16,6))
    rad.rate.mean(axis=1).plot()
    plt.show()
    

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

In [11]:
p.plot_correlation(time_step='1H')
    

Seasonal and diurnal patterns of rain rate.

In [12]:
p.plot_rate(interval='seasonal')
    
In [13]:
p.plot_rate(interval='diurnal')
    

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.title(title)
    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
    rg.plot_boxplots(interval='seasonal')
    

Across hours:

In [14]:
rg.plot_boxplots(interval='diurnal')
    

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
    p.get_wet_rate()
    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))
    rg.rate['1995-08-26':'1995-08-27'].mean(axis=1).plot(ax=axes[0,0])
    rg.rate['1997-07-22':'1997-07-23'].mean(axis=1).plot(ax=axes[0,1])
    rg.rate['2004-09-07':'2004-09-08'].mean(axis=1).plot(ax=axes[1,0])
    rg.rate['2008-08-26':'2008-08-27'].mean(axis=1).plot(ax=axes[1,1])
    plt.show()
    

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

In [54]:
%%capture
    rg.get_storm('1997-07-23 3:00', storm_end='1997-07-23 15:00')
    storm = Event(rg.storm)
    
In [55]:
storm.detrend()
    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)
    fig.colorbar(sc)

    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