A case study implementing the SHyFT modelling framework for the American

(1)

regional calibration methods for a distributed hydrologic

modelling framework

A case study implementing the SHyFT modelling framework for the American

River Basin in Sierra Nevada, California

Marco Westergren

Thesis submitted for the degree of Master in Geosciences

(Hydrology) 60 credits

Department of Geosciences

Faculty of mathematics and natural sciences

(2)

(3)

(4)

regional calibration methods for a distributed hydrologic

modelling framework

A case study implementing the SHyFT modelling framework for the American River Basin in Sierra Nevada, California

Marco Westergren

(5)

hydrologic modelling framework http://www.duo.uio.no/

Printed: Reprosentralen, University of Oslo

(6)

(7)

UNIVERSITY OF OSLO

Abstract

Faculty of Mathematics and Natural Sciences Department of Geosciences

Seasonally snow-covered regions are an important source of water and a crucial component of the yearly water balance for substantial areas of the world. The American River, located in the California Sierra Nevada presents an opportunity for evaluation of a distributed hydrologic framework, spanning a great range of elevation, land covers, climate, and degree of human influence. This project also details the first application of the SHyFT framework in North America.

This study aims to collect and process available, relevant hydrometeoro- logic data for the region, and using currently available methods in SHyFT, assess the performance differences between regional approaches versus local and lumped methods, as well as the effect different snow-routines have on model results.

The experimental design, employing five different region definitions for each of the respective snow-routines, found highly variable Nash- Sutcliffe efficiencies, ranging between -3.249 and 0.611. The poorest results occurring in problematic calibration approaches for local regions, and regulated sub-regions, the Middle Fork in particular. The overall best performing region was the entire Upper American River Watershed, gauged at Folsom lake using af ull−natural−f low adapted time series. Between the three headwater regions, the free-flowing North Fork had the best overall performance between experimental setups, maximum NSE found during validation of 0.523. Between the two snow-routines included, the more advanced gamma-snow outperformed the HBV-snow adapted routine. Results show limited applicability for transfer of locally calibrated parameters to larger regions, for all but one method.

The issues in model performance are attributed to several possible sources, spatial interpolation methods of point-based measurements, data quality issues, terrain capture difficulties, and the presence of regulated basins.

(8)

(9)

Acknowledgements

I’d like to thank my supervisor John F. Burkhart for enabling this project, assistance, and local knowledge of the region. Although it does deserve mention that the final result does deviate from the originally planned project, in more ways than one.

Secondly I’d like to extend thanks to Felix Matt at the Department of Geosciences, for providing vital technical assistance, without it, this project would never have come to fruition.

While the SHyFT framework provides almost unprecedented adaptabil- ity, it comes at a slight cost of accessibility for users with limited knowledge of programming, and a somewhat tall threshold to overcome. While the shyft google group and shyft-doc github repository provide vital support, certain elements in preparation, calibration and running the model has re- lied heavily on immediate access to developers, and personal communica- tion.

In this project I have tried to outline the most relevant structures, parameters and methods, while also including configuration files and python scripts that addresses practical factors in testing and running SHyFT. Perhaps it might serve a secondary purpose in aiding accessibility for future users in situations similar to my own.

The typesetting of this thesis has been inspired by styles provided by Vel¹ through L^ATEXTemplates (http://www.latextemplates.com/) originally developed by Steve Gunn and Sunil Patel.

Finally I would like to thank my permanent roommate, Emilie, for having patience, understanding, and for support along the way.

Marco Westergren

1http://www.vel.nz/

(10)

(11)

List of Figures

1.1 Mapping of drought severity in California . . . 6

2.1 Global Water Cycle . . . 8

2.2 Schematic diagram of a watershed . . . 8

2.3 Flowchart outlining the processes in the SNOW-17 model. From (Anderson, 2006) . . . 12

2.4 General structure of the SAC-SMA model. From NOAA (2016). . . . 14

2.5 Principle of the snow depletion curve. From (Kolberg, Rue, and Gottschalk, 2006) . . . 22

2.6 Illustration of catchment sensitivity function g(Q) and storage- discharge relationshipf(S)as shown in (Kirchner, 2009) . . . 26

3.1 Location of The American River and subregion . . . 32

3.2 Spatial distribution of mean annual precipitation for the Upper American River with a 4km resolution grid from the PRISM dataset (Shamir and Georgakakos, 2006). . . 36

3.3 Mapped locations of temperature (T), precipitation (P), and discharge stations (Q), as well as relative humidity (HUM) which also works as proxy locations for radiation, produced in QGIS using data described herein. . . 38

3.4 Digital elevation model from the "National Elevation Dataset (NED)" supplied by USGS. . . 44

3.5 Land cover types in the research area, from the National Land Cover Database (NLCD 2011) . . . 45

4.1 Catchment ID definitions created for the region model definition covering all catchments with discharge measurement . . . 52

4.2 Spatially distributed precipitation input to SHyFT . . . 54

4.3 Generalized Semivariogram . . . 56

4.4 Spatially distributed temperature input to SHyFT . . . 57

4.5 Date ranges for available aata . . . 59

4.6 Brief description of method for missing data handling. Excerpt from data-parsing script in python. . . 60

5.1 Discharge time series removed from calibration and validation . . . . 69

(16)

5.2 Boxplot of calibrated parameters between the five model set-ups using the gamma-snow module. . . 74 5.3 Boxplot of calibrated parameters between the five model set-ups

using the HBV-snow module. . . 74 5.4 Hydrograph comparing the simulated and observed streamflow for

the entire Upper American Basin using the GS_F calibration. . . 80 5.5 Hydrograph comparing the simulated and observed streamflow for

Dry Creek in the Lower American Basin using the GS_F calibration. . 81 5.6 Hydrograph comparing the simulated and observed streamflow for

the North Fork American River using the HS_NF calibration. . . 82 5.7 Hydrograph comparing the simulated and observed streamflow for

the North Fork American River using the GS_MF calibration. . . 82 5.8 Hydrograph comparing the simulated and observed streamflow for

the Pilot Creek subcatchment using the HS_SF calibration . . . 83 6.1 Scatter plot of observed(x) versus simulated(y) discharge for the

Upper American basin. . . 88 6.2 Comparison observed and simulated yearly total volume at the end

of September for the North Fork and whole Upper American, from both GS and HS methods. . . 89 6.3 Comparison of MODIS imagery and results from Shamir and

Georgakakos (2006) . . . 92 6.4 SHyFT-modelled snow covered area (SCA) from the same dates as

shown in figure 6.3 . . . 93 6.5 Map of dams within the modelled region. From StreamNet Dams,

(StreamNet, 2016) . . . 96

(17)

List of Tables

2.1 Parameters and descriptions for the SNOW-17 model . . . 13

2.2 Parameter names, typical values and units for the SAC-SMA model. From Koren, Smith, and Duan (2013) . . . 15

2.3 Parameter names and descriptions used in gammasnow. From the SHyFT documentation. . . 19

2.4 Parameters with descriptions and units for HBV-snow routine. . . 23

3.1 Climate normals from Folsom Dam, elevation107m. Based on data collected from 10/26/1955 to 04/30/1993 (WRCC, 2016) . . . 36

3.2 Climate normals from Blue Canyon, elevation1610m. Based on data collected from 03/01/1940 to 06/09/2016 (WRCC, 2016) . . . 37

3.3 Precipitation stations from CDEC . . . 39

3.4 Temperature stations from CDEC and NOAA . . . 40

3.5 Relative humidity measuring stations from CDEC . . . 41

3.6 Wind speed measuring stations from CDEC . . . 41

3.7 Incoming radiation measuring stations from CDEC . . . 42

3.8 Discharge measurement stations provided by the CDEC and USGS . 43 4.1 Experimental design framework . . . 57

5.1 Descriptive statistics and trend for precipitation and temperature data 64 5.2 Decriptive statistics and trend for relative humidity and wind speed data . . . 66

5.3 Descriptive statistics and trend for radiation data . . . 68

5.4 Descriptive statistics and trend for discharge data . . . 69

5.5 Invariant parameters in Gamma-snow based model routines . . . 72

5.6 Range and values for calibrated parameters in the gamma-snow routines. . . 73

5.7 Range and values for calibrated parameters in the HBV-snow routines. 73 5.8 Nash-Sutcliffe Efficiency values from calibration of the various regional definitions . . . 73

5.9 Goodness of fit assessment for full-region and lumped gamma-snow models and covered catchments . . . 76

5.10 Goodness of fit assessment for headwater region models and covered catchments . . . 77

(18)

5.11 Goodness of fit assessment for full-region and lumped HBV-snow models and covered catchments . . . 78 5.12 Goodness of fit assessment for headwater region HBV-snow models

and covered catchments . . . 79 5.13 Goodness of fit assessment for subcatchments with limited data

availability. Validation period 2000-2005 . . . 81 6.1 Calibration differences between NSE and KGE calibration . . . 87 6.2 Validation using KGE calibration function for the North Fork American 87

—————————————————————————————-

(19)

Chapter 1 Introduction

1.1 Introduction

Water resources are a necessity for society, and knowledge of the hydrologic cycle, and physical behavior of water flows are needed to be able to provide sufficient quality and quantity of potable water, as well as agricultural and industrial use.

Climate change, population growth and land use affect the relative water demand (RWD) as defined by Vörösmarty et al., 2000 and shows the current level of water stress, with a contemporary level defined as 1985 and predictions toward 2025 found that more than 2 billion will live in areas of high water stress by 2025, an 85% increase from 1985 levels (Vörösmarty et al., 2000).

Although the majority of the world’s population live within 60 km of a coast¹, seasonally snow covered and mountainous regions are an important control factor for water availability, such the western mountains of North America, and vitally the Hindu Kush-Karakorum-Himalaya (HKH) region, as 50-60% of the world’s population live in areas that to various extents depend on water from this region (Barnett, Adam, and Lettenmaier, 2005). There is great interest in future development of water resources with a changing climate, found from trends in snow cover in the HKH region (Immerzeel et al., 2009), as well as coupled land surface hydrology and atmospheric models (Chen and Dudhia, 2001; Fowler, Blenkinsop, and Tebaldi, 2007). Furthermore, impact studies for particular regions associated with climate change are abundant, exemplified for HKH (Akhtar, Ahmad, and Booij, 2008), and smaller regions such as Sweden (Bergstrom et al., 2001), or the Sierra Nevada (Dettinger et al., 2004).

Global fluxes between stocks of water (oceans, rivers and lakes, atmosphere, soil and ground water) is driven by energy supplied by

1UN Atlas of the Oceans: http://www.oceansatlas.org/

(20)

the sun, and dictates the trends and patterns in global distribution of precipitation. For regional distribution of precipitation, local factors such as lakes, proximity to the ocean and orographic effects control the distribution lending to the difficulties in accurately measuring precipitation at a site to represent larger areas (Dingman, 2008).

Adding to the complexity of accurately measuring actual precipitation, is assessing the degree of loss in measurement. Experiments using four pairs of different manufacture rain gauges in controlled, as well as outdoor testing by the USGS found up to 10.64% absolute difference between gauges (Gordon, 2003). Experiments performed by Dyck and Gray using a dense precipitation network in the Canadian prairies between 1968-1977 found a coefficient of variation between gauges ranging from 6.52-100% within a topographically very homogeneous area (Dyck and Donald Maurice Gray, 1979). The difficulties of acquiring correct areal precipitation for a region is further illustrated by Dingman, 2008, p. 13-17 and Bales et al., 2006.

For several regions the seasonal snow cover is a major important control for the yearly water balance. Land masses above 40^◦ N, 42% of the global land surface area, has a seasonal snow cover of significant duration, and snow cover is present in mountains of sufficient altitude for all regions (Dingman, 2008, p. 179-180). The different nature of solid precipitation, compared to liquid precipitation, adds further complexity in modelling the quantity and distribution of water. Snow can store significant volumes over time, and behaves differently influenced by vegetation, both in terms of energy fluxes, transportation of mass and direct interception (Gelfan, JW Pomeroy, and Kuchment, 2004; Hedstrom and JW Pomeroy, 1998). These issues of under estimating measurements as noted by (Gordon, 2003), is even greater for precipitation as snow and the process of redistribution of snow is further addressed by e.g. Essery, L. Li, and John Pomeroy (1999) and Liston and Elder (2006).

Research in finding a surface water supply index (SWSI) by Shafer and Dezman found that snow can account for 65-85% of annual runoff (Shafer and Dezman, 1982), and its importance is further expressed by several reports (Akhtar, Ahmad, and Booij, 2008; Bales et al., 2006; Barnett, Adam, and Lettenmaier, 2005; Diffenbaugh, Scherer, and Ashfaq, 2013; Rees and Collins, 2006). Adding to the complexities mentioned, the topography of mountain regions add to the difficulty of modelling the distribution of precipitation, snow in particular . Additionally due to accessibility, sampling density is often sparse. For these reasons implementation and

(21)

use of remote sensing tools is becoming more common, mentioned by Bales et al. (2006) and exemplified for the HKH region by Immerzeel et al. (2009).

Historical data are needed in order to properly calibrate conceptual models, and to assess and validate the model before the accuracy of the results can be verified. Different approaches to regional modelling addressed herein include:

• Regionalization by means of physical catchment characteristics

• Regionalization of model parameters using lumped parameter hydrological models

• Regionalization of model parameters using rainfall-runoff models Utilizing physical catchment characteristics can be found in many historical publications (Don M Gray, 1961; D. Gray, 1964; Horton, 1932;

Strahler, 1957). This is to some extent owed to technological limitations at the time, but enabled estimation of hydrological properties from available parameters such as catchment area, length, slope of the stream or stream ordering. Even today, similar principles are utilized for instance by the Norwegian Water Resources and Energy Directorate (NVE). Their method for estimating flood magnitude, for given return periods, use a combination of physical characteristics and historical calibration for different regions.

These do however come with some limitations in terms of catchment size, and caution is recommended for catchments smaller than 100km² (Midttømme et al., 2011).

The use of lumped parameter models, such as the HBV (Lindstrom et al., 1997) have been utilized at various scales and in different parts of the world, however results vary. In a study of 56 catchments in the UK by Deckers et al. (2010) found the majority of modelled catchments showed Nash-Sutcliffe efficiency below 0.75, and by their conclusion, did not favor rationalization.

A similar approach in Austria by Merz and Bloschl (2004), with a larger sample size of 306 catchments found similar results, and further indicated spatial proximity as a better approach.

A persistent issue, highlighted by Kirchner (2006) is the issue of over- parametrization, it mentions Getting the right answer for the right reasons, advocating the use of fewer parameter rainfall-runoff models. In Kirchner (2009) a method for calculating streamflow as a generalized function of storage, NSE exceeding 0.9 were obtained, albeit with a limited sample size.

Detailed description of the method is described in 2.3.4. Other methods

(22)

exemplified by Post and Jakeman (1999) and Kokkonen et al. (2003) in Australia and USA, respectively, had trouble reproducing streamflow, even without accounting for snow.

Another rainfall-runoff methodology is shown with the EXP-HYDRO model by Patil and Stieglitz (2014) for 756 catchments in the USA found Nash-Sutcliffe efficiency exceeding 0.6 for the majority of modelled catchments.

1.2 Motivation

The snowpack in the Sierra Nevada has traditionally been an important factor in Californias water balance, as a control and natural water reservoir to balance water supply across the seasons (Chou, 2014). The 1956-2000 average water storage in the Sierra Nevada snow pack amounts to 15 million acre feet, or 18500 million cubic meters of water per year, reported by the California Department of Water Resources² in 2007. Between the State Water Project (SWP) and the Central Valley Project (CVP) the snow in the Sierra Nevada provides water for an estimated 25 million people and 3.6 million acres of farmland (Chou, 2014).

Correctly modelling the behaviour of this vital water resource is of great interest, particularly in the current drought conditions illustrated in 1.1, and is the subject of great interest, also in the media (Fountain, 2016; Stevens, 2016). As stated by the American River Hydrologic Observatory:

With better management, California’s existing water supply could go further to meeting the needs of the state’s urban and agricultural uses. (ARHO, 2014)

An initial motivation for using this particular study region, was the availability of data from regional wireless sensor networks (WSN) as a means to improve regional model performance. However, for reasons beyond our control, it became apparent that availability and quality of the data was not sufficient for our requirements. Therefore, we discarded the intention of evaluating regionalizatoin of WSN data for hydrologic modeling. However, with the release of Shyft, the decision to keep focus on the American River region, despite the rejection of original project intentions, was largely motivated by the interesting nature of the

2http://www.water.ca.gov/climatechange/docs/062807factsheet.pdf

(23)

region and catchment. Particularly in the context of recent drought and overall water availability of the region, it was felt the Sierra Nevada presented a unique opportunity to assess utility of the SHyFT framework for a different topographical and climatological region than for which it was developed. Previous research, using the similar ENKI model, had primary implementations in Norway (Hailegeorgis and Alfredsen, 2016;

Hailegeorgis, Alfredsen, et al., 2016) and Himalayan India (Trine J Hegdahl et al., 2016).

Focus has since been shifted to the aspect of regional modelling, within a distributed framework, and assessment of quality differences compared to local methods.

(24)

FIGURE 1.1: Current drought conditions as of May 2016.

From: https://www.drought.gov/drought/dews/california

1.3 Objectives

The scope of this project will be to evaluate the functionality and performance of the Statkraft Hydrological Forecasting Toolbox for the American River Basin within the available methods and routines.

• Collect and sort necessary data from available sources

• Quality control and analysis

• Conform data to meet existing SHyFT inp repository requirements

• Establish a full regional model for the American River Basin

• Assess the performance on full regional scale, versus lumped and local calibrations using established methods and objective functions

(25)

Chapter 2 Theory

2.1 Hydrology

Hydrology is the geoscience that studies the the processes governing the water cycle. An outline of the water cycle is illustrated in figure 2.1. The totality of the water cycle is an advanced system incorporating water in all phases, and all aspects from atmospheric moisture, liquid water in streams, lakes, oceans and groundwater, as well as solid in snow and ice. For subjects covered herein, the focus will be on the land phase of the hydrologic cycle defined by the U.S Committee on Opportunities in the Hydrological Sciences in Dingman (2008):

The land phase of the hydrologic cycle:

the movement of water substance on and under the earth’s land surfaces, the physical and chemical interactions with earth materials accompanying that movement, and the biological processes that conduct or affect that movement.

2.1.1 Basic Concepts and Water Balance

The basic hydrological unit is the watershed, also referred to as drainage basin, river basin or catchment. A watershed is defined as

...the area that appears on the basis of topography to contribute all the water that passes through a given cross section of a stream (Dingman, 2008).

From figure 2.2, the basic variables for a watershed is illustrated, G_in and Gout denote groundwater flow in and out, P precipitation in, ET

(26)

FIGURE2.1: Summary of the global water cycle. From USGS (http://water.usgs.gov/edu/graphics/watercyclesummary.jpg)

FIGURE 2.2: Schematic diagram of a watershed, from Ding- man (2008).

(27)

evapotranspiration out and Qthe stream outflow. From this it is possible to define the water balance equation 2.1 for any given period∆t

P +Gin−(Q+ET +Gout) = ∆S (2.1)

2.1.2 Hydrologic Modeling

Hydrologic models are a representation of one, or several aspects of the water cycle. A hydrologic model can cover single aspects such as snowmelt or evapotranspiration, or the multiple aspects producing runoff. In Dingman (2008) the following aspects of hydrological models are addressed:

Simulation Basis

Physically based models derives its equations from the physical world such as the conservation of mass, exemplified in equation 2.1.

Conceptual models uses methods not directly linked to physical characteristics, such as linear reservoir model employed by the HBV model (Bergström and Singh, 1995; Lindstrom et al., 1997) in the generalized form q = k·S where q is the streamflow, k is a positive constant and S denotes storage.

Empirical modelsderive their functions from observation and the use of mathematical regression, e.g. temperature index snowmelt.

Stochastic models depends on time-series analysis and periodic components, e.g. autoregressive models (Xu, 2002, p.2-13).

Spatial Representation

Lumped models treats a watershed as a single point, and input variables and parameters are adapted to single values representative for the entire watershed.

Distributed modelsallow for some differentiation within the model, such as subdividing the catchment into elevation bands or by means of surface characteristics such as the distributed HBV model (Lindstrom et al., 1997).

Coordinate system representation links the model to a fixed coordinate system, such as employed with the SHyFT framework.

(28)

Temporal Representation

Steady state models output single values representing the mean across a longer period. Steady state seasonal offers some more variability by sub- division into seasons. Examples include the evapotranspiration routine in the HBV model, using long term mean values of potential evapotranspiration based on the Penman formula (Penman, 1948; Lindstrom et al., 1997).

Single event produces results corresponding to a simulated input, examples of this method include unit hydrographs (NWS, 2005).

Continuous modelsproduces continuous output from a continuous series of inputs across a defined time period and time interval. Interval ranges commonly used are from hourly to monthly.

Solution Methods

0-dimensionalsolutions are not based on a formal coordinate system, such as implemented in lumped models.

Formal-Analytical solutions are differential equations solvable analytically, for instance the Kirchner method described in 2.3.4.

Formal-Numerical solutions are differential equations not solvable by analytical methods, but need to be solved numerically by finite-difference or finite-element methods. Examples of this approach are for instance ground water models.

Hybrid solutions uses different methods for different aspects contained within the same model framework.

2.2 Hydrologic Models

The main purpose for hydrologic models are predictions and forecasts.

Predictions in this case refer to magnitude of some hydrologic quantity, such as design flood of various return periods, floodplain zones used in land- use planning, while forecasts denote a specific expected event, such as the streamflow based on predicted upcoming precipitation for flood warning (Dingman, 2008).

This section will further explain common forecast modelling approaches focused on model approaches in use in the study region, and furthermore present the theoretical basis of the modelling framework used in this study.

(29)

2.2.1 CNRFC

The California Nevada River Forecast Center (CNFRC)¹ is a field office for the National Weather Service (NWS) based in Sacramento, California. The NWS is an agency of the National Oceanic and Atmospheric Administration (NOAA). The operational responsibility of CNRFC is to provide river and flood forecasts on short and long term, as well as technical support in related matters. The short-term aspect include flash flood warnings, and daily scale flood- and reservoir inflow forecasts, as well as longer range Hydrologic Ensemble Forecast Service(HEFS).

HEFS use ensembles of hydrometeorological forecasts into a hydrologic processor containing modules including:

• Snow-17 - Snow ablation model

• SAC-SMA - SACramento Soil Moisture Accounting model

• Unit hydrograph flow routing

• Reservoir models Demargne et al. (2014)

2.2.2 SNOW-17

The SNOW-17 model is an air temperature index snow ablation model first described in 1973 and has since only seen minor changes, and an addition of snow-depth computation. The model is conceptual, and includes most factors in simplified forms (Anderson, 2006). The main processes of the model for a given column of snow are:

• Form of precipitation

• Accumulation of snow cover

• Energy exchange at snow-air interface

• Internal state of the snow cover

• Transmission of water through the snow cover

• Heat transfer at the soil-snow interface.

(30)

FIGURE2.3: Flowchart outlining the processes in the SNOW- 17 model. From (Anderson, 2006)

The general structure is shown in figure 2.3.

SNOW-17 is comprised of 12 parameters, out of which 6 are classified as major, and are calibrated, and 6 are classified as minor. The parameters are described in table 2.1

The complete technical documentation can be found from Anderson (2006).

2.2.3 SAC-SMA

The Sacramento Soil Moisture Accounting (SAC-SMA) model is a conceptual model that is used to account for the distribution of soil moisture characteristics responsible for streamflow generation, modeled as linear reser- voirs. An overview of the model structure is shown in figure 2.4 (NOAA, 2016).

The SAC-SMA generates runoff as six components per time interval:

• Permanently impervious runoff

• Temporary impervious runoff

1http://www.cnrfc.noaa.gov/about/

(31)

TABLE 2.1: Parameters and descriptions for the SNOW-17 model

Parameter Description Unit

MAJOR

SCF Precipitation correction factor -

MFMAX Maximum melt factor during non-

rain

mm·^◦C·6hr⁻¹

MFMIN Minimum melt factor during non-

rain

mm·^◦C·6hr⁻¹

UADJ Average wind function mm·mb⁻¹

SI Mean area water equivalent above

which there is always 100 % snow cover

mm

Areal Depletion Curve

Curve defining areal extent - MINOR

NMF Maximum negative melt factor mm·^◦C·6hr⁻¹

TIPM Antecedent temperature index pa-

rameter

Range0.01−1

PXTEMP Rain/snow treshold temperature ^◦C

MBASE Base temperature for melt during

non-rain

◦C

PLWHC Liquid water capacity %

DAYGM Constant daily melt at base mm·^◦C·6hr⁻¹

• Direct runoff

• Interflow

• Supplemental baseflow

• Primary baseflow

and further stores water in six soil moisture storages per time interval:

• Upper zone tension water content (UZTWC)

• Upper zone free water contents (UZFWC)

• Lower zone tension water content (LZTWC)

• Lower zone free water contens (LZFWC)

• Lower zone free supplemental contents (LZFSC)

• Lower zone free primary contents (LWFPC)

(32)

FIGURE2.4: General structure of the SAC-SMA model. From NOAA (2016).

In order to run the model, time series for precipitation, channel inflow, potential ET, areal extent of snow are used as input, and for each time step all runoff components and soil moisture storages as described are produced (NOAA, 2016). The SAC-SMA model uses 16 parameters in total, which are described in Koren, Smith, and Duan (2013) and shown in table 2.2.

2.2.4 Unit hydrograph

The final soil moisture model output is then fed to a unit hydrograph model described in CNRFC (2016).

The unit hydrograph is a transfer function for hydrologic response relating characteristic response of a given watershed to a given effective water input (Dingman, 2008, p 450).

(33)

TABLE2.2: Parameter names, typical values and units for the SAC-SMA model. From Koren, Smith, and Duan (2013)

Parameter Description Range Unit

UZTWM Upper zone tension water capacity 10-300 mm

UZFWM Upper zone free water capacity 5-150 mm

UZK Interflow rate from upper zone 0.10-0.75 day⁻¹ ZPERC Ratio of max/min percolation rate 5-350 - REXP Shape parameter of percolation

curve

1-5 -

LZTWM Lower zone tension water capacity 10-500 mm LZFSM Lower zone supplemental free wa-

ter capacity

5-400 mm

LZFPM Lower zone primary free water capacity

10-1000 mm

LZSK Depletion rate lower zone supplemental storage

0.01-0.35 day⁻¹ LZPK Depletion rate lower zone primary

storage

0.001-0.05 day⁻¹ PFREE Percolation fraction directly to

lower zone free water storage

0.0-0.8 -

PCTIM Permanent impervious area fraction

- ADIMP Max area of additional impervious

area due to saturation

-

RIVA Riparian vegetation area fraction -

SIDE Ratio of deep percolation from lower zone

- RSERV Fraction of lower zone not transfer-

able to lower zone tension storage

-

(34)

2.3 SHyFT Stacks and Modules

The SHyFT²model framework, Statkraft Hydrological Forecasting Toolbox, is an open source framework for hydrologic modelling made available through github. The core methods, in C++, are made available through a Python wrapper. SHyFT is under development at Statkraft, and the main contributors for the C++ core were Sigbjørn Helset and Ola Skavhaug.

Modules Used in SHyFT

For a complete water balance model, the framework employs various methods for the different aspects of the surface water balance. The Shyft model operates using different stacks, that in conjunction provide a full coverage of the relevant hydrologic water balance. The separate modules presented here cover evapotranspiration, snow-accumulation and ablation, and flow generation. A theoretical basis for the shyft modules implemented in this study is described in the subsequent sections, while the model setup is covered in chapter 4.

2.3.1 Evapotranspiration: Priestley-Taylor

The Priestley-Taylor evapotranspiration (Priestley, 1972) is a simplified form of the Penman-Monteith equation. The Priestley-Taylor approach uses net radiation and ground heat flux, slope of the saturation vapor-pressure curve∆, the psychrometric constantγand a dimensionless correction factor α. ∆is defined in Dingman (2008):

∆≡ de^∗_s

dT = 2508.3

(T + 237.3)² ·exp

17.3·T T + 237.3

(2.2) which is the derivative of thee^∗_s approximation given by:

e^∗_s = 0.611·exp

17.3·T T + 237.3

(2.3) Here∆is given in unitskP aK⁻¹and temperature is given in^◦C.

The psychrometric constant is defined as:

γ ≡ c_a·P

0.622·λ_v (2.4)

2https://github.com/statkraft/shyft

(35)

• c_ais the heat capacity of air. UnitM J kg−1K⁻¹

• λ_v is the latent heat of vaporisation. Unit: M J K⁻¹

• P is pressure. UnitkP a

These values are not strictly constant, but Dingman (2008) gives γ = 0.066kP aK⁻¹

as a typical approximate value using heat capacity values for typical air temperatures and one atmosphere pressure.

The generalized form of the full Penman-Montheith equation is shown in equation 2.5 as described in (Sumner and Jacobs, 2005).

LE = ∆(Rn+G) +ρ_ac_p^(e^s_r^−e)

a

∆ +γ(1 + ^r_r^s

a) (2.5)

• LE Latent heat flux. UnitM J m⁻²day⁻¹

• ∆is the slope of the saturation vapor-pressure curve. Unit: kP a^◦C⁻¹

• R_nis net radiation on the surfaceM J m⁻²day⁻¹

• Gis soil heat flux at the surfaceM J m⁻²day⁻¹

• ρ_amean air density. Unitkgm⁻³

• c_pspecific heat of air. UnitM J kg−1^◦C⁻¹

• e_s−eis the vapor-pressure deficit,e_sis saturation vapor-pressure,eis actual vapor-pressure of the air. Unit: kP a

• γ is the psychrometric constant. Unit:kP a^◦C⁻¹

• r_sis bulk surface resistance. Unit: sm⁻¹

• r_ais aerodynamic resistance, described in equation 2.6. Unit: sm⁻¹ Furthermore the aerodynamic resistance is estimated by:

r_a = ln[^z−d_z

0 ]ln[^z−d_z

0v]

k²u (2.6)

• z is height of wind-speed measurement

(36)

• dis displacement height

• z₀ is roughness height for momentum

• z0v is roughness height for water vapor

• kis von Karman’s constant≡0.4

• uis horizontal wind speed at sensor heightz

Sumner and Jacobs (2005) further presents several approximations for d ranging from0.1zveg from Dingman 1994, to0.67zveg from Brutsaert 1982 wherez_vegis the height of vegetation and0.1z0as an approximation forz_0v. The final form of the Priestley-Taylor approximation, using α as a scaling term to compensate for the aerodynamic term from the full Penman- Monteith approach is:

LE =α ∆

∆ +γ (Rn−G) (2.7)

(Sumner and Jacobs, 2005)

In the SHyft framework, the model uses radiation measurements and an albedo parameter for Priestley-Taylor calculation as well asαapproximated as1.26.

2.3.2 Snow: GammaSnow

In the SHyFT gamma-snow³routine, an energy balance approach for snow ablation, and snow depletion curve with gamma distribution are combined in a single routine. The parameters used are shown in table 2.3 and the methods are described below. The module has 16 input parameters, and produces three output responses: snow-covered area, storage and runoff.

Snow ablation

The energy balance approach is based on methods by DeWalle and Rango (2008) described by Trine J Hegdahl et al. (2016). Represented in terms of energy,W ·m⁻²:

∆E =S·(1−α) +L_in+L_out+H_SE+H_L+E_G (2.8)

3https://github.com/statkraft/shyft/blob/master/core/gamma_snow.h

(37)

TABLE 2.3: Parameter names and descriptions used in gammasnow. From the SHyFT documentation.

Parametername Description Unit

winter_end_day_of_year Last day of accumulation day initial_bare_ground_fraction Initial bare ground fraction at

melt onset

-

snow_cv Spatial coefficient variation of

new snow

-

tx Liquid/solid precipitation

threshold

◦C

wind_scale Slope of turbulent wind func-

tion

ms⁻¹

wind_const Intercept in turbulent wind

function

-

max_water Maximum liquid water con-

tent

surface_magnitude Surface energy scaling factor

max_albedo Snow albedo, maximum

value

-

min_albedo Snow albedo, minimum

value

- fast_albedo_decay_rate Albedo decay rate during

melt

day slow_albedo_decay_rate Albedo decay rate in sub-

zero conditions

day snowfall_reset_depth Depth of new snowfall at

which snow albedo is reset

mm

glacier_albedo Albedo used for areas de-

fined as glacier

- calculate_iso_pot_energy Calculate potential energy

flux

bool snow_cv_forest_factor Modification factor for

snow_cv as a function of forest fraction in cell

-

snow_cv_altitude_factor Modification factor for snow_cv as a function of altitude

-

(38)

∆E is the net energy flux at the surface available for snow ablation, S is the incoming shortwave radiation moderated for albedo α, L_inandL_out represent the incoming and outgoing longwave radiation, HSE and HL

represent sensible and latent heat fluxes and EG is the subsurface energy flux.

For the shortwave energy input, the albedo, α, of the snowpack is needed to assess the energy toward melting. The model hasα_min and α_max as given parameters, and the snow albedo, α_t, is calculated dependent on decay rate, temperature and snowfall for a given day

α_t=







α_min+ (αt−1−α_min·

1 2^{F DR}⁻¹

ifT_a >0^◦C αt−1−(αmax−α_min)· 2·SDR¹ ifT_a ≤0^◦C

(2.9) F DR and SDR denote fast and slow albedo decay rates, respectively.

The rate is interpreted as the time, measured in days, at which the snow surface albedo is reduced by 95% of the difference between maximum and minimum values. If a new snowfall exceeding the parameter snowf all_reset_depthoccurs, the snow albedo is reset to maximum. Glacier albedo, where applicable, is constant.

The longwave radiation is based on surface and air temperatures. The incoming longwave radiation from the atmosphere is calculated by Stefan Boltzmann’s law using the air temperature, while the outgoing longwave radiation is calculated using the snow surface temperature,T_ss, given by:

Tss=







1.16·T_a−2.09 ifT_ss<0^◦C

0 otherwise (2.10)

The sensible and latent heat fluxes are calculated by first establishing the vapour pressure using Bosen (1960)

e_s,water ≈33.864·[(7.38×10⁻³)⁸−1.9×10⁻⁵(1.8T + 48) + 1.316×10⁻³]×RH (2.11) WhereRH denotes the measured relative humidity.

For freezing conditions this is scaled by:

e_s,ice =e_s,water×(1 + 9.72×10⁻³T + 4.2×10⁻⁵T²) (2.12) The energy contribution is then modelled as:

(39)

0.98·σ· e_s

Tk

0.0687

·T_k⁴ (2.13)

whereσ = 5.67×10⁻⁸is the Stefan-Boltzmann constant,e_sis the vapour pressure andT_k is the air temperature in Kelvin.

The turbulent component is modelled by turb = wind_scale · U + wind_const whereU is the wind speed, and the other factors are the wind parameters. The energy contribution is modelled as:







turb×(T + 1.7·(es−6.12))−0.98σ273.15⁴ ifT_ss = 0 turb×(T −T_ss+ 1.7(es−6.132·e^(0.103T^−0.186))−0.98σ(Tss+ 273.15)⁴ ifT_ss 6= 0

(2.14) The case of precipitation event on an existing snowpack, the energy contribution is calculated depending on the type of precipitation, defined bytx, the rain/snow temperature threshold parameter.







rain·T_a·c_water· dt¹ ifT_a>0

snow·T_a·c_ice· _dt¹ ifT_a≤0 (2.15) Whererainandsnoware precipitation types,T_ais air temperature, and c_water and c_ice denotes the specific heat capacity, 4180J kg⁻¹ ^◦C⁻¹ for water, 2108J kg⁻¹ ^◦C⁻¹if snow.

The final energy component,E_Gis calculated by

surf ace_heat=surf ace_magnitude·c_ice·T_ss· 1

2 (2.16)

where surface_magnitude is a calibrated parameter and surface_heat is the snow surface cold contentJ ·m⁻².

All energy factors are summed as they’re calculated to the variable ef f ectand applied for each time step, and for each cell (Helset et al., 2016).

Snow distribution

The distribution of snow in SHyFT uses a snow depletion curve (SDC) approach described by Kolberg, Rue, and Gottschalk (2006). The SDC describes how the snow covered area is depleted over the ablation period, with various approaches to the shape of the curve. A general image of the Gamma approach is shown in figure 2.5.

(40)

FIGURE 2.5: Principle of the snow depletion curve. From (Kolberg, Rue, and Gottschalk, 2006)

The gamma-snow module uses a 3-parameter mixed distribution, and estimates a snow depletion curve for each cell modelled. The bare-ground fractionyof the cell is given as a function of melt depthλas:

y(t) =SDC(λ(t)|m, cv, y0) = y0+ (1−y0)·y1(t) (2.17) y0 is a parameter for the proportional snow-free area,

initial_bare_ground_f raction in SHyFT, and the bare-ground area, y₁, relative to the initially snow-covered area(1−y₀)is calculated by:

y₁(t) = Z λ(t)

0

p(x;m, cv) =γ 1

cv², λ(t) m·cv²

(2.18) where p() is the Gamma density, and γ(., .) is the cumulative Gamma distribution with shape cv⁻², state classalpha in SHyFT, and scalem·cv². Kolberg, Rue, and Gottschalk (2006) notes that:

...m and λ only appear as a ratio, and the observable bare- ground fraction y is sensitive only to the relative magnitude of the two.

Untilwinter_end_day_of_yearis reached, the snow is represented bym andcv. Any input of new snow, applied thesnow_cvparameter, changes the shape and scale of the SDC and updatemax_water. If the melt is calculated to exceed the storage, the SDC returns to initial values.

(41)

In the event of a melt-event in winter, the thinner snow storages will melt fastest by means of reducing alphaand keeping the scalem constant. If no winter melt-events occur theSW E_CV will have the same value assnow_cv at ablation onset, while subsequent snowfalls will adjustSW E_CV towards snow_cvincrementally (Helset et al., 2016).

2.3.3 Snow: HBV-snow

The HBV-snow⁴routine in shyft uses arbitrary quantiles to model snow, and a temperature-index approach for ablation modelling.

The temperature-index snow melt is generally described as:

∆w=







M ·(Ta−T_m) ifT_a≥T_m 0 ifT_a< T_m

(2.19) Where ∆w is the melt, M denotes the melting factor in depth per temperature per time unit, typicallymm·day⁻¹^◦C⁻¹, andT_aandT_m are air temperature and melting onset temperature, respectively (Dingman, 2008, page 210).

The model includes a snowfall redistribution vector(P.s()), an interval containing starting points for the quantiles(P.intervals())and the method uses five calibrated parameters shown in table 2.4

TABLE 2.4: Parameters with descriptions and units for HBV- snow routine.

Parameter Description Unit

cfr Refreezing coefficient -

cx Melting temperature index mm·^◦C⁻¹day⁻¹

lw Max Liquid water content mm

ts Threshold temperature for melt onset

◦C

tx Threshold temperature for

snow/rain differentiation

◦C

Model outputs are two state variables,SW Efor snow water equivalent in mm, and SCA denoting snow covered area of a cell as a fraction. For ablation the model produces outflow from the snow pack inmm(Helset et al., 2016).

4https://github.com/statkraft/shyft/blob/master/core/hbv_snow.h

(42)

2.3.4 Flow: Kirchner

The following routine is based on the paper Catchments as simple dynamical systems by James W. Kirchner in 2009, based on a case study of the rivers Severn and Wye in Wales. In the paper he further noted the possibility of doing hydrology backwards, to infer catchment precipitation directly from gauged streamflow time series. The following section will focus on the method for modelling streamflow.

The basis for this routine begins with the conservation-of-mass equation as stated in equation (2.1), expressed as,

dS

dt =P −E−Q (2.20)

where S is the stored volume of water in the catchment measured in depth (e.g., mm). And P, E and Q represent the rates of precipitation, evapotranspiration and discharge respectively, as depth per time (e.g., mm/hr). Kirchner (2009) points to the fact that among these, only the discharge is measured at an aggregated catchment scale, while precipitation, evaporation, and storage intrinsically represents regions of magni- tudes smaller than the catchment. This holds true for rain gauges for precipitation, and pan evaporation or Penman-Monteith-approach evapotranspiration and soil moisture or piezometer gauging of storage. Therefore em- phasis is placed on the discharge measurement, assuming there is a storage- discharge functionf(S)that produces the discharge

Q=f(s) (2.21)

This approximation is valid in the case where discharge sources not linked to the storage, such as overland flow, and direct precipitation into the channel, are not dominant driving forces for discharge generation. The only assumptions listed is thatdQ/dS > 0for allQandSsuch that there is an inverse function of (2.21):

S=f⁻¹(Q) (2.22)

which, as stated, will generally be a non-linear function.

With Q solely driven by S, discharge will be increasing whenever P − E > Q, and decreasing opposite. The peak discharge occurs when the discharge exactly equalsP −E in which casedQ/dt =dS/dt = 0. Even without accounting for travel time delays, the peak discharge for a storm

(43)

event will occur after the peak precipitation, so that during peak flow the rainfall rate falls below discharge turning the mass balance negative.

By differentiating (2.21) and substituting (2.20) directly yields a differential equation of discharge through time related to the change in storage as shown in (2.23)

dQ

dt = dQ dS

dS

dt = dQ

dS(P −E−Q) (2.23)

Furthermore the expressiondQ/dSis the derivative off(S)in (2.21), and can be expressed as a function ofQ:

dQ

dS =f⁰(S) =f⁰(f⁻¹(Q)) =g(Q) (2.24) Kirchner (2009) names this function, g(Q), the sensitivity function. The use of this form of mathematical relationship is more useful as a function of Qas it is directly measurable, as opposed to a function ofS, which in most cases is not. Kirchner combines (2.23) and (2.24) to estimateg(Q)based on observational data as:

g(Q) = dQ dS =

dQ dt dS dt

=

dQ dt

P −E−Q (2.25)

The implication of (2.25) is that the storage-discharge function, f(S), can be determined from instantaneous measurements of P, E, Q and the derivative of Q with respect to time. Given the nature, and difficulty, of measuring P and Q, (2.25) is best estimated in cases where P and E are small compared toQso that it may be approximated as shown in (2.26):

g(Q) = dQ dS ≈

−dQ dt

Q

PQ,EQ

(2.26) This also implies that it is possible to estimate the sensitivity function g(Q)fromQalone, if one can identify periods where the requirementP <<

Q, E << Q holds true. This does not necessitate exact measurements ofP andE, the only requirement is knowing which time periods the relationship between them behave under the correct conditions.

From the sensitivity function it is then possible to derive the storage- discharge functionf(S)by integration:

Z

dS =

Z dQ

g(Q) (2.27)

(44)

FIGURE 2.6: Illustration of catchment sensitivity function g(Q) and storage-discharge relationship f(S) as shown in

(Kirchner, 2009)

Note that the mathematical form of g(Q) or f(S) is not defined in this approach, and although for some catchments, and some situations it might be solvable analytically, an empirical estimation from streamflow time series, solved by numerical integration is needed. Kirchner (2009) based his further analysis on recession plots as proposed by Brutsaert and Nieber (1977).

From the two rivers analysed, Kirchner (2009) used hourly precipitation data and streamflows to estimate the flow recession, −dQ/dt between successive hours for periods in which Q was at least 10 times larger than the combined P and E, based on rain gauges and the Penman-Monteith evapotranspiration method, as described in equation (2.5), using weather station data.

The order of magnitude difference between streamflow and streamflow recession rates are substantial, thus a logarithmic approach was used.

To account for the large scattering found with low discharges and low recession rates, Kirchner binned individual hourly data points into ranges depending on Q. Since, as previously stated, there was no implicit mathematical form of the function, Kirchner used a quadratic curve to fit with the logarithmic values,ln(dQ/dt)andln(Q). This to achieve a balance

(45)

between smoothness and flexibility, to arrive at the form:

ln(g(Q)) = ln





−dQ dt

Q

PQ,EQ



≈c₁+c₂ln(Q) +c₃(ln(Q))² (2.28) In Kirchner (2009) the following values, and uncertainties were found for the regression parameters:

• c₁ =−2.439±0.017andc₁ =−2.207±0.028

• c₂ = 0.966±0.035andc₂ = 1.099±0.048

• c3 =−0.100±0.016andc3 =−0.002±0.018 in the rivers Severn and Wye, respectively.

To simulate a hydrograph response, this is solvable without explicitly accounting for storage at all, by combining equations (2.23) and (2.24) to get:

dQ

dt = dQ dS

dS

dt =g(Q)(P −E−Q) (2.29) When integrating equation (2.29) using Euler’s method, there is the risk of numerical instabilities due to the non-linearity of g(Q)·Q and the magnitude variations of Q. Kirchner (2009) therefore suggests integrating the log-transform of equation (2.29) as

d(lnQ) dt = 1

Q dQ

dt = g(Q)

Q (P −E−Q) = g(Q)

P −E Q −1

(2.30) For modelling streamflow, Kirchner (2009) used thec₁, c₂, c₃ parameters and calibrated against an adjustable coefficient k_E to scale evapotranspiration estimates.

Using the above mentioned approach, a Nash-Sutcliffe Efficiencies of 0.913 and 0.859 for Severen and Wye were achieved, modelled over the period 1992-1996.

The SHyFT-implemented Kirchner routine⁵ is combined with the two- parameter Priestley-Taylor evapotranspiration routine, as described in section 2.3.1, and thec₁, c₂, c₃parameters are calibrated in the model routine.

5https://github.com/statkraft/shyft/blob/master/core/kirchner.h

(46)

2.4 Model Validation and Objective Functions

2.4.1 Split-Sample Testing

The model was tested using the split-sample method proposed in Klemeš (1986). The available data is split to a calibration and validation period in time, where one segment is used purely for calibration, and the latter for validation purposes.

2.4.2 BOBYQA

The calibration module of SHyFT offers different methods of optimizing calibration performance, the method used herein relies on the BOBYQA algorithm by Michael J. D. Powell with the Nash and Sutcliffe (1970) objective function. The acronym is short for Bound Optimization BY Quadratic Approximationand is distributed under "The GNU Lesser General Public Licence (LGPL) (Powell, 2009).

2.4.3 NSE

The Nash-Sutcliffe efficiency criterion by Nash and Sutcliffe (1970) is used to verify model performance, both in calibration and validation periods. The criterion is calculated as:

N SE ≡1− PN

i=1

Q_i−Qˆ_i2

PN

i=1 Q_i−Q¯2 (2.31) where Qˆ_i is the simulated discharge and Q¯ is the average discharge across the period. Values take on the range (−∞,1], where 1 is a perfect fit, and values below 0 indicates model performance is below substituting the mean (Nash and Sutcliffe, 1970).

2.4.4 Mann-Kendall test for monotonic trend

The Mann-Kendall test for monotonic trend (Gilbert, 1987) in (PNNL, 2016) tests whether to reject or accept the hypothesis of a monotonic trend.

H0 : No trend

H_a: Monotonic trend present

(47)

by computing:

S =

n−1

X

k−1 n

X

j−k+1

sgn(xj−x_k) (2.32) which sums the number of positive differences minus the number of negative differences. In the case of a sample size nlarger than 10, proceed to calculate the variance ofSas:

V AR(S) = 1 18

"

n(n−1)(2n+ 5)−

g

X

p−1

t_p(tp−1)(2tp+ 5)

#

(2.33)

where g is the number of tied groups, and t_p is the number of observations in the p-th group. Then the Mann-Kendall test statistic Z_{M K} is computed as:

Z_{M K} =











√S−1

V AR(S) ifS >0

0 ifS= 0

√S+1

V AR(S) ifS <0

(2.34)

A positiveZ_{M K} indicates an increasing trend.

A negativeZ_{M K} indicates a decreasing trend.

TheH₀ of no increasing trend in rejected ifZ_{M K} ≥ Z1−α, whereZ1−α is the100(1−α)^thpercentile of the standard normal distribution. The opposite case of decreasing trend is true forZM K ≤ −Z1−α.

Implementation of Mann-Kendall test on historical data was performed using a python function made available through github⁶

6https://github.com/mps9506/Mann-Kendall-Trend

(48)

(49)

Chapter 3 Study Area and Data

American River Basin

The American River Basin is part of the larger Yuba, American, and Bear River subregion, situated west of Lake Tahoe in the Sierra Nevada mountain range as shown in figure 3.1. The two major subdivisions are the Upper and Lower American River watersheds, upstream and downstream of Folsom Lake, which in turn is created by the Folsom Dam (Heiman and Knecht, 2010).

The region is located on the western slope of the Sierra Nevada mountain range, a massif stretching more than 400kmfrom the Mojave desert in the south, to the Cascade mountain range in the north. The American River is located where the mountain range is at it’s widest, approximately130km, in the vicinity of Lake Tahoe (Eardley and James, 2010).

3.1 Geography and Geology

Sierra Nevada consists of predominantly granitic rocks with bands of metamorphosed sedimentary layers that merge with the volcanic rock of the Cascades to the north. Its highest peaks range from 3350mto4418mabove sea level, with the highest peaks being located south of the study catchment (Eardley and James, 2010). The highest peak within the modelled area is approximately2540m(Dollison, 2010).

The formations are a result of an uplift of the Earth’s crust, estimated to have begun several millions of years ago, with the most massive changes occurring within the last two million years. This has created an asymmetry with more abrupt topography on the east-facing slopes, and more gentle topography on the west-facing slopes, where the American River flows into the central valley region of California. Above approximately 1500mon the western flank, the Sierra Nevadas are identified by several glacially carved

A case study implementing the SHyFT modelling framework for the American

regional calibration methods for a distributed hydrologic

modelling framework

A case study implementing the SHyFT modelling framework for the American

River Basin in Sierra Nevada, California

Marco Westergren

Thesis submitted for the degree of Master in Geosciences

(Hydrology) 60 credits

Department of Geosciences

Faculty of mathematics and natural sciences

regional calibration methods for a distributed hydrologic

modelling framework

A case study implementing the SHyFT modelling framework for the American River Basin in Sierra Nevada, California

Marco Westergren

Abstract

Acknowledgements

Contents

List of Figures

List of Tables

Chapter 1

Introduction

1.1 Introduction

1.2 Motivation

1.3 Objectives

Chapter 2 Theory

2.1 Hydrology

2.1.1 Basic Concepts and Water Balance

2.1.2 Hydrologic Modeling

2.2 Hydrologic Models

2.2.1 CNRFC

2.2.2 SNOW-17

2.2.3 SAC-SMA

2.2.4 Unit hydrograph

2.3 SHyFT Stacks and Modules

Modules Used in SHyFT

2.3.1 Evapotranspiration: Priestley-Taylor

2.3.2 Snow: GammaSnow

2.3.3 Snow: HBV-snow

2.3.4 Flow: Kirchner

2.4 Model Validation and Objective Functions

2.4.1 Split-Sample Testing

2.4.2 BOBYQA

2.4.3 NSE

2.4.4 Mann-Kendall test for monotonic trend

Chapter 3

Study Area and Data

American River Basin

3.1 Geography and Geology