UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Winter stream temperature in the rain-on-snow zone of the Pacific Northwest Leach, Jason Andrew 2014

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata


24-ubc_2015_february_leach_jason.pdf [ 13.22MB ]
JSON: 24-1.0135613.json
JSON-LD: 24-1.0135613-ld.json
RDF/XML (Pretty): 24-1.0135613-rdf.xml
RDF/JSON: 24-1.0135613-rdf.json
Turtle: 24-1.0135613-turtle.txt
N-Triples: 24-1.0135613-rdf-ntriples.txt
Original Record: 24-1.0135613-source.json
Full Text

Full Text

Winter stream temperature in the rain-on-snow zoneof the Pacific NorthwestbyJason Andrew LeachBSc, University of Guelph, 2005MSc, The University of British Columbia, 2008A THESIS SUBMITTED IN PARTIAL FULFILLMENTOF THE REQUIREMENTS FOR THE DEGREE OFDoctor of PhilosophyinTHE FACULTY OF GRADUATE AND POSTDOCTORALSTUDIES(Geography)The University of British Columbia(Vancouver)December 2014c© Jason Andrew Leach, 2014AbstractStream temperature dynamics during winter are less well studied thansummer thermal regimes, but the winter season thermal regime can becritical for fish growth and development. The winter thermal regimes ofPacific Northwest headwater streams, which provide vital winter habitat forsalmonids and their food sources, may be particularly sensitive to changesin climate because they can remain ice-free throughout the year and areoften located in rain-on-snow zones. This study examined winter streamtemperature patterns and controls in small headwater catchments within therain-on-snow zone at the Malcolm Knapp Research Forest, near Vancouver,British Columbia, Canada.A diagnostic energy budget analysis highlighted that advective fluxesassociated with hillslope throughflow inputs were a dominant control onthe winter stream thermal regime. In addition, stream temperatures dur-ing rain-on-snow events were generally lower than during rain-on-groundevents after controlling for air temperature. Methods for estimating through-flow temperatures embedded in stream temperature models were evaluatedagainst field observations, and were found either not to account for therole of snow or to under- or over-predict throughflow temperatures by upto 5 ◦C. Therefore, a conceptual-parametric hillslope throughflow tempera-ture model that coupled hydrologic and thermal processes, and accountedfor the role of snow was developed and evaluated against field observa-tions of throughflow temperatures. The hillslope throughflow temperaturemodel was linked to stream energy exchange processes in order to predictstream temperature. The stream temperature model accurately predictediistreamflow and winter stream temperature at three study catchments. Themodel also simulated lower throughflow temperatures during rain-on-snowversus rain-on-ground events, although the magnitude of cooling was lessthan suggested by empirical results. A key implication of this research isthat climatic warming may generate higher winter stream temperatures inthe rain-on-snow zone due to both increased temperature of throughflowinputs and reduced cooling effect of snow cover.iiiPrefaceThe following thesis was completed under the supervision of Dr. R. DanMoore. Dr. Moore provided guidance and suggestions on the study designand analyses, as well as editorial revisions. Dr. Moore also helped conceptu-alize the hillslope throughflow temperature model outlined in Chapter 4.I am the lead author and led the conceptualization, design, data collection,data analyses and writing of this thesis.A version of Chapter 3 has been published as ‘Leach, J. A. and R. D.Moore. 2014. Winter stream temperature in the rain-on-snow zone of thePacific Northwest: influences of hillslope runoff and transient snow cover.Hydrology and Earth System Sciences. 18:819–838’.ivTable of ContentsAbstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iiPreface . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ivTable of Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vList of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ixList of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xAcknowledgements . . . . . . . . . . . . . . . . . . . . . . . . . . . . xiv1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.1 Motivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.2 Literature review . . . . . . . . . . . . . . . . . . . . . . . . . 31.2.1 Previous winter stream temperature research . . . . . 31.2.2 Advective fluxes associated with runoff processes . . 61.3 Study objectives . . . . . . . . . . . . . . . . . . . . . . . . . . 91.3.1 Study 1: Diagnostic energy budget of winter streamtemperatures in the rain-on-snow zone . . . . . . . . 101.3.2 Study 2: Hillslope throughflow temperature observa-tions and model comparison . . . . . . . . . . . . . . 101.3.3 Study 3: Development and evaluation of a conceptualhillslope throughflow model . . . . . . . . . . . . . . 10v1.3.4 Study 4: Development and evaluation of a catchmentscale stream temperature model . . . . . . . . . . . . 111.4 Thesis overview . . . . . . . . . . . . . . . . . . . . . . . . . . 112 Study site and field methodology . . . . . . . . . . . . . . . . . . 122.1 Study site . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 122.1.1 Malcolm Knapp Research Forest . . . . . . . . . . . . 122.1.2 Griffith Creek . . . . . . . . . . . . . . . . . . . . . . . 142.1.3 East and South creeks . . . . . . . . . . . . . . . . . . 152.2 Field methodology . . . . . . . . . . . . . . . . . . . . . . . . 152.2.1 Precipitation events and snow cover . . . . . . . . . . 162.2.2 Streamflow and channel geometry . . . . . . . . . . . 182.2.3 Stream temperature . . . . . . . . . . . . . . . . . . . 202.2.4 Above-stream microclimate . . . . . . . . . . . . . . . 202.2.5 Subsurface water levels and temperature . . . . . . . 212.2.6 East and South creeks . . . . . . . . . . . . . . . . . . 233 Diagnostic energy budget of winter stream temperatures in therain-on-snow zone . . . . . . . . . . . . . . . . . . . . . . . . . . . 253.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . 253.2 Energy budget study at Griffith Creek . . . . . . . . . . . . . 273.2.1 Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . 273.2.2 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . 353.3 Historical study at East Creek . . . . . . . . . . . . . . . . . . 453.3.1 Data sources . . . . . . . . . . . . . . . . . . . . . . . . 483.3.2 Modelling snow cover at East Creek . . . . . . . . . . 493.3.3 Stream temperature response to hydroclimatic eventtype . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 523.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 553.4.1 Relative roles of vertical and lateral heat fluxes . . . . 553.4.2 Role of transient snow cover and rain-on-snow events 584 Hillslope throughflow temperature observations and model com-parison . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 60vi4.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . 604.2 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 654.2.1 Spatiotemporal variability of throughflow temperature 654.2.2 Current approaches for modelling throughflow tem-perature . . . . . . . . . . . . . . . . . . . . . . . . . . 664.3 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 714.3.1 Spatiotemporal variability of throughflow temperature 714.3.2 Throughflow temperature model comparison . . . . 734.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 774.4.1 Spatiotemporal variability of throughflow temperature 774.4.2 Utility of heat as a hillslope hydrologic tracer . . . . . 794.4.3 Throughflow temperature model comparison . . . . 795 Development and evaluation of a conceptual hillslope through-flow model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 825.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . 825.2 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 835.2.1 Hydrology component . . . . . . . . . . . . . . . . . . 835.2.2 Thermal component . . . . . . . . . . . . . . . . . . . 905.3 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1005.3.1 Hydrologic calibration and evaluation . . . . . . . . . 1005.3.2 Thermal calibration and evaluation . . . . . . . . . . 1005.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1095.4.1 Evaluation of model performance . . . . . . . . . . . 1095.4.2 Rain-on-snow events . . . . . . . . . . . . . . . . . . . 1126 Development and evaluation of a catchment scale stream tem-perature model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1146.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1146.2 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1166.2.1 Model development . . . . . . . . . . . . . . . . . . . 1166.2.2 Data processing for East and South creeks . . . . . . . 1186.2.3 Calibration . . . . . . . . . . . . . . . . . . . . . . . . . 121vii6.2.4 Evaluation . . . . . . . . . . . . . . . . . . . . . . . . . 1226.3 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1236.3.1 Conditions at Griffith, East and South creeks . . . . . 1236.3.2 Split-sample and proxy-basin test . . . . . . . . . . . 1236.3.3 Pseudo-differential split-sample test . . . . . . . . . . 1306.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1316.4.1 Streamflow . . . . . . . . . . . . . . . . . . . . . . . . 1316.4.2 Split-sample and proxy-basin tests . . . . . . . . . . . 1376.4.3 Pseudo-differential split-sample test . . . . . . . . . . 1397 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1417.1 Summary of important findings . . . . . . . . . . . . . . . . . 1417.2 Implications of climate and land cover changes on winterstream temperature . . . . . . . . . . . . . . . . . . . . . . . . 1437.3 Suggestions for future research . . . . . . . . . . . . . . . . . 145References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 148viiiList of TablesTable 1.1 Summary of research on winter stream temperature. . . . 5Table 2.1 Research catchments used in this study . . . . . . . . . . . 15Table 3.1 Summary of piezometer measurements and hyporheic en-ergy exchange estimates. . . . . . . . . . . . . . . . . . . . 40Table 4.1 Summary of throughflow temperature model comparison. 77Table 5.1 List of hydrologic parameters used in the hillslope through-flow temperature model. . . . . . . . . . . . . . . . . . . . 90Table 5.2 List of thermal parameters used in the hillslope through-flow model . . . . . . . . . . . . . . . . . . . . . . . . . . . 98Table 5.3 Summary of rain-on-snow events during the 2011–2012and 2012–2013 winter periods. . . . . . . . . . . . . . . . . 107Table 6.1 Parameters used in the stream temperature model for thethree study catchments. . . . . . . . . . . . . . . . . . . . . 119ixList of FiguresFigure 2.1 Map of Malcolm Knapp Research Forest . . . . . . . . . . 13Figure 2.2 Open site automated weather station. . . . . . . . . . . . 17Figure 2.3 Bulk precipitation gauge. . . . . . . . . . . . . . . . . . . 17Figure 2.4 Harvest site snow lysimeter . . . . . . . . . . . . . . . . . 18Figure 2.5 Time lapse camera . . . . . . . . . . . . . . . . . . . . . . 19Figure 2.6 Griffith Creek v-notch weir. . . . . . . . . . . . . . . . . . 19Figure 3.1 Boxplots of monthly mean air temperature and total pre-cipitation from 1962 to 2010. . . . . . . . . . . . . . . . . . 36Figure 3.2 Hydroclimatic overview of the 2011–2012 winter period. 37Figure 3.3 Hydroclimatic overview of the 2012–2013 winter period. 38Figure 3.4 Hourly energy budget components for the forest and har-vest reaches for 2011–2012. . . . . . . . . . . . . . . . . . 42Figure 3.5 Hourly energy budget components for the forest and har-vest reaches for 2012–2013. . . . . . . . . . . . . . . . . . 43Figure 3.6 Ratios of the surface energy fluxes to the net heat inputsfor the forest and harvest reaches plotted against hourlydischarge at the lower boundary of the respective reaches. 44Figure 3.7 Comparison of back-calculated effective advective tem-perature and observed near stream subsurface tempera-tures, as well as groundwater table observations for the2011–2012 winter period. . . . . . . . . . . . . . . . . . . 46xFigure 3.8 Comparison of back-calculated effective advective tem-perature and observed near stream subsurface tempera-tures, as well as groundwater table observations for the2012–2013 winter period. . . . . . . . . . . . . . . . . . . 47Figure 3.9 Calculated mean daily effective advective temperaturefor four event types for the forest and harvest reachesduring winter 2011–2012 and 2012–2013. . . . . . . . . . 48Figure 3.10 Comparison of modelled and observed snow cover, andmeasured snow water equivalence from Griffith Creek,and modelled SWE and observed snow cover at the EastCreek, from October to April for 2011–2012 and 2012–2013. 51Figure 3.11 Comparison of modelled historic snowpack water equiv-alent for winters of 1962 to 2013. . . . . . . . . . . . . . . 52Figure 3.12 Daily mean stream temperature at East Creek as a func-tion of air temperature, classified into five event types. . 53Figure 3.13 Box plots of mean daily stream temperature at East Creekduring periods of mean daily air temperature between 0and 5 ◦C for four different event types. . . . . . . . . . . 54Figure 4.1 Water inputs, air temperature, discharge, groundwaterwell water table elevations, and subsurface temperaturesfor the 2011–2012 study period. . . . . . . . . . . . . . . . 72Figure 4.2 Water inputs, air temperature, discharge, groundwaterwell water table elevations, and subsurface temperaturesfor the 2012–2013 study period. . . . . . . . . . . . . . . . 73Figure 4.3 Subsurface temperatures for the four sites where morethan one temperature sensor was installed. . . . . . . . . 74Figure 4.4 Discharge, interquartile range of throughflow tempera-tures, mean throughflow temperature minus air tempera-ture, and regression R2 value using depth below surfaceas predictor. . . . . . . . . . . . . . . . . . . . . . . . . . . 75Figure 4.5 Comparison of different throughflow temperature estima-tion approaches. . . . . . . . . . . . . . . . . . . . . . . . 76xiFigure 5.1 Conceptual hillslope throughflow temperature model. . 85Figure 5.2 Hydrographs of calibration and evaluation periods. . . . 101Figure 5.3 Histograms of Nash-Sutcliffe efficiencies (NSE) for thehydrology calibration and evaluation periods. . . . . . . 101Figure 5.4 Calibrated model parameters and corresponding Nash-Sutcliffe coefficients for hydrology efficiencies. . . . . . . 102Figure 5.5 Hillslope throughflow temperature calibration for 2011–2012. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103Figure 5.6 Hillslope throughflow temperature model evaluation for2012–2013. . . . . . . . . . . . . . . . . . . . . . . . . . . . 104Figure 5.7 Histograms of the Nash-Sutcliffe efficiency and percent-age of time of successful temperature prediction for thecalibration and evaluation periods. . . . . . . . . . . . . . 105Figure 5.8 Comparison of predicted throughflow temperatures forthe forest and harvest snow regimes. . . . . . . . . . . . . 106Figure 5.9 Calibrated model parameters and corresponding percent-age of correct temperature prediction for the thermalmodel calibration. . . . . . . . . . . . . . . . . . . . . . . . 107Figure 5.10 Modelling exercise to evaluate model performance duringindividual rain-on-snow and rain-on-ground events. . . 108Figure 5.11 Rain-on-ground event simulated as a rain-on-snow eventwith three different snow cover scenarios. . . . . . . . . . 110Figure 6.1 Removal of snowmelt component from the Griffith Creekforest lysimeter . . . . . . . . . . . . . . . . . . . . . . . . 122Figure 6.2 Comparison of stream temperature, air temperature, vapourpressure, and runoff for Griffith, East, and South creeksfor 2011–2012 . . . . . . . . . . . . . . . . . . . . . . . . . 124Figure 6.3 Comparison of stream temperature, air temperature, vapourpressure, and runoff for Griffith, East, and South creeksfor 2012–2013 . . . . . . . . . . . . . . . . . . . . . . . . . 125Figure 6.4 Griffith Creek observed and predicted stream tempera-ture and discharge 2011–2012 . . . . . . . . . . . . . . . . 126xiiFigure 6.5 Griffith Creek observed and predicted stream tempera-ture and discharge 2012–2013 . . . . . . . . . . . . . . . . 127Figure 6.6 East Creek observed and predicted stream temperatureand discharge 2011–2012 . . . . . . . . . . . . . . . . . . . 128Figure 6.7 East Creek observed and predicted stream temperatureand discharge 2012–2013 . . . . . . . . . . . . . . . . . . . 129Figure 6.8 South Creek observed and predicted stream temperatureand discharge 2011–2012 . . . . . . . . . . . . . . . . . . . 130Figure 6.9 South Creek observed and predicted stream temperatureand discharge 2012–2013 . . . . . . . . . . . . . . . . . . . 131Figure 6.10 Histograms of goodness of fit statistics for catchment-scale model evaluation at Griffith, East and South creeks. 132Figure 6.11 Modelled throughflow temperatures at Griffith, East andSouth creeks. . . . . . . . . . . . . . . . . . . . . . . . . . 133Figure 6.12 Comparison of observed and modelled relative tempera-ture differences between Griffith and East or South creeks. 134Figure 6.13 Comparison of hourly modelled stream temperature re-sponse to pre-harvest conditions for 2011–2012. . . . . . 135Figure 6.14 Comparison of hourly modelled stream temperature re-sponse to pre-harvest conditions for 2012–2013. . . . . . 136xiiiAcknowledgementsThis thesis, and the enjoyment that has come with working on it, are thedirect and indirect product of countless people. First and foremost, I amindebted to my supervisor, Dan Moore, for providing invaluable guidanceand direction throughout this endeavour. It has been a great pleasureworking with Dan. I would also like to thank my committee members,John Richardson and Mark Johnson, for their input at critical stages in theresearch. Ionut Aron, Cheryl Power, and Paul Lawson at the MalcolmKnapp Research Forest provided important field logistic support. Variousmembers of the UBC Department of Geography, specifically Sandy Lapksy,Brett Eaton, Ian McKendry, and Marwan Hassan, assisted me along theway and created a friendly work environment during my time at UBC. JoelTrubilowicz and Derek van der Kamp were great office mates and werealways willing to act as a sounding board. Help in the field was providedby Eli Heyman, Spencer Pasieka, Saskia Ho¨velmann, and Justin Knudson.Scholarships from the Natural Sciences and Engineering Research Coun-cil, Canadian Geophysical Union, Canadian Water Resources Association,and UBC provided financial support that allowed me to focus on my studies.I would also like to thank the open source software community, as nearlythe entire thesis and its contents were produced using open source software.I would like to thank Mom, Dad, Craig, and my grandparents for theirlove and support during the years. This thesis is dedicated to Kaitlin andClem. Kaitlin has supplied me with all the understanding, love, and laughterI have needed. Ever since he was an itsy bitsy zygote, Clem has providedjust the motivation I needed to finish this thing.xivChapter 1Introduction1.1 MotivationStream temperature is a principal determinant of aquatic ecosystem compo-sition and productivity (Wehrly et al., 2003, 2007; Friberg et al., 2013). Thereare increasing concerns that changes in land cover and climatic conditionscould produce changes in stream thermal regimes that would be deleteriousto existing aquatic communities (Brown et al., 2007; Durance and Ormerod,2007, 2009), particularly for cold- and cool-water species such as salmonids(Battin et al., 2007). These concerns have generated a need for process-basedapproaches for analysing and predicting stream temperature variability ata range of spatial and temporal scales (Webb et al., 2008; Arismendi et al.,2012).Most stream temperature research has been concerned with maximumsummer temperatures, particularly in response to riparian vegetation distur-bance such as harvesting or wildfire (e.g. Johnson and Jones 2000; Bartholow2005; Gaffield et al. 2005; Gomi et al. 2006; Leach and Moore 2010; Imholtet al. 2013). This summer-focused research consistently identifies energyexchanges occurring at the stream surface (primarily net radiation) as akey control on temperature variability on diurnal to seasonal scales andthe response to changes in riparian canopy conditions (Caissie et al., 2007;Hannah et al., 2008; MacDonald et al., 2014b; Garner et al., 2014). In some1reaches, however, hyporheic exchange or groundwater discharge can sig-nificantly moderate temperature variability (e.g. Johnson 2004; Leach andMoore 2011).In contrast to the depth and breadth of research on summer stream tem-perature dynamics, few studies have examined winter stream temperatureprocesses despite its recognized importance for aquatic ecosystems (Beschtaet al., 1987; Holtby, 1988; Ebersole et al., 2006; Brown et al., 2011; Shuter et al.,2012). Salmonids are poikilothermic, so when stream temperatures declineduring winter, so do their metabolic processes and their ability to swim,feed, avoid predators, and defend their locations from other fish (Brownet al., 2011). Therefore, predation by homeothermic predators during winteris believed to be high, particularly when surface ice cover is non-existentor incomplete (Huusko et al., 2007). It is not clear how future changes inair temperature and precipitation regimes will impact winter water temper-ature of coastal mountainous streams, such as those found in the PacificNorthwest (PNW). The winter thermal regimes of PNW headwater streams,which provide vital winter habitat for salmonids, may be particularly sen-sitive to changes in climate because they can remain ice-free throughoutthe year and are often located in the rain-on-snow zone (Jones and Perkins,2010). Currently, lack of research on winter stream temperature processeslimits our ability to predict winter stream temperature response to climateand land cover changes.A detailed review of previous research on winter stream temperatureis provided in the following section. This review highlights that advectivefluxes associated with runoff processes are likely to be important stream tem-perature controls during winter periods. Therefore, the handfull of studiesthat have examined relationships between runoff processes and stream tem-perature, albeit during summer periods, are also reviewed. The literaturereview is followed by a statement of key research objectives to address themajor knowledge gaps in understanding of winter stream thermal regimes.This section ends with an outline of the thesis structure.21.2 Literature review1.2.1 Previous winter stream temperature researchPrevious stream temperature research has primarily focused on the influ-ence of riparian vegetation on summer thermal regimes (e.g. Brown 1969;Gomi et al. 2006; Hannah et al. 2008). Few studies have focused on winterperiods, despite its recognized biological importance (Holtby, 1988; Tetzlaffet al., 2005; Ebersole et al., 2006; Xu et al., 2010; French et al., 2014). Thelack of winter stream temperature research may be due to the assumptionthat winter stream temperatures are relatively invariable and remain at ornear 0 ◦C, particularly in environments where seasonal ice cover forms(Needham and Jones, 1959). Additionally, winter stream temperatures maynot be regarded as biologically important as summer temperatures. How-ever, in coastal temperate environments, where streams remain ice-free andstream discharge variability is greatest, winter stream thermal regimes canbe highly variable and have a considerable influence on resident fish andinsect populations (Holtby, 1988).Studies focused on winter stream temperature are summarized in Ta-ble 1.1. They include an empirical analysis of inter-annual variability ofwinter stream temperature in a forested headwater catchment in the rain-on-snow zone of south coastal British Columbia (Kiffney et al., 2002). At thatsite, winter stream temperature was related to large scale climate processes(El Nin˜o Southern Oscillation and Pacific Decadal Oscillation) that influ-enced stream discharge and air temperature. An increase in discharge ordecrease in air temperature resulted in a decrease in winter stream tempera-ture, while a decrease in discharge or increase in air temperature resultedin an increase in winter stream temperature. Kiffney et al. (2002) foundconsiderable variability between years that could not be explained by themagnitude of large scale climate patterns, which suggests other processesare influencing winter stream temperatures.Seasonal and daily energy budget studies have found that the mag-nitudes of energy exchanges occurring at the stream surface are smaller3during winter than during summer periods (Webb and Zhang, 1997, 1999;Hannah et al., 2008; Leach and Moore, 2010; Garner et al., 2014). However,the relative importance of the surface energy exchanges relative to energyexchanges occurring at the streambed (groundwater advection, hyporheicexchange, and bed conduction) is variable. Some studies have found that,although stream surface energy exchanges are relatively small in magnitudeduring winter compared to summer, they remained the primary streamtemperature control (Webb and Zhang, 1997, 1999; Hannah et al., 2008).However, these studies were conducted at point scales or over short reaches(less than 40 m length) and did not focus on the role of advective fluxesassociated with groundwater and runoff contributions (Story et al., 2003;Leach and Moore, 2011).Studies from winter periods have provided indirect evidence for theimportance of advective fluxes associated with runoff and groundwatercontributions in controlling winter stream temperatures (Webb and Zhang,1999; Malard et al., 2001; Leach, 2008; Danehy et al., 2010). Groundwateradds heat to the stream and, because groundwater typically has a highertemperature than stream water during winter, it has a warming effect. Lan-gan et al. (2001) examined stream temperature dynamics at decadal andevent scales at Girnock Burn, Scotland. They found that mean daily maxi-mum winter and spring stream temperatures had each increased by about2 ◦C between 1968 and 1997. They attributed this finding to an increasein air temperature, which may have also reduced the magnitudes of snowand snowmelt on the catchment during winter and spring seasons. Forprecipitation events, Langan et al. (2001) found that during the winter therewas a positive relation between stream discharge and water temperature.They speculated that this was because precipitation falling on the catchmentwas warmer than the antecedent soil and stream temperatures, particularlywhen a snowpack was present. This speculation was not evaluated but itdoes support the hypothesis that precipitation events and catchment-scaleprocesses (i.e. runoff processes and groundwater dynamics) are importantwinter stream temperature controls.4Table 1.1: Summary of research on winter stream temperature.Reference Location Riparian vege-tationHydrologicregimeStudy period General findingsKiffney et al.(2002)East Ck, BC Dense forest Rain and snow 1981 to 2000 Water temperatureswere higher and precipitationwas lower dur-ing in-phase or warm PDO/El Nin˜o events than in other years.Danehyet al. (2010)4 reaches,Rocky Moun-tainsForested Snowmelt Nov to Apr Differences in winter stream temperatures were associated withamount of groundwater discharge. Warmer streams had greatergroundwater contributions.Langan et al.(2001)Girnock Burn,UKHeather moor-landRain and snow 1968 to 1997 Mean daily maximum temperature increased during winter andspring by about 2 ◦C. Attributed to an increase in air temperaturewhich may have also reduced the snow and snowmelt on thecatchment during the winter and spring seasons.Leach andMoore(2010)Fishtrap Ck,BCWildfire-disturbedSnowmelt May 2007 toMarch 2008Net radiation, sensible and latent heat exchanges were nega-tive and were of roughly similar magnitude from November toMarch.Webb andZhang(1999)River Piddleand Bere, UKOpen field andpasturesRain 1994 Net radiation was the dominant heat gain (∼90% of non-advective fluxes) during both winter and summer, althoughabout a third less in magnitude during the winter. Relative warmwater temperatures during winter enhanced heat losses by evap-oration. Groundwater contributions added heat during winterand summer, but temperature differences between groundwaterinflows and the streammean that the input had a warming effectin winter but a cooling one in summer. Precipitation inputs tothe stream were a negligible heat input.Webb andZhang(1997)11 reaches, UK Moorland todense decidu-ous woodsRain May 1992 toDecember 1993They authors provide a conceptual model of seasonal variationsin energy budget components and suggest vertical energy ex-changes control winter stream temperature dynamics.Malard et al.(2001)Roseg River,Swiss AlpsSparsely vege-tatedGlacier andsnowmeltJune 1997 toJune 1998Groundwater from different sources had different effects onthe seasonal regime of temperature in the hyporheic zone. In-flow of shallow alluvial groundwater had minimal effects onseasonal patterns of hyporheic temperature, whereas upwellingfrom deep alluvial and hillslope aquifers resulted in significanttime lags and differences in seasonal amplitudes between surfaceand hyporheic temperatures.Hannahet al. (2004)Girnock Burn,UKMoorland andwoodlandRain andsnowmeltOct 2001 to Apr2002Energy budget terms and temperatures exhibit seasonal changesin response to meteorologic and hydrologic conditions.Hannahet al. (2008);Garner et al.(2014)Girnock Burn,UKMoorland andwoodlandRain andsnowmeltJan 2003 to Dec2009Net radiation is the dominant heat sink in autumn-winter andmajor heat source in spring-summer.51.2.2 Advective fluxes associated with runoff processesThe preceding literature review has highlighted that winter energy ex-changes at the stream surface are of lower magnitude than during thesummer, and that advective fluxes associated with groundwater and runoffprocesses may be important stream temperature controls during these pe-riods, especially for environments such as coastal PNW where winter isdominated by heavy cloud cover, frequent precipitation events and variablestreamflow. Although not specific to winter or PNW, research from Japanhas identified that the amount of rainfall affects stream temperature bychanging the proportions of different hillslope runoff following flowpaths(surface runoff, subsurface runoff, and groundwater) (Subehi et al., 2010).Differences in topography and subsequent differences in timing and deliv-ery of precipitation will also influence stream temperature because steepslopes produce more rapid surface runoff, resulting in stream temperaturemore closely matching rain temperature (Subehi et al., 2010). In contrast,gentle slopes promote precipitation infiltration, which displaces and movescooler groundwater into the stream (Subehi et al., 2010). Other research fromJapan, which compared stream temperature response for snowmelt and rainevents, found that stream temperature during the snowmelt period waswarmer than the snowmelt temperature of 0 ◦C, suggesting that streamflowwas derived from meltwater that has either warmed in transit through thesoil, or has displaced older, warmer water stored in the soil (Kobayashiet al., 1999). Stream temperature during summer rain events was lower thanthe rain water, suggesting that streamflow was derived from rainfall thathas either cooled in transit through the soil, or has displaced older, coolerwater stored in the soil. The temperature of stream water during stormflowrecession, during both snowmelt and summer rainstorms, were similar tothose of the soil 1.8 m below the land surface (Kobayashi et al., 1999).Work from the alpine region of the French Pyrene´es has identified thatspatial and temporal differences in stream thermal response were controlledby precipitation, discharge, antecedent basin conditions, and runoff flow-paths (Brown and Hannah, 2007). Most storms resulted in decreased stream6temperature. For low-intensity rain storms, precipitation inputs may beequilibrated (in terms of heat exchange) and insufficient in relative contribu-tion to induce a thermal response (Brown and Hannah, 2007). Additionally,Brown and Hannah (2008) found that stream temperature was influencedby stream water source and was associated with basin characteristics, mostnotably altitude, azimuth and stream length.To the author’s knowledge, only one study has coupled streamflowand temperature measurements at the catchment outlet with measurementsof within-catchment processes, such as water table elevations. At PanolaMountain, Georgia, Shanley and Peters (1988) compared streamflow andthermal response for four rain events during 1986 and 1987. They found thatfor dry periods, groundwater provided minimal contributions to the eventflow response, both in terms of volume and heat. However, groundwaterand its thermal signature contributed to event flow during wet periods. Theinfluence of groundwater contributions on stream temperature dependedon the relative temperatures of the groundwater and stream water.Other than the work by Shanley and Peters (1988), all the studies re-viewed here that have focused on stream temperature response to hydrologicevents have focused on water temperature and discharge at the catchmentoutlet and related the observed patterns to catchment or precipitation char-acteristics. The studies did not rigorously evaluate the physical processesbehind the observed patterns, but only speculated on the physical processescontrolling stream temperature. Some studies have applied coupled hy-drologic and stream temperature models and found that advection andrunoff can contribute significantly to the overall heat budget of a stream(St-Hilaire et al., 2003; Brookfield et al., 2009; Birkinshaw and Webb, 2010;Ficklin et al., 2012; Sun et al., 2014); however, the models assume a numberof conditions, such as a constant groundwater temperature, that field studieshave shown not to be valid. In addition, none of these modelling studiesevaluated their estimates of groundwater or hillslope runoff temperaturesagainst observations.The empirical studies reviewed above have found relations betweenstream temperature and runoff at the catchment scale, but do not provide7a basis for conclusively identifying and modelling the physical processesresponsible for the observed responses. In contrast, heat budget analysisis considered the most rigorous approach to understanding the processesgoverning stream thermal regimes (Boughton et al., 2012). A fact high-lighted here, as well as more recently recognized in the literature, is thatadvective fluxes may be a significant control not only on winter stream tem-perature, but on summer stream temperature as well (St-Hilaire et al., 2003;Herb and Stefan, 2011; Hebert et al., 2011; Deitchman and Loheide, 2012;Boughton et al., 2012; Arismendi et al., 2012; Mayer, 2012; Ficklin et al., 2012;MacDonald et al., 2014a). Despite this growing recognition, there has yetto be conducted a detailed heat budget analysis that rigorously quantifiesthe spatiotemporal variability in the advective flux component of a streamthermal regime.The majority of previous heat budget studies have chosen to ignore theadvective flux while focusing on vertical energy exchanges at the streamsurface and streambed (Brown, 1969; Vugts, 1974; Evans et al., 1998; John-son, 2004; Hannah et al., 2008; Hebert et al., 2011; Benyahya et al., 2012).Others have compared streamflow at the upstream and downstream reachboundaries and attributed any net differences to groundwater losses orgains (Pluhowski, 1970; Webb and Zhang, 1997, 1999; Cozzetto et al., 2006;Leach and Moore, 2011), although some studies have complemented thisapproach with within-reach estimates of groundwater – surface water inter-actions (Story et al., 2003; Guenther, 2007). Other heat budget studies haveestimated the advective flux associated with groundwater discharge andhillslope runoff by assuming it was proportional to drainage area (Mooreet al., 2005b) or by inferring it using temperature data (Hannah et al., 2004).Therefore, there is a need for a better understanding of and ability to modeladvective influences on stream thermal regimes.This literature review points to an increasing acknowledgement thatcatchment-scale runoff processes are a dominant control on winter streamtemperature. Understanding catchment scale runoff processes has been amajor focus of recent hydrologic research (Buttle, 2006; McDonnell et al.,2007; Tetzlaff et al., 2010). Although water temperature has been used as a8hydrologic tracer to understand runoff processes (Shanley and Peters, 1988;Kobayashi et al., 1999; Subehi et al., 2010; Birkinshaw and Webb, 2010) andgroundwater - surface water interactions (Briggs et al., 2012; Vogt et al.,2012), water temperature is a non-conservative tracer and identifying watersources based on thermal signatures will be confounded by energy exchangeprocesses. However, none of this process-based hydrologic research hasbeen coupled with an exploration of catchment-scale stream temperature dy-namics. Therefore, it would be valuable to understand interactions betweencatchment runoff processes and stream temperature, not only to furtherunderstand stream temperature controls, but also to further understandingof catchment runoff processes.1.3 Study objectivesAs shown in the preceding section, previous research on winter streamtemperature dynamics highlights uncertainty in the relative importance ofsurface energy exchanges and advective fluxes associated with runoff andgroundwater contributions. In the PNW, winters are dominated by frontalstorm systems resulting in large amounts of precipitation and persistentcloud cover. Therefore, it follows that for PNW streams, advective fluxesmay dominate the stream thermal regime during winter. How these advec-tive fluxes control stream thermal regimes is uncertain as there has beenlack of process-based research examining stream temperature response toprecipitation events at catchment-scales.The overall objective of this research is to enhance our understandingof and ability to model winter stream thermal regimes in the rain-on-snowzone of the Pacific Northwest, and their sensitivity to changes in climateand forest cover. This overall objective is achieved through four studies thatcorrespond to chapters 3, 4, 5, and 6. Outlines of these studies are providedbelow.91.3.1 Study 1: Diagnostic energy budget of winter streamtemperatures in the rain-on-snow zoneThe first study uses an energy budget analysis to diagnose the dominantenergy exchanges controlling winter stream temperature in headwater catch-ments of the PNW. This study addresses the hypothesis that winter streamtemperatures in the rain-on-snow zone are primarily controlled by advectivefluxes associated with hillslope runoff processes. An additional hypothesisis also addressed: that rain-on-snow events will have lower stream temper-atures than rain-on-ground events because precipitation moving throughthe snow pack will cool prior to entering the stream channel. Without atransient snow pack, this cooling process would not occur.1.3.2 Study 2: Hillslope throughflow temperature observationsand model comparisonThe second study focuses on hillslope runoff processes and seeks to under-stand the spatiotemporal variability in throughflow temperatures. Statisticalmodels are used to understand relationships between spatial variability inthroughflow temperatures and site variables that may represent heat ex-change processes, such as vertical heat conduction and lateral heat advectionassociated with throughflow delivered from upslope. In addition, six exist-ing approaches to estimate throughflow temperatures that are embeddedin catchment-scale stream temperature models are evaluated against fieldobservations of throughflow temperatures.1.3.3 Study 3: Development and evaluation of a conceptualhillslope throughflow modelA new conceptual-parametric hillslope throughflow temperature modelis developed that simulates water and heat exchange processes, and canaccount for transient snow cover. The model structure is designed to achievea balance between overly simplistic approaches to estimate throughflowtemperature embedded in stream temperature models (Study 2) and morecomplex groundwater models. The model is calibrated using a generalized10likelihood uncertainty estimate approach (Beven and Binley, 1992). Themodel is evaluated against measured throughflow temperatures, as well asits ability to simulate the influence of transient snow cover during rain-on-snow events.1.3.4 Study 4: Development and evaluation of a catchment scalestream temperature modelThe fourth study builds on the third study and develops a catchment-scalepredictive stream temperature model that links the hillslope throughflowtemperature model to energy exchange processes occurring at the streamsurface. The model is developed and parameterized using detailed fieldinformation available at the primary study catchment (Griffith Creek), thentested at two other catchments in the region (South and East creeks). Themodel is evaluated further by testing its ability to reproduce the temperaturechange between pre- and post-harvest conditions at Griffith Creek.1.4 Thesis overviewDetails on the study site and field methodology that are common to thefour study chapters are found in Chapter 2. Chapters 3, 4, 5 and 6 detailthe four studies that address the overall objective of the thesis. The thesisconcludes with a synthesis of the key findings and suggested future researchin Chapter 7.11Chapter 2Study site and fieldmethodology2.1 Study site2.1.1 Malcolm Knapp Research ForestThis study was conducted at the University of British Columbia’s MalcolmKnapp Research Forest (MKRF), located at 49◦ 16’ N and 122◦ 34’ W, about60 km east of Vancouver (Figure 2.1). The area experiences a maritimeclimate with wet mild winters and warm dry summers. Mean annual pre-cipitation at the Environment Canada climate station located at the researchforest headquarters (Haney UBC RF Admin, station number 1103332), ele-vation 147 masl, is 2184 mm, of which 70% falls, primarily as rain, betweenOctober and April due to Pacific frontal systems. Snowfall comprises only5% of the total annual precipitation at the low elevation headquarters butincreases at higher elevations. Based on data from 1962 to 2013, mean annualair temperature is 9.6 ◦C and mean monthly temperatures for January andJuly during this period are 2.4 and 17.3 ◦C, respectively.Most of MKRF lies in the Coastal Western Hemlock biogeoclimatic zone,with the highest elevation bands lying in the Mountain Hemlock zone12Figure 2.1: Location of the Malcolm Knapp Research Forest, the threestudy catchments (Griffith Creek, East Creek, and South Creek),and the open site meteorological station and snow lysimeter. Q1,Q2, and Q3 indicate the three streamflow gauging sites.13(Krajina, 1965). MKRF has experienced considerable harvesting and forestdisturbance over the last century and therefore comprises a patchwork ofvariable forest stand compositions and ages, dominated by second growthwestern red cedar (Thuja plicata), Douglas-fir (Pseudotsuga menziesii) andwestern hemlock (Tsuga heterophylla).Soils are primarily coarse-textured humo-ferric podzols (Feller and Kim-mins, 1979). Soil depths range up to 2 m, with compacted till or graniticbedrock found on average at 1 m depth. Bedrock outcrops occur in placesalong topographic divides. Average hydraulic conductivities are typically10−4 to 10−3 m s−1 in the soil and 10−7 to 10−6 m s−1 in the underlying till(Utting, 1979; Cheng, 1988). Partially as a result of tree-throw, the forest floorhas a complex topography, with marked changes in convexity/concavityover distances of 2–5 m.Owing to the high permeability of the soils, almost all water reaching theground surface infiltrates the soil and flows downslope in a saturated layerabove the contact between the soil and underlying till or bedrock. Stream-flow typically responds rapidly to rainfall. Runoff generation processes aredominated by subsurface flow and, to a lesser extent, saturation-excess over-land flow (Thompson and Moore, 1996; Donnelly-Makowecki and Moore,1999; Hutchinson and Moore, 2000; Haught and van Meerveld, 2011). Moststormflow occurs in the autumn-winter wet season and many streams dryup during the summer.2.1.2 Griffith CreekThe bulk of field work was conducted at Griffith Creek. At the locationof a hydrometric weir, Griffith Creek drains an area of 10.8 ha, ranging inelevation from from 365 to 572 m. Griffith Creek was previously the siteof research examining the summer stream temperature response to forestharvesting (Guenther, 2007; Guenther et al., 2012, 2014). The upper portionof Griffith Creek is characterized by channel gradients of approximately 20%in an incised channel, with a decrease in steepness in the downstream direc-tion to less than 7% slope. Bed materials change in composition from large14Table 2.1: Research catchments used in this study.Stream Catchment area (ha) Elevation range (masl)East Creek 44.3 295 - 555South Creek 19.7 175 - 320Griffith Creek 10.8 365 - 572cobbles in the headwaters to sand downstream with increasing amountsof organic matter in the lower 150 m of the study reach (Guenther, 2007).Wetted channel widths generally average around 1 m at typical winter flows.In autumn of 2004, the lower section of the catchment was logged un-der a partial retention approach, resulting in 50% of the basal area beingremoved and a 14% reduction in canopy closure (Guenther et al., 2012).Currently, the lower section’s vegetation cover is composed of mature treesand a shrub understory about 1 to 3 m in height. The shrub understoryexperiences considerable dieback during the winter months. The channelis not heavily incised and shading from banks is minimal. However, steepslopes, especially on the west side of the creek provide topographic shading.2.1.3 East and South creeksTwo additional headwater catchments in the MKRF, East and South creeks(Table 2.1), were monitored. East Creek’s catchment is dominantly coveredby mature forest. Forest stands are about 90 years old, having been loggedin the 1920s, and have crown closures of 75–95%. South Creek’s catchmentwas partially logged with a 30 m riparian buffer during the autumn andwinter of 1998 and is thus composed of mature trees and shrub and treeregrowth of about 1 to 3 m in height, outside the riparian buffer.2.2 Field methodologyA two year field campaign was undertaken at Griffith, East and Southcreeks. Griffith Creek was the primary study catchment and most of the15field equipment and effort were focused there. East and South creeks wereused to evaluate the predictive stream temperature model and thereforehad a reduced field component consisting primarily of streamflow, streamtemperature, and above-stream microclimate monitoring. Instrumentationwas installed during summer 2011 and monitoring was continued untilspring 2013. The field program and data collection focused on the autumn-winter winter period (October to April).The primary goal of the field component of the project was to monitorthe stream temperature dynamics at the study creeks and to quantify theenergy and water fluxes controlling those observed stream temperaturepatterns, particularly during periods of snow cover. An overview of thefield research components is provided in the following sections.2.2.1 Precipitation events and snow coverAt a meteorological station located at the open site (Figure 2.2), a tippingbucket rain gauge was used to measure liquid precipitation, a bulk precipi-tation gauge measured rain and snowfall (Figure 2.3), and a snow lysimetermeasured snowmelt and rainfall. The bulk precipitation gauge consistedof a 1.2 m length of PVC pipe (20.32 cm diameter) sealed at the bottom andequipped with a pressure transducer. The pipe contained antifreeze to meltsnow and the water in the pipe was covered by a thin layer of mineral oil tominimize evaporation. The tipping bucket and bulk precipitation gaugeswere logged by Campbell Scientific CR10x loggers and data were storedevery 10 minutes. The snow lysimeter (4 m2 area and 0.0625 mm per tip) wasconstructed following the design by Smith (2011) and was logged with anEm5b Decagon data logger (Decagon Devices, Pullman, Washington). Totaldrainage from the lysimeter (rainfall plus snowmelt) was stored every 1 h.Two additional snow lysimeters with the same specifications as the opensite lysimeter were installed at Griffith Creek, one in the harvested area(Figure 2.4) and one in the unharvested area. Both sites also had tippingbucket rain gauges installed below the forest canopy, which were used toconfirm occurrence of rainfall.16Figure 2.2: Open site automated weather station.Figure 2.3: Bulk precipitation gauge installed at the open site to mea-sure precipitation falling in both rain and snow phases.17Figure 2.4: Harvest site snow lysimeter.Manual surveys of snow depth and density were made during sitevisits at the open site, and at the Griffith Creek harvest and forest sites,when a snow pack was present. Snow depth was measured with a rulerat approximately 2 m intervals along a single transect at each site. Thetransects were 80 to 100 m in length and were located in the vicinity of themeteorological stations. Snow pack density was measured using a Federalsnow sampler. A minimum of five density measurements were made ateach site. Snowpack water equivalent was computed as the product of meandepth and mean density, divided by the density of liquid water.Time-lapse cameras installed at the open site, the Griffith Creek harvestand forest sites, and the outlets of East and South creeks captured imagesdaily at 9:00, 12:00 and 15:00 Pacific Standard Time (Figure 2.5). Images wereused to map snow extent and cover and, in conjuction with meteorologicaldata, help identify occurrences of rain-on-snow events (Floyd and Weiler,2008; Garvelmann et al., 2013).2.2.2 Streamflow and channel geometryStreamflow was monitored at the Griffith Creek catchment outlet (Q1), aswell as at two additional locations along the stream reach (Q2 and Q3).The drainage areas for Q1, Q2 and Q3 are 10.8, 6.6 and 0.8 ha, respectively.Manual streamflow measurements were made using the constant-rate injec-tion salt dilution method (Moore, 2004). Streamflow measurements were18Figure 2.5: Time lapse camera used to determine periods of transientsnow accumulation and melt.Figure 2.6: Griffith Creek v-notch weir.accurate to ±5%, based on replicated gauging. The Q1 station had a v-notchweir and pressure transducer (Figure 2.6), whereas Q2 and Q3 were outfittedwith stilling wells and pressure transducers. Rating curves were developedfor each location in order to estimate continuous streamflow records.Average wetted width and depth of the stream were determined frommeasurements made at 25 locations distributed along Griffith Creek (11 and14 in the harvest and forest reaches, respectively). Sixteen sets of 25 widthand depth surveys were made over a range of stream discharges. Meanwidth and depth for the lower harvested and upper forested reaches wereregressed against discharge in order to fit empirical power-law relationsfor predicting width and depth for the entire study period. Mean widths19ranged between 0.49 and 0.97 m with a relative standard error of ±21% forthe harvest reach, and 0.28 and 0.62 m (±19%) for the forest reach. Meandepths ranged between 0.10 and 0.27 m (±26%) for the harvest reach, and0.08 and 0.21 m (±37%) for the forest reach.2.2.3 Stream temperatureEight submersible temperature loggers (Tidbit v2 Temp, Onset ComputerCorporation, accurate to ±0.2 ◦C) were distributed along Griffith Creek.Three of the locations correspond to streamflow gauging sites Q1, Q2, andQ3. Temperature loggers were installed at sites with sufficient water depthso that loggers were not exposed during low flows. The loggers wereshielded with white PVC pipe with drilled holes to facilitate water exchange(Quilty and Moore, 2007). Sensors were logged at 15 min intervals andaveraged every 1 h. Sensors were calibrated at 0 and 20 ◦C (in an ice bathand at room temperature, respectively), before and after field deployment.Vertical and lateral manual spot measurements were made using a WTW340i handheld conductivity and temperature meter (accurate to ±0.1 ◦C) ateach logger site during site visits to ensure that the logger was placed in a lo-cation with full vertical and lateral mixing. Manual spot measurements werealso used to check logger records for drift throughout the study. Differencesbetween logger and manual measurements were within ±0.2 ◦C.2.2.4 Above-stream microclimateTwo automated weather stations were installed at sites within 2 m of theGriffith Creek channel to characterize the above-stream microclimate in theforested and harvested reaches. The weather stations monitored air temper-ature with a HMP45C-L probe (accurate to ±0.3 ◦C), relative humidity witha HMP45C-L probe (accurate to ±3% for the 0-90% relative humidity rangeand±5% for the 90–100% relative humidity range), incoming solar radiationwith a CMP3 Kipp and Zonen pyranometer, and wind speed with a Met One3-cup anemometer (starting threshold of 0.45 m s−1) at approximately 1.5 mabove the ground surface, and rainfall with a tipping bucket rain gauge20(0.254 mm per tip) located between 0.3–0.4 m above the ground surface.All sensors were scanned every 10 s and averaged (or summed for the raingauges) every 10 min by Campbell Scientific CR10x data loggers.Above-stream net radiation was measured using a roving Kipp andZonen net radiometer to test predictions from a net radiation model. Thenet radiometer was set up at fourteen locations (six and eight locations inthe forest and harvest reaches, respectively) along Griffith Creek during thestudy period in order to sample spatial variability in net radiation due tovariable riparian vegetation structure. At each location, the net radiometerwas positioned approximately 30–40 cm above the stream surface. Thesensor was scanned every second and readings were averaged every 10 minusing a Campbell Scientific CR10x logger.Radiative exchanges were modelled along the stream surface using hemi-spherical canopy images and meteorological data, in order to account forspatial variability in riparian canopy structure and surrounding topography(Moore et al., 2005b; Leach and Moore, 2010). A Nikon fisheye convertedFC-E8 lens and Nikon Coolpix 4500 4.0 megapixel digital camera, set onfisheye mode and highest image quality, were used to capture the images.Eighty-two images (34 at the harvest reach and 48 at the forest reach) weretaken from the centre of the stream at 5 to 10 m intervals along the entirestudy reach. The camera was levelled and oriented to north prior to theimages being captured. Images were taken approximately 10–20 cm abovethe stream surface. Hemispherical images were also taken at the locationswhere net radiation was measured and at the two stream meteorologicalstation pyranometers. Images were taken on an overcast day to ensureuniform sky conditions to facilitate image processing.2.2.5 Subsurface water levels and temperatureShallow groundwater levels were monitored at 50 wells that were installedby hand augering to the soil–till interface (mean depth was approximately0.6 m). Wells were made from PVC pipe (12.7 mm inside diameter) withholes drilled along the lower half and screened with permeable garden21fabric to prevent soil movement into the well. The wells were locatedwithin 5–10 m of the stream in order to characterize throughflow inputs tothe channel. Wells were located to sample a range of hillslope sizes andshapes determined from different topographic indices, including upslopecontributing area and topographic wetness index. Forty wells were outfittedwith Odyssey capacitance water level loggers (Dataflow Systems Pty Ltd,Christchurch, New Zealand) to provide 15 min interval water level records.The remaining wells were monitored manually once a week during thewinter period. The groundwater well network was surveyed using a LeicaFlexline TS06 total station (Leica Geosystems AG, Switzerland).Soil temperatures were recorded within 1 m of 43 of the groundwaterwells. At four sites, soil temperatures were recorded by thermocouplesat three depths (0.05, 0.25 and 0.5 m), which were connected to CampbellScientific CR10x data loggers. Thirty-six of the groundwater well sites werefitted with submersible temperature loggers installed at the soil-till interfaceusing a hand auger and backfilled after installation. Four of these sites hadan additional soil temperature logger installed at a second depth abovethe soil-till interface. In addition, soil temperatures were recorded at twodepths (0.05 and 0.15 m) at each of the snow lysimeter sites located withinthe Griffith Creek catchment. All soil temperature sensors logged data at15 min intervals with the exception of the snow lysimeter sites, where soiltemperatures were logged at 60 min intervals. Temperature loggers weredetermined to be below or above the water table at each time step basedon the known depth of the temperature logger and the corresponding wellwater level.Piezometers were installed in the centre of the streambed at twenty-fivelocations at approximately even spacing along the study reach. At fourstep-pool sequences, piezometers were installed upstream of the step anddownstream in the pool. Piezometers were constructed from 6 mm internaldiameter plexiglass tubing with holes drilled in the bottom 5 cm, and wereinstalled to depths of 10–30 cm below the streambed. Vertical hydraulicgradients were computed as the difference in water level between the insideand outside of the piezometer, divided by the depth of the mid-point of the22perforated section (Scordo and Moore, 2009). Positive hydraulic gradientsindicate upwelling flow. Water levels inside and outside the piezometerswere measured with an electronic beeper and accurate to ±5 mm (Guen-ther, 2007). Piezometers were also used to measure saturated hydraulicconductivity of the bed sediments using a falling head slug test (Freeze andCherry, 1979). Hydraulic conductivity (Ksat) was computed based on anequation derived by Hvorslev (1951) and modified by Baxter et al. (2003) forclosed-bottom perforated piezometers.At each of the twenty-five piezometer locations, thermocouples wereinstalled within the streambed. Twenty-two locations had three thermocou-ples installed at various depths and three locations had only two thermo-couples, due to the shallow depth of bed sediment. Thermocouple depthsvaried by location and ranged between 0.02 and 0.3 m depending on the easeof installation, which was determined by differences in streambed composi-tion and depth to bedrock. Vertical hydraulic gradients, bed temperatures,and hydraulic conductivity were measured once a month during winter sitevisits. Based on previous studies reported in the literature and at sites withinMalcolm Knapp Research Forest (Moore et al., 2005b; Guenther, 2007), it wasanticipated that bed heat conduction and hyporheic energy fluxes would besecondary terms in the energy budget and therefore more resources weredevoted to estimating the other fluxes. Spot measurements were performedto provide a basis for confirming that these fluxes were indeed secondaryterms.2.2.6 East and South creeksBoth East and South creeks were outfitted with v-notch weirs and stagelevel recorders to monitor streamflow. East and South creeks also had sub-mersible temperature sensor loggers installed at four locations extendingfrom the weir upstream at approximately 50 to 100 m intervals. The streamtemperature loggers recorded at 15 minute intervals. For the 2012–2013 win-ter, air temperature and relative humidity sensors (HOBO Pro v2 U23, OnsetComputer Corporation, Bourne, MA) were installed above both streams to23characterize the above-stream microclimate. In addition, time lapse camerawere installed at East and South creeks to determine periods of snow cover.24Chapter 3Diagnostic energy budget ofwinter stream temperatures inthe rain-on-snow zone3.1 IntroductionMost stream temperature research has been concerned with summer tem-peratures, particularly in response to riparian forest disturbance such asharvesting or wildfire (e.g. Brown, 1969; Johnson and Jones, 2000; Bartholow,2005; Gaffield et al., 2005; Gomi et al., 2006; Gravelle and Link, 2007; Leachand Moore, 2010; Janisch et al., 2012; Imholt et al., 2013). This summer-focused research consistently identifies energy exchanges occurring at thestream surface (primarily net radiation) as a key control on temperature vari-ability on diurnal to seasonal scales and the response to changes in ripariancanopy conditions. In some reaches, however, surface water-groundwaterinteractions can significantly moderate temperature variability (e.g. Johnson,2004; Moore et al., 2005b; Leach and Moore, 2011; MacDonald et al., 2014b).In contrast to the depth and breadth of research on summer streamtemperature, few studies have examined winter stream temperature pro-cesses despite its recognized importance for aquatic ecosystems (Beschta25et al., 1987; Holtby, 1988; Ebersole et al., 2006; Brown et al., 2011; Shuteret al., 2012). In many temperate regions, particularly the coastal portionsof the Pacific Northwest of North America (PNW), headwater catchmentsexperience moderate air temperatures and transient snow cover. As a result,headwater streams in these regions typically remain unfrozen during mostof the winter, and their thermal regimes may be particularly sensitive tochanges in winter air temperature and precipitation. Under these condi-tions, even a relatively small (e.g. 1–2 ◦C) but persistent change to streamtemperature could have a significant effect on fish bioenergetics and on ratesof growth and development of invertebrates that dominate salmonid foodsources (Brown et al., 2011; Arismendi et al., 2013).Winter in coastal portions of the PNW region is characterized by frequentcloud cover and precipitation. Therefore, solar radiation should be a lessimportant control on stream temperature than in summer, especially consid-ering the low solar elevation angles. Although there do not appear to be anypublished winter-time stream energy budgets in the coastal portion of thePNW, studies from sites outside this region have found that magnitudes ofenergy exchanges occurring at the stream surface are smaller during winterthan during summer periods (Webb and Zhang, 1997, 1999; Hannah et al.,2008; Leach and Moore, 2010; Garner et al., 2014). An important characteris-tic of the coastal portion of the PNW is that, at low to medium elevations,frequent rain and rain-on-snow events maintain high flows, which couldresult in substantial lateral advection of thermal energy via hillslope runoff.This study addressed two hypotheses: (1) winter stream temperaturesin headwater catchments of the rain and rain-on-snow zones of the PNWregion are primarily controlled by advective fluxes associated with runoffprocesses during storm events, and (2) stream temperatures should be de-pressed during rain-on-snow events, compared to rain-on-bare-ground, dueto the cooling effect of rain passing through the snowpack prior to infiltrat-ing the soil or being delivered to the stream as saturation-excess overlandflow. A diagnostic energy budget analysis was conducted using data col-lected during the winters of 2011–2012 and 2012–2013 from Griffith Creekand supplemented with thirteen years of historical stream temperature data26(1997–2002, 2007–2008, and 2010–2013) and modelled snowpack dynamicsfor East Creek to explore winter thermal regimes and the role of transientsnow cover. The methodology and results of the energy budget study arepresented in Section 3.2, followed by the methodology and results of thehistorical study in Section 3.3. The two complementary studies are discussedtogether in Section Energy budget study at Griffith CreekA reach-scale energy budget analysis was used to assess the relative impor-tance of the various energy fluxes controlling winter stream temperatures atGriffith Creek. The following sections outline the methods used to estimatethe various energy exchanges acting on the stream, followed by a descrip-tion of the reach-scale energy budget equations. Note that the term ‘surfaceenergy fluxes’ is defined as net radiation, sensible and latent heat fluxesoccurring at the stream surface; ’vertical energy fluxes’ is defined as surfacefluxes, bed heat conduction, and stream friction; and ’lateral energy fluxes’is defined as advective fluxes from overland flow and throughflow.3.2.1 AnalysisNet radiationNet radiation was modelled using hemispherical images that were analysedwith Gap Light Analyser (GLA) software (Frazer et al., 1999) to derivegap fractions as a function of zenith angle and azimuth. The model wasevaluated against measured net radiation from the roving net radiometer atGriffith Creek. The net radiation model, adapted from Moore et al. (2005b)and Leach and Moore (2010), computes net radiation, Q∗ (W m−2), as thesum of its short and long wave components:Q∗ = K∗ + L∗ (3.1)27where K∗ is the net short wave radiation at the stream surface and L∗ is thenet long wave radiation (W m−2).The net short wave radiation, K∗, can be expressed asK∗ = (1− α) [D(t)g(t) + S(t) fv] (3.2)where α is the albedo of the stream, D(t) is the direct component of incidentsolar radiation at time t (W m−2), S(t) is the diffuse component of incidentsolar radiation at time t (W m−2) and fv is the view factor. A value of 0.05was used for α (Oke, 1987).Gap fractions as a function of zenith angle (θ) and azimuth (ψ), g∗(θ,ψ),were derived by analysing the hemispherical images with Gap Light Anal-yser (GLA) software (Frazer et al., 1999) using 5◦ increments for both zenithand azimuth angles. The canopy gap at the location of the sun’s disk as afunction of time, g(t), was derived from g∗(θ,ψ) by computing the solarzenith and azimuth angles as a function of time, t, using equations in Iqbal(1983). The view factor was computed asfv =1pi∫ 2pi0∫ pi/20g∗(θ,ψ) cos θ sin θ · dθ · dψ (3.3)where θ is the solar zenith angle (vertical = 0), ψ is the azimuth angle, andg∗(θ,ψ) is the gap fraction at sky position θ, ψ. The double integral wasapproximated by summations using an interval of 5◦ for both zenith andazimuth angles.Solar radiation measured at the open site was used as input to Equation3.2. The measured solar radiation was partitioned into diffuse and directcomponents by calculating the diffuse fraction kd (ratio of the diffuse-to-global solar radiation) using the procedure presented by Erbs et al. (1982):kd =1.0− 0.09kt for kt ≤ 0.220.951− 0.1604kt + 4.388k2t−16.638k3t + 12.336k4t for 0.22 < kt ≤ 0.800.165 for kt > 0.80(3.4)28where kt is the clearness index (ratio of the global-to-extraterrestrial solarradiation).An optimum threshold value for converting the hemispherical photosto binary images was found by applying threshold values from 120 to 240,at 10 unit increments, to the hemispherical images associated with the twoabove-stream pyranometers. For each threshold value the incoming solarradiation was modelled using Equation 3.2 without the albedo component.A threshold value of 160 minimized the root mean square error and meanbias error for observed and modelled incoming shortwave radiation at thetwo pyranometer sites.The long wave radiation that is emitted by the atmosphere, topographyand vegetation, and which reaches the stream surface, was modelled asL ↓= [ fvεa + (1− fv)εvt]σ(Ta + 273.2)4 (3.5)where εa is the emissivity of the atmosphere, εvt is the emissivity of thevegetation and terrain, σ is the Stefan-Boltzmann constant, and Ta is the airtemperature measured at the nearest stream microclimate station (◦C). Avalue of 0.97 was used for εvt (Oke, 1987). The longwave radiation emittedby the water surface isL ↑= εwσ(Tw + 273.2)4 (3.6)where εw is the emissivity of the stream and Tw is the stream temperaturemeasured at the nearest water temperature logger site (◦C). A value of 0.95was used for εw (Oke, 1987). The net longwave radiation exchange is thenL∗ = εwL ↓ −L ↑ (3.7)Day time atmospheric emissivity was calculated by first using the Prata(1996) equation for clear sky conditions:ε0 = 1− (1 + w)exp[−(1.2 + 3.0w)0.5](3.8)29where ε0 is the clear sky atmospheric emissivity and w is the precipitablewater content of the atmosphere (cm), estimated asw = 465[ea/(Ta + 273.2)] (3.9)where ea is the atmospheric vapour pressure (kPa), and Ta is the air temper-ature (◦C), both measured at the nearest stream microclimate station to thesample site.The emissivity was adjusted for cloud cover using a procedure adaptedfrom Arnold et al. (1996) and Brock and Arnold (2000). The procedurerelates the cloud fraction, n, to the ratio K ↓ /K ↓max, where K ↓max is themaximum incoming short wave radiation under a cloud free sky, computedasK ↓max= I0Ψ(P/P0) cos θ (3.10)where I0 is the solar constant (1367 W m−2), Ψ is the clear sky transmissivity(0.75), P is atmospheric pressure calculated as 97 kPa for Griffith Creek,P0 is mean atmospheric pressure at sea level (101.3 kPa) and θ is the solarzenith angle. It is assumed that n = 1.0 for K ↓ /K ↓max ratios ≤ 0.2. ForK ↓ /K ↓max between 0.2 and 1.0, n is assumed to decrease linearly from 1.0to 0.0. The non-clear sky atmospheric emissivity, εa, can then be calculatedasεa = (1 + κn)ε0 (3.11)where κ is a constant depending on cloud type. Following Braithwaite andOlesen (1990), κ was set to 0.26, the mean value for altostratus, altocumulus,stratocumulus, stratus and cumulus cloud types. Equation 3.11 can only becalculated for daytime periods. Therefore, night time values of atmosphericemissivity were calculated as the mean of daytime values for n before andafter the night of interest.Mean hourly net radiation was modelled and evaluated against mea-sured net radiation at each of the fourteen roving net radiometer locations30using the associated hemispherical image for that location. For the energybudget analysis, hourly reach-scale mean net radiation was modelled forboth the forested and harvested reach using hemispherical images from therespective reaches.Sensible and latent heat fluxesThe latent heat flux, Qe (W m−2), was computed using an empirical windfunction fitted to evaporimeter data collected at Griffith Creek for both pre-and post-harvest conditions (Guenther et al., 2012):Qe = 627.8[0.0424 ·U · (ea − ew)] (3.12)where U is wind speed (m s−1) and ea and ew are the vapour pressures ofair and water, respectively (kPa), and 627.8 accounts for the unit conversionfrom mm h−1 to W m−2. Saturation vapour pressure (esat in kPa) was cal-culated as a function of air or water temperature, Ta or Tw (◦C), using thefollowing relation:esat(T) = 0.611 · exp[aTT + b](3.13)where T is in ◦C, and the coefficients a and b are given by (a, b)= (17.27,237.26) for T > 0 ◦C and (a, b)= (21.87, 265.5) for T ≤ 0 ◦C. The vapourpressure at the water surface, ew, was assumed to equal esat(Tw), while theactual vapour pressure of the air (ea) was calculated as:ea =(RH100)esat(Ta) (3.14)where RH is the relative humidity (%) measured at the nearest streammicroclimate station.The sensible heat flux (Qh) was computed as:Qh = β ·Qe (3.15)31where β is the Bowen ratio, calculated as:β = 0.66 · (P/1000) · [(Tw − Ta)/(ew − ea)] (3.16)where P is ambient air pressure, which was assigned the value for a standardatmosphere for the site elevation (97 kPa).Bed heat conduction and hyporheic energy exchangeSpot estimates of bed heat conduction, Qc (W m−2), were calculated as:Qc = Kc · (Tb − Tw)/0.05 (3.17)where Kc is the thermal conductivity of the streambed material (W m−1 K−1),and Tb and Tw are bed temperature at depth of 0.05 m and stream tempera-ture, respectively (◦C). The thermal conductivity was assumed to equal 2.6W m−1 K−1 based on estimates of Lapham (1989) using a porosity of 0.30,which is typical for gravels (Freeze and Cherry, 1979).The heat flux associated with hyporheic exchange, Qhyp (W m−2), for thecombined forest and harvest reaches was estimated as:Qhyp = ρwcwFhyp(Thyp − Tw)/W (3.18)where ρw is the density of water (1000 kg m−3), cw is the specific heat of water(4180 J kg−1 K−1), Fhyp is hyporheic exchange rate per unit length of channel(m3 s−1 m−1), Thyp and Tw are the temperatures of upwelling hyporheicwater and stream water, respectively, and W is the mean stream width (m).The mean of streambed temperatures from upwelling sites was used forThyp. Qhyp was not calculated separately for the forest and harvest reachesbecause of the limited number of piezometer and streambed temperaturesites available to compute Fhyp and Thyp.A reach-scale estimate of Fhyp was computed by assuming that all waterinfiltrating the bed within the reach follows subsurface flow paths thatdischarge within the same reach. It was further assumed that the fraction ofthe total bed area that experiences downwelling is equal to the fraction of32the piezometers with downwelling flow. Following from these assumptions,Fhyp was computed as:Fhyp = (ndw/npiezo) · A · qz/L (3.19)where ndw and npiezo are the number of piezometers indicating downwellingflow and the total number of piezometers, respectively, A is the area (m2)of the streambed, qz is the vertical flux of water infiltrating the streambed (m s−1), and L is the reach length (m). The area of the streambed wascomputed as the product of reach length and mean surface width, whichwas computed from the fitted relation with discharge.Infiltration rates were estimated based on Darcy’s Law:qz = Ksat · |∆h/∆z| (3.20)where Ksat is the saturated hydraulic conductivity of the bed (m s−1), deter-mined from slug tests, and ∆h/∆z is the mean vertical hydraulic gradientfrom piezometers that registered downwelling flow.Stream frictionHeat generated by fluid friction as water flows downstream, Q f (W m−2),was computed as:Q f = ρw · g ·Q · S/W (3.21)where g is the gravitational acceleration (9.81 m s−2), Q is the reach averagedischarge (m3 s−1), which is assumed to equal the mean of the upstreamand downstream discharges, S is the slope of the reach (m m−1) extractedfrom a 5 m resolution digital elevation model of the catchment, and W is theaverage wetted width of the stream reach (m).33Lateral heat fluxes calculated as residual of stream energy budgetEnergy budget analyses were conducted for the 2011–2012 and 2012-2013winter periods (1 October to 1 May) for the harvest reach and the forestreach. The longitudinal heat transfer (Ji) was calculated at ten minute timesteps at each gauging site as:Ji = ρw · cw ·Qi · Ti (3.22)where Qi is the discharge (m3 s−1), and Ti is the stream temperature, atgauging site i.The difference between Ji+1 and Ji is equal to the sum of the net heatinput to the reach and change in storage heat within the reach. This energybudget formulation does not require knowledge of travel times. The storageterm was found to be negligible and is not included in the calculation oflateral advective heat flux, described below. The energy exchanges at thestream surface (net radiation, sensible, and latent heat) and stream frictionwere estimated using the modelling approach described above and weresubtracted from the net heat input to the reach. The residual was attributedto the lateral advective heat input, Jadv (W):Jadv = Ji+1 − Ji − L ·W · (Q∗ + Qh + Qe + Q f ) (3.23)Effective lateral inflow temperatureThe effective lateral inflow temperature (Tadv) was calculated as:Tadv = Jadv/[(Qi+1 −Qi) · ρw · cw] (3.24)The calculated value of Tadv was compared to the 43 measurements ofnear-stream soil temperature. To explore how Tadv responded to differenthydroclimatic event types, daily means of Tadv were calculated and eachday during the study was classified as either (1) rain-on-ground, (2) rain-on-snow, (3) no precipitation and bare ground, or (4) no precipitation and snowcover.34Error analysisA standard approach for error analysis (Bevington and Robinson, 2003) wasapplied to the stream energy budget calculations in order to assess uncer-tainty. The propagated probable error was determined using the followinguncertainties: Q∗ ± 22 W m−2, Qe ± 25%, Qh ± 25%, Q± 5%, Tw ± 0.2 ◦C,L± 10%, S± 20%, W ± 21% (harvest reach) and ±19% (forest reach). Thevalues for Q∗, Q, W, and Tw were determined from uncertainty assessmentsof field measurements, whereas the remaining errors were deemed to bereasonable estimates that erred in the direction of overestimation.3.2.2 ResultsOverview of study periodFigure 3.1 places the 2011–2012 and 2012–2013 field seasons within a broaderclimatic context using long term (1962 to 2013) climate data from the MKRFheadquarters station. Monthly mean air temperatures for October to Aprilfor the two field seasons were generally within the middle 50% of historicvalues. There was more variability in precipitation patterns, as Octoberthrough December of 2011–2012 was dry relative to historic conditions,whereas precipitation during the 2012–2013 season was generally similar toor well above the long term median.Griffith Creek was a net gaining stream during the study. The net gain instreamflow was greater for the forest reach (Q3 to Q2) than for the harvestreach (Q2 to Q1), consistent with the differences in drainage area (Figures 3.2and 3.3). The 2012–2013 winter experienced higher peak flows than 2011–2012 and discharge was generally more responsive to rain events in 2012–2013, particularly for the Q3 site. This difference in discharge responsebetween the two years was similar for other headwater catchments in theresearch forest (data not shown).Mean hourly stream temperatures varied between 1 and 12 ◦C duringthe study period (Figures 3.2 and 3.3). Mean hourly longitudinal streamtemperature differences between Q1, Q2, and Q3 were often less than 0.5 ◦C.35Figure 3.1: Boxplots of monthly (October to April) mean air tempera-ture and monthly total precipitation recorded at the MKRF head-quarters station from 1962 to 2010. Red circles and blue squaresrepresent conditions during the two winter periods (2011–2012and 2012-2013, respectively) during which the detailed field studyat Griffith Creek was conducted.However, in April, temperatures at Q1 tended to be about 1 ◦C greater thanat Q2 and Q3 during the day. Maximum diurnal range in hourly tempera-tures at Q1 were 3.2 and 3.4 ◦C for 2011–2012 and 2012–2013, respectively.Mean diurnal range in hourly temperatures was 0.9 ◦C for both winters.Stream temperature did not have a consistently positive or negative rela-tion with discharge, but seemed to vary with discharge depending on theantecedent air and stream temperatures. Mean hourly air temperatures fellbelow 0 ◦C for only short periods (less than one week) during winter.36Figure 3.2: Hydroclimatic overview for the 2011–2012 winter period(October to May). Water input refers to rainfall and snowmeltmeasured by the snow lysimeter located at the open site.37Figure 3.3: Hydroclimatic overview for the 2012–2013 winter period(October to May). Water input refers to rainfall and snowmeltmeasured by the snow lysimeter located at the open site.38Total water input (rainfall plus snowmelt) measured at the open sitelysimeter for the study periods (1 October to 1 May) was 1736 mm for 2011–2012 and 1950 mm for 2012–2013. The 2011–2012 winter had 52 days ofsnow cover in the forested area and 63 days in the harvested area, spreadover a number of events, whereas the 2012–2013 winter had 71 days of snowcover in both forested and harvested areas, mostly due to a snow pack thatpersisted from 14 December to 25 January.Energy balanceHourly modelled net radiation was compared to measured net radiationat the fourteen net radiometer locations during the study period. The rootmean square error and mean bias error ranged between 17 and 22, and -18and -6 W m−2, respectively, across the fourteen sites. These errors are similarto those found in previous efforts to model above-stream net radiation fromhemispherical photographs (Leach and Moore, 2010; Bulliner and Hubbart,2013).Streambed hydraulic conductivity based on multiple slug tests acrosssites had a mean of 4.2 ×10−5 m s−1 and vertical hydraulic gradients rangedbetween -0.6 and 0.5 cm cm−1. Guenther et al. (2014) measured similardownward gradients during summer, but their upward gradients weregenerally weaker than the winter measurements shown here. Table 3.1shows the number of piezometers with downwelling flow, total number ofpiezometers sampled, vertical hydraulic gradients calculated for piezome-ters with downwelling flow, reach-average differences between Thyp andTw for the upwelling locations and Qhyp estimates, for six dates when allpiezometer sites were sampled. For flows below about 10 L s−1, up to al-most 80% of the piezometers registered downwelling flow, with the fractiondropping to just over 50% at higher flows. Temperature differences betweenThyp and Tw were less than 0.7 ◦C and estimates of the reach average heattransfer associated with hyporheic exchange were less than 50 W m−2. Boththe fraction of piezometers registering downwelling flow and the differencein hyporheic and stream temperatures decreased with increasing discharge,39Table 3.1: Summary of discharge (Q), number of piezometers indicatingdownwelling flow (ndw), total number of piezometers measured(ntotal), reach average vertical hydraulic gradient for ndw (VGHdw),mean temperature difference between Thyp− Tw for upwelling sites,and estimated reach average Qhyp. Note Thyp was not measuredfor the first three dates.Date Q ndw ntotal VGHdw Mean Thyp − Tw Qhyp(L s−1) (cm cm−1) (◦C) (W m−2)30 November 2011 10.4 19 25 −0.26 NA NA6 December 2011 3.5 15 24 −0.18 NA NA1 February 2012 38.1 13 25 −0.36 NA NA9 November 2012 7.1 18 25 −0.56 0.6 40.410 January 2013 24.2 13 25 −0.36 0.3 9.43 March 2013 34.5 13 25 −0.37 0.2 6.4suggesting that the energy transfer associated with hyporheic exchange alsodeclines with increasing discharge. Although separate values of Qhyp for theforest and harvest reaches were not determined due to the small number ofsamples, the limited data suggest that there was not a substantial differencebetween the reaches in terms of Thyp and percentage of piezometers withdownwelling flow. For bed heat conduction calculations, temperature gradi-ents between the stream and bed were small during all spot measurementsand across all locations. Temperature differences between the stream andbed at depths of 15 cm were less than 0.8 ◦C and often less than 0.2 ◦C.Estimates of bed heat conduction were less than 15 W m−2.Figures 3.4 and 3.5 show estimated hourly surface energy fluxes and theheat inputs due to stream friction for the harvest and forest reaches for 2011–2012 and 2012–2013. Net radiation was the dominant surface energy flux.Hourly net radiation ranged between -70 and 150 W m−2 for the harvestreach and -60 and 120 W m−2 for the forest reach, whereas the sensible andlatent heat exchanges at both reaches never exceeded ±30 W m−2. Duringthe months of November to February for both 2011–2012 and 2012–2013,net radiation was mostly an energy loss from the stream. Starting in lateFebruary, net radiation was an energy source to the stream during day anda sink during night. Heat generated by frictional dissipation of potential40energy averaged between 15 and 30 W m−2 over the study periods andranged between 0.2 and 219 W m−2 for the harvest reach and 4.8 and 218W m−2 for the forest reach. Heat inputs due to friction were greatest duringhigh flows.Figures 3.4 and 3.5 summarize the reach-scale energy budget analysisfor the 2011–2012 and 2012–2013 winter periods at the harvest and forestreaches. Except during low flow periods, when lateral inputs of waterbecame relatively small, the vertical energy exchanges were only a smallfraction of the total heat input into either the harvest or forest reaches evenwith accounting for uncertainty in the flux estimates. Figure 3.6 illustratesthat for discharges above approximately 25 L s−1, surface energy fluxesaccount for less than 2% of the heat input to the reach. For discharges below25 L s−1, surface energy fluxes can account for a considerable proportionof the net reach heat input, both as a heat source and sink, as indicated bypositive or negative ratios, respectively, in Figure 3.6. Most of the highlynegative ratios correspond with periods of large relative errors in the heatbudget calculations, particularly due to the small changes in dischargesalong the reach. Surface fluxes were, on average, 4.6% and 4.5% of the netheat input to the forest reach during the 2011–2012 and 2012–2013 winterperiods, respectively, and slightly greater at 9.3% and 5.9% during 2011–2012and 2012–2013, respectively, for the harvest reach. However, during somelow flow periods, when advective inputs and heat generated by frictionwere minimal, the surface fluxes accounted for nearly the entire change inlongitudinal heat flux along the reaches. The pattern of net heat input intothe reaches followed that of discharge, which highlights lateral advection asa dominant heat input to the stream reach. The generally higher magnitudeof net heat inputs for 2012–2013 vs. 2011–2012 reflects the higher magnitudeof streamflow that year. Relative magnitudes of the surface energy inputsbetween the forest and harvest reaches were similar between the two years.41Figure 3.4: Time series of hourly energy budget components (net ra-diation, Q∗; latent heat flux, Qe; sensible heat flux, Qh; and heataddition due to friction, Q f ) for the forest and harvest reaches for2011–2012. Note the difference in scales of the total heat exchangeplots and the vertical heat flux plots.42Figure 3.5: Time series of hourly energy budget components (net ra-diation, Q∗; latent heat flux, Qe; sensible heat flux, Qh; and heataddition due to friction, Q f ) for the forest and harvest reaches for2012–2013. Note the difference in scales of the total heat exchangeplots and the vertical heat flux plots.43Figure 3.6: Ratios of the surface energy fluxes to the net heat inputs forthe forest and harvest reaches plotted against hourly dischargeat the lower boundary of the respective reaches. Red circles in-dicate ratios with relative errors less than or equal to 50%, andblue circles indicate ratios with relative errors greater than 50%.Plots include both 2011–2012 and 2012-2013 data. Positive ratiosindicate that the surface fluxes are adding heat to the reach andnegative ratios indicate they are removing heat from the reach.Note the x-axis scale is logarithmic.Effective lateral inflow temperatureFigures 3.7 and 3.8 show the observed near-stream subsurface temperatures(n = 33 sites for 2011–2012 and n = 43 for 2012–2013) and the proba-ble range of back-calculated effective lateral inflow temperatures (takinginto account the propagated uncertainties), as well as hourly groundwatertable levels. The water table levels exhibited a flashy response to stormevents. The observed subsurface temperatures are temporally variable, dif-fer spatially by up to 5 ◦C and plot within the probable range of effectivelateral inflow temperatures. The few periods when observed subsurface44and estimated lateral inflow temperatures do not agree (e.g. forest reach inmid-January 2012) correspond to periods of low streamflow when water andheat contributions from the hillslope would be relatively small and verticalheat exchanges would be more important. The effective lateral inflow tem-perature has a wider uncertainty range for the harvest reach compared tothe forest reach because of smaller differences in upstream and downstreamdischarge, which produces higher relative uncertainty in the estimates ofnet water input to the reach.Figure 3.9 show daily mean effective lateral inflow temperatures anddaily mean air temperatures for both winter periods classified for four eventtypes: (1) rain-on-ground, (2) no precipitation and bare ground, (3) no pre-cipitation and snow cover, and (4) rain-on-snow. The scatter plots illustratethat for air temperatures in the range of 0–5 ◦C, lateral inflow temperaturesduring rain-on-snow events are dominantly below 5 ◦C in the forest reachand 6 ◦C in the harvest reach, whereas they commonly range up to 8 ◦C forrain-on-ground events in both reaches. Observed subsurface temperaturesfor loggers that were below the water table showed a similar relationship toair temperature when classified into the four event types (data not shown).3.3 Historical study at East CreekThe relationship between hydrometeorologic event type and winter streamtemperature at East Creek was examined using thirteen years of streamtemperature data. The historical study complements the energy budgetstudy by addressing the role of snow cover on stream temperature duringrain events by extending the limited period of study of the detailed energybudget analysis and sampling more precipitation events over variable snowand weather conditions. The following sections outline the data collectedas well as details of a snow accumulation and melt model used to simulatehistoric snow cover at East Creek. The model output was used, in con-junction with air temperature and precipitation data, to classify days intofive categories: (1) rain-on-ground, (2) bare ground and no precipitation,(3) snowing, (4) snow-on-ground and no precipitation, and (5) rain-on-snow.45Figure 3.7: Comparison of back-calculated effective advective tempera-ture and observed near stream subsurface temperatures, as well asgroundwater table observations for the 2011–2012 winter period.The black lines are the calculated effective advective temperatureand the gray bands represent the probable error range. Duringperiods when the subsurface temperature logger was below thewater table the observed temperature line is blue, and red forperiods when the logger was above the water table. The ground-water level plots show the water table elevation below the groundsurface (dotted line) measured at the individual groundwaterwells.46Figure 3.8: Comparison of back-calculated effective advective tempera-ture and observed near stream subsurface temperatures, as well asgroundwater table observations for the 2012–2013 winter period.The black lines are the calculated effective advective temperatureand the gray bands represent the probable error range. Duringperiods when the subsurface temperature logger was below thewater table the observed temperature line is blue, and red forperiods when the logger was above the water table. The ground-water level plots show the water table elevation below the groundsurface (dotted line) measured at the individual groundwaterwells.47Figure 3.9: Calculated mean daily effective advective temperature forfour event types for the forest and harvest reaches during winter2011–2012 and 2012–2013.3.3.1 Data sourcesContinuous stream temperature measurements were made at the East Creekcatchment outlet from 1997 to 2002, 2007 to 2008, and 2010 to 2013, usingsubmersible Onset Tidbit temperature loggers (accurate to ±0.2 ◦C). Theloggers were programmed to record and store temperature readings every20 min. Loggers were fitted with radiation shields made from white PVCpipe with holes drilled through the pipe to facilitate water exchange. Sen-sors were calibrated at 0 and 20 ◦C (in an ice bath and at room temperature,respectively) before and after field deployment. Daily mean winter (Novem-ber to March) stream temperatures were calculated from the continuousrecords. Daily air temperature, rainfall, snowfall, and presence/absence ofsnow were recorded at the Research Forest headquarters weather station.In addition, a time lapse camera was installed near the East Creek outletfrom 2011 to 2013 and was used to determine presence of snow cover on thecatchment.483.3.2 Modelling snow cover at East CreekDue to the elevation difference between the headquarters climate stationand East Creek, as well as data gaps in the snow record at the headquartersclimate station, a simple temperature-based snow accumulation and meltmodel was used to simulate snow cover at East Creek. A temperature-basedsnow model was chosen for two reasons. First, it has a minimal number ofmodel parameters that need to be calibrated compared to more physicallybased models. Second, East Creek has dense canopy cover resulting inlittle solar radiation reaching the ground surface and also low wind speeds;therefore, snow melt is primarily driven by longwave radiation from theforest canopy, which is strongly related to air temperature (Black et al., 1991;Pomeroy et al., 2009). For simplicity, the model did not account for coldcontent and changes in water content. These omissions should not be crit-ical at this site because transient snowpacks are typically wet, isothermal,and shallow. Furthermore, the primary interest is in predicting the pres-ence/absence of snow on the ground, and not the precise timing of snowmelt runoff.The snow model was used to estimate snow cover at East Creek from1962 to 2013. An elevation of 352 m was used in the model, which representsthe mean elevation of the catchment. A threshold temperature (TT) wasused to partition precipitation into rain (temperature > TT) and snow(temperature < TT). For precipitation falling as snow, accumulated SWE(mm) was calculated as:SWEt = SWEt−1 + SFt (3.25)where SF is the amount of snowfall (mm) and t is the day.For days when the lapsed air temperature was greater than the thresholdtemperature and snow was present, snowmelt (SM, mm) was calculated as:SMt = K f · Tt (3.26)49where K f is the snowmelt factor (mm ◦C−1) and T is the air temperature(◦C). The resulting SWE (mm) was calculated as:SWEt = SWEt−1 − SMt (3.27)The snow model used elevation-adjusted air temperature and precipitationfrom the research forest headquarters climate station. The temperature lapserate was assumed to be 0.006 ◦C m−1 (Oke, 1987) and precipitation wasassumed to be 1.16 times the precipitation measured at the Research Forestheadquarters station based on field measurements made by Hetherington(1976) and Donnelly-Makowecki and Moore (1999).Values for TT and K f were calibrated by comparing modelled SWEwith manual snow surveys made in the forested portion of Griffith Creek’scatchment (430 m elevation) and also by comparing predicted snow presencewith records of snow cover determined from time lapse photography at theEast Creek weir (291 m elevation). The comparison resulted in TT = 0 ◦Cand K f = 1.5, which are similar to values determined for forest sites in thesouthern Coast Mountains of British Columbia in other modelling studies(e.g. Moore, 1993). Running the model using these calibrated values for theelevation of the Research Forest headquarters produced modelled snowfallthat generally matched the timing of snowfall recorded at that station.Figure 3.10 shows modelled SWE for 2011 to 2013, observed periods ofsnow cover, and mean snowpack water equivalent with standard errors fortwo snow surveys, for both the forested portion of Griffith Creek and thearea near the East Creek outlet. At both elevations (430 and 291 m), the snowmodel captures the magnitude and duration of the larger, more persistentsnow packs, but tends not to capture the short (< 1 day) duration snowpacks.Figure 3.11 shows modelled SWE for 1962 to 2013. The years used in thisstudy cover a range of snow accumulation and melt patterns, but generallyrepresent low to moderate accumulations compared to the longer referenceperiod.50Figure 3.10: Comparison of modelled and observed snow cover, andmeasured snow water equivalence from the forested area of Grif-fith Creek (430 m elevation), and modelled SWE and observedsnow cover at the East Creek outlet (291 m elevation), from Octo-ber to April for 2011–2012 and 2012–2013. Red horizontal linesindicate periods of snow cover determined from the time lapsecameras. The black points are mean snow water equivalent (withstandard error) from manual snow surveys conducted at theforested site at Griffith Creek.51Figure 3.11: Comparison of modelled historic snowpack water equiva-lent for winters of 1962 to 2013. Black lines represent years usedin the historic East Creek study and grey lines represent yearsfor which stream temperature data are not available.3.3.3 Stream temperature response to hydroclimatic event typeFor the historic winter records (1 November to 30 March for each winter),562 days were classified as days with bare ground and no precipitation,1230 as days with rain-on-ground, 159 as snowing, 228 as no precipitationbut snow on the ground, and 241 as rain-on-snow events. Figure 3.12shows the relation between daily winter stream temperatures and dailywinter air temperatures. On days when a rain-on-snow event occurred,stream temperatures remained below approximately 5 ◦C, regardless of thecorresponding air temperature. In comparison, for the same air temperaturerange during which rain-on-snow events occurred, stream temperatures52Figure 3.12: Daily mean stream temperature at East Creek as a functionof air temperature, classified into five event types: (1) rain-on-bare-ground, (2) no precipitation with bare ground, (3) snowing,(4) no precipitation with snow cover, and (5) rain-on-snow.were as high as 8 ◦C when no snowpack was present on the catchment.Figure 3.13 shows daily winter stream temperature as a function of eventtype for periods when mean daily air temperature was between 0 and 5 ◦C.This figure highlights that, for similar mean daily air temperatures, meandaily stream temperatures are generally lower when snow is present on thecatchment compared to days when snow is absent. Analysis of varianceshows that daily mean stream temperatures during rain-on-snow eventsare about 1.4 ◦C lower than during rain-on-ground events (3.2 versus 4.6◦C, p < 0.001). However, the residuals of the analysis of variance exhibitautocorrelation (lag-one autocorrelation coefficient = 0.81), which meansthat the significance levels may be overstated.In order to address the issue of autocorrelation, a Monte Carlo test of53Figure 3.13: Box plots of mean daily stream temperature at East Creekduring periods of mean daily air temperature between 0 and 5◦C for four different event types: (1) no precipitation with bareground, (2) rain-on-bare-ground (ROG), (3) rain-on-snow (ROS),and (4) no precipitation with snow cover.significance was conducted to test the null hypothesis that, for air tempera-tures between 0 and 5 ◦C, there is no difference in mean stream temperaturefor rain-on-snow and rain-on-bare-ground events. The Monte Carlo ap-proach involved the generation of stochastic realizations of snow cover forthe thirteen winters of historic stream temperature data. The air tempera-ture, precipitation and stream temperature data were kept fixed in order torespect their correlation structures. To account for the inherent autocorrela-tion of snow cover presence/absence, a transition probability matrix wasdeveloped that included the probability of day i having snow dependingon whether there was snow cover on day i− 1. These probabilities werecomputed from the modelled snow cover for the thirteen winters of stream54temperature data.Five thousand realizations of snow cover were generated for the historicstudy period using the transition probability matrix. Days with rain and airtemperatures between 0 and 5 ◦C were extracted to calculate mean streamtemperature for rain-on-snow and rain-on-bare-ground conditions. Themean difference in stream temperature between rain-on-snow and rain-on-bare-ground days for the 5000 realizations had a mean of 0.1 ◦C anda standard deviation of 0.3 ◦C. The mean stream temperature differencebetween the observed rain-on-snow and rain-on-bare-ground days was -1.4◦C. The observed mean difference fell outside the range of differences basedon 5000 realizations, suggesting a p value of < Discussion3.4.1 Relative roles of vertical and lateral heat fluxesSurface heat fluxes estimated for winter at Griffith Creek were of similarmagnitude to those estimated for other studies conducted during winterperiods (Webb and Zhang, 1997, 1999; Hannah et al., 2004, 2008). Netradiation was generally greater in magnitude than the sensible and latentheat fluxes, but was typically less than 100 W m−2, representing a heat lossthrough most of the winter and switching to a heat source in February. Themagnitude of the ratio of surface heat fluxes and stream friction to the totalheat input to the reach decreased with discharge and became negligible fordischarges greater than about 25 L s−1.Spot estimates of bed heat conduction were less than 12 W m−2 and aresimilar in magnitude and general direction as found in previous studiesthat have estimated this term for winter periods (Webb and Zhang, 1999;Hannah et al., 2008). However, Webb and Zhang (1999) did find that bed heatconduction was a large energy sink during the day (as large as -170 W m−2over an averaging period of ten minutes), which was related to diurnalpatterns in stream temperature. Griffith Creek did not exhibit strong diurnalwarming and, therefore, energy exchanges of the magnitudes estimated by55Webb and Zhang (1999) likely did not occur.Although previous studies have quantified the energy flux associatedwith hyporheic exchange across the stream bed during summer (Story et al.,2003; Moore et al., 2005b; Hester et al., 2009; Neilson et al., 2009), no previousstudies have estimated this flux for winter periods. Spot estimates of reachaverage hyporheic energy exchange suggest that the magnitude of thisterm could be up to around 50 W m−2 for lower flows, which makes itcomparable in magnitude to net radiation. Moore et al. (2005b) estimatedthe hyporheic exchange flux to be up to about 120 W m−2 in absolute valueat A Creek in the Malcolm Knapp Research Forest. However, that estimatewas for low flows on a clear summer day within a clearcut, which wouldpromote a greater contrast in temperature between the hyporheic zone andstreamwater, and thus a greater hyporheic exchange flux. Guenther (2007)estimated reach average heat fluxes associated with hyporheic exchangeduring summer both before and after harvesting at Griffith Creek. Underpre-harvest conditions, Qhyp was a heat sink during the day and rangeddown to about -20 W m−2. Following harvesting, the greater contrast intemperature between the hyporheic zone and the stream resulted in Qhypdown to about -100 W m−2. Considering these summer-time estimates andthe limited sampling under winter conditions, it appears reasonable toconclude that hyporheic exchange and bed heat conduction are of similarmagnitude to the surface energy exchanges during winter and thus couldbe important heat fluxes at lower flows, but are negligible at higher flows,despite some recent research that shows that transient storage can be animportant process during high flow periods (Ward et al., 2013). Observedtemperature differences between the stream water and bed were often lessthan 0.3 ◦C, which limits the magnitude of the hyporheic energy flux, despiteuncertainties in the hyporheic water flux. Further research is required toclarify the role of hyporheic heat exchange and its influence on streamtemperature dynamics.Heat generation due to friction is often not included in stream energybudget calculations (Johnson, 2004; Hebert et al., 2011; Benyahya et al., 2012)or is found to be negligible for summer low flow periods (Webb and Zhang,561999; Moore et al., 2005b). Winter energy budget studies from the UnitedKingdom report maximum heat inputs due to friction during winter tobe less than 4 W m−2 for reaches with slopes of 0.0008 and 0.0028 m m−1(Webb and Zhang, 1999) and up to 60 W m−2 for a reach with a slope of0.01 m m−1 (Hannah et al., 2004). Heat generated due to friction for GriffithCreek’s relatively steep reaches (slopes of 0.09 m m−1 and 0.16 m m−1 forthe harvest and forest reach, respectively) contributed on average 15-30W m−2 but ranged up to 200 W m−2 during high flow events. Therefore,heat addition due to friction is of the same magnitude as the surface energyfluxes, but is small compared to the lateral advective flux at higher flows.The reach-scale energy budget analysis highlights that advective fluxesassociated with runoff processes dominate the winter thermal regime forcoastal headwater catchments except during recession periods betweenrain and rain-on-snow events. This result is contrary to a number of pre-vious studies that found that, although stream surface energy exchangesare relatively small in magnitude during winter compared to summer, theyremained the primary control on stream temperature dynamics (Webb andZhang, 1997, 1999; Hannah et al., 2008). However, those studies were con-ducted at point scales or over short reaches (less than 40 m length) and didnot focus on the role of advective fluxes associated with runoff contributionsin a headwater context. The potential importance of advective fluxes forwinter stream temperature has been recognized in previous studies (Webband Zhang, 1999; Malard et al., 2001; Leach, 2008; Danehy et al., 2010); how-ever, this study provides the first energy budget analysis to quantify itsimportance under winter conditions in a headwater catchment. Findingsfrom this study provide process-based support for recent empirical researchsuggesting that winter and summer periods of the PNW are characterizedby distinct thermal regimes and sensitivities to environmental change (Aris-mendi et al., 2013).The significance of lateral advective heat inputs implies that accuratesimulation of runoff generation and the temperature of lateral inflow isrequired to model stream temperature dynamics in regions with frequentwinter rain and rain-on-snow events. Most stream temperature models57have either ignored lateral advection (e.g. Sinokrot and Stefan, 1993) orassumed a spatially uniform subsurface temperature that is a functionof seasonal or annual mean air temperature (e.g. St-Hilaire et al., 2000;Ficklin et al., 2012). Field observations revealed considerable variability insubsurface temperatures both through time and among loggers. Thirty-sixof 43 subsurface temperature sites had a single logger which was installedat the soil-till interface, where lateral flow generation in this system isinitiated (Hutchinson and Moore, 2000). However, analysis of the subsurfacetemperatures at sites with temperature records at multiple depths (Chapter4) suggests that during wet periods when the water table is elevated, verticaltemperature differences within the subsurface are less than 0.3 ◦C.3.4.2 Role of transient snow cover and rain-on-snow eventsTo the author’s knowledge, no previous research has examined streamtemperature response to rain-on-snow events with the exception of Lan-gan et al. (2001), who attributed long term increases in winter maximumdaily temperatures of 2 ◦C over 30 yr at a Scottish catchment to reductionsin snow cover. Findings from both Griffith and East creeks indicate thatstream temperatures are depressed during rain-on-snow events, comparedto rain-on-ground events, by about 1-2 ◦C for a given air temperature. Themagnitude of this temperature depression appeared to be relatively con-sistent among years, between the forest and logged sites at Griffith Creek,and between Griffith and East creeks. It is hypothesized that this temper-ature depression is due to the cooling of rain as it percolates through thesnowpack prior to infiltrating the soil or being delivered to the stream assaturation-excess overland flow.Inferences regarding the historic influence of transient snow cover atEast Creek depend on the accuracy of the snow model for classifying eventtypes. The snow model appears to capture the larger snow events, whichwill have the most pronounced impact on the thermal regime, but appears tomiss the short duration (less than a day) and shallow snowpacks that formsporadically during the winter at these sites. One particular source of error58is that the model did not account for canopy effects, which would causea tendency for the model to overestimate snow accumulation, particularlyduring snow events that contain periods with temperatures slightly abovethe freezing point, during which intercepted snow tends to melt and reachthe ground as canopy drip (Berris and Harr, 1987; Storck et al., 2002). Thenet effect of model errors would be to classify days with snow cover as dayswithout snow cover, and vice versa, which would act to reduce the apparentstream temperature difference between event types rather than introducea spurious effect. Classification of days into different event types is alsoa source of uncertainty, since precipitation can shift phases during a day.However, these phase shifts usually occur over short periods and the longerevent dynamics are likely captured by the modelling approach used here.Indeed, the fact that a clear and consistent signal emerged from the analysesover multiple sites and years supports the robustness of the results, despitethe presence of some misclassification of event types.59Chapter 4Hillslope throughflowtemperature observations andmodel comparison4.1 IntroductionThe dominant role of energy exchanges at the stream surface on streamthermal regimes, primarily net radiation, has been recognized since thepioneering heat budget analysis by Brown (1969). Most predictive streamtemperature models reflect this understanding and their focus has beenon accurately representing the energy fluxes occurring at the stream-airinterface (e.g. Theurer et al. 1984; Boyd and Casper 2003). This focus maybe adequate for wide streams and periods when net radiation and turbulentheat fluxes dominate the thermal regime. However, a growing body ofresearch has highlighted that advection associated with hyporheic exchange,groundwater discharge and hillslope runoff can strongly influence streamtemperature dynamics (Webb and Zhang, 1999; Malard et al., 2001; Storyet al., 2003; Brown and Hannah, 2007; Danehy et al., 2010; Leach and Moore,2011; MacDonald et al., 2014b). Considerable effort has been made to charac-terize the thermal role of hyporheic exchange (White et al., 1987; Story et al.,602003; Moore et al., 2005b; Arrigoni et al., 2008; Hester et al., 2009; Marzadriet al., 2013). A number of modelling studies have simulated heat advectionby hillslope runoff (e.g. St-Hilaire et al. 2000; Ficklin et al. 2012; MacDonaldet al. 2014a) or explicitly modelled temperatures of shallow groundwaterdischarge (Kurylyk et al., 2014); however, none of these studies evaluatedtheir predictions of runoff and shallow groundwater temperatures againstobserved data. To the author’s knowledge, no studies have documentedthe spatial and temporal variability in hillslope throughflow temperatures.Given that throughflow is a dominant runoff process in permeable catch-ments with shallow soils (Anderson and Burt, 1978; Kirkby, 1988; McGlynnet al., 2002), an understanding of the variability in throughflow temperaturesis needed for accurate prediction of stream temperature, particularly duringwinter rain and rain-on-snow events, as seen in the previous chapter.A number of coupled hydrologic and stream temperature models simu-late throughflow and groundwater processes and their influences on streamtemperature (Morin and Couillard, 1990; Bicknell et al., 2001; Boyd andCasper, 2003; Allen et al., 2007; Haag and Luce, 2008; Ficklin et al., 2012;Loinaz et al., 2013; Null et al., 2013; MacDonald et al., 2014a; Sun et al., 2014).The majority represent groundwater temperature by adding 0 to 2 ◦C tomean annual air temperature (e.g. Allen et al. 2007; Ficklin et al. 2012; Nullet al. 2013; Sun et al. 2014; MacDonald et al. 2014a) or require a user specifiedvalue (e.g. Boyd and Casper 2003). Although this term is generally meant torepresent groundwater originating from deep flow pathways, some of themodels also use it to represent shallow groundwater and throughflow in-puts (e.g. Sun et al. 2014). A recent update to the Soil and Water AssessmentTool (SWAT) stream temperature routine (Ficklin et al., 2012) represents thetemperature of lateral inflows using a lagged mean daily air temperaturevalue, where the length of the lag is determined through calibration. Thisversion of SWAT also assumes that water inputs from snowmelt enter thestream at 0.1◦C. The CEQUEAU model (Morin and Couillard, 1990) includesan interflow advection term. The model estimates the interflow tempera-ture as (Tg + Ta)/2, where Tg is groundwater temperature (◦C), assumedconstant throughout the year, and Ta is air temperature (◦C). An updated61interflow temperature equation in CEQUEAU introduced by St-Hilaire et al.(2000) is also primarily driven by air temperature and uses a force-restoremethod (Deardorff, 1978). The Hydrological Simulation Program-Fortran(HSPF) (Bicknell et al., 2001) includes an interflow advection term that iscalculated using a calibrated relationship with air temperature and interflowtemperature from the previous day. The LARSIM model (Bremicker, 2000)was extended to simulate stream temperature (LARSIM-WT) and includesadvection associated with runoff processes (Haag and Luce, 2008). Directrunoff, interflow, and groundwater temperatures are estimated as linearfunctions of air temperature. The MIKE SHE model explicitly simulatesgroundwater flow, but assumes a spatially uniform groundwater tempera-ture for inflows to the stream (Loinaz et al., 2013). MIKE SHE also includes adrainage heat flux, which represents overland flow (Loinaz et al., 2013, 2014).This drainage heat flux uses the equilibrium temperature concept (Boganet al., 2003) to represent the temperature of this input. In contrast to theprevious models, Tung et al. (2014) used measurements of soil temperatureat an unspecified depth to represent the temperature of subsurface waterinputs. They found that representing the temperature of subsurface waterinputs with soil temperature provided better model performance than usinglocal air temperature.A shortcoming of these existing modelling approaches is that, otherthan SWAT, they do not explicitly account for the role of snowmelt or rain-on-snow events, which can have a detectable influence on stream thermalregimes, as seen in the previous chapter. Despite these efforts to simulateshallow groundwater and throughflow temperatures, none of these modelshas compared modelled throughflow and groundwater temperatures to fieldobservations. An evaluation of these simulated temperatures is importantbecause correct representation of the key energy exchange controls on streamtemperature is required for prediction, especially when predicting howfuture climate and land cover changes may alter the magnitude and timingof throughflow to the stream and its subsequent heat advection.In contrast to the relatively simple representation of throughflow advec-tive heat flux in current stream temperature models, a number of models62exist that are sufficiently complex to simulate shallow groundwater flowand heat transport and can represent the influence of snowmelt, includingSHAW (Flerchinger, 2000), COUP (Jansson and Karlberg, 2001), ForHyM2(Yin and Arp, 1993), SUTRA (Voss and Provost, 2002) and HydroGeoSphere(Therrien et al., 2010). SHAW, COUP, and ForHyM2 primarily solve soilthermodynamics at the point scale and cannot be applied to systems wherelateral advection associated with throughflow is an important hydrologicprocess. SUTRA and HydroGeoSphere can be employed for two- and three-dimensional problems and include representation of throughflow. However,these models require substantial computational resources and highly skilledoperators to set up for an individual site. For example, the author set up SU-TRA for a small hillslope in Griffith Creek; because the shallow groundwatersystem is highly dynamic, short time steps and a fine finite-element meshsize were required for stable solutions of the mass and energy exchangeequations. Given these requirements, a single model run to simulate shallowgroundwater and heat flow for a six month period took approximately aweek to complete on a reasonably fast desktop computer processor (IntelCore i7-3930K 3.20GHz CPU).Understanding spatiotemporal variability in throughflow temperaturesis not only helpful to evaluate current catchment scale stream temperaturemodels, but can also be used to address the utility of using heat as a hy-drologic tracer. There is a growing use of heat as a tracer to understandinteractions between surface water and groundwater (Anderson, 2005; Con-stantz, 2008). The usefulness of these approaches relies on sufficiently largetemperature differences between water end-members and is most success-fully applied to systems with relatively stable groundwater temperatures.A few studies have attempted to use variations in shallow subsurface hills-lope water temperatures as tracers to separate event hydrographs and inferprimary runoff flow pathways (Shanley and Peters, 1988; Kobayashi et al.,1999; Birkinshaw and Webb, 2010). In the case of Kobayashi et al. (1999)and Birkinshaw and Webb (2010), subsurface runoff temperatures weremade at multiple depths at one location in the study catchment and theirfindings rely on the assumption that pre-event and event water are spatially63invariant. An evaluation of the assumption that pre-event and event watertemperatures are uniform in space would help to determine the usefulnessof heat as a tracer for hillslope runoff processes.In summary, no previous research has examined the spatiotemporal vari-ability in throughflow temperatures across a catchment. Understanding thisvariability is valuable for developing accurate stream temperature models,as well as evaluating the usefulness of heat as a hillslope runoff tracer. Vari-ous coupled stream temperature and hydrologic models include terms torepresent hillslope throughflow advection; however, none of these modelshas been evaluated against field measurements of throughflow temperature.This study addresses the following objectives:1. To document spatiotemporal variability in throughflow temperaturesfor two winter seasons.2. To evaluate the utility of using heat as a tracer for hydrograph separa-tion.3. To evaluate commonly used methods to estimate throughflow temper-atures against field observations.These objectives are addressed using field data collected during twowinters at Griffith Creek. The focus is on winter conditions since in thecoastal Pacific Northwest heat advected to the stream by hillslope runoff isa major component of the stream heat budget during this season (Chapter3) and is less important during the summer (Moore et al., 2005b; Guenther,2007).The following is a note on terminology used in this chapter. The term‘subsurface temperature‘ is used here to refer to soil temperature measure-ments made both above and below the water table. ‘Throughflow’ refersto hillslope water moving laterally as saturated flow along the soil andbedrock/till interface. ‘Throughflow temperature’ refers to subsurface tem-perature measurements made only during periods when a sensor was belowthe water table. Throughflow temperature was used even during dry pe-riods between rain events, as long as the sensor was below the measured64water table. Because of the steep topography and placement of the welland subsurface temperature sites in relation to the channel, the majority oflocations had a steep topographic gradient from well to channel; however,a few wells likely were installed in local bedrock or till depressions wheresome local hillslope storage occurred even at low flows.4.2 MethodsThe analyses in this chapter draw upon the field measurements described inChapter 2. The following sections present descriptions of the data analysisand modelling.4.2.1 Spatiotemporal variability of throughflow temperatureIt was hypothesized that spatial variability in throughflow temperaturescould be explained by differences in vertical heat conduction, as well asdifferences in lateral advection of throughflow from upslope. Regressionmodels were used to address these hypotheses. Four predictor variableswere used to represent differences in vertical heat conduction: depth belowsurface, forest cover, aspect, and elevation. Observed temperature gradi-ents in soil are typically a function of conduction rates and temperaturesshould vary with depth below the ground surface. Forest cover, aspect, andelevation influence conduction rates by controlling the amount of energyavailable at the soil surface (Kang et al., 2000). Two predictor variableswere used to represent differences in lateral advection of throughflow fromupslope: topographic wetness index and upslope contributing area. Studieshave found that topographic wetness index and upslope contributing areacan be related to the timing and magnitude of shallow hillslope groundwa-ter response and that solutes, such as dissolved organic carbon, also vary inrelation to these hillslope variables (Jencso et al., 2010; Pacific et al., 2010).To the author’s knowledge, no studies have related upslope contributingarea and the topographic wetness index to soil temperature.Statistical analysis focused only on the 2012–2013 winter period, sincethe number of temperature loggers installed and located below the water65table was greater than during the 2011–2012 winter. Each hour, from October1, 2012 to April 30, 2013, linear regression models were fit to the subsurfacetemperatures at the well sites. In order to focus on temperatures that wouldbe representation of lateral subsurface throughflow, only those temperatureloggers that were identified to be below the water table during the hour ofinterest were used in the analysis.Topography-based predictor variables were generated from a 5 m reso-lution digital elevation model for the catchment. Upslope contributing areawas calculated using the deterministic 8 approach (O’Callaghan and Mark,1984). Forest cover was taken from the MKRF data inventory. Spatial datawere processed using the System for Automated Geoscientific AnalysesGeographic Information System (2013). All predictor variables were fit indi-vidually, as well as in combination with the other predictors. For multiplelinear regression models, the interaction terms were also evaluated. Modelresiduals were checked for outliers, non-linearities, and heteroscedasticity.4.2.2 Current approaches for modelling throughflowtemperatureSix published approaches for estimating shallow groundwater and through-flow temperatures for predicting stream temperature were evaluated againstobserved throughflow temperatures measured at Griffith Creek. The dif-ferent methods were evaluated by calculating the percentage of time thatsimulated throughflow temperatures fell within the 10th and 90th percentilesof the observed subsurface temperatures that were below the water table.Most of the approaches include adjustable parameters that require userinput or calibration. In lieu of calibration, parameter values were optimizedto maximize the percentage of time that the estimated throughflow temper-atures fell within the range of observed temperatures. This optimizationrepresents the best possible calibration of the models. Analysis for the modelcomparison focused on the 2012–2013 winter period since a greater numberof temperature sensors were located below the water table than during the2011–2012 winter period. The following section details how the individualapproaches were calculated for Griffith Creek.66Approaches based on mean annual air temperatureMean annual air temperature is used to represent the temperature of ground-water discharging into the stream for a variety of predictive stream tem-perature models (e.g. Allen et al. 2007; Ficklin et al. 2012; Null et al. 2013;Sun et al. 2014; MacDonald et al. 2014a). Mean annual air temperature wascalculated using 1990 to 2013 air temperature monitored at the Haney UBCRF Admin climate station located 2 km from Griffith Creek. Differences inelevation between the climate station (147 m) and Griffith Creek (mean ele-vation = 440 m) were accounted for by applying a lapse rate of 0.006 ◦C m−1(Oke, 1987). This adjustment resulted in a mean annual air temperature forGriffith Creek of 8.3 ◦C. Regression models relating mean daily air tempera-ture between the Griffith Creek forest meteorological station and the MKRFclimate station, based on the two year period that the forest station wasactive, suggest that the mean annual air temperature is 8.1 ◦C. Since somemodels use mean annual air temperature plus one or two degrees Celsius(Ficklin et al., 2012), the range between 8.1 and 10.3 ◦C was considered torepresent groundwater temperatures as represented by mean annual airtemperature.Soil and Water Assessment ToolFicklin et al. (2012) updated the stream temperature model in the Soil andWater Assessment Tool (SWAT). The updated approach includes simula-tion of temperature and water contribution from lateral subsurface flow,snowmelt and groundwater discharge. The temperature of lateral subsur-face flow on day t (Tlat in ◦C) is calculated as:Tlat(t) =λnn∑i=1Ta(t− i) (4.1)where λ is a calibration coefficient, Ta(t− i) is the mean daily air temperature(◦C) on day t − i, and n is the total number of days. Ficklin et al. (2012)recommended determining n by examining the day-to-day variation inobserved stream temperatures. If Tlat is less than or equal to 0.1 ◦C, then67Tlat is set to 0.1 ◦C. Snowmelt temperature is assumed to be 0.1 ◦C. SWATalso includes groundwater temperature, which is assumed to be the meanannual air temperature plus one or two degrees Celsius.For Griffith Creek, the optimization of Tlat was constrained by consider-ing values of λ between 0.5 and 1 and lags between 2 and 14 days. Theserepresent the range of values used by Ficklin et al. (2012) for eight catch-ments for which they calibrated SWAT.CEQUEAUThe hydrology model CEQUEAU (Morin and Couillard, 1990) includes aninterflow term that is used in the model’s stream temperature routine. Theoriginal model calculated the interflow temperature (Tint in ◦C) as:Tint(t) = [Tg(t) + Ta(t)]/2 (4.2)where Tg is the groundwater temperature (◦C) and is assumed constantthroughout the year, and Ta is mean daily air temperature (◦C). St-Hilaireet al. (2000) provided an updated estimate of Tint based on the force-restoremethod for predicting soil temperatures (Deardorff, 1978); however, rea-sonable values of Tint could not be calculated given the equation and infor-mation provided in St-Hilaire et al. (2000). For Griffith Creek, the originalequation for interflow temperature from Morin and Couillard (1990) wasused with a site specific range estimate of mean annual air temperatures(6.1 to 10.3 ◦C) for Tg and computed mean daily air temperatures from theabove-stream meteorologic station.HSPFThe Hydrological Simulation Program-Fortran (HSPF) (Bicknell et al., 2001)computes interflow temperature on day t (Tint(t)) as:Tint(t) = Tint(t− 1) + SHSPF · [Ta(t) + DHSPF − Tint(t− 1)] (4.3)68where Tint(t) is the interflow temperature from the previous time step,SHSPF is a smoothing factor, Ta(t) is air temperature on day t, and DHSPFis the difference between air temperature and the mean temperature ofinterflow, all in ◦C. Both SHSPF and DHSPF are user-adjustable parameters.Optimization of these parameters was contrained by considering rangesof 0.01 to 0.2 and 0 to 5 for SHSPF and DHSPF, respectively. These rangeswere partly informed by the HSPF calibration done by Marce´ and Armengol(2008) for a river in Spain.LARSIM-WTThe LARSIM-WT model (Bremicker, 2000; Haag and Luce, 2008) simulatesthe temperature of three different runoff components (Trc in ◦C) – directrunoff, interflow, and groundwater – using different parametrizations of thefollowing equation:Trc = max{yrc + βrc × (Ta − yrc)0.0(4.4)where yrc is the expected long-term mean water temperature of the runoffcomponent (◦C), βrc is the slope, and Ta is air temperature. The parameter yrcis typically set equal to the long-term mean annual air temperature, whichwas considered to be 8.1 ◦C for Griffith Creek. The slope parameter, βrc, isallowed to vary between 0 and 1. Haag and Luce (2008) set β equal to 0.75,0.30, and 0.05, for direct runoff, interflow, and groundwater, respectively.For Griffith Creek, values for βrc of 0 and 1 were used in the optimization ofTrc.MIKE SHEMIKE SHE models groundwater flow dynamics, but not heat transport, andthus requires an input temperature that represents the temperature of theupper portion of the aquifer and is assumed to be spatially uniform (Loinazet al., 2013, 2014). MIKE SHE does calculate a surface drainage flux that isaimed at representing heat contributions from overland flow. Comparing69the overland flow temperature calculations to the observed throughflowtemperatures is not an equal comparison. However, this comparison wasstill conducted since MIKE SHE does not include throughflow temperatures.The overland flow runoff temperature is calculated using an the equilibriumtemperature approach (Bogan et al., 2003). The equilibrium temperature isthe water temperature for which the net surface energy exchange is zero;that is, for Griffith Creek, the equilibrium temperature was determined bysolving for the surface heat budget:K∗ + L∗ + Qe + Qh = 0 (4.5)where K∗ is net shortwave radiation, L∗ is net longwave radiation, Qe is thelatent heat flux, and Qh is the sensible heat flux, all in W m−2. Equation 4.5can be expanded to the following, based on the equations used in Leach andMoore (2010) and in Chapter 3:(1− α)[D× g + S× fv] + [ fv × εa + (1− fv)εvt]σ(Ta + 273.2)4−εwσ(Tw + 273.2)4 + 627.8[0.0424 ·U · (ea − ew)]+0.66 · (P/1000) · (Tw − Ta)× 627.8(0.0424 ·U) = 0(4.6)where α is the albedo of the water surface, D is the direct component of inci-dent solar radiation, g is the gap fraction derived from digitized hemispheri-cal images for the sun’s position, S is the diffuse component of incident solarradiation (W m−2), fv is the view factor, εa, εvt and εw are the emissivitiesof the atmosphere, vegetation and terrain, and water, respectively, σ is theStefan-Boltzmann constant, Ta and Tw are air and water temperatures (◦C),U is wind speed (m s−1), ea and ew are the vapour pressures of air and runoff,respectively (kPa), P is ambient air pressure (kPa), and 627.8 accounts forunit conversion from mm h−1 to W m−2. For further details on the termsin equation 4.6 and how they were evaluated please see Chapter 3. Equa-tion 4.6 was solved using the Newton-Raphson method as implemented inthe ‘rootSolve’ package (Soetaert, 2009) in the R programming language (RDevelopment Core Team, 2014).704.3 Results4.3.1 Spatiotemporal variability of throughflow temperatureSubsurface temperatures ranged between 2 and 12 ◦C during the studyperiod (Figures 4.1 and 4.2). Stream temperature measured at the catchmentoutlet was generally closer to the lower range of hillslope temperatures,reflecting the cooling influence of energy exchanges acting on the stream(Chapter 3). During the middle of the winter, temperature loggers above thewater table recorded generally lower temperatures than loggers below thewater table. Two locations, GW037 (three depth measurements) and GW018(one depth measurement), which were located in a known groundwaterseepage area, exhibited subsurface temperatures that were more stablethroughout the study period than the remaining sites. Temperatures at thesesites were around 6 to 8 ◦C and are identified as the four traces that plotbelow most temperatures during September and October and above mosttemperatures during December to March. Water from these wells had higherelectrical conductivity (40 - 90 µs cm−1) than all other sites, which typicallyhad electrical conductivities of 8 to 15 µs cm−1 throughout the study period.Figure 4.3 shows subsurface temperatures for well sites with multipletemperature loggers installed at different depths. During rain events, tem-peratures were relatively uniform across depths with differences betweendepths often being less than 1 ◦C. Temperature variations with depth rangedup to 3.5 ◦C in the recession periods between rain events and during periodswhen a snow pack was present on the catchment. This graph also highlightsthe influence of the groundwater seepage zone on GW037, which results inrelatively stable and uniform temperatures around 7 to 8 ◦C for the entirestudy period across the 0.05, 0.25 and 0.5 m depths.Figure 4.4 shows plots of the interquartile range for throughflow tem-peratures, mean throughflow temperature minus air temperature, R2 valuefor the linear model using depth of sensor as the predictor, and dischargeand snow cover for the 2012–2013 winter period. The interquartile range inthroughflow temperatures was as great as 1.4 ◦C during periods with snow71Figure 4.1: Water inputs, air temperature (Ta), discharge (Q), groundwa-ter well water table elevations (WT), and subsurface temperatures(Tsub) for the 2011–2012 study period.on the catchment and during dry periods between rain events. During peri-ods of high runoff the interquartile range was around 0.2 to 0.4 ◦C. Depth ofthe temperature sensor was the only predictor variable of throughflow tem-perature that had some consistency in being statistically significant duringthe study period. Topographic wetness index, upslope contributing area,forest cover, aspect, and elevation were generally weak predictors of spatialvariability in throughflow temperature. The R2 value of the depth modelwas greatest (0.3 to 0.4) during dry periods between events and duringperiods when snow was on the catchment; however, depth was a weakpredictor of throughflow temperature during periods of higher runoff.72Figure 4.2: Water inputs, air temperature (Ta), discharge (Q), groundwa-ter well water table elevations (WT), and subsurface temperatures(Tsub) for the 2012–2013 study period.4.3.2 Throughflow temperature model comparisonFigure 4.5 shows the results of different modelling approaches to estimatingthroughflow temperatures and compares these against observed through-flow temperatures at Griffith Creek. Table 4.1 summarizes the optimizedparameterization of the six models and highlights their prediction success.The commonly used assumption that mean annual air temperature canapproximate throughflow temperatures results in consistent overpredictionfor most of the winter period. Of the range considered, a mean annual airtemperature of 9.2 ◦C provided the highest throughflow prediction rate at9%. Mean annual air temperature consistently overpredicts throughflowtemperatures by 1 to 4 ◦C.732610df13$dateTime.x[df13$siteName == "GW036.5"]df13$Temp[df13$siteName == "GW036.5"]0.05 m0.25 m0.50 m a) GW036124Water inputs (mm hr−1 )2610df13$dateTime.x[df13$siteName == "GW037.5"]df13$Temp[df13$siteName == "GW037.5"]0.05 m0.25 m0.50 m b) GW0372610df13$dateTime.x[df13$siteName == "GW014a"]df13$Temp[df13$siteName == "GW014a"]0.30 m0.42 m c) GW0142610df13$dateTime.x[df13$siteName == "GW022a"]df13$Temp[df13$siteName == "GW022a"]0.10 m0.20 md) GW022Oct Nov Dec Jan Feb Mar Apr MaySubsurface temperature (°C)2012−2013Figure 4.3: Subsurface temperatures for the four sites where more thanone temperature sensor was installed.Throughflow temperatures simulated by SWAT capture some of thetemporal variability in throughflow temperatures; however, overall thisapproach tends to underestimate throughflow temperatures during the win-ter. The model fails most notably during periods with snow and snowmelt,when it underpredicts by up to 5 ◦C. Overall, SWAT successfully predictedthroughflow temperatures 8% of the time.The CEQUEAU equation for throughflow plots within the range ofobserved throughflow temperatures 33% of the time. The CEQUEAU esti-mates are more temporally dynamic than the observed temperatures andmaximum under- and over-predictions are 3.5 and 3.2 ◦C, respectively.Of the models considered here, the HSPF interflow equation was themost successful at predicting throughflow temperatures. HSPF predictions74sptTemp$dateTimelog10(sptTemp$meanQ)0.0110Q (L/s)Snow0.21.0sptTemp$dateTimesptTemp$IQRIQR (°C)−15010sptTemp$dateTimetdiffSnowT tf−T a  (°C)0.00.4lm.run$dateTimelm.run$r2.dptR2  depthOct Nov Dec Jan Feb Mar Apr MayFigure 4.4: Discharge (Q), interquartile range of throughflow tempera-tures (IQR), mean throughflow temperature minus air tempera-ture (Tt f − Ta), and regression R2 value using depth below surfaceas predictor.fell within the range of observed throughflow temperatures 53% of the timeand only overestimated throughflow temperatures by a maximum of 1.3 ◦C.HSPF did not under-predict throughflow temperatures during the studyperiod.The optimized parameter value for the slope, βrc, of the linear LARSIM-WT function was found to be 0.53. This value resulted in a throughflowprediction rate of 31%. LARSIM-WT under- and over-predicted throughflowtemperatures by up to 4.0 and 3.0 ◦C, respectively.The MIKE SHE overland flow temperature using the equilibrium tem-perature approach consistently predicted equilibrium temperatures below75Figure 4.5: Comparison of different throughflow temperature estima-tion approaches. Blue lines represent the best estimated tempera-ture time series based on each approach and grey bands representthe range of the 10th and 90th percentiles of observed subsurfacetemperatures.76Table 4.1: Summary of throughflow temperature model comparison. In-cluded are whether the modelling approach accounts for snowmelt,the optimized parameter values, the percentage of time that theoptimized model predicts throughflow temperatures within theobserved range, and the maximum under- and over-predictionmade by the model using the optimized parameterization duringthe 2012–2013 study period.Includes Optimized Temperature Max under/overApproach snowmelt? parameter(s) prediction (%) prediction (◦C)Mean annual Ta No Ta = 9.2 ◦C 9% 1.8/4.4SWAT Yes λ = 1, n = 6 8% 5.1/2.5CEQUEAU No Tg = 8.6 ◦C 33% 3.5/3.2HSPF No SHSPF = 0.03 53% 0/1.3DHSPF = 2.4LARSIM-WT No βrc = 0.53 31% 4.0/3.0MIKE SHE No - 1% 10.9/6.30 ◦C during the winter, which are set equal to 0 ◦C. For periods with pos-itive equilibrium temperatures, the resulting temperatures overestimatedthroughflow temperatures by up to 6 ◦C. MIKE SHE predictions were suc-cessful 1% of the time.4.4 Discussion4.4.1 Spatiotemporal variability of throughflow temperatureWinter throughflow temperatures at Griffith Creek were variable in time,ranging up to 12 ◦C during autumn and spring and dropping to as low as2 ◦C during December and January. In contrast, throughflow temperatureswere less variable in space, where differences between sites were at timesup to 5 ◦C, but were mostly less than 2.5 ◦C. The magnitude in temperaturedifference across sites varied with hydroclimatic conditions. During dryperiods between rain and melt events, throughflow temperature variabilitywas partly accounted for by differences in depth of measurement below thesoil surface. Between-site throughflow temperature variability was minimalduring rain-on-ground events (0.2 to 0.4 ◦C), but increased to between 0.5 to771.5 ◦C during rain-on-snow events.No relationship was found between spatial patterns in throughflow tem-peratures and hillslope characteristics, such as upslope contributing areaor topographic wetness index, forest cover, aspect, or elevation. The lackof relationship may be partly due to the limited resolution and accuracy ofthe DEM used in this study. In addition, spatial differences in throughflowtemperatures were often modest (<2.5 ◦C); therefore, uncertainties in thevalues of the predictor variables may be too large to detect a signal. Over-all, throughflow temperatures were more variable in time than in space.Therefore, building models to predict throughflow temperatures may beable to treat catchments in a spatially-lumped approach during rain-on-ground events, but need to account for the role of transient snow duringrain-on-snow events.This study highlights that predictability of the spatial patterns of through-flow temperatures can be confounded by hillslopes that have some contri-bution of deeper and more regional groundwater inputs. A known seepagearea in the catchment (Guenther et al., 2014) was characterized by stablewater table elevations, throughflow temperatures, and higher electrical con-ductivities compared to the other hillslope sites. Throughflow contributionsfrom this area had temperatures that were relatively insensitive to changesin air temperature, rain events, or snowpacks. Measured throughflow tem-peratures from this site were within the mean annual air temperature rangeduring October and November, but were up to 2 ◦C lower than the mean an-nual air temperature for the rest of the winter period. While the importanceof this seepage area on the stream heat budget appears to be overwhelmedby the other hillslopes during the winter (Chapter 3), it is thermally im-portant and helps sustain flow in the stream during dry summer periods(Guenther, 2007; Guenther et al., 2014). In terms of being able to predictthe presence of this seepage area from terrain indices, the seepage areadid correspond to the hillslope in the catchment with the largest upslopecontributing area and topographic index values. However, some prelimi-nary groundwater modelling of this hillslope using the SUTRA distributedgroundwater flow model suggests that in order to maintain the stable wa-78ter tables and temperatures observed, the hillslope must be connected todeeper groundwater flowpaths, likely associated with flow in fractureswithin the granitic bedrock underlying the compacted till. Therefore, pre-dicting deeper groundwater discharge areas and their temperature fromonly terrain indices remains a challenge, but an important one for predictingstream temperature.4.4.2 Utility of heat as a hillslope hydrologic tracerThe use of heat as a tracer to determine pre-event and event water contribu-tions to streamflow or flowpath depths (Shanley and Peters, 1988; Kobayashiet al., 1999; Birkinshaw and Webb, 2010) may be confounded by spatial vari-ability in pre-event water temperatures and lack of subsurface temperaturesvariability with depth, respectively. As shown in this study, pre-event sub-surface temperatures can vary significantly across space and with depth,making it difficult to assign a single pre-event water temperature. Duringrain events, throughflow temperatures were relatively invariant with depth,therefore making it difficult to identify the depth-origin of streamflow basedon temperature. This work highlights that heat as a tracer for runoff watermay have limited utility, and the most effective use of temperature tracerswould be for systems dominated by deeper groundwater discharge thathave stable temperatures.4.4.3 Throughflow temperature model comparisonCurrent routines to estimate throughflow and overland flow temperaturesthat are embedded in coupled hydrologic and stream temperature modelsgenerally fail to provide accurate predictions of winter throughflow temper-atures at Griffith Creek. In particular, the effects of transient snow cover andthe associated snowmelt are not well represented.The HSPF model, original CEQUEAU function, and LARSIM-WT werethe most successful approaches at estimating throughflow temperaturesat Griffith Creek with prediction rates of 53%, 33%, and 31%, respectively.However, a major drawback of these approaches is that they do not ac-79count for snow, which limits their use as a predictive tool to understandstream temperature response to changes in climate and land cover, andcorresponding impacts on snowpack dynamics.The majority of stream temperature models use mean annual air tem-perature to represent the temperature of groundwater or throughflow dis-charging into the stream. This approach is generally intended to representgroundwater that has annual or longer residence times. For this situation,using mean annual air temperature is simple and may be a reasonable ap-proach, although decadal scale changes in climate and associated shifts inlongterm mean annual air temperature may confound selecting a represen-tative value (Kurylyk et al., 2014). Some models, such as DHSVM-RBM (Sunet al., 2014), use mean annual air temperature to represent the temperatureof water with much shorter residence times, such as throughflow or inter-flow. The DHSVM-RBM case study model application for a small catchmentnear Seattle overpredicted stream temperatures during high flow winterperiods, likely because the model sets these water inputs to 10 ◦C (the meanannual air temperature for the catchment). DHSVM-RBM also initializesthe headwater upstream temperature boundary using the water and airtemperature regression from Mohseni et al. (1998). That approach does notaccount for the influence of snowmelt. The results from the model compari-son in this study suggest that the winter stream temperature overpredictionby DHSVM-RBM is due, at least in part, to overprediction of throughflowtemperatures or incorrect initializing of the headwater temperatures.The updated SWAT stream temperature model (Ficklin et al., 2012) un-derpredicted runoff temperatures during most of the winter and particularlyduring snow and snowmelt events. SWAT assumes that snowmelt waterenters the stream at 0.1 ◦C based on results presented by Kobayashi (1985).However, Kobayashi (1985) did not assume that snowmelt entered thestream at this temperature, but instead used this information to help withhydrograph separation. In addition, Kobayashi (1985) did not observe soiltemperatures near 0.1 ◦C. Observed soil and throughflow temperaturesduring snowmelt events at Griffith Creek never dropped to 0.1 ◦C. There-fore, it is likely that in areas without frozen soils, SWAT will underpredict80stream temperatures during snowmelt periods because of the assumptionthat snowmelt enters the stream at 0.1 ◦C. Ficklin et al. (2012) acknowledgedthis shortcoming of the model during an evaluation for the Entiat River inthe PNW, which is a river characterized by frequent winter runoff events,similar to Griffith Creek. Similar to the Griffith Creek comparison, SWATunderpredicted Entiat River temperatures during winter runoff events. Thisunderprediction means that for evaluating climate change scenarios, SWATmay overestimate the cooling effect of snowmelt on stream temperature.MIKE SHE drainage runoff temperatures using the equilibrium tem-perature approach generally resulted in negative temperature predictions,thereby resulting in the model assuming overland temperatures were 0 ◦C.This is a known limitation of the equilibrium temperature approach, asBogan et al. (2003) routinely predicted negative temperatures during win-ter periods. The MIKE SHE drainage runoff temperature is intended torepresent overland flow and not throughflow, so the comparison is notentirely appropriate. However, MIKE SHE does not include throughflowtemperature and there is the risk that the drainage runoff temperature maybe used as a surrogate for throughflow temperature in some applications ofthe model to sites where throughflow is a dominant runoff process.The implication of this model comparison is that current models do notadequately simulate throughflow temperatures in catchments like GriffithCreek, which are characterized by high winter runoff and rain-on-snowevents. It is important to develop better models of throughflow temperaturesin order to provide more accurate predictions of climate and land coverchanges on stream temperatures.81Chapter 5Development and evaluation ofa conceptual hillslopethroughflow model5.1 IntroductionChapter 4 highlighted that existing methods to estimate throughflow tem-peratures embedded in catchment scale stream temperature models eitherdo not accurately predict throughflow temperatures observed at GriffithCreek or do not explicitly account for snowmelt contributions to through-flow temperatures. In contrast, process-based models that simulate soilwater and heat transport, such as SHAW (Flerchinger, 2000), COUP (Janssonand Karlberg, 2001), and ForHyM2 (Yin and Arp, 1993), can represent theinfluence of snowmelt; however, these models primarily solve soil ther-modynamics at the point scale and cannot be applied to systems wherelateral advection associated with throughflow is an important hydrologicprocess. Groundwater models, such as SUTRA (Voss and Provost, 2002) andHydroGeoSphere (Therrien et al., 2010), can be employed for two- and three-dimensional problems and include representation of throughflow. However,these models require substantial computational resources and highly skilled82operators to set up for an individual site. Therefore, a conceptual hillslopethroughflow temperature model was developed that was informed by fieldobservations and previous research on the hillslope hydrologic and thermalprocesses that dominate at the MKRF study site (Moore and Thompson,1996; Thompson and Moore, 1996; Donnelly-Makowecki and Moore, 1999;Hutchinson and Moore, 2000; Haught and van Meerveld, 2011). The modelpresented here aims to strike a balance between simple empirical relation-ships with air temperature that do not explicitly represent snow dynamics,as seen in Chapter 4, and complex models that are highly parameterizedand require considerable computing resources.The primary objective of this chapter is to develop and test the newconceptual-parametric hillslope throughflow model. The model is cali-brated using a generalized likelihood uncertainty estimation approach toquantify uncertainty of model predictions and account for issues of equi-finality (Beven and Binley, 1992; Beven and Freer, 2001). The calibratedmodel is first evaluated using a split-sample test based on the throughflowtemperature measurements. The model is also evaluated at the event scaleto test its ability to reproduce the influence of transient snow cover andassociated rain-on-snow events on throughflow temperatures, as revealedby the empirical analyses presented in Chapter 3.5.2 MethodsThe conceptual-parametric throughflow temperature model development,calibration and evaluation procedures are presented in the following sec-tions. The first section focuses on the model’s hydrology component. Thesecond section focuses on the model’s thermal component.5.2.1 Hydrology componentThe following sections outline the structure and governing equations of thehydrologic component of the model, the initial conditions, and the modelcalibration and evaluation procedure.83Structure and governing equationsThe model simulates hillslope throughflow and thermal dynamics of a rep-resentative hillslope (Figure 5.1). The model represents the catchment asthree reservoirs in serial configuration. The first two reservoirs representupper and lower portions of a hillslope and the third reservoir represents thestream channel. This model structure is consistent with the understandingof the dominant runoff processes at the MKRF (Donnelly-Makowecki andMoore, 1999) and those found in other headwater catchments in humidtemperate environments that are characterized by shallow soils overlyingrelatively impermeable compacted till or bedrock (e.g. Seibert and McDon-nell 2002; Seibert et al. 2003). The model operates on a sub-daily time step(hourly for this study). Each hillslope reservoir consists of a duff layer, amineral soil layer and a basal layer representing the underlying bedrock andtill. The mineral soil layer comprises vadose and phreatic zones, separatedby a dynamic water table. Each hillslope reservoir has an independent watertable.The model is formulated as a set of differential equations that representthe water balances of the duff and mineral soil layers for the two hillslopereservoirs and the water balance of the channel reservoir. The rate of changein water storage of the duff layers for reservoirs 1 and 2 (dSd/dt in m s−1) is:dSddt= I − qd (5.1)where Sd is the water storage depth (m) of the duff layer, I is rate of waterinput (rain and snowmelt) entering the duff layer (m s−1), and qd is the rateof water flow from the duff layer to the mineral soil layer (m s−1).Rainfall and snowmelt measured at the forest site lysimeter are used aswater input to the model. Water inputs are retained in the duff layer untilthe water retention capacity of the duff layer is exceeded, after which water84Figure 5.1: Conceptual hillslope throughflow temperature model.inputs are delivered to the mineral soil layer:qd ={0 if Sd < θd(zd − zm)I if Sd ≥ θd(zd − zm)(5.2)where θd is the water retention capacity of the duff layer, and zd and zm arethe heights (m) of the top of the duff layer and top of the mineral soil layer,respectively.The rate of change in water storage of the mineral soil layers in reservoirs1 and 2 is:dSs,1dt= qd − Flat,1/A1 − E1 (5.3)85dSs,2dt= qd + (Flat,1 − Flat,2)/A2 − E2 (5.4)where Ss,i is the soil moisture storage for the mineral soil layers (m) ofreservoir i, Flat,1 is the rate of lateral water flow (m3 s−1) from the upperreservoir to the lower reservoir, Flat,2 is the rate of lateral water flow (m3 s−1)from the lower reservoir to the stream channel, A1 and A2 are the areas (m2)of reservoir 1 and 2, respectively, and E1 and E2 are the evapotranspirationrates (m s−1) from reservoir 1 and 2, respectively. A1 and A2 are calculatedas the total catchment area multiplied by farea and 1− farea, respectively.farea is the area fraction of the upslope reservoir and 1− farea is the areafraction of the downslope reservoir.Water drains laterally from the upper reservoir to the bottom reservoir,and from the bottom reservoir to the stream channel, when the water tableelevation of the mineral soil layer, zwt (m) exceeds the residual storage depth,zdz (m):Flat,i = sin(ζ) cos(ζ)ai(zwt,i − zdz,i)bi (5.5)where ζ is the the mean slope angle (radians), ai is a constant, and bi is therecession exponent of reservoir i.The water table elevation (zwt) for each mineral soil layer is calculatedas:zwt = (Ss − zmθs)/(φs − θs) (5.6)where zm (m) is the depth of the mineral soil layer, φs is the effective porosityof the mineral soil, and θs is the water retention capacity of the mineral soil(m3 m−3).Evapotranspiration losses from the forest were computed using theapproach of Spittlehouse and Black (1981). The energy-limited potentialevapotranspiration rate, Emax (m s−1), is estimated using the Priestley-Taylor86method (Priestley and Taylor, 1972):Emax = αEδδ+ γQ∗Lv(5.7)where αE is the Priestley-Taylor coefficient, δ is the relation between satura-tion vapour pressure and air temperature (kPa ◦C−1), γ is the psychrometricconstant (0.070 kPa ◦C−1 for Griffith Creek), Q∗ is the net radiation at thetop of the forest canopy (W m−2), and Lv is the latent heat of vapourizationof water (J kg−1). Net radiation at the top of the canopy was calculated fromQ∗ = (1− α)K ↓ +L∗ (5.8)where α is the albedo of the forest canopy, assumed to be 0.1 for coniferoustrees (Oke, 1987; McMahon et al., 2013), K ↓ is incoming shortwave radiationmeasured at the open site meteorological station, and L∗ is net longwaveradiation. Net longwave radiation at the top of the forest canopy wascalculated using the equations in Chapter 3 assuming that the temperatureof the canopy was equal to the air temperature measured at the forestmeteorological station and using an emissivity for the vegetation of 0.97(Oke, 1987).A soil limited rate of evapotranspiration, Es (m s−1), was calculated as:Es = βE(Ss/zm)− θwiltφs − θwilt(5.9)where βE is a coefficient determined through calibration, and θwilt is the soilmoisture wilting point (m3 m−3). The actual evapotranspiration (E) is equalto the lesser of Emax and Es.The rate of change in water storage of the stream channel (dVc/dt inm3 s−1) is:dVcdt= Flat2 − Fc (5.10)87where Vc is the water storage volume of the stream channel (m3), and Fc isthe water discharge rate (m3 s−1) from the third reservoir, calculated as:Fc = ac(Vc −V0)bc (5.11)where V0 is residual storage volume in the channel (m3), and ac and bcare the recession constant and exponent, respectively, of reservoir 3. Theparameters ac and bc were estimated using field estimates of channel storagebased on channel width and depth measurements made over a range ofdischarge magnitudes.The set of differential equations is integrated at six-minute intervalsusing a fourth-order Runge-Kutta scheme in the ‘deSolve’ package in R(Soetaert et al., 2010). The six-minute values of q3 were used to calculatehourly means, which were compared with observed stream discharge.Initial conditionsInitial conditions for the hydrologic equations were specified as follows:Sd = 0 (5.12)zwt = zdz (5.13)Ss = zwtφs + θs(zm − zwt) (5.14)Vc = V0 (5.15)It was necessary to set the initial water table depth greater than zerodue to issues with heat content calculations for the phreatic zone. This88assumption is reasonable given that some groundwater wells registered awater table even during dry periods, suggesting that hillslope storage inlocal depressions exists at these sites.The model was run using fifteen years of daily precipitation and air tem-perature data from the MKRF headquarters climate station that precededthe period of interest in order to remove any influence of the initial condi-tions. Precipitation and air temperature were corrected for the elevationdifferences between the MKRF climate station and Griffith Creek.Calibration and evaluationThe hydrologic component of the model was calibrated using a generalizedlikelihood uncertainty estimate (GLUE) approach (Beven and Binley, 1992).The model was run 100 000 times using parameter sets randomly selectedfrom uniform distributions of select parameters. Table 5.1 describes thehydrologic parameters used in the model and the range of values used inthe calibration. Parameter ranges were determined from field observations,previous studies conducted at MKRF (Donnelly-Makowecki and Moore,1999), or typical reference ranges (Campbell and Norman, 1998; Reddinget al., 2005). Single values indicate that the parameter was fixed for thecalibration.The model was calibrated by comparing hourly means of Fc with ob-served discharge at the catchment outlet. The logarithm of discharge wasnot used for calibration because of issues with zero modelled dischargevalues during low flow periods. In addition, heat advection via hillsloperunoff is most critical for stream thermal regimes during high flow periods;therefore, it was deemed more important to successfully represent highflow periods rather than low flow periods. Parameter sets resulting in aNash-Sutcliffe efficiency (NSE) of 0.89 and greater and maximum zwt lessthan 0.6 m (the depth of the mineral soil layer) were retained for use in thetemperature component calibration. For the hydrology calibration, 487 pa-rameter sets met the criteria and were used in the calibration of the thermalcomponent of the hillslope runoff temperature model.89Table 5.1: List of hydrologic parameters used in the hillslope through-flow temperature model.Hydrologic parameter Symbol Units Range in GLUEHeight of duff layer zd m 0.8Height of mineral soil zm m 0.6Duff water retention capacity θd m3 m−3 0.2 - 0.4Mineral soil water retention capacity θs m3 m−3 0.1 - 0.25Duff effective porosity φd - 0.7 - 0.95Mineral soil effective porosity φs - 0.3 - 0.5Mean hillslope slope angle ζ radians 0.2583Outflow constant of reservoirs 1 and 2 a1, a2 s−1 0 - 30Outflow constant of reservoir 3 a3 h−1 1.77Recession exponent of reservoirs 1 and 2 b1, b2 1 - 5Recession exponent of reservoir 3 b3 1.14Priestley-Taylor constant αE 0.6 - 1.26Soil limited evapotranspiration constant βE mm d−1 0 - 20Permanent wilting point θwilt m3 m−3 0.05 - 0.1Residual storage depth of reservoirs 1 & 2 zdz,1, zdz,2 m 0.01Fractional area of upslope reservoir farea 0 - 1Note: All units are converted to m and s in the model.5.2.2 Thermal componentWithin each reservoir, the model simulates vertical conduction betweenthe model surface boundary, duff, vadose, phreatic and basal layers, aswell as vertical advection through the duff and vadose layers. Lateral heatadvection is simulated from the phreatic zone of the upslope reservoir to thedownslope reservoir, and from the phreatic zone of the downslope reservoirto the stream. A key assumption of the model is that temperatures within theduff, vadose and phreatic zones are uniform within each hillslope reservoir.Structure and governing equationsThe thermal component of the model is structured as a set of differentialequations that calculate the heat content and transfer for the duff, vadose,and phreatic layers. The rates of change in heat storage for the duff layer ofreservoirs 1 and 2 are calculated as:dHddt= Qc,s−d + Qc,v−d + ρwcw ITs − ρwcwqdTd (5.16)90where Hd is the heat storage in the duff layer (J m−2), Qc,s−d is heat con-duction between the duff surface and duff layer (W m−2), Qc,v−d is heatconduction between the duff and vadose layer (W m−2), ρw is the densityof water (kg m−3), cw is the specific heat capacity of water (J kg−1 ◦C−1), Tsis the boundary temperature at the surface of the duff layer (◦C), and Td isthe duff temperature (◦C). The last two terms represent vertical advectivefluxes associated with water inputs to the duff layer and water losses to thevadose layer.Heat conduction between the duff surface and duff layer is calculatedas:Qc,s−d = kd(Ts − Td)/[(zd − zm)/2] (5.17)where kd is the bulk thermal conductivity (W m−1 ◦C−1) of the duff. Conduc-tion is calculated using the distance between the duff surface to the middleof the duff layer.Heat conduction between the duff layer and vadose layer is calculatedas:Qc,v−d = (Tv − Td)/[((zd − zm)/2)/kd + ((zm − zwt)/2)/kv] (5.18)where Tv is the vadose temperature (◦C) and kv is the bulk thermal conduc-tivity (W m−1 ◦C−1) of the vadose layer. Conduction is calculated using thedistance from the middle of the duff layer to the middle of the vadose zone.The temperature at the duff surface, Ts, was set equal to the air tem-perature measured by the climate stations, except when snow was presenton the catchment, in which case Ts was set to 0 ◦C. Setting Ts to the localair temperature was assumed to be a reasonable approach, in contrast tosimulating a full energy budget by incorporating radiative and turbulentenergy exchanges, because canopy closure was extensive at the study sitesand heavy cloud cover persisted for much of the study period. Therefore,longwave radiation likely dominated the energy exchanges at the duff layersurface and would be adequately represented by local air temperature. Spot91measurements of the duff surface temperature were made during four dayswith contrasting weather conditions using a thermal imaging camera andfound to be within 1 ◦C of local air temperature, supporting the validity ofthis assumption for winter conditions under a forest canopy.The rate of change in heat storage for the vadose layers is calculated as:dHvdt= Qc,d−v + Qc,p−v + Qwt,v −QE,v + ρwcwqdTd − ρwcwqdTv (5.19)where Hv is the heat storage in the vadose layer (J m−2), Qc,d−v is conductionbetween the duff and vadose layers (W m−2), Qc,p−v is conduction betweenthe phreatic and vadose layers (W m−2), Qwt,v is the change in heat storageassociated with water table changes (W m−2), QE,v is advective heat loss dueto evapotranspiration (W m−2), and the final two terms represent verticaladvective fluxes associated with water inputs to the vadose layer and waterlosses to the phreatic layer.Heat conduction between the duff and vadose layers (Qc,d−v) is calcu-lated as:Qc,d−v = (Td − Tv)/[((zd − zm)/2)/kd + ((zm − zwt)/2)/kv] (5.20)where conduction is calculated using the distance from the middle of theduff layer to the middle of the vadose zone.Heat conduction between the vadose and phreatic layers (Qc,p−v) iscalculated as:Qc,p−v = (Tp − Tv)/[((zm − zwt)/2)/kv + (zwt/2)/kp] (5.21)where Tp is the phreatic layer temperature (◦C) and kp is the bulk thermalconductivity (W m−1 ◦C−1) of the phreatic layer. Conduction is calculatedusing the distance from the middle of the vadose layer and middle of thephreatic layer.The change in heat storage associated with water table changes (Qwt,v)92is calculated as:Qwt,v ={(−dzwt/dt)Tv[(1− φs)ρmcm + θsρwcw] if dzwt/dt > 0(−dzwt/dt)Tp[(1− φs)ρmcm + θsρwcw] if dzwt/dt < 0(5.22)where ρm is the density of the mineral soil (kg m−3), and cm is the specificheat capacity of the mineral soil (J kg−1 ◦C−1).The advective heat loss associated with evapotranspiration is assumedto be drawn from both the vadose and phreatic layers, proportional to therelative depths of the two layers. The evapotranspiration heat loss from thevadose layer (QE,v) is calculated as:QE,v = ρwcwTvE[(zm − zwt)− zdz]/(zm − zdz) (5.23)The rate of change in heat storage of the phreatic zone (dHp/dt) is calcu-lated as:dHpdt= Qc,v−p + Qc,b−p + Qadv,p + Qwt,p −QE,p (5.24)where Hp is the heat content of the phreatic zone (J m−2), Qc,v−p is the heatconduction between the vadose and phreatic layers, Qc,b−p is the heat con-duction between the basal and phreatic layers, Qadv,p is net advective heatflux, Qwt,p is the change in heat storage associated with changes in the watertable elevation, and QE,p is advective heat loss due to evapotranspiration,all in W m−2.Heat conduction between the vadose and phreatic layers (Qc,v−p) iscalculated as:Qc,v−p = (Tv − Tp)/[((zm − zwt)/2)/kv + (zwt/2)/kp] (5.25)where conduction is calculated using the distance between the middle ofthe vadose layer to the middle of the phreatic layer.Heat conduction between the basal and phreatic layers (Qc,b−p) is calcu-93lated as:Qc,b−p = (Tb − Tp)/[(zwt/2)/kp + (dx/2)/kb] (5.26)where Tb is the basal temperature for the uppermost basal layer (◦C), kb isthe bulk thermal conductivity of the basal layer (W m−1 ◦C−1), and dx is thethickness of the uppermost basal layer (m). Conduction is calculated usingthe distance between the middle of the phreatic layer to the middle of theuppermost basal layer.The net advective flux associated with vertical and lateral advection forreservoir 1 is calculated as:Qadv,p,1 = ρwcwqd,1Tv,1 − ρwcw(Flat,1/A1)Tp,1 − ρwcwEp,1Tp,1 (5.27)where the first term represents advective gains from vertical water inputsfrom the vadose zone, the second term represents advective losses due tolateral water losses to reservoir 2, and the final term represents advectivelosses due to evapotranspiration.The net advective flux associated with vertical and lateral advection forreservoir 2 is calculated as:Qadv,p,2 = ρwcwqd,2Tv,2 + (ρwcwFlat,1Tp,1 − ρwcwFlat,2Tp,2)/A2−ρwcwEp,2Tp,2(5.28)where the first term represents advective gains from vertical water inputsfrom the vadose zone, the second term represents lateral advective gainsfrom reservoir 1 and advective losses due to lateral water losses to the chan-nel, and the final term represents advective losses due to evapotranspiration.Change in heat storage associated with water table elevation changesfor the phreatic zone (Qwt,p) is calculated as:Qwt,p ={(dzwt/dt)Tv[(1− φs)ρmcm + θsρwcw] if dzwt/dt > 0(dzwt/dt)Tp[(1− φs)ρmcm + θsρwcw] if dzwt/dt < 0(5.29)The advective heat loss associated with evapotranspiration is assumed94to be drawn from both the vadose and phreatic layers, proportional to therelative depths of the two layers. The evapotranspiration heat loss from thephreatic layer (QE,p) is calculated as:QE,p = ρwcwTpE(zwt − zdz)/(zm − zdz) (5.30)The basal zone is divided into fifty layers extending from the bottom ofthe mineral soil layer to a depth where the temperature is assumed to remainrelatively constant throughout the year (zdeep in m). The value for zdeep wasdetermined by calculating the damping depth (Campbell and Norman, 1998)at which temperature fluctuations are expected to be reduced to 5% of thesurface value, which was around 15 m for the thermal diffusivity rangeconsidered in the calibration. The change in temperature through time(dTb,i/dt) for each basal zone layer (i) was calculated using a differentialequation based on an explicit finite difference solution:dTb,idt= κ(Tb,i−1 − 2 · Tb,i + Tb,i+1)/dx2 (5.31)where κ is the thermal diffusivity (m2 s−1) of the basal layer and dx is thethickness of each basal layer (m). For the uppermost layer, Tb,i−1 was initiallyset equal to Tp and for the lowermost layer, Tb,i+1 was fixed as Tdeep, whichis the mean annual air temperature for the catchment, corrected for thegeothermal gradient (Wood and Hewett, 1982). The stability of the finitedifference scheme for all calibration parameter sets was evaluated using thevon Neumann stability index (Isaacson and Keller, 2012). Unique basal zonetemperatures were calculated for the two hillslope reservoirs. The thermaldiffusivity of the basal zone was calculated as:κ = kb/(ρb · cm) (5.32)where kb is the effective thermal conductivity of the basal layer (W m−1 ◦C−1),and ρb is the density of the basal zone (kg m−3), calculated as:ρb = φbρw + (1− φb) · ρm (5.33)95where φb is the porosity of the basal zone.Bulk thermal conductivities for the duff, vadose, phreatic, and basallayers are calculated as:kd = Φwkw +Φaka +Φoko (5.34)kv,p = Φwkw +Φaka +Φmkm (5.35)kb = Φwkw +Φbedkbed (5.36)where Φw, Φa, Φo, Φm, Φbed are volume fractions and kw, ka, ko, km, andkbed are thermal conductivities (W m−1 ◦C−1) of water, air, organic matter,mineral soil, and bedrock, respectively.Temperatures of the duff (Td), vadose (Tv), and phreatic (Tp) layers (◦C)are calculated as:Td = Hd/[Sdρwcw + (1− φd)ρoco(zd − zm)] (5.37)Tv = Hv/[(Ss − φszwt)ρwcw + (1− φs)ρmcm(zm − zwt)] (5.38)Tp = Hp/(φsρwcwzwt + (1− φs)ρmcmzwt) (5.39)where Hd, Hv, and Hp are heat storage (J m−2) of the duff, vadose, andphreatic layers, ρw, ρo, ρm are the densities (kg m−3) of water, organic matter,and mineral soil, cw, co, cm are the specific heat capacities (J kg−1 ◦C−1) ofwater, organic matter and mineral soil, and the remaining variables aredefined in section conditionsInitial conditions for the thermal component of the model were calculatedas:Hd = Tint[Sdρwcw + (1− φd)ρdcd(zd − zm)] (5.40)Hv = Tint[(Ss − φszwt)ρwcw + (1− φs)ρmcm(zm − zwt)] (5.41)Hp = Tint[φsρwcwzwt + (1− φs)ρmcmzwt] (5.42)where the initial temperature, Tint, was set to 10 ◦C. The temperatures ofeach of the fifty basal layers were linearly interpolated between the phreaticzone temperature at 10 ◦C and Tdeep of 8.5 ◦C. The model was run usingdaily data from the MKRF climate station for fifteen years prior to the periodof interest in order to remove any influence of the initial conditions.CalibrationThe coupled hillslope hydrology and temperature model was calibratedto the 2011–2012 winter data. A GLUE approach was used and the modelwas run 100 000 times. For each run, a set of hydrology parameters wasrandomly selected from the 487 best parameters sets from the hydrologycalibration. The thermal parameters calibrated are listed in Table 5.2. Fourparameters (in addition to the hydrology parameter sets) were allowed tovary for the thermal calibration: the thermal conductivity of organic matter(ko), thermal conductivity of mineral particles (km), thermal conductivity ofthe bedrock material (kbed), and porosity of the basal layer (φb). Thermalconductivity ranges for the organic, mineral particles, and bedrock materialwere based, in part, on values presented by Campbell and Norman (1998).The lower thermal conductivity bound for the soil mineral particles is lower97Table 5.2: List of thermal parameters used in the hillslope throughflowmodelThermal parameter Symbol Units Range in GLUEWater density ρw kg m3 1000Organic matter particle density ρo kg m3 1300Mineral soil particle density ρm kg m3 2650Specific heat of water cw J kg−1 ◦C−1 4180Specific heat of organic matter co J kg−1 ◦C−1 1920Specific heat of mineral particles cm J kg−1 ◦C−1 870Thermal conductivity of water kw W m−1 ◦C−1 0.6Thermal conductivity of air ka W m−1 ◦C−1 0.024Thermal conductivity of organic matter ko W m−1 ◦C−1 0.05 - 0.5Thermal conductivity of mineral particles km W m−1 ◦C−1 1.5 - 3.5Thermal conductivity of bedrock kbed W m−1 ◦C−1 2.5 - 3.5Porosity of basal zone φb - 0.01 - 0.1Temperature of deep layer Tdeep ◦C 8.5Depth of deep zone zdeep m 15Thickness of basal layers dx m 0.3than typical reported values (Campbell and Norman, 1998). This lowerbound is not unreasonable, since it could reflect the substantial presenceof organic matter in this layer (Utting, 1979), which will cause the effectivethermal conductivity to be lower than if no organic matter were present.The temperature at zdeep was set to the mean annual air temperature andadjusted for the geothermal gradient, based on a gradient of 25 ◦C per 1000m (Wood and Hewett, 1982). The best parameter sets were selected basedon the percentage of time that the predicted temperature of the phreaticzone of the second reservoir fell within the 10th and 90th percentiles of theobserved throughflow temperatures. A criterion of 96% or greater was usedfor the temperature calibration, yielding 30 behavioural parameter sets.The model was calibrated using the microclimate and snow cover datafrom the forested portion of the catchment because the majority of thecatchment was forested. It was deemed not possible to classify the observedthroughflow temperatures into forest or harvest hillslopes, because therewas no reliable approach to determine the hillslope drainage area for eachwell based on the resolution of available digital elevation models (5 m). The98snow regimes of the forested and harvested regions of the catchment diddiffer; therefore, the sensitivity of the throughflow temperature predictionswas assessed by running the model for both snow regimes.EvaluationThe model was evaluated by comparing modelled throughflow tempera-tures against observed data for the 2012–2013 winter period. In addition,two sets of modelling exercises were conducted with different snow coverscenarios to assess whether the model could produce lower throughflowtemperatures during rain-on-snow than during rain-on-ground events asdocumented by the empirical analysis in Chapter 3. For each scenario themodel was run with the 30 best fit parameter sets and the resulting through-flow temperatures were compared to the model simulations of current (i.e.baseline) conditions.The first modelling exercise isolated individual rain-on-snow eventsduring both the 2011–2012 and 2012–2013 winters, and the snow cover wasremoved in the model at the onset of the rain event. All the rain-on-snowevents for the 2011–2012 and 2012–2013 winter periods were identified foruse in the analysis. Lysimeter data were used in conjunction with the opensite incoming solar radiation and relative humidity data sets to distinguishrain-on-snow melt events from melt events without rain.The second modelling exercise involved adding a simulated snow coverto a rain-on-ground event that occurred mid-March 2013. Three snow coverscenarios were explored in this second exercise: 1) snow was added onlyduring the period prior to the rain event, 2) snow was added only for theperiod during which rain was falling, and 3) snow was added both duringthe period prior to the rain event and during the rain event. The mid-Marchrain-on-ground event was chosen because it satisfied a number of conditionsthat made it plausible that the event could have been a rain-on-snow eventhad a snowpack been present on the catchment: 1) the event was precededby a dry period with relatively low air temperatures (< 5 ◦C) similar toantecedent conditions for other rain-on-snow events, 2) the incoming rain99temperature during the event (5 to 7 ◦C) was similar to other mid-winterrain-on-snow events, and 3) air temperatures remained relatively low fora period after the rain event. These conditions made the rain-on-groundand rain-on-snow comparison more realistic than using a rain event morecharacteristic of relatively warm autumn or spring rain events. In addition,the event was chosen because it allowed simulation of a larger rain-on-snowevent (198.1 mm of water inputs over a 100 hour period) not otherwiseobserved during the two winters of this study.5.3 Results5.3.1 Hydrologic calibration and evaluationFigure 5.2 shows time series of observed and modelled discharge for thecalibration period (2011–2012) and evaluation period (2012–2013) at GriffithCreek. The model efficiencies for the evaluation period ranged between 0.79and 0.86 (Figure 5.3). During the evaluation period, the model tended tounderpredict observed discharge. The most noticeable period of underpre-diction occurred in late January and February.Figure 5.4 shows ‘dotty plots’ for the hydrology model parameters.The plots highlight the issue of equifinality since many of the parametersperformed equally well across their range.5.3.2 Thermal calibration and evaluationFull winter periodFigure 5.5 shows the results of the calibration of the hillslope model’s ther-mal component. A total of 30 parameter sets resulted in modelled through-flow temperatures (i.e. temperature of the second reservoir phreatic zone)falling within the 10% and 90% quantiles of saturated soil temperatures96% of the time or greater. For the evaluation period, these parameter setsresulted in modelled throughflow temperatures falling within the observedrange between about 75 and 78% of the time (Figures 5.6 and 5.7).100Figure 5.2: Hydrographs for calibration and evaluation periods. Greybands are modelled discharge for the ensemble best parametersets. Black lines are observed discharge.Figure 5.3: Histograms of Nash-Sutcliffe efficiencies for the hydrologycalibration and evaluation periods.101Figure 5.4: ‘Dotty plots’ of calibrated model parameters and corre-sponding Nash-Sutcliffe efficiencies for hydrology calibration.Figure 5.8 compares modelled throughflow temperatures using the forestand harvest snow regimes for the calibration and evaluation periods. Forthe evaluation period, both snow regimes resulted in predicted throughflowtemperatures consistently plotting within the observed range. When theharvest snow regime was used for the evaluation period instead of the forestsnow regime, the percentage of time that the model successfully predictedthroughflow temperatures increased to a range of 81 to 87%, with a mean102Figure 5.5: Hillslope throughflow temperature calibration for 2011–2012. Grey lines are the predicted runoff temperatures, blue isperiods with snow cover, black lines are the 10th and 90th per-centiles of saturated soil temperatures. Observed and predictedhydrographs are provided for reference.of 85%. This increase in performance when using the harvest snow regimeoccurred because the snow was found to have disappeared earlier than inthe forested area for the early March event.Figure 5.9 shows ‘dotty plots’ for the calibrated thermal parameters,as well as the hydrologic parameters that influence the calculation of bulkthermal conductivities. Thermal conductivities of the organic matter and soilmineral particles appeared to be critical parameters, with better throughflowtemperature predictions occurring at the lower range of ko and the higherrange of km. Better temperature predictions appear to be associated with103Figure 5.6: Hillslope throughflow temperature model evaluation for2012–2013. Grey lines are the predicted runoff temperatures, blueis periods with snow cover, black lines are the 10th and 90th per-centiles of saturated soil temperatures. Observed and predictedhydrographs are provided for reference.lower values of the duff layer water retention capacity. This is consistentwith the better model performances being associated with a lower bulkthermal conductivity for the duff layer.Rain-on-snow event evaluationNine rain-on-snow events were identified and are summarized in Table 5.3and shown in Figure 5.10. The maximum temperature difference betweenrain-on-snow and rain-on-ground events varied between 0.01 ◦C for ROSevent 1 to 0.97 ◦C to ROS event 8. The magnitude in throughflow temper-104Figure 5.7: Histograms of the Nash-Sutcliffe efficiency and percentageof time of successful temperature prediction for the calibrationand evaluation periods.105Figure 5.8: Comparison of predicted throughflow temperatures for theforest (blue) and harvest (red) snow regimes during the calibrationperiod (2011–2012) and evaluation period (2012–2013). Horizontallines indicate periods of snow cover.ature difference had a generally positive relationship with the amount ofwater input to the soil during the event.Figure 5.11 shows the throughflow temperature model comparisons ofthe three snow cover scenarios for the rain-on-ground event that was simu-lated as a rain-on-snow event. Throughflow temperatures were up to 0.5 ◦Clower during the rain-on-snow event than during the rain-on-ground eventwhen the snow cover was present only prior to the rain event. In contrast,when snow cover was added only during the rain, throughflow tempera-tures were a maximum of about 2 ◦C lower than the rain-on-ground event.When snow was added both before and during the rain event, throughflow106Figure 5.9: Calibrated model parameters and corresponding percentageof correct temperature prediction for the thermal model calibra-tion.Table 5.3: Summary of rain-on-snow events during the 2011–2012 and2012–2013 winter periods.Event Event start (PST) Duration (h) Water inputs (mm) Max modelled ∆T (◦C)1 2011-11-21 08:00 8 7.4 0.012 2012-01-20 11:00 119 152.9 0.483 2012-01-28 20:00 14 27.0 0.084 2012-02-18 15:00 22 43.5 0.045 2012-03-02 06:00 34 36.2 0.116 2012-03-14 11:00 26 55.3 0.237 2013-01-03 19:00 114 91.6 0.358 2013-02-28 10:00 24 107.8 0.979 2013-03-19 17:00 17 27.8 0.18107Figure 5.10: Modelling exercise to evaluate model performance during individual rain-on-snow and rain-on-ground events. First row of plots show air temperature (Ta), second row show total water inputs(rain and snowmelt) and period of snow cover (horizontal lines), third row shows absolute simulatedtemperatures during rain-on-ground (ROG) and rain-on-snow (ROS) events, and fourth row showsrelative difference between throughflow temperatures simulated as rain-on-ground minus rain-on-snow events. Rain-on-snow events 1 to 6 are from 2011–2012 winter and rain-on-snow events 7 to 9 arefrom 2012–2013 winter.108temperatures were a maximum of 2.3 ◦C lower than during the rain-on-ground event. The results highlight that although presence of a snow packprior to a rain-on-snow event can decrease subsurface and throughflowtemperatures, the majority of the lower throughflow temperatures associ-ated with rain-on-snow events is due to the cold water inputs from snowmelt and cooled rain that has passed through the snow pack. This mod-elling exercise also highlights that the thermal impacts of the snow coverand rain-on-snow event may persist for many weeks after the snow hasdisappeared.5.4 Discussion5.4.1 Evaluation of model performanceThe percentage of time that the model successfully predicted throughflowtemperatures during the calibration period was 96% or greater and droppedto 75 to 78% during the evaluation period. Although prediction rates werelower during the evaluation period, modelled throughflow temperatureswere never greater than 0.7 ◦C outside the 10th and 90th percentiles of ob-served temperatures, and the model performed far better than the existingmodels. The drop in performance during the evaluation period could, inpart, be due to a lower number of throughflow observations used to calibratethe model than were used to evaluate the model. A wider range of through-flow temperatures could have been sampled had more observations beenmade during the calibration period. This could have resulted in selectionof different behavioural parameter sets during calibration that may haveprovided better predictions during the evaluation period.The lower model performance during the evaluation also appears tobe associated with a transition period between a rain-on-snow event torain-on-ground event in early March 2013, as snow cover disappeared. Themodel predicted decreases in throughflow temperatures during this period,whereas observed throughflow temperatures increased. As shown in Figure5.8, when using the harvest snow regime, which had a shorter period of109Figure 5.11: Rain-on-ground event simulated as a rain-on-snow eventwith three different snow cover scenarios. Top plot shows airtemperature (Ta), second plot shows total water inputs (rain andsnowmelt) and periods of snow cover (horizontal lines), thirdplot shows absolute simulated throughflow temperatures duringthe four snow cover scenarios, and fourth plot shows relativedifference between throughflow temperatures for the four snowcover scenarios.110snow cover during the early March rain-on-snow event, the accuracy of themodel improved with 81 to 87% of the predicted throughflow temperaturesfalling between the 10th and 90th percentiles of observed temperatures. Thismodelling error highlights limitations and possible errors in the approachused to classify snow cover. The catchment was assumed to be snow coveredif time lapse images indicated that at least 75% of the ground in the imagewas covered in snow. This was a subjective classification made by visualinspections of the time lapse images. The modelling results suggest thatat the end of the early March event, the patchy snow cover observed fromthe time lapse images had a limited influence on throughflow temperatures.This issue reinforces the non-trivial challenge in defining when a rain-on-snow event occurs, either at the small headwater catchment scale studiedhere, or for larger catchments that span elevation ranges where snow mayonly cover a portion of the catchment or rain is falling at low elevationswhile snow is falling at higher elevations (Jones and Perkins, 2010).There are a number of limits of the hillslope throughflow model, al-though the model does highlight that a relatively simple conceptual modelthat represents dominant processes controlling subsurface temperatures(conduction and advection) can accurately predict throughflow tempera-tures. A key assumption is that the model assumes that hillslope through-flow temperatures are spatially uniform and that temperatures within thephreatic zone are uniform within each reservoir. Field observations suggestthat this is a reasonable assumption, especially during high runoff periods.The model does not resolve hillslopes dominated by deeper groundwatercontributions. In Griffith Creek, a known seep area did not appear to have asignificant influence on the stream thermal regime during winter periods.However, in other catchments that may be dominated by regional ground-water discharge, the model would need to be modified to include a reservoirrepresenting deeper groundwater flow paths.1115.4.2 Rain-on-snow eventsThe model evaluation for the individual rain-on-snow events showed amaximum increase in throughflow temperatures of 0.01 to 0.97 ◦C whenthe events were simulated with the snow cover removed. The modellingexercise that simulated a rain-on-ground event as a rain-on-snow eventsuggests that decreases of up to 2.3 ◦C, mostly associated with cool snowmeltand rain water inputs, could have been produced had snow been presentduring the event. Results from Chapter 3, which determined effectivethroughflow advective temperature calculated using the diagnostic heatbudget, suggested that throughflow temperatures during rain-on-groundevents may be on average 1.1 ◦C, and up to 3 to 4 ◦C, higher than during rain-on-snow events when controlling for air temperature (Figure 3.9). However,there was considerable scatter and much overlap between throughflowtemperatures during rain-on-ground and rain-on-snow events. In addition,there is up to ±4 ◦C uncertainty in some of the calculations of effectivethroughflow temperature from Chapter 3. The results from the modellingexercise fall within those ranges reported in Chapter 3.A possible explanation for the discrepancy between the modelled andempirical rain-on-snow versus rain-on-ground throughflow temperaturedifference is that the empirical analysis overestimated the cooling effectof rain-on-snow events. The overestimation is possible, since throughflowtemperatures are correlated in time. Indeed, the throughflow temperaturemodel suggests that the cooling effect of snow cover and associated rain-on-snow events can persist for a number of weeks after the event. Therefore,throughflow temperature response to a single event will be contingenton the type (rain-on-snow versus rain-on-ground) and characteristics (e.g.magnitude, duration, temperature) of the events that precede it. However,it is also possible that the model underestimated the cooling effect of rain-on-snow events due to errors in model input, model parameterization, ormodel structure.Errors in model parameterizations could also contribute to the modelunderestimating the cooling effect of rain-on-snow events. The model was112calibrated based on whether the throughflow temperature predictions fellwithin the 10th and 90th percentiles of observed throughflow temperaturesfor the full winter period, not specific events. Therefore, the calibrationprocedure may have favoured parameter sets that were less thermally re-sponsive to individual rain and rain-on-snow events, but consistently fellwithin the observed temperature range at the seasonal scale.Finally, errors in model structure could also account for errors in modelprediction. One mechanism not included in the model, that could accountfor cooler throughflow during rain-on-snow events is preferential flowpaths,particularly in the highly porous duff layer. The model assumes all waterentering the ground is warmed or cooled first to the duff layer temperature,then to the vadose zone temperature, before entering the phreatic zone. It ispossible that preferential flowpaths in the duff and soil rapidly deliver coldwater from the snow pack to the phreatic zone.In conclusion, the model does not appear to replicate the magnitude ofthroughflow cooling associated with rain-on-snow events relative to rain-on-ground events suggested by the empirical analysis, although there maybe issues of autocorrelation and up to ±4 ◦C uncertainty in the empiricalresults. However, the model does simulate a cooling influence and themodelling exercise suggests that the magnitude of cooling is related to theduration of the rain event, duration of the snow cover during the rain event,temperature of incoming rain, antecedent soil temperature, and snowpackdepth and spatial extent.113Chapter 6Development and evaluation ofa catchment scale streamtemperature model6.1 IntroductionComputer models provide a framework within which to synthesize ourknowledge of complex hydrologic systems (Blo¨schl, 2006). In addition,critical evaluation of model performance can improve our understanding ofthese systems (Beven, 2001). One approach to assessing our understandingof a system is to evaluate the temporal and geographic transferability ofa model by applying the model to a different time period or location thanwhen or where it was developed and calibrated. In Chapter 4, two primaryissues with current catchment-scale predictive stream temperature modelswere highlighted that limit their ability to properly simulate winter streamthermal conditions at headwater catchments similar to those found in theMalcolm Knapp Research Forest: 1) most models cannot account for therole of transient snow cover on the temperature of lateral throughflow, and2) of the models that do account for transient snow cover, they inaccuratelyestimate the temperature of lateral throughflow by up to 5 ◦C. Therefore, the114primary objective of this chapter is to develop and evaluate a conceptual-parametric catchment-scale stream temperature model that includes the roleof transient snow cover and lateral advection associated with throughflow.The model structure is an extension of the throughflow temperature modeldeveloped in Chapter 5.The catchment-scale stream temperature model is evaluated followingthe hierarchical approach proposed by Klemesˇ (1986), which includes testinga model’s geographic and climatic transposability. Transposability refers toa model’s ability to perform satisfactorily when applied to other catchmentsor climatic conditions for which it was not calibrated. Klemesˇ’ hierarchicalprocedure includes four levels of testing: 1) split-sample test, 2) proxy-basintest, 3) differential split-sample test, and 4) proxy-basin differential split-sample test. For this study, the stream temperature model was evaluatedat the first three levels of testing by calibrating the model to the 2011–2012winter period at Griffith Creek and testing it for the 2012–2013 winter periodat Griffith Creek (split-sample test), testing the model for both 2011–2012and 2012–2013 winter periods at East and South creeks using the calibratedparameters from Griffith Creek (proxy-basin test), as well as testing whetherthe model can reproduce the difference between pre- and post-harveststream temperatures at Griffith Creek (differential split-sample test).The differential split-sample test is based on empirical estimates of thedifference between pre- and post-harvest stream temperatures as deter-mined using a paired-catchment analysis by Guenther (2007); Guentheret al. (2014). He found that for the first year following harvesting, loggingappeared to have minimal impact on daily mean stream temperatures atthe outlet during winter. Post-logging winter stream temperatures werewithin±1 ◦C of the predicted winter temperatures had logging not occurred.However, post-logging daily mean stream temperatures were up to 2 ◦Cgreater during spring. The differential split-sample test conducted in thischapter is considered a pseudo-differential split-sample test and not a truedifferential split-sample test, because although the test is run under differ-ent conditions at Griffith Creek (pre-harvest versus post-harvest), streamtemperature was simulated for seven and eight years post-harvest, whereas115the empirical results against which the model was tested are from the firstyear post-harvest.6.2 MethodsIn this section the catchment-scale stream temperature model structure isdescribed, followed by the methodology and data sources used to calibrateand evaluate the model.6.2.1 Model developmentStructure and governing equationsThe stream temperature model adds stream channel water and heat balancecalculations to the hillslope throughflow temperature model described inChapter 4. First, the equations used to calculate the water balance of thestream channel are outlined, followed by the equations used to calculate theheat balance.The rate of change in water storage in the stream channel is calculatedas:dVcdt= Qlat2 −Q3 (6.1)where Vc is the volume of water stored in the stream channel (m3), Qlat2 isthe discharge from hillslope reservoir 2 (m3 s−1) and is calculated from thewater drainage rate, q2 (m s−1), from Chapter 4, and Q3 is the outflow fromthe catchment (m3 s−1), calculated as:Q3 = a(Vc −V0)b (6.2)where a and b are parameters and V0 is residual water storage in the channel(m3). Parameters a and b were calibrated using a fitted relationship betweenobserved discharge and total channel volumes calculated from the reachlength, mean stream widths and depths determined from field measure-ments. The residual water storage is included in the model to account for116storage in the channel (e.g. in pools) during periods when discharge fromthe catchment is zero, and is calculated as:V0 = L · w · d0 (6.3)where L is the reach length (m), w is the mean stream width during thestudy period (m), and d0 is a calibrated mean residual stream depth (m).The residual water storage was the only parameter calibrated beyond thosecalibrated for the hillslope throughflow temperature model. The calibrationof V0 is described below in Section 6.2.3.Stream temperature, Tw (◦C), is calculated as:Tw = Hc/(Vc · ρw · cw) (6.4)where Hc is the heat content in the channel (J), ρw is the density of water(kg m3), and cw is the specific heat capacity of water (J kg−1 ◦C−1). The rateof change in Hc is calculated as:dHcdt= Lw(Q∗ + Qe + Qh) + Qlat2ρwcwTp2 −Q3ρwcwTw (6.5)where L is the reach length (m), which was assumed constant, w is the streamwidth (m), which varied in relation to Q3 based on an empirical relationshipusing field measurements, Q∗, Qe and Qh are net radiation, latent andsensible heat fluxes (W m−2), respectively, and Tp2 is the phreatic zonetemperature of the hillslope reservoir 2 (◦C) and represents the temperatureof throughflow inputs to the stream. The stream surface energy fluxes (Q∗,Qe and Qh) were estimated using their respective set of equations outlinedin Chapter 3, with the exception that modelled Tw was used in place ofobserved Tw. Griffith Creek was treated as a single reach, although thereach has both forested and harvested riparian conditions with differingmicroclimates and canopy conditions (Chapter 3). Because of these twocontrasting riparian conditions, total incoming shortwave and longwaveradiation for the entire reach was calculated using proportional weightingof the forest and harvest reaches.117Equation 6.5 was integrated at six-minute intervals using a fourth orderRunge-Kutta scheme in the ‘deSolve’ package in R (Soetaert et al., 2010).The six-minute values of Tw were used to calculate hourly means.6.2.2 Data processing for East and South creeksThe model was developed and calibrated for Griffith Creek because of theextensive field data collected at that catchment. In order to transfer themodel to East and South creeks, where less field data were collected, certainassumptions and data processing were required, which are described in thissection.The empirical channel geometry relationships determined for GriffithCreek, which link channel volume and stream width and depth to discharge,were used at East and South creeks. Riparian canopy conditions and hemi-spherical images for modelling below-canopy incoming shortwave andlongwave radiation for East and South creeks were assumed to be similar tothose conditions at the forested section at Griffith Creek, since the riparianzones of both East and South creeks are dominated by mature trees.Above-stream air temperature and relative humidity were measured atEast and South creeks during the 2012–2013 winter periods. These variableswere used in the calculations of the surface energy exchanges. Relativehumidity was converted to vapour pressure. In order to estimate air tem-perature and vapour pressure for the 2011–2012 winter, hourly mean valuesof these variables for 2012–2013 were regressed against air temperature andvapour pressure at Griffith Creek. The statistical relationships establishedfor 2012-2013 were used, in conjunction with the 2011–2012 data from Grif-fith Creek, to predict air temperature and vapour pressure for the 2011–2012winter at East and South creeks.Catchment area, mean slope, and reach length were determined froma 5 m digital elevation model for MKRF (Table 6.1). The temperature atthe depth below the surface where it is assumed to be relatively constantthroughout the year (Tdeep) was determined by regressing mean daily airtemperatures monitored at East and South creek against mean daily air118Table 6.1: Parameters used in the stream temperature model for thethree study catchments.Catchment Elevation (m) Area (ha) Slope (◦) Reach length (m) Tdeep (◦C)Griffith 365-572 10.8 14.8 560 8.5East 295-555 44.3 10.6 1100 7.9South 175-320 19.7 12.6 880 9.5temperature at the MKRF headquarters. The relationship was then appliedto the preceding fifteen years of mean daily air temperatures at the MKRFheadquarters. From these fifteen year records, the long-term mean wascalculated. In addition, a standard geothermal temperature gradient of25 ◦C per 1000 m depth into the earth surface (Wood and Hewett, 1982) wasapplied to the long-term mean, since Tdeep was set at 15 m below the groundsurface.Hourly water inputs (rainfall and snowmelt) were measured with asnow lysimeter at the Griffith Creek forest site. No lysimeters were installedat East and South creeks. In order to use the water input data from GriffithCreek at East and South creeks, it was necessary to account for the differ-ences in snow regime and elevation between the three catchments. Thewater input data measured at Griffith was first corrected by substractingthe snowmelt from the total water inputs in order to estimate the rain-onlycomponent. Snowmelt was estimated using a simple energy budget ap-proach that included net radiation and ground heat fluxes. Turbulent energyexchanges were assumed to be negligible for forested areas (Beaudry andGolding, 1983) due to low observed wind speeds (mean of 0.5 and 0.6 m s−1during the two winter periods for the forest and harvest sites, respectively).Snowmelt (SM, m s−1) was estimated as:SM = Qm/(ρw · L f ) (6.6)where Qm is the energy available for snowmelt (W m−2), ρw is the densityof water (kg m−3), and L f is the latent heat of fusion (J kg−1). The energy119available for snowmelt (Qm) was calculated as:Qm = K∗ + L∗ + Qg (6.7)where K∗ is net shortwave radiation, L∗ is net longwave radiation, and Qgis the ground heat flux, all in W m−2. Net shortwave radiation (K∗) wascalculated as:K∗ = (1− α)K ↓ (6.8)where α is the albedo of snow and is assumed to be equal to 0.8 (Oke,1987), and K ↓ is the incoming solar radiation reaching the snow surfaceand was assumed equal to the modelled incoming solar calculated usingthe hemispherical images for the forested site (Chapter 3). Net longwaveradiation (L∗) was calculated as:L∗ = L ↓ −σεsnowT4snow (6.9)where L ↓ is the incoming longwave radiation and is assumed equal to themodelled incoming longwave radiation using the hemispherical images forthe forested site (Chapter 3), σ is the Stefan-Boltzmann constant, εsnow is theemissivity of snow and was assumed to be equal to 0.97 (Oke, 1987), andTsnow is the temperature of the snow surface (K) and was assumed to be273.15 K.The ground heat flux (Qg) was calculated as:Qg = kd[(Tsoil − Tsnow)/∆z] (6.10)where kd is the thermal conductivity of the soil duff layer (W m−1 ◦C) andwas assumed to be 0.2 based on the calibration of the throughflow tempera-ture model (Chapter 4), Tsoil is the temperature of soil measured at 0.05 mlocated near the forest lysimeter site, and ∆z is the distance between thesnow and the soil temperature measurement (0.05 m).The snowmelt was converted to mm h−1 and subtracted from the total120water inputs measured by the snow lysimeter at the Griffith Creek forest site.The resulting rainfall-only data were applied to East and South creeks (Fig-ure 6.1). The snowmelt contributions at East and South creeks were addedto the rainfall-only data for periods when a snowpack was present on Eastor South catchments as determined by the site-specific time lapse cameras.The snowmelt contributions at East and South creek were estimated usingthe same approach used above for estimating snowmelt at Griffith Creek,but substituting catchment-specific air temperatures for calculation of thenet longwave flux. No adjustment between the rainfall inputs from GriffithCreek were made for East Creek, since these catchments are at similar eleva-tions and have similar forest cover. The rainfall inputs from Griffith Creekwere adjusted using a proportionality factor of 0.81 based on throughfallmeasurements made by Donnelly-Makowecki and Moore (1999) at East andSouth creeks.6.2.3 CalibrationThe 30 behavioural parameter sets determined from the generalized likeli-hood uncertainty estimate (Beven and Binley, 1992) calibration of the hills-lope throughflow temperature model in Chapter 5 were used to simulatestream temperature at Griffith, East and South creeks. In addition, the depthused in the calculation of the residual channel storage was calibrated for the2011–2012 winter period at Griffith Creek. A sequence of depths betweenthe range of 0 to 0.15 m, at 0.01 m intervals, were used in the calculation ofV0. For each value of V0, the model was run 30 times (i.e. for each parameterset). Hourly predicted stream temperatures were compared to observedstream temperatures and five goodness of fit statistics were calculated: 1)Nash-Sutcliffe efficiency (NSE), 2) coefficient of determination (R2), 3) rootmean square error (RMSE), 4) mean absolute error (MAE), and 5) mean biaserror (MBE). A depth of 0.08 m maximized NSE and R2, and minimizedRMSE, MAE, MBE, and was therefore used in all subsequent analyses.121Figure 6.1: Total water inputs (snowmelt plus rainfall) measured at theGriffith Creek forest lysimeter for the 2011–2012 and 2012–2013winter periods plotted against estimated rainfall-only for whichthe snowmelt component was removed. The diagonal line is the1:1 line.6.2.4 EvaluationSplit-sample and proxy-basin testThe stream temperature model split-sample test involved predicting stream-flow and stream temperature for the 2012–2013 winter period at GriffithCreek. The proxy-basin test involved predicting streamflow and streamtemperature at East and South creeks for 2011–2012 and 2012–2013 winterperiods using the 30 behavioural parameter sets determined from the hill-slope throughflow temperature model. Model predictions of streamflowwere compared against observed streamflow using the NSE whereas streamtemperature predictions were evaluated using NSE, R2, RMSE, MAE, andMBE.122Pseudo-differential split-sample testThe ability of the stream temperature model to simulate pre-harvest winterstream temperature conditions at Griffith Creek was evaluated. Forest coverconditions during the study were characterized by partial retention harvestfor the lower portion of the catchment, which was carried out in autumn2004 (Chapter 2). The upper portion of the catchment was covered withmature trees. Pre-harvest conditions were simulated by applying the forestsite net radiation, snow cover, and above-stream microclimate conditions tothe entire catchment. Comparisons were made between simulated streamtemperature under this pre-harvest scenario and the current post-harvestconditions for the 2011–2012 and 2012–2013 winter periods.6.3 Results6.3.1 Conditions at Griffith, East and South creeksFigures 6.2 and 6.3 compare stream temperature, above-stream air temper-ature and vapour pressure, and unit discharge at Griffith, East and Southcreeks for the 2011–2012 and 2012–2013 winter periods. All three streamsexhibited similar temporal patterns in stream and air temperature, vapourpressure and unit discharge. South Creek stream temperature was generallyhigher than Griffith Creek, and Griffith Creek was generally warmer thanEast Creek, although all three creeks were usually within 1 to 2 ◦C of eachother. Similarly, air temperatures were generally highest at South and low-est at East Creek. Timing of streamflow response at the three catchmentswas similar. Unit discharge was similar between South and Griffith creeks,whereas East Creek had peakflows that were higher by about 1 to 3 mm hr−1than South and Griffith creeks.6.3.2 Split-sample and proxy-basin testFigures 6.4 to 6.9 show hourly observed and modelled stream temperatureand discharge at Griffith, East and South creeks for the 2011–2012 and2012–2013 winter periods. Figure 6.10 shows summary histograms of the123Figure 6.2: Stream temperature (Tw), air temperature (Ta), vapour pres-sure (kPa), and runoff (mm/hr) for the Griffith, East and Southcreeks for the 2011–2012 winter period. No above-stream air tem-perature and relative humidity observations were made at Eastand South creeks during 2011–2012.124Figure 6.3: Stream temperature (Tw), air temperature (Ta), vapour pres-sure (kPa), and runoff (mm/hr) for the Griffith, East and Southcreeks for the 2012–2013 winter period.125Figure 6.4: Griffith Creek observed and predicted hourly stream tem-perature and discharge the 2011–2012 calibration period. Blacklines are observed stream temperature or discharge, blue lines aremodelled stream temperature for the harvest snow regime, redlines are modelled stream temperature for the forest snow regime,and grey lines are modelled discharge. Horizontal lines indicateperiods of snow cover.goodness of fit statistics for discharge and stream temperature predictions.NSE for the discharge at South Creek were generally between 0.68 and 0.7.At East Creek, NSE dropped to between 0.41 and 0.63. The model generallyunderpredicted peakflows at East Creek, but captured the recession limbsand low flow conditions. For South Creek, modelled peak magnitudesfor the event hydrographs were similar to observed magnitudes, but themodelled recession limbs were more gradual than those for the observedhydrographs.Stream temperature predictions for all sites and years had RMSE below126Figure 6.5: Griffith Creek observed and predicted hourly stream tem-perature and discharge the 2012–2013 evaluation period. Blacklines are observed stream temperature or discharge, blue lines aremodelled stream temperature for the harvest snow regime, redlines are modelled stream temperature for the forest snow regime,and grey lines are modelled discharge. Horizontal lines indicateperiods of snow cover.1 ◦C. The model tended to overpredict stream temperature at East and Southcreeks, with MBE of around 0.25 ◦C. The stream temperature overpredictionswere more prevalent during the autumn, especially for the 2011–2012 winter.This is reflected in the lower NSE and R2, and higher RMSE, MAE andMBE for the 2011–2012 years (Figure 6.10). The overprediction at East andSouth creeks may be attributable to overprediction of modelled throughflowtemperatures at these two sites. Figure 6.11 shows modelled throughflowtemperatures for Griffith (under forest and harvest snow and microclimateconditions), East and South creeks. East Creek throughflow temperatures127Figure 6.6: East Creek observed and predicted hourly stream tempera-ture and discharge the 2011–2012 evaluation period. The horizon-tal line indicates periods of snow cover.are similar to throughflow temperatures at Griffith Creek, except duringNovember for 2011–2012, when throughflow temperatures at East Creek areabout 1 ◦C higher than at Griffith Creek. Modelled throughflow tempera-tures at South Creek are consistently about 2 ◦C higher than at Griffith andEast creeks.Griffith, East and South creeks have similar thermal regimes and there-fore, simply comparing model predictions to observed stream temperaturesdoes not provide a rigorous test of the model’s transposability. In order toprovide a more rigorous evaluation of the model, the ability of the modelto reproduce the small (1 to 2 ◦C) differences between Griffith, East andSouth creeks was evaluated. Figure 6.12 shows the relative differences be-tween observed and modelled stream temperature between Griffith and128Figure 6.7: East Creek observed and predicted hourly stream tempera-ture and discharge the 2012–2013 evaluation period. Horizontalline indicates periods of snow cover.East (top plots) and South (bottom plots) creeks for the two winter periods.The modelled relative difference between Griffith and East/South creeksgenerally match the observed relative differences. The exception being thefirst few months of the 2011–2012 period for the East Creek comparison,where East Creek was generally cooler than Griffith Creek by 0.5 to 1 ◦C, butthe model predicted that both streams had similar thermal regimes. Thisresult is consistent with the model overprediction for this period at EastCreek (Figure 6.6). For the remaining periods, the model was generallysuccessful at predicting the small relative differences in stream temperatureduring periods with and without snow cover.129Figure 6.8: South Creek observed and predicted hourly stream tem-perature and discharge the 2011–2012 evaluation period. Thehorizontal line indicates periods of snow cover.6.3.3 Pseudo-differential split-sample testFigures 6.13 and 6.14 compare stream temperature predictions for the cur-rent post-harvest conditions and simulated pre-harvest conditions. Thesimulated pre-harvest conditions were generally similar to current condi-tions, although pre-harvest stream temperatures tended to be lower thancurrent temperatures by as much as 0.5 ◦C. This difference is most evidentduring early autumn and spring, when incoming solar radiation is greaterand more of this energy flux will reach the stream surface under the currentpost-harvest conditions compared to pre-harvest conditions.130Figure 6.9: South Creek observed and predicted hourly stream tem-perature and discharge the 2012–2013 evaluation period. Thehorizontal line indicates periods of snow cover.6.4 Discussion6.4.1 StreamflowEvaluation of the discharge component of the stream temperature model forGriffith, East and South creeks resulted in NSE ranging between 0.42 and0.73. Donnelly-Makowecki and Moore (1999) used a similar three reservoirmodel to predict discharge at East and South creeks using the hierarchicaltesting framework proposed by Klemesˇ (1986). Their model evaluationswere similar to those in this study, with NSE in the range of 0.59 to 0.87.There are a number of error sources that could explain the differences inmodelled and observed discharge. Two key potential sources of error areerrors in precipitation data and the discharge rating curves. Spatial vari-131Figure 6.10: Histograms of goodness of fit statistics for model evalua-tion at Griffith (GC), East (EC) and South (SC) creeks for 2011–2012 (11/12) and 2012–2013 (12/13). Q is discharge; Tw is streamtemperature; NSE is Nash-Sutcliffe efficiency; R2 is the correla-tion coefficient; RMSE is root mean square error; MAE is meanabsolute error; and MBE is mean bias error.132Figure 6.11: Modelled throughflow temperatures at Griffith (forest andharvest conditions), East and South creeks for 2011–2012 and2012–2013 winter periods.ability in rainfall was accounted for by applying a correction factor basedon precipitation and throughfall data collected by Hetherington (1976) andDonnelly-Makowecki and Moore (1999). Using a single precipitation inputfor a catchment may not be representative; however, the catchments studiedare relatively small and for the largest, East Creek, Donnelly-Makowecki andMoore (1999) found only minimal differences in throughfall rates betweenthe low and high elevations. Since water inputs (snowmelt and rainfall)were measured only in and near Griffith Creek, differences in the snowregime across the catchments may have led to input data errors. Differencesin snow cover between the three catchments were accounted for by usingthe snowmelt energy budget approach. A potential source of error is that133Figure 6.12: Relative comparison of observed (black lines) and mod-elled (grey bands) stream temperature differences between Grif-fith and East or South creeks for 2011–2012 and 2012–2013. Hori-zontal lines at the top of each plot indicate periods of snow coverfor Griffith (black) and either East or South (blue).134Figure 6.13: Comparison of hourly modelled stream temperature re-sponse to pre-harvest conditions for 2011–2012. Upper plotshows absolute temperatures between baseline (post-harvest)conditions (black lines) and modelled pre-harvest conditions(blue lines). The lower plot shows post-harvest temperaturessubstracted from pre-harvest temperatures. Negative valuesindicate pre-harvest temperatures are lower than post-harvesttemperatures.135Figure 6.14: Comparison of hourly modelled stream temperature re-sponse to pre-harvest conditions for 2012–2013. Upper plotshows absolute temperatures between baseline (post-harvest)conditions (black lines) and modelled pre-harvest conditions(blue lines). The lower plot shows post-harvest temperaturessubstracted from pre-harvest temperatures. Negative valuesindicate pre-harvest temperatures are lower than post-harvesttemperatures.136the turbulent energy exchanges were not included in the snowmelt model.Latent and sensible heat exchanges can be important for snowmelt in clear-ings, particularly during rain-on-snow events (Beaudry and Golding, 1983;Berris and Harr, 1987); however, low wind speeds, such as those measuredunder mature forest canopies, will greatly inhibit the magnitude of theseturbulent fluxes (Marks et al., 1998).A likely source of error for the streamflow predictions is the parameteri-zation of the discharge rating curves for East and South creeks. Establishedrating curves provided for these creeks were supplied by Dr. Michael Feller,who had previously maintained these sites since the 1970s. Replicatedmanual discharge measurements were made at Griffith Creek and an inde-pendent rating curve was built that suggested that the rating curve providedby Dr. Feller for Griffith Creek was in error at times by up to 50%. Becauseof the extensive field campaign at Griffith Creek, there were insufficientresources to conduct additional manual stream flow gauging at East andSouth creeks to establish the error of those rating curves. Therefore, uncer-tainties in the rating curves for those streams may be a cause for differencesbetween modelled and observed discharge.6.4.2 Split-sample and proxy-basin testsAlthough the model generally predicted the pattern and magnitude ofstream temperature accurately, there are a number of sources of model error.The model does not include some known stream energy exchange processes.Heat advected by channel-intercepted precipitation was not included butis generally considered to be a minor heat budget term (Webb and Zhang,1999). However, channel interception of snowfall may be important, partic-ularly during low flow periods, when the energy consumed for the phasechange from solid to liquid could reduce temperatures. This mechanismcould be an explanation for the spiky decreases in stream temperatureduring the snow accumulation period for mid-December 2012 at GriffithCreek. Although snow also fell at East and South creek during this event,the spiky decreases in stream temperature at these two creeks were not137as pronounced. This could be explained by the more open above-streamcanopy conditions for the lower portion of Griffith Creek, which wouldallow greater channel interception than at East and South creeks, whichhave more extensive above-stream canopy closure. Conclusions regardingthe effect of snowfall into the channel on stream temperature are tentativebecause reliable measurements and precise timing of snowfall were notavailable. Other processes not represented in the model include bed heatconduction, hyporheic exchange, stream friction, or groundwater inputsfrom deep, thermally-stable sources. These heat budget terms were found tobe relatively minor for Griffith Creek during winter periods (Chapter 3), andappear to also be minor terms for East and South creeks, since predictionsof stream temperature with RMSE of less than 1 ◦C were made withoutincluding these terms. These processes, in particular groundwater inputsand hyporheic exchange, may be important to model in order to accuratelypredict stream temperature during summer low flow periods (Guenther,2007). Indeed, not including the thermal moderating effect of groundwaterdischarge or hyporheic exchange may be one explanation for overpredictionof stream temperature at East and South Creeks during the autumn andearly winter periods.Two additional sources of error for the stream temperature predictionsconcern model parameterization and input data. It was assumed that param-eterizations of the model for Griffith Creek hold for East and South creeks.In particular, the empirical relations used to determine the channel volumebased on discharge was assumed to be the same at Griffith Creek as Eastand South creeks. Errors in the water input data may have lead to errorsin model predictions. Determining snow cover with the use of time lapsecameras has sources of uncertainty. In particular, only the areas around thecatchment outlets were monitored. In this high relief region, snow may havebeen present at the upper elevations when snow had completely meltedat lower elevations. However, these are small catchments and for largesnowfall events, which will have the strongest influence on subsurface andstream temperatures, the time lapse camera approach is likely adequate. Asdiscussed above, errors in discharge predictions will have consequences for138the stream temperature predictions, since errors in the volume of water inthe channel will result in incorrect simulated stream temperatures even if thesurface energy fluxes and throughflow temperatures are accurately deter-mined. As long as the throughflow temperatures are relatively accurate, thissource of error will be more of a concern during low flow periods. Errorsin the reach length calculated from the DEM will influence the accuracy ofthe stream temperature predictions, by over- and under-representing themagnitude of energy added or lost at the stream surface.Despite these various sources of error, the model appears to capturethe energy exchange processes that dominate the thermal regimes of thesestreams, namely throughflow advection during high flows and surfaceenergy fluxes during low flow periods. Therefore, the improved under-standing of winter stream thermal regimes, gained through the detailedfield work conducted at Griffith Creek, can be generalized to other head-water catchments in the region. What is particularly promising is that themodel is able to replicate the small differences in stream temperature be-tween Griffith, East and South Creeks, during periods with and withoutsnow on the catchments.6.4.3 Pseudo-differential split-sample testGuenther (2007) found that at the outlet for the first year following har-vesting, logging appeared to have minimal impact on daily mean streamtemperatures during winter. However, post-logging daily mean streamtemperatures were up to 2 ◦C greater during spring. Modelled streamtemperature results from this study are generally consistent with Guenther(2007), although pre-harvest simulations were only about 0.5 ◦C lower thanpost-harvest conditions during spring. The smaller modelled response, com-pared to the empirical findings, may be due to some riparian vegetationrecovery. The 2011–2012 and 2012–2013 winters looked at in this studyare seven and eight years, respectively, since harvesting occurred, whichis consistent with other studies that have found some recovery within 5 to10 years post-harvest (Moore et al., 2005a). Another possible explanation139for the differences in response for the spring period is that increases inthroughflow temperature as a result of post-harvest increases in solar radia-tion reaching the soil surface were not considered. Guenther (2007) foundsummer increases in lateral inflow temperature of about 2 ◦C followingharvesting. It is not certain if or what the increases in lateral inflow duringthe winter might be.140Chapter 7ConclusionsThis study investigated the dominant energy exchange processes controllingwinter stream thermal regimes in headwater catchments located in the rain-on-snow zone of the Pacific Northwest. The following sections providea summary of important findings, implications of climate and land coverchanges, and opportunities for future research.7.1 Summary of important findingsIn Chapter 3, a diagnostic energy budget analysis was employed to deter-mine the dominant controls on winter stream temperature for two winters.The energy budget analysis highlighted the role of advective fluxes associ-ated with hillslope throughflow as being the dominant control on winterthermal regimes. During low flow periods energy exchanges occurringat the stream surface were of similar magnitude to the advective fluxesassociated with throughflow. It was found that stream temperatures associ-ated with rain-on-snow events were lower than rain-on-ground events aftercontrolling for air temperature. Therefore, predictive winter stream tempera-ture modelling for these types of catchments needs to consider throughflowadvection, surface energy exchanges and the role of transient snow cover.Chapter 4 focused on hillslope processes and provides the first study todocument spatiotemporal variability of throughflow temperatures. Field141observations suggested that throughflow temperatures were variable overtime, ranging between 2 and 12 ◦C during the winter period, but relativelyinvariable in space, particularly during rain events. Existing approaches tosimulate throughflow temperatures embedded in predictive stream temper-ature models were evaluated against field observations. These approacheseither could not represent the role of transient snow cover or over- or under-estimated throughflow temperatures by up to 5 ◦C. Therefore, these modelsmay not be suitable for predicting stream temperatures for catchmentsthat have stream thermal regimes influenced by throughflow inputs andrain-on-snow events.In Chapter 5, a conceptual-parametric throughflow temperature modelthat is computationally efficient was developed for the study site. Modelcalibration accounted for equifinality and uncertainty, and the model wasevaluated against field observations. The model successfully made pre-dictions of throughflow temperature at the study catchment between 75and 78% of the time, and modelled temperatures were never greater than0.7 ◦C outside the 10th and 90th percentiles of observed temperatures. Themodel was also tested for its ability to reproduce differences in throughflowtemperatures between rain-on-snow and rain-on-ground events. Althoughthe model successfully predicted lower throughflow temperatures duringrain-on-snow events compared to rain-on-ground events, the magnitudewas lower than was suggested by the empirical results from Chapter 3 butwas similar to the differences in the means and medians of the differentevent types.In Chapter 6, a predictive stream temperature model was built by extend-ing the hillslope throughflow model developed in Chapter 5 and includingenergy exchange processes for the stream channel. It was necessary to de-velop a new model since existing models were unable to simulate advectionassociated with throughflow inputs accurately. The stream temperaturemodel was evaluated at the primary study catchment, as well as two otherheadwater catchments in the Malcolm Knapp Research Forest. The modelprovided estimates of stream temperature at the three catchments withRMSE of less than 1 ◦C, supporting the empirical analysis (Chapter 3) which142showed that advection associated with throughflow inputs is an importantstream temperature control for headwater catchments in this region. Streamtemperature predictions made for pre-harvest forest cover conditions wereconsistent with empirical findings for Griffith Creek, in that winter streamtemperatures did not show a strong response to forest harvesting (Guenther,2007).7.2 Implications of climate and land cover changeson winter stream temperatureGiven that climate change projections indicate increased winter air tempera-tures in the Pacific Northwest (Mote and Salathe, 2010) and thus a decreasein transient snow cover, at least at lower elevations, the implication of thisstudy is that climate warming should generate higher winter stream temper-atures due to both the increased rain temperature and the reduced coolingeffect of snow cover. Persistently higher winter stream temperatures couldhave a profound effect on the rate of growth and development of salmonidembryos (e.g. Holtby, 1988) and on the rate and timing of invertebrate emer-gence. Earlier invertebrate emergence could have implications for riparianpredators, such as birds, who depend on this food source (Baxter et al., 2005;Richardson et al., 2005). Studies that have examined harvesting-inducedwinter stream temperature increases suggest that higher winter stream tem-peratures may promote salmonid growth; however, growth is also highlydependent on food availability (Holtby, 1988; Leach et al., 2012). Whereasharvesting-induced stream temperature increases are associated with in-creased solar radiation reaching the stream due to forest canopy removal,stream warming in the context of this study would be driven by warmingdue to advective fluxes without an increase in solar radiation. This differ-ence may be important since increases in solar radiation reaching the streamdue to harvesting could generate increases in primary production, withsubsequent food web linkages to invertebrates and salmonids, which maynot occur under warming generated by alterations in the advective fluxes(warmer rain and decreased snow cover).143There is less certainty in predicted precipitation changes due to climatechange than for air temperature changes; however, the general predictionis an increase in autumn and winter precipitation for the PNW (Rodenhuiset al., 2007; Mote and Salathe, 2010). An increase in precipitation over thewinter months would likely result in lateral advective fluxes associatedwith runoff being an even more dominant control on stream temperature,compared to vertical energy fluxes, than documented here. However, thefrequency of precipitation events may have more influence on the relativeimportance of the lateral advective fluxes than magnitude of precipitationevents. More frequent precipitation events will sustain runoff to the streamfor longer periods of time, thereby reducing the importance of verticalenergy exchanges on stream temperature. In contrast, high magnitude butless frequent precipitation events may result in longer periods of between-event low flows and increased importance of vertical energy exchanges.Land cover changes due to forest harvesting can have a considerableimpact on snow accumulation and melt in rain-on-snow zones (Berris andHarr, 1987; Marks et al., 1998; Jones, 2000; Hudson, 2001). There were dif-ferences in snow accumulation and melt timing between the forested andharvested areas of Griffith Creek; however, these differences were modest,likely a result of the partial retention harvest treatment, which reducedcanopy closure from about 95% pre-harvest to 80% post-harvest. Findingsfrom this study suggest that snow and runoff responses to forest harvest-ing could have greater impacts on winter stream temperature than surfaceenergy exchange responses to riparian forest removal. This study foundthat surface energy fluxes are relatively minor stream temperature controls,except during recession periods, and were relatively insensitive to ripariancanopy structure, as indicated by only slight differences in energy fluxesbetween forest and harvest reaches. In addition, previous studies fromMKRF found relatively small winter stream temperature responses to forestharvesting, compared to responses in spring and summer (Gomi et al., 2006;Guenther et al., 2014). Management activities focused on riparian vegetation,such as buffer strips (Gomi et al., 2006) or afforestation (Webb and Crisp,2006), will be less effective at influencing winter stream temperatures than144for summer temperatures. In order to effectively manage winter streamtemperatures, management activities may need to focus on catchment-scaledecisions that impact snow dynamics and hillslope runoff processes. Onetentative inference drawn from this study that would require further re-search is that forest harvesting treatments that reduce canopy closure morethan observed in this study may result in increased frequency and durationof transient snow cover (Hudson, 2001) with the potential to reduce winterstream temperature.Although this study was conducted in the Pacific Northwest of NorthAmerica, the findings may have broader geographic importance for under-standing winter stream temperature dynamics since rain-on-snow events oc-cur in many temperate regions, including New Zealand (Moore and Owens,1984) and Europe (Garvelmann et al., 2013). Futhermore, widespread shiftsin hydrologic regimes from snow-dominated to rain-dominated and mixedregimes (Barnett et al., 2005; Regonda et al., 2005; Knowles et al., 2006) mayresult in more stream systems experiencing reduced surface ice cover, higherwinter flows and increased frequency of rain-on-snow events. These hydro-logic regime shifts may result in winter thermal regimes being dominatedby advective fluxes associated with runoff for catchments that previouslydid not experience high winter runoff.7.3 Suggestions for future researchThis research has improved the understanding of the dominant energyexchanges controlling stream thermal regimes during winter in the rain-on-snow zone of the Pacific Northwest. Despite the extensive field observationsand modelling conducted herein, the conclusions drawn from this work maybe limited to the geographic context of this research. Therefore, conductingstream energy budget studies during winter in other environments will helpto further understand the dominant controls on winter stream temperaturesand how they might respond to climate and land cover changes. One suchformalized approach would be to test the hillslope throughflow model onother catchments outside of the Malcolm Knapp Research Forest. Some145additional directions for further study are explored in the following section.Although the catchments in this study had some snow cover during thewinter and rain-on-snow events did occur, the majority of the snowpackswere relatively shallow and most of the corresponding rain-on-snow eventswere small compared to rain-on-snow events documented elsewhere (e.g.Beaudry and Golding 1983; Marks et al. 1998). The relatively short durationof snow cover at South Creek, the catchment in this study with the lowest el-evation, suggests that these study catchments may be in the lower portion ofthis region’s rain-on-snow zone. In addition, the impacts of warmer wintersand diminishing snowpacks may already be occurring for the study catch-ments since modelled snowpack depth was relatively shallow for the 1997 to2013 period compared to historical modelled snowpacks since 1962 (Figure3.11). An interesting extension of this work would be to study catchmentsat elevations that are located higher in the rain-on-snow zone, which willhave deeper and longer duration snowpacks, and subsequently larger rain-on-snow events (Brunengo, 2012). The event-based rain-on-snow modellingexercises in Chapter 5 suggest that larger rain-on-snow events could resultin more pronounced stream temperature depressions than examined in thisstudy. Therefore, there may be considerable impacts to aquatic ecosystemsunder a warming climate for higher elevation catchments where histori-cally deeper snowpacks are reduced or removed and stream temperaturesincrease.The lower portion of Griffith Creek was harvested in 2004 under a par-tial retention approach that resulted in a 14% reduction in canopy closure(Guenther et al., 2012). This harvesting approach resulted in modest differ-ences between snow accumulation and melt for the harvested and forestedportions of Griffith Creek. It would be valuable if future research exploredwinter stream temperature response to different harvesting approaches thatresult in greater reduction in canopy closure and therefore, greater changesin snow dynamics between forested and harvested sites (Hudson, 2001).Snow accumulation is generally greater in clearings than under intact forestcanopies (Golding and Swanson, 1986). In addition, snow melt rates canalso be greater in clearings than in forests due to increased melt generated146by turbulent energy exchanges, particularly during rain-on-snow events(Marks et al., 1998). Therefore, a potential question worth addressing iswhether winter stream temperature might actually decrease with more ex-tensive harvesting due to greater snow accumulation and associated largerand more frequent rain-on-snow events.The number of stream energy budget studies conducted for summerperiods far exceeds the number conducted for winter periods. The generalconsensus from these summer studies is that net radiation dominates thestream thermal regime (e.g. Brown 1969; Hannah et al. 2008). However,there is a growing body of research that suggests advection associated withthroughflow inputs may also be an important energy budget componentduring summer periods, particularly during rain events (St-Hilaire et al.,2003; Hebert et al., 2011; Arismendi et al., 2012; Mayer, 2012). Measurementsof throughflow temperatures during summer periods would be valuable tofurther evaluate the ability of stream temperature models to estimate thisadvective heat flux.This study focused on the thermal regimes of small headwater catch-ments. Understanding the thermal regimes of headwater streams is impor-tant because headwater systems can make up almost 80% of total streamlength in many drainage networks (Shreve, 1969; Sidle et al., 2000) and pro-vide critical habitats for various organisms sensitive to stream temperature.An important question for further research is how the influence of lateral ad-vection changes as one moves down the channel network, especially becauselarger streams provide habitat for a broader range of aquatic species. Asdrainage area and stream order increase, streams generally become wider,providing a greater surface area for energy exchange, and the relative changein discharge per unit channel length should decrease. These considerationssuggest that lateral advection might become less important as catchmentscale increases. In addition, temperature signals at downstream locationswould be complicated by the merging of tributaries with varying elevationdistributions in their catchment areas. Further field research and modellingeffort are required to understand the scale-dependence of winter streamtemperature dynamics in coastal catchments.147ReferencesAllen, D., Dietrich, W. E., Baker, P. F., Ligon, F., and Orr, B. 2007.Development of a mechanistically-based, basin-scale stream temperaturemodel: applications to cumulative effects modeling. In Standiford, R. B.,Giusti, G. A., Valachovic, Y., Zielinski, W. J., and Furniss, M. J., editors,Proceedings of the Redwood Region Forest Science Symposium: What Does theFuture Hold?, pages 11–24. Forest Service, US Department of Agriculture,Albany, California.Anderson, M. G. and Burt, T. P. 1978. The role of topography in controllingthroughflow generation. Earth Surface Processes, 3:331–344.Anderson, M. P. 2005. Heat as a ground water tracer. Ground Water,43(6):951–968.Arismendi, I., Johnson, S. L., Dunham, J. B., and Haggerty, R. 2013.Descriptors of natural thermal regimes in streams and theirresponsiveness to change in the Pacific Northwest of North America.Freshwater Biology, 58(5):880–894.Arismendi, I., Johnson, S. L., Dunham, J. B., Haggerty, R., andHockman-Wert, D. 2012. The paradox of cooling streams in a warmingworld: Regional climate trends do not parallel variable local trends instream temperature in the Pacific continental United States. GeophysicalResearch Letters, 39(10):L10401.Arnold, N. S., Willis, I. C., Sharp, M. J., Richards, K. S., and Lawson, W. J.1996. A distributed surface energy-balance model for a small valleyglacier. I. Development and testing for Haut Glacier d’Arolla, Valais,Switzerland. Journal of Glaciology, 42(140):77–89.Arrigoni, A., Poole, G., Mertes, L., O’Daniel, S., Woessner, W., and Thomas,S. 2008. Buffered, lagged, or cooled? Disentangling hyporheic influences148on temperature cycles in stream channels. Water Resources Research,44(9):W09418.Barnett, T. P., Adam, J. C., and Lettenmaier, D. P. 2005. Potential impacts ofa warming climate on water availability in snow-dominated regions.Nature, 438(17):303–309.Bartholow, J. 2005. Recent water temperature trends in the Lower KlamathRiver, California. North American Journal of Fisheries Management,25(1):152–162.Battin, J., Wiley, M. J., Ruckelshaus, M. H., Palmer, R. N., Korb, E., Bartz,K. K., and Imaki, H. 2007. Projected impacts of future climate change onsalmon habitat restoration actions in a Puget Sound river. Proceedings ofthe National Academy of Sciences, 104:6720–6725.Baxter, C., Hauer, F. R., and Woessner, W. W. 2003. Measuringgroundwater-stream water exchange: new techniques for installingminipiezometers and estimating hydraulic conductivity. Transactions ofthe American Fisheries Society, 132(3):493–502.Baxter, C. V., Fausch, K. D., and Saunders, W. C. 2005. Tangled webs:reciprocal flows of invertebrate prey link streams and riparian zones.Freshwater Biology, 50(2):201–220.Beaudry, P. and Golding, D. L. 1983. Snowmelt during rain on snow incoastal British Columbia. In Proceedings of the 51st Annual Meeting of theWestern Snow Conference, pages 55–66. Colorado State University, FortCollins, Colorado.Benyahya, L., Caissie, D., Satish, M. G., and El-Jabi, N. 2012. Long-waveradiation and heat flux estimates within a small tributary in CatamaranBrook (New Brunswick, Canada). Hydrological Processes, 26:475–484.Berris, S. N. and Harr, R. D. 1987. Comparative snow accumulation andmelt during rainfall in forested and clear-cut plots in the westernCascades of Oregon. Water Resources Research, 23(1):135–142.Beschta, R. L., Bilby, R. E., Brown, G., Holtby, L. B., and Hofstra, T. D. 1987.Stream temperature and aquatic habitat: Fisheries and forestryinteractions. In Salo, E. O. and Cundy, T. W., editors, Streamsidemanagement: Forestry and fishery interactions, number 57, pages 191–232.University of Washington, Institute of Forest Resources, Seattle,Washington.149Beven, K. 2001. On hypothesis testing in hydrology. Hydrological Processes,15:1655–1657.Beven, K. and Binley, A. 1992. Future of distributed models: Modelcalibration and uncertainty prediction. Hydrological Processes,6(3):279–298.Beven, K. and Freer, J. 2001. Equifinality, data assimilation, and uncertaintyestimation in mechanistic modelling of complex environmental systemsusing the GLUE methodology. Journal of Hydrology, 249:11–29.Bevington, P. R. and Robinson, K. D. 2003. Data reduction and error analysisfor the physical sciences. Boston: McGraw-Hill, third edition.Bicknell, B. R., Imhoff, J. C., Kittle, J. L., Jobes, T. H., and Donigian, A. S.2001. Hydrological Simulation Program-Fortran (HSPF) user’s manual forrelease 12. US Environmental Protection Agency, National ExposureResearch Laboratory, Athens, GA.Birkinshaw, S. J. and Webb, B. 2010. Flow pathways in the Slapton Woodcatchment using temperature as a tracer. Journal of Hydrology,383(3-4):269–279.Black, T. A., Chen, J. M., Lee, X., and Sagar, R. M. 1991. Characteristics ofshortwave and longwave irradiances under a Douglas-fir forest stand.Canadian Journal of Forest Research, 21(7):1020–1028.Blo¨schl, G. 2006. Hydrologic synthesis: Across processes, places, and scales.Water Resources Research, 42:W03S02.Bogan, T., Mohseni, O., and Stefan, H. 2003. Streamtemperature-equilibrium temperature relationship. Water ResourcesResearch, 39(9):1245.Boughton, D. A., Hatch, C., and Mora, E. 2012. Identifying distinct thermalcomponents of a creek. Water Resources Research, 48:W09506.Boyd, M. and Casper, B. 2003. Analytical Methods of Dynamic Open ChannelHeat and Mass Transfer: Methodology for the Heat Source Model Version 7.0.Oregon Department of Environmental Quality, Portland, Oregon.Braithwaite, R. J. and Olesen, O. B. 1990. A simple energy-balance model tocalculate ice ablation at the margin of the Greenland ice sheet. Journal ofGlaciology, 36(123):222–228.150Bremicker, M. 2000. Das Wasserhaushaltsmodell LARSIM -Modellgrundlagen and Anwendungsbeispiele. Freiburger Schriften zurHydrologie, Band 11: Freiburg, 119.Briggs, M. A., Lautz, L. K., McKenzie, J. M., Gordon, R. P., and Hare, D. K.2012. Using high-resolution distributed temperature sensing to quantifyspatial and temporal variability in vertical hyporheic flux. WaterResources Research, 48:W02527.Brock, B. and Arnold, N. S. 2000. A spreadsheet-based (Microsoft Excel)point surface energy balance model for glacier and snow melt studies.Earth Surface Processes and Landforms, 25:649–658.Brookfield, A. E., Sudicky, E. A., Park, Y. J., and Conant Jr, B. 2009. Thermaltransport modelling in a fully integrated surface/subsurface framework.Hydrological Processes, 23(15):2150–2164.Brown, G. W. 1969. Predicting temperatures of small streams. WaterResources Research, 5(1):68–75.Brown, L. E. and Hannah, D. M. 2007. Alpine stream temperature responseto storm events. Journal of Hydrometeorology, 8(4):952–967.Brown, L. E. and Hannah, D. M. 2008. Spatial heterogeneity of watertemperature across an alpine river basin. Hydrological Processes,22(7):954–967.Brown, L. E., Hannah, D. M., and Milner, A. M. 2007. Vulnerability of alpinestream biodiversity to shrinking glaciers and snowpacks. Global ChangeBiology, 13(5):958–966.Brown, R. S., Hubert, W. A., and Daly, S. F. 2011. A primer on winter, ice,and fish: what fisheries biologists should know about winter iceprocesses and stream-dwelling fish. Fisheries, 36(1):8–26.Brunengo, M. J. 2012. Where is the rain-on-snow zone in the West-CentralWashington Cascades?: Monte Carlo simulation of large storms in theNorthwest. PhD thesis, Portland State University.Bulliner, E. and Hubbart, J. A. 2013. An improved hemisphericalphotography model for stream surface shortwave radiation estimationsin a central U.S. hardwood forest. Hydrological Processes, 27:3885–3895.151Buttle, J. M. 2006. Mapping first-order controls on streamflow fromdrainage basins: the T3 template. Hydrological Processes, 20(15):3415–3422.Caissie, D., Satish, M. G., and El-Jabi, N. 2007. Predicting watertemperatures using a deterministic model: Application on MiramichiRiver catchments (New Brunswick, Canada). Journal of Hydrology,336:303–315.Campbell, G. S. and Norman, J. M. 1998. An introduction to environmentalbiophysics. Springer-Verlag, New York, 2nd edition edition.Cheng, J. D. 1988. Subsurface stormflows in the highly permeable forestedwatersheds of southwestern British Columbia. Journal of ContaminantHydrology, 3:171–191.Constantz, J. 2008. Heat as a tracer to determine streambed waterexchanges. Water Resources Research, 44:W00D10.Cozzetto, K., McKnight, D., Nylen, T., and Fountain, A. 2006. Experimentalinvestigations into processes controlling stream and hyporheictemperatures, Fryxell Basin, Antarctica. Advances in Water Resources,29(2):130–153.Danehy, R. J., Colson, C. G., and Duke, S. D. 2010. Winter longitudinalthermal regime in four mountain streams. Northwest Science,84(1):151–158.Deardorff, J. W. 1978. Efficient prediction of ground surface temperatureand moisture, with inclusion of a layer of vegetation. Journal ofGeophysical Research, 83(C4):1889–1903.Deitchman, R. and Loheide, S. P. 2012. Sensitivity of thermal habitat of atrout stream to potential climate change, Wisconsin, United States.Journal of the American Water Resources, 48:1091–1103.Donnelly-Makowecki, L. M. and Moore, R. D. 1999. Hierarchical testing ofthree rainfall-runoff models in small forested catchments. Journal ofHydrology, 219(3-4):136–152.Durance, I. and Ormerod, S. J. 2007. Climate change effects on uplandstream macroinvertebrates over a 25-year period. Global Change Biology,13(5):942–957.152Durance, I. and Ormerod, S. J. 2009. Trends in water quality and dischargeconfound long-term warming effects on river macroinvertebrates.Freshwater Biology, 54(2):388–405.Ebersole, J., Wigington Jr, P., Baker, J., Cairns, M., Church, M., Hansen, B.,Miller, B., LaVigne, H., Compton, J., and Leibowitz, S. 2006. Juvenile cohosalmon growth and survival across stream network seasonal habitats.Transactions of the American Fisheries Society, 135(6):1681–1697.Erbs, D., Klein, S. A., and Duffie, J. A. 1982. Estimation of the diffuseradiation fraction for hourly, daily and monthly-average global radiation.Solar Energy, 28(4):293–302.Evans, E., McGregor, G., and Petts, G. 1998. River energy budgets withspecial reference to river bed processes. Hydrological Processes,12(4):575–595.Feller, M. C. and Kimmins, J. P. 1979. Chemical characteristics of smallstreams near Haney in southwestern British Columbia. Water ResourcesResearch, 15(2):247–258.Ficklin, D. L., Luo, Y., Stewart, I. T., and Maurer, E. P. 2012. Developmentand application of a hydroclimatological stream temperature modelwithin the Soil and Water Assessment Tool. Water Resources Research,48:W01511.Flerchinger, G. N. 2000. The Simultaneous Heat and Water (SHAW) Model:Technical Documentation, Tech. Rep. NWRC-2000-09. Technical report,USDA-ARS, Northwest Watershed Research Center, Boise, Idaho.Floyd, W. and Weiler, M. 2008. Measuring snow accumulation and ablationdynamics during rain-on-snow events: innovative measurementtechniques. Hydrological Processes, 22(24):4805–4812.Frazer, G. W., Canham, C. D., and Lertzman, K. P. 1999. Gap Light Analyser(GLA), Version 2.0: Imaging Software to Extract Canopy Structure and LightTransmission Indices from True-Colour Fisheye Photographs. User’s Manualand Program Documentation. Simon Fraser University, Burnady, BritishColumbia.Freeze, R. A. and Cherry, J. A. 1979. Groundwater. Prentice-Hall:Englewoods.153French, W. E., Vondracek, B., Ferrington Jr., L. C., Finlay, J. C., andDieterman, D. J. 2014. Winter feeding, growth and condition of browntrout Salmo trutta in a groundwater-dominated stream. Journal ofFreshwater Ecology, 29:187–200.Friberg, N., Bergfur, J., J., R., and Sandin, L. 2013. Changing Northerncatchments: is altered hydrology, temperature or both going to shapefuture stream communities and ecosystem processes? HydrologicalProcesses, 27:734–740.Gaffield, S., Potter, K., and Wang, L. 2005. Predicting the summertemperature of small streams in southwestern Wisconsin. Journal of theAmerican Water Resources Association, 41(1):25–36.Garner, G., Malcolm, I. A., Sadler, J. P., Millar, C. P., and Hannah, D. H. 2014.Inter-annual variability in the effects of riparian woodland onmicro-climate, energy exchanges and water temperature of an uplandScottish stream. Hydrological Processes, DOI: 10.1002/hyp.10223.Garvelmann, J., Pohl, S., and Weiler, M. 2013. From observation to thequantification of snow processes with a time-lapse camera network.Hydrology and Earth System Sciences, 17:1415–1429.Golding, D. L. and Swanson, R. H. 1986. Snow distribution patterns inclearings and adjacent forest. Water Resources Research, 22:1931–1940.Gomi, T., Moore, R. D., and Dhakal, A. S. 2006. Headwater streamtemperature response to clear-cut harvesting with different ripariantreatments, coastal British Columbia, Canada. Water Resources Research,42(8):W08437.Gravelle, J. A. and Link, T. E. 2007. Influence of timber harvesting on watertemperatures in a northern Idaho watershed. Forest Science, 53:189–205.Guenther, S. M. 2007. Impacts of partial-retention harvesting with no bufferon the thermal regime of a headwater stream and its riparian zone.Master’s thesis, University of British Columbia.Guenther, S. M., Gomi, T., and Moore, R. D. 2014. Stream and bedtemperature variability in a coastal headwater catchment: influences ofsurface-subsurface interactions and partial-retention forest harvesting.Hydrological Processes, 28:1238–1249.154Guenther, S. M., Moore, R. D., and Gomi, T. 2012. Riparian microclimateand evaporation from a coastal headwater stream and their response topartial-retention forest harvesting. Agriculture and Forest Meteorology,164:1–9.Haag, I. and Luce, A. 2008. The integrated water balance and watertemperature model LARSIM-WT. Hydrological Processes, 22:1046–1056.Hannah, D. M., Malcolm, I. A., Soulsby, C., and Youngson, A. F. 2004. Heatexchanges and temperatures within a salmon spawning stream in theCairngorms, Scotland: seasonal and sub-seasonal dynamics. RiverResearch and Applications, 20(6):635–652.Hannah, D. M., Malcolm, I. A., Soulsby, C., and Youngson, A. F. 2008. Acomparison of forest and moorland stream microclimate, heat exchangesand thermal dynamics. Hydrological Processes, 22(7):919–940.Haught, D. R. W. and van Meerveld, H. J. 2011. Spatial variation intransient water table responses: differences between an upper and lowerhillslope zone. Hydrological Processes, 25:3866–3877.Hebert, C., Caissie, D., Satish, M. G., and El-Jabi, N. 2011. Study of streamtemperature dynamics and corresponding heat fluxes within MiramichiRiver catchments (New Brunswick, Canada). Hydrological Processes,25:2439–2455.Herb, W. R. and Stefan, H. G. 2011. Modified equilibrium temperaturemodels for cold-water streams. Water Resources Research, 47:W06519.Hester, E., Doyle, M., and Poole, G. 2009. The influence of in-streamstructures on summer water temperatures via induced hyporheicexchange. Limnology and Oceanography, 54(1):355–367.Hetherington, E. D. 1976. Investigation of orographic rainfall in south coastalmountains of British Columbia. PhD thesis, University of British Columbia.Holtby, L. B. 1988. Effects of logging on stream temperatures in CarnationCreek British Columbia, and associated impacts on the coho salmon(Oncorhynchus kisutch). Canadian Journal of Fisheries and Aquatic Sciences,45(3):502–515.Hudson, R. O. 2001. Roberts Creek Study Forest: Preliminary Effects ofPartial Harvesting on Peak Streamflow in Two S6 Creeks. Technical155report, B.C. Forest Service, Vancouver Forest Region, Forest ResearchExtension Note EN-007 (March 2001), Nanaimo, British Columbia,Canada.Hutchinson, D. G. and Moore, R. D. 2000. Throughflow variability on aforested hillslope underlain by compacted glacial till. HydrologicalProcesses, 14(10):1751–1766.Huusko, A., Greenberg, L., Stickler, M., Linnansaari, T., Nyka¨nen, M.,Vehanen, T., Koljonen, S., Louhi, P., and Alfredsen, K. 2007. Life in the icelane: the winter ecology of stream salmonids. River Research andApplications, 23(5):469–491.Hvorslev, M. J. 1951. Time lag and soil permeability in ground-waterobservations. Technical report, Waterways Experimental Station, Corpsof Engineers, U.S. Army, Vicksburg, Mississippi.Imholt, C., Soulsby, C., Malcolm, I. A., and Gibbins, N. C. 2013. Influence ofcontrasting riparian forest cover on stream temperature dynamics insalmonid spawning and nursery streams. Ecohydrology, 6:380–392.Iqbal, M. 1983. An Introduction to Solar Radiation. Academic Press: Toronto.Isaacson, E. and Keller, H. B. 2012. Analysis of Numerical Methods. CourierDover Publications.Janisch, J., Wondzell, S. M., and Ehinger, W. J. 2012. Headwater streamtemperature: Interpreting response after logging, with and withoutriparian buffers, Washington, USA. Forest Ecology and Management,270:302–313.Jansson, P.-E. and Karlberg, L. 2001. Coupled heat and mass transfer model forsoil-plant-atmosphere systems. Department of Land and Water ResourceEngineering, Royal Institute of Technology, Stockholm.Jencso, K. G., McGlynn, B. L., Gooseff, M. N., Bencala, K. E., and Wondzell,S. M. 2010. Hillslope hydrologic connectivity controls ripariangroundwater turnover: Implications of catchment structure for riparianbuffering and stream water sources. Water Resources Research,46(10):W10524.Johnson, S. L. 2004. Factors influencing stream temperatures in smallstreams: substrate effects and a shading experiment. Canadian Journal ofFisheries and Aquatic Sciences, 61(6):913–923.156Johnson, S. L. and Jones, J. A. 2000. Stream temperature responses to forestharvest and debris flows in western Cascades, Oregon. Canadian Journalof Fisheries and Aquatic Sciences, 57(S2):30–39.Jones, J. A. 2000. Hydrologic processes and peak discharge response toforest removal, regrowth, and roads in 10 small experimental basins,western Cascades, Oregon. Water Resources Research, 36(9):2621–2642.Jones, J. A. and Perkins, R. M. 2010. Extreme flood sensitivity to snow andforest harvest, western Cascades, Oregon, United States. Water ResourcesResearch, 46(12):W12512.Kang, S., Kim, S., Oh, S., and Lee, D. 2000. Predicting spatial and temporalpatterns of soil temperature based on topography, surface cover and airtemperature. Forest Ecology and Management, 136(1-3):173–184.Kiffney, P. M., Bull, J. P., and Feller, M. C. 2002. Climatic and hydrologicvariability in a coastal watershed of southwestern British Columbia.Journal of the American Water Resources Association, 38(5):1437–1451.Kirkby, M. 1988. Hillslope runoff processes and models. Journal ofHydrology, 100:315–339.Klemesˇ, V. 1986. Operational testing of hydrological simulation models.Hydrological Sciences Journal, 31(1):13–24.Knowles, N., Dettinger, M. D., and Cayan, D. R. 2006. Trends in snowfallversus rainfall in the western United States. Journal of Climate,19:4545–4559.Kobayashi, D. 1985. Separation of the snowmelt hydrograph by streamtemperatures. Journal of Hydrology, 76(1):155–162.Kobayashi, D., Ishii, Y., and Kodama, Y. 1999. Stream temperature, specificconductance and runoff process in mountain watersheds. HydrologicalProcesses, 13(6):865–876.Krajina, V. 1965. Biogeoclimatic zones and classification of British Columbia.Ecology of Western North America, 1:1–17.Kurylyk, B. L., MacQuarrie, K. T. B., and Voss, C. I. 2014. Climate changeimpacts on the temperature and magnitude of groundwater dischargefrom shallow, unconfined aquifers. Water Resources Research,50:2013WR014588.157Langan, S., Johnston, L., Donaghy, M., Youngson, A., Hay, D., and Soulsby,C. 2001. Variation in river water temperatures in an upland stream over a30-year period. The Science of the Total Environment, 265(1-3):195–207.Lapham, W. W. 1989. Use of temperature profiles beneath streams todetermine rates of vertical ground-water flow and vertical hydraulicconductivity. Technical report, USGS Water-Supply Paper 2337, 35 p.Leach, J. A. 2008. Stream temperature dynamics following riparian wildfire:effects of stream-subsurface interactions and standing dead trees.Master’s thesis, University of British Columbia.Leach, J. A. and Moore, R. D. 2010. Above-stream microclimate and streamsurface energy exchanges in a wildfire-disturbed riparian zone.Hydrological Processes, 24(17):2369–2381.Leach, J. A. and Moore, R. D. 2011. Stream temperature dynamics in twohydrogeomorphically distinct reaches. Hydrological Processes,25(5):679–690.Leach, J. A., Moore, R. D., Hinch, S. G., and Gomi, T. 2012. Estimation offorest harvesting-induced stream temperature changes and bioenergeticconsequences for cutthroat trout in a coastal stream in British Columbia,Canada. Aquatic Sciences, 74(3):427–441.Loinaz, M. C., Davidsen, H. K., Butts, M., and Bauer-Gottwein, P. 2013.Integrated flow and temperature modelling at the catchment scale.Journal of Hydrology, 495:238–251.Loinaz, M. C., Gross, D., Unnasch, R., Butts, M., and Bauer-Gottwein, P.2014. Modeling ecohydrological impacts of land management and wateruse in the Silver Creek basin, Idaho. Journal of Geophysical ResearchBiogeosciences, 119:487–507.MacDonald, R. J., Boon, S., and Byrne, J. M. 2014a. A process-based streamtemperature modelling approach for mountain regions. Journal ofHydrology, 511:920–931.MacDonald, R. J., Boon, S., Byrne, J. M., and Silins, U. 2014b. A comparisonof surface and subsurface controls on summer temperature in aheadwater stream. Hydrological Processes, 28:2338–2347.158Malard, F., Mangin, A., Uehlinger, U., and Ward, J. 2001. Thermalheterogeneity in the hyporheic zone of a glacial floodplain. CanadianJournal of Fisheries and Aquatic Sciences, 58(7):1319–1335.Marce´, R. and Armengol, J. 2008. Modelling river water temperature usingdeterministic, empirical, and hybrid formulations in a Mediterraneanstream. Hydrological Processes, 22:3418–3430.Marks, D., Kimball, J., Tingey, D., and Link, T. 1998. The sensitivity ofsnowmelt processes to climate conditions and forest cover duringrain-on-snow: a case study of the 1996 Pacific Northwest flood.Hydrological Processes, 12(10-11):1569–1587.Marzadri, A., Tonina, D., and Bellin, A. 2013. Effects of streammorphodynamics on hyporheic zone thermal regime. Water ResourcesResearch, 49:2287–2302.Mayer, T. D. 2012. Controls of summer stream temperature in the PacificNorthwest. Journal of Hydrology, 475:323–335.McDonnell, J., Sivapalan, M., Vache´, K., Dunn, S., Grant, G., Haggerty, R.,Hinz, C., Hooper, R., Kirchner, J., Roderick, M., Selker, J., and Weiler, M.2007. Moving beyond heterogeneity and process complexity: a newvision for watershed hydrology. Water Resources Research, 43(7):W07301.McGlynn, B. L., McDonnell, J. J., and Brammer, D. D. 2002. A review of theevolving perceptual model of hillslope flowpaths at the Maimaicatchments, New Zealand. Journal of Hydrology, 257:1–26.McMahon, T. A., Peel, M. C., Lowe, L., Srikanthan, R., and McVicar, T. R.2013. Estimating actual, potential, reference crop and pan evaporationusing standard meteorological data: a pramatic synthesis. Hydrology andEarth System Sciences, 17:1331–1363.Mohseni, O., Stefan, H. G., and Erickson, T. R. 1998. A nonlinear regressionmodel for weekly stream temperatures. Water Resources Research,34:2685–2692.Moore, R. D. 1993. Application of a conceptual streamflow model in aglacierized drainage basin. Journal of Hydrology, 150(1):151–168.Moore, R. D. 2004. Introduction to salt dilution gauging for streamflowmeasurement Part II: Constant-rate injection. Streamline WatershedManagement Bulletin, 8(1):11–15.159Moore, R. D. and Owens, I. F. 1984. A conceptual runoff model for amountainous rain-on-snow environment, Craigieburn Range, NewZealand. Journal of Hydrology (New Zealand), 23:84–99.Moore, R. D., Spittlehouse, D. L., and Story, A. 2005a. Riparian microclimateand stream temperature response to forest harvesting: a review. Journal ofthe American Water Resources Association, 41(4):813–834.Moore, R. D., Sutherland, P., Gomi, T., and Dhakal, A. 2005b. Thermalregime of a headwater stream within a clear-cut, coastal British Columbia,Canada. Hydrological Processes, 19(13):2591–2608.Moore, R. D. and Thompson, J. C. 1996. Are water table variations in ashallow forest soil consistent with the TOPMODEL concept? WaterResources Research, 32(3):663–669.Morin, G. and Couillard, D. 1990. Encyclopedia of fluid mechanic, surface andgroundwater flow phenomena, volume 10, chapter 5: Predicting rivertemperatures with a hydrological model, pages 171–209. Gulf PublishingCompany.Mote, P. W. and Salathe, E. P. 2010. Future climate in the Pacific Northwest.Climatic Change, 102(1):29–50.Needham, P. and Jones, A. 1959. Flow, temperature, solar radiation, and icein relation to activities of fishes in Sagehen Creek, California. Ecology,40:465–474.Neilson, B. T., Stevens, D. K., Chapra, S. C., and Bandaradoga, C. 2009. Datacollection methodology for dynamic temperature model testing andcorroboration. Hydrological Processes, 23(20):2902–2914.Null, S. E., Viers, J. H., Deas, M. L., Tanaka, S. K., and Mount, J. F. 2013.Stream temperature sensitivity to climate warming in California’s SierraNevada: impacts to coldwater habitat. Climatic Change, 116:149–170.O’Callaghan, J. F. and Mark, D. M. 1984. The extraction of drainagenetworks from digital elevation data. Computer Vision, Graphics and ImageProcessing, 28:323–344.Oke, T. 1987. Boundary layer climates (Second Edition). Halsted Press, London,United Kingdom.160Pacific, V. J., Jencso, K. G., and McGlynn, B. L. 2010. Variable flushingmechanisms and landscape structure control stream DOC export duringsnowmelt in a set of nested catchments. Biogeochemistry, 99:193–211.Pluhowski, E. J. 1970. Urbanization and its effect on the temperature of thestreams on Long Island, New York. Technical report, Geological SurveyProfessional Paper, 627-D.Pomeroy, J. W., Marks, D., Link, T., Ellis, C., Hardy, J., Rowlands, A., andGranger, R. 2009. The impact of coniferous forest temperatures onincoming longwave radiation to melting snow. Hydrological Processes,23:2513–2525.Prata, A. J. 1996. A new long-wave formula for estimating downwardclear-sky radiation at the surface. Quaterly Journal of the RoyalMeteorological Society, 122(533):1127–1151.Priestley, C. H. B. and Taylor, R. J. 1972. On the assessment of surface heatflux and evaporation using large-scale parameters. Monthly WeatherReview, 100(2):81–92.Quilty, E. and Moore, R. D. 2007. Measuring stream temperature. StreamlineWatershed Management Bulletin, 10:25–30.R Development Core Team 2014. R: A Language and Environment forStatistical Computing. R Foundation for Statistical Computing, Vienna,Austria. ISBN 3-900051-07-0.Redding, T. E., Hannam, K. D., Quideau, S. A., and Devito, K. J. 2005.Particle density of aspen, spruce, and pine forest floors in Alberta,Canada. Soil Science Society of America Journal, 69:1503–1506.Regonda, S. K., Rajagopalan, B., Clark, M., and Pitlick, J. 2005. Seasonalcycle shifts in hydroclimatology over the western United States. Journal ofClimate, 18:372–384.Richardson, J. S., Naiman, R. J., Swanson, F. J., and Hibbs, D. E. 2005.Riparian communities associated with Pacific Northwest headwaterstreams: assemblages, processes, and uniqueness. Journal of the AmericanWater Resources Association, 41(4):935–947.Rodenhuis, D. R., Bennett, K. E., Werner, A. T., Murdock, T. Q., andBronaugh, D. 2007. Hydro-climatology and future climate impacts in161British Columbia. Technical report, Pacific Climate Impacts Consortium.University of Victoria, Victoria, British Columbia.Scordo, E. B. and Moore, R. D. 2009. Transient storage processes in a steepheadwater stream. Hydrological Processes, 23(18):2671–2685.Seibert, J. and McDonnell, J. J. 2002. On the dialog between experimentalistand modeler in catchment hydrology: Use of soft data for multicriteriamodel calibration. Water Resources Research, 38(11):1241.Seibert, J., Rodhe, A., and Bishop, K. 2003. Simulating interactions betweensaturated and unsaturated storage in a conceptual runoff model.Hydrological Processes, 17:379–390.Shanley, J. and Peters, N. 1988. Preliminary observations of streamflowgeneration during storms in a forested Piedmont watershed usingtemperature as a tracer. Journal of Contaminant Hydrology, 3(2-4):349–365.Shreve, R. W. 1969. Stream lengths and basin areas in topologically randomchannel networks. Journal of Geology, 77:397–414.Shuter, B. J., Finstad, A. G., Helland, I. P., and Zweimu¨ller, Ho¨lker, F. 2012.The role of winter phenology in shaping the ecology of freshwater fishand their sensitivities to climate change. Aquatic Sciences, 74:637–657.Sidle, R., Tsuboyama, Y., Noguchi, S., Hosoda, I., Fujieda, M., and Shimizu,T. 2000. Stormflow generation in steep forested headwaters: a linkedhydrogeomorphic paradigm. Hydrological Processes, 14(3):369–385.Sinokrot, B. and Stefan, H. 1993. Stream temperature dynamics:measurements and modeling. Water Resources Research, 29(7):2299–2312.Smith, R. 2011. Space-time dynamics of runoff generation in asnowmelt-dominated montane catchment. PhD thesis, University of BritishColumbia.Soetaert, K. 2009. rootSolve: Nonlinear root finding, equilibrium and steady-stateanalysis of ordinary differential equations. R package 1.6.Soetaert, K., Petzoldt, T., and Setzer, R. W. 2010. Solving differentialequations in R: Package deSolve. Journal of Statistical Software, 33(9):1–25.Spittlehouse, D. L. and Black, T. A. 1981. A growing season water balancemodel applied to two Douglas fir stands. Water Resources Research,17(6):1651–1656.162St-Hilaire, A., El-Jabi, N., Caissie, D., and Morin, G. 2003. Sensitivityanalysis of a deterministic water temperature model to forest canopy andsoil temperature in Catamaran Brook (New Brunswick, Canada).Hydrological Processes, 17(10):2033–2047.St-Hilaire, A., Morin, G., El-Jabi, N., and Caissie, D. 2000. Watertemperature modelling in a small forested stream: implications of forestcanopy and soil temperatures. Canadian Journal of Civil Engineering,27(6):1095–1108.Storck, P., Lettenmaier, D. P., and Bolton, S. M. 2002. Measurement of snowinterception and canopy effects on snow accumulation and melt in amountainous maritime climate, Oregon, United States. Water ResourcesResearch, 38(11):1223.Story, A., Moore, R. D., and Macdonald, J. S. 2003. Stream temperatures intwo shaded reaches below cutblocks and logging roads: downstreamcooling linked to subsurface hydrology. Canadian Journal of ForestResearch, 33(8):1383–1396.Subehi, L., Fukushima, T., Onda, Y., Mizugaki, S., Gomi, T., Kosugi, K.,Hiramatsu, S., Kitahara, H., Kuraji, K., and Terajima, T. 2010. Analysis ofstream water temperature changes during rainfall events in forestedwatersheds. Limnology, 11:115–124.Sun, N., Yearsley, J., Voisin, N., and Lettenmaier, D. P. 2014. A spatiallydistributed model for the assessment of land use impacts on streamtemperature in small urban watersheds. Hydrological Processes, DOI:10.1002/hyp.10363.System for Automated Geoscientific Analyses Geographic InformationSystem 2013. SAGA GIS, version 2.1.1. Available from:http://www.saga-gis.org.Tetzlaff, D., Carey, S. K., Laudon, H., and McGuire, K. 2010. Catchmentprocesses and heterogeneity at multiple scalesbenchmarkingobservations, conceptualization and prediction. Hydrological Processes,24(16):2203–2208.Tetzlaff, D., Soulsby, C., Youngson, A., Gibbins, C., Bacon, P., Malcolm, I.,and Langan, S. 2005. Variability in stream discharge and temperature: apreliminary assessment of the implications for juvenile and spawningAtlantic salmon. Hydrology and Earth System Sciences, 9(3):193–208.163Therrien, R., McLaren, R. G., Sudicky, E. A., and Panday, S. M. 2010.HydroGeoSphere: A three-dimensional numerical model describingfully-integrated subsurface and surface flow and solute transport.Groundwater Simulations Group, University of Waterloo, Waterloo, ON.Theurer, F. D., Voos, K. A., and Miller, W. J. 1984. Instream watertemperature model. Technical report, Fort Collins, CO: Fish and WildlifeService Instream Flow Information Paper 16, FWS/OBS-84/15. 335 p.Thompson, J. C. and Moore, R. D. 1996. Relations between topography andwater table depth in a shallow forest soil. Hydrological Processes,10(11):1513–1525.Tung, C. P., Lee, T. Y., Huang, J. C., Perng, P. W., Kao, S. J., and Liao, L. Y.2014. The development of stream temperature model in a mountainousriver of Taiwan. Environmental Monitoring and Assessment, 186:7489–7503.Utting, M. G. 1979. The generation of stormflow on a glaciated hillslope incoastal British Columbia. Master’s thesis, University of British Columbia,Vancouver, Canada.Vogt, T., Schirmer, M., and Cirpka, O. A. 2012. Investigating ripariangroundwater flow close to a losing river using diurnal temperatureoscillations at high vertical resolution. Hydrology and Earth SystemSciences, 16:473–487.Voss, C. I. and Provost, A. M. 2002. SUTRA, a model forsaturated-unsaturated variable-density ground-water flow with solute orenergy transport. Technical report, U.S. Geological Survery WaterResources Investigations Report 02-4231, 270p.Vugts, H. 1974. Calculation of temperature variations of small mountainstreams. Journal of Hydrology, 23:267–278.Ward, A. S., Gooseff, M. N., Voltz, T. J., Fitzgerald, M., Singha, K., andZarnetske, J. P. 2013. How does rapidly changing discharge during stormevents affect transient storage and channel water balance in a headwatermountain stream? Water Resources Research, 49:5473–5486.Webb, B. W. and Crisp, D. T. 2006. Afforestation and stream temperature ina temperate maritime environment. Hydrological Processes, 20(1):51–66.164Webb, B. W., Hannah, D. M., Moore, R. D., Brown, L. E., and Nobilis, F. 2008.Recent advances in stream and river temperature research. HydrologicalProcesses, 22(7):902–918.Webb, B. W. and Zhang, Y. 1997. Spatial and seasonal variability in thecomponents of the river heat budget. Hydrological Processes, 11(1):79–101.Webb, B. W. and Zhang, Y. 1999. Water temperatures and heat budgets inDorset chalk water courses. Hydrological Processes, 13(3):309–321.Wehrly, K. E., Wang, L., and Mitro, M. 2007. Field-based estimates ofthermal tolerance limits for trout: incorporating exposure time andtemperature fluctuation. Transactions of the American Fisheries Society,136(2):365–374.Wehrly, K. E., Wiley, M. J., and Seelbach, P. W. 2003. Classifying regionalvariation in thermal regime based on stream fish community patterns.Transactions of the American Fisheries Society, 132(1):18–38.White, D. S., Elzinga, C. H., and Hendricks, S. P. 1987. Temperature patternswithin the hyporheic zone of a northern Michigan river. Journal of theNorth American Benthological Society, 6:85–91.Wood, J. R. and Hewett, T. A. 1982. Fluid convection and mass transfer inporous sandstones – a theoretical model. Geochimica et Cosmochimica Acta,46(10):1707–1713.Xu, C., Letcher, B. H., and Nislow, K. H. 2010. Context-specific influence ofwater temperature on brook trout growth rates in the field. FreshwaterBiology, 55:2253–2264.Yin, X. and Arp, P. A. 1993. Predicting forest soil temperature from monthlyair temperate and precipitation records. Canadian Journal of ForestResearch, 23:2521–2536.165


Citation Scheme:


Citations by CSL (citeproc-js)

Usage Statistics



Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            async >
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:


Related Items