Case Study 6.1 - Extreme Rainfall & Flood Frequency

np.random.seed(101)
years = np.arange(1994, 2024)      # 30 years of record
# synthesise annual maxima by drawing from a Gumbel (extreme-value) distribution
annual_max_1day = np.round(stats.gumbel_r.rvs(loc=185, scale=55, size=len(years), random_state=101), 1)

# .fit() estimates the Gumbel parameters that best describe this record
loc_fit, scale_fit = stats.gumbel_r.fit(annual_max_1day)
print(f"Fitted Gumbel: loc={loc_fit:.1f} mm, scale={scale_fit:.1f}")

# Return period T <-> non-exceedance probability (1 - 1/T); .ppf inverts the CDF
print("\nReturn-period 1-day rainfall:")
rp = {}
for T in (2, 5, 10, 25, 50, 100):
    rp[T] = stats.gumbel_r.ppf(1 - 1 / T, loc=loc_fit, scale=scale_fit)
    print(f"  {T:>3}-year: {rp[T]:6.1f} mm")

fig, axes = plt.subplots(1, 2, figsize=(12, 4.5))
axes[0].bar(years, annual_max_1day, color="steelblue")
axes[0].set_title("Annual Maximum 1-Day Rainfall"); axes[0].set_ylabel("mm")
# smooth curve of design depth vs return period; the fitted rp points sit on it
T_axis = np.linspace(1.01, 120, 200)
axes[1].plot(T_axis, stats.gumbel_r.ppf(1 - 1 / T_axis, loc_fit, scale_fit), color="darkred")
axes[1].scatter(list(rp), list(rp.values()), color="black", zorder=5)   # zorder=5: draw on top
axes[1].set_xscale("log"); axes[1].set_xlabel("Return period (years)"); axes[1].set_ylabel("Rainfall depth (mm)")
axes[1].set_title("Return-Period Curve")
plt.tight_layout(); plt.show()