Skip to content

Rainfall

Risk maps for the 100-year event of the 5-day precipitation.


Annual

Total Precipitation

1961-1990 1990-2019
Screenshot Screenshot

North-East Germany: No trends in annual precipitation between 1961-1990 and 1990-2019 caused by dynamical changes.

Extreme Precipitation

1961-1990 1990-2019
Screenshot Screenshot

North-East Germany: Significant incease in extreme precipitation between 1961-1990 and 1990-2019 caused by dynamical changes.

Ratio to Monthly Mean Sum

1961-1990 1990-2019
Screenshot Screenshot

North-East Germany: A 100 year event exceed a 3 month total precipitation.


Code

Importing

#!/usr/bin/env python
# coding: utf-8
import os
import datetime

import matplotlib  

matplotlib.use('Agg')

import numpy as N
import matplotlib.pyplot as P
from scipy import stats as S
from netCDF4 import Dataset,num2date
from mpl_toolkits.basemap import Basemap,maskoceans
from matplotlib.offsetbox import AnchoredText

P.rcParams["font.family"] = "serif"

Reading

file = '../csv/gen_1961-1990/gen_1961-1990.txt' 

jj=N.genfromtxt(file,usecols=(0),skip_header=1,dtype="I")
mm=N.genfromtxt(file,usecols=(1),skip_header=1,dtype="I")
dd=N.genfromtxt(file,usecols=(2),skip_header=1,dtype="I")
gw=N.genfromtxt(file,usecols=(3),skip_header=1,dtype='unicode')

nd = len(gw)

Processing

jo = N.arange(2200,2301,1);nj = len(jo)

for j in jo:

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

    for d in id:

        date = '%4i-%02i-%02i'%(jj[d],mm[d],dd[d])
        inp = '../csv/cdf/%02i-%s.nc'%(mm[d],gw[d])

        if os.path.isfile(inp):

              inp = '../csv/cdf/%02i-%s.nc'%(mm[d],gw[d])
              t = d

        else:

              inp = '../csv/cdf/%02i-%s.nc'%(mm[t],gw[t])

        out = '../csv/tmp/%s.nc'%(date)

        os.system("cdo -O -settaxis,%s,00:00:00,1day %s %s"%(date,inp,out))

    os.system("cdo -O mergetime ../csv/tmp/*.nc ../csv/gen_1961-1990/tw/%i.nc"%(j))
    os.system("rm ../csv/tmp/*.nc")

Plotting

ko = {'tot':['prcp tot (mm)'],'ext':['prcp ext (mm)'],'ext-tot':['ext = x monthly tot']}


for j in ['1961-1990','1990-2019']:
    for s in ['ja']:

      if(s=='ja'): nt = 360
      if(s=='sh'): nt = 180

      for k in ko:

        if(k=='ext-tot'):

             file = '../csv/gen_%s/dw/prc_ext_%s.nc'%(j,s)

             nc = Dataset(file,'r')

             lons = N.array(nc.variables['longitude'][:])
             lats = N.array(nc.variables['latitude'][:])
             maxs = N.array(nc.variables['pr'][:])
             maxs = maxs*3600.*24.

             nc.close()

             file = '../csv/gen_%s/dw/prc_tot_%s.nc'%(j,s)

             nc = Dataset(file,'r')

             lons = N.array(nc.variables['longitude'][:])
             lats = N.array(nc.variables['latitude'][:])
             avgs = N.array(nc.variables['pr'][:])
             avgs = 30.*avgs*24.*3600./3.
             nc.close()

             dats = (maxs[0,:,:]/avgs[0,:,:])

             lev = [1,2,3,4]
             col = 'cool' 

        if(k=='ext'):

             file = '../csv/gen_%s/dw/prc_ext_%s.nc'%(j,s)

             nc = Dataset(file,'r')

             lons = N.array(nc.variables['longitude'][:])
             lats = N.array(nc.variables['latitude'][:])
             dats = N.array(nc.variables['pr'][:])
             dats = dats*24.*3600.
             nc.close()

             dats = dats[0,:,:]

             #print (dats.max(),fac)

             lev = [50,100,150,200,250,300] 
             col = 'BuPu'

        if(k=='tot'):

             file = '../csv/gen_%s/dw/prc_tot_%s.nc'%(j,s)

             nc = Dataset(file,'r')

             lons = N.array(nc.variables['longitude'][:])
             lats = N.array(nc.variables['latitude'][:])
             dats = N.array(nc.variables['pr'][:])
             dats = dats*24.*3600.*nt/3.

             nc.close()

             dats = dats[0,:,:]

             if(s=='ja'): lev = [500,1000,1500,2000,2500] 
             if(s=='sh'): lev = [250,500,750,1000,1250]

             col = 'YlGnBu'


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

        map = Basemap(projection='cyl',llcrnrlat=35.,urcrnrlat=70.,llcrnrlon=-15,urcrnrlon=30,resolution='l')

        map.drawcoastlines(color='k')
        map.drawcountries(linewidth=1,color='k')
        map.drawlsmask(land_color = "lightgray",ocean_color="lightgray",resolution = 'h')

        lon,lat = map(*N.meshgrid(lons,lats))

        dat = maskoceans(lon,lat,dats)

        CM = map.contourf(lon,lat,dat,levels=lev,cmap=P.get_cmap(col),extend='both')

        cbar = map.colorbar(CM,location='right',extend='max',pad="1%",drawedges=True,ticks=lev)
        cbar.set_label(ko[k][0],fontsize=12,weight='bold')
        cbar.ax.tick_params(labelsize=8)

        at = AnchoredText(r"$\varnothing$ %.1f"%N.nanmean(dat),prop=dict(size=10),frameon=True,loc='upper left') 
        at.patch.set_boxstyle("round,pad=0.,rounding_size=0.2")
        ax.add_artist(at)

        at = AnchoredText(j,prop=dict(size=10),frameon=True,loc='lower left')
        at.patch.set_boxstyle("round,pad=0.,rounding_size=0.2")
        ax.add_artist(at)

        P.savefig('./img/risk_%s_%s_%s.png'%(k,s,j),dpi=240,transparent=False,bbox_inches='tight',pad_inches=0.0)