Interactive comment on “ Precipitation extremes on multiple time scales — Bartlett – Lewis Rectangular Pulse Model and Intensity – Duration – Frequency curves

This paper investigates the ability of the original Bartlett-Lewis model for estimating extreme rainfall at various levels of aggregation. Unfortunately, the paper is not very novel. It is already known for a long period that the Bartlett-Lewis (BL) models have problems in reproducing extremes, especially at shorter aggregation levels. It is not clear why the authors chose for the Original Bartlett-Lewis (OBL) model, while the Modified Bartlett-Lewis (MBL) model or one of the later versions (e.g. Onof and Wheather, 1994) that were further optimized for addressing the problem of the undergeneration of extremes.


Introduction
Precipitation is one of the most important atmospheric variables.Large variations on spatial and temporal scales are observed, i.e. from localized thunderstorms lasting a few tens of minutes up to mesoscale hurricanes lasting for days.Precipitation on every scale affects everyday life: short but intense extreme precipitation events challenge the drainage infrastructure in urban areas or might put agricultural yields at risk; long-lasting extremes can lead to flooding (Merz et al., 2014).Both short intense and long-lasting large-scale rainfall can lead to costly damage, e.g. the floodings in Germany in 2002 and 2013 (Merz et al., 2014), and are therefore the object of much research.
Risk quantification is based on an estimated frequency of occurrence for events of a given intensity and duration.This information is typically summarized in an intensityduration-frequency (IDF) relationship (e.g.Koutsoyiannis et al., 1998), also referred to as IDF curves.These curves are typically estimated from long-observed precipitation time series, mostly with a sub-daily resolution to also include short durations in the IDF relationship.These are indispensable for some hydrological applications -e.g.extreme precipitation characteristics are derived from IDF curves for the planning, design, and operation of drainage systems, reservoirs, and other hydrological structures.One way to obtain IDF curves is modelling block maxima for a fixed duration with the generalized extreme value distribution (GEV; e.g.Coles, 2001).Here, we employ a parametric extension to the GEV which allows a simultaneous modelling of extreme precipitation for all durations (Koutsoyiannis et al., 1998;Soltyk et al., 2014).
Published by Copernicus Publications on behalf of the European Geosciences Union.
C. Ritschel et al.: IDF for precipitation in the BLRPM Due to a limited availability of observed high-resolution records with adequate length, simulations with stochastic precipitation models are used to generate series for subsequent studies (e.g.Khaliq and Cunnane, 1996;Smithers et al., 2002;Vandenberghe et al., 2011).The advantages of stochastic models are their comparatively simple formulation and their low computational cost, enabling them to quickly generate large ensembles of long precipitation time series.A review of these models is given in, e.g., Onof et al. (2000) and Wheater et al. (2005).A canonical class of sub-daily stochastic precipitation models are Poisson cluster models with the original Bartlett-Lewis (OBL) model as a prominent representative (Rodriguez-Iturbe et al., 1987, 1988;Onof and Wheater, 1994b;Wheater et al., 2005).The OBL model has been shown to be able to well reproduce certain characteristics found in precipitation observations (Rodriguez-Iturbe et al., 1987).
Due to the high degree of simplification of the precipitation process, known drawbacks of the OBL model include the inability to reproduce the proportion of dry periods, as reported by Rodriguez-Iturbe et al. (1988) and Onof (1992), and underestimation of extremes, as found by, for example, Verhoest et al. (1997) and Cameron et al. (2000), especially for shorter durations.Furthermore, problems occur for return levels with associated periods longer than the time series used for calibrating the model (Onof and Wheater, 1993).Several extensions and improvements to the model have been made.Rodriguez-Iturbe et al. (1988) introduced the randomized parameter Bartlett-Lewis model, allowing for different types of cells.Improvements in reproducing the probability of zero rainfall and capturing extremes have been shown for this model (Velghe et al., 1994).A gamma-distributed intensity parameter and a jitter were introduced by Onof and Wheater (1994b) for more realistic irregular cell intensities.Nevertheless, problems still remain.Verhoest et al. (2010) discussed the occurrence of infeasible (extremely long lasting) cells, and a too severe clustering of rain events was found by Vandenberghe et al. (2011).Including third-order moments in the parameter estimation showed an improvement in the Neyman-Scott model extremes (Cowpertwait, 1998).For the Bartlett-Lewis variant Kaczmarska et al. (2014) found that a randomized parameter model shows no improvement in fit compared to the OBL model, in which skewness was included in the parameter estimation.Furthermore, an inverse dependence between rainfall intensity and cell duration showed improved performance, especially for extremes at short timescales (Kaczmarska et al., 2014).Here, we focus on the OBL model with and without the thirdorder moment included.This model is still part of a wellestablished class of precipitation models, and the reduced complexity is appealing as it enables it to be used in a nonstationary context (Kaczmarska et al., 2015).
The OBL model and IDF relationships are of particular interest to hydrological modelling and impact assessment.In the following, we address three research questions by means of a case study: 1. Is the OBL model able to reproduce the intensityduration relationship found in observations?
2. How are IDF curves affected by very rare extreme events which are unlikely to be reproduced with the OBL model in a reasonably long simulation?
3. Is the parametric extension to the GEV a valid approach to obtain IDF curves?
For the class of multi-fractal rainfall models, question (1) has been addressed by Langousis and Veneziano (2007).
Section 2 gives a short overview of the OBL model used here, and Sect. 3 briefly explains IDF curves and how to obtain them from precipitation time series.This is followed by a description of the data used (Sect.4).The results section (Sect.5) starts with a description of the estimated OBL model parameters for the case study area (Sect.5.1), then discusses the ability of this model to reproduce IDF curves (Sect.5.2) and the influence of a rare extreme (Sect.5.3), and closes with a comparison of the duration-dependent GEV approach to IDF curves with individual duration GEV quantiles (Sect.5.4).Section 6 discusses results and concludes the paper.

Bartlett-Lewis rectangular pulse model
From early radar-based observations of precipitation, a hierarchy of spatio-temporal structures was suggested by Pattison (1956) and Austin and Houze (1972): intense rain structures (denoted as cells) tend to form in the vicinity of existing cells and thus cluster into larger structures -so-called storms or cell clusters.The occurrence of these cells and storms (cell clusters) can be described using Poisson processes.This makes the Poisson cluster process a natural approach to stochastic precipitation modelling.
The idea of modelling rainfall with stochastic models has existed since Le Cam (1961) modelled rain gauge data with a Poisson cluster process.Later, Waymire et al. (1984) and Rodriguez-Iturbe et al. (1987) continued the development of this type of model.Poisson cluster rainfall models are characterized by a hierarchy of two layers of Poisson processes.
Similarly to various others studies (e.g.Onof and Wheater, 1994a;Kaczmarska, 2011;Kaczmarska et al., 2015), we chose the original Bartlett-Lewis (OBL) model, a popular representative of Poisson cluster models with a set of five physically interpretable parameters (Rodriguez-Iturbe et al., 1987).At the first level cell clusters (storms) are generated according to a Poisson process with a cluster generation rate λ and an exponentially distributed lifetime with expectation 1/γ .Within each cell cluster, cells are generated according to a Poisson process with cell generation rate β and exponentially distributed lifetime with expectation 1/η, hence the Figure 1.Scheme of the OBL model.A similar scheme can be found in Wheater et al. (2005) term Poisson cluster process.Associated with each cell is a precipitation intensity which is constant during cell lifetime, exponentially distributed with mean µ x .The constant precipitation intensity in a cell gave rise to the name rectangular pulse.Model parameters are summarized in the parameter vector θ = {λ, γ , β, η, µ x }.
Simulations with the OBL model are in continuous time on the level of storms and cells.We aggregate the resulting cell rainfall series to hourly time series.Figure 1 shows a sketch of the OBL model illustrating the two levels: the start of cell clusters (storms) are shown as red dots; within each cluster, cells (blue rectangular pulses) are generated during the cluster's lifetime starting from the cluster origin.The cells' lifetime is shown as horizontal and their precipitation intensity as vertical extension in Fig. 1; cells can overlap.This continuous-time model yields a sequence of pulses (cells) with associated intensity; see Fig. 2 (top panel).Adding up the intensity of overlapping pulses yields a continuous-time step function (Fig. 2, middle panel).Although time continuous, this function is not continuously differentiable in time due to the rectangular pulses.Observational time series are typically also not continuously differentiable, as they are discretized in time.Summing up the resulting continuous series for discrete time intervals makes it comparable with observations and renders unimportant the artificial jumps from the rectangular pulses present in the continuous series, Fig. 2 (bottom layer).
An alternative to the Bartlett-Lewis process is the Neyman-Scott process (Neyman and Scott, 1952).The latter is motivated from observations of the distribution of galaxies in space.In the Neyman-Scott process cells are distributed around the centre of a cell cluster.Both are prototypical models for sub-daily rainfall and are discussed in more detail in Wheater et al. (2005).
Due to known drawbacks of the OBL model, several improvements and extensions have been made in the past: Rodriguez-Iturbe et al. (1988) introduced the random param-eter model, allowing for different types of cells, and additionally Onof and Wheater (1994b) used a jitter and a gammadistributed intensity parameter to account for a more realistic irregular shape of the cells.Cowpertwait et al. (2007) improved the representation of sub-hourly timescales by adding a third layer, pulses, to the model.Non-stationarity has been addressed by Salim and Pawitan (2003) and Kaczmarska et al. (2015).Applications of these kind of models include the implementing of copulas to investigate wet and dry extremes (Vandenberghe et al., 2011;Pham et al., 2013), regionalization (Cowpertwait et al., 1996a, b;Kim et al., 2013), and accounting for interannual variability (Kim et al., 2014).
Parameter estimation for the OBL model is by no means trivial.The canonical approach is a method-of-momentbased estimation (Rodriguez-Iturbe et al., 1987) using the objective function This function relates moments of precipitation sums τ i (θ ) derived from the model with parameters θ to empirical moments M i from the time series.The set of k moments M i is typically chosen from the first and second moments obtained for different durations h.Here, we use the mean at 1 h aggregation time, the variances, the lag-1 autocovariance function, and the probability of zero rainfall for h ∈ {1, 3, 12 and 24 h}, similar to Kim et al. (2013), and thus end up with k = 13 moments M i .Their analytic counterparts τ i (θ ) are derived from the model.The weights in the objective function were chosen to be w i=1 = 100 and w i =1 = 1, similar to Cowpertwait et al. (1996a), emphasizing the first moment M 1 (mean) 100 times more than the other moments.
It turns out that Z(θ ; M) has multiple local minima and optimization is not straightforward.To avoid local minima, the optimization is repeated many times with different initial guesses for the parameters sampled from a range of feasible values in parameter space using the Latin hypercube sampling algorithm (McKay et al., 1979).As only positive model parameters are meaningful, optimization is performed on log-transformed parameters.Similarly to Cowpertwait et al. (2007), we use a symmetric objective function A few tests indicate that the symmetric version is robust and faster in the sense that fewer iterations are needed to ensure convergence into the global minimum (not shown).Numerical optimization techniques based on gradient calculations, e.g.Nelder-Mead (Nelder and Mead, 1965) or BFGS (Broyden, 1970;Fletcher, 1970;Goldfarb, 1970;Shanno, 1970), are typically used.For the current study, we use R's optim() function choosing L-BFGS-B as the underlying optimization algorithm (R Core Team, 2016) and 100 differ- ent sets of initial guesses for the parameters sampled using the Latin hypercube algorithm.
Following studies by Cowpertwait (1998) and Kaczmarska et al. (2014), we include the third moment in the parameter estimation using analytical expressions derived by Wheater et al. (2006), replacing the probability of zero rainfall in the objective function.Thus, still 13 moments are used to calibrate the OBL model.Due to comparability with other studies, most of our analyses will not include the third moment though.A comparison between IDF curves of the model calibrated with the third moment and with the probability of zero rainfall will be carried out, to discuss the effect of including the third moment.
Models of this type suffer from parameter nonidentifiability, meaning that qualitatively different sets of parameters lead to minima of the objective function with comparable values (Verhoest et al., 1997).A more detailed view on global optimization techniques and comparisons between different objective functions is given in Vanhaute et al. (2012).
During this work the authors developed and published the R package BLRPM (Ritschel, 2017).The package includes functions for simulation and parameter estimation.
3 Intensity-duration-frequency IDF curves show return levels (intensities) for given return periods (inverse of frequencies) as a function of rainfall duration.Their formulation goes back to Bernard (1932).They are frequently used for supporting infrastructure risk assessment (e.g.Simonovic and Peck, 2009;Cheng and AghaKouchak, 2013).IDF curves are an extension to classical extreme value statistics.The latter aims at better characterizing the tails of a distribution by using parametric models derived from limit theorems (e.g.Embrechts et al., 1997).There are two main approaches: modelling block maxima (e.g.maxima out of monthly or annual blocks) with the GEV distribution or modelling threshold excesses with the generalized Pareto distribution (GPD) (e.g.Coles, 2001;Embrechts et al., 1997).We chose the block maxima approach with the general extreme value distribution as parametric model for the block maxima z.The GEV is characterized by the location parameter µ, the scale parameter σ and the shape parameter ξ .These can be estimated from block maxima using a maximum-likelihood estimator (e.g.Coles, 2001).Here, we use maxima from monthly blocks.To avoid mixing maxima from different seasons, a set of GEV parameters is estimated for all maxima from January, another set for all maxima from February and so on.For a given month, GEV parameters are estimated for various durations, e.g.d ∈ {1, 6, 12, 24, 48 h, . ..}.For a specific return period T = 1/(1 − p), with p denoting the non-exceedance probability, a parametric model can be fitted to the corresponding p quantiles Q p,d from GEV distributions for different durations d (e.g.Koutsoyiannis et al., 1998).
with θ , η, and σ being independent of the duration d.The parameter η quantifies the slope of the IDF curve in the main region and θ controls the deviation of the power-law behaviour for short durations.Furthermore, location is reparametrized by µ = µ/σ d , which is now independent of the durations d as well as of the shape parameter ξ .This leads to the following formulation of a duration-dependent GEV distribution: which allows consistent modelling of rainfall maxima across different durations d using a single distribution at the cost of only two additional parameters.These parameters can be analogously estimated by maximum likelihood (Soltyk et al., 2014).To avoid local minima when optimizing the likelihood, we repeat the optimization with different sets of initial guesses for the parameters, sampled again according to a Latin hypercube scheme.This method of constructing IDF curves is consistent in the sense that curves for different return periods cannot cross.We refer to this approach as the duration-dependent GEV approach (dd-GEV).During this work the authors developed and published the R package IDF (Ritschel et al., 2017a).The package includes functions for estimating IDF parameters based on the dd-GEV approach given a precipitation time series, and plots the resulting IDF curves.However, the data points for different durations are dependent (as they are derived from the same underlying highresolution data set by aggregation), and thus the i.i.d.assumption required for maximum-likelihood estimation is not fulfilled.Consequently, confidence intervals are not readily available from asymptotic theory; however, they can be estimated by bootstrapping.

Data
A precipitation time series from the station Botanical Garden in Berlin-Dahlem, Berlin, Germany, is used as a case study.A tipping-bucket records precipitation amounts at 1 min resolution.For the analysis at hand, a 13-year time series with 1 min resolution from the years 2001-2013 is available.The series is aggregated to durations d ∈ {1, 2, 3, 6, 12, 24, 48, 72, 96 h} yielding nine time series with different temporal resolution.IDF parameters are estimated using annual maxima for each month of the year individually using all nine duration series.

Estimation of OBL model parameters
Minimizing the symmetric objective function (Eq.2) yields OBL model parameter estimates individually for every month of the year, shown in Fig. 3 and explicitly given in Table A2 in Appendix A. The resulting OBL model parameters are reasonable compared to observed precipitation characteristics: during summer months, we observe very intensive cells ( μx between 4 and 8 mm h −1 ).However, in June and August, storm duration is relatively short ( γ between 0.25 and 0.35 h −1 ), which can be interpreted as short but heavy thunderstorms, as are typically observed in this region in summer (Fischer et al., 2018).Conversely, in winter small intensities and long storm durations correspond to stratiform precipitation patterns, typically dominating the winter precipitation in Germany.The storm generation rate λ shows only minor seasonal variation.
With the OBL model parameter estimates (Table A2, Appendix A), 1000 realizations with the same length as the observations (13 years) are generated.From both the original precipitation series and the set of simulated time series, we derive a set of statistics for model validation.The first moment M 1 -the mean -is very well represented (not shown), as it enters the objective function with weight w 1 = 100 compared to weights of 1 for the other statistics.Figure 4 shows the variance for 6-hourly aggregation and the probability of zero rainfall; for all months the 6 h variances of simulated and observed series are in good agreement.This is particularly noteworthy as the 6-hourly aggregation was not used for parameter estimation.Similar to previous studies (e.g.Onof and Wheater, 1994a), the model fails to reproduce the probability of zero rainfall, here for instance shown for the 12hourly aggregation.The model mainly overestimates it, and therefore has shortcomings in the representation of the time distribution of events (Rodriguez-Iturbe et al., 1987;Onof and Wheater, 1994a).
An important aspect for hydrological applications is the model's ability to reproduce extremes on various temporal www.hydrol-earth-syst-sci.net/21/6501/2017/ Hydrol.Earth Syst.Sci., 21, 6501-6517, 2017 6506 C. Ritschel et al.: IDF for precipitation in the BLRPM q q q q q q q q q q q q 0.00 0.02 0.04 Mon λ [1 q q q q q q q q q q q q Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec 0 10 20 30 40 1 γ [h] q q λ 1 γ q q q q q q q q q q q q 0.0 0.4 0.8 Mon β [1 q q q q q q q q q q q q Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec 0 0.5 1 1.5 1 η [h] q q β 1 η q q q q q q q q q q q q 0 2 4 6 8  scales.This behaviour is investigated in the next section with the construction of IDF curves.

IDF curves from OBL model simulations
Monthly block maxima for every month in the year are drawn for various durations (1, 3, 6, 12, 24, 48, 72, 96 h) from the observational time series and 1000 OBL model simulations of the same length.This is the basis for estimating GEV distributions for individual durations, as well as for constructing dd-GEV IDF curves.IDF curves for Berlin-Dahlem obtained from observation are shown as dotted lines in Fig. 5 for January, April, July, and October for the 0.5 quantile (2-year return period, red), 0.9 quantile (10-year return period, green), and the 0.99 quantile (100-year return period, blue).
Analogously, IDF curves are derived from 1000 simulations of the OBL model precipitation series; see Sect.5.1.The coloured shading in Fig. 5 gives the range of variability (5 to 95 %) for these 1000 curves with the median highlighted as a dotted line.Except for January, the curves obtained directly from the observational series can be found within the range of variability of curves derived from the OBL model.The main IDF features from observations are well reproduced by the OBL model: the power-law-like behaviour (straight line in the double-logarithmic representation) in July extending almost across the full range of durations shown, as well as the flattening of the IDF curves for short durations for April and September.The relative differences in IDF curves given in Fig. B1 (Appendix B) suggest a tendency for the OBL model to underestimate extremes, particularly for large return levels and short durations, similar to results found by, for example, Verhoest et al. (1997) and Cameron et al. (2000).Figure 6 shows the relative difference = dd-GEV OBL − dd-GEV obs dd-GEV obs × 100 % (6) between IDF curves (dd-GEV) derived from the OBL model dd-GEV OBL including the third moment in parameter estimation (red lines) or alternatively using the probability of zero rainfall to calibrate the model (blue lines), and directly from the observational time series dd-GEV obs for July and two quantiles: (a) 0.5 and (b) 0.99.Including the third moment in parameter estimation slightly improves the model extremes for July for all durations and both short and long return periods.Nevertheless, those promising results could not be found for all months (not shown), and thus we cannot conclude that including the third moment in parameter estimation improves extremes in the OBL model, in contrast to findings for the Neyman-Scott variant (Cowpertwait, 1998).
We interpret the different behaviour for short durations (flattening versus continuation of the straight line) for summer (July) and the remaining seasons as a result of different mechanisms governing extreme precipitation events: while convective events dominate in summer, frontal and thus more large-scale events dominate in the other seasons.
As an example, we show segments of time series including the maximum observed/simulated rainfall in July for durations 1, 6, and 24 h as observed (RR obs ) and simulated (RR OBL ) in Fig. 7.
Parts of the observed and simulated rainfall time series corresponding to the extreme events for the three different durations are shown in the left and right column, respectively.Additionally, the middle column shows the simulated storms and cells generating this extreme event in the artificial time series.As an example, we only show one single model simulation.Visual inspection of several other simulated series support the main features.For all durations, the extremes are a result of a single long-lasting cell with high intensity.In contrast to an analysis based on the random parameter BL model (Verhoest et al., 2010), these cells are neither unrealistically long nor have an unrealistically high intensity.
For January, IDF curves from observations and OBL model simulations exhibit large discrepancies: for all durations, the 0.99 quantile (100-year return level) is above the range of variability from the OBL model, and the 0.5 quantile (2-year return level) is below for small durations.This implies that the shape of the extreme value distribution characterized by the scale σ and shape parameter ξ differs between the two cases.This is likely due to the winter storm Kyrill hitting Germany and Berlin on 18 and 19 January in 2007 (Fink et al., 2009).We suppose that this rare event is not sufficiently influential to impact OBL model parameter estimation but does affect the extreme value analysis.For the latter only the one maximum value per month is considered.In fact, the shape parameter ξ estimated from the observational time series shows a large value compared to the other months; in contrast, this value is estimated to be around zero from OBL model simulations.The following section investigates this hypothesis by excluding the precipitation events due to Kyrill.
We furthermore find that the OBL model is generally able to reproduce the observed seasonality in IDF parameters; see Fig. 8.For all parameters, the direct estimation (blue) is mostly within the range of variability of the OBL model simulations.For σ , θ , and η, the direct estimation (blue line) features a similar seasonal pattern as the median of the OBL model (red line), whereas for ξ , the direct estimation is a lot more erratic than the median (red).As the GEV shape parameter is typically difficult to estimate (Coles, 2001), this erratic behaviour is not unexpected, and 11 out of 12 months stay within the expected inner 90 % range of variability.

Investigation of the impact of a rare extreme event
The convective cold front passage of Kyrill accounted for a maximum intensity of 24.8 mm rainfall per hour, whereas the next highest value of the remaining Januaries would be 4.9 mm rainfall per hour in 2002, thus being more than . IDF curves obtained via dd-GEV for (a) January, (b) April, (c) July and (d) October: 0.5 (red), 0.9 (green), and 0.99 quantiles (blue) corresponding to 2-, 10-, and 100-year return periods, respectively.Solid lines are derived directly from the Berlin-Dahlem time series.Coloured shadings mark the central 90 % range of variability of IDF curves obtained in the same manner with the same colour code but from 1000 OBL model simulations (Sect 5.1); the dotted lines mark the median of these curves.5 times lower than for Kyrill.We construct another data set without the extreme event due to Kyrill, i.e. without the year 2007.The intention of this experiment is not to motivate removal of an "unsuitable" value.We rather want to show that the OBL model is generally able to reproduce extremes; it is, however, not flexible enough to account for a single event with magnitude far larger than the rest of the time series.Based on the model with parameters estimated from observations with and without the year 2007 (observed), we obtain return periods for the event "Kyrill" for different durations and find this event to be very rare, especially on short timescales (1-3 h); see Table 1.
For this data set, we estimate the OBL model parameters and again simulate 1000 time series with these new parameters.The simulated time series were also reduced in length by 1 year, containing 12 years of rainfall in total.From those precipitation time series, we constructed the dd-GEV IDF curves; see Fig. 9 (right).Without the extreme events due to Kyrill, the OBL model performs in January as well as in the other months with respect to reproducing the IDF relations.In particular, the spread between the 0.5 quantile (2-year return level) and the 0.99 quantile (100-year) return level is reduced, as are the absolute values of extreme quantiles as well; see Fig. 9, left and right panel.Note the different scales for the intensity axes.

Comparing dd-GEV IDF curves to individual duration GEV
In the frame of a model-world study, long time series simulated with the OBL model can be used to investigate the adequacy of the dd-GEV model conditional on the simulated series.To this end, we compare the resulting IDF curves to a GEV distribution obtained for various individual durations.The basis is a set of 1000-year simulations with the OBL model with parameters optimized for Berlin-Dahlem.For a series of this length, we expect to obtain quite accurate (low variance) results for both the dd-GEV IDF curve and the GEV distributions for individual durations.However, sampling uncertainty is quantified by repeatedly estimating the desired quantities from 50 repetitions.The resulting dd-GEV IDF curves are compared to the individual duration GEV distribution in Fig. 10 for January (a) and July (b).
For most durations in January and July, the dd-IDF curves are close to the quantiles of the individual duration GEV distributions.Notable differences appear for small durations and large quantiles (return levels for long return periods); particularly in January, the dd-GEV IDF model overestimates the 10-year and 100-year return levels (duration of 1 h); in July, this effect seems to be present as well, but smaller in size.This is accompanied by a slight underestimation of the dd-GEV IDF for durations of 2 to 6 h in July and 3 to 6 h in January, most visible for the 0.99 quantile (100-year return level).Both effects together suggest that the flattening of the dd-GEV IDF for small durations is not sufficiently well represented.This could be due to deficiencies in the model for the duration-dependent scale parameter (Eq.4) but might also be a consequence of an inadequate sampling of durations (d ∈ {1, 6, 12, 24, 48 h, . . .}) to be used to estimate the dd-GEV IDF parameters.This is a point for further investigation.Additionally, the middle column shows the simulated storms (red rectangles) and cells (blue rectangles) corresponding to the extreme event of the simulated time series.

Discussion and conclusions
The original version of the Bartlett-Lewis rectangular pulse (OBL) model has been optimized for the Berlin-Dahlem precipitation time series.Subsequently, IDF curves have been obtained directly from the original series and from simulation with the OBL model.The basis for the IDF curves has been a parametric model for the duration dependence of the GEV scale parameter which allows a consistent estimation of one single duration-dependent GEV using all duration series simultaneously (dd-GEV IDF curve).Model parameters for the OBL model and the IDF curves have been estimated for all months of the year and seasonality in the parameters is visible.Typical small-scale convective events in summer and large-scale stratiform precipitation patterns in winter are associated with changes in model parameters.
We have shown that the OBL model is able to reproduce empirical statistics used for parameter estimation; mean, variance, and autocovariance of simulated time series are in good agreement with observational values, whereas the probability of zero rainfall is more difficult to capture (see Rodriguez-Iturbe et al., 1987;Onof and Wheater, 1994a).
With respect to the first research question posed in the introduction, we have investigated to what extent the OBL model is able to reproduce the intensity-duration relationship found in observations.We have shown that they do repro- duce the main features of the IDF curves estimated directly from the original time series.However, a tendency to underestimate return levels associated with long return periods has been observed, similar to Onof and Wheater (1993).Including the third moment in parameter estimation has not significantly improved the OBL model's representation of extremes, in contrast to findings for the Neyman-Scott variant (Cowpertwait, 1998).
Furthermore, IDF curves for January show a strong discrepancy between the OBL model simulations and the original series.We have hypothesized and have investigated that this is due to the Berlin-Dahlem precipitation series containing an extreme rainfall event associated with the winter storm Kyrill passing over Berlin on 18-19 January 2007.This event was very rare, in the sense that on short timescales (e.g. 1 and 3 h) such an event has probability to occur only once in a period larger than 1000 years on average.This addresses the second research question: how are IDF curves affected by very rare extreme events which are unlikely to be reproduced with the OBL model for a reasonably long simulation?Having excluded the year 2007 from the analysis, the aforementioned discrepancy in January has disappeared.We conclude that an extreme event which is rare (return period of 23 000 years) with respect to the timescales of simulation (1000 × 13 years) has the potential to influence the dd-GEV IDF curve, as 1 out of 13 values per duration (i.e. one maximum per year out of a 13-year time series) does change the GEV distribution.However, its potential to influence mean and variance statistics used to estimate OBL model parameters is minor.
The third question addresses the validity of the durationdependent parametric model for the GEV scale parameter Shown are the 0.5 (red), 0.9 (green), and 0.99 (blue) quantile from observations at Berlin-Dahlem (solid lines).The shaded areas are the respective 0.05 and 0.95 quantiles for the associated IDF curves obtained from 1000 OBL model simulations.Shown are the 0.5 (red), 0.9 (green), and 0.99 quantile (blue); shaded areas/box plots represent the variability over the 50 repetitions (5-95 %) which allows a consistent estimation of IDF curves.For a set of long simulations (1000 years) with the OBL model, the comparison of IDF curves with the duration-dependent GEV approach with quantiles from a GEV estimated from individual durations suggests a systematic discrepancy associated with the flattening of the IDF curve for short durations.Quantiles from individual durations are smaller for short durations than in the dd-GEV approach IDF curves, which is a challenge for the latter modelling approach.However, instead of altering the duration-dependent formulation of the scale parameter σ d (Eq.4), a different sampling strategy for durations d used in the estimation of the dd-GEV parameters might alleviate the problem.This is a topic for further investigation.
We have not found the OBL model producing unrealistically high precipitation amounts, as discussed for the random-η model (Verhoest et al., 2010).Nevertheless, improvements in reproducing the observed extreme value statistics (especially large return levels) could be made by adding a constraint between intensity and duration parameters in the model equations, as a previous study has shown (Kaczmarska, 2011).
In summary, the OBL model is able to reproduce the general behaviour of extremes across multiple timescales (durations) as represented by IDF curves.Very rare extreme events do not have the potential to change the OBL model parameters but they do affect IDF statistics, and consequently modify the previous conclusion for these cases.A durationdependent GEV is a promising approach to obtain consistent IDF curves; its behaviour at small durations needs further investigation.

Figure 2 .
Figure 2. Example realization of the OBL model.The top panel shows the continuously simulated storms and cells by the model.In the middle panel the cell intensities are combined with a step function.The bottom panel shows the aggregated artificial precipitation time series.Parameters used: λ= 4/120 h −1 , γ = 1/15 h −1 , β = 0.4 h −1 , η = 0.5 h −1 , µ x = 1 mm h −1 .

Figure 4 .
Figure 4. Comparison of statistics derived from the observational record (red dots) and 1000 simulated time series (box plots): (a) variance at 6-hourly aggregation level and (b) probability of zero rainfall at 12-hourly aggregation.

Figure 6 .
Figure 6.Relative differences between observed and simulated return levels obtained with including the third moment (red) and with using the probability of zero rainfall (blue) in parameter estimation for (a) July 0.5 quantile and (b) July 0.99 quantile.Dotted lines show the 0.05 and 0.95 quantile range of 1000 simulations.

Figure 7 .
Figure 7. Visualization of July extremes as observed (RR obs , left column) and simulated by the OBL model (RR OBL , right column).Shown are short segments including the maximum observed/simulated rainfall (red vertical bars) at durations 1 h (a), 6 h (b), and 24 h (c).Additionally, the middle column shows the simulated storms (red rectangles) and cells (blue rectangles) corresponding to the extreme event of the simulated time series.

Figure 8 .
Figure 8. Seasonality of IDF model parameters estimated directly from the Berlin-Dahlem series (blue line), and estimated from 1000 OBL model simulations (red).The red shadings give the range of variability (5 to 95 %) from the 1000 simulations with the median as solid red line.

Figure 10 .
Figure10.dd-GEV IDF curves for (a) January and (b) July and associated quantiles of a GEV distribution estimated for individual durations.Shown are the 0.5 (red), 0.9 (green), and 0.99 quantile (blue); shaded areas/box plots represent the variability over the 50 repetitions (5-95 %) Koutsoyiannis et al. (1998)e IDF T (d).The estimated IDF curve IDF T 1 (d) for return period T 1 is independent of the estimate of another curve IDF T 2 (d) with return period T 2 > T 1 .There is no constraint ensuring IDF T 2 (d) > IDF T 1 (d) for arbitrary durations d.For example, for a given duration d, the 50-year return level can exceed the 100-year return level.Consequently, this approach easily leads to inconsistent (i.e.crossing) IDF curves.To overcome these problems and increase robustness in constructing IDF curves,Koutsoyiannis et al. (1998)suggested a duration-dependent scale parameter σ d :

Table 1 .
Return period for the event Kyrill as estimated from the observational time series with this particular event left out and included for parameter estimation.