UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Morphodynamics and sediment transport in a wandering gravel-bed channel : Fraser River, British Columbia Ham, Darren Gary 2005

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

Item Metadata


831-ubc_2005-10497X.pdf [ 46.41MB ]
JSON: 831-1.0092264.json
JSON-LD: 831-1.0092264-ld.json
RDF/XML (Pretty): 831-1.0092264-rdf.xml
RDF/JSON: 831-1.0092264-rdf.json
Turtle: 831-1.0092264-turtle.txt
N-Triples: 831-1.0092264-rdf-ntriples.txt
Original Record: 831-1.0092264-source.json
Full Text

Full Text

Morphodynamics and sediment transport in a wandering gravel-bed channel: Fraser River, British Columbia by Darren Gary Ham  B.Sc. (Hons.), University of Victoria, 1990 M.Sc, The University of British Columbia, 1996  A thesis submitted in partial fulfillment of the requirements for the degree of Doctor of Philosophy in THE FACULTY OF GRADUATE STUDIES (Geography)  The University of British Columbia July 2005 © Darren Ham, 2005  Abstract This study investigates the relation between sediment transport and channel deformation on a 70 km long wandering reach of lower Fraser River, British Columbia. The reach is characterized by an irregularly sinuous single-thread channel split around large bar and island complexes. An extensive network of secondary channels, produced from channel shifting and abandonment, is found within these complexes and along the adjacent floodplain. Most material accumulates within wide 'sedimentation zones' where deposited sand and gravel create an obstruction to downstream flow conveyance, and in areas of flow expansion where shear stress declines. These sites correspond to the location of prominant diagonal riffles. Changes to reach morphology are dominated by the transfer of coarse alluvial sediments. Small gravel sheets are attached to existing lateral and mid-channel bars, forcing compensating erosion across the channel and propagating instability downstream. Over periods of several decades, entire bars and bar / island complexes migrate downstream in association with the riffles. Sediment budgets constructed for several periods between 1952 and 2003 show that the transport rate into the reach is highly variable over time and does not correlate well with the magnitude or frequency of large floods. Influx into the reach has apparently increased over time, from 180,000 m /yr from 1952 to 1984 to 494,000 m /yr from 1999 to 2003. The present high rates 3  3  are associated with channel alignment upstream of Agassiz-Rosedale Bridge that causes extensive erosion at lower Herrling Island. This sequence is apt to be disrupted within the next several years as the island becomes bissected, and transport rates should subsequently decline. The long-term (1952-1999) rate of 212,000 m /yr appears to be the more conservative figure for applied gravel 3  management strategies. This budget also indicates that channel reaches upstream of the bridge have degraded. It is probable that sedimentation in the gravel reach during most of the 20th century was influenced by placer mining in the Fraser Canyon which introduced significant (though unknown) quantities of bed material to the river. More recent factors such as development of the Coquihalla Highway and numerous pipeline and transmission lines may also have led to higher inputs than would be expected in a natural, undeveloped system (Church et al., 2001). However, it is difficult to confirm that these activities continue to impact the gravel budget, or that recent degradation is in fact a response to an exhaustion of upstream sediment supply. Although the gravel reach remains a relatively pristine environment relative to most large navigable rivers in populated regions, the isolation of several large sloughs and backchannels, lateral confinement through bank hardening, and removal of sand and gravel from the channel have - ii -  cumulatively reduced the complexity of channel morphology and reduced the channel gradient. Despite evidence that recent morphologic evolution is trending towards a more stable single-thread meadering habit, considerable downtream instability is expected to continue for several decades as this pattern continues to develop. The relations between sediment transport and channel deformation developed in this thesis should provide better tools to facilitate sound riverine management during this period, but additional survey is required to evaluate such change.  - iii -  Table of Contents Abstract  •  »  Table of Contents List of Tables  iv vii  List of Figures  i x  Acknowledgements  xii  Chapter 1: Introduction 1.1 1.2 1.3 1.4  1  Problem statement Wandering channels Research Hypotheses Morphodynamic changes from survey data 1.4.1 Planform data and analysis 1.4.2 Limits of topographic data and analysis  1 4 7 9 11 12  Chapter 2: Physical setting  18  2.1 Introduction 2.2 Physiography and hydrology 2.3 Anthropogenic influences on morphology 2.3.1 Hood protection structures 2.3.2 Bank hardening 2.3.3 Sediment removals 2.3.4 Woody debris  18 18 25 27 28 31 33  Chapter 3: Data Compilation  34  3.1 Introduction 3.2 Channel and floodplain mapping 3.2.1 Data acquisition from old maps 3.2.2 Data acquisition from aerial photographs 3.3 Bathymetric surveys 3.3.1 Compilation of available bathymetric data 3.4 Construction of 3-D surfaces 3.4.1 Introduction 3.4.2 Reference surface data 3.4.3 Reference surface derivation 3.4.4 Model comparison Chapter 4: Channel changes  34 34 34 38 42 44 49 49 52 55 60 •  4.1 Introduction 4.2 Active channel width  67 67 68  - iv -  4.2.1 Island evolution 4.2.2 Bank erosion and deposition 4.3 Longitudinal profiles 4.4 Cross-section morphology 4.5 Hypsometry 4.6 Summary  73 75 79 85 95 98  Chapter 5: Bed material transport  101  5.1 Introduction 5.2 Sediment budgets from morphologic change 5.2.1 General approach 5.2.2 Sediment budget construction 5.3 Sediment transport 5.3.1 Summary sediment budgets 5.3.2 Time-integrated bias 5.3.3 Distributed sediment budgets 5.3.4 Accumulated sediment budgets 5.4 Sediment budget errors: Negative transport rates 5.4.1 Violation of the zero transport assumption 5.4.2 Sounding errors 5.4.3 Unreported dredging 5.4.4 Sand and gravel fractions 5.5 Sediment budget errors: uncertainty of the estimates 5.5.1 Introduction 5.5.2 Survey error 5.5.3 Floodplain error 5.5.4 Channel Error 5.6 Comparison with other studies  101 102 102 105 112 112 116 118 120 123 123 123 126 128 129 129 129 130 131 133  Chapter 6: Channel deformation and sediment transport  142  6.1 Introduction 6.2 Channel and floodplain formation on wandering channels 6.3 Sediment transfer length 6.4 Sedimentation and morphology of bar complexes 6.5 GIS-based spatial analysis: channel transitions Chapter 7: Conclusions  142 143 146 152 172 182  7.1 Summary 7.2 Evaluation of research hypotheses 7.2 Summary and evaluation of the sediment budget 7.3 Management implications References  182 184 186 189 191  - v-  Appendix A : Channel cross-sections  210  Appendix B : Calculating the gravel budget  222  B.l Computation of unadjusted volumetric changes and bed level changes B.2 Computation of bed material changes and associated bed level changes B. 3 Some example calculations B.3.1 Simple cell B.3.2 Complex cell B.3.3 Summary difference Appendix C : Atlas of channel changes  222 223 225 226 231 232 239  C. 1 Historic channel changes in the 20th century C.2 Hunter Creek to Spring Bar: Figures C-2 to C-5 C.3 Herrling Island to Carey Point: Figures C-6 to C-9 C.4 Foster Bar to Lower Yaalstrick Bar: Figures C-10 to C-13 C.5 Vedder River to Mission Railway Bridge: Figures C-14 to C-17 C.6 Appendix references  - vi -  239 240 243 248 ;..253 256  List of Tables Table 3-1. Historic airphotos selected for morphologic mapping  40  Table 3-2. Comparison of different interpolation schemes using thinned dataset with reference datase 63 Table 4-1. Average annual bank erosion rates. All values expressed in metres per year. Values that exceed the reach mean are shown in bold 78 Table 4-2. Average annual bank deposition rates. All values expressed in metres per year. Values that exceed the reach mean are shown in bold 78 Table 4-3. Summary statistics for total active cross-section area below the reference datum along the extracted cross-sections 94 Table 4-4. Tests for significant differences between mean (t) and variance (f) of active cross-section area below the reference datum. P is the probability that the test statistic is less than the critical value (shown in bold type) 94 Table 5-1. Gross volumetric changes between Mission and Agassiz  112  Table 5-2. Summary of bed material sediment budget calculations. All values are bulk in millions ofm 113 3  Table 5-3. Average annual transport rates (m / yr) for sand, gravel and total bed material load. Rates are equivalent to tthe influx past Agassiz-Rosedale Bridge 114 3  Table 5-4. Average annual transport rates (m / yr) for sand, gravel and total bed material load incorporating adjustments for bedload transport at Mission and possible sounding errors in cells 11-20. Rates are equivalent to the influx past Agassiz-Rosedale Bridge 126 3  Table 5-5. Comparison of WSC transport estimates with sediment budget estimates for gravel. The budget estimates include the adjustment for bedload transport at Mission and the possible sounding errors in cells 11-20 134 Table 5-6. Sediment budget -1952 to 1984  137  Table 5-7. Sediment budget - 1984 to 1999  138  Table 5-8. Sediment budget - 1952 to 1999  139  Table 5-9. Sediment budget -1999 to 2003  140  Table 5-10. Sediment budget - 1984 to 2003  141  Table 6-1. Summary statistics for delineated meanders, 1949 and 1999  148  Table 6-2. Selected low-water photography  153 - vii -  Table 6-3. Transition area and probability matrices for the period 1943 to 1949. 'FROM' states are given in the first column, while 'TO' states are given in the first row 174  - viii -  List of Figures Figure 2-1. Location map and boundary of Fraser River watershed  19  Figure 2-2. Relief map of Fraser River watershed  22  Figure 2-3. Example annual hydrographs, Fraser River at Hope, 1972 and 2001  24  Figure 2-4. Annual maximum and mean daily discharge, Fraser River at Hope  26  Figure 2-5. Location and summary of bank hardening along lower Fraser River  29  Figure 3-1. Spatial extent of historic bathymetric surveys  43  Figure 3T2. The complete testing data for development of the testing grid and a thinned sub-set of data used to compare different interpolation schemes 53 Figure 3-3. Scatterplots showing difference between actual and interpolated values for the reference dataset using A. TIN and B. Topogrid models, and Semivariograms of C. complete reference dataset and D. residuals from TIN model 57 Figure 3-4. Plot of elevation difference between actual and interpolated values for the reference dataset using TIN and Topogrid models 58 Figure 3-5. Reference surfaces using TIN and Topogrid models with no bankline contours  59  Figure 3-6. Reference surfaces based on TIN (B) and Topogrid (C) models. The mapped channel (A) as interpreted from 1991 aerial photographs is shown for comparison 61 Figure 3-7. Comparison of the reference grid with the modeled grid (Topogrid A) using the thinned dataset 65 Figure 3-8. Histograms of 1991 elevation data including distributions for the original 5352 points, the reference grid and the test grid based on the thinned dataset 66 Figure 4-1 .Variation of active channel zone width for the entire gravel reach over time shown with long-term trends in annual maximum daily discharge 69 Figure 4-2. Site map showing morphologic divisions of lower Fraser River gravel reach  71  Figure 4-3. Changes in active channel width over time for individual study reaches  72  Figure 4-4. Areal changes in island area over time for individual sub-reaches  74  Figure 4-5. Average annual bank erosion and deposition rates for different mapping periods. ...76 Figure 4-6. Thalweg profdes for 1952, 1984 and 1999 produced from bathymetric data  ix  80  Figure 4-7. Longitudinal bed profiles of lower Fraser River as measured along the channel thalweg for different dates 82 Figure 4-8. Plot of superimposed longitudinal profiles showing variability in the location of riffles and pools over time 84 Figure 4-9. Location of extracted (even-numbered) cross-sections at 800 metre average spacing 87 Figure 4-10. Selected cross-section plots to show complexity of morphology and dynamics of scour, fill, and thalweg migration in wandering reaches compared to partly confined straight and meandering reaches 89 Figure 4-11. Relation between surface age and surface elevation for old bars, young islands and old (mature) islands between Mission (km 0) and Laidlaw (km 65) 90 Figure 4-12. Active channel cross-section area (a) below reference line along x-sections with islands, (b) below reference line along unvegetated x-sections, and (c) above the reference line 92 Figure 4-13. Histograms of flow depths (relative to 1999 bankfull discharge) within the active chanel zone 97 Figure 4-14. Hypsogram for a confined and wanderinf section of channel from 1999  99  Figure 5-1. Location of 1-km computing cells and distances upstream from Sand Heads  104  Figure 5-2. Surfaces of difference from individual topographic models between Mission and Agassiz showing complex patterns of scour and fill 106 Figure 5-3. (A) Surface and (B) subsurface grain size versus distance upstream from Sand Heads. 108 Figure 5^1. Subsurface sand fraction (%) plotted against distance upstream from Mission  109  Figure 5-5. Net apparent influx of sand and gravel as a function of the intersurvey period  117  Figure 5-6. Cell-by-cell volumetric changes of bed material (top) and gravel (bottom) computed for the period 1952 to 1999 by direct difference of surveys 119 Figure 5-7. Downstream trends in bed material (top) and gravel (bottom) transport rates estimated from the sediment budget calculations 121 Figure 5-8. Downstream trends in bed material and gravel transport rates, 1984 to 2003  124  Figure 5-9. Downstream trends in bed material (top) and gravel (bottom) transport rates estimated from the adjusted sediment budget calculations 127 Figure 6-1. Location and wavelength of river meanders, 1949 -x-  149  Figure 6-2. Location and wavelength of river meanders, 1999  150  Figure 6-3. Examples of bar edge delineation  154  Figure 6-4. Close range (A) and distal (b) views of gravel sheets prograding on existing bar surfaces at Queens Bar complex. (C) Annotated 1943 photo near lower Carey Bar complex showing example delineation of gravel sheers based on tone and morphology 155 Figure 6-5. 1943 channel morphology illustrating hierarchy of depositional units  156  Figure 6-6. Morphologic changes from (A) 1943 to 1949 and (B) 1949 to 1954  159  Figure 6-7. Morphologic changes from (A) 1954 to 1959 and (B) 1959 to 1964  161  Figure 6-8. Morphologic changes from (A) 1964 to 1971 and (B) 1971 to 1979  164  Figure 6-9. Morphologic changes from (A) 1979 to 1986 and (B) 1986 to 1993  167  Figure 6-10. Morphologic changes from (A) 1993 to 1999. (B) Photo mosaic from low water photography taken December, 2003 with superimposed 1999 outline 170 Figure 6-11. Transition state probabilities for water surface and gravel bar polygons throughout the low-water mapping series 173 Figure 6-12. Similarity matrix values of gravel bar transitions for low-water mapping periods 176 Figure 6-13. Locational probability of channel bars and water surface  178  Figure 6-14. Decay curves showing loss of bar and island surfaces since 1943  180  - xi -  Acknowledgments First and foremost, I wish to express my sincere thanks to my supervisor, Dr. Mike Church, for his invaluable contributions and support during the completion of this thesis. Because of his always thoughtful commentary, thorough review, and stimulating questioning, I have become a better writer and scientist. I should also like to express my gratitude to the other members of my committee: Dr. Ted Hickin, Dr. Brian Klinkenberg and Dr. David McLean for their input and review. Thanks are extended to Mr. Terry Keenhan, formerly of the Ministry of Environment, Lands and Parks, for initiating interest on the gravel management concern in lower Fraser River, and to Mr. Ron Henry, senior hydraulic engineer responsible for lower Fraser River (also with MoELP, later Ministry of Water, Land and Air Protection) for his continued interest in the project. I am appreciative for the efforts of Karen Kranabetter and Jason Barraclough with regard to the monumental task of digitizing so many airphotos, and to Hamish Weatherly for his assistance with the early formulation of the sediment budget. I would be remiss if I did not also thank Scott Paper for access to their archival materials, the Fraser Basin Council for their interest and support of this project and to Kevin and Rosemary in the GIC for always being helpful in accessing materials. Finally, I am grateful to my friends and family, especially my parents, who so strongly supported my decision to undertake a PhD. Financial assistance for this project was provided by the British Columbia Emergency Flood Protection Program, The City of Chilliwack and a National Sciences and Engineering Research Council strategic grant awarded to Dr. Mike Church. The project could not have been completed without this assistance.  - xii -  Chapter 1: Introduction 1.1 Problem statement The floodplains of large alluvial rivers have long been a focus of human settlement, providing ideal locations for development and agriculture, important sources of food and building materials, and vital conduits for transportation and commerce. For as long as such areas have been populated, flooding and channel erosion have prompted river engineering solutions to minimize these impacts, and costly maintenance programs to ensure effectiveness (Schumm and Winkley, 1994). Consequently, many rivers have been altered to some extent to accommodate the needs of people, usually to the detriment of the natural environment. As population growth fuels the need, and technological advances continue to enhance our ability to modify riverine environments, it becomes increasingly important to recognize the processes that govern river behaviour and to understand the resultant morphological response over appropriate spatial and temporal scales. Failure to do so will undermine efforts to effectively manage alluvial channels while preserving ecosystem health and diversity. Attempts to relate form and process relations on alluvial channels have been a dominant theme in fluvial science for the past half-century (Leliavsky, 1955, Leopold, Wolman and Miller, 1964). This avenue of study is complicated by the large number of variables involved, the broad range of possible channel adjustments and the scale of investigation. Over periods of years to decades (the timescale of central interest in this thesis), alluvial channels can exhibit a wide variety of distinct morphologic forms in response to the complex interaction of various governing conditions. Along the length of a channel network, the presence of a given morphology within a homogenous reach is dominantly influenced by the local flow and sediment regimes (cf. Mollard, 1  1973; Kellerhals et al., 1976; Richards, 1982; Church, 1992) and by valley gradient (Carson, 1984; van den Berg, 1995). Further, land use patterns (Dykaar and Wigington, 2000; Kondolf et al., 2002), riparian vegetation (Hickin, 1984; Millar, 2000), bank strength (Schumm, 1963; Millar and Quick, 1993) and anthropogenic modifications (Kellerhals and Church, 1989; Collins and Dunne, 1990) have been identified as secondary factors that locally influence adjustments to channel shape, position and pattern. Together, these variables directly influence the entrainment, transport and deposition of clastic sediments, and hence the morphology of alluvial channels.  1. Reaches are variable lengths of channel within which the dominant factors that influence morphology do not change appreciably such that channel pattern remains reasonably uniform.  -1 -  The direct association between sediment transport and morphology has not been widely pursued until relatively recently. Although the basic relation is well established (following the sediment continuity equation), a complete solution requires a reliable estimate of sediment transport (Hickin, 1983). This problem has itself been a subject of extensive study over the past century, yet has not been satisfactorily resolved (Carson and Griffiths, 1987; Reid and Frostick, 1994). Fluvial scientists have long sought to develop functional relations that correlate transport rates to some property of flow. Attempts to relate sediment transport and morphologic development to hydraulic variables in rivers have been largely unsuccessful, however, due to the non-uniformity and unsteadiness of flow (Gomez, 1991; Ferguson and Ashworth, 1992). Ashmore and Church (1998) point out that functional relations are statements of equilibrium conditions, which are rare in nature. Consequently, predictions may be grossly over-estimated in supplylimited channels, especially where bed structures increase the force required to mobilize clasts (Brayshaw, 1984; Church et al., 1998). Similarly, under-predictions are common where input parameters such as peak instantaneous flow (hence shear stress) and channel dimensions are averaged over larger temporal and spatial scales. In particular, most bedload equations fail to adequately perform on systems for which they were not developed, having been calibrated for a particular grain size distribution and set of flow parameters that may simply not translate to other systems. Furthermore, any approach that stresses equilibrium conditions is not suitable for understanding river behaviour over longer timescales, although this is the period over which management concerns become pronounced. Perhaps more importantly, such approaches to the problem do not emphasize the important links between sediment transport and morphology which are often central to investigations of channel stability (Ashmore and Church, 1998; Hooke, 2003) and which provide the physical basis for aquatic habitat (Lisle, 1989; Kellerhals and Miles, 1996). Nevertheless, their use remains critical in numerical modeling studies where a sediment transport function is specified, though uncertainty in the reliability of such functions may require an independent means of calibration if the models are to be used in predicting bed level changes (cf. Li and Millar, 2004).  2. The morphology of a channel is determined by the bed material load, or sediment that comprises the bed and lower channel banks. Bed material is mainly transported as bedload, where particles travel by rolling or sliding along the bed, but may also include sediments moving in suspension. The bed material load is distinguished from the wash material load, which consists of sediments that mainly travel through a reach in suspension, are deposited only on overbank surfaces during floods, and have minimal influence on channel morphology.  -2-  Numerous attempts have been made to measure bedload transport directly using a variety of sampling devices and traps. This approach requires an extensive field effort because the spatial and temporal variability associated with its movement makes observation and collection of transported sediment exceedingly difficult (see reviews in Hubbell, 1987; Gomez, 1991). It is necessary to integrate a large number of samples over a wide range of flow conditions to adequately characterize a river's sediment transport regime. This is especially true of larger rivers, where adjustments may occur over decades for example, and similarly require long-term monitoring programs (Macklin et al., 1998). Sampling difficulties are amplified by the fact that most sediment transport in gravel bed rivers occurs during high flow events when direct observation and mensuration is difficult, even dangerous. Consequently, most direct sampling programs have been limited to fairly small channels (though measurements on Fraser River provide an unusual exception). Concern over the accuracy of the sampling devices themselves raises an additional problem that requires sampler calibration and adjustment for trapping inefficiencies (Hubbell, 1987; McLean et al., 1999; Sterling and Church, 2002). Although many researchers continue to pursue theoretical, empirical and direct approaches to estimating bedload transport, an alternative strategy has emerged in recent decades which emphasizes the obvious association between sediment transport and morphologic change. Since the morphology displayed by a channel reach is a direct consequence of sediment transport (McLean and Church, 1999), observations of channel deformation can reveal details of sediment transfer in rivers (cf. Ashmore and Church, 1998). The association between morphologic measurement and sediment transport was first recognized by Popov (1962) but the study did not clearly demonstrate how they were related (Neill, 1971). Neill (1971) presented a more formal approach, demonstrating that an estimate of sediment transport along a meandering channel could be made by measuring bank erosion volumes over time. Church et al. (1987) provided a more general statement of Neill's approach and adapted the technique to estimate sediment transfer on rivers with more complex morphologies. More recent studies have modified the technique to describe temporal changes in local bed material transport rates for a variety of river types and sizes (McLean, 1990; Ferguson and Ashworth, 1992; Goff and Ashmore, 1994; Lane et al., 1995; Ham and Church, 2000; Lane et al., 2003). Refinements in measuring techniques have been an integral part of this progress. Although there are circumstances in which the morphologic approach fails to provide a reliable estimate of sediment transfer, the technique reveals details of spatial and temporal channel -3-  change that other approaches do not. Since the underlying principle relating sediment transport to morphologic change involves intersurvey comparison, useful information on rates and locations of morphologic change can be derived as a consequence. McLean (1990), Ham and Church (2000) and Gaeuman et al. (2003) provide examples for bank erosion, channel widening and changes in channel morphology utilizing planimetry-based studies. Other researchers have demonstrated the magnitude and dynamics of channel scour and fill along surveyed cross-sections (Ferguson and Ashworth, 1992; Goff and Ashmore, 1994; Hickin, 1995; Martin and Church, 1995; Paige and Hickin, 2000). Increasingly, three-dimensional surfaces of difference have been generated to illustrate complex patterns of sediment erosion and deposition distributed along entire channel reaches (McLean, 1990; Eaton and Lapointe, 2001; Lindsay and Ashmore, 2001; Lane et al., 2003). These efforts have begun to reveal details of the association between channel form and process rates on braided rivers (Ferguson and Ashworth, 1992; Lane et al., 1995; Brasington et al., 2000; Gaeuman et al., 2003). While such developments remain incomplete, it appears clear that observations of morphodynamic changes may usefully provide a means to better understand form and process relations in alluvial rivers. As such, the approach provides a valuable tool for monitoring and evaluating the sensitivity of riverine systems to environmental changes, hence to predict channel response to disturbance. The purpose of this research is to elucidate the relation between sediment transfer and channel morphodyamics on the 75-km long wandering gravel-bedded reach of lower Fraser River, British Columbia. It is intended that this research will improve our scientific understanding of the style and rate of instability along large wandering channels and provide better tools to facilitate management decisions that impact channel development. At present, our knowledge of the spatial and temporal behaviour of these systems at the channel scale is incomplete. Although it is known that channel morphodynamics is connected to the downstream staging of bed material, the variability of sediment transport has hindered scientific progress towards defining this association.  1.2 Wandering channels Lane (1957) noted that although there is nearly an infinite variety of possible stream forms, similar combinations of the major controlling variables produced similar principal patterns, which he classified as meandering or braided types. Leopold and Wolman (1957) added that channel planform formed a continuum, and introduced straight, meandering and braided channels as quasiequilibrium forms. Both works are notable for introducing standardized definitions such that -4-  different rivers could be classified and compared between regions, and for an attempt to relate channel form to principal controlling factors. Both are largely remembered for the development of discriminating functions distinguishing meandering from braided channels on the basis of slope and discharge (Ferguson, 1987). Although the classification introduced by these authors, particularly that of Leopold and Wolman, has since become firmly entrenched in the literature, it has also been criticized for being too restrictive. Fluvial scientists now recognize a much broader range of planform types, including variations of the original tripartite classification, several transitional forms and a variety of multichannel forms where the flow is separated around semi-permanent vegetated islands. Mollard (1973) presented 17 example channel morphologies illustrating the intergrading between classical end-member types and attempted to link channel and floodplain form with the variables that influence them. Mollard (1973) also gave examples of wandering channels, but reserved the term 'anastomosing' to describe high-energy, low-sinuosity streams split around wooded islands. A similar classification was given by Schumm (1985), although the term 'anabranching' was adopted to characterize channels split around large islands, while anastomosing channels were not considered to be a distinct type because individual channel branches could exhibit straight, meandering or braided patterns. More recently, Church (1992) presented a conceptual classification scheme for large rivers that combined concepts from Mollard (1973) and Schumm (1985) but which emphasized inchannel sediment storage. The style of lateral channel stability (after Kellerhals and Church, 1989) was presented as an additional parameter for classification. Church (1992) distinguished wandering from anastomosing channels on the basis of sediment calibre, with the latter type conveying a much finer load. The importance of sediment load (suspended, mixed, bedload) in influencing channel type was also explicit in Schumm's (1985) work. More recently, Nanson and Knighton (1996) have nominated multiple-channel forms as types of anabranching rivers in an attempt to clarify confusing nomenclature. This categorization includes laterally-active, graveldominated rivers, more commonly referred to as wandering channels. They are largely discriminated from other anabranching channels by greater specific stream power, larger sediment size and less cohesive banks. The term 'wandering' was originally used to describe the process of thalweg migration between channel banks (Leopold and Wolman, 1957) and amongst semi-permanent bedforms (Kopaliani and Romashin, 1970). Neill (1973) first formally nominated the term to characterize a -5-  distinct channel planforrn type in a discussion of the characteristics and behaviour of the gravelbedded Athabasca River, Alberta. Wandering gravel-bed channels are now commonly recognized throughout mountainous and piedmont regions of Western Canada (Desloges and Church, 1989) and have further been identified in China, Europe, New Zealand, Russia and western United States. These rivers typically form part of a downstream continuum of channel planforrn types between braided and meandering reaches and develop in response to environmental controls that vary downstream (Church, 1983a; Brierley, 1989). Wandering channels possess similarities to both braided and anastomosed systems (Knighton and Nanson, 1993) but exhibit several distinguishing characteristics. They commonly flow in irregularly sinuous single-thread channels but are frequently split around large wooded islands, even at peak flows. Reach morphology is dominated by the erosion and mobilization of stored floodplain sediments, though tributary streams may contribute significant loads. Coarse sediments accumulate in locally unstable 'sedimentation zones' (cf. Church, 1983a) which may exhibit low-order braiding, and which are separated by narrower, stable 'transport zones.' Aggradation along the channel is spatially and temporally discontinuous - both the location and number of sedimentation zones are apt to evolve over time. Sedimentation zones locally increase water-surface elevations, causing overbank flooding, and direct flows around them resulting in bank instability and channel shifting. Over time, these processes typically produce an extensive network of seasonal, perennial and abandoned sidechannels along the adjacent floodplain. These pathways delineate positions where the active channel flowed in the past, and are subject to avulsions and cutoffs during high flows (Desloges and Church, 1987). Bedrock or erosion-resistant bank materials commonly abut channel margins along transport zones, creating deep scour holes. Floodplain modification and channel instability along wandering rivers are directly related to the downstream transfer and storage of coarse alluvial sediments, although the association is not well understood. Along lower Fraser River, British Columbia, these natural processes create significant management concerns for infrastructure developments, silvicultural investments, aquatic habitat diversity, navigation and can even limit recreational opportunities. Insofar as applied management strategies (i.e. dyking, placement of riprap and gravel removals) may conflict with the maintenance of aquatic habitat for the spawning, rearing and migration of numerous fish species, including great runs of Pacific salmon, mitigating potential environmental, economic and engineering conflicts requires an increased understanding of the links between channel and floodplain morphodynamics and the contemporary sediment cascade. This knowledge has -6-  practical implications for successful land-use and riverine management along other wandering channels where similar conflicts are shared or expected. 1.3 Research Hypotheses Although bed material generally represents only a small fraction of the total load carried by most alluvial rivers (Church 1985; Meade et al., 1990), the staging of coarse sediment along a river channel directly alters channel morphology. As a consequence of spatial and temporal displacements in channel form, both natural and anthropogenic systems may be directly impacted. Excessive scour or fill can have a detrimental influence on the amount and quality of spawning gravels (Reiser, 1998) and can reduce the diversity of available habitat types by changing the distribution of flow depth and velocity (sensu Hogan and Church, 1989). Erosional processes may also undermine riprap revetments, flood protection structures, bridge piers or pipelines and lead to the loss of developed lands, while depositional processes can fdl reservoirs, increase the flood risk, and create navigational problems. Providing reliable estimates of rates of bed material erosion, transport and deposition along river channels, therefore, has a number of practical environmental and engineering applications. Along unconfined alluvial channels, rivers typically modify their pattern through the downstream progression of bank erosion and sediment deposition in channel bars (Church, 1983a). Neill (1971, 1987) recognized this association on rivers displaying a regular downstream progression of meanders and used it to estimate the rate of bed material transport on a large braided river. Along meandering rivers, quantitative theories have been sufficiently developed to predict facies patterns and planform evolution (Allen, 1965, Ferguson and Ashworth, 1992). Bend migration has been found to relate to bank erosivity, channel width and bend curvature (Hooke, 1980; Hickin and Nanson, 1984) and has been partly described by kinematic wave models (Beck, 1984; Mosselman et al., 1995; Gilvear et al., 2000). The relations between braided channel morphology and causative variables has received comparatively less attention (Ferguson and Ashworth, 1992), although the conditions necessary for braiding and mechanisms of braid development appear fairly well known (Miall, 1977; Ashmore, 1991; Knighton, 1998). Recent studies (Ferguson and Ashworth, 1992; Goff and Ashmore, 1994; Lane et al., 1995; Murray and Paola, 1994) have begun to reveal details of the spatial patterns of sediment transport and morphologic development in small braided streams and attempts at predicting channel evolution in braided rivers have been made (Lane et al., 1995; Paola, 2001; EGIS, 2002).  In comparison to the classic forms of meandering and braided channels, the style of instability along wandering rivers remains poorly described. Quantitative models of deformation along wandering channels do not currently exist. These rivers are not easily amenable to scale model studies (e.g. compared to flume studies of straight or braided channels) given the difficulties of adequately representing a vegetated channel-floodplain system (although recent work by Gran and Paola, 2001 demonstrates progress). As a result, the physical processes which determine channel form and development have not been well addressed. The most suitable approach to observe and measure patterns and rates of sediment transfer and morphologic development along large rivers involves mensuration of channel and floodplain processes at appropriate temporal scales. This appears to be years to decades along lower Fraser River, which corresponds to the persistence of major sedimentation zones (McLean, 1990; McLean and Church, 1999). This extendedtimescaleoccurs because the volume of bed material that is mobilized along wandering channels is modest relative to that stored within the channel and floodplain (Desloges and Church, 1989; Ham and Church, 2000). As a consequence, many years of persistent erosion or deposition are required to cause significant changes to the morphology. Similarly, large floods may complete or accelerate channel changes currently in progress (McLean and Church, 1999) but are not necessarily associated with major, distinct morphologic developments. Insofar as the patterns of sediment transfer and morphologic evolution occur relatively slowly, archival data records provide the appropriate temporal resolution for investigating this behaviour. The matter is discussed in further detail at the end of the chapter. The main objective of this study is to elucidate the association between the staging of sediment and the style of morphologic change along a large, wandering gravel-bedded channel. The lower Fraser River is an ideal site for this investigation because historic documentation of the river is extensive. Planimetric maps and aerial photographs date back more than a century, while topographic surveys of the channel bed and banks are available since 1952. As well, the river has been continuously gauged since 1912, while prior information on sediment flux is available from a period of direct measurement (Tywoniuk, 1973; McLean et al., 1999) and from a comparison of computational formulae with morphologic estimates (McLean, 1990). The specific research hypotheses to be examined are: 1. that the spatial and temporal variability of sediment transport on lower Fraser River can be directly related to observable sequences of channel morphodynamics;  2. that historic patterns of morphodynamic change can be used to predict locations and rates of future instability; and 3. that past gravel removals, construction of flood dykes and riprap bank protection, and a reduction of naturally occurring woody debris deposits have significantly altered the sediment transport and morphology regimes of the river over the past century, and may have caused degradation of the channel  The remainder of this chapter reviews the potential for inferring channel morphodynamic changes from survey data, including historic maps, aerial photographs and hydrographic surveys. The following chapter provides an overview of the physiographic setting of the study area and includes a detailed description of major anthropogenic changes that have occurred since early European settlement of lower Fraser Valley. Chapter 3 gives details on the data available for this study and explains how these data are analyzed for subsequent use. The interpretation of this information is given in Chapter 4, where quantitative analysis of channel change over time is provided, and in Chapter 5, where the association between the processes which transport bed material and the deformation of the channel is used to provide a detailed bed material sediment budget. Chapter 6 investigates how sediment is staged through the reach and relates this movement to spatial and temporal patterns of morphologic change. The final chapter discusses probable future evolution of the river based on observed patterns and trends and presents implications for successful management of the system in the future. 1.4 M o r p h o d y n a m i c changes f r o m survey data Rivers of different size and morphology exhibit a wide variety of depositional forms in response to the calibre of sediment being transported. Accordingly, any successful attempt to model these forms and record their displacement requires that the resolution of channel surveys corresponds to the deformation style of a given river. Application of the morphologic approach, therefore, requires some a priori knowledge on the details of sediment transport and storage. By definition, efforts to relate sediment transport to channel deformation are limited to the movement of bed material only, that is sediments which comprise the bed and lower banks of alluvial channels (cf. McLean et al., 1999). Along gravel-bedded channels, this material typically consists of coarse sand and larger particles which move relatively short distances during flood events. Finer sediments may be introduced to the system from erosion of channel banks and islands, where such -9-  materials may constitute a significant proportion of overbank flood deposits. Once entrained, however, these sediments are mainly transported long distances in suspension, though some fraction is deposited interstitially within the bed. Suspended sediments do not appreciably influence channel morphology except where they enhance upper bank strength (Martin and Church, 1995). Since sediment erosion, transfer, and deposition result in deformation of the channel, a measurement technique must provide a means to capture this deformation in sufficient detail to prevent, or at least minimize, unmeasured changes that otherwise occur between successive surveys. While methods to continuously monitor the movement of bedload at fixed locations have been developed (cf. Rennie et al., 2002; Rennie and Millar, 2004) no method currently exists for measuring real-time deformation at the reach scale. Bed material exchanges alongriversprimarily include scour and fill along the channel bed, erosion and deposition along islands and channel banks, and sediment accumulations and movements amongst bed and bar forms. As movement of this material is known to exhibit high spatial and temporal variability (Meade, 1985; Hubbell, 1987; Gomez, 1991; Nicholas et al., 1995), morphologic estimates provide a reasonable means to average this variability over the reach scale provided sufficient survey detail is available (Knighton, 1998; Fuller et al., 2003). Along gravel-bed rivers, most transported bed material is stored in large-scale bars (Neill, 1987) which scale to channel width or greater and have thicknesses comparable to mean bankfull depth (Church and Jones, 1982). The frequency with which bar migration occurs varies enormously by channel type and scale. Along small meandering and braided channels, changes in morphology have been observed daily (Meade, 1985; Goff and Ashmore, 1994; Lane et al., 1995). Conversely, detectable change along large wandering channels may take several years if there are no large floods, although some minimum volume of transport usually occurs seasonally. Ideally, survey data would be available at a resolution appropriate to record this displacement. However, most studies must make use of commonly available survey data, including available maps, aerial photographs, satellite imagery and cross-section surveys, though application of emerging remote3  3. Satellite imagery can provide information similar to that derived from conventional aerial photography, although multiband sensors allow individual spectral signals to be combined for increased interpretive capabilities. Few riverine studies make use of this technology, however, because commercially available images date back only to the 1970s, and coarse spatial resolution limits application to very largeriversonly. Given recent improvements in spatial resolution, their use may become increasingly common for geomorphologic studies. However, since the spatial and temporal limits of airphotos can be extended to include satellite images, the latter are not specifically discussed herein.  - 10-  sensing technologies is becoming more common. The technical and practical limits of these different techniques for measuring changes in sediment storage, hence estimating bed material transport rates, are discussed below. 1.4.1 P l a n f o r m data a n d analysis  The use of historic sequences of maps and aerial photography for quantifying channel changes is well-established in the geomorphic literature. Maps have particular utility for channel change mapping because this type of record extends further in time than any other type of riverine record (see Chapter 3). Aerial photographs have a much more limited temporal coverage, but their use eliminates cartographic biases often associated with maps. As well, photographs frequently allow measurements to be made over much greater spatial scales than is generally practical in the field. Along rivers, their use is largely restricted, however, to phenomena that can be observed in planform above the subaqueous zone. Typically, the main application has been to measure changes in channel width and lateral bank erosion over some length of channel (Hooke, 1980; Hickin and Nanson, 1984). Neill (1971; 1987) first formally demonstrated that an estimate of sediment transfer could be made on meandering rivers by measuring average bank recession rates from maps and airphotos. Assuming that the height of the bank could be measured, the only limiting assumption is that an estimate of average travel length can be made (given by Neill as half the meander wavelength on a meandering channel). Neill (1987) cautioned, however, that the method would underestimate actual transport volumes when the assumed transfer length was too short; similarly, the technique will overestimate actual short-term transport volumes if material moves only locally (Ashmore and Church, 1998). McLean (1990) tested the approach on an eroding island on Fraser River that exhibited similar progressive meandering over a 20-year period and compared eroded volumes over time to estimates of bed material transport at an adjacent gauging station where bedload samples were available. McLean found the approach underestimated actual transport values because migrating gravel sheets passed undetected through the study area, suggesting the technique is further limited along rivers that exhibit multiple modes of sediment transport. This bias of undetected sediment transfer obviously is magnified as the intersurvey period increases. Carson and Griffiths (1989) applied the travel distance approach to a 5-km section of the braided Waimakariri River, New Zealand. Travel distance was assumed equivalent to the average distance between zones of scour and fill along the reach. Erosional and depositional areas were measured from overlays of photos taken at different dates over 16 months, and converted to -11 -  volumes using representative scour and fill depths from repeat cross-sections. The authors noted that, by acquiring photographs following each significant flow event, undetected throughput should be minimized (assuming throughput is actually proportional to discharge), but added that the morphologic approach can provide only a lower bound estimate. Ham and Church (2000) demonstrated the magnitude of unmeasured flux due to compensating scour and fill along a 49-km length of the wandering Chilliwack River, British Columbia. Estimates of bed material transport measured at intervals of 7 to 14 years were compared to shorter-term estimates along a single downstream reach. The study revealed that long-term average transport rates clearly underestimate the magnitude of short-term transport because sites of distinct erosion and deposition are masked by compensating scour and fill. This circumstance is a consequence of the timing of individual floods. Photographs would ideally be available after each major flood to reduce unmeasured scour and fill over longer periods. This would also ensure that channel morphology was sufficiently modified to be detected, hence measured. Assuming photos are available at suitable dates, reliable estimates of sediment transport are mainly affected by errors in storage change estimates. This largely results from uncertainty in measuring the thickness of deposited bed material (Ham and Church, 2000; Gaeuman et al., 2003) since planimetric (mapping) errors are comparatively small. The uncertainty can be reduced where repeat cross-sections are available, but otherwise requires sufficient field measurements to establish reach-averaged or denser estimates of bar and bank heights. However, Ham (1996) recognized that planimetric errors could dominate between successive dates if water levels (discharges) were significandy different. A condition of low water on the first map, and high water on the second results in the interpretation of erroneously small deposition zones, and correspondingly large erosional zones, leading to an apparent [false] impression of net scour, even if no sediment was actually moved. Church et al. (1987) demonstrated a technique for resolving this issue by relating areal exposure to stage, but recognized the approach was limited by several general assumptions and a lack of sufficient data. Their technique was modified by Ham and Church (2000) for wide, shallow reaches by replacing flow depth with flow 'area'. While this adjustment resolves most of the limiting assumptions, the procedure has not been tested on different morphologies. 1.4.2 Limits of topographic data and analysis Most morphologically-based studies have been based on between-survey comparisons of topographic data collected through regularly-spaced cross-sections, ground survey, or remote - 12-  sensing, though tracers and scour chains have been used to measure deformation on very small channels (Hassan and Church, 1992; Haschenberger and Church, 1998). Through the comparison of interpolated surfaces, a distributed, three-dimensional view of erosion and deposition is provided, which can reveal details on the links between flow hydraulics, channel topography and sediment transport (Ferguson and Ashworth, 1992; Goff and Ashmore, 1994, Lane et al., 1995). The obvious advantage of this pursuit is that it captures details of bed deformation that are not detectable by planimetric means. Nevertheless, there are a number of concerns relating to the spatial and temporal resolution of the surveys that can limit applicability, though recent advances in data acquisition techniques have begun to address these. It is also important to recognize that historic data sources may be sparse, even non-existent, so longer-term trends in channel deformation can not always be determined. Neill (1987) suggested that the approach developed for channels with systematic meander shifting could be adapted for channels characterized by large-scale bar migration. Given that this process results in thalweg shifting, repeat cross-section surveys could be used to measure scour and fill volumes across the channel. Though not explicitly stated, the approach also assumes that crosssections are taken at the same fixed locations, limiting use to smaller channels. Large rivers are conventionally surveyed by depth sounding, and sampling positions are seldom replicated over time. Church et al. (1987) were able to roughly compare sounding surveys over two periods by subtracting averaged depth values to calculate net storage change. McLean (1990) followed this procedure along the gravel-bed reach of lower Fraser River using more sophisticated computerized modeling to interpolate and compare surfaces generated from two bathymetry surveys and was able to both illustrate and calculate spatially distributed volumes of scour and fill. More recent studies have taken advantage of improved data collection and analysis techniques to more fully explore the spatial and temporal limits of data acquisition that 'constrain practical application of the morphologic approach and limit the precision of results' (cf. Ashmore and Church, 1998). Ferguson and Ashworth (1992) and Goff and Ashmore (1994) used closely spaced cross-sections to interpolate distinct sites of erosion and deposition, but it was found that intervening scour and fill or channel migration can mask the distinctiveness of erosion and deposition sites on small braided channels even within a few days. Lane (1997) was able to quantify the effects of changing temporal resolution, demonstrating an apparent time-dependent bias where the magnitude of net volumetric change decreased as the volume of water through a reach increased. While unmeasured scour and fill represent a negative bias in volume change - 13 -  estimates (cf. Ashmore and Church, 1998), event-based surveys may not always be practical or necessary. On large rivers dominated by snowmelt floods, annual surveys may be adequate (Ashmore and Church, 1998). Correspondingly, ideal boundary conditions (i.e. zero or minimal downstream transport) can eliminate, or at least greatly reduce, temporal bias. Several studies have identified the magnitude of loss of information that results from decreased survey spacing (McLean; 1990, Lane et al., 1994; Martin and Church, 1995) although Lane et al. (1995) noted that larger cross-section spacing may be appropriate in larger rivers. It is therefore reasonable to assume that increased survey density should obviously produce the opposite effect. Lane (1998) suggested that this issue had not been fully explored because of the legacy of at-a-station hydraulic geometry, but there were more practical constraints: until recently, there was simply no method available by which spatially distributed data could be economically or conveniently collected. Existing approaches, including manual photogrammetry and handlogged field survey were simply too time-consuming to collect distributed topographic data over large spatial scales or short temporal intervals. Lane et al. (1994) attempted to resolve this problem by collecting topographic data using a combination of close-range photogrammetry for exposed regions with tacheometric survey of the subaqueous zone, recognizing that widely spaced crosssections are insufficient to reveal details of within-reach erosion and deposition. This approach overcame the traditional trade-off between collecting data at higher spatial densities over large areas and the frequency of repeat survey (Lane et al., 1995). Several researchers have since followed this work using electronic total stations (Milne and Sear, 1997; Eaton and Lapointe, 2001; Fuller et al., 2003) and kinematic GPS surveys (Higgitt and Warburton, 1999; Brasington et al., 2000). Not only do these techniques help overcome spatial and temporal limitations of data collection, but because they are explicitly field-based, collection procedures can be modified on-the-fly to account for complex terrain or features of particular interest such as banklines. There is of course a lower limit for temporal spacing set by the precision with which any detectable change can be made for a particular morphology (Ashmore and Church, 1998). Surveying may be limited by fast or deep water, limiting use to small, wadeable systems. Chandler (1999) demonstrated the use of automated digital elevation extraction for quantifying landform topography using vertical and oblique photographs, an approach followed by Lindsay and Ashmore (2002) to measure scour and fill on a model braided river over periods as short as ten minutes. Accurate topographic reproductions of the landscape may be limited by vegetative cover (though analytical techniques can largely resolve this issue) so the use of photos  - 14-  remains limited to exposed surfaces. Applications of through-water photogrammetry have been explored on shallow braided channels (Westaway et al., 2000) including shallow, modestly turbid waters (Westaway et al., 2003) but deep water, surface turbulence and high turbidity levels limit practical adoption along rivers where such conditions may be common. Water also represents a measurement barrier to conventional L i D A R (light detection and ranging), an emerging technology that uses pulsed laser to measure the distance from an aircraft to the ground. L i D A R has demonstrated potential to rapidly survey large areas at high spatial resolution (Charlton et al., 2003; French, 2003). Further, the technique is unaffected by shadows, interpolation schemes can estimate ground elevations beneath the vegetation canopy, and data collection appears more costeffective than other approaches, an advantage which increases with areal coverage. A n alternative type of blue-green L i D A R known as SHOALS (scanning hydrographic operational airborne lidar survey) uses multiple survey frequencies to survey up to 60 metres depth in clear water, but is currently restricted to a maximum of 2-3X the Secchi depth (http://shoals.sam.usace.army.mil; see also Guenther et al., 2000). While such advances show promise, large turbid rivers (including Fraser River) still require survey of subaqueous regions using more conventional approaches. As such, there will be an obvious mis-match in survey frequency above and below the water line that will affect the overall quality of individual digital elevation models (DEMs). Current research has also begun to explore the accuracy by which irregularly spaced topographic data can be converted to a digital representation of a geomorphic surface, either in the form of a triangulated irregular network (TIN) or regular array of interpolated grid cells. There are several issues to consider, including the accuracy of individual survey points (Lane et al., 1994; Chandler, 1999) and the density and distribution of points (Brasington et al., 2000). The size of interpolated grids is also important (Milne and Sear, 1997; Fuller et al., 2003), especially where TIN facets are converted to a regular grid array. Studies have also examined the validation of individual DEMs, a process which involves comparing modeled elevations to independently surveyed checkpoints (Westaway et al., 2000; Brasington and Smart, 2003) and error propagation between DEMs of difference (Lane et al., 2003). Such comparisons are not only important for producing reliable estimates of volumetric changes, but they identify the lower limit from which morphologic change can be detected (cf. Lane et al., 2003). Somewhat surprisingly, the choice of an appropriate modeling strategy has been poorly investigated. Although the effects of interpolation accuracy are related to data density (and can largely be ignored for spatially dense datasets), there can be significant differences in apparent storage volumes when the different - 15-  models are applied to well-spaced points, such as along individual cross-sections. This concern is analyzed in greater detail within Chapter 3. Extending measures of net volumetric changes to estimates of sediment transport is also complicated by assumptions about the transport lengths of particles during transport or by the need to define an upstream or downstream transport rate (Lane, 1997; Ashmore and Church, 1998). Neill (1987) recognized that reliable application of the approach required knowledge of an appropriate sediment step length, with surveys taken at closely spaced spatial and temporal intervals. Church et al. (1987) presented a more generalized method by which information on net storage changes could be combined with a local estimate of the transport rate to extend transport estimates upstream and downstream, termed a sediment budget approach. For some length of channel, the transport rate out of a reach was determined as the transport volume entering, minus the change in volume stored. This procedure eliminates the requirement that the transport steplength be known, provided the surveys are sufficiently dense to define loci of erosion and sedimentation (Lane et al., 1995) and the length of the study reach exceeds typical transport distances (Eaton and Lapointe, 2001). Although Church et al. (1987) lacked sufficient data to construct a sediment budget, the approach has since been widely pursued. Ferguson and Ashworth (1992) and Goff and Ashmore (1994) prepared within-reach sediment budgets using a bounding transport estimate. Their efforts were confounded, however, by the requirement to identify a reference transport, which both studies identified as a weakness in the approach that biases the estimates, especially along active braided channels that change quickly. A flux rate must be determined at some location along the channel to propagate transport rate estimates. Martin and Church (1995) also used prismatic approximation (cf. Goff and Ashmore, 1994) to estimate volumetric bed changes between cross-sections along a larger wandering channel. Since the study river had an assumed zero transport boundary, the authors were able to extend reach-scale transport beyond the bar scale estimates following a sediment budget approach. However, the calculation of negative transport rates for some dates demonstrated that even under more ideal circumstances, the approach still provides only a lower-bound estimate of actual transport. A violation of the zero transport assumption, imprecision in the volume calculations or unreported dredging were suggested as possible explanations for the unexpected results. Since the downstream end of the gravel reach of lower Fraser River is also characterized by zero gravel transport (transition to a sand-bed river), McLean (1990) was able to extend - 16-  transport calculations upstream from this reference point by 2-km long cells (reaches) and illustrate spatial variations in the gravel transport rate over a 32 year period. Although the survey interval is large, there is no practical upper limit on rivers with no downstream transport since all material becomes trapped. Although this eliminates the problem of unmeasured bed material throughput that otherwise occurs, compensating scour and fill between surveys masks shorter-term details of spatial and temporal trends in bed material transport within the reach.  The active pursuit of topographic data collection at increasing spatial and temporal scales has provided unprecedented views of form and process feedback relations along alluvial rivers. This advance, however, has been largely restricted to relatively short lengths of fairly small channels. To the extent that these approaches can be modified for, or extended to large rivers remains unclear. Although Lane et al. (2003) and Westaway et al. (2003) have demonstrated progress for capturing fully-distributed three-dimensional topography on a large braided river, the techniques have not been tested on rivers with greater vertical relief and greater flow depths. It is also important to recognize that survey comparisons are limited by the accuracy, density and availability of the poorest quality data set. On very large channels like lower Fraser River, limitations include aerial photographs taken with non-metric cameras and bathymetric data surveyed at low-density along widely spaced cross-sections. However, we can not ignore these data because large channels are modified over timescales that necessitate the inclusion of such information if we are to arrive at a meaningful understanding of channel morphodynamics at appropriate scales of investigation. This perspective is developed in the following chapters.  - 17-  Chapter 2: Physical setting 2.1 Introduction The Fraser River drainage basin encompasses an area of 233,000 km , nearly one-fourth the land surface area of British Columbia, and is the largest watershed located essentially within provincial boundaries (see Figure 2-1). The headwaters are located near Yellowhead Pass along the British Columbia / Alberta border, from where the river flows northwest along the Rocky Mountain trench to Prince George, south to Hope townsite in a deeply incised channel through the Interior Plateau and west through a broad alluvial plain to its delta on the Pacific Ocean, a total distance of nearly 14X10 km. Along much of its length upstream of Hope, the course of the river is structurally controlled and lateral movement is commonly restricted by bedrock or Pleistocene terraces. Below Hope, however, the channel becomes largely unconfined, coarse sediments are deposited as gradient declines and the river has formed a laterally-active wandering channel that extends downstream to Mission. Sediment transfer and deposition in this gravel-bedded reach of Fraser River dominates the stability, hence morphology, of the channel and adjacent floodplain, and presents a variety of management concerns. This chapter provides an overview of basin physiography and hydrology. Flow and sediment regimes within the gravel reach are a reflection of upstream controls established [chiefly] by the Quaternary history of the province, and are manifested in the contemporary channel morphology. The remainder of the chapter examines anthropogenic influences on channel morphology within the gravel reach - activities that have altered the natural riverine environment over the past century.  2.2 Physiography and hydrology The Fraser River basin lies wholly within the Canadian Cordilleran Region, one of five major physiographic subdivisions of Canada (Holland, 1964) and has been further sub-divided into major units including the Rocky and Columbia Mountains, the Interior Plateau, the Coast and Cascade Mountains and the Fraser Lowland. The basin was repeatedly covered by the Cordilleran Ice sheet during the Pleistocene Epoch and thick deposits of Quaternary sediments were deposited in main valleys and lowlands during both glacial and nonglacial periods (Clague, 1989). The most recent (late Wisconsinan) Fraser glaciation terminated about 13,000 years ago (BP) (Armstrong, - 18-  Figure 2-1. Location map and boundary of Fraser River watershed. Principal hydrometric stations are shown. The gravel reach extends from station 7 to 8. Base maps from B.C. GIS warehouse (Ministry of Sustainable Resource Management).  -19-  1981). Fraser glaciation drift consists mainly of till overlain and underlain by glacial fluvial and glacial lacustrine deposits, while glacial marine sediments are common in some coastal areas (Ryder and Clague, 1989). Repeat cycles of sedimentation followed by partial or complete excavation have resulted in a complex reworking of Quaternary sediments in most lowlands and valleys, making it difficult to correlate discontinuous statigraphic units in the absence of dating techniques (Ryder and Clague, 1989). Following deglaciation, mass wasting and fluvial processes began to rapidly redistribute the vast volumes of unstable, unvegetated drift that were left behind. Initially, this resulted in the formation of large alluvial and colluvial fans and aggrading valleys, while finer sediments were commonly transported and deposited in lacustrine or marine deltas and offshore sediment basins (Ryder and Clague, 1989). As slopes began to stabilize, vegetation was established and the supply of drift declined. Many streams, including Fraser River and its major tributaries, incised into their former floodplains creating deep valleys cut into Quaternary fills (Clague, 1993) and producing relict terraces and dissected fans (Ryder and Clague, 1989). The rate of incision along lower stream courses may have been accelerated by the effects of isostatic rebound following final ice retreat (Holland, 1964). An observed downstream increase in sediment yield within the Fraser River watershed (Church and Slaymaker, 1989; Church et al., 1999) demonstrates that degradation of Quaternary valley fill in the Interior Plateau (cf. Church, 1990) continues today. These sources produce mainly fine grained sediments as banks are eroded during high flows (Church and McLean, 1994). Annual total suspended yields are strongly correlated to the size of the freshet, and average roughly 17 million tonnes below Fraser Canyon (McLean et al., 1999). In comparison, upland regions are not prolific sediment sources, as these areas are poorly connected to the channel network (Church, 1990). The major alluvial reaches of Fraser River are located within the Fraser Lowland, a triangular region of roughly 3000 km that extends westward from Hope and forms the southwest 2  corner of the province, including part of northwestern Washington State (Holland, 1964). The valley preserves evidence of three intense late Pleistocene glaciations that left deposits of unconsoliated Quaternary sediments up to 300 metres thick above Tertiary and Cretaceous bedrock. Fraser River forms the dominant geomorphologic feature of the lowland, occupying a late glacial and post-glacial valley up to 5 km wide and 225 metres deep fronted by an extensive, growing delta (Armstrong, 1984). East of Mission, the valley is comprised mainly of Salish and Fraser River (Holocene) sediments, including sand, gravel and overbank fines. West of Mission, -20-  sand and silt become the dominant floodplain material. Aggradation on the floodplain continues today, probably at a lower rate than during the early Holocene (Ryder and Clague, 1989; Church, 1990) but there are almost no contemporary records of sedimentation for comparison (Church and Ryder, 1972). Further, since the floodplain is almost entirely isolated from the river, fine sediment deposition is largely restricted to vegetated islands. From Hope downstream to Laidlaw, the river is nearly continuously confined by bedrock, landslides or Pleistocene terraces (McLean et al., 1999). Below Laidlaw, the river exhibits a wandering channel planforrn within a partly confined floodplain that extends 55 km downstream to Sumas Mountain (Church and McLean, 1994). The deposition of gravel in this reach is a consequence of a reduction in channel gradient below the steep confines of the Fraser Canyon. The terminus of gravel deposition is marked by an additional (order-of-magnitude) decrease in gradient and transition to a single-thread sand-bed river. The typical thickness of the active gravel layer is estimated as 5 to 15 metres along much of the reach (Armstrong, 1981). GPR profiles of large bar complexes within the reach (cf. Wooldridge, 2002) show alluvial sediments 8-24 metres thick, underlain by steeply dipping deposits interpreted as former positions of a prograding delta front following deglaciation. Similar measurements on the wandering Squamish River (Wooldridge, 2002) provide additional evidence to support Desloges and Church's (1987) assertion that the wandering planforrn appears vertically stable. Although the gravel reach of Fraser River is known to be slowly aggrading at present, channel behaviour, hence morphology, appears to be dominated by lateral migration within the floodplain. Therefore, anthropogenic modifications such as bank hardening may have important implications for both recent and future channel development along the reach. Despite extremes in climatic diversity that arise from the watershed's large and varied geography, the hydrologic regime is strongly characterized by nival flows. Elevations in the watershed range from near sea level at the delta and lower Fraser Valley to 3600 metres in the headwaters. Hypsometric calculations reveal that half of the watershed lies at elevations above 1150 metres, while 1.2% of the watershed remains glaciated with the greatest concentration in the Coast and Cariboo Mountains (Figure 2-2). Although autumn and winter maritime air masses produce heavy rains at lower elevations, snowfall is common at higher levels and interior locations where these masses encounter colder arctic air. Two-thirds of annual precipitation in the basin falls in the form of snow (FBMP, 1994) so widespread, prolonged melting of the accumulated snowpack dominates runoff patterns (Moore, 1991). Large rainstorms can increase runoff in the -21 -  Figure 2-2. Relief map of Fraser River watershed. The location of existing glaciers and a hypsographic curve are shown for reference. Data are derived from Provincial 1:250,000 base mapping.  -22-  summer to autumn period and even prolong the freshet (i.e. as occurred in 2001) but rarely generate flows at or near bankfull discharge. The annual freshet occurs from May to July (median date is June 13) and commonly takes the form of a single-event hydrograph (Church and McLean, 1994; see example, Figure 2-3). Most gravel transport, hence morphologic change, occurs during this time. The mean annual flood at Hope, site of the longest continuously operating gauge in the watershed (since 1912) is 8750 m /s. For comparison, the mean minimum discharge is 600 m /s 3  3  the ratio between them (14:1) is typical of large rivers and reflects the attenuation of snowmelt and rainfall inputs over large areas (cf. Dunne and Leopold, 1978). Mean annual flow is 2700 m /s, 3  while major tributaties, including Coquihalla, Chilliwack and Harrison Rivers, increase discharge to 3350 m /s at Mission, located at the distal limit of the gravel reach . The known flood of record, 3  1  in 1894, was estimated from local high water marks as 17,200 m /s (at Mission) and the measured 3  flood of record was 15,200 m /s in 1948 (at Hope). Other major floods that pre-date the gauging 3  record occurred in 1876, 1882, 1896, and 1903 (Johnston, 1921). The historic sequence of annual peak flows is shown in Figure 24A. The largest floods occur when cool spring temperatures give way suddenly to extended warm periods, such that melt occurs quickly rather than over a protracted period. These conditions were the principal cause of the 1894, 1948 and 1972 floods, with significant rain amplifying the 1948 and 1972 events. A large snowpack with high water equivalent, above-freezing overnight temperatures at high elevation and heavy rainfall have been identified as additional factors which increase the rate and volume of runoff during the melt phase (FBMB, 1994). Discharge in the watershed is partly regulated by a number of storage facilities for hydroelectric power generation. The largest of these is a diversion located on the Nechako River (completed 1952), which effectively reduced the size of the basin by 14,000 km (5.6%) and mean 2  annual flows by 3% (Moore, 1991). Smaller hydro reservoirs are found on the Bridge River system, lower Fraser Valley tributaries and Mabel Lake. During the 1972 flood, it is estimated that that utilized storage capacity at the Nechako and Bridge River facilities reduced the peak flow at Hope by 1130 m /s and peak flood crest elevation at Mission by 0.3 metres (FBMB, 1994). The 3  attenuation effects of water storage have contributed to a reduction of the mean annual flood by  1. Annual maximum discharge at Hope ( W S C gauge 08MF005) is calculated for the period 1912 to 2002, while mean annual flows are availabile since 1913. The Mission gauge (08MH024) has been in continuous operation since 1965 but the reported figure is based on published data available for the period 1966-1992 inclusive.  -23-  14000  Figure 2-3. Example annual hydrographs, Fraser River at Hope, 1972 and 2001.  -24-  13% over the 50 years since 1952, though natural long-term climatic variability confounds such a simple conclusion. Over the complete record, flood flows have declined nearly 5% while, surprisingly, the mean annual discharge has actually increased by the same proportion. This trend may in fact be a signal of glacial retreat throughout the watershed (D. Moore, pers comm., 2004). Individual floods are not as important in influencing channel morphology on large rivers as are longer-term trends in flood flows. McLean (1990) found that individual floods may accelerate local changes such as bank retreat, but do not directly influence large-scale morphologic patterns, because the volume of material stored in the river is so large relative to the capacity of the flow to erode and transport the stored sediment during a single freshet. This effect can be observed qualitatively by comparing low-water photo series in 1943 and 1954 (McLean, 1990) and from the early 1940s to 1949 (Church and Ham, 2004) — these show a surprising lack of major channel shifting and erosion compared to that observed since 1928 despite the magnitude of the 1948 flood (examples are provided in Appendix B). In order to better understand long-term patterns of channel behavior, it is neccessary to relate observed changes to long-term trends in flood flows. Plots of cumulative departures for the annual maximum and mean flow series are given in Figure 2-4 B. Mean flows were generally above the long term average to 1922, below average to 1952, and above average to 1977. Flows have generally been below average since, but there is limited evidence to suggest that they have been increasing since the mid-1990s. These dates reveal shifts in climatic regime, and are thought to coincide with the Pacific Decadal Oscillation (Mantua et al., 1997, in McLean and Church, 1999). Flood flows follow a broadly similar pattern, but dates of major regime changes do not conform with mean flows and intervening flow trends exhibits greater variability. Nevertheless, it is expected that major changes in bed material transport, and hence changes in channel morphology, should correspond with trends in flood flows. 2.3 Anthropogenic influences on morphology The lower Fraser River remains today a relatively pristine environment compared to many large rivers that flow through heavily urbanized regions. Nevertheless, efforts have continued since the late 1800s to protect developed areas and infrastructure along the river from flooding and erosion through a system of dykes, bank protection structures, channel dredging and removal of woody debris. The known history of these developments is summarized in the following sections. It is probable that these efforts have influenced morphologic development over this same period. The sensitivity of channel response to anthropogenic disturbance is examined in Chapter 4. -25-  16000  14000  12000  10000  e?  8000  CJ  6000  4000  2000  y = 1.422x- 67.23 0 1910  1920  1930  1940  1950  1960  1970  1980  1990  2000  2010  Year 15000  10000 3  ts OS  a, u  T3  5000  o  53  3  B  6  -5000  u > -10000 3  -15000 1910  1920  1930  1940  1950  1960  1970  1980  1990  2000  2010  Year Figure 2-4. Annual maximum and mean daily discharge, Fraser River at Hope (08MF0O5). A. Annual flow series with best-fit trend lines. B. Accumulated departure series showing major trends in mean and flood flows over the period of gauging. An ascending plot indicates periods of persistently above-average flows, a descending plot indicates periods of persistentiy belowaverage flows, and a level plot indicates flows near the long-term mean. See text for discussion. -26-  2.3.1 Flood protection structures  Much of the initial European settlement of lower Fraser Valley was established on lowlying areas along the Fraser River floodplain, largely because of the high demand for suitable agricultural lands. These areas were also at risk of flooding during freshet in most years. North and Teversham (1984) note that early farmers recognized the need to lower the water table and protect their planted crops and grazing lands from damaging floods and many responded by building small dykes along their property (Sewell, 1965). Following a major flood in 1876, an act was passed [1878] to build privately-financed dykes at Chilliwack, Sumas and Matsqui in exchange for land, but the project was discontinued (Winter, 1966). Also in 1878, a separate plan was initiated to drain the seasonally-flooded Sumas Lake in an attempt to reduce the area of floodplain that was inundated, create additional arable land, and open east-west travel corridors. However, the project was not completed until 1923 (Siemens, 1966). In 1892, earthen dams were built to close several small sloughs and the upstream end of the major slough at Nicomen Island (Public Works Canada, 1949). Additional dams were built at Hope Creek and Camp Slough near Chilliwack (Sinclair, 1961). Despite these efforts, the dyking system was discontinuous, often built without engineering advice, and many sections had been abandoned during construction (Sewell, 1965). The devastating 1894 flood emphasized these problems, and spurred efforts to protect newly developing areas from future flooding by establishing an integrated system of dykes. Through substantial government assistance, most flood-prone areas had been protected within ten years (Siemens, 1966). There were no significant improvements to the dyking system in following decades, although additional flood dykes were built to support the draining of Sumas Lake. This inaction was largely due to financial difficulties, as most local municipal districts did not have sufficient revenues to both service accumulated debts and finance dyke construction and maintenance. As a result, many of the dykes were never completed to recommended engineering standards and commonly fell into disrepair (Sewell, 1965). Fortuitously, there was a lack of damaging floods during this same period, and this circumstance likely reinforced the status quo. A major program to repair and upgrade existing dykes was, however, initiated following the disastrous 1948 flood, which breached and damaged much of the existing flood-protection framework, resulting in considerable economic loss and disruption. The upgrades were sufficient to withstand a second, less severe flood in 1950, but engineering tests revealed that the desired design elevation was not being met (Fraser River flood control program, 1968). As a result, another major program of dyke -27-  upgrades and flood control works was initiated under the plan of protection outlined in the 10-year Federal-Provincial Flood Control Agreement of 1968. Much of the work in the gravel reach involved new or improved bank protection works where lateral erosion threatened to undermine dykes or cause loss of developed lands. The program also provided funds for improvements and additions to the pumping / ditching system designed to provide internal drainage of lands enclosed within the dykes (Fraser River flood control program, 1968). Notably, there were no planned extensions to the existing 130 km of flood dykes within the gravel reach and, since then, the only major addition to the system has been a 3.4 km branch dyke along Island 22 near Chilliwack. The current network of flood dykes is shown in Figure 2-5. Along much of the floodplain, the dykes are set back a sufficient distance from the margins of the active channel zone such that they do not directly interfere with morphologic development of the river. One of the main impacts has been the separation of many prominent sloughs and backchannels from the main channel. McLean (1990) notes that several of these have since filled with finer sediments and vegetation has encroached. This has resulted in a loss of aquatic habitat, and water once conveyed through these secondary channels at high flows has since been forced back into the main channel. The primary impact of the dykes on channel morphology has been indirect. As the river has long been recognized as being laterally unstable, bank protection works have been used to restrict erosion where the dykes are vulnerable. 2.3.2 Bank hardening  The construction of bank protection structures within the gravel reach dates back to 1892 when crude buttresses of brush and rock were sunk along wood pilings to limit erosion at Miller's Landing near Sumas (Public Works Canada, 1949). By 1910, additional similar structures had also been placed at Sumas and Matsqui, while a wing dam was built at Chilliwack. In 1914, severe bank erosion along Nicomen Island resulted in several attempts to limit erosion using brush, rock, and pilings, but all measures were unsuccessful. The erosion was not stopped until a channel was dredged to divert flow away from the sharp developing bend, although the diversion channel was not actually adopted by the river until 1922 (Public Works Canada, 1949). Other major projects included placed rock at Rosedale (1928 to 1939) and installing pile groins followed by placed rock near Agassiz (1928 to 1948). As McLean (1990) points out, there is little specific documentation of these early protection measures. Much of the work completed to protect the dykes and settled land from bank erosion likely took place after the late 1920s, when the federal government made its first financial contribution to flood control (Winter, 1966). By the date of signing of the 1968 -28-  S u m m a r y of bank hardening between Mission and Hope  1. 2. 3. 4. 5.  Reach  Total bank length (m)  Railway  Sumas  31 925  4715  Chilliwack  35 739  Rosedale  26 346  Cheam Hope  1  2  Dykes  Riprap  3  Bedrock  4  Total protected  Hope reach 5  % protected  13 760  4 873  23 348  73.1  14 898  6 525  21 596  60.4  4 091  14 202  833  19 126  72.6  42 501  18 358  11 739  30 097  69.2  32 202  6 207  8 692  17213  53.5  173  2314  Hope  outer banks of main channel only; does not include island shoreline railway protects many banks that otherwise would be classified as bedrock riprap includes 2907 m of rock berms, mainly in Cheam reach bedrock includes non-alluvial, non-erodible banklines (eg. Mission bend) all categories are exclusive (i.e. any length is counted in one category only); total protected (and % protected) includes bedrock  Bank Hardening Railway  (\J r\j  Riprap Dykes bedrock  Mining Volumes (m ) + < 50,000 • 50,000 - 100,000 3  Hassiz I Rosedale bridge  •  •  •  •  100,000 - 200,000 200,000 - 400,000 400,000 - 800,000 > 800,000  Channel Feature (1999) water surface floodplain  Gravel site numbers correspond to locations for which removal volumes are known (see Weatherly and Church, 1999)  gravel bar reach boundary  0 Figure 2-5. Location and sv  2  6 km  raening along lower Fraser River. The location of known gravel mining sites is also shown -29  Federal-Provincial Flood Control Agreement, 36 km of shoreline were protected by riprap (McLean, 1990), of which 13 km was slated for upgrading and an additional 9.4 km was proposed (Fraser River flood control program, 1968). Roughly 55 km of bankline are currently protected, including 1.9 km along Peters Island (Figure 2-5). Bank protection structures have also been used throughout the gravel reach to protect railway lines from bank erosion. Construction on the transcontinental Canadian Pacific Railway began in 1880 near Yale, then the head of navigation on lower Fraser River, and by 1885 the line extended downstream to Port Moody (Roy, 1966). Additional railways established several new lines in the Fraser Valley at the turn of the century to provide service to communities along the southern shore of the river, but these generally extended north-south, or ran close to the U.S. border. A second transcontinental railway (Canadian Northern Railway, now Canadian National) was completed in 1914 and connected Hope with New Westminster along the southern side of the river (Roy, 1966). Both the CP and CN railway lines parallel the river through lower Fraser Valley and run adjacent to the active channel in several locations (see Figure 2-5). Since the railways are vital transportation and economic corridors, the rail bed has been heavily engineered with riprap to protect against undercutting or lateral erosion by the river. It is estimated that the railways currently protect 33 km of floodplain from erosion. Upstream of Agassiz, the railway also functions as a flood dyke since the rail bed has been built to sufficient height along the floodplain to maintain rail transport during high-water events. Lateral channel movement is further restricted at several locations by natural phenomena. Between Hope and Laidlaw, the river is almost continuously confined by bedrock, landslides and Pleistocene terraces (McLean et al., 1999) and there is little adjacent floodplain. Downstream of Laidlaw, bedrock is found adjacent to the active channel margins along the Skagit Range near Herding Island and at Hopyard Hill, Mount Woodside, Harrison Knob, Chilliwack Mountain and Sumas Mountain. Terraces are also found at Mission, near the downtream limit of the gravel reach. The shape and extent of the modern floodplain, therefore, is strongly affected by the presence of these features (it is effectively a confined alluvial fan). This circumstance is not unusual along wandering gravel-bed channels (cf. chapter 1). The location of natural constrictions (for mapping convenience, all are plotted as bedrock) is also shown in Figure 2-5. In total, 14.5 km of channel are protected by natural constriction or confinement. A summary of bank hardening, including the percentage of outer banks that are currently protected, is summarized as an inset table in Figure 25. As it is known that the net effect of this bank hardening has been to narrow the width of the -30-  channel zone (the average width between outer banks) along all reaches (cf. Ham and Church, 2002) major adjustments in active channel width must result from erosion or construction of channel islands. Forced channel narrowing may also lead to deeper, faster currents in order to maintain flow conveyance, which has implications for sediment transport capacity, hence overall channel stability. This possibility is reviewed in Chapter 5. 2.3.3 Sediment removals  Sediment has been dredged from the channel since at least the start of the 20th century when steam-operated hydraulic suction dredges became available. The purpose of these removals was to maintain a safe navigable channel for shallow-draft vessels (that once extended upstream to Yale). Earlier removals may have taken place, but the volumes would have been very modest given the limited available means that would have existed to physically dredge and transport the material. Public Works Canada (1949) reports that 3.77 million cubic yards (2.88 million m ) of material 3  (bulk volume) were dredged between the delta and Chilliwack during the period 1901-1910. Given that these removals extended well into the gravel reach, and the diameter of the suction pipe (20 inches, or 508 mm) exceeds that of the largest mobile clasts in the channel, some gravel is likely to have been extracted along with sand andfinermaterials. The records are nonetheless inadequate to estimate the proportion of the total volume removed that was taken upstream of Mission. The only other documented early extraction is the 393,000 m dredged between 1914 and 1922 for the 3  diversion channel at Nicomen Island, but this material was replaced within the channel. Gravel was also dredged from the river during the construction of railways and roadways where this infrastructure was located in close proximity to the active channel (several examples are documented in Weatherly and Church, 1999). Some gravel may also have been obtained for dyke construction, but borrow pits were usually located on the landward side of the dykes (Sewell, 1965). Large scale channel removals certainly took place for dyke upgrades following the 1948 flood, however, since extraction pits and heavy trucks are directly visible on 1949 aerial photographs (eg. near Minto Channel). However, as there was no documentation of these removals, it is not possible to reliably estimate the total volume taken. The earliest summary inventory of gravel extraction in the river was completed by Kellerhals Engineering Service Ltd. (1987) for the period 1973-1986. The report gathered information given by local gravel operators and representatives from various government agencies. Weatherly and Church (1999) note that, as gravel extraction has been regulated since 1974 by British Columbia Assets and Land Corporation, the Kellerhals inventory appears fairly accurate. -31 -  However, the authors caution that the records may be incomplete because of (1) inaccurate selfreporting of amounts removed by gravel operators, (2) the lack of documentation of gravel removals by the Ministry of Transportation and Highways, which had established gravel reserves but was not subject to regulation, and (3) unregulated removals on private lands. The latter includes removals on Indian Reserves, for which records appear non-existent. The completeness of extraction records has improved since 1980, when the Department of Fisheries and Oceans first required permits for any removal because of fisheries concerns (McLean, 1990; Weatherly and Church, 1999). As well, pre- and post-extraction surveys were sometimes requested to verify removal volumes, adding credibility to reported figures. Gravel removal volumes presented herein are primarily based upon the 1999 Weatherly and Church report, which provides an update to the existing Kellerhals inventory, and extends the reporting period from 1964 to 1998. As the authors were able to contact a number of individuals, additional information was obtained for pre-regulation commercial operations and removals related to road and dyke construction. A second report (Kerr Wood Leidal, 2001) was used to extend the inventory to 1999 and amend reported [1998] values at Hamilton Bar and Powerline Island. A total of 4.2 million m of gravel (bulk volume) is documented to have been removed 3  between Mission and Laidlaw during the period 1964-1999, an average of 120,000 m per year. 3  The total extraction increases to 4.8 million m if extended upstream to Silverhope Creek near 3  Hope, or 138,000 m per annum. Although far more gravel was removed in 1995 (428,000 m ) than 3  3  in any other year of record, maximum removal volumes generally took place from the early 1970s to the early 1980s. The dominant location of these removals has been along Minto Channel, though large removals have also occurred at Foster Bar and Powerline Island. A nearly complete record of known gravel extractions between Mission and Hope is summarized in Weatherly and Church (1999) and is also reported in the sediment budget tables (Chapter 5). Corresponding site locations are shown in Figure 2-5. Of particular interest in this thesis are those removals which took place between the years 1952 and 1999, as this corresponds to the period for which estimates of sediment transport are given (see Chapter 5). The incomplete information on gravel mining prior to industry regulation, and missing information on removals between 1952 and 1964 represents a negative bias for estimates of total bed material transport in the gravel reach. The magnitude of this bias remains unknown, but may be small, given that individual recollections of these early extractions typically report modest volumes (cf. Weatherly and Church, 1999). -32-  While scuffle dredging along shallow riffle zones was until recently practiced to maintain navigational safety, most current removals are conducted as commercial operations. Sand and gravel is taken from targeted locations where deposition is thought to cause local increases in water-surface elevation. The sediment is sorted to eliminate smaller grain sizes, and the gravel is sold for profit, both offsetting removal costs and providing crown royalties. The success of this practice in actually reducing water levels for a given discharge on Fraser River has been examined by UMA (2002, 2004). It was found that reported sand and gravel removals from 1952 to 1999 lowered the flood profile up to 20 cm between Chilliwack and the Agassiz-Rosedale Bridge, but that channel alignment appears equally important in determining the flood profile. 2.3.4 Woody debris The removal of woody debris (snags) from the gravel reach dates to the early 1870s when the river was a principal conduit for transportation and commerce. A total of 2300 snags were removed by 1894 to improve navigational safety along the lower 135 kilometres of channel (approximately to Agassiz). An additional 35,400 snags had been removed by 1949 between the delta and Harrsion River confluence (Public Works Canada, 1949). Despite these efforts, woody debris jams still acumulated along the upstream margins and distributary channels of major bar/island units, as visible in historic aerial photographs. By 1979, a boom to trap floating debris was established near Hope, replacing the boats that had previously been used. Hales (2000) reported that the trap has caught 90,000 m of wood annually since that time, roughly 90-93% of all floating debris in the 3  river. The continuous operation and greater efficiency of the trap has since resulted in a steady decline of deposited wood in the gravel reach. At present, both the gravel reach and the tidal flats (cf. Hales, 2000) are largely free of woody debris. Since woody debris is known to be an important nucleus for the establishment, growth and maintenance of channel islands and side-channels, there may be observable changes in both over the past two decades. The morphologic and ecologic implications of this practice are explored in Chapters 4 and 6.  -33 -  Chapter 3: Data Compilation 3.1 Introduction This chapter provides an overview of the primary data used for mapping and measuring morphologic changes and estimating sediment transport rates within the gravel reach of lower Fraser River. The pattern and rate of morphologic changes over time is measured from a sequence of historic maps and aerial photographs available for various dates since the early 1900s. Patterns of aggradation and degradation within the gravel reach are determined by comparing topographic surfaces of the channel bed and banks from hydrographic surveys completed in 1952, 1984, and 1999. The quantitative interpretation of this information is given in Chapter 4, while the links between observed changes in channel deformation and the processes that produce them are discussed in Chapter 5. A complete description of data collection and analysis procedures is provided in the following sections. 3.2 Channel a n d floodplain mapping 3.2.1 Data acquisition from old maps The availability of old maps and charts represents a valuable resource from which to investigate the morphologic evolution of river channels. Their use is particularly relevant for large river systems, where significant changes may occur over timescales that exceed the temporal limits of available aerial photography (Lawler, 1993). By incorporating available maps, studies of channel migration and development can be extended by decades or even centuries beyond the photographic record. Braga and Gervasoni (1989) provide an exceptional example of chronologic evolution on the Po River, Italy, where migration of the channel within its floodplain was dated back to the Dth century at several locations. Similarly, Decamps et al. (1989) were able to track changes on the Garonne River at the city of Toulouse, France from the 17th century, while Sundborg (1956) was able to infer local details of bank erosion on farmland on River Klaralven, Sweden from local tax records over a similar timeframe. These studies demonstrate that archival records provide a means for evaluating the long-term stability of the magnitude and frequency of dominant channel forming processes (cf. Hooke, 1980) and can illustrate the nature of channel adjustment in response to natural (Vandenberghe and Maddy, 2000) or anthropogenic (Hooke and Redmond, 1989) influences. The availability of old maps is, however, commonly restricted to settled areas, limiting applicability at the reach or larger scale.  -34-  Quantitative assessments of channel changes can reasonably be made since the development of cadastral surveys in the 18th century (cf. Church, 1983b). However, the quantitative use of old maps is complicated both by the interpretive biases of the surveyors and cartographers who prepared them, and by the accuracy of the instruments used in map construction (cf. Leys and Werrity, 1999). Depicted channel features such as gravel bars and islands may also change over time due to changes in surveying practices, a concern further complicated by varying standards adopted by different mapping organizations (Hales, 2000). In practice, this means that maps may not only be inaccurate, but may not show all features that are of interest, or may show particular features in insufficient detail where such features are incidental to primary mapping themes. The level of detail also varies with map scale, since some generalization of information is neccessary when abstracting reality to a cartographic representation (cf. Vitek et al., 1996). Veregin (1988) adds that this process can further introduce conceptual (fuzzy) errors which affect the accuracy of spatial overlays. The use of old maps in morphological studies can also be compromised by an inaccurate scale, a lack of referencing grids, or scale changes between different dates. These problems can be resolved by georeferencing all maps through the identification of common ground control points, a procedure which relates all data to a common frame of reference (Burrough and McDonnell, 1998). Positionally stable (fixed) points are assigned a set of cartesian coordinates which can be obtained, for example, from field surveys, planimetric maps, orthophotos, stereo aerial photographs, or differential GPS. The location of the digitized points is subsequently compared to the mathematical position of the assigned coordinates and a least squares fit is calculated. Large errors are usually caused by the inclusion of one or more incorrectly located points (or incorrectly recorded coordinates) but can be reduced by iteratively removing individual points until the 'fit' improves. Since older maps commonly suffer from unequal paper stretch (Burrough and McDonnell, 1998) an additional transform is generally required before any measurements can be collected. This step is known as "rubbersheeting", and is the mathematical stretching, rotating and translation of points along vectors linking the observed location of fixed points to their true location. Despite potential difficulties, even qualitative assessments of gross planform changes such as channel shifting or the permanence of islands and major deposition zones can provide information about process changes on large rivers. A recent study on the Brahmaputra River (EGIS, 1997) compared the position of the modern channel to that depicted on a British admiralty -35-  chart, revealing tens of kilometres of migration in response to tectonic effects nearly three centuries previous. Strahler and Strahler (1991) provide an example of bend migration and chute cutoff on the Mississippi River from 1765 to 1930, while Winkley (1994) notes that details shown on these same maps have been used to infer the effects of tectonics on river meandering and sedimentation. Even in cases where inferences of process changes can not be directly qualified, historic planform studies may provide valuable corraborative evidence to support other lines of inquiry. Archival records in Western Canada are not nearly as extensive as in Europe. The first complete map of lower Fraser River dates to 1859, and merely provides a coarse snapshot of channel and floodplain conditions at that time. Although the detail shown is insufficient to characterize any channel changes that have since occured, the map usefully confirms that several of the major island groups within the river have been present for at least 150 years. More detailed maps are available through five early mapping projects in the Fraser Valley, completed using established surveying procedures (North etal., 1977; North and Teversham, 1984). These include the Crown Colony and Royal Engineers town surveys from 1858-1863, the Provincial surveys from 1873-1877, the Dominion surveys from 1873 onwards, and the surveys of Indian lands. However, these surveys did not generally extend into the channel beyond the main river banks, so the condition of major deposition zones can not be assessed. An exception is the Indian land surveys, which were almost exclusively of channel islands (North et al., 1977), but these surveys do not encompass the entire lower gravel reach. North and Teversham (1984) compiled a vegetation map of the islands and floodplain of lower Fraser River from the original land surveyor's ground coordinates, notes and sketches. Whilst the map represents the period 1859 to 1890, the authors note that most of the floodplain was surveyed between 1872 and 1878, a reasonably short timeframe on a large river from which a general morphologic assessment can be made. The reliability of this map must be questioned though, because the number of different surveyors who originally worked on the project introduces problems of interpretive bias and consistency in data collection practices which may have been translated into the final production map. The first available map showing the complete gravel reach in sufficient detail to allow morphologic comparisons with more recent mapping was first published in 1905, and a modified version was reprinted as the "1913 Map of New Westminster District." The map was revised by land surveyors in 1912 from 'official sources' which appear to be the original township surveys, given the similarities to several maps shown in McLean (1990). The map provides a detailed -36-  cadastral survey of sectional boundaries and property lots within the lower Fraser Valley. Since land title extended to active channel margins and included channel islands, it can be used fairly reliably to measure active channel width or to compare with more recent maps and photographs to calculate rates of bank advance or recession. Exposed gravel bars are also shown in locations where deposition would typically occur, such as where channel width increases or along the margins of larger islands. Additionally, the size, shape, orientation and variety of bars are consistent with established sedimentation patterns on wandering rivers, and on Fraser River in particular. Nonetheless, several areas within the channel, including lower Herrling Island and present day Queens Bar, are shown as having no accumulation of gravel at all, which is improbable. It is unknown whether the inconsistent mapping of channel bars is due to interpretive differences amongst surveyors, or simply reflects river conditions during high water periods when no bars would be visible. Although the map appears to reliably depict channel banks and island positions, the representation of channel bars is therefore ignored. The 1912 map lacked both a coordinate system and a suitable number of well-spaced identifiable ground control points to allow rubbersheeting. However, a two-mile cadastral grid was drawn over most of the map surface area, providing a means to scale the map. Since the map was too large to be digitized in one piece, a set of 6 reference coordinates were arbitrarily established along both sides and the map centre (hence it was registered twice using common central coordinates). The location of outer channel banks and all channel islands was digitized along with the cadastral grid, then imported into the GIS for further processing. The map was scaled to realworld coordinates by converting the length of the 2-mile grid in millimetres to the representative length in metres. An artificial grid (2 miles by 2 miles) covering the same area as the map was generated and the scaled map was translated to overlap the dummy grid. Observed differences between the two grids clearly showed locations of distortion through paper shrinkage or expansion. To correct this, vectors were digitized at each grid intersection from the digitized map to the dummy grid, and a rubbersheeting transform was applied. The final step involved rotation and further translation to align the digitized map with a more recent (and positionally accurate) photogrammetric map, using the railway, Harrison River, backchannels and bedrock outcrops as a guide. Despite these efforts, correspondence between the 1912 map and the modern channel could not be completely resolved for all areas. In particular, that section of the map between Vedder Canal and Agassiz-Rosedale bridge did not align as closely as the remaining sections. Although the 2-mile grid was coincidentally absent over this same area, differences are attributed to -37-  cartographic errors introduced during map construction since remaining displacement usually exceeded that associated with paper shrinkage observed elsewhere. The logistical difficulties of accurately surveying this section of river likely reduced the accuracy of cartographic representation relative to other mapped areas. Consequently, the 1912 map is primarily used for qualitative assessment of major channel changes only (Appendix B). 3.2.2 Data acquisition from aerial photographs  When available, aerial photographs should be used in preference to maps to minimize intrepretive or omission errors. They are available much more frequently than maps, an important consideration on dynamic rivers. Though primarily a two-dimensional resource, a wide variety of historical information can be derived from photographs depending on operator knowlege and available measuring devices. The qualitative use of aerial photographs for interpreting fluvial landforms is well-established, and a number of classification schemes have been developed for associating channel pattern with some combination of different formative variables (cf. Mollard, 1973; Kellerhals et al., 1976; Schumm, 1985). Provided photo scale is known or can be determined, quantitative information such as channel width and sinuosity can also be estimated, but distortion (e.g. due to camera tilt, changes in topographic relief) limits accuracy. By orienting photographs relative to one another through common scaling (cf. Leys and Werrity, 1999), more precise measurements of change can be observed and measured, including the rate of change over time. Church et al. (1987) demonstrate a method for overlaying successive images using an optical projection system. Relative scaling was achieved by matching control points on successive images. This approach has the advantage that external survey control is not required, allowing large numbers of photographs to be analysed rapidly, but small angular deviations introduced during successive photo registration introduces displacement errors relative to 'true' map position. Nanson and Hickin (1986) and Piegay et al. (1997) follow a similar approach, but projected successive images to a reference map. However, even where control points are matched to topographic maps, image displacement due to camera tilting, high topographic relief and radial distortion can introduce positional errors that are not easily rectified. More recently, computerized techniques have become widely available for georeferencing individual photos, either printed images on a digitizing tablet or scanned images on a computer screen, using simple image processing tools. The procedure involves registering control points by relating their location to known coordinates and computing a best-fit (least-squares) adjustment. -38-  An additional processing step known as polynomial rectification (effectively a rubbersheeting procedure) rotates and positions imagery relative to the known coordinate space, making direct overlays possible. Further, image distortion due to camera tilt or lens distortion can be reduced if sufficient ground control points are identified to re-establish an orthogonal plane (Livingstone et al., 1999). Known coordinates can be determined from ground survey using conventional bearing and distance or GPS survey (Chandler, 1999), the latter technique being suitable for large spatial areas. Topographic maps are more commonly used, but the accuracy of absolute orientation may be compromised both by map scale and the ability to record points at fine resolution using rulers or digitizers. To fully remove the effects of image displacement caused by terrain relief requires the use of orthorectified imagery (geometrically correct image maps). While orthophotos permit more accurate mensuration, the additional cost of having them prepared relative to more simple georeferencing approaches may be difficult to justify, especially for mapping of relatively flat areas like rivers and floodplains. Although there apppears to be no standard approach for calculating georeferencing errors, they are commonly reported using the RMSE statistic. Stated errors range from ±5-10 metres using recent large scale maps and photos (Winterbottom and Gilvear, 2000; Gilvear et al., 2000; Micheli and Kirchner, 2002; Gaeuman et al., 2003), but errors quickly rise as scales decrease and either map or photograph sources become dated. The accuracy with which photographs can be georeferenced is an important consideration in historic analysis because it establishes the resolution with which detectable change can be reliably determined. This limits effective measurement on small channels, over short time intervals, or on systems with modest rates of deformation. Use of the RMSE statistic has been criticized, however, as it is assumed to apply uniformly over the entire modeled surface (Burrough and McDonnell, 1998), masking the difference between the accuracy and precision of individual points (Lane et al., 2003). As the photographic records on lower Fraser River are extensive, they are used as the primary source of information on channel morphodynamic changes in this thesis (see Table 3-1). The earliest photographs of the entire gravel reach were taken in July, 1928. However, as these photographs were taken at a large scale (1:10,000), individual frames often did not contain sufficient fixed points to permit georeferencing. As well, the lack of fiducial marks, high radial distortion in older photographs (the displacement of image points from their true positions due to lens distortion of unknown magnitude; Slater, 1975) and inconsistent stereo overlap makes it difficult to establish relative and absolute orientation parameters. The 1928 photographs were -39-  therefore mapped by mono-restitution on a standard digitizing tablet. Individual frames were scanned and re-scaled to 50% size, then oriented (best-fit by eye) to produce photo-mosaic mapsheets. Stable locations, such as road intersections that could be identified on more recent mapping, were used to georeference the mosaics (registration error averaged roughly ±50 metres RMS). Since many individual frames in the mosaic did not contain any control points, rubbersheeting was used to reduce displacement errors. The 1938 to 1943 photo sets — the combination of which forms a complete record of channel change between the 1928 and 1949 mapping dates — were registered to 1949 mapping (see below) using control point matching tools in the GIS. Registration accuracy was less than ±10 metres RMS for all points, while the calculated least-squaresfitbetween digitized and true coordinate positions was used to eliminate poor quality control points to improve overall results. These errors are small relative to the size of the channel (average width of 900 metres - see Chapter 4) but nevertheless affect the accuracy by which quantitative comparisons of channel position over time (i.e. bank erosion rates) can be made. For comparison, channel parameters which can be measured from individual years of mapping (i.e. channel width, length) are more reliable, since the photos are correctly scaled.  Table 3-1. Historic airphotos selected for morphologic mapping Discharge  Year  Date  Source  Scale  1999  March 20  BCB99001  1:40,000  700  1991  Sept 4  BCB91079  1:50,000  4340  1983  July 19 July 22 July 30  BC83007 BC83008,12 BC83017/18,20  1:20,000  6110 5380 4820  1971  March 19  BC5406/7  1:32,000  800  1962  May 7  BC5042  1:32,000  2940  1949  March 23 March 31 April 6  BC717-20 BC720-23 BC730/31  1:12,000  730 730 980  1943  Dec 5  A7077/78  1:15,000  930  1940  July 15  BC204/5  1:32,000  4870  1938  April 7  A5869/70  1:22,500  750  1928  July 15  A288/89,296  1:10,000  5780  1  1. Scale is approximate only. 2. Mean daily discharge (m /s) measured at Hope gauge. 3  -40-  2  All remaining photo sets were digitized using an analytic stereoplotter, a device which mathematically relates two-dimensional positions on photographs to their real-world threedimensional equivalents using standard equations for interior, relative and absolute orientation (Slama, 1980). Following these procedures, distortion and displacement errors are minimized, while common scaling is achieved through a statistical fit of digitized control points on the stereo image to coordinates derived from a defined reference system. The stereoplotter also has the advantage of providing a magnified, three-dimensional view of channel topography which aids in the accurate interpretation of channel features. Control points for orientation were obtained from 1:70,000 TRIM (Terrain Resource Information Mapping) diapositive airphotos taken in 1987. The TRIM photos have been used to produce 1:20,000 topographic maps of British Columbia, and are marked with up to seven control points, the locations of which were determined through geodetic field survey and aerotriangulation. The stated accuracy of individual pugged points is ± 0.5 m RMS. Control points were bridged from the TRIM photos to the 1999 photos, yielding absolute orientation errors for individual points of 1-3 metres (calculated as approximately ± 2 m RMS). Images from other years were bridged from the 1999 photography or from successively registered years of photography where common, identifiable control points became less distinct because of development (i.e. there are significant changes in land use patterns from 1949 to 1999, making it difficult to identify common control points). Between successive years, plotting errors are estimated as ± 2.8 metres RMS. Absolute accuracy is important in this study since mapped imagery is combined with 1952 and 1984 survey data (cf. McLean, 1990) and more recent bathymetric and altimetry data from 1999 to produce complete floodplain-to-floodplain topographic surfaces of the channel (see Section 3.4). Improper alignment of these data will produce erroneous models, hence spurious measures of bed material erosion and deposition volumes along the gravel reach. Absolute mapping errors for the planimetric mapping are consistent with GPS-based positioning errors reported for both (sounding and altimetry) 1999 surveys, and exceed the horizontal accuracy of the 1984 surveys (see McLean, 1990). The horizontal uncertainty of the 1952 surveys is unknown, but is estimated to be 5-10 metres maximum based on the juxtaposition of digitized soundings and contour lines with photogrammetrically mapped banklines. Channel features on each photo pair were digitized as a series of lines coded by feature type (e.g. banklines). A complete set of mapped channel features includes outer channel banks, islands, exposed channel bars, partly submerged channel bars and backchannels. All digitized data were -41 -  subsequently imported into the GIS for further processing, coding and analysis. The initial step involved converting all mapped lines to polygons (enclosed features having the same start and end coordinates) sharing contiguous boundaries. Since mapping did not extend beyond the main outer channel banks, complete channel and floodplain maps were produced by importing the location of bedrock outcrops, Harrison Bay / River and the outer extent of the floodplain. These data were provided by the Provincial Ministry of Environment from proprietary digital files of lower Fraser River gravel reach (data layers for railway lines, flood dykes and bank hardening were also supplied but are not used for producing floodplain maps). Once complete maps of the channel and floodplain were finalized, all polygons were classified by type of morphologic feature and by channel reach, and sequential maps were superimposed. The quantification and interpretation of historic channel characteristics based upon these maps is presented in the following chapter. 3.3 Bathymetric surveys Repeat topographic surveys of the channel bed and banks are used to measure the spatial distribution of bed material erosion and deposition volumes within common survey boundaries. Individual surveys are analysed to produce topographic models of the channel bed and banks and consecutive models are differenced to determine elevation increases (aggradation) or decreases (degradation) over time. On large, wandering gravel-bed rivers, channel deformation occurs over years to decades with most erosion and deposition occuring in distinct zones. Therefore, bathymetric surveys can be used to adequately describe these patterns provided the temporal and spatial resolution of individual surveys is sufficiently dense. Concerns related to estimating bed material transfer from observations of morphologic change have been reviewed in Chapter 1, while the relation between bed elevation changes and sediment transport is formally established within the context of a sediment budget framework in Chapter 5. The first complete survey of the river was undertaken in 1898. The survey extended upstream to Hope and included lower reaches of all major tributaries. The final dataset consisted of triangulations, topography, soundings and calculations of discharge, water-surface slope, and flood contours for the 1894 flood (Public Works Canada, 1949). Unfortunately, the original records were stored in New Westminster and subsequently destroyed in a major fire that occured later that same year. Other early records of development in the navigable sections of Fraser River were also lost in the fire (Public Works Canada, 1949). The gravel reach was not resurveyed until 1952. Subsequent bathymetric surveys have been completed in 1984, 1991, 1999 and 2003. -42-  Figure 3-1. Spatial extent of historic bathymetric surveys (black dots).  -43-  Various sounding transects have also been taken by government and private consultants for sitespecific projects within the gravel reach over the past several decades (including the full channel width between Hopyard Hill and Agassiz in 1962), but the spatial extent of these data are not sufficient for repeat topographic modeling at the reach scale. Full hydrographic surveys of all study reaches (Mission to Laidlaw) are available for only the years 1952 and 1999. The additional datasets are spatially constrained by comparison — the 1984 survey extends from Mission to Agassiz, the 1991 survey terminates at the confluence with Harrison River, and the 2003 survey ends at lower Herrling Island (Figure 3-1). The 1952 and 1984 data are re-analysed to contrast the methodological approach adopted in this study with results given in McLean (1990) using the same base data. The 1991 data are not used to calculate bed material storage changes, but are incorporated to determine an optimal method of spatial interpolation for bathymetric surveys collected along transect lines (see Section 2.4). The 1999 survey is used to measure the spatial and temporal variability of storage and transport zones within the river more reliably than was possible using a single intersurvey period (cf. McLean, 1990) and allows sediment transport estimates to be extended upstream of Agassiz. The 2003 survey was completed during a smaller than average freshet and no data were collected along elevated bar surfaces, or along major sidechannels where access was limited. However, more recent (March, 2004) acquisition of spatially dense LiDAR data allows these 'data gaps' to be filled in and an accurate 2003 surface constructed. Comparison of 1999 and 2003 modelled surfaces is used to construct a sediment budget in which sites of sediment erosion and deposition are apt to remain distinct. 3.3.1 C o m p i l a t i o n of available bathymetric data  The 1952 survey of the gravel reach was undertaken by Public Works Canada and extended from Barnston Island to Yale. River bathymetry was surveyed by fathometer, while exposed (above waterline) surfaces were mapped photogrammetrically (McLean and Church, 1999). Soundings were collected along cross-sections with an average spacing of 140 metres between Mission and Laidlaw. The final dataset was compiled as a series of 1:4800 charts with bed depths given to the nearest 0.1 foot (3 cm) and contours presented at 5 foot (1.52 m) intervals. The data were recovered to digital form by manually digitizing all charts covering the study reaches (approximately 12,000 individual sounding points). In addition, all contour lines from bar and island surfaces within the active channel zone and adjacent floodplain margins were digitized to provide nearly complete spatial coverage of the channel. The data were imported into Arc/Info GIS -44-  (Version 8.3, ©ESRI, 2002) where elevations were converted to metres and positions were adjusted from latitude/longitude to UTM coordinates (NAD83 datum). As depths were given to chart datum (geodetic mean sea level) no additional adjustment for water level was required. A TIN model (surface of irregular triangles) was created from the data to visually examine the output for coding errors (i.e. points where no elevation value was assigned, or the assigned value was clearly wrong based on neighbouring values) and obvious discrepancies were corrected by hand. The accuracy of contour line elevations was similarly verified, but small positional errors could not be detected (hence corrected) because the dense spacing and thickness of the drawn lines on the charts hindered interpretation in some regions. The TIN model also revealed that contour spacing was commonly too sparse on most of the larger exposed bar surfaces, resulting in an interpolated surface with little or no topographic variability (i.e. many of the modeled bars appeared unrealistically 'flat'). Similarly, contour density was found to be inadequate within several prominant sloughs and most secondary channels, resulting in a modeled surface in which these features appeared either unrealistic in shape, or the model could not distinguish the channels from surrounding topography. To increase the accuracy with which the 1952 surface could be interpolated, the original survey was amended with photogrammetrically derived elevation data. Large-scale (1:10,000) aerial photographs of the channel from 1949 and small-scale (1:40,000) photographs from 1952 were used to digitize a coarse array of points in regions of insufficient spatial density. Additional points were recorded separately on adjacent island and floodplain surfaces such that relative elevation differences could be calculated. This step was necessary because of the difficulty in establishing absolute vertical accuracy of stereo-model orientation (because the vertical precision of the TRIM points is commonly more than 2 metres) hence the accuracy of recorded elevations. The data were imported into the GIS and relative elevations were calculated locally by comparing surveyed with digitized elevations on the islands and floodplain. Measured differences of at least 0.5 metres were subsequently used to adjust the recorded elevations of nearby points digitized on bars and secondary channels. A total of 1037 photogrammetrically-derived points were added to the original dataset. The second major survey was completed in 1984 using a combination of automated hydrographic survey, conventional cross-sectional survey, and terrestrial ground mapping (McLean, 1990). The hydrographic data were collected along transects (120 metre average spacing) over 38.5 km of main channel using Environment Canada's HYDAC survey system -45-  (Zrymiak, 1984). A copy of the reduced dataset was obtained on CD-ROM from Environment Canada, Ottawa (UTM coordinates, NAD83 horizontal datum, geodetic vertical datum). The conventional cross-section sounding surveys (average spacing 250 metres) were completed by a UBC crew along 18 km of sidechannels between Chilliwack Mountain and Agassiz and 10 km of main channel between lower Sumas Mountain and Yaalstrick Island. The UBC crew also surveyed exposed bars, islands and channel banks using a terrestrial-based transit level. The UBC data were supplied [from D. McLean] in the same coordinate system as the HYDAC data. The complete 1984 channel survey consisted of 66,600 individual elevation points. The data were imported to the GIS and overlaid on the the 1983 base map to check positional accuracy (i.e. proper alignment with respect to islands and the outer banks). Several sounding lines near Calamity bar were oriented using a different (NAD27) datum, hence were adjusted. The 1984 survey was found to be more spatially limited than either the 1952 or 1999 surveys between Mission and Agassiz, as relatively low water conditions during the 1984 freshet precluded access to several areas within the channel. This effectively limits the precision with which bed deformation can be detected, as no surface can be interpolated in regions without sufficient topographic data (hence no surface comparison can be made). Rather than exclude these regions (cf. McLean, 1990; McLean and Church, 1999), it was decided to amend the original survey with photogrammetrically derived points, as completed for the 1952 survey. Low-water aerial photographs from 1979 were selected as the large scale (1:10,000) allowed bar and island topography to be clearly discernable and photos were available for the entire gravel reach. Ideally, the selected photos would be nearer to the date of the survey because intervening erosion or accretion of sediment on mapped bar surfaces introduces an unknown bias in elevation measurements. McLean (1990) shows photo mosaics of low-water aerial photographs for 1982 and 1984 that covered much of the gravel reach, but flightline records for those dates could not be recovered. However, a visual comparison of the 1979 photographs with these mosaics shows modest morphologic change, coincident with the lack of major floods during this period. It is therefore assumed that island and floodplain surface elevations would have remained constant, whilst bar surface elevations likely changed by no more than several tens of centimetres on average, based on long-term mean bed-level changes downstream of Agassiz (Church etal., 2001). Elevation changes of this magnitude are equivalent to the precision with which elevations can be measured from aerial photographs using the stereoplotter.  -46-  Following initial model orientation of each stereo pair, a point fde of coordinates from the 1999 laser altimetry survey was imported into the stereoplotter digitizing program, and positions were displayed on-screen. A replacement set of control points was recorded at these locations (all of which were chosen on stable island / floodplain surfaces), and absolute orientation parameters were re-calculated based on the more accurate altimetry elevations (-15 cm, compared with the 15 metre vertical error associated with the TRIM data). This procedure also obviated further vertical datum adjustments in the GIS, as required for the 1952 data. A text fde containing coordinates for all 1984 bathymetry data and altimetry data was then imported and displayed to provide simple visual reference on the photographs for identifying locations where elevation data were lacking (and similarly avoided the collection of redundant data within surveyed regions). A total of 13,038 points were digitized from exposed bars and islands and added to the 1984 sounding file. The 1999 bathymetric survey data were delivered by Public Works and Government Services, Canada in October, 1999. The complete data set consisted of 301,500 sounding points in ASCII format. Although survey lines were spaced roughly 200 metres apart, elevations were obtained at sub-metre spacing along each line, thus accounting for the large total number of points. To reduce the density of the sounding points, a subset of the original data was created using a Fortran program to retain points at 20 metre spacing between the start and end points along each individual survey line (additional data were deemed redundant for the intended bed surface model construction). The edited dataset consists of roughly 12,000 points, comparable with the 1952 survey. To complete the 1999 dataset, elevation points from exposed bars and floodplain surfaces were incorporated from a laser profding topographic survey completed in March, 1999. The data were collected using a single-track pulsed laser mounted on a helicopter. Measurements were continuously recorded along 200 metre flightlines across the floodplain (roughly between the flood dykes) along the length of the gravel reach. The raw data were processed by the contractor to calculate ground positions and eliminate outliers, recorded where a signal is returned from buildings or vegetation rather than the true ground surface. The completed dataset consisted of roughly 60,000 points at a spacing of 1 to 10 metres. Initial modeling attempts revealed the presence of remaining outiiers (usually isolated single points) that had not been detected. A Fortran program was written to search for adjacent points at least 3 metres higher than their immediate 2 neighbours along each transect line, an experimentally derived value that successfully 'flagged' the suspected individual outliers while preserving real topographic variation. The fdtered dataset was -47-  then overlaid with the 1999 channel map, so each topographic point could be assigned a channel map classification code (i.e. water, gravel bar, floodplain). Within the GIS, points falling on either the river surface, backchannels, or bar edges were deleted since these points obviously do not represent the actual bed of the river. Following this additional editing, the laser altimetry dataset was reduced to 55,700 points. The 2003 bathymetric survey data were collected by Public Works and Government Services, Canada during the summer freshet (late June to early July) and delivered in October, 2003. The complete data set consisted of 332,000 sounding points in ASCII format from Mission to lower Herrling Island. The data were collected following the same procedures as in 1999. including near replication of sounding line positions. However, additional lines were collected between the 1999 transects, giving an average line spacing of 100 metres. The modest increase in the total number of recorded points despite the increase in survey density partly reflects the absence of data collection from lower Herrling Island upstream to Laidlaw. More importantly, low water conditions during the data collection period precluded surveying on major bars and bar complexes which were either exposed, or covered by water that was too shallow for navigation. Consequently, the 2003 survey was spatially limited in extent, similarto the 1984 survey. The 2003 dataset was not completed until June, 2004 when the availability of additional LiDAR data fortuitously became available as part of an independent hydraulic modeling study to update the design flood profile. LiDAR data were collected in April, 2004 across the channel zone extending to the limits of the floodplain (i.e. beyond the flood dykes) from Chilliwack Mountain to lower Herrling Island, and at Matsqui Prairie (the left bank floodplain near the downtream end of the study reach). No data were acquired on theriveror floodplain between these locations. Data were collected using a swath-beam system that provided near continuous ground coverage, and were post-processed at 1-metre ground resolution. The complete dataset consisted of approximately 250 million individual 'bare-earth' ground returns. A copy of this survey was obtained in ASCII format (78 separate files) and further edited in the GIS. Given the physical size of the datafiles, it was first necessary to thin the data to a more manageable survey spacing due to computer memory limitations. A TIN model was created for each of the ASCII files and was converted to a regular grid array using a user-defined spacing (chosen as 5 metres). A second subroutine was used to convert each grid array to points, where each point location is equal to the centroid of a grid cell. This file format is used to store all of the historic bathymetric and terrestrial survey data. Since the LiDAR survey recorded ground returns -48-  on the water surface, these points required deletion or else would conflict with the sounding elevations. The water surface was delimited from 2004 aerial photography (coincidentally taken near the date of the LiDAR survey) and depicts channel conditions at lower flow conditions than during the sounding survey. All points that intersected the mapped water-surface were removed from the point fdes. The remaining data gap between Matsqui Prairie and Chilliwack Mountain was closed by amending the 2004 LiDAR survey with points from the 1999 LiDAR survey, as it is assumed that the vertical elevation of the floodplain and islands did not change over this period. Despite these efforts, there remained a paucity of elevation data at lower Yaalstrick bar. However, a comparison of 1999 elevations on the exposed bar surface with several randomly collected photogrammetrically-derived points revealed little or no apparent change in bar elevation. Therefore, 1999 survey data were assumed to adequately represent the true channel topography at this location, though of course no net change will be computed at this location between 1999 and 2003. The final dataset consisted of all sounding points, 4.73 million points from the 2004 LiDAR survey and 7700 points from the 1999 LiDAR survey. The bathymetry data were not thinned, because the 2003 surface is compared only to the 1999 surface, which was also remodeled using the complete set of available data points. Procedures for converting these different sources of elevation data into topographic surfaces that can be directly compared are presented in the remaining sections of this chapter.  3.4 Construction of 3-D surfaces 3.4.1 Introduction  The distribution of elevation behaves reasonably predictably because elevation is not influenced by complex, or unknown, interactions with other variables. Consequently, models that represent topography can usually be constructed with greater fidelity than can models for most other spatial phenomena. Nevertheless, several factors may complicate the accurate portrayal of the earth's surface, including biased or inadequate distributions of data, elevations available only along contour lines, sharp breaks in slope, local low (sink) and high (peak) elevations and the presence of vertical features such as cliffs (ESRI, 1991; Burrough and McDonnell, 1998). Hutchinson (1996) adds that spatially varying surface complexity and anistropic behaviour also preclude a simple analytic approach to interpolation.  -49-  Although different interpolation schemes tend to produce similar results given sufficiently dense data, they can produce widely varying surface representations given the presence of one or more complicating factors (Burrough and McDonnell, 1998). For example, a recent study (EGIS, 1997) estimated erosion and deposition volumes for channel islands along a 17 km reach of Brahmaputra River and found differences ranging from -13 million m (net erosion) to +5 million 3  m (net deposition) using two different mathematical functions to interpolate sparse elevation data. 3  Choosing the most appropriate technique for the available data is therefore a key consideration for interpolating sparse data (Isaaks and Srivastava, 1989; Englund, 1990; Bishop and McBratney, 2002) and reducing unwanted artifacts (noise). Despite this, the quantitative comparison of different modeling techniques is rarely considered and does not appear to have been undertaken in the geomorphic literature. There in no known universally accepted (correct) technique by which to topographically model the bed and banks of a river channel. Indeed, even experts will opt to apply different interpolation approaches to the same dataset based on partly subjective judgements (Englund, 1990). Fluvial researchers have presented diverse approaches including Laplacian and spline interpolation (McLean, 1990), inverse distance weighting (Eaton and Lapointe, 2001) and kriging (Fuller et al., 2003), all of which may be appropriate when the collected data are semi- or fullydistributed. However, most recent studies have suggested that TIN models (specifically Delaunay triangulation) present the most appropriate technique for representing complex river bed topography (Lane et al, 1994; Milne and Sear, 1997; Brasington et al, 2000). Triangulated irregular network models (TINs) are often used to represent a continuous surface because the method preserves known spot elevations and the density of the triangulated mesh can be adjusted during data collection to the complexity of the surface (Burrough and McDonnell, 1998). As well, interpolation can be restricted across breaks of slope (i.e. river banks) to limit topographic distortion (Lane et al, 1994). Delaunay triangulation consistently appears to provide a realistic representation of a channel surface where the data are well distributed (as in the examples cited above). However, the bathymetry data available for lower Fraser River (as is common with available conventional cross-sections on many systems) are heavily biased in the cross-stream direction. Although studies have recognized the weakness of cross-sections for representing distributed channel topography (Lane et al., 1994; Westaway et al., 2000) they have not examined the potential weakness of the TIN model itself for interpolating between cross-sections. Nevertheless, cross-section data continue to be collected because they are convenient for many -50-  other applications (e.g. ID hydraulic modeling) and because cross-section lines simplify data collection (i.e. river navigation). One of the major difficulties in choosing the most appropriate interpolation model is that there rarely is an independently acquired true surface with which direct comparisons can be made. As a result, qualitative judgements are commonly applied to reject or accept particular modeling strategies based on the perceived visual representativeness of the output. Because minimum curvature techniques (eg. splines) tend to produce smoother, more realistic appearing surfaces, researchers may be biased in choosing these models without subjecting them to error analysis or verification. A more rigorous approach for model evaluation involves cross-validation, where an estimate is made for a point at a known location by excluding that point and re-estimating its value from remaining data (Isaaks and Srivastava, 1989; Carr, 1995). The difference between actual and estimated values can subsequently be evaluated for the entire datatset to calibrate estimation accuracy (Carr, 1995). Deutsch and Journel (1992) add that different modeling strategies can be compared by examining the distribution of differences (errors) where the best approach produces a symmetrical error distribution centered on zero, minimum variance and no spatial correlation (the errors are independent). Violation of these properties can be used as a guide to adjusting model parameters (eg. sample weight, search radius) until the distribution improves. Cross-validation is appropriate for evaluating different interpolation schemes, including exact interpolation techniques (i.e. TIN, kriging) where estimates at known data locations are the same as actual data values. However, this procedure is subject to a number of limitations, especially the inability to evaluate model performance where actual data are sparse or absent (Isaaks and Srivastava, 1989). An alternative technique known as jackknifing can be used to overcome this common problem. The term applies to sampling without replacement (Deutsch and Journel, 1992) and refers to separating a dataset into estimation and validation components. Since the 1991 bathymetry data are spatially more dense than in the other years for which soundings were collected, this dataset can best be used for jackknifing, wherein the estimation data can be sampled from sounding lines spaced similarly to the other surveys and the remaining points can represent the validation data. A modified jackknifing strategy is employed herein, whereby the estimation data are used to produce an estimated surface to compare against a true 'reference' surface. The reference surface is derived from both estimation and validation datasets because removal of points between cross-sections may affect the quality of interpolation at retained data locations. As well, the goal -51 -  of selecting the most appropriate model is to reproduce the entire channel bed as accurately as possible using the thinned data, rather than simply minimizing estimation errors in regions lacking data. An additional benefit of this approach is that an estimate of model error over the entire channel bed can be provided. The 'best' modeling strategy, therefore, is that which reproduces a continuous surface whose spatial characteristics are most similar to that of the reference surface and which minimizes the absolute difference between the two surfaces. A true reference surface is derived along a morphologically complex 7-km long test reach extending from Queen's Bar to Harrison Bar using the 1991 survey data. The development of this surface, and a comparison of alternate modeling approaches available in the GIS, is presented below. 3.4.2 Reference surface data  Lower Fraser River was surveyed by echo-sounder in 1991 by the Canadian Hydrographic Service (Institute of Ocean Sciences). Soundings were taken along cross-sections spaced 50 metres apart from New Westminster to Foster Bar, and included lower Harrison River. Both instream and off-channel planimetric features (ie. gravel bars, islands, banklines, backchannels) were mapped photogrammetrically from 1:30,000 aerial photographs taken in 1990, but elevations were not recorded. The data were published as a series of 1:10,000 nautical charts, with soundings reduced to water datum (related to geodetic datum by the water-surface slope). A digital copy of the sounding data (15,650 points) was obtained within the study area, with all depths adjusted to geodetic datum, and positions converted to the UTM (NAD83) coordinate system. As the sounding lines were not always perpendicular to the shoreline (the survey also contains longitudinal, oblique and curvilinear cross-sections), the complete survey represents the channel bed as a welldistributed set of sample points, each spaced roughly 30-50 metres apart. This regular distribution is ideal for interpolating spatial phenomena (cf. Burrough and McDonnell, 1998) and is sufficiently dense for capturing the major variability of bed topography within the wetted channel. Since the purpose of the 1991 survey was to provide updated hydrographic charts within navigable sections of channel, no attempt was made to survey smaller, secondary channels or exposed (above waterline) regions. Upstream of Sumas Mountain (where the channel bifurcates around islands and emergent gravel bars) this results in a significant paucity of data with which to produce a complete topographic surface of the entire channel zone (see Figure 3-2 A). This effectively limits comparison of different interpolation techniques to the morphologically simple channel between Mission and Sumas Mountain. In an attempt to reduce the effect of this restriction, a GIS-based procedure was developed to interpolate elevations in regions of no data, -52-  Figure 3-2. (A) The complete testing data for development of the reference grid. A thinned sub-set of data (B) is used to to compare different interpolation schemes (see text).  -53-  the boundaries of which were identified in the GIS (there are no low-water photographs near this date from which elevation could be directly measured). The interpolation is based on a difference grid (25 m resolution) calculated from 1984 and 1999 surveys . Since the river is so large and 1  change occurs relatively slowly, it is assumed that there is little compensating scour or fill within this 15 year intersurvey period (i.e. if there was aggradation at a given location between 1984 and 1999, it is likely that there was aggradation in both 1984-91 and 1991-99 intersurvey periods). Absolute errors do not really matter in any case, since this is merely a test dataset for the interpolation methods. This procedure ensures that modeled surfaces appear morphologically (visibly) realistic, which aids in evaluation of interpolator comparison. The change in elevation between these dates was adjusted by 0.5 (eg. half the total bed level change from 1984 to 1999 was assumed to have taken place by 1991), then added to the 1984 surface to produce an estimate of the 1991 elevation. These grid-cell elevations were then converted to a set of equally-spaced points, where each point corresponds to the centre of a grid cell. To make the data spatially comparable to the 1991 bathymetry data, the grid was sub-sampled at a spacing of 50 metres. These equally spaced points were then appended to the bathymetry data. To complete the point dataset, the 1999 laser altimetry transects were also incorporated into the test dataset, but first were parsed to eliminate points falling on water or bar surfaces and points falling outside the test reach boundary. The complete point dataset consists of 5353 sample points; 2528 from soundings, 1007 from interpolation between 1984 and 1999, and 1818 from laser altimetry. The point records were supplemented with available contour data from the 1952 survey where the contours were coincident with island and floodplain surfaces stable within the 19521991 intersurvey period, since contours are known to improve the continuity of a modeled surface. Additional contours were manually interpolated at 5 metre intervals along the channel bed for completeness. A series of breaklines (linear features that define abrupt breaks in surface continuity) mapped as island and floodplain edges from 1991 aerial photography completed the input parameters for model development and comparison. The complete spatial distribution of the input data is shown in Figure 3-2 A.  1. Construction of 3-dimensional surfaces, differencing of surfaces, and the estimation of bed material transport rates have been previously reported using the 1952,1984 and 1999 survey data in Church et al. (2001). A copy of the section of the manuscript detailing computational procedures is given in Appendix B of the thesis.  -54-  3 . 4 . 3 Reference surface derivation  Although there are a variety of local and global models within the GIS from which interpolated surfaces can be generated, the reference surface was initially constructed using only the TIN and Topogrid modules in Arc/Info since these models have been specifically designed for and tested on elevation data and both models can incorporate the available point and contour data. TIN algorithms are particularly well suited for modeling regular arrays of elevation data provided the point spacing is sufficient to characterize topography in regions of complex terrain (Burrough and McDonnell, 1998). In addition, the inclusion of contours and breaklines reduces the number of invalid triangles (flat triangles across breaks of slope or topographically unrealistic, long, thin triangles). The Topogrid model is based on a "discretized thin plate spline," (i.e., one in which an exact spline surface is replaced by a locally smoothed average). The model replaces the minimum curvature interpolation criterion inherent in the spline model with a minimum profile method that better matches landform processes and preserves drainage structure (Hutchinson, 1996). Effectively, this allows the model to follow abrupt changes in terrain, such as along channel banks, which minimum curvature techniques can not (John E. Hughes Clarke, pers comm, 2000). Although elevation contours can be assimilated into the interpolation process, natural breaks in slope can not be explicitly defined except where coincident with a contour or sufficient point data. The accuracy of both the TIN and Topogrid interpolations was first evaluated by superimposing the original data points on each interpolated surface and comparing the two datasets using the RMSE statistic. For this comparison, the TIN model was first converted to a 10 metre regular raster array (coincident with the resolution of the Topogrid model). This step is neccessary because points can not be directly superimposed onto the TIN data structure in the GIS (similarly, consecutive TINs can not be directly subtracted). The best results (RMSE ± 0.345 m) were produced using a quintic (5th degree polynomial) interpolator to convert the TIN to a raster. For comparison, the best results using the Topogrid model produced a calculated RMSE of ± 0.547 m. The TIN models produce smaller errors, in general, because the interpolation method is exact - the location and value of the original data points is preserved in the model. However, deviations (RMSE) occur during conversion to a regular grid since each grid cell may spatially overlap, hence average, more than one original data value. As the 10 metre grid cell size is small relative to the size of triangle facets (which typically connect 50 metre sounding points) this averaging does not generally lead to large differences between original and modeled elevations.  -55-  A scatterplot of the data reveals that only a small number of points display significant deviation between actual and interpolated values in either model (differences are as large as 9 metres near vertical banks - see Figure 3-3 A, B). However, because the RMSE statistic is highly sensitive to outliers, these points account for a large proportion of the total sum of squared deviations. Analysis of these residuals demonstrates that they are slightly skewed for both interpolated grids (more frequently underestimated) while variograms show a pure nugget effect for the TIN, an indication that errors are randomly distributed over the study area (Figure 3-3 C). Errors are similarly random for the Topogrid model, and both variograms exhibits a dip (hole effect) at a distance of 5 km (cf. Figure 3-3 D), indicating noise or a periodic trend in the data (Isaacs and Srivistava, 1989; Burrough and McDonnell, 1998). This trend likely corresponds to channel gradient. However, a simple plot of the residuals in the GIS clearly showed that the largest errors were found along island and bank perimeters, not within the main channel (Figure 3-4). The variogram model likely fails to detect this pattern because island and bank features are welldistributed throughout the study area. A visual plot of the actual interpolated surfaces illustrates poor interpolation of the TIN model between channel banks and the end of sounding lines and along channel margins between sounding lines, despite the density of available data (Figure 3-5). In comparison, surface continuity in these regions appears more realistic using the Topogrid models (Figure 3-5) even though, statistically, they perform worse overall (larger RMSE). Because statistical measures do not necessarily capture these types of errors (especially when the interpolation is accurate near known data values in the TIN), visual checks are an essential part of the interpolation process (Isaaks and Srivistava, 1989). This problem probably occurs because breaklines were not explicitly defined with elevation values and contours encompass only a fraction of the total bankline length (cf. ESRI, 1991). Since gradient changes are greatest between the channel bed and floodplain, interpolation errors along 400 km of channel margins (i.e. over the entire gravel reach) can cumulatively result in significant total errors that could dominate erosion and deposition measurements, especially given the larger transect spacing of the full surveys. Reducing interpolation error along the bank edges is therefore essential to minimize bias in the entire sediment budget. To reduce these errors, a method was developed to estimate the elevation along channel margins in an attempt to provide a more reasonable representation of actual topographic variation. The laser altimetry data were isolated and used to create a surface grid model that included island and floodplain surfaces only. Because these surfaces are reasonably flat, the survey line spacing -56-  Figure 3-3. Scatterplots showing difference between actual and interpolated values for the reference dataset using A. TIN and B. Topogrid models, and Semivariograms of C. complete reference dataset and residuals from TIN model. Plots for the Topogrid model (not shown) exhibit a similar pattern. Arrow in C indicates a dip in the variogram (see text for discussion).  Figure 3-4. Plot of elevation difference between actual and interpolated values for the reference dataset using TIN and Topogrid models.  does not bias modeled surface continuity. The banklines mapped from airphotos were subsequently overlaid on the floodplain models to produce a set of banktop contour lines (i.e. breaklines with elevation values) that could be incorporated into both models. A second set of parallel contours (with 10 metre offset) was derived for the bottom of the banks using only sounding points. Several TIN and Topogrid models were re-estimated using the test dataset and the bankline contours and altering available input parameters (Topogrid) or trying different combinations of banklines (for example, by eliminating bank bottom contours in TIN). In general, the addition of bankline contours results in visually smoother, more continuous surfaces with a corresponding reduction in RMSE to 0.468 m for the Topogrid model using recommended default input parameters for tolerances that influence data smoothing. By comparison, the TIN models generally showed improved interpolation along the banks, but many invalid triangles remained and obvious interpolation errors were occasionally observed on channel bar and island surfaces. The RMSE actually increased to 0.420 m because the bankline contours frequently overlaid actual survey data, resulting in averaging at those coincident locations. Although the increased RMSE for the TIN model remained smaller than that of the Topogrid interpolator^ the Topogrid model is adopted as the reference surface because, visually, it appears to provide a slightly more faithful reproduction of an actual river channel (see Figure 3-6). 3 . 4 . 4 Model comparison In an attempt to test the stability of surface models with reduced information, different interpolation techniques were applied to a thinned subset of the idealized test dataset. The reference dataset was parsed by removing data along sounding transects until the remaining point spacing and density resembled those of the 1952,1984 and 1999 surveys (i.e., retained points formed a set of parallel sounding lines roughly 200 metres apart). A total of 2864 points (55%) were deleted however, since all altimetry points were retained, the points removed represent 81% of those within the channel bed (see Figure 3-2 B). A series of tests was then performed by subtracting models based on the thinned dataset (with and without bank contours) from the reference surface. This allows the suite of interpolation tools available in the GIS to be directly compared, and allows the benefits of the bank contouring method to be demonstrated. Since only the TIN and Topogrid modules can incorporate contour lines directly into the modeling routines, all contours were converted to a series of nodes for all other models to provide a comparable set of test data. The optimal modeling strategy is that which minimizes not only the net difference between the model surface and the reference surface, but also the magnitude of the apparent scour and fill, since -60-  Figure 3-6. Reference surfaces based on T I N (B) and Topogrid (C) models. The mapped channel ( A ) as interpreted from 1991 aerial photos is shown for comparison.  significant values of either represent interpolation error between survey lines (i.e. both values should be zero). Although there will be no net difference between model and reference surfaces if scour and fill volumes are equal, large values of either measure represent deviation from the reference surface and will produce commensorate changes in the magnitude of the RMSE. Beyond those models already presented, there exists a wide variety of tools for modeling point data within a GIS. The following presentation is, however, limited to those models commonly applied to interpolating a continuous surface representation of elevation, and includes inverse distance weighting (IDW), splines, and kriging models. While it is recognized that other GIS or modeling programs may have different tools (and may therefore produce very different results) the purpose of this analysis is to compare the models available in a widely-available GIS program. The use of any model is complicated by the variety of available options (i.e. for sampling weights, search distances, mathematical 'fit' functions) that necessitate some operator knowledge of spatial statistics and the distribution of the phenomenon under study. Each of the tools was used to create a 10-m resolution grid starting with recommended default parameters and this grid was subtracted from the reference surface to calculate apparent scour, fill, and net volumetric change. Default options were then modified within reasonable limits, and each 'modified' 10-m surface was also compared to the reference surface, a process that continued either until no further reduction in apparent scour or fill volumes was observed or until the net volumetric difference became increasingly large. The best results for each model (from a total of 30 trials) are presented in Table 3.2. The widely varying range of results reflects the appropriateness of the individual interpolators for modeling topographic data. Inverse distance weighting determines cell elevations using a weighted combination of sample points, where weight is a function of distance. As such, it is a local method of interpolation, where it is assumed that the value at unmeasured locations is similar to nearby measured values (Burroughs and McDonnell, 1998). Since the interpolated value for any cell is an average of surrounding values, the method is not able to preserve local maxima or minima, hence can not reproduce ridges and valleys unless these features have been sampled (ESRI, 1992). While the technique may reasonably model densely sampled surfaces, it clearly is not appropriate for complex river topography that is irregulary sampled. Splines are also a local interpolator, but attempt to fit a continuous surface that passes exactly through known data points  -62-  while minimizing surface curvature. While this approach results in the relatively small net change observed, smoothing between data points produces unreliable results in general. Since an underlying trend is the bed surface was previously identified (see Figure 3-4), universal kriging was tested using a linear drift, where drift is defined as a systematic change in elevation over a particular direction (ESRI, 1992). Kriging is a geostatistical interpolator which recognizes that the spatial variation of a continuous attribute may be too irregular to be modeled with a smoth mathematical function and alternatively may be better described by a stochastic surface (Burroughs and McDonnell, 1998). While the technique has been commonly used to create digital elevation models, abrupt local changes in elevation can create discontinuities in surface reproduction that vary with the size of the interpolation neighborhood (Meyer, 2004). In general, the technique is not appropriate for modeling complex river topography because anomalous maxima and minima produce a high degree of variability, which the technique tries to reproduce, resulting in high values of apparent scour and fill. Although these effects could be reduced by stratifying the data into more homogenous regions (ESRI, 1992), this was not done.  Table 3-2. Comparison of different interpolation schemes using thinned dataset with reference dataset Net change (m )  RMSE  1,651,979  -405,426  0.444  1,508,029  1,438,153  -69,876  0.901  Topogrid - B  3,425,290  3,977,451  +552,161  1.085  IDW  4,844,037  8,991,107  +4,147,070  2.030  Spline  8,245,762  7,986,886  -258,876  2.120  24,773,392  22,548,868  -2,224,523  3.522  Apparent Ero-  Apparent Deposi-  sion (m )  tion (m )  TIN  2,057,405  Topogrid - A  Method  3  Kriging  3  3  By comparison, the TIN and Topogrid models produce a much more realistic surface from the thinned data. In particular, these results demonstrate that the Topogrid model [A] replicates the reference surface much more accurately than any of the other available tools. Bishop and McBratney (2002) reached a similar conclusion based on a comparison of elevations modeled using the Topogrid tool with interpolation techniques available in different statistical programs. The 'best' Topogrid results in Table 3.2 incorporate the simulated contours along the bank top and bottom and include a 10 cm vertical error term (equivalent to the reported accuracy of the sounding data). The model was also run without the simulated contours [B] to demonstrate their impacts on -63-  model accuracy. These results show that cut and fill volume errors more than double, and the net volumetric difference increases by a factor of 8. Although the net volumetric difference is small relative to that produced by other interpolation schemes, the individual cut and fill volumes appear large. However, much of the difference occurs in locations where the thinning process resulted in large data gaps - the thinned data are actually less spatially dense than the full bathymetric surveys at several locations. Adding additional contours by hand along the bed would likely reduce the magnitude of cut and fill volumes, but not necessarily the net difference. A visual comparison of the 'thinned' model with the reference model is given in Figure 3-7 and shows that the interpolation appears reasonable given the significant reduction in data. The most obvious modeling discrepancies appear along the channel bed, where the continuity of depth contours is sometimes interupted between sounding lines, and along the transition from the wetted low-flow channel to shallow bars where the modeled surface appears slightly wavy. These differences are more clearly illustrated in a histogram (Figure 3-8). This plot shows a bimodal distribution of elevations, corresponding to the channel bed and bar and island deposits. The jagged interpolation along shallow bar edges is seen in the 0 and 1-metre elevation bins. In order to test whether the thinned and reference grids are statistically similar, the calculated histograms were subjected to the Kolmogorov-Smirnov goodness-of-fit test for continuous data (Zar, 1984). The K-S test is a non-parametric (distribution-free) method for comparing two independent datasets measured at the ordinal or higher scales of measurement and there are no limitations of the size of the sample or frequencies in categories (Norcliffe, 1982). The test statistic, D, measures the maximum absolute difference between cumulative frequency distributions. The null hypothesis states that the two samples are drawn from populations with the same distribution. Critical values for the K-S test were determined from tables in Zar (1984, after Smirnov, 1939) for large n, as: r\ («, D  x  =  n) -  /-In (a/2) T  n  0.16693 —  The calculated test statistic, 0.0183, is nearly equivalent to the critical value at a = 0.05 (0.0185) but at a = 0.01, the critical value increases to 0.0222, so it can be concluded that the two populations are statistically similar. The net volume difference between scour and fill (70,000 m ) represents net interpolation 3  error and can be used as a measure of the accuracy in sediment budget calculations provided survey density is similar (Chapter 5). The net difference from the reference surface could actually be -64-  700  elevation bin (x-1 < value < x)  Figure 3-8. Histograms of 1991 elevation data including distributions for the original 5352 points, the reference dataset and the thinned dataset. The reference and test point datasets are extracted from the respective interpolated grid surfaces.  -66-  Chapter 4: Channel changes 4.1 Introduction In this chapter, adjustments to channel size, shape and position over the past century on lower Fraser River are investigated. A complete description of channel form must include planforrn geometry, channel gradient, cross-sectional shape and bed topography. By measuring these characteristics from available planimetric and topographic surveys for different dates, observed spatial and temporal changes can be related to changes in the flow and sediment regimes (the dominant governing factors) and to secondary anthropogenic influences. The magnitude of channel adjustments may provide an indication of the sensitivity of the river to environmental changes that have occured over the past century, and may therefore be used to predict future impacts on the development of channel morphology and complexity. Lateral migration is a natural, persistent process along alluvial channels that occurs in response to sediment exchange along the system. Channels may widen during periods of persistently above-average flood flows through bank erosion, vegetation removal (stripping) on depositional surfaces, or activation of secondary channels. However, the increase in sediment influx generally associated with higher flows, combined with inputs from bank erosion, will correspondingly cause aggradation on the channel bed and higher water levels if this material is not removed. Consequently, widening may be associated with a change in the shape of the wetted cross-section. In comparison, channel narrowing occurs during periods of persistently belowaverage flood flows from the establishment and maturation of vegetation on elevated bar deposits, but also because of bar and bank deposition. It is apt to be associated with comparatively modest exchanges of bed material or re-shaping of the cross-section. Channel width along wandering channels is typically variable over time and irregular along the channel length in response to variations in peak discharge, sediment influx and bank composition. Wide reaches are associated with laterally unstable 'sedimentation zones' whilst narrow reaches are associated with quasi-stable 'transport' zones (cf. Church and Jones, 1982). Although sedimentation zones are known to migrate over time, narrow reaches of channel may persist where the channel is bounded by erosion resistent materials. While lateral bank shifting and channel width are obviously related, these processes are examined separately since channel width will remain stable if bank erosion and deposition rates match, regardless of magnitude.  -67-  4.2 Active channel width Active channel width is defined as the average width of the water surface and exposed bars between established edges of perennial, terrestrial vegetation. It is calculated in the GIS by summing the total area of digitized water and gravel bar polygons, then dividing by channel thalweg length (measured for each date of mapping to account for thalweg migration). Active channel width represents that part of the channel which experiences active exchange of mobile sediments. Changes in width are related to long-term rates of sediment exchange where width generally increases in response to a greater sediment load (Schumm, 1977, Lyons and Beschta, 1983; Kellerhals and Church, 1989), although the relation is poorly understood (Warburton, 1996). Conversely, narrowing may be indicative of decreased sediment supply, although other factors must be considered (Friedman et al., 1996). Since the active width of the channel further represents the area of potential habitat substrate for fishes and aquatic insects, observed changes in width over time can provide a useful, albeit broad, indicator of changes to potential ecosystem productivity. Historic trends in active channel width between Laidlaw and Mission are shown in Figure 4-1. Trends in active channel width broadly parallel those of flows, but exhibit less variability. This is partly attributed to inertia in the system since peak flows can be highly disparate from year to year, but changes in channel width occur more gradually on large rivers. For example, many years of below average peak flows are neccessary to allow vegetation to establish and become sufficently dense to be classified as island or floodplain surfaces (i.e. there is a delayed response). Conversely, increases in width can occur more rapidly if channel alignment is conducive to extensive erosion of island or floodplain deposits (so the two trends do not neccessarily balance). This pattern is also influenced by the timing and dates of the photo mapping. The area classified as vegetated may be underestimated during leaf-free periods, while the analysis of only eight mapping dates masks any variability in width trends that might occur over periods shorter than an individual mapping interval, and further, does not always correspond with timing of dates associated with major trends in the flow. Nevertheless, the figure reveals two significant patterns in the long-term development of the river. Overall, the active channel zone has narrowed 22% from 1912 to 1999 betwen Mission and Hope even though flood flows have declined only 5% over the same timeframe (Section 2.2). Since width is generally proportional to the square root of discharge for a wide range of channels (Knighton, 1998), the expected decrease in width can be estimated as (0.95) , which is only about 3%. While this relation is approximate only (being subject to 05  variations in boundary conditions, the magnitude of sediment transport and the dominant channel-68-  Figure 4-1. Variation of active channel zone width for the entire gravel reach over time shown with long-term trends in annual maximum daily discharge. Major recent floods (1948,1972) are indicated. The flow trend is calculated as the cumulative departure from the mean annual flood.  -69-  forming discharge) the relation does indicate a major loss of potential aquatic habitat. Further, the volume of water through the reach must be flowing at greater depth and / or velocity to compensate for the reduced capacity to maintain conveyance. As well, the rate of channel narrowing (i.e. 19121949; 1971-1999) is obviously greater than the rate of channel widening (1949-1971) implying a response to additional forcing unrelated to discharge. Variations in width changes along the gravel reach were examined by first partitioning the river into 13 individual sub-reaches (Figure 4-2). The sub-reaches delineate the major, wide sedimentation zones which are separated by narrower, stable transport zones. Additional reach breaks were added at channel constrictions (i.e. bedrock outcrops) in an attempt to delineate the river into roughly equal lengths. Changes in width for individual reaches (Figure 4-3) generally follow the trends observed over the entire channel, except that confined reaches (i.e. 1,2,13) have changed little over the period of record. These reaches are characterized by modest sand and gravel storage (largely transport zones), no channel islands, and are bounded by erosion resistant terraces. In contrast, alluvial sections of channel have experienced a slightly greater (24%) decrease in average active channel width compared to the entire channel, but there is considerable variability. Channel width actually increased 11% overall since 1912 in reach 5 (adjacent to Chilliwack Mountain) because of the gradual erosion of mid-channel [Yaalstrick] islands, though more recent narrowing is consistent with flow trends. Apart from the rather extraordinary loss of 1.25 km of mean width (54%) in reach 6 from 1912 to 1999, the greatest decrease in width since 1928 (the earliest verifiable mapping) is found along reach 7, where the channel narrowed by 560 m, or 39% by 1999. In fact, significant narrowing has occurred along the entire channel from Chilliwack Mountain to Herrling Island, corresponding to an average decline of roughly one-quarter the active channel zone since 1912, and 19% since 1928. The observation that increases and decreases in width over time are reasonably uniform amongst the reaches qualitatively suggests that a single process is responsible (cf. Lyons and Beschta, 1983). Significant changes in patterns of sedimentation, for example, would result in greater discontinuity in width trends along a wandering river, since aggradation (or degradation) tends to occur in distinct zones rather than being evenly distributed (sedimentation patterns are discussed in Chapter 5). Further, the mapped area of the water surface for 1949 and 1999 (when flows are similar) is nearly identical, a circumstance that would not occur if width changes were mainly influenced by sedimentation. Instead, most of the loss is directly attributed to the isolation of major sidechannels around island and floodplain deposits through a combination of riprap, flood -70-  island / floodplain gravel bar water surface  Laidlaw  sub-reach boundary  Figure 4-2. Site map showing morphologic divisions of lower Fraser River gravel reach (based on 1999 photography). The river is divided into 13 sub-reaches which correspond to major deposition zones and narrower, stable transport zones.  5 km  2 <g o S a- 3 3 f>  1600  1600  1400  1400  3  "g  1200  1000  •a '%  1000  c c  800  c c  800  o > u  600  2  400  <  1200 to  p  >-• OQ n  S n 3  cn —  o  <*  •a  o  o 3. 3* t» O 3 CO  3  U  co  03  <  >  200  3  51  o  3-  1910  o <  9  600  ^  400 200  1920  1930  1940  1950  1960  1970  1980  1990  1910  2000  1920  1930  1940  3  1960  1970  1980  1990  2000  1970  1980  1990  2000  Year of mapping  Year of mapping  CH  1950  (2333 m)  CT>  3  a. c  1600  1600  1400  1400 £  1200  1200  C  CL *<  o sr en CO  CO  .e? OQ CH  1000  1000  c c  C  800  o u > u  600  3  600  13  400  400  200  200  & 3"  800  CS  0 1910  1920  1930  1940  1950  1960  Year of mapping  1970  1980  1990  2000  1910  1920  1930  1940  1950  1960  Year of mapping  gate construction, dyking and, eventually, siltation and vegetation growth. These activities are further responsible for the discrepancy between flow changes and width changes over time. 4.2.1 Island evolution Most of the bank hardening along lower Fraser River has been limited to outer channel banks. In comparison, channel islands remain largely unprotected with the main exception being populated Peters Island. Major adjustments in channel width, therefore, should largely occur as a result of growth or destruction of islands. Changes in island area over time are shown in Figure 44. The complete gravel reach including all alluvial reaches shows a dramatic decline in the total area of forested islands during the period 1912 to 1928, with a continued, though reduced decline to 1949 (note that 'missing' 1940 photos do not allow total island area to be computed for that date). This pattern runs counter to the observed decrease in active channel width during a period of mainly below average flood flows. This circumstance could occur if the rate of natural floodplain reconstruction was greater than the loss of island area through erosion, but this explanation seems improbable from a morphologic perspective. In fact, many islands identified on the 1912 map became increasingly isolated from the active channel zone by 1949 because of dyking, flood gate construction and subsequent infilling (siltation), hence re-classification of the islands to floodplain. In effect, artifically high rates of new floodplain generation were created. This process effectively increased the total area of off-channel floodplain, while the dominant cause of channel narrowing was the corresponding loss of side and backchannels around and through the islands. Island area continued to decline to 1971 (and probably to the late 1970s; see Figure 4.1) but has since increased, consistent with both long-term flow trends and declining active channel width. By 1999, the total surface area of the islands had actually recovered to the 1928 total area, and was greater than during any other year since the 1948 flood, after which the dyking program was expanded. This recovery is particularly noteworthy given that woody debris accumulations — a major nucleus for island formation (Hickin, 1984; Abbe and Montgomery, 1996) — have been mostly eliminated since the debris trap was established in 1979 (see Chapter 2). However, this statement must also be tempered by the possibility that woody debris has never played an important role in island formation on the river. That overall (total area) changes are minor reflects the fact that changes on the river occur slowly, with major sites of erosion and deposition persisting for years to decades because the magnitude of bed material transport is modest relative to the magnitude stored within the channel. -73-  2 3  10000  10000  1000  1000  o ere  c3  3*  gff  4^  8, 3  s 5 00  ere  JS" 3 a  3*  g3* 2  n  a ioo  -o £  1  3  510  10  W  5  3  a  /  Z B 85 rt  1  s»  4^  1  1910  1920  1930  _0  1940  1950  1960  1970  1980  1990  2000  1910  1920  1930  1940  1950  1960  1970  1980  1990  2000  1980  1990  2000  Year of mapping  Year of mapping  § 9 * Z 3 o  100  „ (11,617 ha)  o i-i c  10000  10000  Total (missing data) 1000  1000  a ioo  a ioo  rt 0  a O.  o 3 rt  1d  3  ed  13  •a  M3  10  10  1  12  l 1910  1920  1930  1940  1950  1960  Year of mapping  1970  1980  1990  2000  1910  1920  1930  1940  1950  1960  1970  Year of mapping  Although most individual reaches share a markedly similar pattern of island development with the complete channel, others reveal somewhat more dynamic change. Within a reach, the loss of island/floodplain deposits can occur during a single freshet, but re-establishment takes many years so a direct correlation with long-term flood flows is complicated by temporal discontinuities in mapping resolution. The pattern can be particularly erratic where total island area is modest (cf. reach 4) because avulsions or individual floods can influence a high proportion of the surface area (including complete removal of all islands). Reach 3 shows a continuous increase in vegetated area since 1928, though the total area is similarly small and confined to a single island (Strawberry) where mostly sands and fine to medium gravels are deposited. The main exception is found in reach 9, which was characterized by dramatic loss before 1983 as changes in channel alignment related to bar growth caused persistent erosion of large islands up- and downstream of AgassizRosedale Bridge. There are signs of recovery since as new islands have established and grown. 4.2.2 Bank erosion and deposition  Bank erosion and deposition are ubiquitous characteristics of alluvial channels, and can occur even along rivers which exhibit stable equilibrium form (Schumm, 1977). These processes are intimately linked with sediment exchange, with floodplain deposits acting as both a sink for material entrained upstream, and a potential source of erodible material (Meade, 1985; Kelsey et al., 1987). Floodplain surfaces are distinguished from gravel bars in that they are higher in relative elevation and usually well-vegetated, so that sand and silt sized materials (overbank sediments) become included in the process of sediment exchange. Although the rate of bank erosion depends upon the complex interaction of many variables (Hooke, 1980; Hickin and Nanson, 1986), most bank erosion results from high shear stresses, particularly along the outer apex of bends, that directly entrain material and undercut the banks causing gravitational slumping. The magnitude of bank erosion should therefore be related to the magnitude of flood flows, hence erosional (and depositional) rates should broadly conform to longterm trends in flood flows, though individual large floods may skew averaged short-term rates (Hooke, 1980). Bank erosion was calculated in the GIS by overlaying successive maps and summing the total area of island and floodplain polygons on the early date converted to water on the later date (the reverse situation was used to calculate deposition). These definitions are unaffected by stage since vegetated surfaces are exposed at all flows below bankfull discharge. Normalized bank changes for the gravel reach are given in Figure 4-5. Lateral bank deposition rates exceed erosion rates from 1928 to 1949 and since 1983, consistent with channel narrowing during 1  -75-  40,000  1928-49  1949-62  1962-71  1971-83  1983-91  1991-99  mapping period Figure 4-5. Average annual bank erosion and deposition rates for different mapping periods. A rate of 200,000 m per year is equivalent to 1.4 m per year (for each bank) averaged over the 72.5 km length of the gravel reach, but this rate effectively increases to 4.7 m per year if averaged over the total length of unprotected bank (see Figure 2-5). 2  -76-  these periods of mostly below-average floods. The 1928-1949 depositional bias would actually be greater if the potential erosive impact of the large 1948 flood is discounted, but the mapping sequence is not sufficient to isolate its effect on channel morphologic change. Conversely, bank erosion rates were much larger than deposition rates from 1949 to 1971 when channel width increased and flood flows were mainly above-average. However, the modest increase in channel width from 1949 to 1962 is surprising given the rate of bank erosion was much greater than the rate of bank deposition (19:1). This circumstance occured because the major expansion and upgrading of the dyking program narrowed the active channel zone along floodplain margins, countering the effect of erosion along channel islands and unprotected banks. The period 1971 to 1983 also displays an anomalous pattern given that channel width decreased despite the erosional bias. This can be explained by vegetation growth on gravel bars. Although there was a very large flood in 1972, this was a period of predominantly smaller floods conducive to vegetation establishment and expansion. Bank erosion and deposition rates along individual reaches can obviously be much larger or smaller than the averaged values presented above. The maximum recorded erosion rate was 8.8 metres per year in reach 8 (1949-1962) while maximum deposition averaged 15.2 metres per year in reach 4 from 1971 to 1983. Tables 4-1 and 4-2 summarize these data by reach for the period of mapping. Notwithstanding the complications involved in relating these changes to temporal flow trends, erosion and deposition rates are locally influenced by the amount of bank hardening (cf. inset table, Figure 2-5), the commercial removal of trees, and especially by local alignments in flow caused by bar migration and deposition. Consequently, a detailed examination of temporal trends is not provided. For comparison with the trends described for the planform mapping, volumetric calculations of bank change were calculated from surfaces of difference modeled from the bathymetric data. Volumes are equivalent to multiplying areal changes by the mean of differences in elevation over the same mapping locations. These calculations reveal that the total volume of sediments lost to bank erosion from 1952 to 1999 was about 10% greater than the total volume deposited, although the total areal gain is actually larger (1.65:1 based on a direct overlay of 1949 1. Bank erosion and deposition rates from 1912 to 1928 are not given because uncertainty in the absolute map orientation for 1912 may inflate rates of actual temporal change. However, since the map was correctly scaled, relative accuracy is preserved, allowing parameters such as channel width and island area to be calculated and presented. Nevertheless, bank shifting appeared to be far more vigorous early in the century compared to the present, with several sites along the gravel reach experiencing apparent lateral erosion up to 500 metres, with a maximum of 1000 metres recorded near Seabird Island (refer to discussion in Appendix C).  -77-  Table 4-1. Average annual bank erosionrates.All values expressed in metres per year. Values that exceed the reach mean are shown in bold. Reach  1928-49  1949-62  1962-71  1971-83  1983-91  1991-99  Average  1  0.88  1.15  0.11  0.35  0.32  0.57  0.64  2  1.12  2.61  0.51  0.64  0.36  0.89  1.12  3  1.97  1.54  2.52  1.57  2.17  2.50  1.97  4  3.10  1.97  5.24  3.49  3.54  2.71  3.24  5  4.46  7.58  4.00  3.00  2.85  1.71  4.23  6  3.09  2.20  1.93  1.74  2.41  1.48  2.29  7  3.50  6.89  3.24  3.76  1.97  0.76  3.65  8  3.59  8.78  4.72  4.76  4.69  3.70  5.02  9  0.95  3.49  4.01  2.85  1.54  0.66  2.16  10  4.68  7.61  2.50  5.20  639  2.21  4.94  11  5.68  1.94  1.35  2.21  0.70  0.46  2.71  12  2.09  3.92  1.74  2.09  1.46  0.90  2.17  13  0.25  0.96  0.00  039  0.39  0.00  0.36  Table 4-2. Average annual bank deposition rates. All values expressed in metres per year. Values that exceed the reach mean are shown in bold. Reach  1928-49  1949-62  1962-71  1971-83  1983-91  1991-99  Average  1  0.52  0.01  0.82  0.13  0.93  0.52  0.44  2  0.87  0.06  0.32  0.14  0.98  0.65  0.52  3  7.68  0.56  0.36  11.69  3.05  1.94  4.96  4  10.98  0.16  1.09  15.15  1.09  0.72  6.18  5  4.68  0.16  0.46  5.71  2.34  2.31  2.96  6  5.00  0.29  2.40  1.04  4.24  3.77  2.91  7  4.55  0.46  4.48  0.32  6.26  5.18  3.34  8  4.24  0.30  0.97  0.64  2.17  3.09  2.13  9  2.68  0.07  2.01  0.48  0.95  1.03  1.36  10  5.41  0.33  2.19  0.42  4.57  3.77  2.95  11  1.53  0.11  0.61  0.02  2.70  2.96  1.19  12  1.73  0.30  0.26  0.12  0.83  1.38  0.84  13  1.64  0.00  0.27  0.09  1.19  0.82  0.76  -78-  and 1999 channel maps). This circumstance occurs since eroded areas include established (mature) deposits that are higher in absolute elevation (i.e. they are thicker sediment bodies) relative to sites of more recent deposition because of vertical aggradation of fines over time. The positive relation between overbank thickness and vegetation age along lower Fraser River has also been demonstrated by Boniface (1985) and McLean (1990). If overbank sediments are removed from the calculations, the deposited volume of gravel and coarse sand is roughly double the eroded volume, a bias that is much greater upstream of Agassiz-Rosedale bridge (3.2 X) than downstream (1.8 X). These results confirm the growth of channel islands over the past half-century, but the matureriparianforests on old islands are increasingly being replaced by young vegetation on new islands, a pattern qualitatively confirmed by visually comparing pre-1950's with recent aerial photography. The transition to lower elevation, immature islands implies a change in the topographic characteristics of the channel. This possibility is investigated by examing changes in the longitudinal profile, cross-sectional area, and distributed water depths over the available period of record. 4.3 Longitudinal profiles The longitudinal profile of an alluvial channel represents one dimension of adjustment to the prevailing discharge and sediment regimes. The long profile is commonly represented as the elevation of the high water-surface at distances downstream. The water surface parallels the average bed surface at the reach scale, but also masks local variations in bed topography (Richards, 1982; Rice and Church, 2001). These local variations correspond to riffle and pool morphology in alluvial gravel channels. The shape of the bed profile changes over time in correspondence with patterns of bed material erosion and deposition, as the spatial distribution of riffles and pools changes. Persistent bed material aggradation or degradation will further influence the shape of the bed profiles by filling or scouring riffles and pools, and may alter the average slope of the longitudinal profile. Longitudinal profiles of the channel bed were extracted from 1952, 1984, and 1999 modeled surfaces (see Chapter 3) along the respective channel thalwegs. Thalwegs were defined by digitizing the 'most probable' flow path along the bed between exposed bars and islands. Although well-defined within confined reaches, avulsions that produce isolated pools, the lack of a distinct dominant channel where the flow bifurcates and where flow migrates across shallow riffles neccessitates some qualitative judgements with respect to determining thalweg position. -79-  Figure 4-6. Thalweg profiles for 1952, 1984 and 1999 produced from bathymetric data. -80-  However, this approach produces a smoother, more realistic flow path than one which simply connects the lowest-elevation points along successive sounding lines. The plotted profdes (see Figure 4-6) display the 'irregularly sinuous' pattern characteristic of wandering channels (cf. Church, 1983a). Both avulsions and lateral migration are evident between Sumas Mtn. and Laidlaw, where most of the bed material load is deposited. Beyond these limits, modest volumes of mobilisable sediment are stored, channel banks are comparatively stable, and the thalweg has remained remarkably fixed in position. Bed profiles for each date were sampled in the GIS at regular 100 metre intervals from the interpolated grids between Mission Railway Bridge and the upstream extent of each surface model. However, as the extraction procedure automatically records additional points along the vertices that define each digitized thalweg profile, the average sampling interval decreases to roughly 50 metres. Although this strategy depends on the accuracy of the grid interpolation, it is computationally efficient, produces a much denser profile than could be directly obtained from the sounding points, and avoids sampling difficulties where sounding data are inadequate to represent channel curvature. The geographic coordinates and elevation of each sampled point were written to a spreadsheet to calculate the elevation at cumulative distances upstream. A plot of these data clearly shows a pattern of alternating riffles and pools, several of which exceed 30 metres depth at high discharge (Figure 4-7). The deepest pools are usually found near the outside boundary of sharp bends with erosion-resistant banks. Neill (1973) similarly observed natural deep scour holes at sharp bends on the wandering Athabasca River, Alberta. At such locations, major eddies form and descend, increasing shear stress along the bed (Richards, 1982) and can result in extreme localized bed scour (Mlynarcyzk and Rotnicki, 1989). A water surface profile collected at the time of the 1999 survey is shown for reference, and exhibits a distinct flattening near Sumas Mountain (about 16 km upstream) demarkating the transition zone where sand becomes the dominant bed material. Although it is clear that several riffles and pools have remained remarkably stable over the period of study, much of the bed is unstable. Most unstable pools are found along the outside edge of sharp bends within sedimentation zones. Pool instability results from migration of the thalweg in response to the downstream staging of bed material stored within transient bars (cf. Church and Jones, 1982). This process fills existing pools and creates new ones through erosional processes that maintain flow conveyance. The stability of the complete pool-bar-riffle sequence, therefore, is related to the supply and mobility of sediments within the active channel zone. Bars develop on the -81 -  l  O  O  ^  O  ^  O  v  -  i  O  1  ^  (JSBUI) UOlJBA3ia  Figure 4-7. Longitudinal bed profiles of lower Fraser River as measured along the channel thalweg for different dates. The profiles clearly show a pattern of alternating pools and riffles. Differences in thalweg length may affect the relative position of features (horizontal distance) providing the illusion of modest up- or downstream movement.  -82-  upstream edge of stable riffles, producing an apparently stable morphology that may persist at the temporal scale for river regime stability (Church and Jones, 1982). On Fraser River, this scale corresponds to decades or more at the spatial scale of large bar/island complexes, although local (unit bar scale) changes are common during the annual freshet. Major avulsions can also occur within a freshet, although they are a culminating response to many years of morphologic 'preparation'. A direct overlay of channel thalwegs (Figure 4-8) suggests an apparent decrease in the number of large pools, an increase in the elevation of several major riffles, and a decline in topographic variability (the range between minimum and maximum channel depths) between 1952 and 1999. Superficially, these characteristics imply that channel morphology has become simplified. However, this conclusion ignores the degree of variability along each profile between major riffles and pools, and that at unsampled flow bifurcations within the channel. Since each profile is a sample of bed topography, the total number of deep pools and riffles on each date remains unknown. Riffles and pools have been found to exhibit regular downstream spacing in alluvial channels that scales with meander wavelength and channel width (Leopold and Wolman, 1957, Keller and Melhorn, 1978, Church and Jones, 1982). Autoregressive models of bed topography have confirmed cyclic patterns of bedform development (Richards, 1976; Church and Jones, 1982). The well known (and widely accepted) association is that pool spacing averages 5 to 7 channel widths in self-formed channels, although both width and spacing can be highly variable over short distances (Knighton, 1998). An objective and consistent definition of riffles and pools is prerequisite to quantitative description (Madej, 1999). That a wide variety of criteria have been presented to define pools (cf. Montgomery et al., 1995; Madej, 1999) illustrates that this task can be problematic, but two contrasting techniques appear useful. Richards (1976) presents a conceptually simply approach for accomplishing this task whereby a linear trend line is fitted to the bed profile so that riffles (elevations above the trend) become easily distinguished from pools (elevations below the trend). Lisle (1987a) presents an attractive method for defining pools independent of discharge based upon the concept of residual depth. Assuming zero discharge, a pool is defined as the area below a constant line of elevation extending upstream from each riffle crest to the next riffle and residual depth is the maximum depth of each defined pool. Neither technique seems entirely appropriate along lower Fraser River. Defined 'pools' vary greatly in size, with some sufficiently small to seemingly ignore, while some larger pools -83-  o  ©  o  O  (iSBUl) UOpBA3{3  Figure 4-8. Plot of superimposed longitudinal profiles showing variability in the location of riffles and pools over time. The 1952 profile has been 'stretched' to be equivalent with the 1999 profile length to facilitate direct comparison. Note the apparent decrease in the size and frequency of very large pools since 1952. The arrows indicate several prominent riffles that have developed over the past 47 years (approximate locations are noted).  -84-  appear to consist of two or more 'sub-pools' having intermediate riffle crest elevations similar to, but less than the downstream riffle crest elevation. Further, it is clear that spatially dense surveys along the long profile will produce higher riffle and pool counts than a more generalized survey because very small riffles and pools are resolved. Therefore, each approach was modified using qualitative judgements to include only reasonable large, distinct pools. Results show that each technique produces a similar count, averaging 27 in both 1952 and 1999, and 20 in the shorter surveyed reach in 1984. Over the length of the gravel reach, pools are spaced every three channel widths for each date, about twice the frequency generally accepted for gravel bed rivers. However, if an equvalent single channel regime width (~ 525 m at Mission and Agassiz gauges) is applied, pool spacing converges towards modal values. Similar results are realized if riffle/pool spacing is normalized by the braid index (roughly 2.0 for lower Fraser River), taken as the average number of flow channels per unit length. Linear regression lines drawn through the longitudinal profiles also reveal changes in the mean slope of the gravel reach. The slope decreased from 4.59 m / km in 1952 to 4.41 m / km in 1999, corresponding to an increase in thalweg length (by roughly 1 km), hence channel sinuosity (from 1.140 to 1.155). Increasing sinuosity has been associated with a reduced supply of mobile sediments (cf. Schumm, 1977) but such statements must be tempered unless corroborating evidence can be provided. Ham and Church (2002) compared the area beneath each profile (analagous to cross-section plots above a minimum base level) and found a decline in area upstream of Agassiz (signifying degradation) and a downstream increase (signifying aggradation) which is consistent with a lowering of the average profile. More robust investigation of these trends requires an examination of bed level changes that better represents the entire active channel zone. The following sections of this chapter (cross-section and hypsometric analysis) and the following chapter (distributed sediment budget) provide appropriate tools for this inquiry.  4.4 Cross-section morphology In addition to adjusting channel geometry through changes in channel width, alluvial channels may adjust their depth and velocity to maintain flow conveyance. For any given flow, a reduction in active channel width will generally result in increased flow depth, since this parameter increases proportionately greater than velocity . Along rivers with mobile beds, depth can increase 2  2. This can be demonstrated, for example, from the Manning and Chezy equations where depth varies with the 1.5th and 2nd power of velocity respectively. Furthermore, the dominant factors influencing flow velocity (flow resistance, bed slope and channel pattern) would remain fairly stable over years to decades on Fraser River.  -85-  through degradation of the bed, or by raising water levels for any given flow. As the lower Fraser River gravel reach is known to be aggrading, the latter circumstance must be occuring, although the bed may be degrading locally. Simulations of design flood profdes confirm this general pattern over the past several decades (UMA Engineering, 2001). Specific gauge analysis has also been used to demonstrate persistent bed level changes since gauges are typically sited where channel banks are stable (hence channel width remains fixed). Along wandering channels, however, such sites are typically associated with sediment transport zones, rather than zones of significant erosion or deposition, hence this avenue of study is unlikely to reveal meaningful interpretations. Adjustments to the shape of the channel can be determined by comparing historic crosssection profiles of the river. These profiles were derived from the bathymetric and topographic surveys completed in 1952, 1984 and 1999. Conventional cross-section comparisons require survey lines that match between survey dates — in practice, this result is achieved by surveying the channel banks and bed between fixed endpoints. However, matching profile lines could not simply be derived from the available data because each survey encompasses a different spatial area of channel, the density of the elevations varies between the surveys, and the transverse survey lines are rarely coincident. To overcome this limitation imposed by the data, a GIS-based technique to derive cross-sections was developed to simulate fixed sections. A triangulated irregular network (TIN) model of the bed and banks was constructed for each survey date using available bathymetric and altimetric elevations, photogrammetricallyderived elevations, hand-interpolated bed contours and mapped floodplain contours. The models are based on the Delaunay method, wherein each point is connected to its nearest two neighbours to form a series of continuous facets that represent a surface. The main advantages of the TIN versus regular gridded interpolations are that the data structure is based on irregularly spaced point, line and polygon data (ESRI, 1991) and actual elevation values are preserved in the output model (Burrough and McDonnell, 1998). A series of lines representing the location of each desired crosssection was then digitized at roughly 800 metres (four times the spacing of the 1999 survey) from the right-to-left channel banks (see Figure 4-9) to roughly correspond with active channel width. The endpoints of each transect were placed at stable regions beyond the outer extent of active channel margins (between 1949 and 1999) such that areas of intervening erosion or deposition would be captured. Cross-section profiles were extracted from each TIN by interpolating elevation values at regular (20 metre) intervals along each cross-section. The elevation at each sampling distance is -86-  -87-  estimated from linear interpolation between survey points where the data are not coincident. Since the location of the transects was based on the 1999 tracklines (taking advantage of the regular spacing and spatial completeness of the survey) cross-section profiles extracted from the 1999 TIN more faithfully preserve actual bed elevations than those of the other years. A small subset of crosssections profiles is plotted in Figure 4-10 to illustrate the magnitude of vertical and lateral instability that occurs along the river, especially within major sedimentation zones. The complete set of extracted profdes is given in Appendix A. The profiles are each shown with a reference line that approximates the bankfull waterlevel surface. The region below each line, therefore, represents the active channel zone crosssection area, while the region above represents island or floodplain cross-section area. The reference line corresponds to the average elevation of young island surfaces along the channel. Although the vegetation trim line can be identified for each survey from airphoto mapping, the young island elevations are a convenient surrogate because the high water mark along outer channel banks is not always easily defined, especially along confined sections bounded by high terraces or bedrock (e.g. RB of XS 4, Figure 4-10). High bar surfaces recently stripped of vegetation, steep unvegetated banks, and natural variation in the age (hence elevation) of vegetated surfaces also introduces considerable variability in the elevation of the trim line across most crosssections. This makes it difficult to establish a mean bankfull elevation for an individual section, to compare areas above and below the water surface between survey dates, and creates an unrealistic downstream elevation profile. The reference elevations were derived from the laser-altimetry survey data collected in 1999. Channel maps produced from 1949, 1971 and 1999 low-water aerial photography were superimposed to identify three surface types: old bars (polygons mapped as gravel bars on all dates); mature islands (polygons mapped as island or floodplain on the all dates); and young islands (polygons mapped as bars in 1971 and island or floodplain in 1999). The overlay map was subsequently divided into 1-km computing cells and overlaid with the altimetry data such that elevations could be averaged for each surface type and plotted by distance along the channel (Figure 4-11). A best-fit polynomial line was computed for each surface class to smooth scatter or anomalies due to deposit age differences (some island surfaces are older than 50 years) or to insufficient data (some 1-km cells have no islands or bars). The difference in elevation between old bar and young island surfaces represents the thickness of recent overbank (fine sand and silt) deposition, whilst the difference in surface elevation between young and mature islands represents -88-  Wandering reach  Straight / meandering reach  xs 82 0  200 400 600 800 1000  -25-|—i—|—i—|—i—|—i—|—i—|  0  0  O-j  200 400 600 800 1000  0  500  1  1  1000 1500 2000  1  1  400 600 800  1000  1  1  1  1000 1500 2000  500  1000 1500 2000  Figure 4-10. Selected cross-section plots to show complexity of morphology and dynamics of scour, fill, and thalweg migration in wandering reaches compared to partly confined straight and meandering reaches.  -89-  1  500  1  200  1  1  1  1  1  2500  3000  1  1  2500  1  1  2500  ]  3000  r  3000  I I 200 400 600 800 1000 Distance from RB (m) I  1  1  the thickness of erodible overbank deposits (adopted for sediment budget calculations - see Chapter 5). These thicknesses (1 and 3 metres respectively) closely correspond with limited field measurements made by McLean (1990, p. 74). Active cross-section area was calculated for all 78 profiles, and subdivided into forested (island) and unvegetated (bar) groups. A further division was made for profiles located up- and downstream of Agassiz-Rosedale Bridge for comparison with the observed decrease in channel slope (Section 4.3). Since it is known that vegetation is associated with greater bank stability (cf. Millar, 2000) and influences channel morphology (cf. Hickin, 1984; Murray and Paola, 2003), differences in area (and shape) are expected between the two groups. Forested profiles were identified from areal photography because defining islands based strictly on the reference elevation was found to exclude some cross-sections where the elevation of immature islands was equal to or less than the interpolated bankfull elevation. It is assumed that a gain in active channel area over time indicates bank erosion if width at that location also increases, and a deepening of the channel if width is stable or declines. A gain in area above the reference water-level is associated with overbank fine deposition and island growth. The data are summarized in Figure 4-12 and Table 43. The significance of changes in cross-sectional area over time was examined through a statistical comparison of section mean and variance values. Clearly, any channel with a labile bed and banks will exhibit sufficient variability in cross-sectional form that half the sections will be larger than the mean area, on average, with the remainder smaller. Along morphologically complex channels, however, it is expected that the range of size deviations from the mean is much greater than in a channel with a comparatively simple morphology. Therefore, if the complexity of channel morphology is changing over time, both the mean and variance of cross-section area should be altered. Changes in these properties are tested using the t-test for means, and the F-test for variance. Results for cross-section area below the bankfull reference level are summarized in Table 4-4. Mean active channel area between Mission and Agassiz increased from 1952 to 1984 for both forested and unvegetated sections, although the change is not statistically significant. However, as the active channel width was slightly narrower in 1984 compared to 1952, there must have been some deepening of the channel to maintain flow conveyance. The area above the reference surface also increased despite an overall (albeit modest) decline in island area. While the change is not statistically significant, growth is consistent with the known trend of generally aboveaverage flood flows, a neccessary prerequisite for vertical aggradation. From 1984 to 1999, mean -91 -  8000 7000 6000 5000 E. 4000 3000 2000 1000 0  10000 9000 8000 7000 CvP  6000  10  5000  E  <u  ra  4000 3000 2000 1000 0  3000  2500  _  g ra  2000  CD  m  1500  1000  500  Figure 4-12. A c t i v e channel cross-section area (a) below reference line along x-sections with islands, (b) below reference line along unvegetated x-sections, and (c) above the reference line (solid area). Refer to Figure 4-9 for location (ordering is sequential, but not continuous). _92-  cross-section area declined amongst vegetated sections. This change is consistent with the observed increase in island area within the same period on most reaches. However, since the area above the reference water level did not also increase, island and floodplain growth must be associated with younger, low elevation surfaces. Lack of growth on mature, higher elevation surfaces is coincident with the general lack of large floods since the late 1970s. During the same period, there was an areal increase amongst unvegetated sections, the magnitude of which was sufficient to produce an overall balance in the mean area of all cross-sections. There is also weak evidence to suggest that the variability of vegetated cross-section area has been in decline downstream of Agassiz. This finding implies that individual islands have become larger in total area, but that they have correspondingly declined in number. In-filling and vegetation growth along sidechannels separating adjacent islands may be partly responsible for this pattern. Overall, channel width has declined by 82 metres (8%) since 1949, but channel area has not. This provides clear evidence that the channel has become deeper, on average, over the past half-century. Dividing mean active channel area by mean active channel width illustrates these changes. Average depth is estimated as 4.84 metres in 1952,5.07 metres in 1984, and increased to 5.50 metres by 1999. Upstream of Agassiz, mean active channel area below the bankfull datum increased along both vegetated and non-vegetated cross-sections between 1952 and 1999. Although the change is not statistically significant, an increasing area and corresponding narrowing along most upper reaches (18% overall) indicates a deepening channel overall. As calculated above, estimated mean channel depth increased from 3.53 to 4.55 metres. That the channel appears to be degrading more rapidly upstream than downtream of Agassiz is consistent with the overall decline in channel profile (Section 4.2). There is also a distinct pattern of decreased area variability amongst upstream reaches that is statistically significant overall. The decline is also statistically significant for vegetated cross-sections. The overall increase in island area implies a morphologic development toward fewer, larger islands, paralleling the pattern downstream. A significant decrease in variability above the high water mark provides additional evidence to support this claim. At present, the channel is narrower and deeper with reduced ecotone habitat as smaller island clusters have coalesced into larger island complexes.  -93-  Table 4-3. Summary statistics for total active cross-section area below the reference datum along the extracted crosssections.  n = 50 x = 4730 m s.d. = 1298 m  All cross-sections  n = 50 x = 4924 m s.d. = 1390 m  2  n = 34 x = 4387 m s.d. = 1143 m  n = 34 x = 4623 m s.d. = 1293 m  2  n= 16 x = 5460 m s.d. = 1340 m  Unvegetated sections  n= 16 x - 5563 m s.d. = 1412 m  2  1952 n = 28 x = 3485 m s.d. = 1310 m  All cross-sections  n = 22 x - 3608 m s.d. = 1437 m n =6 x = 3033 m s.d. = 528 m  n= 16 x = 6085 m s.d. - 1482 m  1984  1999  no data  n = 28 x = 3695 m s.d. = 894 m  2  2  2  no data  n = 22 x = 3858 m s.d. = 928 m  2  2  2  Unvegetated sections  2  2  2  2  Vegetated sections  n = 34 x = 4376 m s.d. = 969 m  2  2  2  2  2  AGASSIZ-LAIDLAW  2  2  2  n = 50 x = 4923 m s.d. = 1398 m 2  2  2  Vegetated sections  1999  1984  1952  MISSION-AGASSIZ  2  no data  n=6 x = 3100 m s.d. = 392 m  2  2  2  2  2  Table 4-4. Tests for significant differences between mean (t) and variance (f) of active cross-section area below the reference datum. P is the probability that the test statistic is less than the critical value (shown in bold type). MISSION-AGASSIZ  1984-1999  7952 -1984  1952 -1999  All x-sections  F= 1.15 p = 0.32  t = 0.72 p = 0.47  F= 1.01 p = 0.48  t = 0.002 p = 0.999  F= 1.16 p = 0.30  t = 0.72 p = 0.48  Vegetated sections  F= 1.28 p = 0.24  t = 0.80 p = 0.43  F= 1.78 p = 0.05  t = 0.89 p = 0.38  F= 1.39 p = 0.17  t = 0.04 p = 0.97  Unvegetated sections  F= 1.11 p = 0.42  t = 0.21 p = 0.83  F= 1.10 p = 0.43  t= 1.02 p = 0.32  F= 1.22 p = 0.35  t= 1.25 p = 0.22  AGASSIZ-LAIDLAW 1952 - 1999  All x-sections F = 2.15 p = 0.03  t = 0.70 p = 0.49  Vegetated F = 2.40 p = 0.03  -94-  t = 0.68 p = 0.50  Unvegetated F= 1.82 p = 0.26  t = 0.25 p = 0.81  4.5 Hypsometry The ecological productivity and diversity of river channels is directly related to the extent and variety of available physical habitats (Kondolf and Wilcock, 1996). Indeed, most channel restoration efforts emphasize increasing morphologic complexity by excavating pools and adding large roughness elements to create aquatic habitats (Kellerhals and Miles, 1996; Kondolf, 1998). The wandering planform in particular provides excellent habitat for fishes and aquatic insects because the complexity of characteristic channel forms provides a variety of depth and velocity conditions over the full range of discharge (Payne and Lapointe, 1997; Church, 2002). Long-term mensuration of changes to the physical habitat can provide a surrogate means for evaluating potential habitat loss. While the isolation of backchannels and ecotone habitats from the river corridor is apparent, changes within the active channel zone are not as easily recorded because they are three-dimensional in nature. Hypsometry provides a convenient means for summarizing the cumulative distribution of elevations for a defined landscape above a reference plane. A hypsometric plot can reveal details of the evolutionary [morphologic] history of a geomorphic surface (cf. Strahler, 1952) and is particularly valuable for comparing landform development over time. In recent years, landscape evolution modelers have begun to adopt the technique for comparing the evolution of simulated to real landscapes (Willgoose and Hancock, 1998), a task made increasingly facile by the availability of elevation data in digital form. The distribution of active channel area for different depth intervals can reveal information about the nature of morphologic complexity within a given river, while temporal comparisons should provide a means for illustrating the response of a river to changes in natural (scour and fill) or anthropogenic (gravel removals and bank hardening) impacts. In the simplest case of a rectangular canal with fixed banks and no bed material transport, the entire bed would be covered by water of constant depth for any given discharge. The presence of mobile sediments and deformable banks provides increasingly complex channel morphologies with a broader distribution of area-depth classes for a given discharge. It is therefore supposed that channels with highly complex morphologies (i.e. braided and wandering reaches) exhibit a higher degree of spatial complexity than simpler morphologies (straight or meandering reaches). For this investigation, the total area for different depth bins is calculated for bankfull discharge, a procedure that is analogous to the representation of ocean depths (cf. Harrison, 1998). A GIS coverage containing water surface elevations was constructed by adding a new database -95-  field, surface elevation, to an existing coverage depicting 1-km long computing cells. For each 1km cell, the water-surface elevation, Zs, was estimated using the polynomial equation fitted to the young island elevations. The 1952 and 1999 interpolated grids were then converted to coverages and combined with the water surface coverage, and morphologic maps produced from the airphotos. The composite GIS layer, therefore, gives the bed surface elevation for each computed grid cell, and a single water-surface elevation for all grid cells within each 1-km computing cell. Bankfull depths were then calculated as Zs - Zbed for all grid cells defining the active channel zone, including channel islands because some immature islands may be inundated at bankfull flow. A histogram showing changes in the area-depth distribution of the gravel reach from 1952 to 1999 is shown in Figure 4-13. This plot demonstrates a loss of shallow-water habitats at high flow (less than 4 metres total depth) and a corresponding increase in deep-water habitats, a pattern that is consistent with evidence of channel degradation presented in section 4-3. It should be realized, however, that this assertion is complicated by a lack of knowlege with respect to actual water levels in 1952. Water levels would have been somewhat lower, on average, in 1952 because of aggradation that has since occurred. Changes in flow alignment will also have affected water levels. Therefore, the loss of shallow-water area may be exaggerated. Correspondingly, the observed increase in deeper water environments over this same period may be over-estimated, but it provides further evidence of a trend towards a narrower, deeper, channel. To examine the differences in hypsometry that are expected for distinct morphologies, confined and wandering sub-reaches were partitioned from the complete dataset. The single-thread meandering section includes 10 km of channel upstream from Mission, plus 4 km downstream from Laidlaw. A continuous 25 km long stretch from lower Yaalstrick Island to Big Bar near the Agassiz-Rosedale Bridge was used to represent the wandering planforrn. Most artificial training works and gravel mining activities have been undertaken in this reach of channel (see Figure 2-5). The confined reaches actually exhibit a wider range of depth classes than the wandering reaches because they include the deepest pools along the lower river. Roughly 80-85% of the active channel lies at depths between 4 and 16 metres in both periods, while the area of shallow (0-4 m) areas increased by half from 11% in 1952 to 16% by 1999. Overall, there is no clear trend in channel depth changes in these reaches, reflecting their relatively small total area, and the dominance of sandy bed-material transport and deposition, especially near Mission. Although the available data do not permit direct mensuration, changes in the frequency distibution of some channel depth intervals possibly occur every few years in response to sediment transfer. -96-  30000000  25000000  20000000  15000000  >  10000000  5000000  1800000  Confined section  1600000 1400000  n_  1200000  ns  o  •S  n  1000000 800000 600000 400000 200000  14000000  It  .fl.rl.rl.m. Wandering section  12000000  10000000  ID  c e  8000000  cS  o to >  6000000  4000000  2000000  <-2  <2  <4  <6  1™ <8  <10 <12  <14 <16  <18  <20 <22 <24 <26 <28 <30  water depth (2 m intervals)  Figure 4-13. Histograms of flow depths (relative to 1999 bankfull discharge) within the active channel zone. Confined and wandering sections do not represent a complete inventory of total active channel area but illustrate comparative differences due to morphology.  The wandering sub-section displays a similar overall pattern to that of the entire gravel reach because this is the dominant morphology. Roughly 96% of the active channel lies at depths less than 8 metres in both periods. Very deep pools (>8 metres depth) occupy the remaining 4% of the active channel, and appear to be a minor component of morphologic complexity in wandering reaches. Shallow (0-4 m) areas show a similar decline in proportion compared to the entire gravel reach, but changes are mainly pronounced in the 0 to 2 metre depth class (which is the most sensitive to bar aggradation). Similarly, the apparent loss of area from 2 to 6 metres may be influenced by bed aggradation since the actual change is fairly modest. The differences in morphologic complexity between the two channel typologies is further emphasized using a hypsogram plot (Figure 4-14). The dominance of shallow depths as a percentage of total active channel area on the wandering planforrn is clearly shown. 4.6 S u m m a r y There are several lines of evidence which indicate that the morphologic complexity of the gravel reach has become reduced over the past century in response to the cumulative effects of anthropogenic influences. The channel has narrowed more quickly than can be explained by a corresponding decline in bankfull discharge over the period of record, hence must be responding to other forces. Since changes in width are fairly uniform amongst the different study reaches, a single process dominates, though the impacts of individual factors are difficult to isolate. However, it is known that much of the loss in active channel area early in the century is the direct consequence of dams and flood gate construction that isolated prominent sloughs and secondary channels, which have since silted up over time and become vegetated. As many of the secondary channels encompassed large islands, this isolation also led to a dramatic loss of island area early in the century, as these deposits were effectively converted to floodplain. The loss of island area continued to the 1970s, consistent with generally above-average flood flows, but their total area has since recovered to levels not seen since before the dyking program was expanded. Superficially, this pattern implies that dyking, bank hardening and woody debris removal have had no significant impact on the morphologic development of channel islands, but available evidence indicates that mature, riparian forests are becoming replaced by young vegetation on newly established islands. This has been accompanied by a trend towards fewer, but larger islands, resulting in a loss of ecotone habitat.  -98-  1  normalized cumulative area  Figure 4-14. Hypsogram for a confined and wandering section of channel from 1999. Confined sections have depthsrangingfrom 0 to 26 metres, while wandering sections range from 0 to 22 metres.  -99-  Channel narrowing and changes in island morphology have also modified the topographic characteristics of the channel. The gravel reach has become less steep, with degradation upstream of Agassiz-Rosedale and corresponding aggradation in downtream reaches. Analysis of extracted cross-section profiles demonstrates that the channel has actually degraded in these reaches, but the rate is more modest than observed upstream. Hypsometric analysis reveals that sedimentation in lower reaches has become concentrated on higher bar surfaces, reducing the area of shallow habitats at bankfull discharge. However, given the absence of a bankfull water-surface profile for 1952, it is not possible to definitively state that aggradation since has reduced channel capacity, or that gravel mining and bank hardening have reduced morphologic complexity. The seemingly remarkable resilience of the river to respond markedly to environmental and anthropogenic changes since 1952 may primarily reflect inertia given its large size. If the volume of material stored within the active channel zone is far greater than the average annual influx of bed material, or of sand and gravel removals, the morphology may alternatively be dominated by internal modifications over periods of decades. A sediment budget is presented in the following chapter and is used to investigate this supposition more fully.  -100-  Chapter 5: Bed material transport 5.1 Introduction In this chapter, an estimate of the bed material transport rate along the lower Fraser River is provided over the past half-century. Sediment transfer and deposition dominate the morphologic development of the gravel reach and present a number of engineering and management challenges because of development on the floodplain. The primary response to combat river morphologic development has been to harden outer channel banks and selectively remove gravels. Establishing the relation between the nature and magnitude of these channel modifications and the sediment transport regime provides a means for evaluating the sensitivity of the channel to environmental changes, and provides a tool for predicting the response to additional forcing. Given the large size of the river, many decades may be required to resolve the impacts, if observable, of these activities on the rates and volumes of sediment input and transfer, and hence channel stability. Assuming these activities will continue, it is important to understand the past effectiveness (or failure) of these applied management strategies for protecting developed areas while maintaining the quality of available aquatic habitat. Establishing a detailed bed material sediment budget provides a means to pursue this line of inquiry and establish guidelines for future management decisions. Measured morphologic changes are related directly to sediment transport rates within the sediment budget framework. This technique is suitable for large gravel-bedded channels in which significant channel changes occur over years to decades. Estimates of bed material transport for the lower Fraser River are provided for 1952 to 1984,1984 to 1999, and 1999 to 2003, and are used to illustrate the dynamics of sedimentation over periods commensurate with timescales for management practices. The availability of recent (1999 and 2003) bathymetry and altimetry data further provides a means for investigating the linkages between channel deformation and sediment transport at a resolution that may be unique for such a large river system. The sediment budget estimates are subject to a number of limitations, including the reliability of historic channel and floodplain topographic surveys, the accuracy with which these data can be converted to interpolated surfaces that reflect actual topographic variability and the reliability of assumptions required to parse complete topographic models of the bed and banks into channel and bank sediments, and hence to estimates of the total load distributed into washload, coarse sand and graven;! components. These problems are addressed in a discussion of sediment  - 101 -  budget errors at the end of the chapter. The validity of the estimates is further examined by comparing results presented herein with those of past studies. 5.2 Sediment budgets f r o m morphologic change 5.2.1 G e n e r a l approach  The basis for relating changes in channel morphology to sediment transport is the morphologic sediment budget, essentially an accounting procedure for tabulating net sediment exchanges within a defined reach of channel. The approach is most suitable for describing the movement of bed material, since this is the material that determines channel morphology and governs stability. An estimate of washload can be made under certain conditions (Church et al., 1987; Ashmore and Church, 1998) although none is provided in this thesis. On the gravel reach of lower Fraser River, exchanges of bed material include sediments entrained upstream of Laidlaw and transported downstream, channel scour and fill, erosion and construction of island and floodplain deposits and output at the downstream limit of the study area. To complete the budget, sediment evacuations through dredging (documented in Chapter 2) are added back to the reach from which they were removed as a depositional volume. The incomplete knowledge of all past gravel removal volumes represents a negative bias in the sediment budget (meaning transport estimates represent a minimum) although the magnitude of this bias remains unknown. Nominally, the contribution from major tributaries must also be accounted for. However, as neither of the two major tributaries that enter the gravel reach (Chilliwack-Vedder and Harrison Rivers) contributes material larger than sand which mainly travels through the reach in suspension (McLean et al., 1999), this volume is ignored. The sediment budget can be summarily expressed as: V = Vi-(l-p)AV-V 0  d  where V is the volumetric sediment output, V; is volumetric sediment input to the reach during 0  some specified time period, and V is the volume removed from dredging. The storage term, AV, d  is measured as the net difference between scour and fill of the channel bed and banks, and is estimated from differences in modeled bathymetry surfaces over the period between surveys. The equation can be reduced to a mean transport rate by dividing all terms by the time between successive surveys. An additional adjustment for sediment porosity (1-p) or unit weight (t/m ) is 3  conventionally included to express all terms as mineral volumes. Finally, total volumetric changes can be converted to bed material volumes by adjusting for <I>, the proportion of material finer than - 102-  medium [wash material] sand (<0.177 mm). Within the active channel zone, the adjustment is negligible, but considerable volumes of fines can be stored overbank on island and floodplain surfaces where they become trapped or settle at high discharge due to lower flow velocities. McLean (1990) assumed that the transport of gravel past Mission was negligible, noting that material >2 mm comprises less than 5% of channel bed material below Mission, and coarse gravels are nearly absent downstream. Given this assumption, the transport rate of bed material at Mission was assigned a value of zero, and calculations of sediment influx were propagated upstream. However, it is not clear that the assumption is entirely valid because some gravel has been found in sand dune troughs downstream of Mission. Although this material could be derived locally (i.e. from erosion of terrace deposits), the possibility remains that gravel is actually transported past Mission. Bed material particle-size analysis data published by the Water Survey of Canada will be examined to revisit this question. Nevertheless, the zero transport assumption remains desirable because it implies that all material that enters the gravel reach is trapped, hence there is no practical upper limit for the temporal interval between successive surveys. Further, the requirement to estimate a reference transport rate is eliminated. In channels with no defined 'zero' boundary condition, the possibility for unmeasured throughput increases with the length of the survey interval, and only a lower bound estimate of transport can be made (Carson and Griffiths, 1989; Ham and Church, 2000). Nevertheless, details on sediment transfer and channel deformation within lengthy intersurvey intervals are, of course, limited by compensating scour and fill. The complete sediment budget expressed as a mean transport rate is: Qi = Qo + Qd +  AV / At + (1-pXl-O) A V / At C  f  where Q represents the mean flux rate, while V and V define channel and floodplain changes c  f  respectively. If the value of Qj exceeds Q at any location, the channel is locally aggrading. The 0  sediment budget can further be defined for any arbitrary length of channel smaller than the total study reach length. In this thesis, the channel is divided into 65 computing cells of approximately 1-km in length, allowing the transport rate to be estimated at any location within the study reach (Figure 5-1). This approach is useful for illustrating the variability in transport that can occur downstream and provides a means to relate the displacement and storage of sediment to channel stability, details of which are provided in the following chapter.  - 103-  -104-  5.2.2 Sediment budget construction  The construction of period sediment budgets for lower Fraser River is based on a sequence of computing steps that have been refined and adjusted to minimize error margins within the constraints of available information (cf. Ham and Church, 2003). Topographic models of the bed and adjacent floodplain were constructed for 1952, 1984 and 1999 following the conventions described in Chapter 3. As the 1999 data were collected at higher spatial density along individual sounding lines than in previous surveys, they were thinned to be compatible with the 20-metre spacing of the earlier surveys. Individual datasets should have similar spatial density to ensure selfconsistent comparisons between modeled surfaces over time. However, as the 2003 soundings were collected at roughly the same density as in 1999 (although sounding line spacing was reduced from 200 to 100 metres) no data thinning was required for the 1999-2003 intersurvey comparison. In addition, the availability of 1-metre ground resolution LiDAR data on stable islands and floodplains enabled construction of topographic surfaces for both 1999 and 2003 that are more accurate than for the older datasets. Therefore, the sediment budget results presented herein are based on modeled surfaces computed at higher (10-metre) spatial resolution. This change preserves more of the actual topographic variability that exists on bars and islands, such asridgesand chute channels, that is otherwise'smoothed out' at lower resolution. The net of bed level changes for each intersurvey period is calculated by subtracting consecutive modeled surfaces using GIS tools. These 'surfaces of difference' illustrate the complex patterns of channel scour and fill that occur in response to the downstream staging of bed material (Figure 5-2). Each was then superimposed with a masking layer to replace interpolated differences with a 'no-data' value in regions where the interpolation was known to be weak due to insufficient data, outside the margins of known channel change (i.e. on large islands and floodplain) or where the modeled region extended beyond the area of principal interest (i.e. Vedder Canal, major sloughs). The same mask was applied to the 1952-84, 1984-99 and 1952-99 surfaces to ensure all computing areas were identical and contain the actual active channel zone over the past 50 years. The masked surfaces were also overlaid with morphologic maps produced from available airphotos nearest to the date of the soundings. The morphologic maps permit each 20metre grid cell to be associated with a particular type of channel transition, including: • channel scour and fill (polygons coded as water or gravel bar on both mapping dates) • bank erosion (island or floodplain on early date, water or gravel bar on later date • bank deposition (water or gravel bar on early date, island or floodplain on later date) - 105 -  Figure 5-2. Surfaces of difference from individual topographic models between Mission and Agassiz showing complex patterns of scour and fill.  • no change (island orfloodplainon both dates; assumed that these are stable surfaces)  It is assumed that a 'no change' polygon has a constant surface elevation over the study period, so these may be ignored (although some surface compaction or overbank sedimentation could technically occur). The remaining cells, therefore, define common regions of channel scour and fill, bank erosion, and bank deposition: the constituent components of morphological change. Cells divided by more than one transition type are assigned separate transition codes for each cell 'piece'. Since the dates of the airphotos and bathymetry surveys do not correspond exactly, some codification error is possible. Two additional processes, vegetation stripping (island or floodplain to gravel bar) and recovery (gravel bar to island or floodplain) were also identified, but are collapsed into the bank erosion and deposition categories respectively since the data are reduced in an identical manner. Data reduction refers to the process of converting gross volumetric changes into estimates of bed material (sand and gravel) and wash material (which is eliminated from the calculations). Procedural details vary somewhat by transition type and are described below. The volumetric calculations for channel scour and fill are reasonably straightforward. The total volume change for an individual grid cell is simply the product of the elevation change and the cell (or cell piece) area. For each of the 65 Lkm reaches, individual channel grid cell volumes are summed, giving the net of scour or fill over the intersurvey interval. This total volume consists of bed material only — material finer than medium sand is not deposited on the active bed. Although silt lenses have been observed on some exposed bars at low flows, the small size, and accordingly volume, of these deposits can be considered negligible compared to the volume of mobilized bed material in any computing reach. Grain size data were collected in 1983 along bar heads and flanks throughout the entire gravel reach (McLean, 1990) while suites of 9 20-kg bulk samples were collected from upper-, mid-, and lower-bar sites along the reach in 2000, then combined to yield 3 measurements for each bar. Church et al. (2001) found no systematic difference between the data collected in the different years, allowing all 70 samples to be pooled and plotted as grain size versus distance upstream (Figure 5-3). An additional division of bed material into sand and gravel fractions is also presented based on the estimated proportion of sand along the gravel reach as determined from each field sample. Afittedtrend line is used to estimate the sand fraction for different reaches because of considerable scatter, which increases linearly from 17% at Laidlaw to 27% at lower Sumas Mountain (Figure 5-4). The given equation is used to interpolate the sand fraction for each 1-km computing reach between these locations (computing - 107-  Figure 5-3. (A) Surface and (B) subsurface grain size versus distance upstream from Sand Heads. Data are based on pooled samples collected in 1983/84 and 2000.  -108-  cells 9 to 65). Information on grain-size distributions below Sumas Mountain is nearly absent because of sampling limitations. The sand fraction was arbitrarily increased to 30% and 40% in the following two cells downtream. This ratio of sand to gravel is apt to completely fill all pore spaces and minimal gravel transport would continue downstream. The remaining cells were all assigned a value of 95% following values reported by McLean (1990). The adopted values are also given in the sediment budget tables (Tables 5-5 to 5-8). The calculations for bank erosion and deposition volumes are considerably more complex because of a 1 to 3 metre layer of fine sand and silt that overlies the basal layer sand and gravels. Variations in the thickness of this material have been correlated to the age of overlying vegetation at different sites (Boniface, 1985; McLean, 1990). McLean (1990) estimated bank erosion and deposition volumes by multiplying eroded areas by the thickness of the basal gravels based on direct field measurement at each site, although extrapolation was necessary where former islands were completely eroded. An alternative GIS-based approach is adopted that relates site surface elevation to the deposit age (rather than the vegetation age) and the estimated overbank volume is subtracted from the total volumetric change at eroded and deposited island and floodplain sites. Combinations of morphologic maps from low-water airphotos taken in 1949,1971 and 1999 were used to define young island, old bar and old island surfaces, given as: • young island: bar surface in 1971, vegetated by 1999 • old bar: bar surface in both 1949 and 1999 • old island: island (orfloodplain)in both 1949 and 1999 The polygons corresponding to these surface types were intersected with the 1999 LiDAR data, and the elevations of individual spot heights were spatially averaged by each computing reach. These data are plotted as average surface elevations relative to the distance upstream from Mission (cf. Figure 4-11). Best-fit exponential lines are fitted to average scatter or anomalies that may exist because of insufficient elevation data (i.e. several reaches have few or no spot heights on young islands). It is assumed that the difference in elevation between old bar and young island surfaces represents the thickness of recent overbank deposition, while the difference between old bar and old islands represents the thickness of eroded overbank deposits. A value of 0.84 metres, or half the maximum deposited thickness, is adopted as these surfaces will be in various stages of development. The maximum thickness of the overbank sediments on eroded surfaces is 3 metres, considerably larger than the deposited value, but consistent with the previously cited thickness / -110-  vegetation age relation. The same values are applied to all 1-km computing reaches since the bestfit lines are very nearly parallel. Although these conventions' appear straightforward, they do not always lead to a simple adjustment for erosion and deposition volumes. The possibility exists for apparentiy mature surfaces to actually be former (1949) islands that have been eroded and replaced with initially bare bar surfaces that have re-vegetated within the following 50 years. Similarly, some eroded islands may have considerably less than 3 metres of overbank deposits. Following these examples, removing 3 metres of sand along an eroded polygon, multiplied by its erosional area, may result in a volume of overbank material to be subtracted that is larger than the total eroded sand and gravel volume calculated by the GIS. In such cases, the total (GIS) erosional volume is simply assigned to sand. Where the product of the stripping area and the 3 m overbank thickness was found to be smaller than the total volume recorded by the GIS, gravel was also assumed to have been eroded. Similarly, measured revegetation volumes in excess of 0.84 m times the revegetated area were assumed to have also resulted from additional gravel deposition. These gravel volumes are separated to gravel and coarse sand volumes using the same proportions described for channel scour and fill. While it is recognized that assigning different adjustments for erosional and depositional volumes based on the relative age of different surfaces (which can be estimated from available mapping) is theoretically feasible, the technical details of managing such a large database make such an exercise practically onerous. A final calculation is required to adjust the volume of removed overbank fines by  the  proportion of material finer than 0.177 mm. From direct sampling, this is estimated as 70% of the overbank deposit volume. This proportion is considered wash material and is eliminated from the computations. The remaining 30% is combined with the sand fraction associated with erosion and deposition within the basal gravel layer and with the volume of sand removed by mining. The total sand sum is reported in the sediment budget tables. The gravel sum includes the total volume of gravel eroded or deposited within the channel bed and lower banks and the volume removed by mining. A more detailed description of the calculations for sand and gravel volumes for individual computing cells is given in Church et al. (2001) and is repeated in Appendix B. The total sand and gravel sum deposited within each computing cell is also given in the tables. This quantity, divided by the computing cell area, gives the relative change in the average elevation of the bed over the intersurvey period. The final budget is as complete as possible given imperfect knowledge of sand distribution throughout the study reach. This includes not only the thickness of overbank sands, but - Ill -  also the proportion found within the bed and banks. It should be recognized that this limits the precision of both sand and gravel budget estimates.  5.3 Sediment transport 5.3.1 S u m m a r y sediment budgets  A strict summation of all positive and negative elevation changes for modeled surfaces demonstrates the magnitude of total sediment flux (including overbank deposits) over time within the gravel reach (see Table 5-1). Up to 53 million m of sediment was eroded (with 49 million m deposited) between 1952 and 1984, giving a net loss of 4.3 million m . An additional 33 million 3  m was eroded between 1984 and 1999, but a larger volume (39 million m ) was deposited, 3  3  resulting in 5.6 million m net aggradation. Adding the net changes from the intersurvey periods 3  (1952-1984 + 1984-1999) produces the same result as the direct difference from 1952 to 1999, consistent with the basic continuity equation for mass flux. In contrast, the sum of cut and fill volumes for the intersurvey periods do not equate to the direct difference from 1952 to 1999. This discrepancy occurs because compensating scour and fdl masks the full range of channel changes, the magnitude of which increases with the survey interval. Conversely, shorter survey intervals should display increasingly large volumes of sediment scour and fill (although the net difference would remain constant). This is demonstrated by the 1999-2003 surveys that show up to 22 million m eroded and 30 million m material deposited within the gravel reach over 4 years, a far greater 3  3  average annual rate than is observed between 1952 and 1999, since scour and fill areas remain largely distinct. These numbers illustrate the magnitude of sediment exchange that occurs as bars and islands are eroded and migrate downstream. The actual volume of 'new' material that is transported into the reach (i.e. net volume) is, however, comparatively modest. Table 5-1. Gross volumetric changes between Mission and Agassiz Cut Volume (m )  Fill Volume (m )  Net Volume (m )  1952-1984  53,150,000  48,809,000  - 4,341,000  1984-1999  32,991,000  38,598,000  + 5,607,000  1952-1999  56,277,000  57,543,000  + 1,266,000  1999-2003  22,289,000  29,933,000  + 7,644,000  Survey Period  3  3  3  Any discussion of total sediment flux is complicated by both the entrainment and deposition of fine sediments on islands and floodplain, and by the suspended fraction of the sand load. The movement of this material inflates the magnitude of apparent morphologic change in the reach. By removing these components of the total load according to the conventions described -112-  earlier, the bed material fraction can be isolated. A summary of bed material exchanges in the gravel reach is given in Table 5-2. All values correspond to the influx past Agassiz-Rosedale Bridge. Reported totals including values for individual 1-km reaches are given in the sediment budget tables. Table 5-2. Summary of bed material sediment budget calculations. AH values are bulk in millions of m . Period 1952- 1984  1984-1999  1952- 1999  1999 - 2003  Total  Sediment  Channel bed  Channel banks  Mining  gravel  4.918  -0.999  1.554  5.473  sand  2.208  -3.519  0.463  -0.848  Total  7.126  -4.518  2.017  4.625  gravel  3.576  -0.118  0.951  4.409  sand  1.079  -0.005  0.288  1.362  Total  4.655  -0.123  1.239  5.771  gravel  2.831  3.075  2.505  8.412  sand  1.806  -1.723  0.751  0.834  Total  4.638  1.352  3.256  9.246  gravel  3.404  -0.473  0.135  3.067  sand  1.734  -0.585  0.038  1.188  Total  5.139  -1.058  0.174  4.254  The summary sediment budget shows that the net influx of gravel exceeds the input of sand in all periods, with most of the material deposited on the channel bed and bars. Channel banks represent a negative net flux (erosion exceeds deposition) for all periods except for the full 1952-99 intersurvey comparison. This occurs because banks and islands can be significantly eroded by a few floods depending on local channel alignment, but many decades are required to establish large areas of new floodplain. From 1952 to 1984, there was a net gain of 145,000 m /yr of bed material 3  to the gravel reach. The influx of gravels was actually larger, 171,000 m /yr, because of an 3  apparent sand deficit. From 1984 to 1999, the total bed material influx rate increased to 385,000 m /yr, 76% or 294,000 m /yr of which was gravel, and to 1,292,000 m /yr from 1999 to 2003 (59% 3  3  3  or 767,000 m /yr gravel). Average annual transport volumes are summarized in Table 5-3. The 3  bulk volumes can be converted to weight using a bulk density of 1770 kg m" ±110 m" . This 3  3  conversion is based on six 100 kg field measurements along the gravel reach (with a 2 standard deviation range).  - 113-  Table 5-3. Average annual transport rates (m / yr) for sand, gravel and total bed material load. Rates are equivalent to tthe influx past Agassiz-Rosedale Bridge. 3  Period  Gravel  Sand  Total bed material  1952- 1984  171,000  -27,000  145,000  1984- 1999  294,000  91,000  385,000  1999 - 2003  767,000  297,000  1,292,000  1952 - 1999 (Direct)  179,000  18,000  197,000  1952- 1999 (Sum)  210,000  11,000  221,000  The direct difference of the 1952 and 1999 surveys yields a total bed material influx rate of 197,000 m /yr, or 179,000 m /yr of gravel (91%). For comparison, the direct sum of the 3  3  intersurvey budgets gives 221,000 m /yr (total bed material) and 210,000 m /yr (gravel). These 3  3  figures are respectively larger by a factor of 1.12 and 1.17 times and represent a negative bias associated with the direct 1952-1999 budget. While reasonably close to unity, the sediment budgets should express mass continuity since the intersurvey budgets are simply intermediate subcomponents of the full computing period. Although sand certainly is transported past Mission, it is assumed that all gravels are 'trapped' and there should be no negative bias in the gravel budgets. Although there is a small possibility that a large volume of gravel actually does move past Mission which could account for the bias, the problem likely relates to uncertainty in the actual sand to gravel proportion downstream of Sumas Mountain. Since large volumes of mainly sand are mobilized in this region, small errors in the sand fraction could produce large errors in the gravel fraction. Amplifying this problem is the possibility that gravels may be prograding downstream over a former sand bed and the gravel estimates may be negatively biased as a consequence. As there are no bed samples available from 1952, and few since, there is insufficient information to fully resolve this problem. Upstream of Agassiz-Rosedale Bridge, direct difference of the 1952 and 1999 modeled surfaces reveals a modest net loss of bed material (3300 m /yr) because of a sand deficit. In 3  contrast, there has been a minor gain of 2400 m /yr of gravel (see Table 5-7). Adding these figures 3  to the transport rate at the bridge gives the bed material influx at Laidlaw, estimated as 193,000 m  3  (94% gravel). The positive net flux of sand and gravel along channel banks is consistent with the pattern observed downstream over this same period, but the channel bed and bars show substantial degradation. Although interpretation of these figures must be viewed cautiously because of limited  - 114-  spatial resolution associated with the 1952 survey data between the bridge and Laidlaw, a degrading channel zone is consistent with a lowering of the thalweg (Chapter 4). The summary sediment budget data imply that average annual transport rates have been increasing over time. In particular, the most recent (1999-2003) transport rate is considerably greater than during the longer term from 1952 to 1999 (Table 5-2). The total influx of bed material between 1984 and 1999 exceeds that recorded from 1952 to 1984 despite the shorter computing period and generally smaller flood flows. The estimate of gravel from 1952 to 1984 is, however, reasonably consistent with a value of 4.012 million m reported by McLean and Church (1999) 3  using a similar morphologic study approach (no estimate of sand was given). Further, those authors estimated the gravel influx past Agassiz to be 4.54 million m from Water Survey of Canada 3  measurements. Nevertheless, the limited spatial extent of the 1984 survey means that potential errors associated with the surface model could account for the discrepancy between the subperiods. For example, a modeled 1984 surface that was too low, on average, would inflate the 1984-1999 transport rate (while the 1952-1984 estimate would correspondingly decline). Although the 1999 to 2003 figure appears anomalously large, the topographic data upon which the estimates are based produces a more accurate surface of difference than is possible for the earlier intersurvey periods. The high influx of bed material derives in part from 1.5 kilometres of bank erosion (up to 130 m lateral shifting) on mid-Herrling Island across from Tranmer Bar and scouring of the adjacent bed. This eroded volume does not, however, account for the total sediment influx. Lateral erosion at this site has been on-going since the 1940s, when Tranmer Bar was initially established (Church and Ham, 2004). High rates of sediment influx past Agassiz-Rosedale Bridge are likely to continue for several more years until the island is bissected, whence rates are likely to decline significantly. Additional discussion of the progression of these events is given in Appendix C. Although the estimated influx rates at Agassiz-Rosedale Bridge for the different computing periods appear reasonable, the above discussion also acknowledges that summary sediment budgets may be sensitive to a variety of factors that influence reliability of the results. These include compensating scour and fill that masks the magnitude of channel changes when converting gross bed changes into component values of sand and gravel and uncertainty in the ratio of sand to gravel, especially in lower computing cells. These factors are related to, or are amplified by, the time interval between successive surveys, suggesting the potential for systematic time-dependent bias that negatively influences the results. This potential is reviewed below. - 115-  5.3.2 Time-integrated bias  A major assumption of the sediment budget approach is that there is no compensating scour and fill between successive surveys. While major sites of erosion and deposition remain distinct over short intersurvey periods, the possibility for undetected changes increases over time. For example, a cell may experience 2 metres aggradation from 1952 to 1984, followed by degradation of 2 metres to 1999. In such a case, there is no apparent change in elevation, hence storage, over the full computing period. While any sequence over the three dates that is persistent (erosion, deposition or no change) will not produce any bias, any succession from scour to fill will introduce this complication. The morphologic coding of individual 20 metre grid cells permits such circumstances to be identified. Over the period 1952-1999,23.2% of the total surface of difference area was found to exhibit non-persistent transitions. Nonetheless, continuity suggests that the gravel component of the bed material load should still be located within the gravel reach regardless of potential displacement between 1952 and 1999 (because sand may truly be transported past Mission, the same assertion does not necessarily hold). There are, however, transition sequences where mass continuity fails because wash material is discarded from the calculations. A simple case involves the transition island > bar > island which, for an individual grid cell, is equivalent to a net loss of 2.16 metres elevation (overbank erosion minus deposition thickness) over the intersurvey calculations, but the same cell is regarded as stable (no net change) between the endmember dates. Since many transitions involve either non-persistent bank erosion or deposition, the appearance of bias is introduced. A more detailed review of this problem is given in Appendix B. The sediment budgets also present an apparently remarkable circumstance whereby the transport rate of gravel decreases in proportion to the length of the intersurvey interval. The sand fraction (hence total bed material load) displays a similar pattern (Figure 5-5). This occurrence could be the outcome of some time-dependent bias, of which undetected transport (compensating scour and fill) is the most likely source. If this negative bias is real, then extrapolating back to an even shorter period than 4 years may be a necessary means for estimating an unbiased transport rate (i.e. when there is little or no compensating scour and fill). However, it is not known what this interval actually is. Using the equation from a best-fit logarithmic trend-line for this adjustment, the annual influx of gravel is estimated as 1,064,000 m per annum, a far greater influx rate than 3  is estimated from the intersurvey budgets or from direct measurements. There are, however, several reasons not to accept this result as being credible. First, the inverse relation between the transport rate and intersurvey interval may be spurious given that there - 116-  15  20  25  30  35  Years between surveys  Figure 5-5. Net apparent influx of sand and gravel as a function of the intersurvey period. A best-fit (logarithmic) trend line is given for gravel.  -117-  are only 4 data points upon which the analysis is based. Second, the 1999-2003 surface of difference shows a number of very distinct erosional sites and it is not likely that large volumes of intervening compensating scour and fill could have occurred within only 4 years. Finally, extrapolating the unbiased estimates for gravel transport over 47 years (assuming this material does not exit the study reach) yields a total gravel deposition of 50 million m between Mission and 3  Agassiz. This volume works out to an average deposited depth of 1.25 metres over the entire active channel zone, or 2.5 metres over the main 500 metre wide thalweg zone. The total volume of dredged bed material (3.3 million m ) is insufficient to compensate for this deposition (although 3  this is lower bound estimate given undocumented mining). As well, the actual survey data (either gross volumetric changes or extracted cross-sections - see Chapter 4) simply do not support this quantity of aggradation. Such a large quantity of deposition would also create morphologic changes including a transition to a dominantly braided channel, yet there is no evidence of this. Nevertheless, the possibility for unmeasured sediment flux remains and is further examined in the discussion of sediment budget errors. 5.3.3 D i s t r i b u t e d sediment budgets  Details on the distributed pattern or magnitude of erosion and deposition within the gravel reach are desirable for many direct management concerns, including where erosion might threaten local infrastructure or where deposition increases the local flood risk. Stored volume changes for the bed material load and for gravels are summarized on a cell-by-cell basis for the period 195299 in Figure 5-6. The volume changes summarize the constituent components of the sediment budget and display an alternating pattern of erosional and depositional lengths of channel. However, since displayed values are net changes, local details of complex scour and fill patterns within each computing cell are not revealed (refer to Figure 5-2). Eroding areas are typically associated with partly confined lengths of channel but are also common where the channel is shifting laterally, while depositing areas correspond to major bar/island complexes. Over shorter comparative intervals, the association between erosional and downstream depositional zones becomes more distinct and can be linked with changes to channel morphology (see Chapter 6). These relations are largely obscured over the full study interval. Upstream of the Agassiz-Rosedale Bridge, most computing cells degraded from 1952 to 1999, although there is significant storage 5 to 6 km upstream on Tranmer Bar. In comparison, downstream cells are dominantly aggradational with major storage zones near Greyell, Harrison, and Yaalstrick Islands. Most lateral instability occurs in these areas because deposited sediments -118-  2500  • sand + gravel removal • bank volume changes • channel volume changes  2000  o o o  x  m  BED  MATERIAL  1500  1000  E  CD  CT) c  500  CO CD  E _g  o >  CD  -500  -1000  -1500 1  5  9  13  17  21  25  29  33  37  41  45  49  53  41  45  49  53  computing cell 2500  2000  o o o  • gravel removal • bank volume changes • channel volume changes  GRAVEL  1500  1000  O)  c  CO  500  o <D  E _3  o >  CD  -500  -1000  -1500 13  17  21  25  29  33  37  57  computing cell Figure 5-6. Cell-by-cell volumetric changes of bed material (top) and gravel (bottom) computed for the period 1952 to 1999 by direct difference of surveys.  force the channel to flow around them to maintain conveyance. Upstream of Vedder River, net stored gravel volumes are a scaled (74-80%) fraction of the total bed material load. Below Vedder River, the rapidly increasing sand fraction begins to dominate channel changes and there is a major zone of degradation with a corresponding (albeit more modest) downstream zone of deposition at Matsqui bend. In contrast, there is a major zone of gravel degradation, but no downstream storage of this material. It is reasonable that the sand fraction in the computing cells immediately downstream of Vedder may be underestimated and the [apparent] degradation of gravel is actually a loss of sand. This circumstance could also occur if the assumption of zero gravel transport past Mission is violated. The magnitude of this potential problem can be better illustrated by assuming the sand and gravel input to a cell is equal to the output from the next cell upstream and extending calculations on a cell-by-cell basis. 5.3.4 Accumulated sediment budgets  By accumulating stored volumes of coarse sand and gravel along each 1-km computing cell upstream from Mission, transport rates can be estimated at different locations within the gravel reach. The accumulated transport plots more clearly delineate the major storage and transport zones along the gravel-bed reach. Figure 5-7 plots the higher transport rate of bed material and gravel for 1999-2003 on a separate scale than for the other periods so trends in transport are more clearly discernible than when all data are plotted at equivalent scaling. The 1952-99 plot shows that the transport rate increases between Laidlaw and mid-Herrling Island because the channel is degrading, while all plots show that the transport rate generally declines downstream to Mission as sand and gravels become deposited. In the 1952 to 1984 period, there are three major and 2 minor deposition zones that match reasonably in location (and magnitude) with those identified by McLean and Church (1999). From 1984 to 1999 the transport rate declines fairly steadily from Agassiz-Rosedale bridge to Queens Bar and a second major deposition zone extends downstream to Sumas Mountain. The direct 1952 to 1999 comparison shows a fairly consistent decline in bed material transport that extends from Agassiz-Rosedale to Vedder River though several distinct depositional zones can be identified. Downstream of Vedder River, the transport rate generally increases to Mission, reflecting a net loss of bed material sands. The 1999-2003 sediment budget extends upstream to lower Herrling Island and reveals a consistent decline in the transport rate to Chilliwack Mountain (cell 21) though there is only modest storage between Agassiz-Rosedale and Gill Island (cell 37) and downstream of Harrison River. Below Chilliwack Mountain, the transport rate increases to cell 18, then declines to zero at Mission, paralleling the pattern in 1952-84. - 120-  (5  500,000  1,500,000  400,000  1,200,000  300,000  900,000  200,000  600,000  100,000  300,000  -100,000  -200,000  -300,000 9  13  17  21  25  29  33  37  41  45  49  53  57  61  65  500,000  1,500,000  400,000  1,200,000  300,000  900,000  200,000  600,000  100,000  300,000  -100,000  -200,000  GRAVEL -300,000 1  5  9  13  17  21  25  29  33  37  41  45  49  53  57  61  65  computing cell  Figure 5-7. Downstream trends in bed material (top) and gravel (bottom) transport rates estimated from the sediment budget calculations. A falling plot (downstream direction) indicates a deposition zone, while a rising plot indicates an erosional zone.  With the exception of the most recent (1999-2003) results, the period gravel budgets present a troubling aspect whereby the calculated transport rate is negative for some computing cells (Figure 5-8). Assuming no gravel transport past Mission, negative transport implies that gravel is moving in an upstream direction, which is not physically reasonable. Since sands are mobilized past Mission, the bed material load may legitimately be negative if the magnitude of sand erosion exceeds that of deposition. Martin and Church (1995) encountered the problem of negative gravel transport on Vedder River and attributed the result to a violation of the zero downstream transport assumption. From 1952 to 1984, the rate of gravel 'loss' is equivalent to 10% of the gravel load entering the reach past Agassiz-Rosedale. The actual volume is relatively modest and could remain undetected if it remained along the channel bed downstream of Mission. However, the apparent losses in 1984-99 and 1952-99 represent respective losses of 64% and 25% of the total gravel influx. The magnitude of these volumes over many years should be manifested in the formation of gravelly barforms downstream, though none exist, at least on the surface. There exists the possibility that gravel lenses are transported along the bed in the active thalweg (where they can not be seen) or that gravels are buried by sand as the freshet wanes, but there is no simple means to sample the bed in this location to validate this supposition. A second troubling aspect of the budgets is that the large negative transport rate near Vedder River (1984-1999 and 1952-1999) corresponds with a large positive transport rate at the same location (1999-2003). The sediment budget tables reveal mainly large negative bed changes in downstream cells for the early computing periods, and mainly large positive bed changes in the most recent computing period. Such a circumstance could occur if the 1999 bed elevations for these cells were too low, on average. The 1999 soundings were collected over a period of 11 days, with soundings collected on the second day of the survey covering cells 11 through 20 - the cells with the potential bed error. Although there is no reason to suspect any systematic bias in the 1999 survey, some human or computing error introduced during data collection, or post-processing, may have affected the accuracy of those survey points. There are additional potential problems that could also lead to, or at least amplify, the magnitude of the negative transport rates. The major factors include unreported sand and gravel mining and uncertainty in the ratio of sand to gravel, especially among computing cells below Vedder River. The magnitude of each source of error, and an adjustment to correct for its effects, is examined below in a discussion of sediment budget errors.  - 122-  5.4 Sediment budget errors: Negative transport rates 5.4.1 V i o l a t i o n o f the zero transport assumption  Although previous attempts to calculate the sediment budget of lower Fraser River have assumed the gravel transport rate at Mission is zero, or at least negligible (McLean, 1990; McLean and Church, 1999; Church et al., 2001), sediment budget results presented herein show a potential conflict with this assertion. While there are no contemporary records, the WSC analysed particle size for 74 bedload and 124 bed material samples collected at Mission between 1966 and 1986. A brief review of these data shows clasts up to 32 mm collected in some samples, so it is reasonable that the percentage of smaller sized gravels moving along the bed may be significant. Analysis of all bed material samples reveals that, over a range of flows up to 13,000 m /s, a consistent 15% of the bed material sampled at Mission is greater than 2 mm. For flows above 5000 m /s, the data reveal that 10% of the bedload is greater than 2 mm. This is the preferred figure 3  for adjustment because bed material moves in intermittent suspension near Mission (McLean et al., 1999). The mean annual bedload transport rate at Mission (1966 -1986) was reported by McLean et al. (1999) to be 154,000 tonnes / year, so the adjustment to the budget becomes 15,400 tonnes or 8800 m per annum. This value, multiplied by the length of each survey interval, is added to the 3  first cell of each accumulated sediment budget. As the unaccounted sum is modest in relation to the size of the gravel and bed material deficits for the periods 1984-1999 and 1952-1999, this adjustment has a minimal impact on eliminating the negative transport rates, but it is sufficiently important to reduce the bed material deficit for 1952-1984 by 25% and the gravel deficit by half. To fully eliminate negative transport rates from 1952 to 1984,40% of the mean bedload transport measured at Mission would have to be added, but even this amount would not eliminate negative transport rates for the other periods. Consequently, there must be additional sources of error remaining. 5.4.2 S o u n d i n g errors  As reported, the possibility of data collection or processing errors for computing cells 11 through 20 may explain the anomalously large negative transport rates near Vedder River. A simple means to explore this issue is to ignore the 1999 survey and compute a sediment budget directly from the 1984 and 2003 surveys. This comparison yields a total bed material influx rate of 470,000 m /yr, or 350,000 m /yr of gravel (74%) - see Figure 5-8; Table 5-10. The transport rate 3  3  declines rapidly to Harrison Bar as sediment is deposited, while a second, smaller, sedimentation zone extends from Wellington Bar to lower Yaalstrick Bar. Negative transport rates — the - 123 -  600  500  -1001  0  5  10  15  20  25  30  35  40  45  Computing cell  Figure 5-8. Downstream trends in bed material and gravel transport rates, 1984 to 2003. The plots include the correction factor for bedload transport past Mission.  -124-  1  50  consequence of accumulated sand and gravel deficits on the bed within cells 11 through 20 — remain, but are strongly abated. The net (total bed material) change in these cells is -279,000 m , 3  compared with -1.11 million m calculated for the sum of the 1984-1999 and 1999-2003 3  intersurvey budgets. These findings indicate that an error in the 1999 survey remains likely, but additional sources of error must also be considered since the comparison does not satisfy mass continuity. These are reviewed in the following sections of the chapter. If the 1999 survey is in error, then it is necessary to estimate the magnitude of the error by increasing (adding sediment to) the 1984-1999 budget and decreasing (removing sediment from) the1999-2003 budget by some volume. Sediment must also be added to the 1952-1999 budget for consistency. This approach is considerably more efficient than changing the sounding elevations and recomputing topographic models of the bed. Since there is no obvious way to estimate this adjustment, an iterative procedure was developed whereby a fixed volume was added or subtracted from cells 11 to 20, and the effects were reviewed by examining the updated sediment budget plots (cf. Figure 5-7). To adjust the gravel budgets, the fixed volume was multiplied by 0.75 to correct for the sand fraction. The largest accumulated negative transport rate for all budgets is roughly 4.5 million m between 1984 and 1999. Distributing this volume over the 10 cells (11 to 20) reduces the magnitude of the annual bed material transport rate from -300,000 to -87,000 m /yr, and gravel 3  transport from -220,000 to -63,000 m /yr. The running sum remains negative because there was 3  more sand and gravel erosion than deposition in the first 10 computing cells. The adjustment also eliminates negative transport rates for both bed material and gravel for the 1952-1999 computing period. However, subtracting the same volume from cells 11 to 20 for the 1999-2003 budget creates such a large sediment 'deficit' that the accumulated transport rate for bed material and gravel attains a maximum of -564,000 and -520,000 m /yr respectively in cell 22, and the transport 3  rate remains negative for all cells upstream to Agassiz-Rosedale Bridge. This adjustment is obviously too great, so was incrementally reduced by 50,000 m to 150,000 m per cell for the bed 3  3  material budgets (and correspondingly, 112,500 m for gravel budgets). Dividing 150,000 m by 3  3  the average cell area yields an estimated elevation adjustment of 15 cm. This value provides the best compromise between reducing the size of the negative transport rates and maintaining positive transport rates for the 1999-2003 computing period. The budget modifications also result in a minor reduction of the negative bias (from 1.12 to 1.10 times) between the direct sum of the intersurvey budgets and the direct difference of the 1952 and 1999 surveys. Updated transport rates - 125-  for sand and gravel are summarized in Table 5-4, and accumulated transport plots are shown in Figure 5-9.  Table 5-4. Average annual transport rates (m /yr) for sand, gravel and total bed material load incorporating adjustments for bedload transport at Mission and possible sounding errors in cells 11-20. Rates are equivalent to the influx past Agassiz-Rosedale Bridge. 3  Period  Gravel  Sand  Total bed material  1952 - 1984  180,000  -27,000  153,000 ± 16,000  1984- 1999  378,000  116,000  494,000 ± 33,000  1999 - 2003  494,000  203,000  697,000 + 41,000  1952 - 1999 (Direct)  212,000  25,000  237,000 ± 11,000  1952- 1999 (Sum)  243,000  19,000  262,000 ± 15,000  The persistence of negative transport rates in lower computing cells remains a concern, especially from 1984 to 1999. The possibility of systematic sampling bias in 1984 must also be considered. Cells 10 to 17 were surveyed by a UBC crew using a much smaller, lighter boat than was used elsewhere along the main channel in 1984 (McLean, 1990) or in 1999 and 2003. This circumstance could yield measurements of the bed elevation that were too high on average relative to the mapping datum. Terrestrial surveying was additionally used near the mouth of Sumas River (cell 15) which is coincident with the largest sand and gravel deficit of any single computing cell. Although there is no available means to confirm that topographic points measured in 1984 may have been consistently biased, the precision of individual points was certainly lower relative to more recent measurements given contemporary survey techniques and available ground control. Indeed, McLean (1990, Table A3) shows that vertical closure errors in survey traverses were much larger between Sumas Mtn. and Chilliwack Mtn. than for the remainder of the channel. Given the very real likelihood that the large bed material deficit in this region from 1984 to 1999 is merely an artefact of survey imprecision, the local transport rate must be considered unreliable relative to that estimated between Agassiz-Rosedale and Chilliwack Mtn. Likely errors associated with the 1999 survey in this same region only compound the uncertainty, but there are also remaining possible sources of error. 5.4.3 Unreported dredging  An additional explanation for negative transport rates is that dredged gravel volumes have been underestimated. Removal volumes are known to be negatively biased before 1974 because of undocumented removals (Chapter 2). Adding an additional 800,000 cubic metres of unreported - 126-  Figure 5-9. Downstream trends in bed material (top) and gravel (bottom) transport rates estimated from the adjusted sediment budget calculations. A falling plot indicates a deposition zone, while a rising plot indicates an erosional zone. Note the change of ordinate scaling.  gravel to the 1952-84 budget would eliminate the negative transport rate. This volume is comparable with that removed from a single removal site near lower Minto channel, hence appears credible. A similar adjustment would also correct the 1952-99 budget, but the required volume (2 million m ) is considerably larger and more difficult to justify given available evidence. The 3  remaining deficit is sufficiently large from 1984 to 1999 (3.4 million m ) that unreported volumes 3  of this magnitude are not believable given that permitting requirements were in place during this period. Therefore, an alternative consideration that could negatively bias the transport rate must also be examined. 5.4.4 S a n d a n d gravel fractions  Perhaps the main factor that could produce negative transport rates is uncertainty in the actual sand and gravel fractions for individual computing cells, since this could inflate the volume change of sand or gravel, while suppressing the other. From 1952 to 1984, there is a loss of 1.8 million m of gravel in cells 19-22 (near modern Wellington Bar) that likely accumulated in 3  downstream cells (mainly 14 to 18). If this was assumed the true zero boundary of gravel transport — given that the proportion of sand in the deposits become increasingly large downstream — the accumulated transport rate upstream would remain positive and the influx at Agassiz-Rosedale would increase to 186,000 m /yr. While some modest amount of gravel probably is deposited below cell 14, the accumulation of a large gravel deficit in these cells produces a negative upstream transport rate. If most of this deficit actually reflected the erosion of bed material sands, the gravel budget would again remain positive. Conversely, an underestimate of the gravel fraction in cells which show apparent large stored volumes of sand would also correct budget problems. As discussed earlier, there has been minimal bed material sampling in the lower end of the reach because there are few emergent bars and no islands so the assigned sand / gravel fractions are imprecise. Similar arguments obviously apply to the 1984-99 and 1952-99 budgets, but problem cells are more difficulty to identify. In 1952-99, there is a gravel deficit of 700,000 m in cells 15 to 17, 3  yet no obvious downstream storage of this material since there is a further sand and gravel deficit of 2.2 million m in cells 9 to 14 with no corresponding zone of accumulation. A combination of 3  uncertainty in the actual sand to gravel fraction and unrecorded gravel removals may explain most of the deficit. There is a further possibility that gravels have been prograding downstream over a former sand bed which would result in a negative bias of the gravel estimates. However, there are no bed samples for any date near 1952 that could resolve this question. From 1984 to 1999, a very - 128-  large (3.8 million m ) deficit of both sands and gravel in cells 7 to 17 produces correspondingly large negative transport rates. While these cells are consistent in location with the transition to a sand-dominated bed, the possibility'of sounding errors remains a significant concern.  5.5 Sediment budget errors: uncertainty of the estimates 5.5.1 Introduction The development of the sediment budget approach is also subject to a number of errors and assumptions that produce uncertainty in the estimates. The main source of error relates to imprecision associated with the channel bed and banks modeling from available bathymetry and topographic data, since all following computations are based upon manipulations of these modeled surfaces. Errors in individual models can result in apparent erosion or deposition when none has actually occurred, inflating the magnitude of actual change. The budgeting approach to parse gross volumetric changes into estimates of sand, gravel and wash material introduces additional error since there is insufficient field sampling to confirm the actual proportions of these sediments as they vary along the gravel reach. While field samples collected in 1983 and 2000 show no significant difference in the sand fraction, there is considerable local scatter that could affect the precision by which the sand / gravel ratio is proportioned for any 1-km cell. The possibility that these ratios change over time must be recognized, but there is no physical means to consider this problem. The following discussion of different sources of error provides estimates of the error margins associated with the budgets, hence the probable range of sand and gravel transport influx for each computing period.  5.5.2 Survey error Surveying and mapping errors limit the accuracy by which individual modeled surfaces can be represented (Chandler, 1999) and can be introduced at any stage from initial observation to presentation of results (Burrough and McDonnell, 1998). Lane et al. (1994) point out that individual points may be subject to random, gross, or systematic errors. Random errors refer to the precision within which repeat measurements yield the same result, gross errors relate to equipment or operator blunders and systematic errors (bias) refer to a persistent measurement problem, such as an incorrect datum. Since conventional approaches such as ground truthing for checking the quality of individual points are neither practical over such a large area nor feasible for historic data, it is difficult to quantify these errors. However, the vertical precision of the 1999 and 2003 sounding data is reported by PWGSC to be roughly ±10 cm, reflecting the accuracy with which the sounder can measure the bed depth relative to the water surface, and to relate the water surface - 129-  elevation to the geodetic datum. Horizontal positioning errors are larger (up to 2 metres) and are caused by imprecision in the GPS measurements. The altimetry data collected in 1999 by Terra Surveys Ltd. are reported to have the same vertical and horizontal precision as the sounding data. For comparison, McLean (1990) estimated maximum vertical errors for the 1984 survey at several tens of centimetres, while horizontal errors ranged from 2 to 5 metres at most locations. The precision of the 1952 soundings is unknown, but the error of photogrammetrically derived contour elevations is conventionally estimated as half the contour interval, or ±75 cm. Photogrammetric points appended to both the 1952 and 1984 surveys have similar errors, estimated as ±50 cm based on stereoplotter parameters and operator experience. While the potential error associated with the surveys is substantial, it is reasonable to assume that the errors are random. Therefore, the effect of individual points is reduced as these errors become compensating or are averaged out by the Topogrid modeling procedure (Chapter 3). The possibility for gross errors was examined by producing TIN models from the bathymetry data and visually checking for anomalous values. The utility of this exercise is limited for the 1984 and 1999 datasets, but several digitizing and coding blunders were discovered, then corrected, for the 1952 dataset. Additional discussion of survey errors is pursued indirectly by examining the accuracy of the individual surface models into which the survey errors are incorporated, which may introduce bias into the final sediment budget results. Several different approaches for estimating modeling error are considered, following results presented in Ham and Church (2003). 5.5.3 Floodplain error  An empirical test of modeling bias was made by comparing interpolation errors on stable floodplain areas between survey dates. The surface models extend beyond the active channel zone (which includes channel scour and fill, bank erosion and bank deposition) to the floodplain. If it is assumed that the surface elevation of these areas has remained constant (although some amount of compaction or overbank deposition could technically occur) then any computed volumetric differences on these areas should correspond to model interpolation error. Extending this error over the entire active channel zone can therefore provide an indication of the volumetric error included in the sediment budget calculations. These errors are calculated by dividing computed volumetric changes over the stable floodplain by the corresponding area over which the measurements are made. These differences are -3 cm from 1952 to 1984, +12 cm from 1984 to 1999 and +3 cm between 1952 and 1999 (positive numbers indicate net aggradation). The 1999-2004 surface of difference shares common 2004 Lidar data on much of the stable floodplain so the floodplain error - 130-  is reduced overall to -0.4 cm. Most of the error occurs along modeled surface boundaries or where 1999 Lidar data overlap — elevation differences between the 2 floodplain surfaces is otherwise typically less than 1 mm (the maximum resolution of the models). The 1984-99 errors are significantly larger than for the remaining computing periods, but most of the error is limited to three adjacent computing cells. These cells are characterized by a paucity of elevation data near the edge of the floodplain and the channel bed, so model interpolation is poor in this region. Eliminating these three cells decreases the average floodplain error to 4 cm. An approximation of the potential effect of these errors on the sediment budget calculations can be made by pooling the volumetric precision in bed elevation changes for individual computing cells. Each cell has an average area of the order 1,000,000 m , so the corresponding 2  volumetric precision is roughly 30,000 - 40,000 m . From 1952 to 1984 and 1952 to 1999, the 3  pooled error is estimated as 197,000 m , rising to 262,000 m from 1984 to 1999. These figures 3  3  represent roughly 5% of a sediment budget. The pooled error from 1999 to 2003 is reduced by a factor of lOx from the 1952-84 and 1952-99 periods. In general, floodplain modeling errors are probably smaller than bed modeling errors because a relatively flat surface can be modeled by a GIS with greater fidelity than a more complex topography. This is especially relevant for the 19992003 comparative period since the spatial density of the 2004 LiDAR data (even thinned from 1 m to 5 m) produces a much more accurate model of the floodplain than is possible for any previous comparative period, and is certainly more accurate than bed topography interpolated from spaced sounding lines. The reliability of individual bed models may therefore provide a better indication of survey error. 5.5.4 Channel Error The comparison of different interpolation schemes in Chapter 3 provides a direct means to empirically examine modeling bias in the channel. It was demonstrated that the adopted Topogrid model yielded the smallest net difference between a reference surface and one derived from a thinned subset of the data used to create the reference surface. The calculated difference was 70,000 m over a 7-km long test reach between Queen's Bar and the mouth of Harrison River 3  (indicating an erosional bias). This test reach was selected because it is morphologically complex, hence should provide an indication of the maximum expected modeling error. By pooling this area by the number (9) of morphologic reaches between Mission and Agassiz, a total pooled error of 210,000 m is estimated for an individual survey. Between successive surveys, the combined error 3  becomes 297,000 m , similar to the pooled floodplain error. These results show that the volume 3  - 131 -  bias associated with the selected modeling procedure is actually not large relative to the volume of sediment that is mobilized within the gravel reach, hence does not appear to be a major factor limiting the reliability of the sediment budgets. A further analysis of channel modeling error can be pursued by validating the precision of the actual channel models, or the difference in elevation between the modeled surface elevations and the survey elevations at corresponding locations. The bathymetry, altimetry and photogrammetric spot elevations were overlaid on the grids interpolated from these same data and the root mean square error was calculated for each survey. This measure actually quantifies accuracy in addition to precision unless the mean error of the data is zero (Lane et al., 2003). However, to calculate accuracy independently requires comparison to a true model, but none exists. The magnitude of RMS errors is partly dependent on the size of the modeled grid cells. As cell size increases, the possibility for multiple survey elevations to overlap the cell also increases. Consequently, large RMS values can occur where these points have disparate elevations, since the grid cell represents a single, spatially averaged elevation. Lowering cell size reduces this possibility, and also increases the total number of computing cells, further suppressing RMS values. However, cell sizes smaller than the average spatial density of the survey data can not reliably reproduce any real variability in channel morphologic complexity at the same spatial resolution. The RMS error is calculated as ±0.97 m in 1952, ±0.85 m in 1984, ±0.79 m in 1999 for the 20 m grids. Errors are reduced to ±0.63 m in 1999 and ±0.56 m in 2003 using 10 m grids. The reduction in error occurs despite the use of original (non-thinned) bathymetry data so individual grid cells may spatially overlap more than 20 bathymetry points. The estimated precision between surveys is E  = (Ej + E 2 2  d i f f  2  ) ' 0  5  where E is the RMS of the successive surveys. This produces an  estimate of ±1.29 m from 1952 to 1984, ±1.16 m from 1984 to 1999, ±1.25 m from 1952 to 1999 and ±0.84 m from 1999 to 2003. The RMS error of the mean difference in elevation is given as E ;ff / ( n ) d  05  where n is the number of independently determined points in each surface of  difference. Since interpolated cell elevations are estimated from the value of surrounding points centered around each grid cell, n must be divided by 9 (the assumed number of points in a computing array). From 1952 to 1984, there are 13,410 independent cells, yielding a pooled error of ± 517,000 m , or ± 4.1 cm over the entire reach. This result is slightly better than the RMS error 3  of ± 4.7 cm estimated by McLean (1990) and probably reflects the inclusion of bed and bank contours in the modeling procedure (see Chapter 3). Corresponding values for the period 1984 to - 132-  1999 are ± 491,000 m (± 3.9 cm), ± 524,000 m (± 4.1 cm) from 1952 to 1999 and ± 163,000 m 3  3  3  (± 1.3 cm) from 1999 to 2003. These values are also reported as annual values in Table 5-4. Although these errors are small relative to the total bed material influx for each period (5 to 10%) it is recognized that not all potential sources of error have been considered. Consequently, the total actual uncertainty is almost certainly larger than these reported values, but there is no obvious (and perhaps possible) manner by which all error terms can be calculated and propagated.  5.6 Comparison with other studies Several previous attempts have been made to establish bedload influx to the lower gravel reach. Most of these efforts have been associated with the program of direct sediment sampling undertaken by the Water Survey of Canada from 1965 to 1987 in support of various river management concerns. Initial results were presented by Tywoniuk (1973) who found that coarse sediments aggraded in the Mission to Agassiz reach, but there was a much greater loss of sand and finer sediments. Tamburi (1979) made no attempt to compute annual loads for bedload transport, but demonstrated a highly non-linear relation between discharge and transport for material coarser than 6 mm at Agassiz. The first comprehensive review of the entire sampling program was given by McLean and Church (1986), while McLean and Tassone (1987) examined the precision of the WSC sampling procedures. McLean et al. (1999) also provided an overall assessment of the program in an attempt to make the results more accessible. The authors reported that 62 bedload measurements were collected at the Agassiz-Rosedale Bridge during freshet from 1968 to 1976. Adjustments were made to correct for sampler inefficiency and the proportion of fine (<6 mm) gravel that was not trapped, while annual loads were estimated from a rating curve based on daily measurements. McLean et al. (1999) estimated the average annual influx of bedload past Agassiz to be 227,000 tonnes per year, or 130,000 m (adopting a gravel unit weight of 1.75 t/m ). The 3  3  bedload and bed material load can be considered equivalent at Agassiz, since few particles < 2 mm are trapped within the bed and these move in suspension once entrained (McLean et al., 1999). Church (2001a) related estimated annual loads (reported in McLean et al., 1999) to the annual maximum daily flow, assumed to be the most relevant flow index given the highly nonlinear nature of sediment transport on Fraser River. The results yielded a smaller standard error of estimate because the annual rating curve exhibits less scatter. Church's back-transformed equation is given as: - 133-  G = 2.231xl0- Q 19  a  where G is given in tonnes per year, and Q a  m a  x'  s  6037 max  given in m/s, measured at Hope gauge. Using  this relation, annual loads at Agassiz can conveniently be calculated for different periods, and thus can be compared with the sediment budget estimates (see Table 5-5).  Table 5-5. Comparison of W S C transport estimates with sediment budget estimates for gravel. The budget estimates include the adjustment for bedload transport at Mission and the possible sounding errors in cells 11-20. Period  W S C (total m )  W S C ( m / yr)  Budget (total m )  Budget (m / yr)  4,334,000  134,000  5,473,000  180,000  1,668,000  97,000  4,409,000  378,000  6,002,000  123,000  8,412,000  212,000  670,000  168,000  3,067,000  494,000  2,338,000  186,000  6,604,000  350,000  3  1952 - 1984  1  1984- 1999  1  1952- 1999  1  1999 - 2003 1984 - 2003  2  3  1. Values reported in Church et al., (2001) using a unit weight of 1.75 t / m  3  3  3  2. Calculated from Hope flow records using relation given in Church (2002) and unit weight of 1.75 t / m  3  These results show that sediment budget estimates of gravel transport are consistently greater than the sediment transport estimates, although there is reasonably good correspondence from 1952 to 1984. Since sediment traps are known to underestimate total load because of sampling difficulties and trap inefficiencies, McLean et al. (1999) adjusted annual loads based on sampler calibrations and assigned an error bound of 40% around the annual estimates. Church et al. (2001) calculated pooled errors from these annual estimates to be 1.77 million m over the 1952-1999 period. Incorporating this result, the upper bound estimate from the sediment transport data becomes 7.77 million m , which is close to the sediment budget estimate for the same period and almost certainly 3  overlaps the lower bound estimate of the sediment budget error margins (these margins are presented in the following section of this chapter). However, the annual transport rate does not necessarily correlate with the size of the peak annual flow since the river is often transporting less material than its hydraulic capacity. Consequently, the transport rate can vary considerably for a given flow (McLean and Church, 1999). Certain channel configurations may therefore yield large volumes of available sediment and produce correspondingly high transport rates, at least for short durations. The present channel alignment whereby Tranmer Bar forces the current directly towards an eroding Herrling Island provides an example of this, and appears to explain the relatively high,  -134-  recent (1999-2003) transport rate. Regardless, the 1984 to 1999 sediment budget estimate remains suspiciously large. The most comprehensive attempt to estimate sediment transport on the river was by McLean (1990) who contrasted the WSC measurements with a number of bedload formulae and morphologic estimates based on the 1952 and 1984 bathymetric surveys. Estimates for gravel influx at Agassiz-Rosedale were found to vary from a low of 70,000 m /yr using the Meyer-Peter and Muller equation, to a high of 150,000 m /yr following Neill's approach based on bank erosion. 3  Some of the differences were attributed to the unequal time periods considered by each method. McLean and Church (1999) presented a re-examination of these methods and concluded that the long-term average annual inflow of gravel is 200,000 t/yr, or 114,000 m /yr (equivalent to the bulk 3  volumetric rate of gravel input from the sediment budget). The figure of 180,000 m /yr presented 3  herein is not unreasonable given this range of estimates. Although this number is larger than that given by McLean and Church (1999) using the same bathymetric data, the discrepancy can be attributed to several factors incorporated into this study, including a less restrictive computing mask, additional topographic data on bar surfaces, a larger adjustment for eroded overbank deposits and a different modeling algorithm. A distributed gravel budget has also been computed by Li and Millar (2004) based on a depth-averaged 2-D numerical model. Five runs at different discharges were executed until equilibrium transport conditions were achieved. An annual transport rate based on 1999 river topography was calculated by integrating a discharge-transport rating curve with the flow duration histogram at boundaries that discriminate aggradational and degradational zones. The analysis yielded an estimate of 206,000 m /yr at Agassiz-Rosedale, which the authors suggested would be 3  most directly comparable with the 1999-2003 morphologic budget estimates because of similar channel configuration. The underestimate of actual transport over this same period (compared to the sediment budget) occurs because the transport rate was effectively calculated as a function of discharge. However, the transport rate does not necessarily correlate with discharge, and is largely dependent on sediment availability, hence estimates of bed material transport may be accordingly over- or under-estimated during computing periods of several years or less. It is expected that these aberrations are averaged over periods of decades.  Given the issue of short-term variability, and uncertainty in the accuracy of the 1984 survey, the long-term sum (1952-1984 + 1984-1999) rate of 243,000 m /yr appears to be the best, 3  - 135-  or at least most conservative, estimate of gravel influx for river management purposes. Construction of period sediment budgets is a complex task, and the numbers should not be considered as final. There are a variety of assumptions and errors that affect the reliability of the estimates, and there is no means at present to fully resolve all of the issues presented in this chapter. It does appear, however, that sediment influx has been much larger since 1984 than before that date, with corresponding implications for gravel and flood management. Additional surveys at fairly regular intervals will be required in the future to appraise this observation.  - 136-  Table 5-6: Sediment budget -1952 to 1984  Cell 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 R1-43I  bed area (m ) 2  channel change -86,440 -195,763 -142,856 286,127  105,543 483,895 437,985 372,880 544,760 629,784 721,734 800,881 631,219 760,305 722,515 607,124 722,274 545,692 856,885 982,603 996,985 1,217,050 1,344,434 1,424,141 1,140,624 1,146,263 1,248,533 1,528,863 1,221,580 1,036,854 1,174,899 1,187,079 1,119,470 1,498,734 1,463,718 1,117,674 844,961 643,393 665,710 1,428,674 1,127,393 1,118,303 931,040 1,037,152 797,303 546,391 537,182  -338,891 -415,525 751,176 869,463 273,600 345,853 1,145,592 -353,017 -569,560 -364,005 -607,697 -27,135 794,661 654,204 403,807 -166,094 -61,051 126,516 565,829 -293,513 -97,404  39,470,480!  7,126,250  445,859 419,511 377,450 372,310 -584,136 -324,778 121,628  218,083 171,437 -242,543 -354,184 -435,767 1,260,010 346,121 1,404,870 417,917 964,026 50,556  bed changes % channel gravel sand 95 95 95 95 95 95 40 30 27 27 26 26 26 26 26 26 25 25 25 25 25 24. 24 24 24 24 24 23 23 23 23 23 23 22 22 22 22 22 22 21 21 21 2.1  -4,322 -9,788 -7,143 14,306 22,293 20,976 226,470 260,617 -428,104 -238,580 89,556 -250,106 -307,374 556,948 646,137 203,793 258,202 857,218 -264,757 -428,135 -274,243 -458,881 -20,536 602,778 497,356 307,683 -126,840 -46,727 97,049 435,006 -226,153 -75,217 168,780 132,973 -188,540 -275,929 -340,231 985,926 271,423 1,104,080 329,154 760,920 39,991 4,918,0251  channel sand -82,118 -185,975 -135,713 271,820 423,566 398,535 150,980 111,693 -156,031 -86,198 32,073 -88,784 -108,151 194,228 223,326 69,808 87,651 288,374 -88,259 -141,425 -89,762 -148,816 -6,599 191,883 156,849 96,124 -39,254 -14,324 29,468 130,822 -67,360 -22,187 49,303 38,464 -54,003 -78,255 -95,535 274,084 74,698 300,790 88,764 203,106 10,565 2,208,224  erosion (sub 3m) 0 0 0 0 -293,103 -64,199 0 0 0 0 -62,691 -64,738 0 -75,744 -308,990 -48,014 -100,725 -270,619 -95,730 0 -486,929 -114,689 0 -34,266 -10,342 0 0 0 0 0 -110,569 -102,583 -267,230 0 -156,736 0 0 0 0 0 -17,588 -417,213 0 -3,102,6991  deposition (sub 0.84m) 48 1,792 0 75 47,011 1,187 0 0 0 0 0 0 134,272 22,793 3,618 0 0 0 0 60,385 38,258 0 98,790 168,602 7,954 146 54,609 60,500 43,566 113,526 0 48,247 45,469 92,866 0 0 0 0 289,629 102,922 34,643 18,865 0  bank changes bank total 48 1,792 0 75 -246,092 -63,012 0 0 0 . 0 -62,691 -64,738 134,272 -52,952 -305,372 -48,014 -100,725 -270,619 -95,730 60,385 -448,671 -114,689 98,790 134,336 -2,389 146 54,609 60,500 43,566 113,526 -110,569 -54,336 -221,762 92,866 -156,736 0 0 0 289,629 102,922 1.7,055 -398,348 0  1,489,7751 -1,612,924|  bank gravel  bank sand  2 90 0 4 -12,305 -3,151 0 0 0 0 -46,160 -47,778 99,325 -39,260 -226,936 -35,763 -75,197 -202,498 -71,796 45,391 -338,031 -86,604 74,767 101,898 -1,816 112 41,703 46,305 33,419 87,278 -85,194 -41,959 -171,627 72,030 -121,838 0 0 0 227,123 80,886 13,433 -314,422 0  46 1,702 0 72 -233,787 -59,861 0 0 0 0 -16,531 -16,961 34,948 -13,691 -78,436 -12,250 -25,527 -68,122 -23,934 14,994 -110,640 -28,086 24,023 32,437 -573 35 12,906 14,195 10,147 26,248 -25,375 -12,377 -50,135 20,836 -34,898 0 0 0 62,506 22,036 3,622 -83,926 0  -998,567|  -614,357  ;  O/B sand (>0.177 mm) -526 -816 467 -318 -50,632 -46,883 -9,435 -21,098 -23,659 -33,927 -20,948 35,131 93,377 -40,806 -88,685 -40,160 -113,206 -133,464 -128,480 -92,360 -108,529 -127,688 29,460 -81,406 -25,353 -62,021 13,472 6,195 42,086 4,281 -29,857 -377,856 -258,049 -6,732 -53,061 -48,563 -204,953 -393,468 -51,488 -60,991 -139,679 -193,422 -60,585  gravel removal 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 597 0 1,725 4,585 57,259 0 0 0 0 0 0 499,790 415,453 35,365 0 0 0 237,344 0 0 0 173,709 72,145 2,986 0 0 52,999  2,904,634| I 1,553,9581  sand removal 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 203 0 575 1,515 18,741 0 0 0 0 0 0 153,210 126,147 10,635 0 0 0 68,656 0 0 0 48,291 19,855 814 0 0 14,001  462,642  gravel sum (m ) 3  -4,320 -9,699 -7,143 14,310 9,988 17,825 226,470 260,617 -428,104 -238,580 43,396 -297,884 -208,049 517,688 419,202 168,029 183,602 654,720 -334,828 -378,159 -555,015 -545,485 54,231 704,677 495,539 307,794 -85,137 499,368 545,921 557,649 -311,347  sand sum (m ) 3  -82,599 -185,088 -135,246 271,574 139,148 291,791 141,545 90,595 -179,691 -120,124 -5,407 -70,613 20,174  -117,176 -2,847 442,347 -310,377 -275,929 -340,231 1,159,636 570,691 1,187,953 342,587 446,498 92,989  139,730 56,205 17,398 -50,879 86,789 -240,098 -217,276 -290,190 -304,589 46,885 142,914 130,923 34,138 -12,876 159,276 207,848 171,986 -122,592 -412,420 -258,881 121,224 -141,962 -126,818 -300,488 -71,094 105,571 262,648 -47,293 -74,243 -36,019  5,473,4161  -848,1241  total s+g (m ) 3  -86,918 -194,787 -142,389 285,884 149,136 309,616 368,015 351,212 -607,795 -358,704 37,989 -368,498 -187,875 657,419 475,407 185,427 132,722 741,509 -574,926 -595,435 -845,205 -850,074 101,116 847,591 626,462 341,932 -98,013 658,644 753,769 729,636 -433,939 -529,597 -261,728 563,571 -452,339 -402,747 -640,719 1,088,542 676,262 1,450,601 295,293 372,256 56,971  Cell 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43  4,625,292| R1-43  -137-  T a b l e 5-7: S e d i m e n t b u d g e t - 1 9 8 4 to 1999  area Cell 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 R1-43  (m ) 2  channel change  106,741 494,455 447,283 375,795 604,410 677,044 747,773 829,248 658,920 806,125 740,098 655,217 715,199 589,780 954,076 1,040,531 1,117,195 1,329,174 1,297,617 1,586,769 1,289,589 1,256,391 1,173,212 1,344,348 1,163,611 957,233 1,040,860 965,073 731,208 1,037,948 973,838 1,285,608 1,084,795 648,575 614,565 1,187,823 1,215,904 1,380,573 1,299,809 1,033,797 871,156 631,587 594,330  -48,635 82,881 -328 5,144 52,096 124,390 -105,554 -314,575 -69,325 -315,853 -449,745 -188,703 -395,300 -656,839 -1,222,480 -35,765 -481,197 294,607 451,885 228,709 655,965 726,596 450,744 363,750 -6,992 -377,085 -118,270 -60,949 4,898 -50,577 1,349,104 1,012,155 370,701 12,637 587,388 386,391 907,498 432,158 663,833 176,544 317,871 28,691 -133,434  39,555,284  4,655,031  bed changes % channel sand gravel 95 95 95 95 95 95 40 30 27 27 26 26 26 26 26 26 25 25 25 25 25 24 24 24 24 24 24 23 23 23 23 23 23 22 22 22 22 22 22 21 21 21 21  channel sand  erosion (sub 3m)  deposition (sub 0.84m)  bank changes bank total  -2,432 4,144 -16 257 2,605 6,220 -63,333 -220,202 -50,807 -232,024 -331,150 -139,266 -292,413 -487,003 -908,480 -26,640 -359,245 220,447 338,907 171,920 494,207 548,664 341,135 275,917 -5,315 -287,322 -90,319 -46,649 3,757 -38,884 1,039,491 781,601 286,895 9,802 456,603 301,020 708,543 338,153 520,568 138,745 250,357 22,646 -105,549  -46,203 78,736 -312 4,886 49,491 118,171 -42,222 -94,372 -18,518 -83,829 -118,596 -49,437 -102,887 -169,836 -314,000 -9,125 -121,952 74,160 112,978 56,790 161,758 177,933 109,610 87,833 -1,676 -89,763 -27,951 -14,300 1,141 -11,694 309,613 230,554 83,806 2,835 130,785 85,371 198,955 94,005 143,265 37,799 67,514 6,045 -27,884  0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -394,081 -60,835 0 0 0 -1,241 -8,051 0 0 0 0 0 0 0 -41,540 0 0 0 0 -392,494 0 0 0 -3,839 0 0 0 0  0 0 0 0 19 216 0 1,275 628 0 0 0 177 0 0 0 108 30,596 225,284 166 23,015 6,910 72,693 0 0 0 0 0 0 0 228,058 0 0 0 104,001 50,283 0 0 0 0 0 0 1,494  0 0 0 0 19 216 0 1,275 628 0 0 0 177 0 0 -394,081 -60,726 30,596 225,284 166 21,774 -1,140 72,693 0 0 0 0 0 0 -41,540 228,058 0 0 0 -288,493 50,283 0 0 -3,839 0 0 0 1,494  3,575,554  1,079,477  -902,081  744,922  -157,159  bank sand  bank gravel 0 0 0 0 ' 1 11 " 0 892 460 . 0 • ~0 0 131 0 0 . -293,533 -45,336 22,894 168,959 125 16,405 -861 55,016 0 0 0 0 0 0 -31,936 175,720 0 0 0 -224,259 39,173 0 0 -3,010 0 0 0 1,182 -117,967  O/B sand (>0.177 mm)  0 0 0 0 18 205 0 382 168 0 0 0 46 0 0 -100,548 -15,390 7,702 56,324 41 5,369 -279 17,677 0 0 0 0 0 0 -9,604 52,338 0 0 0 -64,234 11,110 0 0 -828 0 0 0 312  147 779 638 982 147 1,225 -23 -1,972 -259 643 -2,049 -4,301 -763 -3,854 4,397 -82,736 -32,325 9,094 45,200 -3,087 -140,869 -25,560 17,581 26,447 -4,885 -49,465 22,919 50,195 92,803 109,714 136,171 69,581 -16,870 -1,608 -50,300 80,212 -3,437 -117,066 -93,159 -11,020 15,335 9,400 -14,153  -39,192  33,849  gravel removal  sand removal  0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 5,123 35,868 0 15,171 0 0 0 491,371 223,644 0 0 0 0 135,736 0 0 0 0 0 44,167 0 0 0  0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1,677 11,632 0 4,829 0 0 0 150,629 67,906 0 0 0 0 39,264 0 0 0 0 0 12,033 0 0 0  951,0801  287,970  gravel sum (m ) 3  sand sum (m ) 3  -2,432 4,144 -16 257 2,606 6,230 -63,333 -219,310 -50,347 -232,024 -331,150 -139,266 -292,282 -487,003 -908,480 -320,173 -404,581 243,341 507,867 172,044 515,735 583,671 396,150 291,088 -5,315 -287,322 -90,319 444,722 227,401 -70,820 1,215,211 781,601 286,895 145,538 232,345 340,193 708,543 338,153 517,558 182,912 250,357 22,646 -104,367  -46,056 79,515 327 5,868 49,656 119,600 -42,244 -95,962 -18,609 -83,186 -120,644 -53,738 -103,603 -173,689 -309,603 -192,409 -169,668 90,956 214,503 53,744 27,935 163,725 144,868 119,110 -6,561 -139,228 -5,033 186,524 161,850 88,416 498,122 300,134 66,936 40,491 16,251 176,692 195,518 -23,061 49,278 38,812 82,850 15,445 -41,725  4,408,6671  1,362,1041  total s+g (m ) 3  -48,487 83,659 310 6,125 52,262 125,831 -105,577 -315,272 -68,955 -315,210 -451,794 -193,004 -395,886 -660,693 -1,218,083 -512,582 -574,249 334,297 722,370 225,789 543,669 747,396 541,018 410,198 -11,877 -426,550 -95,352 631,246 389,251 17,597 1,713,333 1,081,736 353,831 186,029 248,596 516,885 904,061 315,092 566,836 221,724 333,207 38,091 -146,092  Cell 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43  5,770,7711 R1-43  -138-  Table 5-8: Sediment budget -1952 to 1999 area Cell 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 1-43 44-65 1-65  (nft .  channel change  105,179 -138,245 482,079 -133,395 437,314 -152,327 372,511 288,660 549,284 653,942 626,766 615,558 721,655 267,010 800,690 80,037 631,089 -667,285 -676,335 762,666 -250,859 728,274 -398,307 631,705 -811,601 742,715 80,091 545,030 854,257 -365,873 1,026,199 -59,995 995,226 -67,842 1,188,927 1,321,780 1,158,355 286,897 1,427,735 -578,662 237,332 1,179,466 35,236 1,129,061 326,782 1,157,396 485,807 1,248,368 1,133,537 553,252 888,176 -450,143 1,030,108 -458,164 933,607 -269,375 712,461 -71,976 999,054 168,167 908,136 -34,189 835,266 341,712 800,935 452,340 67,822 607,604 230,080 559,180 -178,097 1,128,482 318,709 1,016,240 976,715 1,229,069 840,067 941,589 849,411 868,125 718,940 363,007 414,306 350,150 528,159 -162,769 615,997 -328,943 623,688 -72,173 -455,036 1,248,249 -365,144 1,568,237 1,238,489 1,614,182 1,012,703 672,998 68,228 706,241 -82,218 509,400 1,021,036 -714,512 939,529 -603,867 -510,143 1,031,268 896,256 151,451 951,185 -263,733 931,581 -1,166,805 740,250 -487,900 1,213,784 -417,360 -21,851 847,341 681,076 110,517 484,233 37,618 409,336 168,848 -138,399 478,153 -161,714 549,080 35,382,330 4,637,713 18,697,111 -2,965,956 54,079,441 1,671,757  bed changes channel % sand gravel 95 95 95 95 95 95 40 30 27 27 26 26 26 26 26 26 25 25 25 25 25 24 24 24 24 24 24 23 23 23 23 23 23 22 22 22 22 22 22 21 21 21 21 21 21 20 20 20 20 20 20 19 19 19 19 19 19 18 18 18 18 18 17 17 17  channel sand  erosion (sub 3m)  deposition (sub 0.84m)  bank changes bank bank gravel total  bank sand  O/B sand (>0.177 mm)  gravel removal  sand removal  •gravel sum (m) 3  sand sum (nf)  -6,912 -6,670 -7,616 14,433 32,697 30,778 160,206 56,026 -489,043 -496,832 -184,709 -293,957 -600,362 59,383 -271,897 -44,688 -50,648 989,055 215,168 -434,978 178,807 26,607 247,317 368,502 420,607 -342,988 -349,884 -206,173 -55,212 129,286 -26,343 263,875 350,077 52,605 178,851 -138,747 248,837 961,716 738,380 682,255 285,906 276,379 -128,755 -260,765 -57,337 -362,279 -291,335 1,290,659 539,263 54,787 -66,162 -576,193 -488,000 -413,132 122,909 -214,482 -950,905 -398,455 -341,561 -17,920 90,824 30,979 139,338 -114,447 -134,004  -131,333 -126,725 -144,710 274,227 621,245 584,780 106,804 24,011 -178,242 -179,503 -66,150 -104,351 -211,240 20,709 -93,976 -15,308 -17,193 332,725 71,728 -143,685 58,525 8,629 79,465 117,305 132,645 -107,154 -108,280 -63,202 -16,764 38,881 -7,846 77,837 102,263 15,217 51,228 -39,350 69,872 267,353 203,209 185,870 77,101 73,771 -34,015 -68,178 -14,835 -92,757 -73,808 323,522 133,735 13,441 -16,057 -138,319 -115,867 -97,011 28,542 -49,251 -215,900 -89,444 -75,799 -3,931 19,694 6,639 29,510 -23,952 -27,711  0 0 0 0 -437,146 -138,008 0 0 0 0 -114,824 -178,084 0 -89,838 -273,476 -327,978 -257,748 -170,964 -172,074 0 -614,685 -84,805 0 0 0 0 0 0 0 -27,921 -37,459 0 -257,198 0 -690,342 0 0 0 0 0 0 -266,617 0 0 -399,468 0 0 0 0 0 0 0 0 0 0 0 -173,299 -157,520 -289,219 0 0 0 0 0 -15,530  362 5,033 0 0 38,676 12,586 0 0 0 0 0 0 402,615 42,748 7,598 0 334 51,362 133,496 152,844 66,853 0 270,708 543,646 46,689 136,274 142,569 111,441 202,241 483,400 1,250,371 421,512 119,982 97,186 229,582 139,705 93,815 300,319 621,529 701,568 290,024 436,348 29,395 124,287 4,389 353,111 35,122 715,238 891,794 248,302 98,617 0 0 0 98,604 200,793 75,134 0 101,404 279,595 15,025 446 4,508 11,012 76,710  362 5,033 0 0 -398,471 -125,422 0 0 0 0 -114,824 -178,084 402,615 -47,090 -265,878 -327,978 -257,414 -119,601 -38,578 152,844 -547,833 -84,805 270,708 543,646 46,689 136,274 142,569 111,441 202,241 455,479 1,212,912 421,512 -137,216 97,186 -460,760 139,705 93,815 300,319 621,529 701,568 290,024 169,730 29,395 124,287 -395,079 353,111 35,122 715,238 891,794 248,302 98,617 0 0 0 98,604 200,793 -98,165 -157,520 -187,815 279,595 15,025 446 4,508 11,012 61,180  18 252 0 0 -19,924 -6,271 0 0 0 0 -84,545 -131,428 297,824 -34,914 -197,586 -244,296 -192,176 -89,495 -28,933 114,892 -412,740 -64,038 204,879 412,374 35,495 103,835 108,875 85,294 155,136 350,170 934,555 325,498 -106,195 75,381 -358,170 108,838 73,248 234,992 487,394 551,359 228,424 133,971 23,252 98,527 -313,869 281,131 28,022 571,886 714,581 199,386 79,358 0 0 0 80,022 163,296 -80,001 -128,642 -153,705 229,294 12,348 367 3,720 9,107 50,696  344 4,782 0 0 -378,547 -119,151 0 0 0 0 -30,279 -46,655 104,791 -12,176 -68,292 -83,682 -65,238 -30,107 -9,645 37,952 -135,093 -20,768 65,829 131,271 11,194 32,439 33,694 26,147 47,105 105,309 278,357 96,014 -31,021 21,805 -102,591 30,867 20,568 65,327 134,135 150,209 61,600 35,760 6,143 25,760 -81,210 71,980 7,099 143,352 177,213 48,917 19,259 0 0 0 18,582 37,497 -18,164 -28,877 -34,110 50,300 2,677 79 788 1,906 10,483  357 2,997 3,460 604 -51,877 -47,189 -6,169 -26,482 -19,264 -23,936 -24,663 74,279 88,196 -42,283 -92,728 -91,902 -143,026 -121,921 -83,801 -28,992 -202,085 -.121,170 62,544 9,514 -12,477 -10,166 50,044 76,867 146,111 98,877 106,593 -269,286 -271,073 15,886 -95,859 69,308 -113,993 -457,680 -173,830 -44,032 -101,746 -158,053 -60,948 -1,591 -221,394 61,213 20,753 -1,201 50,526 59,129 65,822 44,796 -178,383 -14,535 29,866 53,428 -9,632 -100,507 -201,362 32,430 -28,009 -6,476 -715 4,635 -5,653  0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 597 0 1,725 4,585 62,382 35,868 0 15,171 0 0 0 991,161 639,097 35,365 0 0 0 373,080 0 0 0 173,709 72,145 47,154 0 0 52,999 237,939 214,500 0 27,127 98,348 43,269 0 0 0 0 0 0 0 0 0 0 0 0 0 0 62,847 0  0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 203 0 575 1,515 20,418 11,632 0 4,829 0 0 0 303,839 194,053 10,635 0 0 0 107,920 0 0 0 48,291 19,855 12,846 0 0 14,001 62,211 55,500 0 6,873 24,652 10,731 0 0 0 0 0 0 0 0 0 0 0 0 0 0 13,153 0  -6,894 -6,418 -7,616 14,433 12,774 24,507 160,206 56,026 -489,043 -496,832 -269,254 -425,385 -302,537 24,468 -469,483 -288,984 -242,227 899,561 187,960 -315,500 -171,551 -1,563 452,195 796,047 456,101 -239,154 -241,009 870,282 739,021 514,820 908,212 589,373 243,882 501,067 -179,318 -29,910 322,084 1,370,418 1,297,919 1,280,768 514,330 410,350 -52,503 75,701 -156,706 -81,148 -236,185 1,960,893 1,297,113 254,172 13,196 -576,193 -488,000 -413,132 202,931 -51,186 -1,030,906 -527,098 -495,266 211,374 103,171 31,346 143,058 -42,493 -83,307  -130,632 -118,947 -141,251 274,831 190,821 418,440 100,634 -2,471 -197,506 -203,439 -121,092 -76,727 -18,253 -33,750 -254,996 -190,892 -225,254 180,698 -21,142 -133,210 -258,235 -121,677 207,838 262,921 131,362 -84,880 -24,542 343,651 370,505 253,703 377,104 -95,436 -199,831 160,828 -147,221 60,825 -23,553 -76,709 183,369 304,894 36,955 -48,522 -74,819 18,202 -261,941 40,436 -39,083 490,325 372,204 121,488 69,025 -93,523 -294,250 -111,547 76,990 41,675 -243,696 -218,829 -311,271 78,800 -5,638 242 29,583 -4,258 -22,880  2,831,340 -2,418,219 413,121  1,806,373 -547,737 1,258,636  -4,139,168 -1,035,036 -5,174,203  7,582,810 3,334,090 10,916,900  3,443,643 2,299,054 5,742,697  3,075,244 1,845,523 4,920,767  368,398 453,532 821,930  -2,090,994 -346,860 -2,437,854  2,505,037 684,031 3,189,069  750,613 173,119 923,731  8,411,621 111,335 8,522,956  834,390 -267,946 566,444  total s+g (m) . Cell 3  -137,526 -125,365 -148,867 289,264 203,594 442,946 260,840 53,555 -686,549 -700,271 -390,346 -502,112 -320,790 -9,282 -724,479 -479,876 -467,481 1,080,258 166,818 -448,711 -429,786 -123,240 660,034 1,058,968 587,463 -324,034 -265,551 1,213,932 1,109,526 768,523 1,285,316 493,937 44,051 661,894 -326,540 30,915 298,531 1,293,708 1,481,288 1,585,662 551,285 361,828 -127,322 93,903 -418,646 -40,712 -275,269 2,451,218 1,669,318 375,660 82,221 -669,716 -782,250 -524,678 279,920 -9,511 -1,274,602 -745,926 -806,537 290,174 97,534 31,588 172,641 -46,751 -106,187  1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65  9,246,012 1-43 -156,612 44-65 9,089,400 1-65  -139-  Table 5-9: Sediment budget -1999 to 2003  Cell  area (m ) 2  channel change  1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48  491,616 446,169 378,771 614,312 681,529 753,194 839,439 663,649 816,715 764,404 632,324 706,619 598,314 963,466 1,144,653 1,157,859 1,329,963 1,291,600 1,604,463 1,445,963 1,256,428 1,151,379 1,310,159 1,190,207 898,045 1,025,285 941,042 719,021 952,625 958,261 1,282,512 1,113,235 640,218 674,793 1,160,741 1,290,229 1,546,355 1,360,952 1,003,890 874,886 625,426 636,681 613,024 562,427 1,302,627 1,658,056 1,389,520  -47042 -71264 -40512 586443 254728 59711 20336 196513 294716 311400 210799 95971 127275 501633 520425 387661 63429 -135751 -164659 -199022 -102606 16364 76255 68563 -199884 208647 11689 -118451 24915 -13697 42005 208650 313052 195437 705856 508391 75451 75738 -8805 -16168 -1014 95532 70609 -18864 488956 409337 803477  R1-43 R1-48  39,937,390 45,463,045  5,138,710 6,892,224  bed changes channel % gravel sand 95 95 95 95 95 95 40 30 27 27 26 26 26 26 26 26 25 25 25 25 25 24 24 24 24 24 24 23 23 23 23 23 23 22 22 22 22 22 22 21 21 21 21 21 21 20 20 20  channel sand  erosion (sub 3m)  deposition (sub 0.84m)  bank changes bank bank total gravel  bank sand  O/B sand (>0.177 mm)  gravel removal  sand removal  gravel sum sand sum (m ) (m) 3  3  0 -2,352 -3,563 -2,026 29,322 12,736 35,827 14,235 144,022 216,497 229,285 155,573 70,992 94,366 372,786 387,641 289,414 47,462 -101,812 -123,774 -149,944 -77,479 12,385 57,842 52,125 -152,303 159,336 8,947 -90,862 19,154 -10,554 32,437 161,480 242,814 151,922 549,901 396,934 59,038 59,393 -6,920 -12,734 -801 75,568 55,974 -14,987 389,285 326,596 642,440  0 -44,690 -67,701 -38,486 557,121 241,991 23,884 6,101 52,492 78,219 82,115 55,226 24,979 32,909 128,847 132,784 98,247 15,967 -33,940 -40,886 -49,078 -25,127 3,979 18,413 16,438 -47,581 49,311 2,743 -27,589 5,760 -3,143 9,568 47,171 70,238 43,515 155,955 111,457 16,412 16,345 -1,885 -3,434 -214 19,964 14,635 -3,878 99,671 82,741 161,037  0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -114,065 0 0 0 -127,421 -258,849 0 0 0 0 0 0 0 -8,457 -70,512 0 0 0 0 0 -97,947 0 0 0 0 0 0 0 0 0 0 0 0  0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 17,986 0 0 34,516 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 5,670 0 0 0 0  0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -114,065 0 0 0 -127,421 -258,849 0 17,986 0 0 34,516 0 0 -8,457 -70,512 0 0 0 0 0 -97,947 0 0 0 0 0 0 0 5,670 0 0 0 0  0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -84,962 0 0 0 -95,782 -195,018 0 13,612 0 0 26,300 0 0 -6,487 -54,209 0 0 0 0 0 -76,306 0 0 0 0 0 0 0 4,495 0 0 0 0  0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -29,103 0 0 0 -31,639 -63,831 0 4,374 0 0 8,216 0 0 -1,970 -16,303 0 0 0 0 0 -21,641 0 0 0 0 0 0 0 1,175 0 0 0 0  0 1,326 -1,538 0 112 -419 0 0 0 0 64 59 633 0 0 -24,414 -4,729 -615 -2,844 -29,819 -88,930 -10,981 1,747 5,364 3,113 19,028 -5,428 -22,327 -22,511 -60,755 852 1,304 -9,952 5,077 10,856 -61,701 -32,068 -41,957 -38,813 -3,008 -1,103 464 -19,397 10,353 148 10,251 -63,284 -79,457  0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 54,055 0 0 0 0 39,038 7,825 0 34,422 0 0 0 0 0 0 0 0  0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 15,945 0 0 0 0 10,962 2,175 0 9,378 0 0 0 0 0 0 0 0  0 -2,352 -3,563 -2,026 29,322 12,736 35,827 14,235 144,022 216,497 229,285 155,573 70,992 94,366 372,786 302,679 289,414 47,462 -101,812 -219,556 -344,962 -77,479 25,996 57,842 52,125 -126,003 159,336 8,947 -97,349 -35,055 -10,554 86,492 161,480 242,814 151,922 473,595 435,972 66,863 59,393 27,503 -12,734 -801 75,568 60,468 -14,987 389,285 326,596 642,440  0 -43,364 -69,238 -38,486 557,233 241,572 23,884 6,101 52,492 78,219 82,179 55,285 25,612 32,909 128,847 79,267 93,518 15,351 -36,783 -102,344 -201,838 -36,108 10,099 23,777 19,551 -20,337 43,882 -19,584 -52,070 -71,297 -2,291 26,817 37,219 75,315 54,371 72,614 90,351 -23,369 -22,468 4,484 -4,537 250 566 26,163 -3,729 109,923 19,457 81,579  3,404,312 4,803,619  1,734,398 2,088,604  -677,251 -677,251  52,502 58,171  -624,749 -619,079  -472,852 -468,358  -151,897 -150,721  -433,310 -555,299  135,340 135,340  38,460 38,460  3,066,800 4,470,602  1,187,651 1,421,044  total s+g (m ) 3  0 -45,716 -72,801 -40,512 586,555 254,308 59,711 20,336 196,513 294,716 311,465 210,858 96,605 127,275 501,633 381,946 382,932 62,814 -138,595 -321,900 -546,800 -113,587 36,096 81,619 71,676 -146,340 203,218 -10,638 -149,420 -106,352 -12,845 113,309 198,698 318,129 206,293 546,209 526,323 43,494 36,925 31,987 -17,272 -550 76,135 86,631 -18,716 499,207 346,053 724,019  Cell 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48  4,254,451 R1-43 5,891,645 R1-48  -140-  Table 5-10: S e d i m e n t budget - 1 9 8 4 to 2003  Cell  area (m ) 2  channel change  1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44  495,016 447,149 378,188 601,015 675,940 749,909 829,823 660,646 807,174 871,884 649,052 716,355 591,070 932,700 1,025,857 1,119,585 1,327,000 1,284,070 1,583,828 1,300,340 1,228,679 1,149,185 1,284,732 1,173,560 896,570 1,047,038 961,557 728,103 969,842 954,302 1,285,752 1,092,650 670,750 626,578 1,227,774 1,226,747 1,403,832 1,275,311 1,005,632 874,516 629,436 621,895 622,971  -122,019 -99,298 -48,484 640,624 229,344 -67,701 -396,496 67,836 -67,105 275,579 -37,432 -216,284 -589,955 -637,458 499,096 -65,407 532,135 424,101 92,988 521,957 490,486 509,914 237,273 72,851 -622,933 159,590 -216,235 -380,656 -116,562 1,336,949 889,377 569,422 347,087 744,567 863,156 1,381,635 502,022 480,729 95,991 206,996 -103,536 -105,472 -541,758  R1-48  39,381,038  7,736,914  bed changes % channel sand gravel 95 95 95 95 95 95 40 30 27 27 26 26 26 26 26 26 25 25 25 25 25 24 24 24 24 24 24 23 23 23 23 23 23 22 22 22 22 22 22 21 21 21 21 21  channel sand  erosion (sub 3m)  deposition (sub 0.84m)  -6,101 -4,965 -2,424 32,031 11,467 -40,621 -277,547 49,716 -49,295 202,910 -27,626 -159,991 -437,413 -473,724 371,754 -48,831 398,183 318,069 69,899 393,245 370,373 385,915 179,980 55,385 -474,647 121,874 -165,501 -291,996 -89,612 1,030,126 686,790 440,690 269,213 578,785 672,446 1,078,732 392,820 376,980 75,439 163,031 -81,723 -83,431 -429,470  -115,918 -94,333 -46,060 608,593 217,877 -27,080 -118,949 18,120 -17,810 72,669 -9,807 -56,293 -152,542 -163,734 127,342 -16,576 133,952 106,032 23,089 128,712 120,113 123,998 57,293 17,466 -148,286 37,717 -50,734 -88,661 -26,950 306,823 202,587 128,732 77,874 165,781 190,710 302,903 109,202 103,748 20,552 43,965 -21,813 -22,041 -112,287  0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -687,766 -111,593 0 0 -138,344 -561,580 -84,135 -11,541 0 -640 0 0 0 0 -227,078 0 0 0 0 -287,109 0 0 -47,758 -58,208 0 0 0 0  7,711 0 11,736 25,773 3,282 9,537 22,375 5,489 1,594 6,145 1,731 5,145 0 0 0 637 43,234 231,748 41,044 9,451 0 169,302 25,137 9,629 0 0 71,320 0 193,708 412,144 164,300 2,519 4,794 233,864 476,437 241,049 147,869 29,555 0 0 12,548 4,718 0  5,580,938  2,155,975  -2,215,7501  2,625,5281  bank changes bank total  bank gravel  bank sand  7,711 0 11,736 25,773 3,282 9,537 22,375 5,489 1,594 6,145 1,731 5,145 0 0 0 -687,129 -68,359 231,748 41,044 -128,893 -561,580 85,167 13,596 9,629 -640 0 71,320 0 193,708 185,066 164,300 2,519 4,794 233,864 189,328 241,049 147,869 -18,203 -58,208 0 12,548 4,718 0  386 0 587 1,289 164 5,722 15,663 4,023 1,171 4,525 1,277 3,806 0 0 0 -512,987 -51,151 173,807 30,853 -97,109 -424,057 64,457 10,313 7,321 -487 0 54,587 0 148,922 142,594 126,875 1,949 3,719 181,793 147,497 188,203 115,704 -14,275 -45,745 0 9,905 3,732 0  7,326 0 11,149 24,485 3,118 3,815 6,713 1,466 423 1,620 453 1,339 0 0 0 -174,143 -17,208 57,940 10,191 -31,784 -137,522 20,711 3,283 2,309 -152 0 16,734 0 44,786 42,472 37,425 569 1,076 52,071 41,831 52,846 32,165 -3,929 -12,463 0 2,644 986 0  3,072 811 1,836 1,441 7,669 5,090 3,281 -4,261 806 2,226 2,989 -4,746 -2,208 -2,675 1,310 -113,621 -32,600 52,912 9,628 -50,559 -200,969 -4,213 74,049 -13,182 -38,730 28,364 68,122 80,428 138,328 99,371 64,051 -586 -36,516 31,811 -8,726 56,957 -29,195 -180,951 -141,088 -10,071 36,916 3,453 -36,489  0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 49,000 0 0 0 0 35,000 7,000 0 30,660 0 0 0 0  0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 21,000 0 0 0 0 15,000 3,000 0 13,140 0 0 0 0  -5,715 -4,965 -1,837 33,320 11,631 -34,898 -261,884 53,739 -48,124 207,435 -26,348 -156,185 -437,413 -473,724 371,754 -561,817 347,032 491,876 100,752 296,137 -53,684 450,372 190,293 62,705 -475,135 121,874 -110,914 -291,996 59,310 1,172,721 862,665 442,639 272,932 760,578 819,943 1,301,935 515,524 362,705 60,354 163,031 -71,818 -79,698 -429,470  -105,520 -93,523 -33,075 634,518 228,664 -18,176 -108,955 15,326 -16,581 76,515 -6,364 -59,700 -154,749 -166,410 128,652 -304,341 84,144 216,884 42,909 46,369 -218,379 140,496 134,625 6,593 -187,169 66,081 34,121 -8,233 156,164 448,666 325,062 128,716 42,434 249,663 223,815 427,707 115,173 -81,132 -119,859 33,894 17,747 -17,602 -148,776  409,777|  305,031  104,746  -136,466  121,6601  52,140  6,007,630!  2,176,3951  O/B sand (>0.177 mm)  gravel removal  sand removal  gravel sum (m ) 3  sand sum (m ) 3  total s+g (m ) 3  Cell  -111,235 -98,487 -34,912 667,838 240,295 -53,074 -370,840 69,065 -64,705 283,950 -32,713 -215,885 -592,163 -640,133 500,407 -866,158 431,176 708,760 143,660 342,506 -272,063 590,868 324,918 69,299 -662,303 187,954 -76,793 -300,228 215,474 1,621,387 1,187,727 571,355 315,366 1,010,242 1,043,758 1,729,642 630,697 281,574 -59,505 196,925 -54,071 -97,300 -578,246 8,184,025| R1-48  -141-  Chapter 6: Channel deformation and sediment transport 6.1 Introduction In this chapter, spatial and temporal patterns of channel development are examined, and a direct connection is made to the underlying processes that influence the erosion, transport and deposition of coarse sediments. Over time, the morphology of the gravel reach has evolved in accordance with the annual volume of water and sediment supplied from upstream but, over the past century, has also been conditioned by bank hardening, blockage of secondary channels and sand and gravel removals. The most geomorphically active locations are associated with large accumulations of gravel and large channel islands. In these 'sedimentation zones', the deposition of gravel bars alters local flow alignment, eroding adjacent deposits and propagating instability downstream where entrained material becomes redistributed. In the long-term, fluctuations in the supply of bed material may be associated with a variety of channel evolution processes acting at various timescales. Between Vedder River and Agassiz-Rosedale Bridge, variation in the position and number of sedimentation zones is related to changes in channel morphology over time. Major accumulations occur where channel alignment reduces current velocity, encouraging deposition (and compensating erosion). These zones may persist for periods of years to decades in accordance with the magnitude of local transport rates. Indeed, historic airphoto sequences reveal that major bar and island groups have remained within the same computing cells over the period of record, although the morphology of the deposits has evolved as bars migrate and island and floodplain sediments are eroded. These developments indicate that transport rates are not consistent over time. Significant changes may occur within relatively short periods when bed material accumulations force the main channel to migrate. McLean and Church (1999) identified a 1972 avulsion near Harrison mouth as an example of this process. Change in sediment transport rate, both downstream and over time, is analogous to the passage of sediment in waves or slugs (Meade, 1985; Nicholas et al., 1995). This pattern of sediment transport appears to be characteristic of supply-limited rivers, although other processes can also produce pulses (Hoey, 1992). The next section of the chapter briefly reviews common examples of channel and floodplain formation on wandering channels to provide context for the types of processes and corresponding channel deformation observed on lower Fraser River. This is followed by a discussion on the how sediment is displaced on the river between major morphologic features to - 142-  investigate spatial (length) scales of sediment transfer, then by a more detailed analysis of the direct relation between channel bar growth and channel deformation for a major deposition zone extending from Agassiz-Rosedale Bridge to Carey Point using a series of low water photographs taken between 1943 and 1999. By elucidating the association between the spatial and temporal variability of sediment transfer and morphodynamic change, locations and rates of future instability may be predicted. A model of channel evolutionary development is proposed, and validated by comparing expected changes with actual changes documented by additional low water photography collected in 2003. Model results from the qualitative analyses will be contrasted to a GIS-based probabilistic approach developed from transition-state analysis.  6.2 Channel andfloodplainformation on wandering channels The characteristic morphology of wandering channels is a product of the underlying mechanisms that establish and modify channel bar, island, and floodplain elements as sediment is transferred through the system. We can therefore improve our understanding of channel morphodynamics by examining the erosional and depositional processes that maintain channel pattern (Dykaar and Wigington, 2000). Wolman and Leopold (1957) presented a generic model of floodplain development consisting of lateral accretion by point bars on convex bends with corresponding bank erosion on the opposite concave bend. Although typically associated with meandering rivers, the authors stated that the process is common to all channel typologies, but only rivers that transport significant quantities of bed material sediment actually conform with this style of lateral instability. Rivers that transport mainly wash material or fine to medium sand build their floodplains vertically, and are relatively less active laterally although meander migration remains important. More recently, Nanson and Knighton (1996) have suggested that avulsion and accretion-based process are fundamental formative mechanisms on all types of anabranching channels. Sediments deposited as point or complex medial bar forms dominate the depositional environment. In addition to forcing compensating erosion, growth of these bars also creates instream and floodplain channel avulsions (cf. discussion in Chapter 1). Neill (1973) observed that channel shifting and deposition of point bars and islands were dominant channel processes on the wandering Athabasca River, Alberta, although the progression was more erratic than described by Wolman and Leopold (1957). A cursory examination of maps and photographs of the planform channel reveals a more complex channel and floodplain morphology than could be created by the simple point bar model, indicating that additional -143-  processes must be involved. Church and Jones (1982) discussed the accumulation of extensive bar assemblages in wandering channels, which they termed sedimentation zones. These preferred zones of deposition are found at discrete locations along the channel as pulses of bed material propagate downstream. Although they may be spaced at regular intervals, tributary inputs, alluvial fans, and erosion resistant boundaries can introduce deviations. Each assemblage comprises a complex variety of barforms produced by channel migration and cutoffs, and its presence is associated with local channel instability. However, no details of the association between channel and floodplain deposits were given. Church (1983a) examined historic maps and airphotos of Bella Coola River, British Columbia, allowing a descriptive analysis of erosional and depositional development to be constructed. It was found that most sedimentation zones generally persisted over time because they were anchored by alluvial fan constrictions, but extensive local modifications occurred in response to changes in sediment supply. For some sedimentation zones, continued bar accumulations forced avulsions and reoccupation of old flood channels. Correspondingly, others became stabilized where the upstream supply was reduced and much of the deposit was subsequently incorporated into the existing floodplain. In contrast to the dominant island and medial bar complexes found in the sedimentation zones, stable reaches contained mainly bank-attached lateral bars (Desloges and Church, 1987). Ferguson and Werrity (1983) examined bar and channel development over a period of 5 years on the River Feshie, Scotland, concluding that progradation of alternate diagonal bars was the dominant process driving morphologic evolution. Flow divergence and convergence around migrating bars was found to cause erosion on the outer banks opposite accreting bar margins. However, avulsions, initiated from chute channels across elongated bars or into the floodplain, were found to interrupt the process. Based on these observations, the authors concluded that channel evolution on wandering rivers is neither simple nor follows regular cycles of development. A similar review of morphologic development over a period of 60 years on Willamette River, Oregon is given by Dykaar and Wigington (2000). It was found that a combination of erosional and depositional processes was responsible for maintaining the multi-channel planforrn. The establishment of point and medial bars produced proximal lateral bank erosion and scalloping, driving overall channel migration. Accretionary processes which constrict flows in secondary channels and promote infilling, or which force migration of the channel thalweg, provide the nucleus for incipient floodplain formation and maturation. Similarly, bar deposition maintains - 144-  multiple channels by splitting the flow, forcing new cut channels across point bars and enlarging secondary channels where accumulations block the main channel flow (Dykaar and Wigington, 2000). Woody debris has also been recognized as providing an important feedback mechanism linking channel and floodplain changes on wandering rivers. Channel migration promotes a significant influx of woody debris that affects channel process and form (Desloges and Church, 1989; Gottesfeld and Gottesfeld, 1990). Wood plugs at the head of sidechannels can force their eventual abandonment and establishment of new floodplain, while large jams can cause main channel avulsions through former channels and through stable floodplain, generating further wood inputs (Gottesfeld and Gottesfeld, 1990). Ice buildup can play a similar role in northern climes. Woody debris accumulations have also been associated with island and floodplain formation. Jams that establish on bar tops may create conditions favorable for the establishment of vegetation, hence growth into islands, while jam development on meander bends can lead to the creation of new floodplain elements (Hickin, 1984; Abbe and Montgomery, 1996). The importance of wood in maintaining multi-channel morphologies was demonstrated by Sedell and Froggatt (1984) and by Collins et al. (2002) who found that woody debris removals since the mid-1800s have reduced both instream channel complexity and interaction with the floodplain on large Pacific Northwest rivers. However, a paucity of wood on Fraser River bars relative to these examples suggests that large jams are not a ubiquitous formative mechanism of the wandering planforrn - indeed, wandering channels have been identified in the Arctic tundra (cf. Forbes, 1983). Nevertheless, the possibility remains that the role of wood in influencing island and floodplain formation may have been greater in the past when it was not removed for navigational safety. Several studies have directly examined erosional and depositional processes on lower Fraser River. Morningstar (1987) reviewed mechanisms of floodplain formation near Chilliwack, summarized as the consequence of island formation in the main channel followed by attachment of islands to the floodplain by chute infilling. He found that vertical accretion of sands in the lee of grounded woody debris on a deposited gravel platform creates suitable conditions for vegetation establishment, while depressions in the bar become chute channels as flow becomes concentrated with increasing discharge. The chute channels separate islands from each other or from the adjacent floodplain (Roberts and Morningstar, 1989). Sediment blockage at a chute head decreases flow intensity and frequency of inundation downstream, allowing fines to be deposited in the chute tail because of the backwater effect (Morningstar, 1987). The infilling and abandonment of chutes -145-  separating islands from the main channel banks was argued to be the main mechanism for floodplain accretion. Wooldridge (2002) used a combination of bathymetric soundings, historic aerial photography and ground penetrating radar signatures to examine the evolution of bank-attached, partially bank-attached and mid-channel macroforms. The prominence of sub-horizontal reflections in the GPR surveys was interpreted as the vertical and lateral stratification of migrating bedload sheets 2 to 3 grains thick, while the widespread occurrence of slipface accretions, formed where the current is aligned directly with the emergent bar, was found to distinguish wandering rivers from braided and meandering deposits. Wooldridge (2002) also noted that large floods initiate channel changes that are completed by successive floods, but do not cause immediate channel changes, suggesting that this is a characteristic feature of wandering channels. McLean and Church (1999) make a similar claim, but submit that in Fraser River, the rate of change is tempered by the modest transport rate relative to the size of the channel and stored sediment volume. McLean (1990) was the first to describe channel morphodynamics on lower Fraser River, associating bar and island developments at major storage zones to morphologic changes immediately downstream. Eroded sediments are transferred downstream, creating or adding to existing bar deposits which affect local flow alignment. McLean and Church (1999) suggested that these local disturbances propagate down the river as a wave-like phenomenon that travels faster than the rate of sediment transport. McLean (1990) and McLean and Church (1999) noted that erosion and deposition sites change over roughly a decade as channel realignment alters the pattern of erosional attack, while migrating zones of instability can be traced for several decades. The next section of the chapter examines the processes by which sediment is routed downstream and is followed by a review of the genesis of bar and island formation in a major sedimentation zone to examine whether distinct repeating cycles of depositional and erosional patterns exist. 6.3 Sediment transfer length While the movement of individual clasts has been demonstrated to be fundamentally stochastic, particles tend to aggregate into major accumulations where bars act as sediment traps (Hassan and Church, 1992). Despite this association, most efforts to measure particle movement have focussed on channel hydraulics, while neglecting potential morphologic controls (Pyrce and - 146-  Ashmore, 2003). On many gravel channels, bars establish on stable riffles which have a somewhat regular sequence that scales to channel dimensions (Richards, 1976; Church and Jones, 1982). On regularly meandering channels, sediment transfer proceeds in a reasonably orderly manner, with erosion on concave bends and deposition on point bars at downstream convex bends. Neill (1971) recognized that the downvalley migration of a channel meander belt associated with this development could be used to provide an estimate of sediment transport at an arbitrary crosssection if the average travel distance of the sediment could be estimated. Neill (1987) stated that the travel distance during flood was equal to one-half the meander wavelength on rivers with regular meanders, but recognized that this distance may need to be modified for different morphologies. Church et al. (1987) presented a generalization of this approach applicable for more complex channel patterns, where sediment step lengths were assumed equivalent to riffle spacing. Investigations on braided channels have demonstrated similar results (cf. Goff and Ashmore, 1994) but other studies have adopted the distance between distinct erosion and deposition zones (Carson and Griffiths, 1989; Ferguson and Ashworth, 1992), although these can be difficult to identify (Goff and Ashmore, 1994). On the wandering gravel reach of lower Fraser River, the links between major erosion and deposition zones are more clearly distinct because of the large size of the channel (Desloges and Church, 1989). Since changes to channel morphology are largely governed by the rate and magnitude with which bars develop and migrate, the distribution of primary morphologic features may provide a means to investigate the length scales of sediment transfer. Several approaches are presented below. The relation between meander wavelength and sediment transport over the past halfcentury on Fraser River is first investigated by defining successive meander bends, and measuring their respective lengths from Laidlaw to Mission Bend, roughly 4 km upstream of Mission railway bridge. These locations do not correspond exactly with the linear extent of the gravel reach, but clearly demarcate start and end points of the meandering sequence . The channel has developed a 1  more regular sinuous planforrn since 1949 that has been propagating downstream (Church and Ham, 2004; also see Appendix B), but meanders on the river are not well-developed. Consequently, it is difficult to identify meanders definitively and hence the interpretation remains partly subjective. A total of 13 meanders were delineated in both 1949 and 1999 (Figures 6-1, 61. A single meander is defined as a complete sine wave, or twice the distance between successive points of inflection on the meandering wavelength path (after Leopold and Wolman, 1957).  - 147-  2). Mean meander wavelength for both dates is 5.1 km since the thalweg lengths are nearly equal. Summary statistics for defined meanders are given in Table 6-1. Following Neill's (1987) assumption, the characteristic step length on the river (assumed during the annual freshet) should average roughly 2.5 km, and ranges from 1.6 to 3.4 km where variations in meander wavelength are associated with changes in local morphology. Table 6-1. Summary statistics for delineated meanders, 1949 and 1999. Statistic  1949  1999  n  13  13  mean (m)  5123  5115  S.D. (m)  799  1097  range (m)  3900 - 6400  3200 - 6700  total length (m)  66,600  66,500  Although wavelength variability has increased, the change is not signficant (F=1.88, p=0.14). In 1949, the shortest wavelengths were all associated with major bar / island complexes upstream of Agassiz-Rosedale Bridge, but this pattern was partly interupted by 1999 with the development of a more regularly sinuous thalweg near Herrling Island and re-organization of downstream deposits at Greyell complex and Queens Bar. Although mean meander wavelength remained unchanged, the change in wavelength variability may be accomodated by changes in other morphologic parameters. Meander wavelength has been shown to scale with the square root of dominant discharge, hence channel width, as this parameter displays a similar geometric association (Leopold et al., 1964). Since the active channel zone has been narrowing over the past half-century, the meander wavelength, hence sediment step length, should be declining. Although there were 3 meanders in 1999 shorter than the smallest 1949 meander, this pattern does not hold for the entire gravel reach, but the ratio of meander wavelength to main channel width has remained within the expected range of 10 to 14. Neill (1987) added that smaller values of meander wavelength are associated with greater values of bank material thickness or bank erosion rates to maintain sediment continuity. However, the known pattern of decreased bank erosion from 1949 to 1999 (Chapter 4) is not consistent with Neill's statement. Indeed, it appears possible that the combined effects of channel narrowing and decreased bank erosion are compensating such that meander wavelength has not changed over the past half-century.  -148-  island / floodplain  • i  gravel bar  1  1  5600  water surface channel thalweg  •  meander 'end point'  4900 m  Q,  J  Laidlaw  3500 m  5800 m  3600J 5600 m  5500 m 6300 m  3200 m  5800 m 'railway bridge Matsqui Island  Figure 6-2. Location and wavelength of river meanders, 1999.  Greyell complex  Herrling complex  The spacing of major riffles and pools provides an alternative means to investigate average step lengths. Leopold and Wolman (1957) found that the distance between successive riffles scales to one-half the meander wavelength. It is expected that a similar relation holds for pools, which are sequentially opposed forms along the thalweg and represent discrete erosion zones. In Chapter 4, it was shown that there were 27 major pools in both 1952 and 1999, giving an average pool spacing of 2.7 km for each date. This result is close to 3 bankfull widths, but increases to 6 channel widths if scaled to the average width of the dominant main channel, or normalized by the braiding index. This is within the expected range (5 to 7 widths) for self-formed channels and matches the findings of Church and Jones (1982) for the wandering Bella Cool a and Peace rivers. It is also consistent with the expected spacing based on meander wavelength. The spacing of major deposition zones is another obvious candidate for estimating path lengths, since these features are clearly preferential sites of clast accumulation because of increased flow resistence. Hassan and Church (1992) note that bar spacing is potentially a significant transport length scale on gravel channels, especially in large floods. In 1949, there were 26 major bars within the gravel reach, giving an average bar spacing of 2.7 km, increasing to 3.1 km by 1999 when there were 3 fewer bars. Although river meander lengths did not correspondingly change, the increase in bar spacing is assuredly modest. Similarly, McLean (1990) measured the distance between major zones of deposition along the river (identified as accretions on existing bars and the creation of new unit bars) and found mean separations of 2.75 km from AgassizRosedale to Carey Point, increasing to 5 km downstream to Mission. The seemingly close correspondence between bar and riffle/pool frequencies provides an indication that bar location is associated with, and likely conditioned by, riffle location on wandering rivers. This relation is explored in further detail in the following section of the chapter. Based on the preceding results and the assumptions about clast movement, typical sediment step lengths on the river would fall consistently within the range of 2.5 to 3.1 km, with modal values near 2.8 km. These lengths are equivalent to the spacing of major morphologic features on the river. Pyrce and Ashmore (2003) point out that path length is also controlled by flow hydraulics, but add that there is significant evidence linking this distance to the length scale of channel morphology on low-gradient rivers with large (or long duration) floods. Sediment transport rates might therefore be estimated by multiplying this value by the volume of sediment ,-\Y  eroded, or deposited, along the channel. McLean and Church (1999) investigated this possibility at the downstream end of Herrling Island, and found reasonable correspondence between bank - 151 -  erosion and sediment transport rates (compared to sediment budget estimates at Agassiz-Rosedale Bridge) except during periods when channel curvature was reduced. Consequently, the technique may under-predict the true transport rate over longer time periods, and may not be adequate for situations where there are multiple modes of transport (McLean, 1990). Indeed, since many studies have indicated that most bed material transport occurs in the form of unit bar construction and migration, mensuration of the size and spacing of these features may provide a more reliable means for estimating spatial variations in the local transport rate. Such an approach also recognizes that sediment transport may be the product of different volumes, distances and rates of transfer (Nicholas et al., 1985; Gaeuman et al., 2003). This possibility is explored further in the following section of the chapter. 6.4 Sedimentation and morphology of b a r complexes The availability of a sequence of low-water aerial photographs at fairly regular frequency (5 to 8 years - see Table 6-2) provides an opportunity to examine in detail the direct association between styles of bar deposition and channel deformation in a major sedimentation zone extending from Agassiz-Rosedale bridge downstream to Carey point, a distance of 10 km. The study area is comprised of a large island complex (Greyell Island) that has been developing since at least the mid 1800s, and several large bar complexes that have undergone extensive reorganization over the period of study, although the nucleus of contemporary Gill and Carey bars had established by the date (1943) of the earliest photographs in the sequence. Before discussing bar and channel development, it is necessary to clarify terminology related to morphologic classification, since a wide array of terms has been presented in the literature. Crowley (1983) pointed out that fluvial scientists recognize a large number of distinct bar types based on external morphology or position in the channel but there is no definitive classification of these features. Different classification schemes have also been criticized because planforrn morphology is largely stage dependent (Jackson, 1978; Smith, 1978; Bridge, 1993). This can lead to interpretive confusion since bedform size, frequency and geometry will vary with discharge. Most definitions of islands (vegetated bars) are less ambiguous because these deposits remain exposed except at extreme flows, although as Bridge (1993) points out, it is difficult to determine exactly when bars become islands. With the exception of 1986 photography, flows on all dates fall within a reasonably small range, so the area, and morphology, of exposed bar forms is sufficiently uniform to permit - 152-  Table 6-2. Selected low-water photography Flow at Hope (m /s) Date of photography 3  Water surface area (m ) 2  Mar 20, 1999  701  4,182,379  Feb 24, 1993  581  3,268,080  Sept 4/5, 1986  2500  6,317,756  Mar 22, 1979  1010  4,679,129  Mar 19,1971  799  4,914,108  Apr 6, 1964  1110  4,854,138  ?, 1959  1172*  4,450,243  May 7, 1954  1170  4,462,750  Mar 23, 1949  733  3,809,258  Dec 5,1943  929  4,510,484  is not known  consistent classification and discussion. Nevertheless, to account for small variations in bar exposure produced by differences in discharge, a morphologic feature termed 'bar edge' was additionally mapped from the photography (see Figure 6-3). Bar edge defines submerged bars, which are either partly or fully inundated by shallow water and are visible because water clarity in the river is highest during winter and early spring periods (when most of the low-water photo sets were flown) as suspended sediment concentrations are lowest. Incorporating the bar edge and exposed bar polygons makes the extent of total bar area similar for all dates, thus minimizing complications of channel and bar classifications which change with discharge. In the following discussion, a hierarchical classification of bed forms is adopted, where 3rd order unit bars refer to individual bedload sheets, crescents, or plugs deposited by a single flow event onto an existing (2nd order) bar that itself forms a major sedimentary body. Following the nomenclature of Church and Jones (1982), medial, point, and diagonal bars appear most prominent in the wandering reach. A bar complex (1st order) is a collection or group of individual bars and is separated from other bar complexes by open water at low discharge. Bar complexes may include immature vegetation, but mature islands are considered a separate depositional unit for purposes of discussion. Examples of migrating gravel sheets are provided in Figure 6-4, while Figure 6-5 illustrates the hierarchical identification and coding of the different depositional units. These features correspond to macroforms (3rd order) and megaforms (1st and 2nd order) which scale to the size of the channel on gravel bed rivers (Church and Jones, 1982). Similar - 153 -  Figure 6-4. Close range ( A ) and distal (B) views of gravel sheets prograding on existing bar surfaces at Queens Bar complex. (C) Annotated 1943 photo near lower Carey Bar complex showing example delineation of gravel sheets based on tone and morphology.  ordering schemes have been presented by Williams and Rust (1969) and Bristow (1987). They have also been associated with temporal scales of development, with macroforms having a duration of formation from days to years^and megaforms from weeks to millennia (after Jackson, 1975). Hierarchical schemes may therefore provide a template for explaining form-process relations at a range of scales (cf. Brierley, 1996). Primary bedform accumulations, termed microforms and mesoforms (cf. Church and Jones, 1982) are not discussed as these features are rarely visible on aerial photographs (Kellerhals et al., 1976) and their longevity is too brief to be examined over periods of several years. The relation between bar and channel morphologic development is specifically concerned with the following major questions, the scope of which is limited by the information that can be derived from planform analysis. The major questions include: • what are the dominant unit bar types characteristic of each date; • can their movement and development into bars and bar complexes be tracked; • what is the magnitude of a unit bar deposition event; • is their rate of transit related to the development of islands, avulsions and cutoffs; and • are repeating patterns of morphologic evolution observed which could be used to predict shortterm channel changes It should also be realized that the following discussion is based upon an interpretation from the historic photo sequence and there is no direct confirmation that any of the erosion / deposition sequences discussed are definitively correct. The 1943 map (Figure 6-5) shows extensive, partially bank-attached sediment deposits at Gill and Carey bar complexes. A third, smaller, bar complex is found upstream of Hopyard Hill while several medial bars are found near the downstream end of the study reach. Prominent channel islands with mature trees are found at all three sites, indicating that these surfaces have remained stable for decades. At least part of the large islands near the left bank at the Gill complex are more than a century and a-half old given their appearance on the 1859 sketch map. Secondary channels between large islands and major sloughs adjacent to outer channel banks preserve evidence of small migrating sheets and plugs. Minor accumulations of sandy mesoforms (dunes) are observed at a few sites. Smaller, immature islands are also common, growing on elevated bar deposits. The persistence of these features throughout the historic photo record (see Appendix C photomaps)  - 157-  demonstrates that sediment exchange within the reach does not encompass the entire active channel zone. There is evidence of recent unit bar deposits along wetted channel margins at all sites. Here, migrating gravel sheets a few grains thick are stacked upon existing major bar assemblages, and appear to dominate the depositional environment (refer to Figures 6-4, 6-5). These elongated features are oriented both parallel and obliquely to the main flow, and have typical lengths of several hundred metres, and stored volumes in the range of 10 to 10 m . Some of these bar units 3  4  3  have a detached tail that extends into the main channel, giving the distinct appearance of a 'claw'. Miall (1977) suggests that gravel sheets are the primary depositional units of medial bars, deposited where the current wanes. Church and Jones (1982) add that all common bar types in actively changing gravel-bed channels appear to be formed from sheets. Research by Wooldridge (2002) has confirmed the prominence of stratified bedload sheets on lower Fraser River macroforms. Unit bar typology away from active channel margins becomes less distinct. Although remnant sheets, plugs and crescents can be partly distinguished, repeat cycles of washing, erosion and deposition produce a more homogenous bar surface. Dissection of these deposits by numerous chute channels and scour holes during waning flood flows compounds the difficulty in identifying individual bar units. By 1949, there was lateral accretion on existing bar margins at all three bar complexes, but the basic morphology of the reach was maintained despite the magnitude of the 1948 flood (Figure 6-6 A). There was only a single additional flood (in 1946) larger than the mean annual flood during this intersurvey period. Near Hopyard Hill, two new gravel sheets (1: see figure for site numbers) extended the existing bar both laterally and downstream, directing the main flow towards the upstream end of Gill bar complex which experienced compensating erosion. Sedimentation on a large diagonal riffle (amplified by the lower flow stage) connected this feature to Gill Bar across the low flow channel (2). This sediment probably derived from a unit bar that became fully detached from the upstream bar and was subsequently modified (dissected) by the flow. There was more extensive lateral and downstream extension by new unit bar deposits on the convex bend at Gill Bar (3) which forced the flow into the right bank, eroding up to 100 metres of floodplain and trimming bar and island deposits at the head of the Carey complex. However, minor expansion of existing islands here provides evidence that this deposition did not significantly increase flows through the minor chute channels separating the islands (4). At Carey Bar, sediment deposition in existing chute channels (5) coalesced existing smaller bars into a single bar complex, while at least -158-  maps depict channel conditions on the early date; outlines show changes by later date. Numbers denote features in the text. Dashed lines indicate diagonal riffles. -159-  one unit bar was deposited downstream, extending the north side of the bar (6) and reducing the size of the secondary channel. The reduction of flow through the chute channels and lateral extension of lower Carey Bar complex forced the thalweg south, eroding the tail of a medial bar (7) to maintain conveyance, though minor unit bars were also deposited on the east and south sides of this feature. By 1954, lateral extension of a bar at Ferry Island forced the flow north, partly eroding two small islands and removing 14 ha (est. 290,000 m ) from the Hopyard bar complex (1; Figure 6-6 3  B). Most of this erosion likely occurred during the large 1950 flood; freshets for all other years were at or below the long-term average. This development must have increased flow through the prominent right bank secondary channel, because there was significant erosion at the bar tail near Hopyard Hill (2) where the current is forced to take a hard turn left. The diagonal riffle (3) appears reduced in size by the main current (although the feature may simply be underwater) which deflects north off the hardened left bank (4). These eroded sediments were transported downstream and deposited as new unit bars that further extended existing bars at Gill complex (5,6). Deposition at lower Gill Bar (6) forced further compensating erosion on the right bank and increased flow through the head of the Carey complex, trimming the edge of a large island (7). A new unit bar at Carey (8) was formed from the material eroded just upstream, and directed flows towards the left bank. Downstream, the medial bar at (9) was further eroded. This alignment also reduced flows across lower Carey Bar, providing suitable conditions for vegetation extension (10). As the location of the new unit bars (5, 6, 8) is coincident with the location of older bar tails, these sites appear conducive for trapping and accumulating migrating bedload sheets. Wooldridge (2002) found that sheet-like deposition dominates topographically low areas that grade into the main channel at Queens Bar, suggesting that these features form at falling stage. It is expected that further accumulations and compensating erosion will be observed by 1959 at these same sites. There were in fact a number of significant modifications to reach morphology by 1959 (Figure 6-7 A) with annual flood flows greater than the mean annual flood in all years. The Ferry Island bar migrated downstream, and a distinct crescentic-shaped unit bar (1) developed in mid channel. This event marks the origin of contemporary Big Bar (see Appendix C, Figure 1). This caused extensive compensating erosion at Hopyard Bar / Island, but the eroded material rebuilt the point bar that had been removed in the previous period (2), while a small medial bar (3) was deposited at the site of an existing bar tail on the diagonal riffle. Consequently, the flow was increasingly forced toward the left bank, removing unit bars on the upstream end of Gill Bar. - 160 -  Figure 6-7. Morphologic changes from ( A ) 1954 to 1959 and (B) 1959 to 1964. Shaded maps depict channel conditions on the early date; outlines show changes by later date.  -161-  Hardening of the right bank to protect the flood dyke helped direct the current away from the right bank. The material from these unit bars became trapped at the Gill Bar tail, also the upstream end of the major diagonal riffle (4) and more directly aligned the current towards the head of Carey Bar complex, which was extensively eroded (5). This resulted in a considerable widening of the lowflow channel, creating an area of flow expansion where a new cresentric unit bar formed on the riffle. Sediment eroded from upper Carey was deposited at mid Carey Bar (6), extending the bar more than 200 metres laterally and causing erosion at Gill Island. Unit bars were also attached to existing deposits at the mid-channel bar downstream (7) and the tail of Carey Bar (8). Overall, the pattern of increased sedimentation between Carey and Gill bar complexes indicates that flow conveyance through the main channel had become reduced. As a result, a partial avulsion is in prospect in order to maintain flow conveyance. The morphologic developments described to this point present an emerging picture of a reasonably orderly pattern of channel change directly associated with the location of persistent diagonal riffles. These features appear to be stable 'framework' riffles (sensu Church and Jones, 1982) rather than transient ones associated with migrating bar fronts, although they do migrate and develop over decades. Unit bars develop on the convex edge of existing bars at the downstream limit of the riffles, where their progression becomes stalled. This deposition forces compensating erosion in pools on the opposite bank, providing source material to the next downstream riffle. Less common variations of this progression include medial deposition on a riffle crest (e.g. '2' on 194349 map) and unit bar deposition on the upstream limit of the riffle (e.g. '6' on 1949-54 map). The features probably become trapped at falling stage when shear stress falls below critical values. The major riffles are established at a fairly regular interval of 1500 metres, or one-half the average step-length estimated in Section 6.3. Therefore, the average annual volume of sediment influx to the reach is apt to be represented by the volume stored in 'recent' unit bars at four consecutive riffles (which essentially define a complete meander unit). For the periods discussed, this volume falls within the range of 100,000 to 300,000 m per year, consistent with estimates 3  from the long-term sediment budget. However, because sedimentation in the reach is focussed at the same loci for many years, the channel must eventually become too constricted to maintain an orderly development. Growth of the medial bar between Carey and Gill complexes (1954-59) portends such a change. By 1964, additional sedimentation at developing Big Bar (1) caused further compensating erosion at Hopyard and extension of the existing point bar (2), continuing the pattern observed in - 162-  the previous mapping period (Figure 6-7 B). This growth, combined with the downstream migration of the riffle, forced sufficient flow into Gill/Greyell complex to create a small erosional channel separating them at low flow (3). The new current channel alignment favoured additional unit bar accumulations at the riffle tail and increased flow, hence erosion, through the island complexes. Downstream, the previous pattern of unit bar accretion on the convex bank of Gill Bar effectively stopped, but a new unit bar was deposited at the tail of Gill Bar (4). The sequence of morphologic developments near the head of Carey Bar is fairly complex. Sediments continued to accumulate on the existing medial bar in the area of flow expansion (5) which was maintained by erosion of upper Carey Bar. New unit bars, constructed from this material, were deposited at the downstream end of the migrating riffle (6) and the upstream end of the next riffle (7). Compensating erosion occurred at both locations. However, the gradual westward erosion at the head of Carey Bar also resulted in the expected partial avulsion by 1964. This created a prominent secondary channel along the right bank extending the length of Carey Bar complex (8). These events further developed a tortuous 90° bend that reduced conveyance through the main channel. As a result, amalgamation of the medial bar to Carey bar is a probable future development because of the reduced flows, while the right bank channel along Carey Bar is therefore likely to increase in prominence. Overall, the changes in reach morphology were fairly dramatic given the absence of significant floods. This provides further evidence that sediment transport is not strongly correlated with the size of the freshet, rather that the alignment of major morphologic features influences the rate and magnitude of future channel development. Changes to 1971 were less dramatic overall than in the previous reporting period despite two major floods, but there were several important developments (Figure 6-8 A). Accretion of additional unit bars became amalgamated on the downtream end of Big Bar (1) causing extensive compensating erosion at Hopyard complex (est. 600,000 m ). The resulting channel re-alignment 3  also created a new secondary channel along the right bank (2) which uncoupled the complex from the floodplain. Material eroded from Hopyard was deposited on the riffle tail (3), reducing downstream (northward) flow conveyance, but increasing the flow between Gill and Greyell complexes. This altered the shape of the existing secondary channel and excavated a new one (4), providing sediment for unit bars at the channel egress, and downtream on the riffle tail at Carey Bar (5). Compensating erosion of a bar across the channel provided material for several unit bars deposited along the existing riffle downstream (6). Material was also deposited on the north and west banks of the large medial bar at the upstream end of Carey, respectively causing compensating - 163-  Figure 6-8. Morphologic changes from ( A ) 1964 to 1971 and (B) 1971 to 1979. Shaded maps depict channel conditions on the early date; outlines show changes by later date. -164-  right bank floodplain erosion and channel narrowing, given the tortuous entrance approach (7). Eroded floodplain sediments formed new unit bars at the lower end of Carey Bar complex, effectively blocking the secondary channel and partly attaching the existing bar to the right bank floodplain (8). This appears to have created sufficient backwater effect to direct the secondary channel's flow directly south, re-aligning the existing chute channel (9). The period 1971 to 1979 reflects changes following the second largest flood (1972) within the period of the photo record and two subsequent freshets larger than the mean annual flood (Figure 6-8 B). Changes to reach morphology are substantial in comparison to the fairly modest modifications observed after the 1948 flood. Several new unit bars were deposited on Big Bar, increasing flow towards the right bank. This trimmed existing bar material and removed remnants of Hopyard Island that was well-established by 1943 (1), effectively splitting the flow and forming an elongated medial bar (contemporary Hopyard Bar/Island). Eroded material was deposited as a new unit bar on the downstream riffle tail, conjoining the upstream bar (2), and on a large point bar — an early stage of contemporary Hamilton Bar (3). Removal of Hopyard Island further aligned the right bank channel towards a bedrock knob (Hopyard Hill), forcing a hard turn west through the middle of the bar complex, then rejoining the former main channel. Overall, these developments created a more regularly sinuous thalweg along the right bank, but it was no longer the dominant channel. Sedimentation at the riffle tail (2) caused compensating erosion at upper Greyell Island (4) but, surprisingly, did not force additional flow through the existing chute channel (5), possibly because the channel entrance was blocked by wood or a sediment plug. Instead, the thalweg was directed to the northwest, eroding a large (250 metre wide) channel directly through Gill Bar (6) before re-connecting with the main channel near the head of Carey complex. This erosion initiated a number of significant developments downstream. Most of the eroded material was deposited in the former main channel that separated Carey from Gill and Greyell complexes, nearly reducing it to a summer channel only, while the right bank channel at Carey complex was enlarged and had become dominant. A prominent chute channel was also infilled, creating a larger, single bar complex (7). Downstream, continued sedimentation in the minor right bank channel (8) directed flow from the right bank channel through the second chute channel, which had become an extension of the main channel thalweg by 1979 (9). This alignment promoted extensive unit bar deposition on the convex edge of lower Carey Bar, but also caused considerable erosion of established bar deposits at lower Greyell Island (10). Overall, the re- 165 -  organization of channel elements appears to be related to the migration of major stabilizing riffles as persistent sedimentation eventually inundates them. This sediment congestion forces flow into, hence enlarges, existing erosional channels or initiates a channel avulsion. Consequently, entire channel bars can be displaced downstream as these features become separated from outer channel banks, or from one another. Despite the rather chaotic changes that occurred during this period, however, there is an emerging appearance of a return to a more regularly sinuous thalweg. If this assertion is correct, changes in the following period should be comparatively modest, with consolidation of major bars as secondary channels become sedimented, similar to the reach configuration of previous decades. There were three floods in excess of the mean annual flood between 1979 and 1986, the largest of which had a return period of roughly 5 years. Two chute channels developed across Big Bar (1), but they may be an artefact of the relatively high discharge when the 1986 photography was taken (Figure 6-9 A). Similarly, the apparent removal of small unit bars and minor trimming of larger medial bars throughout the reach may not assuredly be real (e.g. site 2). There was neither sedimentation nor erosion near Hamilton Bar (3) and the establishment of vegetation is an indication that this area was becoming increasingly stable. However, lowering of the riffle connecting Hopyard Bar with the downstream bar (4) indicates that the right bank channel was becoming more active. Eroded sediments were deposited on outer Hamilton Bar (5) causing compensating erosion with downstream deposition at the head of Gill Bar (6). Headward extension of Hopyard Bar drove sufficient flow towards Greyell Island that the old chute channel was reactivated (7) and enough material was removed to separate Greyell from Gill Bar complex. The large unit bar at Gill (6) caused compensating erosion at lower Hamilton Bar and directed the thalweg towards the right bank. Consequently, flow through the chute channel that once separated upper Carey Bar from Gill Bar was reduced, while sedimentation increasingly consolidated both features. Minor extension of vegetation on Gill bar complex illustrates that this site, much like Hamilton Bar, was largely inactive. Overall, these patterns of morphologic development imply that the reach was becoming increasingly stable, likely the consequence of reduced sediment influx. The thalweg remained hard against the right bank towards the head of Carey Bar, partly reactivating the minor secondary channel (8) but most of the flow was maintained in the existing channel. However, unit bar deposition at the outlet of the eroded chute channel between Greyell and Gill (9) forced the thalweg to migrate 300 metres westward, removing nearly one-third the  - 166-  Figure 6-9. Morphologic changes from (A) 1979 to 1986 and (B) 1986 to 1993. Shaded maps depict channel conditions on the early date; outlines show changes by later date. -167-  surface area of Carey Bar and directing flow toward Carey Point (left bank) where major erosion occurred (10). There were several small freshets from 1986 to 1993, with only the 1990 freshet exceeding the mean annual flood. Not surprisingly, there was considerable extension and consolidation of establishing vegetation throughout the reach as a result (Figure 6-9 B). There is also the appearance of substantial sedimentation, but the actual magnitude is tempered by a decrease in discharge from 2500 m /s to 580 m /s on the respective dates of photography. There was extensive deposition on 3  3  the right margin of Big Bar (1), but this did not cause compensating erosion because of extensive right bank armouring. Similar sedimentation appears at the head and left margin of Hopyard Island (2), but the absence of compensating erosion at the unprotected Greyell complex suggests that these sediments were more likely exposed than deposited as new unit bars. In fact, much of the flow had become focussed in the right bank channel, as a riffle between Hopyard and Big Bars reduced flow into the left bank channel. This development promoted siltation within the channel between Greyell and Gill complex, effectively rejoining these two features. However, because the chute channel was only partly infilled (a smaller, but topographically distinct, channel remained along much of the right bank of Greyell Island), re-activation during a large flood remains probable. Vegetation growth also continued at Hamilton B