GSMaP

import matplotlib.pyplot as P
import numpy as N
from matplotlib.offsetbox import AnchoredText
from scipy import signal,stats

params = {'legend.fontsize': 8,\
          'font.family': 'serif'}
P.rcParams.update(params)
````

### Defining

```python
def ma(a,n=3):
    ret=N.cumsum(a,dtype=float)
    ret[n:]=ret[n:]-ret[:-n]
    return ret[n-1:]

file = '../../data/gsmap/gsmap_amman.csv'
#file = '../../data/gsmap/gsmap_wadimusa.csv'
#file = '../../data/gsmap/gsmap_irbid.csv'
#file = '../../data/gsmap/gsmap_potsdam.csv'
#file = '../../data/gsmap/gsmap_aqaba.csv'

tt = N.genfromtxt(file,usecols=(0),delimiter=',',skip_header=1,dtype='str')
r0 = N.genfromtxt(file,usecols=(1),delimiter=',',skip_header=1,dtype='f')
r1 = N.genfromtxt(file,usecols=(2),delimiter=',',skip_header=1,dtype='f')
r2 = N.genfromtxt(file,usecols=(3),delimiter=',',skip_header=1,dtype='f')

#r0=r2

nh = len(tt)

r0 = r0[:nh-3]

hh = []
dd = []
mm = []
jj = []

for t in tt[3:]:

    jj.append(int(t[0:4]))
    mm.append(int(t[5:7]))
    dd.append(int(t[8:10]))
    hh.append(int(t[11:13]))

jj = N.array(jj)
mm = N.array(mm)
dd = N.array(dd)
hh = N.array(hh)

ri = N.array([1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24])
nr = len(ri)

i0 = N.where((jj>=2001)&(jj<=2010))
i1 = N.where((jj>=2011)&(jj<=2021))
ii = N.where((jj>=2001)&(jj<=2021))


rint0 = {}
rint1 = {}
rint = {}


for i in ri:

    rint0[str(i)] = []

    tmp = ma(r0[i0],i)

    for t in tmp:

        rint0[str(i)].append(t)

    rint1[str(i)] = []

    tmp = ma(r0[i1],i)

    for t in tmp:

        rint1[str(i)].append(t)

    rint[str(i)] = []

    tmp = ma(r0[ii],i)

    for t in tmp:

        rint[str(i)].append(t)

import numpy as np
from scipy.stats import genextreme

jo = N.arange(2001,2022,1)
nj = len(jo)
ho = N.arange(1,25,1)
nh = len(ho)

tmp = N.zeros((nj,nh),float)

for j in range(nj):

    id = N.where(jj==jo[j])[0]

    for h in range(nh):

        tmp[j,h] = N.max(ma(r0[id],ho[h]))

        print (jo[j],ho[h],tmp[j,h]) 

return_periods = np.array([2,5,10,20,50,100])
nr = len(return_periods)

lev = N.zeros((nr,nh),float)

for h in range(nh):

    data = tmp[:,h] 

    shape, loc, scale = genextreme.fit(data)
    param = genextreme.fit(data)

    P.figure(figsize=(6,4))
    ax = P.subplot(111)

    xx = N.linspace(0,80, num=1000)
    yy = genextreme.pdf(xx, shape, loc, scale)

    hist, bins = N.histogram(data, bins=20, range=(0,80), density=True)

    P.bar(bins[:-1], hist, width = 3, align='edge',fc='dodgerblue',ec='k')

    P.plot(xx,yy,'r',lw=2)
    P.xlim(0,80)
    P.ylim(0,0.12)

    P.xlabel('Precipitation (mm)',fontsize=12,weight='bold')
    P.ylabel('Density',fontsize=12,weight='bold')

    at = AnchoredText('%02ih: shape=%.1f|loc=%.1f|scale=%.1f'%(ho[h],shape,loc,scale),prop=dict(size=10,weight='bold'),frameon=True,loc='upper right')
    at.patch.set_boxstyle("round,pad=0.,rounding_size=0.1")
    ax.add_artist(at)

    P.savefig('./img/gev%02i.png'%ho[h],dpi=80,transparent=True,bbox_inches='tight',pad_inches=0.0)

    return_levels = genextreme.isf(1/return_periods, shape, loc, scale)

    lev[:,h] = return_levels

f = open('./includes/return_Amman.md','w')

f.write('Amman|%i yr|%i yr|%i yr|%i yr|%iyr|%i yr|GEV\n'%(2,5,10,20,50,100))
f.write('---|---|---|---|---|---|---|---\n')

#lev = N.median(lev,2)

for h in range(nh):

    img = '![](../img/gev%02i.png){: style="height:40px;width:100px"}'%ho[h]

    f.write('**%02i h**|%.1f|%.1f|%.1f|%.1f|%.1f|%.1f|%s\n'%(ho[h],lev[0,h],lev[1,h],lev[2,h],lev[3,h],lev[4,h],lev[5,h],img))

f.close()

fig = P.figure(figsize=(12,7))
ax = P.subplot(221)

blue = dict(marker='o', markeredgecolor='b', markersize=3)
reds = dict(marker='o', markeredgecolor='r', markersize=3)


x = []
y = []
z = []

for i in ri:

    tmp = N.array(rint0[str(i)])
    id = N.where(tmp>i)[0]

    box = P.boxplot(tmp[id],positions=[i-0.15],showfliers=True,zorder=50,flierprops=blue)
    P.setp(box['boxes'], color='b')
    P.setp(box['whiskers'], color='b')

    tmp = N.array(rint1[str(i)])
    id = N.where(tmp>i)[0]

    box = P.boxplot(tmp[id],positions=[i+0.15],showfliers=True,zorder=50,flierprops=reds)
    P.setp(box['boxes'], color='r')
    P.setp(box['whiskers'], color='r')

P.text(1,110,'2001-2010',color='b')
P.text(1,100,'2011-2021',color='r')

P.grid(lw=0.5,ls=':')

P.xlim(0,25)
P.ylim(0,120)
P.xticks(ri,N.array(ri,str))
P.ylabel('Accumulated Rainfall [mm]',fontsize=12,weight='bold')
P.xlabel('Duration [h]',fontsize=12,weight='bold')

at = AnchoredText('GSMaP: Irbid 2001-2021',prop=dict(size=8,weight='bold'),frameon=True,loc='lower right')
at.patch.set_boxstyle("round,pad=0.,rounding_size=0.1")
ax.add_artist(at)

ax = P.subplot(222)

for i in ['1','6','12','18','24']:

    y = N.array(rint[str(i)])
    #id = N.where(y>0)[0]

    y = N.sort(y)
    x = N.arange(1,len(y)+1,1)[::-1]

    P.plot(x,y,'k',lw=0.5,label='%s h'%i)

P.legend(loc=2,shadow=True)

P.ylim(1,120)
P.xlim(10000,10)

P.semilogx()

P.grid(which='major',lw=1.0,ls=':')
P.grid(which='minor',lw=0.5,ls=':')

P.subplot(212)

xnew = N.linspace(-1,24,num=100)

#ii = N.where((r0>1)&(jj>=2001)&(jj<=2010))[0]
#density = stats.kde.gaussian_kde(hh[ii])
#P.plot(xnew,density(xnew),'b',lw=2,label='2001-2010')

#ii = N.where((r0>1)&(jj>=2011)&(jj<=2021))[0]
#density = stats.kde.gaussian_kde(hh[ii])
#P.plot(xnew,density(xnew),'r',lw=2,label='2011-2021')

#ii = N.where((r0>1)&(jj>=2003)&(jj<=2021))[0]
#density = stats.kde.gaussian_kde(hh[ii])
#P.plot(xnew,100.*density(xnew),'r',lw=2,label='2003-2021')

nn = []
tt = N.arange(0,25,1)

for t in tt:

    ii = N.where((r0>2)&(jj>=2003)&(jj<=2021)&(hh==t))[0]

    nn.append(len(ii))

P.bar(tt,nn,0.8,color='dodgerblue',ec='k')

P.xlim(-1,24)
P.ylim(0,61)
P.xticks([0,3,6,9,12,15,18,21])

#P.legend(loc=2,shadow=True)
P.grid(which='major',lw=1.0,ls=':')

P.xlabel('Time of Day [h]',fontsize=12,weight='bold')
P.ylabel('Number of Events > 2 mm/h',fontsize=12,weight='bold')

P.tight_layout()

P.savefig('./img/gsmap.png',dpi=240,transparent=False,bbox_inches='tight',pad_inches=0.0)