Rare snowmelt estimation (GB)

I read Hough and Hollis’ 1997 paper recently which uses Met Office synoptic stations to estimate a magnitude – recurrence relationship for snowmelt in the UK. i.e. how often do we get how large a snowmelt event? I’m particularly interested in their work as I’m aiming to build on their work by running snowmelt models over long time period (30+ yrs) and using the Snow Survey of Great Britain for calibration. Keep an eye out for this in the future!

Hough and Hollis focused their effort on stations with long hourly data and supplemented these with higher altitude stations with a daily record. The stations they used are shown below.

Table 1 from Hough and Hollis (1997), Rare snowmelt estimation in the United Kingdom.

The purpose of their work was to call into question the melt rate recommended in the Flood Studies Report (FSR) in 1975 which is applied particularly for the estimation of Probable Maximum Flood in dam design. The figure recommended is 42 mm/day, applied country wide. This is generally taken as a 100 year snowmelt event. This value is still recommended today by the 1999 Flood Estimation Handbook. It seems glaringly obvious that the melt potential of Hampshire is much less than the Cairngorms, which is what Hough and Hollis demonstrated. If you’d like to know more I recommend reading their paper.

The focus of this post is that at the end of the paper Hough and Hollis present four regression equations for estimating five year return period 24 hour snowmelt rates. These are below:

Table 7 from Hough and Hollis (1997), Rare snowmelt estimation in the United Kingdom.

In 1997 the parameters in table 7 were not easily available on a grid, and so making estimates from these equations was likely only done for individual Met Office stations. In 2015 these data are available from Ordnance Survey and the Met Office. The Met Office data are available on a 5 km grid with long term averages available from 1961 for temperature and 1971 for wind. I used the 1971-2000 long term average period. The Ordnance Survey does not provide Northern Irish data, so this country is not included – sorry!

The following R code reads the rasters and performs the calculations to estimate melt:

# ---------------------------------------
# Calculate melt grids from Hough & Hollis 1997
# ---------------------------------------
# Taken from table 7 (p136) in
# 'Rare snowmelt estimation in the United Kingdom'

# ---------------------------------------
# Set up
# ---------------------------------------
library(raster)

# ---------------------------------------
# read rasters
# ---------------------------------------
DEM = raster("~/dir/OS/Terrain50")
# max Jan mean daily ave temp
MAX = raster("~/dir/data/MaxTemp_Jan_1971-2000_LTA.txt")
# mean Jan windspeed (knots at 10 m above ground)
Wind = raster("~/dir/data/MeanWindSpeed_Jan_1971-2000_LTA.txt")

# Resample DEM to 5km
DEM = resample(DEM, MAX, "bilinear")
DEM = mask(DEM, MAX)
DEM = mask(DEM, DEM, maskvalue = 0)

MAX = mask(MAX, DEM)
Wind = mask(Wind, DEM)

ex = extent(0, 7e+05, 0, 1250000)
DEM = crop(DEM, ex)
MAX = crop(MAX, ex)
Wind = crop(Wind, ex)

# ---------------------------------------
# Single non-weather factor
# ---------------------------------------
melt.1.nw = 5.09 + .085 * DEM

# ---------------------------------------
# Double non-weather factor
# ---------------------------------------
North = coordinates(DEM)
North = cbind(North, z = round(North[,2]/100))
North = rasterFromXYZ(North)
North = mask(North, DEM)

melt.2.nw = -3.8 + .083 * DEM + .00187 * North

# ---------------------------------------
# Single weather factor
# ---------------------------------------
melt.1.w = 71.16 - 9.45 * MAX

# ---------------------------------------
# Double weather factor
# ---------------------------------------
melt.2.w = 52.52 - 9.08 * MAX + 1.46 * Wind

The results from these are summarised in the table below to two significant figures:

As you can see from table 7 above, the fittings for these models are quite good. However, when you look at the results the obvious problem is the negative predicted melt. This is discussed in the paper with the solution being to set negative melt to zero. Possibly the regression solution should have been forced through zero.

Below you can see the output from the four regression models, and below the code used to generate this plot.

5 yr return period melt, as estimated by Hough (1997) table 7. Contains Ordnance Survey and Met Office data.

# ---------------------------------------
# Set up
# ---------------------------------------
library(raster)
library(RColorBrewer)
library(rasterVis)
library(lattice)
library(grid)

png("results/Hough_melt.png", width=550, height=750,
    res=100)
x = stack(melt.1.nw, melt.2.nw, melt.1.w, melt.2.w)
myTheme=rasterTheme(region=brewer.pal('Blues', n=9))

levelplot(x, names.attr=c("Single non weather", "Double non weather",
    "Single weather", "Double weather"),
    par.settings=myTheme, layout=c(2,2))
trellis.focus("legend", side="right", clipp.off=T, highlight=F)
grid.text("Melt\n(mm/day)", 0.25, 0, hjust=0.5, vjust=1.2)
trellis.unfocus()
dev.off()