Non-invasive measurement of fluvial bedload transport velocity using an acoustic Doppler current profiler by Colin D. Rennie B.Sc.(Eng.), The University of Guelph, 1995 M.A.Sc, The University of British Columbia, 1998 A thesis submitted in partial fulfillment of the requirements for the degree of Doctor of Philosophy in The Faculty of Graduate Studies (Civil Engineering) We accept this thq£i) as qonfoi/n)r$ to the required standard The University of British Columbia October 2002 © Colin D. Rennie, 2002 In presenting this thesis in partial fulfilment of the requirements for an advanced degree at the University of British Columbia, I agree that the Library shall make it freely available for reference and study. I further agree that permission for extensive copying of this thesis for scholarly purposes may be granted by the head of my department or by his or her representatives. It is understood that copying or publication of this thesis for financial gain shall not be allowed without my written permission. Department of cC/'">/ ^>\? /Vi e

30 min ADP sample). Error bars represent precision (+ standard error), with standard error of bedload based on variability of VuV samples 57 Figure 1.9: Bedload transport rate versus mean primary bed speed for individual concurrent samples without sand transport 58 Figure 1.10: Mean primary bed speed versus shear stress (O, included in regression; • , outlier excluded from regression (see text)). Error bars represent precision (± standard error), with standard error of shear stress based on standard error of u* from log-law fit 59 Figure 1.11: Influence of measurement depth on variance of va 60 Figure 2.1. Single ping measurements of apparent bed velocity (cm/s) for a trial with expected bed speed of 6.5 cm/s and expected direction of 180° (Trial 12). Speed radii for 23, 47, and 70 cm/s are shown 81 Figure 2.2. Measured versus expected apparent bed velocity: a) mean velocity component in direction of transport, b) mean magnitude of velocity vector. Precision bars indicate ± one standard error. Dashed line is perfect agreement. Solid line is best fit weighted regression. Symbols indicate orientation of transport with respect to the ADP x axis (see text), and if the artificial bed was present: • no transport, • 180° with bed, • 0° with bed, • 198° with bed, O 18° with bed, • 214° no bed, A34° no bed 82 Figure 2.3. Measured versus expected direction of transport. Precision bars indicate ± one standard error. Dashed line is perfect agreement. Solid line is best fit weighted regression. Symbols: + with artificial bed, X without artificial bed 83 Figure 2.4. Increasing standard deviation with increasing transport. Symbols: • measured standard deviation, O estimated instrument noise. Solid line is best fit regression to measured: ( y = 1.28x - 0.19; r2 =0.93). Dashed line is best fit regression to instrument noise: (y = 1.26x - 0.33; r = 0.91) 83 Figure 2.5: Running mean and coefficient of variation of the magnitude of the bed velocity vector for Trial 12 84 Figure 2.6: Running mean and coefficient of variation of the direction of the bed velocity vector for Trial 12 84 Figure 2.7: Probability density function and lognormal fit (smooth line) for magnitude of the bedload velocity vector for Trial 12 85 Figure 2.8: Probability density function for bedload velocity resolved in the direction of travel for Trial 11 85 Figure 2.9: Beam 1 velocity for a portion of Trial 11/12. Expected beam velocity based on meanfm. 86 Figure 2.10: Beam 2 velocity for a portion of Trial 11/12. Expected beam velocity based on meanfm 86 Figure 2.11: Beam 3 velocity for a portion of Trial 11/12. Expected beam velocity based on meanfm 87 Figure 3.1: Model procedure and results for optimized compound Poisson-gamma fit to Fraser River data set 120 Figure 3.2: Deconvolved distribution of within beam spatial average bedspeed for the Fraser River and Norrish Creek data sets, using both an optimized compound Poisson-gamma and a gamma distribution. Both linear and log scales shown 121 Figure 3.3: Norrish Creek measured data and fitted model using optimized compound Poisson-gamma distribution 122 Figure 3.4: Comparison of measured distribution to gamma distribution for run L2 bedload transport rate of Kuhnle and Southard (1988) -; 123 Figure 3.5: Estimated gamma distributions of entrained particle velocities from cPg model for Fraser River and Norrish Creek data sets. Both linear and log scales shown 124 Figure 3.6: Spectral densities of primary bedload velocity and near bed water velocity components. Mean spectral density indicated by dashed line, and 95% confidence interval indicated by gray lines 125 Figure 4.1a: Sand bed reach depth (m), based on flow depth at beginning of data collection. Bed sediment sample locations marked with +, and labelled with d50.. 164 Figure 4.1b: Sand bed reach spatially block averaged depth average velocity in a) the first hour of measurement, and b) the second hour of measurement. Tail of vector arrow is centre of block 165 Figure 4.1c: Sand bed reach bed shear velocity (m/s) distribution in first hour of data collection. Location of data collection points shown with small + symbols 166 Figure 4. Id: Sand bed reach raw bedload velocity (m/s) from first hour of data collection, prior to stationary bedload sampling 167 xi Figure 4.1e: Sand bed reach block average bedload velocity (m/s) in first hour of data collection overlain on depth (m) contours. Centre of vector arrow is centre of block. 168 Figure 4.1: f) Sand bed reach kriged bedload velocity distribution in first hour of data collection. Helley-Smith bedload sample locations marked with +, and bedload samples from second hour of data collection marked by * with value labeled, g) Sand bed reach kriged near-bed water velocity distribution for first hour of data collection 169 Figure 4.1.h: Sand bed reach (Sea Reach, Fraser River) bedload transport rate versus ADP bedload velocity calibration curve. Two outliers (denoted by O) were not included in the regression: it was presumed that these large bedload samples were due to bottom dragging of the Helley-Smith sampler 170 Figure 4.2a: Gravel bed reach depth (m) 171 Figure 4.2b: Gravel bed reach shear velocity distribution (m/s). Location of data collection points shown with small + symbols 172 Figure 4.2c: Gravel bed reach spatially block averaged depth average water velocity (m/s). Tail of vector arrow is centre of block 173 Figure 4.2d: Gravel bed reach raw bedload velocities 174 Figure 4.2e: Gravel bed reach block averaged bedload velocities (m/s) overlain on depth (m) contours. Centre of vector arrow is centre of block 175 Figure 4.2: f) Gravel bed reach kriged bedload velocity distribution. Location of bedload sampling from 2000 freshet marked by + 176 g) Gravel bed reach kriged near-bed water velocity distribution 176 Figure 4.3a: Sand-bed reach June 22nd raw bedload velocities 177 Figure 4.3b: Sand-bed reach June 22nd kriged bedload velocities 177 Figure 4.3c: Sand-bed reach June 22nd kriged bedload velocity, with boat track shown by red + symbols: a) east component, b) north component 178 Figure 4.3d: Sand-bed reach June 22nd bedload velocity vector angle versus boat trajectory angle 179 Figure 4.3.e: Bedload velocity angle (counter clockwise angles with East at zero) versus ADP heading angle (clockwise angles with North at zero) 179 X l l Figure 4.3f: Sand-bed reach June 22nd bedload velocity vector angle versus boat trajectory angle after adjusting heading by -10° 180 Figure 4.3.g: Bedload velocity angle versus ADP heading angle after adjusting heading by-10° 180 Figure 4.3h: Sand-bed reach June 22nd kriged bedload velocities following -10° rotation of raw velocities for ADP heading error 181 Figure 4.4.a: Minto side channel bedload velocity vector angle versus boat trajectory angle 182 Figure 4.4.b: Minto side channel bedload velocity vector angle versus boat trajectory angle after adjusting heading by -8.5° 182 Figure 5.1: Mean magnitude of bedload velocity vector versus mean bedload transport rate for Fraser River gravel reach data of Chapter 1 191 Figure 5.2: Standard deviation versus mean of primary bedload velocity for Fraser River gravel reach data of Chapter 1 192 Figure 5.3: Probability density function and beam gamma error model fit for bedload velocity resolved in the direction of travel for Trial 11 of Chapter 2 laboratory data. 193 Figure C. 1: Fraser River at Agassiz-Rosedale Bridge, facing upstream. Former WSC gauge site 08MF035. Width is 510 m. Described as a wandering gravel-bed river (D50,surf = 42 mm), with a slope of 4.8x10-4 (McLean et al. 1999). If bridge spans are numbered from left to right (facing upstream), the thalweg flows through span 5, and a shallow bar occurs at span 2. All 2000 freshet bedload sampling stations were downstream of the bridge at spans 3 and 4 (Chapters 1,3, and 4) 211 Figure C.2: a) ADP, b) Vi size VuV bedload sampler, and c) Helley-Smith bedload sampler 211 Figure C.3: Boat for Fraser River gravel-bed data collection (Chapters 1, 3, and 4).... 212 Figure C.4: Laboratory set-up for synthetic bedload. Note: an improved artifical bed was implemented for Trials 7, 8, 11, 12, 17, and 18. (Chapter 2) 212 Figure C.5: Norrish River at Hawkins Pickle Road Bridge, facing downstream from bridge, February 22, 2002. WSC gauge site 08MH058 (Chapter 3). River width at bridge is 28 m.. 213 Figure C.6: Norrish Creek on Feb 22, 2002. ADP securely mounted to bridge deck (Chapter 3) 213 X l l l Figure C.7: Boat for Fraser River sand bed Sea Reach data collection (Chapter 4) 214 Figure C.8: Fraser River sand bed Sea Reach, facing upstream (Chapter 4) 214 Figure D.l: Model procedure and results for estimated compound Poisson-gamma fit to Fraser River data set 216 Figure D.2: Model procedure and results for optimized compound Poisson-gamma fit to Fraser River data set 217 Figure D.3: Model procedure and results for beam gamma (spatially averaged bedload velocity within a beam area is gamma distributed) fit to Fraser River data set 218 Figure D.4: Model procedure and results for total gamma (linear combination of spatially averaged bedload velocity from each beam is gamma distributed) fit to Fraser River data set 219 Figure D.5: Model procedure and results for estimated compound Poisson-gamma fit to Norrish Creek data set 220 Figure D.6: Model procedure and results for optimized compound Poisson-gamma fit to Norrish Creek data set 221 Figure D.7: Model procedure and results for beam gamma (spatially averaged bedload velocity within a beam area is gamma distributed) fit to Norrish Creek data set.... 222 Figure D.8: Model procedure and results for total gamma (linear combination of spatially averaged bedload velocity from each beam is gamma distributed) fit to Norrish Creek data set 223 L i s t o f S y m b o l s The following symbols are used in this thesis: ADP = Acoustic Doppler Profiler; Ai = projected planar area of the z'th mobile grain; at = coefficient in error model; aDcp = acoustic Doppler current profiler; Bt = random variable for actual beam velocity for beam i; bi(t) = actual beam velocity for beam i; C = concentration; c = speed of sound; Ci = coefficient in error model; cPg = compound Poisson gamma distribution; cv = coefficient of variation; D = random variable for measured data; D = particle diameter; Di = particle diameter; DCh = characteristic particle diameter; DGPS = differential global positioning system; D# = #th percentile particle diameter by mass; 2>e,az = major, azimuthal diameter of individual beam sampling area; jbeilr = secondary, transverse diameter of individual beam sampling area; 14 = sampling area diameter; d = vertical distance between the bed elevation and the aDcp; dio = 50th percentile particle size of subsurface or of supply; da = depth of the active transport layer; dai - depth above the bed of the ith mobile grain; dr = distance from centre of rotation to transducer; F = acoustic operating (carrier) frequency; Fd = Doppler frequency shift; FR = Fraser River; \f\ = scattering form function; f# = probability density function for variable #; fm = percentage of the bed mobile within a beam; fmi = percentage of bed surface occupied by particle i; fR = received frequency; fs = sent frequency; {/}"* = n-fold convolution of random variables; Ga(u,a2) = gamma distribution with mean u and variance a 2 ; g - gravitational acceleration; gb = local bedload transport rate per unit width; HS = Helley-Smith bedload sampler; h = depth above the mean bed elevation; i = individual beam or individual particle; KS = Kolmogorov-Smirnov test of difference between distributions; ks = bed roughness; lp = bottom track pulse length; mi = mass per unit bed area of particle /; N = number of pings; 2 2' N(p,,a ) = normal distribution with mean \i and variance a ' NC = Norrish Creek; Nbi = random variable for random noise error in beam i; Np = number of particles within view of each aDcp beam; Ny = random variable for noise in the direction of transport; n = number of samples; n0i(t) = random noise error in beam i; 0(#) = of order #; PDT = Pacific daylight savings time; p = Tweedie parameter; p = probability r2= coefficient of determination; raDcP = horizontal distance from center of aDcp head to center of a transducer; n = particle radius; rt = transducer radius; SE = standard error of estimate; Ss = specific gravity of the sediment; std = standard deviation; t = time; T = pulse length; u = mean depth averaged primary water velocity; u5 = estimated mean primary water velocity at 5 cm above the bed; u* = shear velocity; V = random variable for spatially averaged bedload velocity; V = random variable for model of measured data; Vi - random variable for spatially averaged bedload velocity in beam i; Vt = volume of the bedload layer above the bed at the location of the z"th mobile th (including the i particle volume); Vsi = volume of the z'th mobile grain; v = spatially averaged bedload velocity; va = apparent velocity of bedload; vad = mean bed speed resolved in the direction of particle movement; vm = the mean magnitude of the apparent bed velocity vector; vBT = boat velocity determined by bottom tracking; voi = beam velocity along beam i; v0i(t) = measurement of bedload velocity along beam i; VDGPS ~ boat velocity determined by use of DGPS data; VE = east velocity; ve = expected average bed velocity; v,- = actual spatial average bedload transport velocity within the ith beam; VAT = north velocity; vp = actual average bedload velocity; vPi = velocity of bedload particle i; vr = rotational velocity; = bottom track velocity error due to changing tilt; vv = virtual velocity of particle; vxi = velocity in x direction through beam i; vyi = velocity in y direction through beam i; vx = forward velocity; v = transverse velocity; XYZj = coordinates of impingement point for beam x = coordinate in forward direction; x = boat displacement; y = coordinate in transverse direction; z = coordinate in vertical direction; Az = depth cell size; a = variable in gamma distribution; a = significance level, probability of Type I error; a = roll angle; P = variable in gamma distribution; P = pitch angle; X = transducer pointing angle from the vertical; X2 = chai squared distribution; (j) = half beam width angle of central lobe; <|) = beam deflection angle from vertical; xix n = parameter in a Tweedie model; 9 = counter-clockwise angle from the x-axis to the direction of bedload transport (Chapter 3); cp = heading error counterclockwise to the true coordinate system when facing down (Chapter 1); K = von Karman constant; X = dune length; X = variable in Poisson distribution; Xa = porosity of the active transport layer; Xai = porosity of the bedload layer for the z'th mobile grain; (j. = mean; 6 = half-intensity (-3 dB) half beam width; Gexp = actual direction of transport; Qobs = observed mean direction of transport; p = correlation coefficient; p = density of water; ps = density of sediment particles; pv2 = correlation coefficient of vector fields of 2 dimensional vectors; Ey = 2x2 covariance or cross-covariance matrix of vectors / and j; a 2 - variance; CJD = standard deviation of measured va data; GINST = standard deviation of vfl due to instrument noise; XX from the central axis. The angle 2 defines the beam width of the central lobe. Similarly, the beam intensity diminishes to half the central intensity (-3 dB) at an angle 0from the central axis, and the 8 angle 2 # defines the half-intensity beam width. Small magnitude sound intensity, which is defined as side lobe intensity, is also observed at angles greater than from the central axis. For our 1.5 MHz ADP, the half-intensity beam width is about 1.5° (V. Polonichko, 2000, SonTek Inc., personal communication). The reported beam width compares favourably with theoretical-numerical estimates of 1.0° and 1.2° for half intensity and main lobe beam widths for a 1.5 MHz, 0.03 m radius circular piston transducer (Medwin and Clay, 1998, p. 140). The reported beam width may be slightly larger than expected due to minimization of side lobe intensity. Typically, minimizing side-lobe intensity will increase the main-lobe beam width (R. Pawlowicz, 2000, Earth and Ocean Science, UBC, personal communication). SonTek claims to have minimized side lobe intensity through careful design of their transducers (SonTek 1998), and our ADP did not appear to generate errors related to beam side-lobes. Perhaps, the relatively weak side lobes, if present, were attenuated by the turbid Fraser River water. Significant backscatter from the side lobes will obscure the sampling area, as areas of the bed outside of the beams will contribute to the signal. However, bottom track side lobes will be received at a different delay time than the main lobes, thus signal processing should be able to remove side lobe effects while bottom tracking. Further, side lobe intensity is much less than main lobe intensity and, because both lobes are scattering off the same boundary while bottom tracking, the side lobe backscatter, if present, should be negligible. The non-linearity of sound intensity within the main lobe of the beam complicates assessment of the sampling area. The center of each beam contributes greater intensity than the outer portions of the beam. For simplicity, we have assumed the sound intensity to be uniform within the half-intensity beam width, i.e. a rectangular pulse, and negligible elsewhere. With this assumption, trigonometry determines the sampling area to consist of three discrete, approximately elliptical areas of uniform intensity. Assuming that the ADP is not tilted and the bed surface is perfectly flat, the quasi-ellipses are centered at horizontal distances of [d tan x + raDcP], which equals [0.47c? + 0.065] in our case, along relative azimuth angles of 120° from the ADP. The parameter raDCp is the horizontal distance between the center of the aDcp head and the center of a transducer (0.065 m). The length of the major, azimuthal axis of each quasi-elliptical area is found to be: ^ d + r.sinx d — r, sin y jbeaz = - - — +2rtcosx eflZ tan(90-x-0) tan(9O~x + 0) = 0.032c? + 0.066 (m) The length of the minor axis of each quasi-elliptical area is: (1.2) 2) =2 tan 9 + r. cosx = 0.029c/ + 0.060 (m) (1.3) Finally, the sample space can be defined as a circle that encloses the insonified sampling areas, which has a diameter: 2> d + rt sin % tan(90-%-Q) = 0.96d + 0.197 (m) -rrtcosx + raDcp (1.4) It should be noted that beam width is usually specified assuming a point source, whereas (1.2)-(1.4) include consideration of the transducer diameter. 10 The sampling area is actually a sampling volume, as suspended particles can act as acoustic scatterers if reflected and forward propagating portions of a pulse are coincident. The sampling volume depends on the sampling area and the bottom track pulse length (lp), where lp equals the duration of the acoustic pulse times the speed of sound. The size of the sampling volume varies dynamically while the pulse reflects off the bed (Figure 1.2). Importantly, the height of the sample volume varies during an acoustic pulse (ping). The sampling volume for a three-beam aDcp consists of three independent volumes located at the impingement points of the three beams. In order to determine the height above the bed of the sample volume, we have assumed that the backscatter receive window occurs only for returns which contain bottom reflections from both the leading and trailing edge of the acoustic pulse (Figure 1.2); this provides the strongest and most unambiguous return signal (Gordon 1996). A further assumption is that at any moment during the receive window, the sample volume occurs where the reflected and forward propagating portions of the pulse are coincident, as scatter of the forward propagating beam from suspended particles will be received at the same time as reflections from the bed. Figure 1.2 shows the dynamic change in the sample volume while a ping is reflecting off the bed. The maximum vertical extent of the sampling volume for a beam is [ 0.5 lp sin(90 ] above the bed. We used a short pulse length of 0.20 m in the present measurements, which minimized the sample volume to within 9 cm of the bed. However, the short pulse length may have increased variance of the velocity estimates (Brumley et al. 1991), although perhaps not very much (see Section 5.2.2). 11 Thus, each measurement of apparent bed velocity is a spatial average based on a sampling of three distinct areas/volumes, which are spread apart within a circle that has a diameter approximately equal to the distance between the bed and the ADP. Only the velocity component parallel to the beam can be measured within each area. The size of the individual quasi-elliptical areas increases linearly with depth due to beam spreading, but the individual areas never overlap regardless of the depth. The bottom track pulse length affects the degree to which suspended scatterers can contribute to the signal. 1.3.2 Preferential particle size It may be that the ADP preferentially registers particles of a particular size or velocity. Acoustic theory suggests that our ADP may preferentially register particles that have sizes greater than or equal to coarse sand (> 0.6 mm). The sound backscattered by a particle is proportional to its form function (| /1). The form function is defined relative to the particle size, thus higher values of | /1 indicate that a particle displays greater backscatter than its size would suggest. Values for | /1 depend on kr{, where k is the wavenumber of the acoustic pulse (6300 m"1 for our 1.5 MHz instrument) and r,- is the particle radius. Theoretical values of | /1 for a rigid sphere are constant at unity for values of kr(» 1 (specular scatter due to reflection), and less than unity with a dependency on {kri)2 for krt« 1 (Rayleigh scatter due to diffraction) (Medwin and Clay 1998, Figure 7.5.4). For 1 < krt < 20, | /1 oscillates about unity due to interferences between reflected and diffracted wavelets. However, this theoretical derivation of | /1 assumes an isolated rigid sphere, whereas bedload particles are irregularly shaped and 12 display some acoustic elasticity, and bottom tracking backscatter will be complicated by interactions with the boundary. Theoretical predictions have been made for backscatter from a hard spherical object near a hard flat bottom (Gaunard and Huang 1996, Bishop and Smith 1999). For a 25° angle of incidence of the acoustic pulse, particles situated within a particle diameter of the bed display increased amplitude oscillations of | /1 in the interference region (2 < krt < 20) due to interactions with the bottom. There has also been some experimental study of the form factor for isolated irregularly shaped quartz particles (Hay 1991, He and Hay 1993, Thorne et al. 1995). The Rayleigh scattering region was observed for krt < 2. Increased oscillations of | /1 in the interference region were observed, but these were highly dependent on the angle of incidence of the acoustic pulse. An average | /1 for several angles of incidence showed little oscillation (Thorne et al. 1995). Thus, the form factor for an aggregate of irregular particles can be assumed, on average, to be a uniform value of unity for kr{ > 2 and less than unity for krt < 2. For our 1.5 MHz ADP, krt = 2 corresponds to a particle diameter of 0.6 mm. Thus, bottom track backscatter from medium and fine sands (A < 0.6 mm) may have been de-emphasized with respect to their size, and negligible from fine silts and clays (D,- < 0.02 mm). It remains unclear, however, if average | /1 will approach uniformity for irregular particles that are near the bed. 13 1.3.3 Instrument Error Error due to the instrument for bottom tracking over solid, stable substrates is an order of magnitude less than error for water column velocity estimates (V. Polonichko, 2000, SonTek Inc., personal communication) because a strong backscatter signal is received from the solid boundary. The standard deviation of the horizontal water velocity measurements is reported to be (Theriault 1986a, SonTek 1998): 140 c 0 \ v = — — 7 = (1-5) F Az V/V where c is the speed of sound (nominal value 1500 m/s), F is the acoustic operating frequency (1,500,000 Hz), Az is the depth cell size (0.25 m), and N is the number of pings. The number of pings equals the averaging interval (5 s) times the pinging rate (9 pings per second for the 1.5 MHz ADP). Substitution of these values yields crw equal to 8 cm/s, which suggests that the standard deviation of bottom tracking velocities due to instrument error for a population of five second ensembles should be on the order of 1 cm/s. This estimated value compares well to laboratory towing tank tests of the bottom tracking capabilities of our instrument (Table 1.1). These tests were performed in the towing tank described in Chapter 2, with the instrument level, positioned at the centre-line of the tank, and attached to the towing carriage. In all but two of the tests, the instrument x-axis (beam 1) was aligned with the direction of travel (following a 2° rotation correction in post-processing for instrument misalignment). The data were collected in xyz coordinates, and thus the internal compass was not used. It is noteworthy that the measured bottom tracking, while remarkably accurate, displays a slight bias towards 14 underprediction (Figure 1.3). This will translate into a positive bias in va (see Equation 1.1), although the bias is of 0(mm/s) and thus is insignificant. However, bottom tracking error is greater when the bed is mobile. There are two sources of instrument error in a three-beam aDcp: 1) instrument noise, i.e. inability of the aDcp to determine the average Doppler shift within a beam during a single ping, 2) error inherent to the averaging process for a three-beam aDcp if the transport is heterogeneous amongst the three beams. In our field testing on Fraser River, gravel bedload transport was thought to be patchy or sporadic across the bed. Bed particle velocities within the sampling area would have been variable, including zero velocities. The nonuniformity of bedload transport velocities would have caused instrument error due to both instrument noise and heterogeneous transport. Instrument noise is increased because the internal signal processing is not presently designed to resolve the velocity of a mobile bed. If a mono frequency sound pulse insonifies a bed that has nonuniform velocity, a spectrum of frequencies will be reflected as different portions of the pulse undergo varying amounts of Doppler shift. The SonTek bottom tracking algorithm is proprietary, but involves a pulse-to-pulse incoherent ("narrowband") technique (V. Polonichko, 2000, SonTek Inc., personal communication), thus it seems likely that a single, strong amplitude frequency is found and used to determine the Doppler shift (see Brumley et al. (1991)). The wideband echo from a mobile bed would not be processed well by such a technique. Signal processing that considers the entire spectrum may produce more accurate velocity estimates (e.g. Hansen 15 1986). Bedload velocity error would also be reduced by an increased bottom track pinging rate, which could be achieved if water column pinging was reduced or eliminated. Optimization of the bottom track pulse length may also reduce the velocity variance. Finally, use of broad-band Doppler sonar techniques would reduce the instrument noise (Brumley et al. 1991), at the expense of the introduction of a maximum unambiguously measurable velocity. Due to the geometry of a three-beam aDcp, some measurement error inherently arises when bedload varies between the three beams. Only the velocity component parallel to the beam can be determined in each insonified area. Theriault (1986b) presented equations to calculate the forward (vx) and transverse (v ) velocities given the velocity determined along each beam, assuming homogenous vertical velocity in all three beams (see Equation 2.6). Note that towing tank tests of ADP bottom tracking with data collected in beam coordinates confirmed Equation 2.6. Further consideration of the 3-beam geometry (see Equation 2.8) yields predictors based on the forward (v*,) and transverse (v,,,-) velocities through each beam: v, = [2 v„ + 0.5 (v,2 + v„ ) + 0.866 (v„ - vy2)] / 3 (1.6a) v, = [0.866 (v„ + v „ ) + 0.5 (v 3 - v J ] / V3 (1.6b) Beams are numbered counterclockwise when the ADP is facing down, with the horizontal component of beam 1 parallel to the x-axis (Figure 1.1). With a three-beam ADP, due to the unequal weighting of vxt in the estimated \>x, heterogeneous vxi will necessarily produce an incorrect vx. Furthermore, vx and vy are contaminated by heterogeneous velocities in the opposite component (Theriault 1986b). During sparse 16 transport when particle velocity differs between the three insonified areas, the direction and magnitude of the measured velocity vector will necessarily be wrong. Even in the ideal case with uniform particle size, particle velocity, and percentage of the bed mobile (fm) in all three of the beams, partial transport can produce heterogeneous average particle velocity between the three beams, which will also result in incorrect resolution of the bedload velocity vector. This can arise when individual particles are large with respect to the beam width, and occupy a substantial portion of a beam. The percentage of the bed mobile within a beam will vary as particles move into and out of the view of the beam. Heterogeneous transport can create error variance in the forward and transverse directions, with maximum possible direction errors ranging from 60°, which can occur if transport is purely along the x axis, to 90°, which can occur for transport along the y axis. In theory, it may be possible to account partially for this error by specifying the expected direction of transport. As we do not know the actual spatial variability of transport during our field tests, we can not estimate the instrument error due to heterogeneous transport. We can, however, provide an estimate of instrument noise. Preliminary laboratory tests where particle speed was controlled and variability infm was accounted for suggest that instrument noise with a mobile bed has a standard deviation of about: where vov is the mean magnitude of the apparent bed velocity vector and N is the number of bottom tracking pings in an ensemble average (see Chapter 2). However, vm may (1.7) 17 itself be a measure of the variance of the data. Refer to Section 5.2.1 for further discussion of this point. Regardless of the uncertainty in Equation 1.7, there is a great deal of variance due to instrument noise, and instrument noise is the dominant source of the measured variance (see Section 1.6.3). 1.3.4 DGPS error The DGPS was of sub-meter precision. However, only differences in sequential position estimates are relevant to velocity estimates. A 0.2 Hz, 15 minute time series for a stationary position had average and maximum velocity vector errors of only 0.7 cm/s and 2 cm/s, respectively, with no directional bias. 1.3.5 Heading error A heading direction from a compass is used to resolve aDcp measurements into earth coordinates. Systematic differences between VBT and VDGPS have been noted previously when there is misalignment between the orientation of the transducers and an external compass (Joyce 1989), which causes incorrect resolution of the east and north velocity components. We used an internal compass/tilt sensor (Precision Navigation TCM2), thus we should not have any misalignment errors. However, the accuracy of the compass was reported to be ± 2°. If this is also a measure of the compass precision, then heading error may have contributed random noise to our bottom tracking measurements. A stationary test of the ADP compass using 5-second ensemble averaging yielded better heading precision: the range of data was only ± 0.6° and the standard error was 0.01° (n = 67). 18 Error in bottom tracking velocity components due to heading error can be predicted (Joyce 1989). Velocity errors due to compass inaccuracy in the east (AVE) and north (AVM) directions are: AvE = (coscp - \ )vE + sin

vE (1.8b) where cp is the heading error measured counterclockwise to the true coordinate system when facing down, and VE and v# are the true east and north velocities A compass error of 2° can produce a maximum bedload velocity component error of 3.5% of the magnitude of the bedload velocity vector. 1.3.6 Error due to dynamic instrument tilt In Rennie et al. (2002), we estimated an error due to rocking of the instrument during beam transmission or reception. We now believe that we were mistaken, as the instrument can measure instantaneous velocities only in the radial beam direction. The original description of error due to dynamic instrument tilt is given below, but further discussion is based on our new understanding. Rocking of the instrument during beam transmission or reception will impart an angular velocity to the beam that could be interpreted as a bottom track velocity. Assuming random fluctuations in the pitch and roll of the instrument, this will create random bottom track velocity errors that will depend on the angular velocity of the instrument and the range of the beam (i.e. the depth to the bed). The magnitude of these errors can be 19 approximated by determining the change in the XY coordinate of the impingement point of each beam due to the change in instrument tilt during an ensemble of measurements. Changes in the values of roll (a) and pitch (P) measured by the ADP from one ensemble to the next were used to estimate bottom track error due to dynamic instrument tilt. The ADP coordinate system follows a right hand rule, with the z-axis positive upwards and the x-axis positive in the horizontal direction of beam number one (Fig. 1.1). For a downward looking ADP, positive roll is defined as a clockwise rotation about the x-axis, and positive pitch is a counterclockwise rotation about the y-axis, when looking out from the origin. All the roll values were made negative to compensate for this irregularity. For each ensemble, the ADP provides estimates of the vertical depth for each beam (dj). Coordinates of the points where the beams impinged on the bed were found to be: d\ tan(x) -d2tan(x) cos (60°) '- d3tan(x)cos(60°) XYZ] = 0 XYZ2 - d2tan(%)sin(60°) XYZ3 = - d3 tan(x)sin(60°,) ~d2 -d3 (1.9) Coordinates of the impingement points following the change in instrument tilt were found using a three-dimensional coordinate transform (Weisstein 1999, p.1580): XYZ, = 1 0 0 0 cos Act - sin Act 0 sin Act cos Act cos A P 0 sin A p 0 1 0 - sin A p 0 cos A P XYZ; (1.10) Values oiXi' and 7,-' were scaled by Z, /Z,-assuming a flat horizontal bed surface. The change in X and Y divided by the ensemble time increment provides a first order estimate of the bottom track velocity error due to changing tilt (v,xi, v^ ,). However, each beam will 20 only sense the error component parallel to the beam. Equation 1.6 yields the sensed error in x and y due to tilt. In our measurements, estimated errors in x and y for individual sampling stations had standard deviations of about 1 cm/s, with a maximum of 4.2 cm/s. As expected, larger error is incurred for stations with greater depth. The estimate described above assumed that a ping acts as a single coherent pencil beam, and thus will be smeared through space as it is transmitted, reflected, and received from a dynamically tilting instrument. The assumption was that curving the beam due to dynamic tilt will result in an observable Doppler shift. In reality, the ping is a series of acoustic wavelets, each of which will act independently. If the instrument is dynamically tilting, each individual wavelet will have a rotational velocity with respect to the instrument. A simple analogy is a pitcher throwing a series of balls from a rotating pitching mound. In the frame of reference of a stationary observer, each ball will fly straight in a different direction than the preceding ball. In the frame of reference of the pitcher, each ball will curve due to a rotational velocity. For stationary scatterers on the river bed, the wavelet will reflect with zero Doppler shift. However, when the wavelet is received by the instrument in the reference frame of the instrument there is a rotational velocity of the backscattered wavelet. Considering rotation about the instantaneous axis of rotation, the rotational velocity (vr, in m/s) will be the product of the angular rotation rate (co) and the distance from the centre of rotation to the transducer face (dr). The centre of rotation should be in the centre of gravity of the boat, assuming a rigid mount attaching the ADP to the boat. The instantaneous direction of vr will be perpendicular to dr. The magnitude of vr should be of 0(cm/s): assuming co equals l°/s i 21 and dr equals 2 m, vr will be 3.5 cm/s. In effect, the magnitude of rotational velocity reflects the degree of translation of the transducer due to rocking of the boat. The question is whether this rotational velocity is parallel or orthogonal to the instantaneous beam axis. If it is orthogonal, then no Doppler shift will be recorded, and dynamic instrument tilt will not result in a bottom tracking error. Clearly, the direction of vr depends on dr, and thus, the relative locations of the ADP and the instantaneous centre of rotation of the boat. If the transducer and rotational centre are in line the rotational velocity will be perpendicular to the transducer face and there will be no Doppler shift recorded. We can not make precise estimates of the magnitude of potential errors due to dynamic instrument tilt, as we do not have information on the relative locations of the centre of rotation of the boat and the ADP. In future measurements, it would be worthwhile to note these locations. 1.3.7 Bedload transport rate Our conceptual model of bedload transport consists of an active transport layer of moving particles which is immediately above the static bed surface. Finer material travels in suspension in the water column above the bedload transport layer. Particles within the active transport layer move either by saltating (short hops) or rolling. Larger particles tend to roll and smaller particles saltate, depending on the dimensionless shear stress (Abbott and Francis 1977, Andrews and Smith 1992, Hu and Hui 1996a). Dimensionless shear stress is defined as: 22 where x is the bed shear stress, p is the density of water, g is gravitational acceleration, Ss is the specific gravity of the sediment, D is the particle diameter, and u* is the shear velocity (x = pw,2). See Equation 1.20 for a method to estimate u*. Andrews and Smith (1992) estimated theoretically that saltation should begin to occur once x* exceeds 0.06. Flume experiments by Hu and Hui (1996a) showed that 80% of particles roll when x* < 0.08, while they tend to saltate or suspend when x* > 0.2, and travel almost exclusively in suspension when x* > 2.8. Suspension occurs if the ratio of particle fall velocity to bed shear velocity is less than about 0.8 (Bagnold 1973). The Hu and Hui (1996a) transport mode criteria were based on measurements of shear velocity and transport mode of individual grains. Thus, effects of particle packing and bed morphology were not considered. Particle packing can increase the shear stress required to initiate motion, but may not effect the mode of transport once motion is initiated. However, increased bed roughness due to dunes or particle clusters will increase the total shear stress without a comparable increase in shear stress imposed on bed particles. Drake et al. (1988) used motion picture photography in a natural fine gravel-bed river to observe particle movements. Bed shear velocity at the time of measurement was about 0.078 m/s. Particles smaller than 3 mm saltated (x* > 0.13). Particles greater than 3 mm usually rolled and particles greater than 7 mm always rolled (x* < 0.054). These observations are consistent with those of Hu and Hui (1996a), thus we will adopt the Hu and Hui criteria. We have observed u* in gravel bed rivers during bed mobilizing flows ranging from 0.1 to 0.2 m/s. Thus, sands will saltate in the active transport layer of gravel bed rivers. Gravels will usually roll (D > 8 mm), but gravels may saltate in the highest imposed 23 flows (D < 30 mm). Materials finer than fine sands (D < 0.2 mm) will travel in suspension. Typical shear velocities we have observed in sand-bed channels during flood flows range from 0.04 to 0.08 m/s. Kostaschuk and Villard (1996) reported u* values as high as 0.19 m/s in the Main Channel of Fraser River. However, these values are for total bed shear velocity, including form resistance due to dunes. Segmentation of their velocity profiles suggests that grain shear velocity may have been about half the total shear velocity (Villard and Kostaschuk 1998). Using a range of u* from 0.04 m/s to 0.08 m/s and the criteria of Hu and Hui (1996a), particles finer than silts (D < 0.03 mm) will always, and fine sands may sometimes, travel in suspension. Sands will either saltate or travel in suspension, but will not roll. For bed slopes typical of natural streams, fluid momentum does not transfer below the bed surface, thus only surface particles are transported and the active transport layer is thin (Bagnold 1973). The depth of the active layer for rolling particles is one particle diameter. For bedload by saltation, the depth of the active layer is limited by the maximum height of saltation. Saltation trajectories have been observed to have heights ranging from 1.2 to 9.2 particle diameters, with a typical value of about 3 particle diameters (Abbott and Francis 1977, Lee and Hsu 1994, Hu and Hui 1996a, Lee et al. 2000). 24 In a sand-bed channel where transport is by saltation, the active layer depth should be about 3Dgo thus the active layer depth is of O(mm). D90 is defined such that 90% of a sample mass is comprised of particles smaller than the D90 particle size. During high transport rates in a sand-bed it is conceivable that particles below the initial surface are entrained following entrainment of a surface particle and prior to deposition of another particle. Thus, the active transport layer may consist of more than surface particles. Definition of the active bedload transport layer is complicated by the presence of a steep gradient in suspended sediment concentration within the water column. In practice, conventional measurements of bedload in sand-bed rivers have been conducted with samplers that have orifices of O(cm) depth . For example, a standard Helley-Smith sampler has a 7.6 cm by 7.6 cm square orifice. Thus, conventional sampling captures both bedload proper and near-bed suspended load. Drake et al. (1988) observed that saltating particle heights never exceeded the highest points of the poorly sorted gravel-bed surface. For a gravel-bed, it appears that the maximum depth of the active transport layer should approach the size of the largest particles on the bed (Andrews and Parker 1987), as the active layer thickness is determined by the size of the largest rolling particles. Thus the active layer depth is of O(cm). However, assessment of the active layer depth is complicated by spatial and temporal heterogeneity of transport. In a poorly sorted gravel-bed, the rate and calibre of bedload transport is non-linearly related to shear stress (Parker et al. 1982). For typical shear stresses in gravel bed channels only partial transport occurs (Wilcock and McArdell 1997), wherein only a portion of the bed surface is active at any one time. At low shear 25 stresses only the smallest particles are mobilized (Jackson and Beschta 1982). As shear stress increases larger particles are mobilized and a greater proportion of the bed surface is active instantaneously. At the highest flows the largest alluvial particles are mobilized, and the active layer depth approaches the depth of the coarse surface layer (about D9o), but the entire bed surface is still not simultaneously mobilized (Haschenburger 1999, Rennie and Millar 2000). Thus, the spatially averaged instantaneous active layer depth should always be less than D90. In gravel-bed channels bedload transport is further complicated by the degree of sediment supply. For low supply channels the bed surface will armour, wherein the bed surface is composed of coarser particles than the subsurface and the supply (Kellerhals 1967). Transport tends to occur as entrainment and transport of individual clasts. For higher sediment supply channels transport has been observed to occur as bedload sheets, which are low relief bedforms (amplitudes of 1 to 2 grain diameters) consisting of alternate congested (coarse) and smooth (fine) zones (Dietrich et al. 1989). The aDcp measures the spatially averaged velocity of the bed surface. For mobile grains, this will include both translational velocity and rotational velocity. For rolling particles, the surface velocity is equal to both the translational velocity and the rotational velocity, which are the same. Saltating grains, however, may have both a translational velocity and an independent rotational velocity. Drake et al. (1988) observed that saltating particles had a strong rotational velocity immediately following initiation of motion or impacts with the bed, but the rotational velocity quickly decayed and diminished to weak wobbling about randomly oriented axes. Thus, it appears that saltating grains do not have significant additional surface velocity due to rotation, and rotational velocity will not be considered further in this thesis. See Appendix A for further analysis of this issue. The aDcp provides an apparent bedload velocity (va), whereas a bedload transport rate often is required. If we assume that va is an unbiased measure of the actual spatially averaged bedload velocity (vp), then the local bedload transport rate per unit width (gb) can be calculated kinematically if the depth (da) and porosity (Aa) of the active transport layer are known (Haschenburger and Church 1998): * * = v , « / . ( l - A . ) p f (1.12) where ps is the density of the sediment particles. It is important to recognize that vp is a spatial average velocity, which depends on both particle velocities and the percentage of the bed surface that is mobile (see Equation 1.14). Application of (1.12) to sand-bed channels is relatively straightforward, wherein da is the depth of the bedload saltation layer that can be assumed to be about 3Z)po and (l - Xa ) is the sediment concentration within the bedload layer (see van Rijn 1984a, Villard et al. in press). Similarly, in gravel-bed channels (1.12) can be applied to time-integrated assessment of bedload transport (Haschenburger and Church 1998), wherein vp is the virtual velocity of the bedload (total travel distance divided by total time including rest periods), da is the time integrated mobilization depth of the bed sediments, and Xa is porosity of the static bed. However, (1.12) is problematic for instantaneous measurements of spatially averaged gravel-bed bedload transport. In this case the active 27 layer is the instantaneous bedload transport layer above the static bed. If bedload occurs as partial transport wherein much of the bed surface is inactive at any one time, it is difficult to estimate a depth or porosity of the active layer. Recall that the percentage of the bed surface that is instantaneously mobile is accounted for in vp. Thus, both da and Xa must be averaged by areal weighting over only the mobile grains: I dai A, da=- (1.13) i Z Ki Ai Xa=-i (1.13b) a V-V • Ki=-±TJL (i-i3c) i where dai is the depth above the bed of the z'th mobile grain, At is the projected planar area of the z'th mobile grain, Xai is the porosity of the bedload layer for the z'th mobile grain, Vsi is the volume of the z'th mobile grain, and Vt is the volume above the bed at the location of the z'th mobile grain (including the z'th particle volume). In practice, neither the depth nor the porosity of the active layer is known explicitly, but it should be possible to estimate da and Xa within acceptable bounds. The instantaneous da should range from about 3D 5 (saltating smallest grains) to Do0 (rolling largest grains). Thus da will range from a few mm to several cm, and can be assessed based on the proportion of each grain size in the bedload. The porosity of a layer with individual grain movements will range from Xa = 7/^ = 0.78 (saltating spherical single grain with a saltation height of 3D) to 0.33 (rolling spherical grain). 28 Application of the kinematic model in this case is complicated. Alternatively, as will be shown, a calibration curve can be developed to relate gb to the measured va, provided that the transport rate has been independently sampled. Also, an alternative to the kinematic model can be developed to predict bedload transport rate from vp. Actual velocities of individual bedload particles (vpi) will vary, depending on particle size and shape, local bed roughness, and stochastic variability in fluid force. Furthermore, particularly in a gravel bed, some of the bed surface within the sampling area will be immobile. The actual average bedload velocity (vp) is a function of the actual particle velocities and the percentage of the bed surface occupied by each particle (fmi): v„=2>„./mi (1.14) Furthermore, a simple predictor of bedload transport rate is: gb =zZsbi = Z vPi m, (1.15) /' i where m,- is the mass per unit bed area of particles size A moving with velocity vpi. Assuming spherical particles and that transported particles do not overlap vertically within the active transport layer, it can be shown that: mi = ri Ps fmi (1-16) where n is the particle radius (A/2). Substituting (1.16) into (1.15) yields: Sb=%Ps }ZvPi fmi ri (1-17) i If the sediment is of uniform size such that r,- is constant, then from (1.14) and (1.17) gb=y3PsDvp (1.18) 29 where D is the uniform particle diameter. Thus, bedload transport rate can be estimated directly from the average bedload velocity i f the bedload particles are of uniform size. It may also be possible to define a "characteristic" particle size (DCh) for a poorly sorted bed. Comparison of (1.12) and (1.18) shows that, for transport of particles that do not overlap vertically in the active layer, Dch=l.5(l-Aa)da (1.19) 1.4 Methods 1.4.1 Study site We field-tested this technique during the 2000 freshet at the Agassiz-Rosedale bridge site of Fraser River in British Columbia, Canada (49.21° N , 121.78° W, Water Survey of Canada former gauge station 08MF035). The study site location is shown in Figure 1.4. The hydrology, bed material, and sediment transport characteristics of the site were described by McLean et al. (1999) on the basis of a 20-year program of sediment transport measurements carried out by the Water Survey of Canada. The river width is 510 m, and the channel gradient is 4.8X10"4. The D50 and D90 of the surface sediment were reported to be 42 mm and 80 mm, respectively. The D50 and D90 of the subsurface sediment were 25 and 80 mm. The bed material is bimodal, mostly gravel-cobble with some sand. The thalweg flows along the left bank, and a shallow mid-channel bar occurs towards the right side of the channel. This bar is the tail of an upstream mid-channel island. Photographs of the site and sampling methods are provided in Appendix B and in Rennie and Millar (2001). 30 1,4.2 Apparent bedspeed Measurements were taken at flows ranging from 5600 m /s to 6800 m /s, which were less than the mean annual flood of 8800 m /s, but greater than the 5000 m /s threshold for gravel mobilization (McLean et al. 1999). The 1.5 MHz ADP, with internal compass and tilt meter, was deployed from a boat. The ADP was interfaced with the DGPS. Refer to Appendix B for the ADP operating parameters. Bottom track measurements occurred at 1 Hz, and ensemble averages were collected at 0.2 Hz. ADP and DGPS data were collected at 31 nearly-fixed locations for time periods ranging from 2 to 112 minutes. Positions were held either by motoring or by tying to the bridge, although there was some boat motion within a restricted area. Flow was steady during sampling at a station. Poor bottom track data quality eliminated eleven of the ADP stations: we were incapable of bottom tracking in depths less than 2 m or greater than 8 m, which limited good stations to the center of the channel. Presumably, the bottom track acoustic signal was excessively attenuated in the turbid Fraser River for depths greater than 8 m. Two bottom track data quality criteria were utilized. First, ensemble averages were accepted only if they had greater than 33% good pings (which guaranteed at least two good pings in an ensemble average) and if the recorded depth was reasonable (see below). Second, greater than 50% good ensemble averages were required in a time series. A bedload velocity vector was calculated for each ensemble, and resolved in the direction of the mean primary water velocity (Bathurst et al. 1977). Finally, an average apparent primary bedload speed was determined for each station. 31 1.4.3 Bedload Concurrent bedload samples were collected at 10 of the locations using a half-size VuV sampler (Novak 1957) for the coarse fraction (> 4.75 mm), and a Helley-Smith (HS) sampler (Helley and Smith 1971) for fines (>0.147 mm, < 4.75 mm). The aperture of the VuV sampler was 255 mm wide by 115 mm high, and the wire mesh gap of the sampler was 4 mm. The aperture of the HS sampler was 76 mm by 76 mm, with a collection bag mesh size of 0.2 mm. Typically, five 5-minute VuV samples and three 3-minute HS samples were collected over the course of about an hour at each station. HS samples were not collected at two stations, so HS samples from a similar station were used to estimate the fine fraction. The samplers were deployed from the bridge using a rope and pulley system. The shipboard ADP was positioned within a few meters downstream of the sampler. The samplers could not be deployed in the thalweg due to excessive depths and high velocities. Sampling efficiencies were assumed to be 0.33 for the VuV sampler (McLean et al. 1999) and 1.5 for the fine fraction from the Helley-Smith (Glysson 1993). Average fractional bedload transport rates were determined for each station by averaging the multi-sample transport within each size fraction, from which average cumulative particle size distribution curves were developed (Figure 1.5). We also examined individual VuV bedload transport rate samples and corresponding concurrent 5 minute ADP bedload velocity averages. 32 1.4.4 Shear stress A local mean bed shear stress was also calculated for each station by a log-law fit to the mean primary water velocity profile. It was necessary to correct the depth recorded by the ADP. It appeared that the estimated depth for each beam was not corrected for tilt of the instrument. The counter-clockwise values of roll (a) and pitch (P) measured by the ADP for each ensemble average were used to resolve a corrected depth for each ensemble. Coordinates of the points where the beams impinged on the bed in the tilted coordinate system were found using (1.9). Coordinates of the impingement points in the tilt-corrected coordinate system were found using (1.10). If an individual beam depth diverged from the mean depth for the beam by greater than 0.5 m, or a depth was not recorded, the beam depths for that ensemble were omitted. A mean depth for the ensemble was found by averaging Z, from the three beams. Finally, a mean depth for the station was found by averaging all the ensemble depths. Water column velocities were collected along three beams in bins vertically spaced 25 cm apart, with velocity in each bin based on backscatter from a 50 cm deep triangularly weighted window (SonTek 1998). The first good bin was centered > 0.25 m from the bed. An average water velocity was determined for each bin above the mean depth for the station. The average velocities for each bin were resolved in the primary direction for the mean profile by minimizing the depth averaged secondary mean velocity (Bathurst et al. 1977). Finally, the shear velocity («») and the bed roughness (ks) were estimated using the log-law with the slope and intercept of the least squares linear regression to the entire profile: 33 u = —— ln{h) + — In K K \ k s J (1.20) where K is the von Karman constant (0.41), and h is the depth of the bin above the mean bed elevation. The mean local bed shear stress (x) is simply pu,2, where p is the density of water. A standard error of each u* estimate was calculated based on the regression fit (Zar 1996,p.330). 1.5 Results The boat velocities determined by bottom tracking and by DGPS were used to estimate an apparent bedload velocity for each ensemble using Equation 1.1 (Figure 1.6a). The bottom track boat velocities were integrated to yield a boat trajectory (Figure 1.6b), which could be compared to the DGPS position fixes. A time series of va resolved in the direction of the primary water velocity was also constructed (Figure 1.7a). The bedspeed was resolved in the direction of the primary water velocity under the assumption that the bedload samplers aligned in this direction. In .most cases there was little difference between mean primary bedspeed and the mean bedspeed vector, as indicated by the small magnitude of the mean secondary bedspeed (Table 1.2). Figures 1.6 and 1.7 provide data from one, typical sampling station (sample 07132). Wide scatter was evident in estimates of va from the individual 5-second ensemble averages (Figure 1.6a). The coefficient of variation (standard deviation divided by the long-term mean) for va resolved in the direction of the average bedload speed vector 34 ranged from 0.9 to 6.4 for the 20 sampling stations. However, the long-term average for a station was consistently in the expected downstream direction. Similarly, the boat trajectory by bottom tracking ran upstream, in contrast to the essentially fixed position of the boat by DGPS. For sample 07132, bottom tracking falsely indicated that the boat moved upstream 290 m over the course of 1.5 hours, whereas the boat was essentially stationary and DGPS indicated that the boat remained within 2.4 m of the starting position (Figure 1.6b). The apparent boat trajectory indicated by bottom tracking reflects bias introduced by bed mobility. It is apparent that a large sample is required to estimate reliably the mean va. Both a running coefficient of variation (cv) (Kuhnle and Southard 1988) and a running average of the apparent primary bedspeed sequence were calculated to determine the duration of sampling required (Figure 1.7b). For each sampling station, a Monte Carlo process of 400 randomized sequences yielded an estimate of the number of 5-second samples required to achieve a reliable estimate within ±5% of the long-term mean or the long-term cv. A randomized sequence was used to ensure that any temporal trend in the data did not influence the result. For all of the stations with greater than 30 minutes of data, about 25 minutes of sampling was required to achieve stable estimates of the mean and the cv. Unfortunately, it appears that short-term sampling is not presently dependable due to excessive noise in the data. Similarly, assessment of temporal variability of at-a-station bedload for scales of seconds to minutes is not presently possible (however, see Chapter 3 for one such assessment). Quasi-periodic pulsing of bedload in plane-bed gravel rivers at periods ranging from 5 to 30 minutes is widely reported (see Gomez et al. 35 (1989) for review). The bedspeed time series of Figure 1.7a may display an irregular pulsing with a period of between 4 and 5 minutes. However, the trend is obscured by noise in the data. Spectral analysis of the bedspeed time series did not reveal any significant periodic or quasi-periodic pulsing. Again, noise in the data may have obscured temporal pulsing, but it could also be that the pulsing was not sufficiently periodic to emerge spectrally. Despite wide scatter in individual 5-second ensembles, the estimates of mean v0 appeared to be coherent when compared with bedload sampler data. The average apparent primary bedspeed was well correlated with mean bedload transport rate determined from-the bedload sampler data (Figure 1.8): Yb = 1.4 v~a- 0.046 (r2 =0.93, SE = ±0.0032, p< 0.0005 ) (1.21) Overbars indicate mean quantities, and units are SI. The linear regression was weighted (Montgomery and Peck 1982, pp. 362-363) by the reciprocal of the standard error of each v0 in order to account for the variable precision in va estimates, which was due to variable durations of sampling and increasing variance with increasing transport rate. SE is the standard error of estimate based on a weighted average va and a weighted sum of i residuals, andp is the probability that the slope is zero (i.e. zero correlation). We have plotted va on the abscissa in order that (1.21) is comparable to (1.12). It is apparent from (1.21) that a calibration curve can be developed for this site to relate va measured using an ADP to the mean bedload transport rate. The non-zero intercept of (1.21) will be discussed in Section 1.6.2. 36 The individual VuV bedload transport rate samples and corresponding concurrent 5 minute ADP bedload velocity averages displayed moderate correlation (Figure 1.9). It appears that relatively short sampling times of 5 minutes can produce useful data, although this correlation was heavily influenced by one high transport rate sample. One datum was deleted for which the sample was observed to fall out of the VuV sampler during retrieval. In this analysis we used only stations where sand was not measured in transport by the HS sampler, as the VuV sampler could not measure the sand fraction. For three of these stations we also compared the cv of the concurrent 5 minute ADP bedload velocity samples to the cv of the VuV bedload transport rate. The cv of va was consistently less than the cv of ga (Table 1.3), which suggests that the ADP bedload velocity is a more reliable measurement of bedload transport than conventional sampling. A significant correlation was also obtained for va versus mean bed shear stress (x) (Figure 1.10): ~a= 0.00128 X 1 ' 4 5 (r 2 =0.44, p< 0.0005 ) (1.22) The standard error of the exponent is 0.40, thus the shear stress exponent of 1.45 is the same as in the typical relationship between bedload transport rate and x 1 ' 5 (Yalin 1972). A power relation was achieved using a log transform of the data, which also helped reduce heteroscedasticity. One outlier with a low shear stress was neglected in the regression, which was justified as this station also had an unrealistic value of ks that cast doubt on the x estimate. 3 7 There is more scatter evident in Figure 1.10 than Figure 1.8, but this is expected as mean shear stress is merely a measure of mean flow competence, whereas bedload transport is also dependent on sediment supply, entrainability of the bed, and fluctuations in fluid force. Similarly, gb was not significantly correlated with x despite the fact that the mean velocity profiles displayed semi-log linearity, which suggests that x should have provided a reasonable estimate of mean local bed shear stress (cf. horizontal standard error bars in Figure 1.10). Tentatively, it appears that partial, patchy mean bedload transport rate can not be tightly predicted even with good measurements of mean local bed shear. 1.6 Discussion 1.6.1 Significance of a linear relation between bedload velocity and bedload transport rate Our results (1.21) indicated that apparent bedload velocity varied linearly with bedload transport rate, in accordance with the kinematic relation (1.12) if da (l - Xa ) was constant. This suggests that the depth and porosity of the active layer do not vary greatly, and the velocity of transport is the most important parameter in the kinematic equation. However, actual particle velocities do not necessarily increase linearly with transport rate, as the percentage of the bed surface that is mobile may also increase (1.14). 38 1.6.2 Influence of suspended sediment in the sample volume We are presently uncertain of the site specificity of (1.21). Possibly, variations in the bedload particle size distribution from site to site will alter the curve. It is noteworthy that (1.21) does not have a zero intercept. It appears that fine sand moving over the essentially stable bed created a velocity response at low mass transport rates. We suspect that bottom tracking was positively biased by suspended sediment near the bed ("water bias"), as well as by bedload proper. Fraser River is very turbid during freshet, with near-bed (15 cm above the bed) total fractional suspended sediment concentrations of 0(1000 mg/L) at the location of our measurements (Environment Canada 1996). If water bias does influence the intercept of the calibration curve, then it is likely that variability from site to site in near-bed suspended load will make the calibration site specific. We roughly estimated the degree of backscatter from the average sampling volume that could be expected due to an estimated near-bed fractional suspended sediment concentration (Table 1.4). The suspended sediment concentration was estimated from previous observations (Environment Canada 1996) between 1968 and 1983 at a vertical 274 m from the right bank (i.e. the centre of the channel, near the location of our measurements). The vertical distribution of suspended sediment concentration for each size fraction was estimated from all available data for flows between 5000 and 8000 m3/s. The fractional concentrations at 9 cm above the bed were estimated by extrapolating these curves (for each profile the measurement nearest the bed was 15 cm from the bed). The number of particles in the average sample volume during a ping was estimated for each size fraction using the concentration estimate at 9 cm above the bed. The average 39 sample volume depended on the depth and the pulse length (see Figure 1.2). The fractional mass was calculated by multiplying the fractional concentration by the average sample volume. The number of particles in a size fraction was determined by converting the fractional mass to a volume and dividing by the volume of a spherical grain with a diameter equal to the geometric mean of the size fraction. The total projected area from these grains was estimated, again assuming spherical particles. Literature values for form factors (|/|) (Table 1.4) were used to modify the expected backscatter from the total projected area of each particle size. It appears that suspended sediment could have effectively covered up to 3 percent of the sampling area (depending on the flow depth). With an estimated typical near-bed (9 cm depth) velocity of 1.2 m/s, near-bed suspended sediment may have contributed 3 cm/s to the va signal, which agrees well with the intercept of (1.21). However, the Helley-Smith sampler should have collected most of the suspended sediment within 7.6 cm of the bed. In future experiments, concurrent measurements of near-bed fractional suspended sediment concentration would be useful. Assuming that suspended sediment introduced the same degree of bias for all sampling stations, then it may be reasonable to remove this bias as estimated by the regression intercept: gb =1.4v„ (1.23) where va =v.- 0.033 (1.24) is the corrected apparent velocity of bedload. 40 We should note that the form factors used in the analysis above depend on the acoustic wavenumber (k), and thus the operating frequency of the aDcp. We used a 1.5 MHz aDcp. The form factors would have been smaller for a lower frequency instrument, and thus suspended scatterers would have less of an impact on the bottom track velocity for an equivalent pulse length. Thus the calibration expressed in (1.21) is likely dependent on the acoustic operating frequency. This leads to an important question: where in the active bedload layer is the bottom track ping reflecting? It is well established that acoustic signal attenuation is inversely related to operating frequency. It may be that lower frequency instruments will penetrate through the mobile layer entirely to reflect off the solid, immobile boundary, with the result that no bedload velocity will be recorded (see, for example, Section 4.4.1). This will also depend on the size of the mobile particles, as particles for which krt>2 should influence the bottom tracking. Further, with a high frequency instrument and a high concentration of suspended scatterers, the bottom track pulse may be attenuated before it reaches the bedload layer. Presumably, in this case, the bottom would not be observable and the bottom tracking processing would not calculate a velocity. This likely occurred in our deep water measurements on Fraser River, where bottom tracking was not possible. 1.6.3 Measurement error The scatter in individual estimates of va was probably due to both measurement error and real temporal variability of bedload transport. It is difficult to quantify the relative importance of measurement error and real variability as contributors to the observed 41 variance in va (Figure 1.6a). However, real transport should only be in the downstream direction. The variance was nearly isotropic: the standard deviation of va in the transverse direction was consistently about 0.9 times the standard deviation in the vector direction (r =0.83). This variance in the transverse direction may have been a result of incorrect resolution of the bedload velocity vector due to sparse bedload transport, as outlined in Section 1.3.3. Furthermore, a large number of physically unrealistic negative apparent bedload velocities were recorded. Thus, it appears that the variance was largely due to measurement error. It is interesting to note that the variance of vfl appeared to depend on the measurement depth, with very high standard deviations observed at shallow depths of 2 m, and minimum standard deviation at about 3 m depth (Figure 1.11). Again, it is difficult to separate real variability from measurement error (see Chapter 3 for a method to do so), but it appears that the ADP worked best in depths of about 3 m. Real temporal variability would have been associated with stochastic variability in entrainment of individual or groups of particles and/or quasi-periodic transport of bedload sheets (Gomez et al. 1989). The passage of large bedforms was not likely a factor: dune migration was not apparent in time series of individual beam depths, nor did spectral analysis of the bed speed time series reveal significant regular periodicity that would be expected with dune migration. Real temporal variability should be high for 5-second samples of bedload, even under steady flow and plane-bed conditions (Einstein 1937, McLean and Tassone 1987, Kuhnle and Southard 1988). A flume study by Kuhnle and Southard (1988) of statistically steady equilibrium bedload transport of a fine gravel mixture moving as bedload sheets provides a lower-bound estimate of the expected real 42 variability. The coefficient of variation decreased logarithmically with increasing sampling duration (their Figure 22). Extrapolating their result to a 5-second sampling duration produces a cv of 0.65 (for their run L2). Our measurements were collected during steady flow, but the steadiness of sediment supply was unknown. Furthermore, while the mean transport rate for run L2 (0.041 kg/m/s) was similar to our field results, the sediment was finer (D50 = 3 mm), thus transport in our measurements was probably more sporadic. It is likely that the real temporal variability of bedload transport would have produced a cv for our 5-second samples on the order of a value of one. The rest of • the variability in va can be attributed to measurement error. There are several methods to estimate the standard deviation of va due to measurement error (crf). First, a lower bound on ae can be crudely approximated by fitting a normal distribution, with a mean of zero, to the negative values of va. This assumes that the negative values are purely erroneous and do not include any positive signal. Second, if it is true that the coefficient of variation due to real variability equals one, then as should equal the measured standard deviation minus the mean. Third, a theoretical cr£ can be approximated by pooling (root sum of squares) the estimated errors due to instrument noise, dynamic tilt (assumed to be of 0(1 cm/s)), heading error, and DGPS error described in Section 2.3. For our station data, all three estimates of ae equal each other on average. It is noteworthy that virtually all of the pooled error is due to instrument noise. This issue is explored further in Chapter 3. 43 1.6.4 Apparent bedload velocity The ADP yields a single estimate of bedload velocity, which we have termed the apparent bedload velocity (vfl). At the present time we do not know if va is a true measure of the actual average bedload velocity (vp), but it appears that v„ may have been an overprediction of vp. Use of (1.18) with observed va and gb consistently yielded a small "characteristic" bedload particle size (< D5 of the bedload, see Table 1.5). These values were an order of magnitude less than expected for DCh calculated using (1.19) and (1.13) from the observed fractional bedload transport rates (Table 1.5). The discrepancy may have been due to high values of va . Similarly, an attempt was unsuccessful to use observed va with the kinematic model (1.12) to predict the bedload transport rates measured with the samplers. Assuming an active layer porosity of 0.4, unrealistic values of da (about 0.3 mm, see Table 1.5) were required to match predicted to measured gb . Assuming a porosity of 0.78, which is typical of individual saltations, yields da equal to about 0.9 mm. These values for da are consistently less than estimated using Equation 1.13 and the observed fractional bedload transport rates. As mentioned above, va may have been positively biased by near-bed suspended sediment transport. If the corrected va is used, estimates of DCH and da are both about 0.8 mm (Xa = 0.4) or da of 2.4 mm (ka = 0.78). The estimate of da equal to 2.4 mm is comparable to the estimated da using Equation 1.13, and may be reasonable for saltating sand particles. It may be that saltating sand dominated the bedload velocity signal, or that va was an overprediction of vp. Controlled experiments were conducted in a laboratory towing tank to investigate vfl 44 systematically (Chapter 2). However, sand transport was not considered in the laboratory. 1.7 Conclusions A new technique has been presented for remote measurement of bedload transport. An acoustic Doppler current profiler was used to measure apparent bedload velocity. Mean apparent bedload velocity correlated well (r^ O.93) with mean bedload transport rates measured using conventional samplers. Thus, a calibration curve can be developed to relate apparent bedload velocity to bedload transport rate. A long sampling duration, on the order of 25 minutes, was required to achieve a reliable estimate of the mean apparent bedload velocity. Remote measurement of bedload transport using acoustic Doppler technology holds great promise, as measurements can be taken with relative ease and safety at channel-forming discharges throughout a study reach. This makes the technique especially favorable for measurement of bed movement in large rivers. Improvements to the technology may allow for collection of useful data with greater temporal resolution. 45 Table 1.1: Towing tank tests of single ping ADP bottom tracking tow velocity (m/s) pulse length (cm) n X (m/s) % error (%) std x (m/s) cv y (m/s) stdy (m/s) -2.00 60 40 -1.985 -0\8 0.083 0.04 0.010 0.050 -1.50 60 30 -1.488 -0.8 0.039 0.03 0.008 0.051 -1.02 60 56 -1.018 -0.2 0.029 0.03 0.004 0.031 -1.00 60 71 -0.992 -0.8 0.034 0.03 0.005 0.025 -0.50 60 96 -0.499 -0.2 0.015 0.03 <|0.001| 0.012 -0.30 60 122 -0.300 0 0.017 0.06 -0.001 0.014 -0.10 a 60 125 -0.098 -2.0 0.019 0.19 <|0.001| 0.005 -0.05 a 60 199 -0.049 -2.0 0.015 0.31 < 0.001 0.004 0 60 751 0 0.001 0.05 a 60 232 0.047 -6.0 0.016 0.34 <|0.00T| 0.005 0.10 a 60 122 0.096 -4.0 0.024 0.25 <|0.001| 0.006 0.30 60 121 0.298 -0.6 0.019 0.03 -0.001 0.013 0.50 60 93 0.495 -1.0 0.017 0.03 <|0.001| 0.015 1.00 60 67 1.002 0.2 0.022 0.02 0.004 0.023 1.02 60 55 1.012 -0.7 0.059 0.06 <|0.001| 0.031 1.50 60 35 1.500 0 0.049 0.03 0.021 0.055 2.00 60 49 1.981 -1.0 0.070 0.04 -0.006 0.070 -0.50 20 77 -0.501 0.2 0.018 0.04 0.002 0.016 0.50 20 76 0.496 -0.8 0.020 0.04 0.002 0.018 -0.50 b 60 108 -0.300 -0.6 0.017 0.03 0.497 0.016 0.50 b 60 110 0.003 -0.8 0.018 0.04 -0.496 0.018 a Tow carriage had trouble maintaining steady velocity in these runs (tow velocity appeared to be periodic about the expected velocity, with maximum error of about 0.01 m/s). bADP was rotated 90° such that motion was in y direction, % error and cv based on y component 06261 06301 06302 06303 06304 06305 06306 06307 07131 07132 07172 07173 07174 07181 07182 07183 07184 07261 07262 07263 Sample (1) W N) S O) CO O) O) 0 ) 0 - ' l\5 A O M ( D W J i > I M M O ( O M - < 0 ) 0 ) 0 ) S U M ( D 0 ) - k J ( 0 - i U 0 ) ( f l 0 1 O 0 l > M O M - i M ^ < » N > C » 2.14 1.99 7.07 3.60 3.62 2.91 5.81 5.80 3.30. 3.41 2.81 2.56 2.18 2.08 2.32 3.21 2.60 3.02 4.50 5.11 Depth (m) (3) 0.114 0.134 0.033 0.080 0.066 0.133 0.129 0.136 0.104 0.120 0.120 0.121 0.172 0.129 0.129 0.112 0.117 0.107 0.106 0.099 0.005 0.026 0.003 0.006 0.008 0.009 0.005 0.008 0.005 0.005 0.006 0.015 0.030 0.002 0.006 0.004 0.007 0.004 0.004 0.004 0.034 0.075 2.5e-7 0.384 0.004 0.195 0.021 0.038 0.008 0.027 0.038 0.043 0.345 0.089 0.083 0.022 0.037 0.037 0.012 0.005 0.052 0.154 0.015 0.008 0.011 0.018 0.216 0.276 0.070 0.052 0.040 0.069 0.110 0.046 0.055 0.041 0.042 0.032 0.078 0.079 Prin (m mean (7) 0.239 0.429 i 0.135 0.125 0.123 0.081 0.224 0.254 0.174 0.135 0.136 0.160 0.291 0.305 0.207 0.120 0.141 0.152 0.188 0.158 Appare nary Is) st.dev (8) 0.009 -0.038 0.003 0.021 -0.020 0.006 -0.009 -0.006 -0.005 0.008 0.004 -0.010 0.040 -0.017 -0.004 0.000 0.004 0.004 0.007 -0.008 ;nt Bedsf Secoi (m mean (9) 0.224 0.444 0.138 0.097 0.115 0.098 0.152 0.181 0.154 0.126 0.123 0.112 0.222 0.282 0.187 0.120 0.136 0.146 0.165 0.168 )eed ndary Is) st.dev. (10) ro o s cs co o i A ^ K , cn o -»• _ ,„ ,„ _ , , ro U - ' O N O M U ^ U l O O l M M M O l - ' - ' U - 1 * — 0.0344 0.0385 0.0863 0.0206 0.0125 0.0104 0.0067 0.0027 0.0735 Transpt (kg/: mean (12) — 0.0334 0.0103 0.0853 0.0147 0.0033 0.0101 0.0035 0.0001 0.0895 Drt Rate s/m) st.dev. (13) 0.028 0.016 0.015 0.020 0.00045 0.009 0.00047 0.00034 0.050 3edload D50 (m) (14) ^ | u o i c n m s o i o i l col I 1 1 1 1 1 1 1 n VuV (15) w | W U U W C O M O J o | | | | | | | | | n HS (16) 3 CT\ Table 1.3 Comparison of coefficient of variation of concurrent samples of ADP velocity (va) and VuV bedload transport rate (gb). Station cv of gb cv of V a n 07132 0.98 0.24 6 07183 1.62 0.46 5 07263 0.54 0.073 2 Table 1.4 Estimate of near-bed (9 cm above bed) fractional suspended sediment concentrations and form factors (k = 6300 m"1) for calculation of proportion of va due near-bed suspended scatterers. D, mm C, mg/L kri | / | 0.0020 60 0.0063 0.00014 a 0.0028 28 0.0089 0.00028 a 0.0057 51 0.018 0.0011 a 0.011 62 0.036 0.0045 a 0.022 66 0.070 0.017 a 0.044 77 0.14 0.068 a 0.088 127 0.28 0 .10 b 0.18 231 0.56 0 .18 b 0.35 212 1.1 0.50 b 0.71 58 2.2 1.0 b a From Medwin and Clay (1998), Rayleigh scattering, | / | = 24n{krif b From Thorne et al. (1995), Figure 9. 4 9 Table 1.5 Estimated characteristic particle diameter, active layer depth, and active layer porosity for bedload samples. All units in mm. Sample D5 D50 da* Aa* Dch b DChc da dae dal d g 07132 7.0 28 17 410 10 0.37 0.4 1.1 1.1 3.0 07173 0.31 16 3.1 710 3.3 0.32 0.4 0.7 1.0 1.8 07174 0.44 15 5.0 640 4.9 0.43 0.5 0.7 1.4 2.0 07181 0.29 20 3.3 710 3.6 0.15 0.2 0.8 0.4 2.3 07182 0.18 0.45 1.4 770 1.6 0.09 0.1 0.5 0.3 1.5 07183 0.25 9 1.8 750 2.0 0.27 0.3 1.9 0.9 5.2 07184 0.18 0.47 1.4 770 1.6 0.13 0.1 0.3 0.4 0.9 07261 0.16 3.4 1.1 780 1.3 0.05 0.05 - 0.1 -07263 8.5 50 24 380 13 0.47 0.5 0.8 1.4 2.3 a Calculated using Equation (1.13) with observed fractional gbi. This required estimation of the number of particles in transport in each size class and the projected planar area of a particle for each size class. Spherical particles with specific gravity of 2.6 were assumed. It was also assumed that particles smaller than 2 mm saltated with dt = 3* Dj and Xai =7/9, and particles larger than 2 mm rolled with dai = A and A,ai = 1 /3. b Calculated using Equation (1.19) and estimated da and Xa. c Calculated using Equation (1.18) with observed va and gb . d Calculated using Equation (1.12) with observed va and gb and assuming A„ = 0.40. e Calculated using Equation (1.12) with observed va and gb and assuming Xa = 0.40. f Calculated using Equation (1.12) with observed va and gb and assuming Aa - 0.78. 8 Calculated using Equation (1.12) with observed va and gb and assuming Aa = 0.78. O o ° p * CT CT p figs ) « )§Hs o o o o O O O cr CD Q 3 CO o Q o I 2.5 -2.5 -i 1 1 1 1 1 ; 1 1 1 1 1 -2.5 -2 -1.5 -1 -0.5 0 0.5 1 1.5 2 2.5 actual velocity, m/s Figure 1.3: Towing tank test of ADP bottom tracking for 60 cm pulse length. 53 Figure 1.4: Study site locations. The study site for chapter 1 was the gravel study site on Fraser River. Also shown are the study sites on Norrish Creek (Chapter 3) and the sand bed Sea Reach of Fraser River (Chapter 4). 54 Sieve Size, mm Figure 1.5: Sampling station composite bedload cumulative grain size distributions. 55 East Bedspeed, m/s -1 -0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 _50 -il 1 I i I i 0 50 100 150 200 250 East, m Figure 1.6: Sample 07132 (91 minute sample), a) Apparent bed velocity (+, individual five second ensemble average; 0, mean). Mean and [st. dev.] of East and North bed speeds are -0.053 [0.136] m/s and -0.007 [0.126] m/s, respectively, b) Boat trajectory by bottom tracking. Trajectory by DGPS nearly fixed at location (0,0). 56 Time, s 0 1000 2000 3000 4000 5000 i 1 1 1 1 1 — -0.05 J 1 1 1 1 1 — L 5 0 200 400 600 800 1000 K. Number of Five Second Samples Figure 1.7: Sample 07132. a) Time series of primary apparent bed speed, b) Running average (lower dashed line) and running coefficient of variation (upper solid line) of primary bed speed time series. Monte Carlo simulation (see text) suggests that 25 minutes of sampling is required for a reliable estimate of the mean. 57 Figure 1.8: Mean bedload transport rate versus mean primary bed speed ( • , 2 min ADP sample; O , > 30 min ADP sample). Error bars represent precision (± standard error), with standard error of bedload based on variability of VuV samples. 58 0.25 Primary bedspeed from ADP, m/s Figure 1.9: Bedload transport rate versus mean primary bed speed for individual concurrent samples without sand transport. 59 10 20 T , N/m2 30 Figure 1.10: Mean primary bed speed versus shear stress (O, included in regression; • , outlier excluded from regression (see text)). Error bars represent precision (± standard error), with standard error of shear stress based on standard error of u* from log-law fit. w 0.40 i E 0.0 2.0 4.0 Depth, m 6.0 8.( Figure 1.11: Influence of measurement depth on variance of va. 61 Chapter 2: Laboratory measurements of bedload transport velocity using an acoustic Doppler current profiler 2.1 Abstract A laboratory test was conducted to assess the accuracy and precision of bedload transport velocity measurement using acoustic Doppler current profiler (aDcp) bottom tracking. A synthetic bedload o f known velocity was created by dragging strings of stones over an artificial river bed in a towing tank. Trials were conducted at various particle speeds, percent mobile bed area, and aDcp orientation with respect to the direction of transport. Estimated apparent bed velocity (v a) vectors showed a large degree of scatter within a trial, and the bed velocity resolved in the direction of transport was biased toward zero. However, mean values of the magnitude and direction of v a predicted reasonably well the average bed velocity. The variance within a trial was due to three factors: 1) instrument noise, 2 ) variability in the percentage mobile bed area within a beam due to the large size of mobile particles with respect to the beam area, and 3 ) the inability of a three beam aDcp to resolve heterogeneous transport. We argue that instrument noise was the primary source of variance. 62 2.2 Introduction We are exploring the potential of acoustic Doppler current profiler (aDcp) bottom tracking for measurement of bedload velocity, with the goal of developing a non-invasive technique for gauging bedload transport. Refer to Section 1.3 for a review of the method. We have tested the method in the field, and found the measured mean apparent bed velocity to correlate with the mean bedload transport rate determined using conventional samplers (Chapter 1). However, in the field-testing we could not ascertain whether or not the measured va was an accurate measure of the actual average bed velocity. In this chapter we present laboratory tests in which we created a synthetic bedload of known velocity and compared the measured apparent bed velocity to the expected average bed velocity. 2.3 Methods We created a synthetic bedload by dragging strings of stones over an artificial river bed beneath an aDcp. We used a commercially available 1.5 MHz three beam Acoustic Doppler Profiler (ADP™) made by SonTek to measure the apparent bed velocity. The tests were performed in a towing tank of painted concrete that is 60 m long, 12ft (3.66 m) wide, by 8 ft (2.44 m) deep. The width and depth appeared to be sufficiently large for the ADP bottom tracking to operate effectively, at least for static bed measurements (see Table 1.1 and Figure 1.3), although the bottom tracking pulse length had to be adjusted to 60 cm in order for the ADP to measure the depth accurately. Adjacent to the tank is a 63 track that runs two trolleys, which are pulled by a 25 hp motor via a cable system. As both trolleys are clipped to the same cable, they both move at the same velocity. Beams were cantilevered out over the tank from each trolley, and strings of stones were strung between the two trolleys. The trolleys were pulled back and forth, so that the stones were pulled back and forth through the sample area beneath the aDcp at an essentially steady, uniform velocity (see Figure C.4 in Appendix C). The velocity was varied between trials, ranging from 0 to 50 cm/s. The aDcp was set at a depth of either 1.9 m or 2.0 m above the artificial river bed, which consisted of small river cobbles (b-axis diameters of 2 cm to 5 cm, rounded to sub-rounded) epoxied to 8 ft x 8 ft of V" plywood sheeting. The artificial bed was used to ensure that the backscattering received by the aDcp was similar to the backscattering from the towed particles, and to improve simulation of a field condition. The area of the bed insonified (or "seen") by the aDcp bottom tracking depends on the distance between the aDcp and the bottom, the beam geometry, and the spread of each beam (Section 1.3.1). We used an aDcp with beams oriented at 120° relative azimuth angles and projecting down at 25° from vertical. We have assumed that the half-intensity beam width represents the beam spread. For a sampling depth of 2 m, the sample area consists of three quasi-elliptical areas, each with a major, azimuthal axis of 0.13 m, an area of 0.012 m2, and centred at 0.94 m along 120° relative azimuth angles from the aDcp. The orientation of the aDcp beams with respect to the direction of particle motion was varied between trials. 64 The stone lines consisted of the same small river cobbles epoxied to 70 lb test nylon fishing line, with buttons strung along the strings providing the contact surface between the stone and the string. The rocks were spaced 10 cm apart, and each string had 50 stones, for a total of 5 m. The length of artificial river bed was 2.4 m, thus 2.6 m of sampling length was available with the artificial bed fully covered by a stone line. Strings were passed through the aDcp sampling area at 15 cm spacing for some trials and 10 cm spacing for other trials. String spacing was maintained by attaching the ends of each string with swivels to eyescrews spaced along wood boards that were in turn roped to the trolley beams. The strings were also passed through eyescrews spaced along wood boards weighted down on the tank bottom, to maintain line spacing and to ensure that the stones were dragged over the artificial bed rather than through the water column. Time series were collected for several particle speeds at each line spacing. We used an aDcp that sends and receives a bottom tracking signal once per second, and each traverse of the sample area produced only a few seconds of data (depending on the particle speed); thus time series were collected as aggregates of several traverses. Only data in which the stringed stones covered the length of the artificial bed were used. Time series were aggregated separately for traverses in the forward and backward directions. Three trials were also conducted without the artificial river bed. The expected average bed velocity (ve) was the product of the particle velocity (vp) and the percentage of the insonified bed surface that was mobile (fm): ve=vpfm (2.1) 65 The particle speed was well controlled, but the percent mobile bed area varied during a trial. In order to determine/„, a photograph of a sample of 72 stones was digitized: the mean and standard deviation of the projected surface area of the stones were 2.2x10"3 m2 and 7.4x10"4 m2, respectively. Similarly, a pebble count of 100 stones (a and b axes) yielded a mean and standard deviation of projected surface area (assuming an elliptical shape) of 1.9xl0"3 m2 and 5.9X10"4 m2, respectively. Using the estimate from digitization, an average stone occupied about 19% of the insonified area from a single beam. The small size of the insonified area may have been problematic: the percentage of the insonified area occupied by moving particles would have varied during a trial as particles moved into and out of the view of the beams. A Monte Carlo simulation was performed to estimate the variability of coverage within a beam, based on possible positions of the stones and the measured distribution of projected area of the stones. The stone areas were modelled as circles. The instantaneous mobile coverage for a single beam could have ranged from 0 to 48% with a mean and standard deviation of 14.5% and 8.5% for the 15 cm line spacing, and from 6% to 48% with a mean and standard deviation of 21.7%o and 6.1% for the 10 cm line spacing. In this chapter it has been assumed that fm equalled 14.5% for the 15 cm spacing and 21.7% for the 10 cm spacing when calculating the expected mean velocities. The effect of variable/„ on the results will be discussed further below. 66 2.4 Results For each trial the va data were plotted in polar coordinates (for example from Trial 12 see Figure 2.1). Wide scatter was evident, with measurements from individual pings being completely unreliable. However, the measured apparent bed velocities tended to be in the direction of bed movement. The data sets for all the trials were summarized (Table 2.1). Comparisons between observed va and expected ve were plotted for both the mean bed speed resolved in the direction of particle movement (vad ) (Figure 2.2a) and the mean magnitude of the velocity vector (vav) (Figure 2.2b). The observed mean direction of transport (Qobs) was also compared to the actual direction of transport (Qexp) (Figure 2.3), with angles defined as counter-clockwise rotations (facing down) about the aDcp x axis (see Figure 1.1). The following regression equations were derived for Figures 2.2a, 2.2b, and 2.3: v~d = 0.21v -0.0060; ad e ( 2 2 ) SE = ±0.048 cm/s, r2 =0.67 = 1.06 ve+0.15; (2.3) SE = ±0.1 lcm/ s, r2=0.9\ Qobs = 1.02 9^ + 4.6; SE = ±3.5°, rl = 0.97 (2.4) SE is the standard error of the estimate. The linear regressions of Equations 2.2, 2.3, and 2.4 were weighted (Montgomery and Peck 1982, pp. 362-363) by the reciprocal of the standard error of each mean in order to account for the variable precision in mean 67 estimates, which was due to variable durations of sampling and increasing variance with increasing va (Figure 2.4). It is encouraging that vad shows an increasing trend with ve (Figure 2.2a); however, vad tends to underpredict ve by 79%. Alternatively, vav displays much closer agreement to the expected ve (Figure 2.2b); and the mean direction of the va vector also tends to be reasonably well predicted (Figure 2.3). In other words, the measured magnitude and direction of va both tended to be distributed around their expected means. It appears that the aDcp predicted the mean magnitude of the bed velocity, and the mean direction of the bed velocity, but not the mean bed velocity component in the direction of transport. Thus, it may be necessary to estimate the mean ve using vav and the mean direction of va. There is a great deal of scatter in the measurements. About 2 minutes of data were required for stable estimates of the mean and cv of the magnitude of the bed velocity vector (Figure 2.5) and of the bed velocity vector direction (Figure 2.6). The coefficient of variation of the va vector magnitude was about 1.28 on average (Figure 2.4). Probability distribution functions of the va vector magnitudes were not significantly different from lognormal (Kolmogorov-Smirnov test) with means near the expected ve (Figure 2.7). Lognormal distributions indicate a combination of random terms by a multiplicative process, which seems reasonable for combination of error terms that produce instrument noise. The distribution of va(/ for the same data set is shown in Figure 2.8, and will be discussed in Chapter 5 in the context of the error modelling presented in 68 Chapter 3. The thrust of the discussion will be that the lab data differ from the field data in their error structure, and the noise observed in the lab may have biased vact toward zero. The aDcp appears to measure smaller values of va when the transport is in the negative x direction than with positive x transport. The x direction is defined as the horizontal direction of the first beam. An aDcp is typically deployed with the flow, and thus the transport, in the negative x direction. The negative x direction corresponds to a transport direction of 180° (solid symbols, Figure 2.2). Orienting the aDcp with transport in the negative x direction yielded a better estimate of ve for vav, but the converse was true for It appears that trials without the artificial bed yielded estimates of va in closer agreement to ve (triangles in Figure 2.2). Owing to the 25° beam elevation angle from vertical, the backscattering from a smooth, flat bottom is probably less than from a rough bottom, as most of the acoustic pulse reflects away from the aDcp. Thus, it is reasonable to expect superior registering of mobile particles in the absence of the artificial bed, as backscattering from the immobile bottom is reduced. 2.5 Discussion 2.5.1 Measurement error The laboratory tests suggest that an aDcp can, on average, correctly resolve the magnitude and direction of bed velocity. However, it appears that vaa was biased toward 69 zero. The underprediction in the direction of travel was due to inability of the ADP to correctly resolve the direction of transport for individual pings. There are three conceivable sources of error that resulted in the scatter. First, due to imprecision in the system hardware or signal processing, the instrument may have incorrectly assessed the Doppler frequency shift in one or more of the beams. We will describe this as instrument noise. Second, variable percentage mobile area within a beam due to the large particle size and spacing would have produced temporal variability in va, despite the fact that ve was steady if measured over the larger scale of the artificial bed. Third, the same variability in percentage mobile area would have caused heterogeneous transport between the three aDcp beams, which, as will be explained below, would have produced erroneous va vectors. We will argue that instrument noise was the predominant source of scatter in the data. In particular, we will show that the inability to correctly resolve the direction of transport for individual pings was the result of large erroneous velocity measurements in individual beams. Error due to instrument noise for bottom tracking over solid, stable substrates is an order of magnitude less than error for water column velocity estimates (V. Polonichko, 2000, SonTek Inc., personal communication) because a strong backscatter signal is received from the solid boundary. The standard deviation of the horizontal water velocity measurements is reported to be (Theriault 1986a, SonTek 1998) 140 c ow = = (2.5) F Az V A where c is the speed of sound (nominal value 1500 m/s), F is the acoustic operating frequency (1.5 MHz), Az is the depth cell size (0.25 m), and N is the number of pings. 70 The number of pings equals the averaging interval (Is) times the pinging rate (9 pings per second for the 1.5 MHz aDcp). Substitution of these values yields crw equal to 19 cm/s, which suggests that the standard deviation of bottom tracking velocities due to instrument noise for a population of one-second ensembles should be on the order of 2 cm/s. Bottom tracking tests in the towing tank with a mobile instrument showed that the standard deviation ranged from 2 to 8 cm/s, with higher standard deviations observed for higher towing velocities (Table 1.1). However, bottom tracking error is greater when the bed is mobile, and the internal signal processing is not presently designed to resolve the velocity of a mobile bed. As a model of gravel transport in rivers, in our laboratory tests only a portion of the bed within each beam was mobile, while the majority of the area was stationary. The monofrequency sound pulse was likely reflected with two different frequencies: the original frequency for the portion insonifying the stationary bed, and a Doppler shifted frequency from the portion reflecting off mobile particles. The bottom tracking algorithm is proprietary, but it involves a pulse-to-pulse incoherent ("narrowband") technique (V. Polonichko, 2000, SonTek Inc., personal communication), thus probably a single, strong amplitude frequency is found and used to determine the Doppler shift (see Brumley et al. (1991)). It seems likely that the instrument would have more difficulty interpreting the received frequency spectrum of a mobile bed. Instrument noise can be attributed to two broad factors: 1) Doppler noise, which is uncorrelated from sample to sample; and 2) signal dwell time in the sample volume, 7 1 which is related to the pulse length (Lemmin and Lhermitte 1999, see Section 3.6.2). We suspect that the Doppler noise was high in the towing tank tests, although the ADP does not output signal to noise ratio for bottom tracking, which makes assessment of Doppler noise difficult. Signal to noise ratio is a measure of the power of the signal in the acoustic return versus the power of the background noise in the return. Doppler noise is proportional to the standard deviation of the backscattered Doppler spectrum, thus spectral broadening increases Doppler noise. Spectral broadening can occur when different scatterers in the sample volume have different velocities, as discussed above. This cause of random error tends to increase with higher particle velocities as a broader signal spectrum is required to account for velocities ranging from zero to the particle velocity. This may explain why we observed higher variance with higher mean velocities (Figure 2.4). Spectral broadening also occurs with receiver or backscatter noise, which may have occurred in the lab trials. It is useful to examine individual beam velocities when evaluating the noise. Theriault (1986b) presented equations to calculate the forward (vx) and transverse (v ) velocities given the velocity determined along each beam, assuming homogenous velocity in all three beams: 72 In this notation, beams are numbered counterclockwise when the aDcp is facing down, with the horizontal component of beam 1 parallel to the x-axis. The angle <|> is the beam deflection angle from vertical (denoted % in Chapter 1). Beam velocities were extracted from the data in xyz coordinates by solving (2.6): v M = vx sin(§) - v z cosfa) (2.7a) - vx sin(§) + V3 vy sinfa) - 2 v z cos(<\>) vb2 = — (2.7b) - vx sin(<\>) - S v y sinfa) - 2 vz cos(§) vb3 = (2.7c) In Figure 2.9-2.11 we present a short time series of the measured and expected vol, Vb2, and Vb3 beam velocities for Trials 11 and 12 (i.e. both forward and reverse directions). In general, along each beam the ADP was able to register zero velocity when the bed was immobile and register a velocity when the particles were moving. The direction of travel was aligned with the ADP x-axis (parallel to horizontal component of Vbi), and vb] shows the most coherence between observed and expected bedload velocity. We should acknowledge, however, that vw was underpredicted (Table 2.2). On the other hand, the Vb3, and the in particular, displayed large excursions from the expected velocity. Interestingly, vb2 was closer to the expected value than vbX or vb3 (Table 2.2). The expected beam velocities for Vbi and were only 1.4 cm/s, which approaches the lower limit of velocity resolution for the ADP (measurements are recorded to the nearest mm/s). This might explain the inability of the ADP to resolve the velocities in beams 2 and 3. We should note that the expected beam velocities in this case were similar to what might 73 be expected in a field condition of partial transport (Wilcock and McArdell 1997) of coarse gravel without fine material. Due to the three-beam geometry of the ADP, some measurement error inherently arises when bedload varies between the three beams. Consideration of the 3-beam geometry of the ADP, with beams separated by 120° relative azimuth angles and slanted at fj> (25°) from the vertical, yields the measured beam velocities in terms of the x, y, and z bottom track velocity components in each beam: vb\ = v\x sin(§) - v u cosfo) (2.8a) vb2 = - v2x sin(30°)sin(§) + v2y cos(30°)sin(§)-v2z cos(§) (2.8b) vbi = - v3x sin(30°)sin(§)- v3y cos(30°)sin($)-v2z cos(§) (2.8c) Combination of (2.6) and (2.8), and assuming that all three viz are equal, yields predictors based on the forward and transverse velocities through each beam: v , = [2 vxl +0.5 (vx2 +VJ+0.866 ( v , 3 - v y 2 ) \ / 3 (2.9a) v, = [0.5 (v , 3 -v x 2 )+0M6 (vy2 +vy3)]/S (2.9b) When the percentage mobile bed area differs between the three beams, the direction and magnitude of the measured velocity vector will necessarily be wrong, due to the unequal weighting of vxi in the estimated vx. As discussed above, this may have arisen in our tests because individual particles were large with respect to the beam width, and occupied a substantial portion of a beam. The percentage of the bed mobile within a beam varied as particles moved through the beam. This error mechanism could have created variance in 74 the forward and transverse directions. In fact, if only one of beams 2 or 3 registered a velocity, then the direction of observed va would have been 60° offset from the expected direction, which corresponds to several of the velocity vectors in Figure 2.1. The errors due to a large beam 2 or 3 velocity are as follows: However, based on Equation 2.9a if all vyi are zero, the magnitude of the observed va vector can not exceed the magnitude of the largest actual within-beam v&. Furthermore, direction errors can not exceed 60°. In Figure 2.1, some va vectors have magnitudes that exceed possible values of v0, even if the areal coverage within a beam was 48%. Also, some direction errors exceed 60°: the direction of the va vector had a standard deviation of about 97° (Table 2.1). The observed data can only be explained if vy and/or negative vx bedload velocity components were recorded. Examination of beam velocities showed that beam 2, in particular, had large positive and negative excursions in some trials (Figure 2.10). In our laboratory tests, there were no actual negative vx components, and only minimal actual vy components were caused by random jostling of the particles as they were dragged over the artificial bed. It appears then, that most of the variance was caused by instrument noise rather than real variability infm. Large value Direction of va vector Positive beam 2 Negative beam 2 Positive beam 3 Negative beam 3 -60° 120° 60° -120° Another possible explanation for the scatter in beams 2 and 3 is interference from side lobes or multiple scattering effects. Beams 2 and 3 were pointing partially toward the tank walls, and thus beam side lobes may have reflected off the tank walls. While side lobe interference was not observed for field tests or for mobile instrument towing tank tests (see Chapter 1), it is conceivable that side lobes reflecting off the tank walls or tank bottom contaminated the mobile bed towing tank tests. Side lobe interference, if present, would have been related to processing the signal from a partially mobile bed within the confines of the towing tank, and would have biased the beam velocities toward zero. In fact, beams 1 and 3 displayed more bias than beam 2 (Table 2.2). It would be desirable to evaluate the instrument noise independent of error due to unsteady percent mobile area within a beam. To this end, it would be useful to conduct further trials with a steady percentage of mobile bed area within beams. Alternatively, the expected variance in v„ due to unsteady percent mobile bed area and the corresponding heterogeneous transport can be determined using our estimates of the probability distribution of percent mobile bed area. Synthetic time series of values of vxj, vX2, and vX3 were created based on the estimated percent mobile area probability distributions. The distribution of actual ve was estimated by averaging instantaneous values of vxi from the three beams. The coefficient of variation (cv) of actual ve was 0.33 for 15 cm string spacing and 0.16 for 10 cm string spacing. Similarly, the distribution of expected measured va was derived using Equation 2.9, assuming vyi were zero. As expected, the maximum deviation from the expected direction was ±60°. The cv of expected measured vav was 0.41 for 15 cm string spacing and 0.19 for 10 cm string spacing, as compared to the measured cv of about 1.28 for all the trials (Figure 2.4). The 76 standard deviation due to instrument noise (Ginstrument) can be estimated assuming that the variance is pooled: ® measured ~ ® instrument ^®expected (2-10) Equation 2.9 yields a cv of va due to instrument noise of about 1.26 (Figure 2.4), i.e. virtually all the measured variance was due to instrument noise. We considered the possibility that frequency tracking biased the measurements. Frequency tracking utilizes feedback from previous pings to narrow the passband of the instrument's noise filters (Chereskin et al. 1989, Chereskin and Harding 1993). This technique increases the signal to noise ratio by narrowing the search range for the Doppler signal, but can introduce bias if the filter passband is too wide or if it is not centered on the signal. However, SonTek ADPs do not use frequency tracking (C. Ward, 2002, SonTek Inc., personal communication). 2.5.2 Future Tests Further lab trials to test the precision and accuracy of the method may be unwarranted with the present instrument. Development of an improved bedload Doppler sonar prior to rigorous lab testing may be more fruitful. A full gamut of trials could be performed at various particle sizes, particle speeds, percent mobile bed area, and aDcp orientation with respect to the direction of transport. It would be possible to conduct multivariate analyses to determine if the aDcp is more capable of determining bedload transport velocity for certain particle sizes, transport rates, and aDcp orientations. It would also be possible to determine if the aDcp is more sensitive to increasing particle velocity or 77 increasing percent mobile bed area. The procedure would have to be modified in order to make the materials robust for such a large number of tests. The testing procedure was eventually abandoned due to breakage of strings. There were two major problems with the existing procedure: 1) the stones were large with respect to the beam area, which limited our ability to control the percent mobile bed area, and 2) the strings were easily tangled if lines were run close together within a beam area. It would also have been useful to use stones of perfectly uniform size, in order to control better the percent mobile bed area. A new method of conveying the stones should be used in any future trials. The problem is to ensure that not all the area insonified by a beam is mobile. Possible methods include attaching small stones to a net and dragging the entire net, and epoxying granules to thick ropes and dragging the ropes through the sample area. Alternatively, a calibration could be conducted in a mobile bed canal or flume of sufficient depth using underwater video to monitor particle velocities through the sample area. 2.6 Conclusions A laboratory test was conducted to evaluate the applicability of acoustic Doppler current profiler (aDcp) bottom tracking for accurate and precise determination of bedload transport velocity. A synthetic bedload of known velocity was created by dragging strings of stones over an artificial river bed in a towing tank. Trials were conducted at various particle speeds, percent mobile bed area within the aDcp beams, and aDcp orientation with respect to the direction of transport. Estimated apparent bed velocity (vfl) vectors showed a large degree of scatter within a trial, and va resolved in the direction of travel was biased toward zero, but mean values of the magnitude and direction of va predicted reasonably well the average bed velocity. The variance within a trial was due to three factors: 1) instrument noise, 2) variability in the percentage mobile bed area within a beam due to the large size and spacing of mobile particles with respect to the beam area, and 3) the inability of a three beam aDcp to resolve correctly heterogeneous transport (i.e., transport that varies between the three beams). Based on the estimated variability of percent mobile bed area within a beam during a trial, virtually all of the measured variance was due to instrument noise. Ul O l O l Ol o o o o 4^ O J w ro ro - i o o o o o o o o o ro ro —*. —x ro ro V i ro ro J>. J > J > '-j b i cn b i cn Ol b i b i b i b i b i 01 ro o o o o < < z •z.<<-z.-z.<<<<<<-z.< teg 2 " • C D CD CD fflfflOifflOiOiKSoiOlSHu-'M^n^ ( o c o - ' r o M o i S S f f l r o g g u r o i o o g S I: cT 00 c 3 3 •3 o & o >-t o CO 51 1-1 & o T3 CD P S CO n a 3 n> »-+• CT* CT & < CT, o o o - N i - ~ j o i o i O ) o > r o r o r o r o o o ( J J Q J I O M O O C O O I O I I D I O M M ^ ^ S W o o oo o o -»• £ oo J>. J>. o - ^ J > . O 0 0 C » 0 0 C » C U CD o' 3 m X $2. O O S" DO CD CD CO ( D J i O l O l S C O M O W t O l M O l O O O ^ ^ b i c » a 5 o ' ^ c o ' ^ b i r o ^ 4 > c o ^ - - > -Cf 3 ; M 3 . CO 03 CD a. ro oo o ^ CO CO CO CO w oo ro Ol oo ro -»• — - v l o —J Ol co — ^ oo - J . ro ~ ro ro Ol £ 8 w O) CD CO oo T3 CD CD — ,o r~ ^ < 13 Q . ® t/T CD O O < O < CD O O C D C O O O C O C O C O C D n n O l ( J ) 0 0 O ) O ) O l - ' 5 ° CD CD o CD CD 5 - J CD g 00 S 3 CD CD |CQ CO W a . CD • < ro bo J> ro o ro ->• o -»• ro CD ro J>. co ho ro p o o o o r o o o o CO co - s i o o CD 00 '->• q s 3 CD < . CO " 3 CO o r* 3 a. E < CD O O O => o o q - i "O co o •S § o ~ * 5-Table 2.2: Beam velocity statistics for Trials 11 and 12. All units cm/s. Trial Beam Mean St. dev. Expected 11 V ; 0.1 1.7 2.8 11 v2 -1.4 7.3 -1.4 11 V3 0.1 2.6 -1.4 12 Vy -0.1 2.4 -2.8 12 0.8 4.5 1.4 12 VJ 0.1 4.2 1.4 Figure 2.1. Single ping measurements of apparent bed velocity (cm/s) for a trial with expected bed speed of 6.5 cm/s and expected direction of 180° (Trial 12). Speed radii for 23, 47, and 70 cm/s are shown. 1% CL L > 8. o 3 O cr CD CL TO' cr a 6. O g CD cr CD CL p •a CD cn CD a a o >T3 O a cn g. N> CD CD " cn < % 2_• CD 3 CD era' a* CD CL I-I CD CTQ >-t CD Cn cn S' 3 a cn CD X r 1 o w m y-En* 3' 9*. o' oo o o cr CD CL • o o O 3. CD a *—t-P «-t-o' a o a* cr CD CL oo o cn O a 3 § o CD CL O pa cn' * 0 §1 cr g cn c r CD a CL S*. < O CD £ o CD O H~ o • • a P CD W p CD a B CL a P <-C L £L CD g a 2. CD cn CD O o a o ° 1 P cn a-CD CL a* cr CD CL O oo o a 4 cr CD CL • o a o a* CD p T3 CD P a CL ".' CD O r-f-o' a o 6 CL CD H 1 B P cn CD ~ a _ CD CD cn M e a s u r e d m e a n a p p a r e n t b e d ve loc i ty , t r an spo r t d i rec t i on ( v a d ) , cm/s i - o m X "O CD Q_ CD C L cu CD s cr CD C L < cp_ o o < CO o M e a s u r e d m e a n m a g n i t u d e o f a p p a r e n t b e d ve l oc i t y v e c t o r ( v a v ) , cm/s o co ro cn T X 3 CD O n, 5 8 era' O ~ fil P £5 °> J " g t o t r oo CD ^ U) CL o • e a X s 3 s CL p cr CO CL 00 "3 CD 1-1 c v o rt-p era i - ! CD CD CD o oo CD CL o a cr P a 3 + a* o p_ cr CD CL < CD !-| 00 . a 00 00 a - cr U J CD » . 00 " CD ocra CD l — i 00 oo O a a oo 3 CD a a o 3 O CD oo oo rt-o" 3 a a o CL 3 5' CD P Cf S 3 2 « Fr a ^ a v; o II ss' CD 00 o i p.' O • ~» a co "VJ cr II S 31 era' S 3 g * ° o • a BH o 3 00 o- 5' O 0Q oT cn ° £ CL 3 % CD CL CD P C . ft 2 p 3 af- CL ^ . n ™^ • ~ - <=L ~ . 2. o (-K cr. 22 CL 3 5 2 5. CD P P rt- 00 2 3 3 era Standard deviation of bed velocity vector magnitude, cm/s 4*. 00 O 00 84 12 TO 10 0 c o ro ro 6 > .2 4 o M— M— CD O 2 O mean coefficient of variation - i i _ 100 150 200 250 300 350 400 N u m b e r pings Figure 2.5: Running mean and coefficient of variation of the magnitude of the bed velocity vector for Trial 12. ro 5 0 )] K 2 + + 3 ' 3 [- 0.5 cos(- cp) + 0.866 sin (cp)] [p. 5 cos (- cp) + 0.866 sin (cp)] V2)cos(y) V3 \sin(q>) + 2Nbl -Nb2 -Nb2 „ „ „ ^ _ , _ i v 62 _ J V 6 3 3sz'n((|)) cos\q> I + —?i sin *J3 sin{§) ( 9 ) Simplifying yields the final model: V = -cos2(q>)Vl + .. ^ 4- ras((p) - Vo. 5 sin (cp) 1 3 + 0.5 cos((p) + VoTJ sm((p) F 3 +A K^(o(3.733a^) (3.12) (3.12) Note that the distribution of the noise in the direction of transport (Ny) is the same as the noise distributions in the x or y directions, and that Ny is independent of the direction of transport (cp). Thus, the model noise is isotropic. We assumed that the Vt were independent and identically distributed, which is reasonable for sporadic gravel transport in three separate but proximate areas of river bed. Note that for the coefficients (c,) on the Vt terms in (3.12), Ydcl: = 1 and ^c, = 0.5. Now, since /=i 1=1 2 c,- = 1, the mean of the model (V) will match the data mean (v) if Vt = v. The noise 1=1 term in (3.12) accounts for instrument noise from each of the three beams, but does not include measurement errors due to averaging heterogeneous transport between beams 100 (see (3.8) and Chapter 2). These averaging errors are included in the model as a result of combinations of random realizations of the F,-, which simulates heterogeneous transport. The measured data can be described as a random variable D with probability density function fa- The probability density functions of the V signal and the noise are, respectively, fy andfy. Assuming/^ and O N , the model distribution of V (i.e.,/V) was generated using (3.12) by brute-force computational combination of-20,000 random realizations of F, and N. The random realizations were generated proportionally to increments of the integrated density function. Some care was required in generating V from fy, as the required/^ were highly positively skewed, and log-scale increments were required for Fnear zero where the probability asymptotically approached infinity. The mean and standard deviation of the generated Fwere within 2% of the expected values. The unknown parameters CTN and fy, which included the standard deviation of the signal (rjv), were optimized by heuristic Monte Carlo simulation and with a genetic algorithm (Houck et al. 1995), wherein the optimum model distribution (fv) was chosen by minimizing the mean absolute difference between/v and fo. The two optimization procedures produced similar results. The final fit statistics were calculated from a run of ~100,000 random realizations of Vi and N. This procedure was conceptually simple but computationally expensive. Further, the optimization was limited by the strong influence of/// on fv; in essence, many possible fy fit reasonably well as long as the mean of the signal equaled the mean of the data. Efficient algorithms for maximum likelihood deconvolution remain an area of active 101 research. Most algorithms require a priori knowledge of the error distribution (aN) (e.g. Cordy and Thomas 1997, van Es and Kok 1998), are nonparametric and thus lend little insight to the actual physical process involved, and are also prone to difficulties with large a N (Wand 1998). 2 „2 2 It should be noted that we could have set = — (Harris 1966), where c, are N 3.733 the coefficients in (3.12), but the resulting/v would have fit poorly as the higher moments would not have been considered. By ignoring this constraint, the variance of our model does not necessarily match the variance of the data. In essence, the optimum fit is a tradeoff between the second and higher moments, and the shape of the distribution is fit. 3.4.3 Model distributions of spatially averaged bedload velocity Two possible signal distributions were fit to the data. A semi-theoretical additive compound Poisson-gamma distribution was developed based on the assumption that particle entrainment followed a Poisson process and particle velocities followed a gamma distribution. Also, a gamma distribution was empirically fit. If we correctly surmised that we measured partial bedload transport, it is reasonable to assume that the number of particles entrained was random and followed a Poisson process. A Poisson process occurs if events are randomly distributed in a space-time domain that can be divided into a large number of equal intervals such that the probability of more than one event occurring in a space-time interval is negligible, and 102 the probability of occurrence of an event in an interval is independent of events in all other intervals. The probability density function for a Poisson random variable is f(n) = e-^-, n = 0,l,... (3.13) nl with mean A, and standard deviation , where n is the number of particles in motion in the interval. It is worth noting that the Poisson distribution is additive, i.e. the sum of two Poisson variables with parameters A.i and X2 is distributed as Poisson(A,i + X2). Thus, a sampling regime is not required to resolve the Poisson space-time interval in order for the Poisson distribution to describe the measured data of a random phenomenon. The Poisson distribution describes the number of particles in motion. In order for the distribution to describe the spatially averaged velocity, the distribution of the particle velocities must be considered. A compound Poisson process occurs when the events in a Poisson process display values which are themselves random variables (Vpj). The density function of a compound Poisson process is (Feller 1957, p.270) f*(VJ = \ X e-X2Z^-{ff forVp>0 n=\ nl e~x forVp=0 (3.14) where {/}"* is the n-fold convolution of the random variables (i.e. the (additive) density function of Tn where T n = Vpi+VP2+.. .+Vpn). We develop our theoretical model based on the assumption that the particle velocities are gamma distributed. Although laboratory measurements show normal distributions for sand saltation velocity (Hu and Hui 1996b, 103 Lee et al. 2000), field measurements of gravel velocity (Drake et al. 1988), which include both rolling and saltating particles, show left-skewed distributions that can be modeled as gamma distributions. The gamma Ga((^ ,a2) density function is: f(x) = £^xa-Ae-*x ,x>0 (3.15) r ( a ) with mean a/p and standard deviation Va/p (Larsen and Marx 1986). In the case where vpi are independent and identically distributed (i.i.d) gamma variables, the distribution is additive compound Poisson-gamma, with density: f*{vP)-P e n=\ r(cc n) n! Vp>0 (3.16) Vp=0 with mean = XNp °-/r. and var = XNp (a + j- (Jorgensen 1997), where Np is the / P P number of particles within view of each aDcp beam, X Np is the mean number of particles in motion, vp = ^ is the mean particle velocity, and avp = is the standard deviation of the particle velocities. The distribution of spatial average velocity can be obtained from the resulting additive distribution if the velocities are divided by, and the probabilities multiplied by, Np. As the mean of the additive distribution is XNpvp, the spatial average velocity signal mean will match the mean of the data if X = = (3.17) VP 104 Finally, it can be shown that the standard deviation of the spatial average velocity distribution is avbs = (3.18) N v p p To limit the search space for the optimization, we constrained the particle velocity gamma distribution by ensuring the maximum particle velocity (mean + one standard deviation) did not exceed the maximum near-bed water velocity. Only one standard deviation was used to ensure that the constraint did not excessively restrict the search of the solution space. The near-bed water velocity was measured by the ADP as a 50 cm deep spatial average (25 cm bin), centered >25 cm above the bed (41 cm for Fraser River, 27 cm for Norrish Creek). Single ping data were collected at Norrish Creek, but the averaging time was 5 seconds (5 pings) for the Fraser River data set. The additive compound Poisson-gamma (cPg) distribution developed above is similar to the additive cPg distribution for total particle displacement suggested by Yang and Sayre (1971), with a total sampling area Np substituted for the total sampling time t. The difference is that (3.16) yields a spatial average of velocity for several particles, whereas the model of Yang and Sayre (1971), when divided by t, yields a temporal average "virtual" velocity for individual particles. As mentioned previously, the two should be equivalent for bedload transport that is statistically steady and uniform. The unknown parameters are ON, N P , vp , and o\p- We adopted two approaches for parameterization of the model. First, N P , vp , and rjvpwere estimated and ON was 105 optimized. Second, all four unknown parameters were optimized. For the latter approach we recognize that an infinite number of combinations of Np, vp , and rjyp yield distributions with the same mean and variance, thus a unique solution may not be possible. However, the optimization procedure fit the higher moments of the distribution. The number of particles insonified in each beam (Np) was estimated. Each measurement of bedload velocity is a spatial average over the acoustic footprint of the three-beam aDcp bottom echo (Section 1.3.1). The vertical range from the ADP to the bed for the Fraser River sample was 3.05 m, which resulted in an acoustic footprint of 0.019 m2 for each beam (Equations 1.2 and 1.3). For a bedload D50 of 0.028 m, N p equals 31 particles. Similarly, the vertical from the ADP to the bed for the Norrish Creek sample was 2.42 m, for a beam footprint of 0.014 m . The NC bedload D50 was crudely estimated by grid-by-number photo-sieving of an ortho-rectified photograph of a post-flood downstream unarmoured gravel deposit. The estimated bedload D50 was 0.028 m, which corresponds to an N p of 23 particles. We assumed vp = 3u *, where u* is the shear velocity, which is a typical value for saltating gravel (Nino et al. 1994, Drake et al. 1988). The measured non-dimensional shear stress for the bedload D50 was 0.097 for Norrish Creek, thus bedload particles were probably in the transition region between rolling and saltating (Drake et al. 1988, Hu and Hui 1996a). The measured non-dimensional shear stress was 0.032 for Fraser River, which is near the critical value for entrainment, thus it is probable that the gravel particles were rolling instead of saltating. Drake et al. (1988) observed larger, rolling particles to have average velocities 30% less than smaller, saltating particles. Regardless, the estimated vp of 0.36 m/s (FR) and 0.62 m/s (NC) are reasonable when compared to field measurements of mean instantaneous velocities of individual gravels and small cobbles, which range from 0.15 to 0.60 m/s (Drake et al. 1988, Schmidt and Ergenzinger 1992, Chacho et al. 1994, Habersack 2001). Furthermore, we assumed CTvp equals 0.6 vp (i.e. cvvp is 0.6, where cv denotes coefficient of variation), as observed in field measurements by Drake et al. (1988). Note that Habersack (2001) observed a maximum particle velocity of 2.0 m/s when average particle velocity was 0.60 m/s, which implies a larger deviation than cvvp equal to 0.6. The gamma distribution (3.15) was also empirically fit. Two methods were employed using (3.12). The first method (denoted by beam gamma) used independent but identically distributed gamma random variables for Vi. This assumes that the spatially averaged bedload velocity in each beam was gamma distributed. We also fit a single gamma distribution which accounted for the first three terms in the simplified version of (3.12) (denoted by total gamma). This assumes that the bedload velocity signal following transformation through the three beams was gamma distributed. Unfortunately, this does not provide information about the spatially averaged bedload velocity, as it is not possible to estimate the distribution for bedload velocity within a beam. 107 3.5 Results We were able to generate model distributions that fit the data distributions very well (Figure 3.1 and Appendix D). A Kolmogorov-Smirnov test (a = 0.05) did not find any of the model distributions to be significantly different from their respective data distribution (Table 3.2). The compound Poisson-gamma distribution was better fit when input parameters were optimized rather than estimated. The optimized cPg and beam gamma distributions were very similar (Figure 3.2). Interestingly, a marginally better fit was achieved using the total gamma distribution, where the linearly combined bedload velocity signal was modelled as gamma distributed. The Norrish Creek model did not yield quite as good a fit as the Fraser River model (Figure 3.3); the NC data set had a peak at 20 cm/s, but this peak could not be adequately modelled by either the cPg or gamma without losing the fit to the tail. For both the FR and NC data sets, optimized bedload speed distributions (fv) were highly positively skewed indicating that the bedload speed averaged within a beam area (v) tended to be mostly near zero with a few high values. Positively skewed distributions were expected for partial bedload transport of gravel, where the particles remain at rest most of the time (Einstein 1937). The optimum a v equalled 0.14 m/s (i.e. a coefficient of variation of 2.6) for the FR data, and 0.54 m/s (cvv of 2.1) for NC. The NC data consisted of single pings, whereas the FR data was collected as 5 ping averages. Standard deviations for FR need to be standardized by multiplying byV5 for comparison with the NC data. Standardized FR a v is 0.31 m/s (cvv of 5.9). The highest estimated v (calculated 108 as v+2av ) equalled u5 for NC, but was substantially less than u5 for the FR data, even when calculated using standardized ov The greater value of standardized rjv in the higher v NC data suggests variability of transport increases with greater transport rate although the coefficient of variation is less. Kuhnle and Southard (1988) observed the same trend. The estimated linearly combined three-beam total noise distributions had cjN.totai of 0.094 m/s for the FR data and 0.31 m/s for the NC data. Multiplying CTN.totai for FR byV? , for a value of 0.21 m/s, standardizes the error estimates. It appears that greater error was observed at NC than FR. This is despite the fact that the two data collections used the same pulse length in similar depths, and the FR data were collected from a boat, which added error sources related to instrument motion and use of the compass and DGPS (Chapter 1). This supports our previous assessment (Chapter 1) that instrument noise dominates the error. It also suggests that instrument noise may be partially dependent on the bedload velocity, as somewhat greater noise was observed in the much higher v NC data. The values of 0.21 m/s and 0.31 m/s are comparable to the expected standard deviation of water column velocity due to instrument noise (SonTek 1998, see Equation 1.5), which is 0.23 m/s for a one second average (9 pings) and a 20 cm pulse. 3.6 Discussion 3.6.1 Bedload velocity distribution The deconvolution procedure yielded highly positively skewed gamma or cPg distributions for spatially averaged bedload velocity. We are unaware of published data 109 that can directly validate our estimated bedload velocity distributions. However, previously measured distributions of bedload transport rate, particle velocities, and particle displacements corroborate our findings. In particular, the good fit achieved using a gamma distribution for bedload velocity compares favourably to measurements of particle displacements (Stelczer 1981, Hassan et al. 1991). Similarly, the good fit with a compound Poisson-gamma distribution supports the cPg as a semi-theoretical distribution for spatially averaged bedload velocity, as has been found for bedload virtual velocity (Einstein 1937, Yang and Sayre 1971, McLean and Tassone 1987). Finally, bedload transport rate observed by previous investigators has followed a gamma distribution. For example, Kuhnle and Southard (1988) measured fine gravel bedload transport rate across the full-width of a flume every 30 s for 150 minutes. We plotted (Figure 3.4) the distribution of transport rate for their low transport rate run (L2) from the time series data (their Figure 7). The distribution closely follows a gamma distribution, with parameters determined by the method of moments to be a = 4.0 and P = 2.3. The coefficient of variation of run L2 was estimated to be 0.6 for 5 s samples (see Section 1.6.3). Similarly, Gibbs (1973), as part of a laboratory calibration of a basket sampler, measured bedload transport using a central slot trap for three different sampling durations for each of three different flow rates. The sediment size distribution was a 1:5 scale model of the gravel-bed Elbow River, thus the model distribution ranged from 1 mm to 10 mm with a d5o of 4.7 mm. In all of his runs he observed that the square root of bedload transport rate followed a normal distribution. We found that the bedload transport rate distributions can 110 also be fit by gamma distributions. Gibbs (1973) observed that the coefficient of variation of bedload transport rate increased as sample duration decreased (Table 3.3). The coefficients of variation observed by Gibbs (1973) and Kuhnle and Southard (1988) compare to our estimated cv of within-beam spatially averaged bedload velocity for five second samples of 2.6 for FR and 0.9 for NC. Flow condition 2 of Gibbs data approximately corresponds to the bedload transport rate for our Fraser River data (0.0344 kg/s/m), but the bed material in his flume was finer. The cv for 60 s samples was 0.79, and the cv for 5 s samples may have been about 1.4 based on a logarithmic extrapolation of the data in Table 3.3. This extrapolated coefficient of variation is comparable to our estimated cv, but still less than we observed at FR. The model Froude scaling suggests that the Fraser River data would be modelled by a laboratory transport rate of 0.00307 kg/s/m. Clearly, very low transport rates were observed on Fraser River, which suggests the partial transport was very sporadic, and very high cv can be expected. Also, we were observing bedload within the relatively small space-time increment of 5 seconds in a beam area of between 0.014 m2 (NC) and 0.019 m2 (FR). We hypothesize that in such a small space-time increment, most often the bed would be stable, but occasionally the entire bed would be mobilized. Averaging over a greater spatial or temporal increment should reduce the skew of the measured distribution. It is apparent from Figure 3.2 that optimized cPg and beam gamma distributions were very similar. The cPg has some theoretical justification in that it is derived assuming random particle entrainment. The good fit with the beam gamma distribution, on the I l l other hand, is difficult to justify theoretically. The gamma distribution describes the distribution of space-time periods between every a t h event in a Poisson process. As the gamma distribution describes the distribution of the periods between events, rather than the distribution of the events themselves, there is little theoretical grounding for the use of this function. To understand why the gamma distribution provides a good empirical fit to the data, it may be helpful to consider Tweedie models. Jorgensen (1997) describes the cPg distribution and the gamma distribution as special cases of a Tweedie model. Several conventional distributions can be described as Tweedie models depending on the Tweedie parameter (p): the normal (p=0), Poisson (p=l), cPg (l £ H 3. b <° c/> vo 4 C O i ro CD 3 CO Z J CO 0) o ' ZJ Probability density Probability density II Probability density 3 |Z1 «P CD •'••a 2 n 0) cr c/> o_ 0 C L 31 —1 CD 3 n o b cn ro -vi o 3 co" C L CD < II O CO 3 5r -o CO a. o CD < CD O S •D X, CT II II P oo ro o> P W ^ CO CU o ^ <° 8 -i cn r*-co _. ^ IO 3 J » CT ^ "D CD • CU CU o a-. 3 -* < S 3 o o 3 3 CD CD CU CU 3 3 =1 II 2. o CO • «s II ro o 0 0 3 w o Q_ CD "0 0) - i Qi 3 CD I—I-CD —\ IP 3 O C/) CD to o 121 40 c/> 30 c CD X 5 ~ 20 ~\ r T 1 1 r • * ' FR cpg opt. •*• FR beam gamma NC cpg opt. NC beam gamma _] T 1 1 1 1 1 r 0.04 0.08 0.12 0.16 Within beam spatial average bedspeed, m/s 0.2 1E+052 *Z 1E+042 C _g 1E+032 >^ i = 1E+022 1Q CD 1E+012 O °- 100 a 1E-008 FR cpg opt FR beam gamma — NC cpg opt. — NC beam gamma o c? Uj o 8 UJ LU Uj Within beam spatial average bedspeed, m/s Figure 3.2: Deconvolved distribution of within beam spatial average bedspeed for the Fraser River and Norrish Creek data sets, using both an optimized compound Poisson-gamma and a gamma distribution. Both linear and log scales shown. 122 Figure 3.3: Norrish Creek measured data and fitted model using optimized compound Poisson-gamma distribution. 123 Figure 3.4: Comparison of measured distribution to gamma distribution for run L2 bedload transport rate of Kuhnle and Southard (1988) 124 0 0.1 0.2 0.3 0.4 0.5 Particle velocity, m/s 0.001 0.01 0.1 1 Particle velocity, m/s Figure 3.5: Estimated gamma distributions of entrained particle velocities from cPg model for Fraser River and Norrish Creek data sets. Both linear and log scales shown. 125 Period, s 200 100 50 20 10 2 -pi—i 1 M i i i i—i i L_ 1.0 0.7 H 0.2 co c 0 T3 "cn o CO 0.7 H 0.2 Primary bedload velocity i I 11111 ' ' 1 1 I 1 i 1 i I 11111— Primary water velocity f-r-rrrf i I i n i l 1 1 1 1 I 1 i 1 I Secondary water velocity I i m l 1 1 1 ' I r i I i i 11 r i i i 111 1 1—i i i 1111 r 0.01 0.1 Frequency, Hz 1.0 0.7 Figure 3.6: Spectral densities of primary bedload velocity and near bed water velocity components. Mean spectral density indicated by dashed line, and 95% confidence interval indicated by 0.2 gray l i n e s-0.2 126 Chapter 4: Spatial distribution of bedload transport velocity 4.1 Abstract Maps are presented of the spatial distribution of two dimensional bedload transport velocity vectors. Bedload transport determines the morphodynamics of fluvial and coastal environments. However, bedload remains poorly understood, in part due to inadequate sampling techniques that can not measure the spatial and temporal variability of bedload. We have been developing a new technique to measure bedload transport velocity in the field, using the bottom tracking feature of acoustic Doppler current profilers (aDcps) (Chapter 1). If a differential global positioning system (DGPS) is used concurrently, then estimates of the apparent velocity of bedload can be gained at a nominal rate of 1 Hz from a moving boat. The apparent velocity should be a measure of the mean velocity of the bed surface. However, the data are noisy, and smoothing is required. Data were collected using an aDcp within a selected area near the Agassiz-Rosedale bridge in the gravel-bed reach of the Fraser River, as well as in the sand-bed Sea Reach of the Fraser River delta. Two dimensional vector maps were generated of the spatial distribution of measured bedload transport velocity. Bedload velocity vectors interpolated onto a uniform grid revealed coherent patterns in the bedload velocity distribution. Concurrent Helley-Smith bedload sampling in the sand-bed reach corroborated the trends observed in the bedload velocity map. Contemporaneous 2D vector maps of near-bed velocity (velocity in bins centered between 25 cm and 50 cm from the bottom) and depth 127 averaged water velocity were also generated from the aDcp data. Using a vector correlation coefficient, which is independent of the choice of coordinate system, the bedload velocity distribution was significantly correlated to the near-bed and depth averaged water velocity distributions. The bedload velocity distribution also compared favorably with variations in depth and estimates of the spatial distribution of shear stress. 4.2 Introduction At all spatial scales, river morphology is determined by the spatial pattern of sediment transfer along a river reach through local erosion and deposition. Understanding of the spatial distribution of bed material transport is thus required for many aspects of river management, including prediction of channel change, analysis of stability of engineered structures such as bridges, and assessment of aquatic habitat. However, relatively little is known about the spatial distribution of bed material transport. Measurement of bedload with conventional samplers is time consuming, difficult, and inaccurate, thus measurement programs are rarely undertaken. Typically, river scientists rely on either scaled mobile boundary physical models or multidimensional morphodynamic numerical models to predict the spatial distribution of bedload and associated channel dynamics. Besides difficulties with appropriate scaling, physical models are tedious and expensive to construct. Given sufficiently detailed boundary topography, which is now readily available through digital surveying and photogrammetry, numerical models can reasonably predict 2D and 3D flow patterns in rivers (e.g., Yusuf 2001). However, morphodynamic models must also predict the spatial distribution of bedload transport, 128 and thus rely on analytic expressions for bedload transport rate. Bedload transport formulae are notoriously poor predictors of bedload (Gomez and Church 1989, Habersack and Laronne 2002), generally because the formulae are derived from empirical observations of reach average conditions, whereas bedload transport is a complex spatially and temporally heterogeneous phenomenon that depends on local interactions between the turbulently varying flow field, bed structure, and sediment supply. Calibration of morphodynamic models usually relies on their ability to predict water levels. The spatial distribution of bedload transport rate, which drives the morphodynamics yet is unreliably predicted, is not calibrated, due to the lack of available data for calibration. Thus, it would be extremely useful to develop a procedure to measure conveniently the spatial distribution of bedload. We have been developing a new technique to measure bedload transport velocity non-invasively using the bottom-tracking feature of commercially available acoustic Doppler current profilers (aDcps) (see Chapter 1). In this chapter we apply the technique to moving boat applications, such that bedload velocity measurements are collected throughout a reach of river. The goal is to generate maps of the spatial distribution of bedload velocity. We begin with a review of previous assessments of bedload spatial variability. We then describe the field measurements, which were collected in both sand and gravel bedded reaches of Fraser River. The new technique yields noisy data, thus data smoothing is required, and we present our methods for doing so. Finally, maps of bedload velocity are produced and compared to measured distributions of bedload transport rate, shear stress, and near-bed water velocity. 129 4.3 Background on spatial distribution of bedload Particle entrainment is thought to be driven by the bursting process of the macroturbulent flow. Bedload tends to occur during sweeps with excess instantaneous streamwise and downward (i.e. downstream slanting) velocities (Drake et al. 1988, Nelson et al. 1995, Williams 1996). Assuming that mean local boundary shear stress correlates with the intensity of the macroturbulent flow, one may expect the spatial distributions of bedload transport rate and size distribution to follow the spatial distribution of bed shear stress in the channel (Dietrich and Whiting 1989, Ashworth et al. 1992b). However, the spatial distributions of bedload rate and caliber are also influenced by: upstream availability of sediment, which determines the volume and distribution of sediment delivered to the reach; sorting of bed sediments in vertical, lateral, and downstream directions, which modifies local particle entrainment; local bed topography, which determines gravitational direction, as well as the shear stress distribution; and near-bed secondary currents, which can cause cross-stream stresses and lateral bedload transport (Dietrich and Smith 1983). In this section we discuss available techniques to measure bedload spatial distribution in the field, and then review previous studies in sand-bed and gravel-bed rivers. 4.3.1 Measurement Two basic methods have previously been adopted for field assessment of the spatial distribution of bedload: direct sampling and morphological approaches. Direct sampling involves physically collecting bedload at points throughout an area of interest. The 130 extreme spatial variability of bedload requires collection of many, closely spaced samples (Rennie and Millar 2000). Replication at each point is also required due to sampler inaccuracy and temporal variability of transport. Further, the time required to collect sufficient samples often exceeds the time of steady flow and bed conditions. Thus, point sampling is usually impractical for assessment of the spatial distribution of bedload transport rate throughout a reach, and only a few studies have been conducted (Davoren and Mosley 1986, Dietrich and Whiting 1989, Bridge and Gabel 1992, Ferguson et al. 1992). More typically, sampling occurs at single cross-sections to assess total bedload transport rate. The measured transport rate is valid only at the section if there is streamwise variation in sediment availability. However, cross-channel variability in transport can be assessed. The morphological method uses a serial record of river bed topography to quantify time-integrated transport rates and patterns of erosion and deposition by assessing sequential changes in channel form (e.g. Simons et al. 1965, Neill 1987, Kostaschuk et al. 1989, Ferguson and Ashworth 1992, Ferguson et al. 1992, Goff and Ashmore 1994, Lane et al. 1995, McLean and Church 1999, Ashworth et al. 2000). There are three basic techniques: 1) morphometric sediment transfer between an identified zone of erosion and associated deposition, 2) a within-reach sediment budget, which yields streamwise variation in transport from two sequential river bed topographies and knowledge of transport rate at a section, and 3) dune tracking in sand-bed rivers, which involves measurement of dune migration rates. The simultaneous use of particle tracers can yield 131 information on transport pathways between sites of erosion and deposition (e.g. Laronne and Duncan 1992, Habersack 2001). Morphological methods can yield good estimates of spatial patterns of bed material transfer in dynamic rivers. However, in order to apply morphological methods, the topography time series must be of sufficient temporal and spatial resolution to capture relevant changes in channel morphology. If the channel changes slowly, required historical data may be unavailable. If the channel changes rapidly, data collection requirements may be intensive. Further, assessment of the spatial distribution of bed material transport requires that river morphology changes. Morphological approaches can not measure bed material transport that occurs as throughput or in suspension. Finally, time and space integrated assessments of bed material transport may be difficult to relate to the fluctuating forcing flow, which limits their utility for developing predictive models. 4.3.2 Sand-bed rivers Dietrich and Whiting (1989) reported spatially intensive Helley-Smith sampler surveys of bedload transport in sand-bed and gravel-bed meander bends. In the sand-bed meander, maximum bedload transport of coarser fractions followed the path of maximum boundary shear stress, which began at the upstream inner bend but was forced out by the point bar to the downstream outer bend pool. On the other hand, fine fractions exhibited maximum transport from the upstream outer bend to the downstream inner bend, with the result that finer sands were transported up onto the bar tail. 132 Bridge and Gabel (1992) and Gabel (1993) performed spatially distributed bedload transport measurements by dune tracking and Helley-Smith point sampling at three locations in a small sand-bed river divided by a mid-stream island. They found that local values of flow velocity, bed shear stress, and sediment transport were highly variable, both within and between dunes. The measurements did not yield coherent maps of bedload transport rate for individual floods, even when spatially averaged. However, averaging measurements from similar floods yielded estimates of the spatial distribution of bedload transport for different flows (Bridge and Gabel 1992). Patterns of transport were similar to those observed by Dietrich and Whiting (1989) in the simple meander bend, with the locus of transport following the zones of highest shear stress, which crossed from the inner bank to outer bank as each channel progressed around the island. This pattern allows for bar/island growth by downstream and lateral accretion, as deposition occurs in the zone of low shear stress and sediment transport. Ten Brinke et al. (1999) used daily echo soundings to estimate bed material transport by dune tracking throughout a 1 km reach of the large sand-bed River Rhine-Waal. The soundings were spatially intensive. Streamwise tracks of both 10 m spaced single beam sonar and approximately 19 m spaced multibeam sonar were taken. The multibeam tracks overlapped allowing for topographical mapping of the entire bed with about 6-7 points per m . Large dunes spanned the channel but were more pronounced on the slightly inner bending right bank. Small dunes were superimposed on the large dunes. The persistence and relatively slow migration of the large dunes allowed for effective 133 dune tracking, which was not the case at a site further upstream. However, bedload and suspended load measurements from further upstream suggested that bedload by dune tracking may have comprised <40% of the total bed material load. They found a regular pattern of dune migration parallel to the river banks, with erosion on stoss faces and deposition on lee faces. Bedload transport was fairly consistent across the channel, but was somewhat greater on the right side associated with the large dunes. They did not plot a spatial distribution of bed material transport rate. Kostaschuk and his colleagues have evaluated the hydrodynamics of individual large dunes in the Main Channel of the Fraser River Estuary (Kostaschuk and Illersich 1995; Kostaschuk and Villard 1996). Transport occurs during spring freshet at ebb tide. Bed material transport was measured using a pump sampler for suspended load and a Helley-Smith sampler for bedload. Sediment concentration and sediment transport rates increased up the stoss side of dunes and decreased down the lee side. Thus, maximum sediment transport was observed near the crest of dunes. However, bedload measured by the Helley-Smith sampler was not found on average to be different between stoss and lee sides of dunes. This was attributed to sampling error. Most sediment was transported in suspension. Bedload on average comprised only 0.6% of the measured bed material load, and 5% of the bed material transport estimated by dune migration. However, the Main Channel during high transport is an exceptional environment with Reynolds numbers of 0(107) and mean velocities ranging from 1.1 to 1.7 m/s. They observed that some of their dunes were symmetric, which they attributed to high suspended transport and associated lee-side deposition. 134 Carling et al. (2000) measured the hydrodynamics of individual large sand dunes over a stable gravel-bed in the River Rhine. Sediment supply in this case is relatively limited, although sediment supply within an individual dune is not. Bedload over the stoss side of dunes was spatially variable, but tended to increase with the local shear stress towards the crest. During dune development when dune height was increasing, the highest transport rate was observed on the lower stoss. They found that suspension played a relatively minor role in bed material transport, and measured bedload rates were comparable to bed material transport estimated by dune tracking. 4.3.3 Gravel-bed rivers Sear (1996) summarized previous studies of sediment transport through pool-riffle sequences in small streams. The conventional model is for transport of fine material from high shear riffles into low shear pools during low flows. As stage rises the locus of maximum shear stress shifts to the pools, and fine and coarse material are scoured from the pool. However, this pattern is heavily influenced by antecedent flow history and associated sediment availability (Sidle 1988). In the gravel-bed meander bend studied by Dietrich and Whiting (1989), the pattern of boundary shear stress was similar to the sand-bed meander, with maximum shear stress observed at the upstream inner bend to the downstream outer bend pool. The coarsest bedload followed this pattern. However, most bedload transport was focused in a narrow band in the centre of the channel, where excess dimensionless shear stress was maximum. 135 Presumably, grain sorting altered the pattern of bedload, and the coarsest surface sediments, which were in the zone of maximum shear, were not mobilized. Despite bankfull conditions the largest size fractions were immobile, and the relatively fine bedload was likely mostly throughput from upstream. Similarly, using a basket sampler at a cross-section on the River Inn, Mulhoffer (1933, reported in Church 1987 and Gomez 1991) observed that bedload transport tended to be concentrated in only a portion of the streambed. Further, he observed a shift in the locus of bedload transport away from the thalweg toward the bar top as discharge increased (corroborated by Hassan and Church 2001 using pit traps on Harris Creek). It appears that the bedload followed the zone of excess dimensionless shear, as observed by Dietrich and Whiting (1989). Low sediment supply, poorly sorted gravel-bed channels adjust vertical and lateral grain sorting such that bed mobility is minimized (Kellerhals 1967, Parker and Klingeman 1982, Brayshaw 1984, Dietrich and Whiting 1989, Dietrich et al. 1989, Church et al. 1998). In general, the bed surface is coarsest in the zones of high flow maximum shear. However, the shear stress distribution changes with discharge, and the relative immobility of gravel beds results in a lag in channel bed structure adjustment to the changing flow. The result is that zones of excess dimensionless shear, which drive the transport of bed material, shift location as discharge changes. It should be noted that it is also possible that the shift of locus of maximum transport is a result of a change in the trajectory of throughput material (Church 1987), which could be related to near-bed secondary currents driving sediment up the slope of the bar. 136 As already mentioned, in gravel-bed rivers there is the added complication that measured bedload may simply be throughput transported through the reach from upstream. Throughput does not interact with the bed locally, and thus does not constitute local bed material transport. The quantity of throughput is heavily dependent on upstream sediment supply. Many studies have emphasized the influence of vagaries of upstream sediment supply on bedload transport in gravel channels (e.g. Mosley 1981, Ashworth and Ferguson 1986, Davoren and Mosley 1986, Sidle 1988, Adenlof and Wohl 1994, de Jong and Ergenzinger 1998, Hassan and Church 2001). This confounds attempts to develop predictive relationships for bedload transport rate from flow variables. Under low excess shear conditions most sediment in motion is derived from upstream sources and relatively fine bedload over a stable bed is transported along the path of maximum shear (Phase I transport, Jackson and Beschta (1982)). As shear stress increases, more of the bed sediment itself is mobilized, and it appears that the locus of transport dynamically shifts to the location where maximum excess dimensionless shear exists. As shear stress rises further, there is an increase in the proportion in transport of the coarsest bed material (Wilcock and McArdell 1997), and transport may be observed throughout the reach (de Jong and Ergenzinger 1998). In divided wandering and braided gravel-bed rivers, the locations of maximum shear are associated with zones of flow convergence such as confluences and thalwegs. Minimum shear occurs in areas of flow divergence, where mid-stream bars and islands tend to form as a result of deposition. Davoren and Mosley (1986) used a basket sampler deployed from a jet boat to measure bedload transport across four cross-sections of an upstream chute and downstream mid-channel bar complex in the large braided gravel-bed Ohau River. Bedload transport was restricted to zones where flows were deep and fast, with largest rates in the thalweg of the chute. They inferred that sediment was transported through and scoured out of the chute to be deposited on the bar head. Thus, sediment transport in a braided gravel-bed river is related to bar construction and modification. They also noted that similar flow conditions in a supply limited reach resulted in no bed material transport due to armouring of the river bed. Ferguson and Ashworth (1992) used both morphological and point-sampling approaches to study several braided and near-braided gravel-bed rivers. They observed that just-competent flows tend to move sediment from zones of flow convergence in pools and confluences to zones of flow divergence on bars and in channel bifurcations. Cutoffs and avulsions create widespread channel change during high flows (e.g. Laronne and Duncan 1992). The most intensive sampling occurred in a Y shaped convergent chute and downstream divergent central bar complex on the Sunwapta River (Ferguson et al. 1992, Ashworth et al. 1992a). Shear stress from velocity profiles and 10 minute bedload transport samples were measured at between 26 and 39 locations in the 35 m long by 20 m wide sampling area on 4 different days. Over the study period the upstream chute initially aggraded (with bed surface fining), but then the chute degraded and a new unit bar lobe formed in the downstream right channel and accreted to the bar head. The sediment scoured out of the chute deposited to form the unit bar. This was reflected in the bedload transport measurements with highest rates observed in the locus of highest shear in the chute during the period of degradation. We note that the fining of the chute 138 during the initial aggradation phase prevented thalweg surface coarsening from limiting transport in the zone of maximum shear. The median bedload size (D50) was similar throughout the reach, which suggests the majority of bedload was not sorted during transport. However, the coarsest bedload followed the trend for maximum shear stress and bedload transport rate in the thalweg of the chute. Their measured rates compared favourably with morphometric estimation of sediment transfer from the chute to the unit bar. However, width averaged transport rates by downstream sediment budgeting between sampling dates differed from width-averaged Helley-Smith point samples. Lane et al. (1995) used intensive digital surveying and photogrammetric methods to map with high precision the patterns of erosion and deposition during one diurnal flood in a channel of a braided proglacial stream. Additionally, they were able to estimate the spatial distribution of bed material transport rate using a two-dimensional sediment budget. This required an estimate of upstream sediment supply from bedload sampling, and use of a 2D flow model to predict the direction of transport. Their results correlated to bedload transport rates sampled throughout the reach. They found that the cross-channel variation of sediment supplied to the reach had a strong influence on the spatial distribution of bed material transport within the reach. Powell and Ashworth (1995) measured bedload caught in pit traps across a 21 m wide section of the River Wharfe that included a low-relief mid-channel bar. The bed was strongly imbricated in the right channel, but the left channel was loosely packed. Bedload transport tended to be recorded at lower flow stages in the left hand channel and 139 over the bar tail. Transport was observed on the bar head and the imbricated right channel only at relatively high flow stages. They inferred that particle imbrication limited transport in the right channel. However, sediment transported in the left hand channel at lower flow stages was relatively fine, and may have been supplied from upstream. McLean and Church (1999) used morphometric and sediment budget methods to study several decades of channel change in a 45 km length of the wandering gravel-bed reach of Fraser River. They were able to quantify the transfer of sediments from zones of erosion to zones of deposition, and resulting changes in channel alignment that initiated new points of erosion. Erosion was associated with bank attack, and deposition with a downstream zone of flow divergence. The points of erosion propagated along the river faster than the sediment. 4.4 Methods 4.4.1 Data Collection Bottom tracking data were collected with a 1.5 MHz SonTek ADP™ linked to a DGPS. The ADP also measured water column velocity profiles and depths for each ping. Water column velocities were collected along three beams in bins vertically spaced 25 cm apart, with velocity in each bin based on backscatter from a 50 cm deep triangularly weighted window (SonTek 1998). The boat was driven somewhat randomly in transverse lines across the channel, with the intention of providing even spatial coverage, and 140 measurements were collected as single pings at a nominal rate of 0.5 Hz. The bottom track pulse length was set to a relatively long 60 cm in an attempt to improve single ping accuracy, but as a result the depth above the bed where suspended scatterers could influence the bottom tracking was 27 cm (see Section 1.3.1). Measurements were taken in both a gravel-bed and a sand-bed location of Fraser River, British Columbia (Figure 1.4) during the 2001 spring freshet. The gravel-bed location was immediately downstream of the Agassiz-Rosedale bridge (49.21° N, 121.78° W, Water Survey of Canada former gauge station 08MF035) in the wandering reach (McLean et al. 1999, Rennie et al. 2002). This is the same location as the 2000 freshet measurements described in Chapter 1. Measurements were taken on June 27 when flow was 5400 m3/s, which was less than the mean annual flood of 8800 m3/s, but marginally greater than the 5000 m3/s threshold for gravel entrainment (McLean et al. 1999). The ADP could bottom track only in depths ranging from 2 m to 8 m, and thus the measurements were confined to the centre of the channel. Measurements were also collected in the thalweg using a 300 kHz RDI ADCP, with some overlap between ADCP and ADP measurements. However, mean observed bedload velocity was zero with the ADCP, perhaps because the lower frequency 300 kHz sound pulse could penetrate through the mobile layer to the stable bed. It is well established that signal attenuation in sonars is inversely related to acoustic frequency. The sand-bed location was in Sea Reach (49.10° N, 123.16° W), which is an estuarine distributary of Fraser River (Figure 1.4). Flow in Fraser Estuary is tidally influenced, 141 with a maximum tide range approaching 5 m (Kostaschuk et al. 1989). Maximum depth in Sea Reach at low tide during moderate freshet is about 5 m, thus the ADP was capable of bottom tracking throughout the reach, except near the banks. Maximum transport is observed in Fraser River Main Channel during high river discharge, ebb tide, and unstratified flow (Kostaschuk et al. 1989). Measurement of bed material transport in an estuarine environment is particularly challenging, as flow conditions dynamically change in space and time. In Sea Reach, bedload transport diminished to near zero at slack tide. We saw no evidence of salt-wedge intrusion in the relatively shallow Sea Reach, nor is a salt-wedge expected in the Main Channel upstream of Steveston during freshet (Thomson 1981). Also, flow reversal due to the tide does not occur during freshet (Thomson 1981). Measurements were taken in a 300 m wide straight reach between Westham Island and Woodward Island (Figure C.8 in Appendix C). Bed material was sampled at nine locations between June 18th and June 22nd (Figure 4.1a). The bed material varied across the channel, being somewhat coarser in the centre. On channel margins and towards the left side bed material dso was about 0.20 mm, whereas dso was about 0.25 mm in the centre of the channel. Bed material towards the left side of the channel included 5% by weight silts and/or clays. Unfortunately, there were no bed material samples collected in the centre-right portion of the channel. The data presented here were collected on June 21st between 10:12 a.m. and 12:13 p.m. PDT during ebb tide. Low tide occurred at about 1 p.m. in the channel, assuming a 45 142 minute delay for low tide at the measurement location with respect to the Straight of Georgia (Steveston Marine and Hardware 2001, Thomson 1981). This tidal lag matches our observations on June 20th. Flow at the non-tidal Hope gauge 151 river km upstream (WSC gauge 08MF005) on July 21st was 5100 m3/s. The tidal range was 4.1 m, with low tide only 0.3 m above the lowest normal tide (Steveston Marine and Hardware 2001), thus flow through the reach was vigorous during the period of measurement. Water level in the channel dropped 0.9 m during the measurements. Depth averaged water velocity spatially averaged across the section was estimated to diminish from about 1.1 m/s in the first hour of measurement to about 0.8 m/s in the second hour. The direction of flow within a specific location also varied with time, presumably in response to evolving bedforms and changing river stage (Figure 4.1b). At 10:34 a.m. boils were observed at the water surface towards the right side of the channel, indicating that bedforms and bed material transport were active. The largest water velocities (>1.5 m/s) were observed in the downstream portion of the right side of the channel. Large velocities were also observed in the upstream left side. The smallest velocities (<0.5 m/s) were observed in the downstream left portion of the channel. Six point measurements of bedload transport rate were collected across the sand-bed section from the right bank to the left bank in the second hour (11:12 a.m. to 12:13 p.m). Five minute samples were collected using a standard Helley-Smith (HS) sampler, with a 76 mm by 76 mm aperture. A further 64 concurrent ADP and bedload samples were collected in the study reach between June 18 and June 22. Many of these samples were collected at two stationary positions. 143 4.4.2 Data Processing Vector maps of raw bedload velocity were plotted. It was immediately evident that the data were noisy and some sort of smoothing filter would be required. Two approaches were adopted: 1) the vectors were block averaged, and 2) evenly spaced grids of both the primary and secondary bedload velocity components were interpolated using ordinary kriging. Kriging was accomplished using Surfer™, with model variograms determined by fitting to empirical variograms derived from the data. A 2.4 m grid spacing was used for kriging. The primary direction was chosen to be the direction that minimized the reach average secondary bedload velocity component, which was essentially parallel to the river banks. In the sand-bed reach, it was necessary to krige the primary and secondary components, as opposed to the east and north components, because the east and north components were correlated as the streamwise direction did not follow a cardinal direction. A bedload velocity map was plotted from the interpolated primary and secondary components. Beam depths for each ensemble were corrected for instrument tilt and averaged (see Section 1.4.4). Depths were interpolated by kriging onto a 5 m spaced grid and contour plotted. In the sand-bed estuarine reach, prior to contouring, the depths were corrected for tidal fall measured in the reach during data collection, such that all depths were based on river stage at the beginning of data collection. The bed topography was qualitatively compared to the estimated spatial distribution of bedload velocity. In the sand-bed reach, the ADP depth measurements were not of sufficient precision to map the dune field, even 144 if individual beam depths and a fine grid were used. This was because the dunes were small. An ADP track up the centre of the channel immediately following data collection showed small dunes with amplitudes of about 10 cm, wavelengths of between 8 and 14 m, in water depth of around 4 m. Dune heights may have been greater during the period of measurement while the tide was ebbing, but maximum dune heights observed in the reach between June 18 and June 21 did not exceed 20 cm. Shear velocity throughout the reach was estimated using the law of the wall (see Equation 1.20). In the sand-bed reach, only the first hour of data collection was used. A moving window for velocity and depth was utilized to ensemble average 15 adjacent profiles. This should have reduced the standard deviation of velocity estimates to 6 cm/s (SonTek 1998, see Equation 1.5). This averaging occurred in both space and time, as the boat had an average velocity of 0.64 m/s in the sand-bed reach in the first hour of measurement (0.47 m/s for gravel-bed reach), and would have moved an average of 19 m in the 30 s required to collect 15 profiles. Averaging a greater number of profiles would have increased the likelihood that the depth changed significantly during the averaging interval, particularly near channel margins. Further filtering was required to remove profiles that produced unreasonable roughness (ks) values. In the gravel-bed reach, ks should have been around D50 of the surface material (0.042 m) (Yusuf 2001). We accepted profiles that produced ks values within two orders of magnitude of 0.04 m; larger and smaller values were rejected. For the sand-bed reach, an empirical predictor for ks in sand-bed channels with dunes was used (Van Rijn 1984b): (4.1) where A is the average dune height and X is the average dune length. Conservatively assuming A equaled 0.10 to 0.20 m, X was 10 to 20 m, and D90 was 0.3 to 0.4 mm, ks should have been between 0.01 and 0.09 m. Further, ks for flat beds is reported to range from \25D35 to 5ADs4 (van Rijn 1982). Conservatively assuming D35 is equivalent to our finer bed material samples at 0.17 mm, the lower bound of ks should be 0.0002 m. Thus, we ignored shear velocity estimates if ks exceeded 0.2 m or was less than 0.0001 m. Shear velocity showed higher variance in shallow depths of 2 to 3 metres, presumably due to the relative lack of points in the vertical, and increased instrument noise near the operating limit of the instrument. The near-bed velocity for each profile was measured in a bin centered between 25 and 50 cm from the bed. A vector plot of interpolated near-bed velocity was determined using kriging, using the same grid points as the bedload map. The bedload and near-bed velocity maps were statistically compared using a vector correlation coefficient, which is independent of the choice of coordinate system (Crosby et al. (1993)): P v - T V ^ n ) " 1 ^ ^ ) " 1 ^ , ] (4.2) where E,y is the 2x2 covariance or cross-covariance matrix of vectors i and j. Crosby et al. (1993) provide an expansion of (4.2) that is convenient for computation. For two-dimensional vectors, the value of pv ranges between 0.0 and 2.0. Use of this correlation coefficient requires that each vector pair is independent of all others. As the data were spatially autocorrelated, it was necessary to use only vector pairs with sufficient spatial lag to eliminate potential for spatial autocorrelation. The variogram model used in the kriging provided information on the spatial lag required for adjacent vectors to be 146 uncorrelated. For large sample sizes (n>64) and most realistic distributions for vectors i and j, n*pv2 approaches %2 with four degrees of freedom (thus n*pv2 = 9.488 at a = 0.05), which allows for significance testing (Crosby et al. 1993). The sample volumes of near-bed water velocity and bedload velocity overlapped, due to the relatively long bottom track pulse length of 60 cm. Accordingly, the vector fields of depth averaged water velocity and bedload velocity were also statistically compared using the above procedure. 4.5 Results 4.5.1 Sand-Bed Reach The channel topography in the sand-bed reach is generally triangular in section, with the thalweg closer to the left bank (Figure 4.1a). There is an anomalous shallow depth in the upstream left portion of the channel, around (488445, 5438840). This may be the crest of a large bedform. The vector plot of the raw single ping bedload velocity data for the sand-bed reach from the first hour of data collection is presented in Figure 4. Id. The blue line represents the path of the boat, with the starting position at 488410 m Easting, 5438828 m Northing. It is encouraging that most of the bedload velocity vectors point in the downstream direction. It is also evident that larger bedload velocities were observed in the upstream left portion and the downstream right portion of the channel. Block averaging the vectors into a 5x5 grid of 43 m by 48 m spatial blocks yields Figure 4.1e. Block size was selected to ensure that sufficient spatial averaging occurred to incorporate several dunes in each block. This precluded the possibility of differential transport between positions on a dune obscuring larger scale spatial trends. Primary bedload 147 velocities were determined by resolving the east and north components of each ping into the downstream direction, which was NW at 310°. Significant differences were observed between block averaged mean primary bedload velocities, despite the rather large variance of each spatial average (see Appendix E for statistical comparisons). It appears that the ADP was able to resolve differences in transport between different locations in the channel. Interpolating bedload velocities onto a 2.4 m by 2.4 m grid using kriging gives the vector plot in Figure 4.If. The same trends are apparent, with high bedload velocities occurring in a band from the upstream left portion to the downstream right portion of the channel. Near zero bedload velocities are observed in the downstream left. Estimated velocity magnitudes range from 0.008 to 0.906 m/s. The kriging standard deviations of the velocity vectors range from 0.28 m/s to 0.38 m/s, with lower standard deviations at locations that had raw data. The standard deviations are less than the range of the estimates, thus spatial differences in the reach are marginally significant. The bedload transport rate measured at six distinct points across the channel during the second hour of data collection are plotted as stars with associated values in Figure 4.If. The bedload transport rate measurements tend to corroborate the observed pattern of bedload velocity, although comparison is awkward because the HS was usually located in transition zones between high and low velocity. A calibration curve (Figure 4.1h) for measured bedload velocity to bedload transport rate from the 70 HS samples collected in the reach suggests that: v = 0.02gb, where v is 148 primary bedload velocity [m/s] and gb is bedload transport rate [g/s/m]. This relation has been used to scale the vectors in Figure 4.If. We should acknowledge, however, that the largest mean v observed during the 5 minute HS samples was only 0.17 m/s. This compares to maximum vector magnitudes of 0.90 m/s for kriged estimates (Figure 4. If), 0.48 m/s for block averages (Figure 4.1e) and 1.6 m/s for raw data (Figure 4.Id). There are two possible explanations for the larger bedload velocities observed with a moving boat: 1) no stationary spot measurements were collected in a zone of maximum transport, 2) bedload transport was more intense at the time of mobile boat measurements, or 3) values measured with a mobile boat are higher than with a stationary instrument. In Figure 4. If we show the locations of the 70 HS samples (plus symbols). It is apparent that some HS sampling occurred near high bedload velocity zones, although bedload velocity may have been higher at these locations during the time of mobile boat data collection. It is tempting to assume that higher bedload velocities are observed with the moving boat because the measurements are instantaneous rather than time averaged. However, the ADP consistently measured high values in certain locations during several sweeps of the reach (Figures 4. Id), and block averaging and kriging (Figures 4.1e and 4.If) represent averages in both time and space. We compare the bedload velocities of Figure 4. If to the near-bed water velocities in Figure 4.1g. Encouragingly, the patches of highest bedload velocity occurred in the zones of high near-bed water velocity. The bedload velocity and near-bed water velocity variograms indicated spatial autocorrelation up to a spatial lag of about 10 m. The empirical pv 2 for a subsample of the raw vector pairs spaced 10 m apart was 0.23, which exceeds the a = 0.05 confidence level (0.06 for n = 168). Thus, the spatial distributions of bedload velocity vectors and near-bed velocity vectors were significantly correlated. Similarly, the vector fields of raw depth averaged water velocity and bedload velocity were also significantly correlated (pv2 = 0.58 > 0.06, n = 168). Finally, the bed shear stress was estimated throughout the reach (Figure 4.1c). Zones of high shear velocity correspond to the zones of high bedload velocity, with a swath of high shear velocity observed from the upstream left side to the downstream right side of the channel. The fine structure in the shear velocity contour plot is the result of differing values in shear velocity determined from adjacent boat tracks taken at different times in the unsteady flow. 4.5.2 Gravel-bed reach We present the same series of plots for the gravel bed reach in Figure 4.2. The measurements were taken in the centre of the channel, with increasing depth towards the south (Figure 4.2a). Recall that the thalweg ran along the southerly left bank. The measurements were taken immediately downstream of a bridge, with bridge piers located at (589090,5451040) and (589060,5450920). The bank lines suggest the river flows due west, but a spatial average of the depth averaged flow indicates mean flow in the measurement area deviated by 8° to the north. This is likely the effect of a secondary current as the river curves around the mid-channel island/bar complex. The raw bedload velocities (Figure 4.2d) tend to the downstream direction, with larger values observed in the areas of deeper flow. The mean and standard deviation of the raw bedload velocity 150 resolved in the downstream direction were 9.4 cm/s and 26.3 cm/s, respectively. The boat starting position was (589059,5450996). The 5x5 block size was 22 m x 25 m (Figure 4.2e). Only the bedload velocity in the most southerly block was significantly different (a = 0.05) from bedload velocity in other blocks (Figure 4.2e). Kriging the bedload velocity and near-bed velocity spatial distributions (Figures 4.2f and 4.2g) was less successful in the gravel-bed case, as neither distribution was spatially autocorrelated. Accordingly, the kriging standard deviation magnitudes of the bedload velocity vectors were relatively large, ranging from 31 cm/s to 33 cm/s. The bedload and near-bed velocities were interpolated onto a 1.25 m x 1.25 m grid to reveal more spatial structure, but the lack of confidence in individual bedload velocity vectors limits interpretation. Subsampling was not required when comparing the raw bedload and near-bed velocity distributions due to the lack of spatial autocorrelation. The two distributions were weakly but significantly correlated (empirical pv2 = 0.062 compared to significant pv 2 of 0.024 for n=402, a = 0.05). Similarly, the raw vector fields of depth averaged water velocity and bedload velocity were significantly correlated (pv2 = 0.165 > 0.024, n = 402). As in the sand-bed reach, the bedload velocity vectors in Figure 4.2e have been scaled for bedload transport rate. The calibration curve for the gravel-bed reach was presented in Figure 1.7. The maximum observed time averaged bedload velocity was 0.11 m/s. Again, the time averaged bedload velocity is less than the maximum bedload velocity observed with a mobile boat, but in this case no bedload samples were collected near the zone of high bedload velocity (Figure 4.2f). 151 The shear stress distribution (Figure 4.2b) is the most uncertain. In order to ensure reasonable depth estimates for all profiles, some erroneous depths were eliminated prior to moving window averaging. The moving window averaged water velocities produced reasonably semilog-linear profiles with estimated shear velocities ranging from 0.06 m/s to 0.25 m/s with a mean and standard deviation of 0.12 m/s and 0.03 m/s, and standard errors of shear velocity estimates between 0.007 m/s and 0.08 m/s with a mean of 0.02 m/s. Interpolating shear velocities onto a grid yielded a fairly uniform map with most locations showing values between 0.10 m/s and 0.12 m/s. The standard deviation of kriged estimates ranged from 0.024 m/s to 0.041 m/s, thus there is little certainty in the spatial trend. Higher values of estimated shear velocity (0.15 m/s) were observed in the deepest water to the extreme south. This was also the location of highest observed bedload velocity (Figures 4.2c,d,e). However, another anomalous location was the shallowest extreme north location, which displayed shear velocity estimates of 0.18 m/s. Despite smaller average velocity, it appears that velocity gradients were high in this shallow location. This may have been related to flow convectively accelerating as it passed up onto the bar, although we expected the locus of maximum shear stress to shift towards the outside bend thalweg as flow passed around the island/bar (Dietrich and Whiting 1989, Bridge and Gabel 1993). High bedload velocities were not observed in this zone of apparent high shear stress. 152 4.6 Discussion We have demonstrated that the ADP is capable of measuring bedload velocity from a moving boat with sufficient reliability to allow coherent mapping of bedload velocity fields. Such velocity fields are useful for characterization of sediment transport processes and as calibration data for morphodynamic river models. Estimation of the velocity fields required smoothing by either block averaging, which limits spatial resolution, or kriging, which introduces a quantifiable degree of uncertainty in the estimated spatial distribution. Based on the results presented here, it appears that bedload velocity mapping is more reliable in sand-bedded than gravel-bedded rivers. This was expected due to the temporally and spatially sporadic nature of bedload transport in gravel-bed rivers, with a corresponding need for intensive sampling to characterize the bedload transport field. 4.6.1 Spatial averaging in a dune field Typically, higher transport rates are observed near the crest on the stoss side of dunes (Kostaschuk and Villard 1996, Carling et al. 2000). Based on ADP measurements taken while motoring straight up the channel in the sand-bed reach between June 18th and June 22nd 2001, we found that measured bedload velocity tended to be greater on the stoss side of dunes, but high variability precluded finding significant differences. In sand-bed applications, the assumption of spatial homogeneity required to spatially average information from three beams may be violated by spatial differences in transport over a dune. If individual beams impinge on different parts of a dune, then the assumption of 153 spatial homogeneity will not hold. This depends on the relative values of the depth of bed elevation below the ADP (d) and dune length (A.). The three beams impinge within a circular area with a diameter equal to 0.96d+0.197 (Equation 1.4), thus the sampling area has diameter of about d. In our case, the maximum d was about 5 m, for a sampling area diameter of 5 m. The small dunes in this study had wavelengths of about 14 m, thus it is possible that different parts of the dune were averaged together during a single measurement. Empirical relations suggest that equilibrium dune length should be 5 to 7 times depth (Allen 1984, Gabel 1993), which is supported by a theoretical prediction of \=2ndtotai (Yalin 1977). Thus, in an equilibrium dune environment, the sampling diameter should be less than the dune length and some spatial coherence should be possible. We have already mentioned that block size was chosen to ensure that several dunes were included in each block, in order that large scale trends would be observed and not spatial differences in transport over a dune. The patterns observed by kriging are more likely to incorporate both large scale trends and local differences. Kriging searches over a specified area to combine neighbouring measurements when estimating a value for a location. Thus, a search area large enough to incorporate several dunes can be specified. However, kriging gives greater weight to nearby samples, according to the pattern of spatial autocorrelation modeled in the variogram. Thus, local estimates will more likely reflect adjacent measurements. If the adjacent measurements are all located on one part of a dune, say the dune stoss, then the kriged value will be local in nature, and the spatial distribution of bedload transport will be patchy. If the adjacent measurements represent 154 an entire dune or more, then a value reflecting the average for the dune will be estimated, and the spatial distribution will smoothly reflect large scale trends. In our case, spatial autocorrelation up to a lag of 10 m was observed over a dune field with wavelength of about 14 m, thus most a of dune will have been included in each kriged value. 4.6.2 Heading error bias It is encouraging that measurements collected on June 22nd in Sea Reach show high bedload velocities in the same locations as the June 21st data (Figure 4.3). The highest transport velocities observed on June 22nd were further upstream than the measurement area of June 21st. However, the June 22nd data revealed an unexpected bias in measured bedload velocity, which contaminated the estimated bedload velocity distribution (Figure 4.3b). The east bedload velocity component was dependent upon the direction of travel of the boat (Figure 4.3c), with bedload velocities tending towards the west if the boat was travelling downstream (NW) and tending towards east when travelling upstream (SE). The result was that when the boat was travelling NW, the bedload velocity had a counterclockwise bias towards south, and when the boat was travelling SE the bedload velocity had a counterclockwise bias towards north (Figure 4.3d). The depth averaged water velocities also showed this bias. This does not appear to be a problem related to instrument tilt. We corrected the bedload velocity for instrument tilt using the method of Alderson and Cunningham (1999) modified for a three-beam system, with negligible impact on results. However, the bias may have been associated with an error in ADP heading (measured by an internal compass). The boat was facing more-or-less upstream at all times, with the ADP x-axis (beam 1) pointing due east (90°, uncorrected heading) 155 on average (Figure 4.3.e). Heading is measured clockwise from true north (accounting for magnetic declination of 19.6°). There is a slight bias in estimated bedload velocity direction with measured heading angle (note the slight negative slope for the majority of the data in Figure 4.3.e), which was not immediately apparent. We optimized a simple coordinate rotation algorithm for VBT that corrected the measured bedload velocity for an unknown constant heading error. The algorithm first converted the bottom track data to xyz coordinates using the measured heading and the tilt correction algorithm, and then estimated the correct velocities in ENU coordinates using the corrected heading (measured heading plus the heading error). Optimization was achieved by assuming that, according to (1.1), the slope of a regression between VBT and VDGPS should equal one (with the intercept equal to v). For a given heading error, we regressed the new estimated VBT ' versus VDGPS for both the east and north components. The optimum heading error minimized the sum of the absolute differences between the regression slopes and a value of one. We found that the optimum heading error adjustment was -10.0°, i.e. the ADP internal compass appeared to overestimate by 10°. Adjustment of the bottom tracking data using this heading error reduced the regression slope for the north component from 1.12 to 1.00. However, the regression slope for the east component was 0.86, which was only a slight improvement from the original value of 0.79. Importantly, the apparent bias related to boat trajectory angle is reduced (compare Figures 4.3d and 4.3f), and the bias related to ADP heading angle is eliminated (compare Figures 4.3.e and 4.3.g). Kriging the corrected data yields a more reasonable spatial distribution of bedload velocity, although some distortion is still apparent in the vector field (Figure 4.3h). 156 We have noticed a similar heading error bias in another data set when the boat travelled in upstream and downstream lines. The second data set was collected in a nearly quiescent backwater area of Minto side channel in the gravel-bed reach of Fraser River on April 24 2002. The bed was immobile, but the river was turbid with suspended load. The data were collected with a mean boat velocity of 0.74 m/s in depths ranging from 1.9 to 5.2 m (mean depth of 4.0 m). The bed surface was silty fine material. The mean ADP heading in this case was 85°. The Minto data set showed a similar directional bias (Figure 4.4a), but the directional bias was eliminated by an ADP heading adjustment of -8.5° (Figure 4.4b). In this case the regression slopes changed from 1.0 to 0.98 for the east component and from 0.72 to 1.0 for the north component. Following the heading error correction, the mean measured bedspeed was zero, with standard deviations in the east and north directions of. 0.099 m/s and 0.071 m/s, respectively. We did not observe similar heading error bias in the gravel-bed reach nor the June 21st sand-bed data presented earlier in this chapter. Nor was a heading error bias observed in data collected on the same day in the cross-stream direction at Minto. Thus, the problem does not appear to be associated with poor compass calibration. We do not know if this heading error was related to collecting data in upstream-downstream lines, or if this was coincidental. 4.6.3 Bottom track error The Minto estimated standard deviations of between 0.07 and 0.1 m/s give an estimate for mobile boat bottom tracking error in Fraser River in the absence of bedload transport. These standard deviations compare to stationary boat standard deviations in the east and 157 north directions of 0.05 and 0.04 m/s for data collected immediately following the mobile run at a 4.8 m deep location within the mobile test area. We tested mobile ADP bottom tracking over a stable bottom in the towing tank described in Chapter 2, with better results (see Table 1.1 and Figure 1.3). For a tow velocity of 0.74 m/s, the lab data suggest the standard deviation should have been about 0.02 m/s (Table 1.1) for the mobile boat data. Presumably, the irregular bottom, turbid water, random ADP motions, and DGPS error experienced in the field increased measurement error. A mobile bottom greatly increases bottom track error (see Chapters 1, 2, 3). 4.6.4 Mobile versus stationary measurements We noted in Section 4.5 the possibility that higher bedload velocities were recorded during moving boat applications than while the ADP was stationary, particularly for the sand-bed reach. This would be problematic, as it will be very difficult to calibrate the method without concurrent sampling of bedload by some other means such as physical sampling or video. Encouragingly, the mobile boat data collected in the gravel-bed reach in 2001 during a flow of 5400 m3/s did not show a large increase in v compared to stationary measurements collected in 2000 (see Figure 4.2f for locations). The most southerly stationary sampling point 07263 (at 589055,5450989) had a v magnitude of 0.079 m/s in 2000 during a flow of 5600 m3/s (see Table 2.2), whereas the mean bedload velocity vector from the three adjacent blocks (Figure 4.2e) had a magnitude of 0.146 m/s in 2001. The 2000 sample 07183 (at 589072,5451049) had a v magnitude of 0.041 m/s when flow was 6200 m3/s, which corresponds to the extreme north-east 2001 block that had a magnitude of 0.042 m/s. The 2001 spatial data were collected with a 60 cm pulse 158 length, but the 2000 time averaged data were collected with a 20 cm pulse length. An increased measured bedload velocity is expected for the spatial data simply due to the increased pulse length (see Section 5.2.2). Thus, similar values were observed during mobile and stationary measurements at the same locations during comparable flows. A better test would be to compare contemporaneous stationary measurements and nearby mobile boat data. We collected limited transect data during the 2000 freshet, and one transect from July 18 intersected station 07183, which was collected on the same day. For these mobile boat data, the 20 pings nearest to station 07183 were within 16 m of the station. These 20 pings had a mean (and standard deviation) primary bedload velocity of 0.116 m/s (0.176 m/s), which compares to 0.041 m/s (0.120 m/s) for station 07183. The means are significantly different by ANOVA (a = 0.05, p = 0.007) despite the small sample of mobile boat data. It appears that larger bedload velocities were observed with the mobile boat. We are uncertain why v collected from a mobile boat may be higher than from a stationary point. Still, it is possible that the ADP is biased toward zero for stationary measurements, or that mobile boat measurements are positively biased. For stationary measurements, the beam velocity will be close to zero, as it will depend on only the bedload velocity resolved along the axis of the beam. For mobile boat applications, the magnitude of each beam velocity will be much greater than zero, as it will include the velocity of the boat. It may be that the ADP accuracy is better for higher velocities than lower velocities, although our laboratory tow tank tests did not find this to be true (see 159 Table 1.1 and Figure 1.3). Further, the error modelling of Chapter 3 suggests the stationary measurements are unbiased. It would be worthwhile to conduct tests specifically to ascertain the difference, if any, between bedload velocity measured from a moving boat versus stationary measurements. The mobile boat data should be collected as a long time series within a small area centered on the stationary point. 4.6.5 Kriging assumptions We used ordinary kriging to interpolate the vector fields. Ordinary kriging (Matheron 1971, see, for example, Isaaks and Srivastava 1989) produces the best linear unbiased estimator (BLUE) of a spatial distribution. The estimate is linear as it is a linear combination of available data, unbiased in that the mean of residuals (prediction errors) is zero, and best according to the Gauss-Markov theorem that the variance of the residuals is minimized (Isaaks and Srivastava 1989). In ordinary kriging, the variogram is used to calculate weights to be assigned to values at adjacent points when determining the estimated value at a point. The variogram is simply a form of the spatial autocorrelation function. However, ordinary kriging requires the assumption of data stationarity. Data stationarity means that the joint distribution between two points depends only on their separation distance, not their location (Isaaks and Srivastava 1989). In other words, the mean and variance of the data are uniform throughout the sample area. In the derivation of the ordinary kriging equations, the assumption of data stationarity is made three times: 1) when deriving the unbiasedness constraint that the sum of the kriging weights for the linear combination of adjacent values within the search area of a particular location should equal one, 2) when deriving the form of the error variance to be minimized, and 3) 160 when estimating the variogram, which is used to define terms in the minimized error variance (Isaaks and Srivastava 1989). Data stationarity is rarely true in most applications, yet stationarity is assumed. In our case, the variogram was derived from the entire data set, yet we have shown that the data deviated from stationarity with statistically significant differences between the means in subpopulations (see Figures 4.Id and 4.2d). Thus, the kriging estimates may not be BLUE. Data with trend in the mean can be detrended prior to calculation of the variogram (universal kriging). This was not helpful in our case: neither the near bed water velocity nor the bedload velocity distributions had sufficient trend to affect the variograms. In the sand-bed reach the flow was unsteady, with the result that closely spaced data taken at different times were dissimilar. Some investigators have interpolated unsteady phenomena by developing spatiotemporal kriging equations, which account for variability in both space and time (e.g. Rouhani and Myers 1990, Bentamy et al. 1996). The variogram is specified for both a spatial and a temporal lag. We adopted a simpler approach and considered data from only the first hour of measurement, which allowed for a more reasonable assumption of steady flow. However, flow unsteadiness during the first hour had the effect of increasing the variance and reducing the spatial autocorrelation of the data. We did not lose a great deal of spatial coverage by ignoring the second hour of measurement, as we were measuring at stationary points during the second hour. 161 An advantage of kriging over other interpolation methods is that, because the error variance is minimized in the algorithm, the standard deviation of the estimated value (CTR) is calculated. It is important to recognize that CTR is hot the standard deviation of model estimates from the true values, but rather is the square root of the difference between the model variance (the variogram sill, which is usually the data variance) and the model sum of weighted covariance between adjacent points and the estimate location. In effect, OR is an estimate of the amount of information the model was able to derive from nearby locations when assigning a value. If there is little spatial autocorrelation in the data, then CTR will be close to the standard deviation of the data, and we can surmise that there is little certainty in the estimated spatial pattern. This was the case for the gravel-bed bedload velocity and shear stress distributions. In this chapter we have used kriging to interpolate velocities, which are two-dimensional vector quantities, by decomposing the velocities into scalars and kriging the scalar quantities. Bentamy et al. (1996) adopted an approach similar to ours by kriging wind velocity magnitudes and then separately kriging the east and north components to interpolate direction. However, it is possible that the vector information is lost in the decomposition into scalar components. In particular, the choice of coordinate system may affect the results. One possible method for vector kriging is cokriging, which utilizes the cross-covariance between multiple variables to yield improved estimates of the spatial distribution of each variable (e.g. Cressie 1993). However, our bedload velocity and near-bed velocity data showed little correlation between the decomposed scalars, thus cokriging would not have improved the estimates. Furthermore, in the sand-bed reach we 162 resolved the velocities into a downstream and cross-stream component, which eliminated cross-correlation. Interestingly, in the gravel-bed reach there was less cross-correlation (p = 0.05) using an east-north decomposition (which followed the full reach downstream-cross-stream directions) than with a downstream-cross-stream coordinate system determined by the direction of transport in the measurement area (which had p = 0.11). Feliks et al. (1996) proposed a vector interpolation scheme utilizing the Gauss-Markov theorem that is similar to cokriging. In their formulation, it is necessary to detrend the spatial distribution by removing the local mean, but this should be unnecessary if stationarity is assumed. Feliks et al. (1996) proved that their method is independent of the choice of coordinate system. Young (1987) proposed a vector kriging procedure that is independent of the coordinate system. He reduced the vector problem to a scalar by using the magnitude of the difference vector as the kriging variable for estimation of the variogram and minimization of the residual variance. In effect, a single variogram, and thus a single set of interpolation weights, is determined for both velocity components. Inggs and Lord (1995) applied the method to interpolation of wind field data. A comparison of the various methods would be useful. There are other details in the kriging procedure that could be optimized to improve the results. In particular, it may be useful to look for outlier pairs in each spatial lag prior to fmalization of the variogram (Isaaks and Srivastava 1989). It may be that outliers precluded finding spatial autocorrelation in the gravel-bed reach. Also no optimization was performed for the kriging search area that specified which neighbouring points were 163 included in the weighted linear summation. It may be that the search area could be optimized by evaluating the effect on O"R. 4.7 Conclusions Maps of bedload transport velocity have been produced for sand-bed and gravel-bed reaches of Fraser River. Measured bedload velocity spatial distribution was significantly correlated with near-bed water velocity and depth-averaged water velocity, and compared favorably to spatial distributions of bedload transport rate, depth, and bed shear stress. The bedload transport rate distribution was based on physical point sampling. The sand-bed reach showed a distinct swath of higher bedload velocities. The data in the gravel-bed reach were more erratic, as expected for relatively infrequent partial gravel transport, but showed significantly greater bedload velocity in the deepest location. Developing vector plots of bedload transport velocity with kriging was complicated in the sand-bed estuarine reach by unsteady conditions and in the gravel-bed reach by a lack of observed spatial autocorrelation, thus kriging assumptions were violated in both cases. 164 Figure 4. la: Sand bed reach depth (m), based on flow depth at beginning of data collection. Bed sediment sample locations marked with +, and labelled with d50. £1 era' O cn H * &> 5* CD p a °* CD >-t a CD o H ^ B . ^ ^ p < c= 8 •< o S I 3 o o a 5 5 r £ era a °* % s-° CJ. >-+> t r a * p O i-t ?r P • era co < 2-o* o rt' 5' cr CD tfl 1-1 V) rt cr o e •1 o 3 CO p cn a CO 3 CO a cr cr co cn CO o o a Northing, m cn Ul cn cn 4*. 4s. 4S. •N CO CO CO CO CO GO CO GO CD CD 00 CD cn O a CJl O cn O cn O o a O a Northing, m cn cn cn cn cn .u 4*. 4s. -ts. 4s. CO CO CO CO CO 00 CO 00 CO CO ao CO CD o a cn o cn o cn o a o a 4s. CO CO 4s. O O 4=* CO m g 0) cn =3 4S. CO CO cn o o 4s. CO CO cn cn o 1 / ••! —r — / / / / cr ' . / ' / / Scale 1 m/s 1 1 / 166 488400 488450 488500 488550 Easting, m Figure 4. lc: Sand bed reach bed shear velocity (m/s) distribution in first hour of data collection. Location of data collection points shown with small + symbols. 167 4 8 8 4 0 0 4 8 8 4 5 0 4 8 8 5 0 0 4 8 8 5 5 0 Easting, m Figure 4. Id: Sand bed reach raw bedload velocity (m/s) from first hour of data collection, prior to stationary bedload sampling. 168 i o «o i o »o o o o d o o o o o 5439050 5439000 5438950 O) c o 5438900 5438850 \ 4.2 \ \ It Reference Vectors 0.07 m/s 0.54 m/s 488400 488450 488500 488550 Easting, m Figure 4. le: Sand bed reach block average bedload velocity (m/s) in first hour of data collection overlain on depth (m) contours. Centre of vector arrow is centre of block. « ~ 1 n» o 3 p o " Northing, m 2T§ LO 3 CJQ _, — o 3 « -* -t ~t s 8 - a S? 8. cr CO o p + o a. ELS: 2. pr cr s. & 8'E 1 § I P> 2 o 2 as 2. ^ 3 ^ ^ 3 & o- 2, a 8 I = o g n o " 1 Q B O o s a § P 3 8 p BS 2- n> ET 2 a. §• cr 3 6 ^ < P , g c » n> c r 5* cr cr o & S* g C L P • C L CJQ co § U p CD a CL A 00 O l o co o: CO a to 03 CD 01 o m 3 ' CQ *T~* —t 1 T—7" 1: % "t * t f * f / / o o CO CO o cn o t t t f t f t / fst f t A AAAAA A +/ / /+ fj f / % I / / / ? t t / / S A A. A As A A ///// o o 3 * § 8 » O cn cn I* A A.A A A AA Ag^*> s // / A A A / A ? \*SA/*„ A / A S A / f A*s s /'fy A s s / s s 7 S S 7 >r /* -* J> A A ALA A A Afff m ? ^ A . AAA. CO CO CO cn o m CO co as 03 tm o o "7 r—7—7 cn -e^ 00 CO CO cn ? r /1 CO co o o o - ' 'f ft t / / / / f * <*AA ////////> A A A S s / » f - 1 -cn 00 C D o cn o CQ / / / , A S S A / Ass/*/ A A A f ^ A A A* / / -# / Af L I A A LA / / f A ft A t i l 170 14 0.18 Mean bedload velocity in transport direction, m/s Figure 4.1.h: Sand bed reach (Sea Reach, Fraser River) bedload transport rate versus ADP bedload velocity calibration curve. Two outliers (denoted by O) were not included in the regression: it was presumed that these large bedload samples were due to bottom dragging of the Helley-Smith sampler. 171 5451060 5451040 5451020 O) .E 5451000 sz tr o 5450980 5450960H 5450940 i 589000 589020 589040 589060 Easting, m Figure 4.2a: Gravel bed reach depth (m). 589080 589000 589020 589040 589060 Easting, m Figure 4.2b: Gravel bed reach shear velocity distribution (m/s). Location of data collection points shown with small + symbols. 173 5451060 5451040 h 5 4 5 1 0 2 0 h Ui o 5451000 5450980 h 5450960 r 5450940 r 589020 589040 589060 589080 589000 Easting, m Figure 4.2c: Gravel bed reach spatially block averaged depth average water velocity (m/s). Tail of vector arrow is centre of block. 174 1 1 ' U-»l{ , , , J 5 8 8 9 8 0 5 8 9 0 0 0 5 8 9 0 2 0 5 8 9 0 4 0 5 8 9 0 6 0 5 8 9 0 8 0 5 8 9 1 0 0 Easting, m Figure 4.2d: Gravel bed reach raw bedload velocities. 175 5451040 5451020-5451000^ 5450980 5 4 5 0 9 6 0 ; 5450940-6 Reference Vectors 0.009 m/s 0.495 m/s 589000 589020 589040 589060 589080 Figure 4.2e: Gravel bed reach block averaged bedload velocities (m/s) overlain on depth (m) contours. Centre of vector arrow is centre of block. era hfl r-y Op' -1 rt g ft! (A 3 as 3 a 1 CP CD S s. CD **' C L "7 rc 3 p CD O g tr cr ™ CD ^ ° " ffl CT < o CD j a o ?• o < » n ^ o | ' ^ ' r- 68 § EF ? C rt-5' 3 — o o e rt 5' B O —. CD Sr O so c cn B rtj 3 C o c h* -i CD CO — 3 p -i O zr + p p p o p p p p p p o o o o o o o o o o o o o o o w co K>s±a)a ro o> co B i S S iStmrn Northing, m on cn o CO & cn o CO O) o CR o co CO o cn A . cn o o o Ul O KJ O cn •u cn o *. o cm _^ o cn o 111 m as 09 rt cn 3 ' CO CO CO o 3 n* o cn GO CO o 05 O Ul CO to o~ co o 4 * I f 4 f k If f / '•* t ( Tl 5E 4 4 4 4 4 4 4 4 4 4 4 * '4 4 ! 4 / 4 f 4 4 4 "fl 4 ~ ' / '4 4 4 4 1 f 4 / 4 4 4 ,* .» f 4 4 4 f fl fl .4 fl 4 4 4 / » j» / 4 IA 4 / 4 f P 4 4 4 4 f 4 '* f 4 1 < 4 / 4 4 4 4 4 f 4 4 4 4 • t I 4 ' 4. * > 4 4 1 4 4 1 F 1 4 4 4 1 4 ' 1* ' 4 4 4 4 4 1 1 4 t * ^. * a fl 4 ' 4 4J 4 4% 4 - 1 t / / 4) t t .* A 4 * / f * * .4 * ,4 .» f M f 4 4 / 4 4J ,4 * t f 4 f 4 4 < 4 4 fl 4 ,< * t t fl 4 4* 4 * * t t ^ c 4 4 4 ? .< 4 4 4 ^ \ \ f f 4 1 4 4 f ^ \ "t ^ 1 * 4 4 / 4 * * • ^, "c M * 4 4 * If 4 H • + 4 S f » 4 , 4 cn cn O CO o —r~ cn cn O co 03 cn cn o CO CO o o o o cn CR O o m c# CQ cn co CO o CO o 4 ~k\ sr 4 4 4 * I i f ¥ f ¥ f I f t 3D CO 1 rn/s rn/s ence N) 1 < CO CO CD & 3 1 O 1 CO 4 t f f f f J f 4 ! 4 * t f f 4 4, : « 4 4 : cn Ul o S o 1 A k * 4 i 4 — r ~ cn on o 01 o M M M M f t i 4 4 * * 4 * k 4-'• ' 1 "' 4 4 4 4 4 4 4 f 4 7 f i * f * | 1 i 1 4 I 4 * I 4 4 J T 4 4 4 f * | f 4 1 4 4 * ! 5439050 5439000 E 5438950 c 2 5438900 5438850 5438800 1 ' 1 . i_J 488450 488500 488550 488600 488400 Easting, m Figure 4.3a: Sand-bed reach June 22nd raw bedload velocities 5439050 5439000 i . ,. » . * * f 4 » *' —a *—* E 5438950 CD C I o Z 5438900 5438850 5438800-*. V. ft f> v . * • it \ »* *— **~5 4 * • < » t t »>•»:». ft •' t f ? - - - -t f « * f » ' .« ? .»f * *. Reference Vectors * • 0.01 m/s 0.79 m/s 488400 488450 480500 468550 488600 Easting, m Figure 4.3b: Sand-bed reach June 22nd kriged bedload velocities 178 5439050 5439000 E 5438950 5438900 488400 488450 488500 488550 488600 Easting, m Figure 4.3c: Sand-bed reach June 22nd kriged bedload velocity, with boat track shown by red + symbols: a) east component, b) north component. 179 250 200 ® 150 _oT 50 Q . W " O 0) o - Q -50 ++ r + + + + •4 + + + + + + t + + + + + -150 -100 -50 300 0 50 100 150 200 250 boat trajectory, deg Figure 4.3d: Sand-bed reach June 22nd bedload velocity vector angle versus boat trajectory angle. 200 250 0 50 100 150 heading, deg (East is 90) Figure 4.3.e: Bedload velocity angle (counter clockwise angles with East at zero) versus ADP heading angle (clockwise angles with North at zero). 180 _g 150 -50 200 250 300 0 50 100 150 boat trajectory, deg Figure 4.3f: Sand-bed reach June 22nd bedload velocity vector angle versus boat trajectory angle after adjusting heading by -10°. 300 -50 0 50 100 150 200 250 heading, deg (East is 90) Figure 4.3.g: Bedload velocity angle versus ADP heading angle after adjusting heading by-10°. 181 o 5439051 5439000-E 5438950-| ra C i o Z 5438900; 5438850 • . \ \ \—r^ r* 1 V*. \ \ \ t * , . . . » » » • , * t ^ t « « t 1 1 *. «. « fc *%*• * * * * 1 \ f t » * < > • . , . » » » . * » * *. * \» . . »., t ", \ \ \ "• »• -\ *• V «. * * * * * » ». » » ».. ».s > v 1 * V , \ " " «- »-• * *• *• \ M 1 U w v ^ ^ , . . ». ». * * * «. « » * » * ' « * ».. * 5438800 Reference Vectors + * 0.02 m/s 0.91 m/s 488400 488450 488500 488550 488600 Easting, m Figure 4.3h: Sand-bed reach June 22nd kriged bedload velocities following -10° rotation of raw velocities for ADP heading error. 182 -150 -100 -50 0 50 100 150 boat trajectory, deg Figure 4.4.a: Minto side channel bedload velocity vector angle versus boat trajectory angle. 250 F 200 - g 150 r £?100 CO "D a> 0) 50 Q_ V) TD -50 ++ -150 -100 -50 0 50 100 150 boat trajectory, deg Figure 4.4.b: Minto side channel bedload velocity vector angle versus boat trajectory angle after adjusting heading by -8.5°. 183 Chapter 5: Concluding Remarks 5.1 Summary of key results This thesis has demonstrated that an aDcp can yield useful information about the spatial and temporal distribution of bedload transport velocity in rivers. In Chapter 1 it was shown that, for stationary sampling locations in a gravel-bed reach of Fraser River, the mean bedload transport rate measured with conventional samplers correlated with the mean bedload velocity measured by the ADP. The regression fit serves as a calibration for ADP bedload velocity measurements, although the calibration appears to be site specific, presumably related to the particle size distribution of the bedload. In Chapter 4 we briefly presented a similar calibration curve for the sand-bed Sea Reach, which is a distributary of Fraser River. The data from an individual sampling station were extremely noisy. An error model was presented in Chapter 3 to separate (deconvolve) the bedload velocity probability density function from the observed data. The model was successfully applied to stationary gravel-bed data sets from Fraser River and Norrish Creek. The estimated bedload velocity distributions were extremely positively skewed, as expected for partial transport of gravel, and were well fit by both a gamma distribution and a compound Poisson-gamma distribution. The estimated instrument noise was of a similar magnitude to noise for aDcp water column velocity measurements. 184 The technique was applied using a mobile boat in both the gravel-bed and sand-bed Fraser River reaches (Chapter 4). Despite noise in the data, coherent patterns were observed in the bedload velocity spatial distributions. The bedload velocity distributions were significantly correlated to the near-bed and depth averaged water velocity distributions, and were reasonable when compared to the estimated distribution of shear stress in the reach. A calibration was also attempted in a laboratory towing tank (Chapter 2). It appeared that the bedload velocity in the direction of transport was biased toward zero. Instead, the mean magnitude and the mean direction of the bedload velocity vector provided good predictions of the expected bedload velocity. It appeared that noise dominated the laboratory measurements. 5.2 Outstanding Issues 5.2.1 Measurement bias The primary question that remains unresolved is whether or not the measurements are biased. The laboratory data (Chapter 2) tentatively suggest that the velocity resolved in the direction of transport is biased toward zero. Unfortunately, the va in the direction of transport tended to be underpredicted by about 79% (Chapter 2). This may be inconsistent with the gravel-bed field results (Chapter 1): in the field it seemed that va in the direction of transport may have been overpredicted, as the expected transported particle size for the observed va was very small (< 1 mm). It may be that, in the field, near-bed suspended sediment transport positively biased va in the direction of transport. In the laboratory, the expected beam velocity was very small, approaching the lower limit of resolution of the instrument, especially for beams 2 and 3. In the field, the contribution of small near-bed suspended scatterers would have provided a more uniform and higher magnitude signal. Another difference between the field and the laboratory data was that in the field there was little correlation between the va mean vector magnitude and mean bedload transport (Figure 5.1). This makes Equation 1.7 a bit awkward: in the field data the va mean vector magnitude may actually be a measure of variance rather than a measure of mean transport. The coefficient in Equation 1.7 is 1.3 based on the laboratory results of Chapter 2. Based on the results of the error modelling presented in Chapter 3, the coefficient in Equation 1.7 is also 1.3 for the Fraser River data set (CTN = 0.093 m/s, vav = 0.158 m/s, and N = 5), but is only 0.58 for the Norrish Creek data set (ON = 0.31 m/s, vm = 0.533 m/s, and N = 1). The Fraser River data set was collected from a moving boat, and thus the estimated O"N included all the errors discussed in Chapter 1. The Norrish Creek data was collected with the ADP securely mounted to a bridge, and thus the estimated ON was generated only from instrument noise. It appears that the laboratory and field data may have differed in their noise structure. A constant coefficient of variation was observed in the laboratory data, which implies that error scaled with the signal. This is consistent with the notion that the laboratory results were dominated by low signal to noise ratio, which should produce errors that scale with measured velocity (Lhermitte and Lemmin 1999). The field data also moderately suggested an increase in scatter with increasing bed velocity (Figure 5.2, 186 compare to Figure 2.4). However, this could have been due to increased real variability (which is expected, e.g. Kuhnle and Southard 1988)) as well as increased noise with increasing bed velocity. In order to test the assertion that the noise structure in the laboratory differed from that in the field data, the error model of Chapter 3 was applied to the laboratory data (Figure 5.3). A reasonable fit was obtained for observed to modelled data distribution, and the two distributions were not significantly different. However, the negative tail was not well fit. In other words, the convolution of a Gaussian noise distribution with a positive signal distribution could not model the observed negative values. This suggests the errors in the laboratory data differed from those in the field data. It appears that there was noise in the laboratory that biased vfl in the direction of transport toward zero. This was likely the result of excessive noise along individual beams, which rendered the ADP incapable of resolving the direction of transport. This did not appear to be the case in the field data, thus it appears reasonable to assume that the field data were unbiased. In fact, we had to assume unbiased field data in order to obtain the error model fits of Chapter 3. 5.2.2 Effect of pulse length The effect of pulse length on the influence of suspended scatterers on the observed bedload velocity was presented in Chapter 1. Assuming that suspended scatterers will positively bias the measured bedload velocity due to their elevated velocity, it can be hypothesized that a greater bedload velocity will be measured if a longer pulse length is used. In the laboratory it was observed that in the shallow depth (2 m) of the towing 187 tank, the pulse length had to be increased to 60 cm in order for the bottom tracking to accurately measure the tank depth (Chapter 2). This is somewhat inconvenient, as the long pulse length required for accuracy also increases the depth above the bed to which suspended scatterers can influence the bedload velocity. A true measure of bedload should measure only particles moving close to the bed. In the field measurements at Norrish Creek (study site described in Chapter 3), successive stationary time series were collected in a steady flow using a 20 cm pulse length and a 60 cm pulse length. As expected, the mean of the bedload velocity resolved in the direction of transport was greater for the long pulse length of 60 cm (0.403 m/s) than for the shorter 20 cm pulse length (0.241 m/s). The standard deviations of the bedload velocity were similar (0.450 m/s for 60 cm pulse, 0.494 m/s for 20 cm pulse), although, as expected, the longer pulse length had slightly lower variance. The similarity of the standard deviations suggests that use of a short pulse length is preferred. 5.3 Recommendations for Further Research The investigations presented in this thesis suggest that measurement of bedload velocity using aDcps is feasible. This represents an exciting new opportunity for river scientists to elucidate river processes, both through direct measurements and by means of better calibrated morphodynamic river models. However, a great deal more research is required before the method can be conveniently applied by practitioners. In particular, calibration of the technique in various fluvial environments is required. Further, it would be helpful if a controlled mobile bed testing facility (canal or flume) of sufficient depth could be 188 used to obtain synchronous estimates of aDcp bedload velocity and particle velocities by some independent means such as video. Finally, the spatial and temporal resolution of the data will not be improved without development of an improved instrument specifically designed to measure bedload. A short list of recommendations for further research in the field, the laboratory, and in instrumentation development is offered below. Field • Collect a long stationary time series (multiple hours to days) for assessment of temporal trends. In a gravel-bed reach the object of study would be pulsing of bedload. In a sand-bed reach the temporal trend in transport velocity would be due to the passage of dunes. For such a long time series, steady flow conditions could likely only be collected in a large river or in a snowmelt freshet stream. • Collect a spatially intensive survey throughout a long reach of river (several channel units) in order to reveal trends in the spatial distribution of bedload transport that are geomorphically significant. This is probably most suited to a gravel-bed reach with reasonably fixed bathymetry, unless the dynamics of the dune field in a sand-bed reach could also be characterized. • Conduct tests specifically to ascertain the difference, if any, between bedload velocity measured from a moving boat versus stationary measurements. The mobile boat data should be collected as a long time series within a small area centered on the stationary points. • Collect more concurrent bedload transport and aDcp data, for calibration of the method in different fluvial environments. 189 Laboratory • Perform a full gamut of trials for various particle sizes, particle speeds, percent mobile bed area, and aDcp orientation with respect to the direction of transport. Multivariate analyses could determine if the aDcp is more capable of determining bedload transport velocity for certain particle sizes, transport rates, and orientations. • Perform tests with different high frequency aDcps to see if different aDcps yield different results. • Utilize a sufficiently deep mobile bed testing facility (canal or flume), which allows for synchronous measurements of aDcp bedload velocity and particle velocities by some independent means such as video. Improvements to the Method • Investigate the possibility of using the estimated noise distribution, as determined by the deconvolution technique of Chapter 3, to filter error from the bedload velocity time series. At present, the deconvolution technique yields the probability density function of the bedload velocity, but the chronology of the data is lost. Extraction of the time series would allow for improved spectral analysis of periodicity in the data. • Collaborate with an instrument manufacturer to improve the bottom track signal processing for the purpose of measuring bedload transport velocity. • Develop a new instrument designed specifically to measure bedload. In particular, signal processing that considers the entire spectrum may produce more accurate velocity estimates for heterogeneous bedload particle velocities. Bedload velocity 190 error would also be reduced by an increased bottom track pinging rate. Use of backscatter signal strength could be used to estimate particle concentration, which, combined with velocity estimates, will give bedload transport rates. 191 0.00 0.05 0.10 Mean bedload transport rate, kg/m/s Figure 5.1: Mean magnitude of bedload velocity vector versus mean bedload transport rate for Fraser River gravel reach data of Chapter 1. 192 0.50 0.00 0.10 0.20 0.30 Mean primary bedload velocity, m/s Figure 5.2: Standard deviation versus mean of primary bedload velocity for Fraser River gravel reach data of Chapter 1. 193 J I I I I I L -0.3 -0.2 -0.1 0 0.1 0.2 0.3 0.4 0.5 Bed speed in transport direction, m/s Figure 5.3: Probability density function and beam gamma error model fit for bedload velocity resolved in the direction of travel for Trial 11 of Chapter 2 laboratory data. 194 References Abbott, J.E., and Francis, J.R.D. (1977). "Saltation and suspension trajectories of solid grains in a water stream." Phil. Trans. Royal Soc, 284A, 225-254. Adenlof, K., and Whol, E., E. (1994). "Controls on bed movement in a subalpine stream on the Colorado Rocky Mountains, USA." Arctic and Alpine Research, 26(1), 77-85. Alderson, S. G., and Cunningham, S. A. (1999). "Velocity errors in acoustic Doppler current profiler measurements due to platform attitude variations and their effect on volume transport measurements."/. Atmosp. Oceanic Technol, 16(1), 96-106. Allen, J. R. L. (1984). Sedimentary Structures: Their Character and Physical Basis, Elsevier, Amsterdam. Andrews, E. D., and Parker, G. (1987). "Formation of a coarse surface layer as a response to gravel mobility." Sediment Transport in Gravel-Bed Rivers, C. R. Thorne, J. C. Bathurst, and R. D. Hey, eds., John Wiley & Sons Ltd., 269-300. Andrews, E. D., and Smith, J. D. (1992). "A theoretical model for calculating marginal bedload transport rates for gravel." Dynamics of Gravel-bed Rivers, P. Billi, R. D. Hey, C. R. Thorne, and P. Tacconi, eds., John Wiley & Sons Ltd, Chichester, England, 41-48. Ashworth, P. J., Best, J. L., Roden, J. E., Bristow, C. S., and Klaassen, G. J. (2000). "Morphological evolution and dynamics of a large, sand braid-bar, Jamuna River, Bangladesh." Sediment., 47, 533-555. Ashworth, P. J., and Ferguson, R. I. (1986). "Interrelationships of channel processes, changes and sediments in a proglacial river." Geograf. Ann., 68A, 361-371. Ashworth, P. J., Ferguson, R. I., Ashmore, P. E., Paola, C , Powell, D. M., and Prestegaard, K. L. (1992a). "Measurements in a braided river chute and lobe 2. Sorting of bed load during entrainment, transport, and deposition." Water Resour. Res., 28(7), 1887-1896. Ashworth, P. J., Ferguson, R. I., and Powell, M. D. (1992b). "Bedload transport and sorting in braided channels." Dynamics of Gravel-bed Rivers, P. Billi, R. D. Hey, C. R. Thorne, and P. Tacconi, eds., John Wiley & Sons Ltd, Chichester, England, 497-515. Bagnold, R. A. (1973). "The nature of saltation and of "bed-load" transport in water." Proc. Royal Soc. London, A332, 473-504. 195 Bathurst, J. C , Thorne, C. R., and Hey, R. D. (1977). "Direct measurements of secondary currents in river bends." Nature, 269(6 October), 504-506. Bentamy, A., Quilfen, Y., Gohin, F., Grima, N., Lenaour, M., and Servain, J. (1996). "Determination and validation of average wind fields from ERS-1 scatterometer measurements." Global Atmosphere Ocean System, 4, 1-29. Bishop, G. C., and Smith, J. (1999). "Scattering from rigid and soft targets near a planar boundary: Numerical results."/. Acoust. Soc. Am., 105(1), 130-143. Brayshaw, A. C. (1984). "Characteristics and origin of cluster bedforms in coarse-grained alluvial channels." Sedimentology of Gravels and Conglomerates. Canadian Society of Petroleum Geologists, Memoir 10, E. H. Koster and R. J. Steel, eds., 77-85. Bridge, J. S., and Gabel, S. L. (1992). "Flow and sediment dynamics in a low sinuosity, braided river: Calamus River, Nebraska Sandhills." Sediment., 39, 125-142. Brumley, B. H., Cabrera, K. L. D., and Terray, E. A. (1991). "Performance of a Broad-Band Acoustic Doppler Current Profiler." IEEE J. Oceanic Eng., 16(4), 402-407. Busskamp, R. (1994). "The influence of channel steps on coarse bed load transport in mountain torrents: case study using the radio tracker technique 'PETSY'." Dynamics and Geomorphology of Mountain Rivers, P. Ergenzinger and K.-H. Schmidt, eds., Springer-Verlag, Berlin. Callede, J., Kosuth, P., Guyot, J.-L., and Guimaraes, V. S. (2000). "Discharge determination by Acoustic Doppler Current Profilers (ADCP): a moving bottom error correction method and its application on the River Amazon at Obidos." Hydrological Sciences Journal, 45(6), 911-924. Campo, L. (2002). Environment Canada, Pacific and Yukon Region. Carling, P. A., Williams, J. J., Golz, E., and Kelsey, A. D. (2000). "The morphodynamics of fluvial sand dunes in the River Rhine, near Mainz, Germany. II. Hydrodynamics and sediment transport." . Cordy, C. B., and Thomas, D. R. (1997). "Deconvolution of a distribution function." J. American Statistical Association, 92(440), 1459-1465. Chacho, E. F., Emmett, W. W., and Burrows, R. L. "Monitoring gravel movement using radio transmitters." Hydraulic Engineering '94, Buffalo, New York, 785-789. Chacho, E. F. J., Burrows, R. L., and Emmett, W. W. "Detection of coarse sediment movement using radio transmitters." International Association for Hydraulic Research, XXIII Congress, Proceedings of Technical Session B: Fluvial Hydraulics, Ottawa, Canada, August 21-25,1989. 196 Chereskin, T. K., Firing, E., and Gast, J. A. (1989). "Identifying and screening filter skew and noise bias in Acoustic Doppler Current Profiler measurements." J. Atmos. Oceanic Technol, 6, 1040-1054. Chereskin, T. K., and Harding, A. J. (1993). "Modeling the performance of an Acoustic Doppler Current Profiler." J. Atmos. Oceanic Technol, 10, 41-63. Church, M. A. (1987). "Discussion of'Formation of a coarse surface layer as the response to gravel mobility'." Sediment Transport in Gravel-bed Rivers, C. R. Thorne, J. Bathurst, C , and R. Hey, D., eds., John Whiley and Sons Ltd., 314-322. Church, M., Hassan, M., A., and Wolcott, J., F. (1998). "Stabilizing self organized structures in gravel-bed stream channels: Field and experimental observations." Water Resour. Res., 34(11), 3169-3179. Cressie, N. A. C. (1993). Statistics for Spatial Data, John Wiley & Sons, Inc., New York. Crosby, D. S., Breaker, L. C , and Gemmill, W. H. (1993). "A proposed definition for vector correlation in geophysics: theory and application."/. Atmos. Ocean. Technol., 10, 355-367. Davoren, A., and Mosley, M. P. (1986). "Observations of bedload movement, bar development and sediment supply in the braided Ohau River." Earth Surf. Proc. Landf, 11,643-652. de Jong, C , and Ergenzinger, P. J. (1998). "On the nature of spatial and temporal variations of bedload transport in a gravel-bed stream." Z. Geomorph. N.F., Suppl.-Bd. 112, 29-54 DeVries, P. E. (2000). "Scour in Low Gradient Gravel Bed Streams: Patterns, Processes, and Implications for the Survival of Salmonid Embryos," Ph.D., University of Washington, Seattle. Dietrich, W. E., Kirchner, J. W., Ikeda, N., and Iseya, F. (1989). "Sediment supply and the development of the coarse surface layer in gravel-bedded rivers." Nature, 340(20 July 1989), 215-217. Dietrich, W. E., and Whiting, P. (1989). "Boundary shear stress and sediment transport in river meanders of sand and gravel." River Meandering, S. Ikeda and G. Parker, eds., American Geophysical Union, 1-50. Dietrich, W. E., and Smith, J. D. (1983). "Influence of the point bar on flow through curved channels." Water Res. Res., 19(5), 1173-1192. 197 Drake, T. G., Shreve, R. L., Dietrich, W. E., Whiting, P. J., and Leopold, L. B. (1988). "Bedload transport of fine gravel observed by motion-picture photography." J. Fluid Mech., 192, 193-217. Einstein, H. A. (1937). "Bed load transport as a probability problem," D.Sc, Federal Institute of Technology, Zurich, Switzerland, translated by W.W. Sayre in Sedimentation: Symposium to Honor Professor H. A. Einstein, 1971, H.W. Shen, ed., Colorado State Univ., C1-C105. Environment Canada. (1996). "HYDAT Version 96-1.05.8, Surface Water and Sediment Data." , E. C. Atmospheric Environment Service. Ergenzinger, P., and Schmidt, K.-H. (1995). "Single particle bed load transport: First results from new approaches." Z. Geomorph. N.F., Suppl.-Bd. 100, 193-203. Ergenzinger, P., Schmidt, K.-H., and Busskamp, R. (1989). "The Pebble Transmitter System (PETS): first results of a technique for studying coarse material erosion, transport and deposition." Z. Geomorph. N.F., 33(4), 503-508. Feliks, Y., Gavze, E., and Givanti, R. (1996). "Optimal vector interpolation of wind fields."/. Appl. Meteor., 35, 1153-1165. Feller, W. (1957). An introduction to probability theory and its applications, Vol. 1, Wiley, New York. Ferguson, R. I., Ashmore, P. E., Ashworth, P. J., Paola, C , and Prestegaard, K. L. (1992). "Measurements in a braided river chute and lobe 1. Flow pattern, sediment transport, and channel change." Water Resour. Res., 28(7), 1877-1886. Ferguson, R. I., and Ashworth, P. J. (1992). "Spatial patterns of bedload transport and channel change in braided and near-braided rivers." Dynamics of Gravel-bed Rivers, P. Billi, R. D. Hey, C. R. Thorne, and P. Tacconi, eds., John Wiley & Sons Ltd, Chichester, England, 477-492. Ferguson, R. I., and Wathen, S. J. (1998). "Tracer-pebble movement along a concave river profile: virtual velocity in relation to grain size and shear stress." Water Resour. ites., 34(8), 2031-2038. Gabel, S. L. (1993). "Geometry and kinematics of dunes during steady and unsteady flows in the Calamus River, Nebraska, USA." Sediment., 40, 237-269. Gaunard, G. C , and Huang, H. (1996). "Sound scattering by a spherical object near a hard flat bottom." IEEE Trans. Ultrasonics, Ferroelectrics, Frequency Control, 43(4), 690-700. 198 Gibbs, C.J. (1973). "Model study of the basket-type bed-load sampler." M.Sc. Dept. Civil Engineering, University of Alberta, Canada, 152 pp. Gintz, D., Hassan, M. A., and Schmidt, K.-H. (1996). "Frequency and magnitude of bedload transport in a mountain river." Earth Surf. Proc. Landf, 21, 433-445. Goff, J. R., and Ashmore, P. (1994). "Gravel transport and morphological change in braided Sunwapta River, Alberta, Canada." Earth Surf. Proc. Landf, 19(3), 195-212. Gomez, B. (1991). "Bedload transport." Earth-Science Reviews, 31, 89-132. Gomez, B., and Church, M. (1989). "An assessment of bed load sediment transport formulae for gravel bed rivers." Water Resour. Res., 25(6), 1161-1186. Gomez, B., Naff, R. L., and Hubbell, D. W. (1989). "Temporal variations in bedload transport rates associated with the migration of bedforms." Earth Surf. Proc. Landf, 14, 135-156. Gordon, R. L. (1996). Acoustic Doppler Current Profiler Principles of Operation: A Practical Primer, RD Instruments, San Diego, CA. Glysson, G. D. G. (1993). "U.S. Geological Survey bedload sampling policy." Hydraulic Engineering '93, Part 1, 701-706. Habersack, H. M. (2001). "Radio-tracking gravel particles in a large braided river in New Zealand: a field test of the stochastic theory of bed load transport proposed by Einstein." Hydrol. Process., 15, 377-391. Habersack, H. M., and Laronne, J. B. (2002). "Evaluation and improvement of bed load discharge formulas based on Helley-Smith sampling in an alpine gravel bed river." J. Hydraulic Eng. (ASCE), 128(5), 484-499. Hansen, D.S. (1986). "Oceanic incoherent Doppler sonar spectral analysis by conventional and finite-parameter modeling methods." IEEEJ. Oceanic Eng., OE-11(1), 26-41. Harris, B. (1966). Theory of Probability, Addison-Wesley, Reading, Massachusetts. Haschenburger, J. K., and Church, M. (1998). "Bed material transport estimated from the virtual velocity of sediment." Earth Surf. Proc. Landf, 23, 791-808. Haschenburger, J. K. (1999). "A probability model of scour and fill depths in gravel-bed channels." Water Resour. Res., 35, 2857-2869. 199 Hassan, M. A., and Church, M. (2001). "Sensitivity of bed load transport in Harris Creek: Seasonal and spatial variation over a cobble-gravel bar." Water Resour. Res., 37(3), 813-825. Hassan, M. A., Church, M., and Ashworth, P. J. (1992). "Virtual rate and mean distance of travel of individual clasts in gravel-bed channels." Earth Surf. Proc. Landf, 17, 617-627. Hassan, M. A., Church, M., and Schick, A.P. (1991). "Distance of movement of coarse particles in gravel bed streams." Water Resources Res., 27(4), 503-511. Hay, A. E. (1991). "Sound scattering from a particle-laden turbulent jet." J. Acoust. Soc. Am., 90(4 (Pt. 1)), 2055-2074. He, C , and Hay, A. E. (1993). "Broadband measurements of the acoustic backscatter cross section of sand particles in suspension."/. Acoust. Soc. Am., 94(4), 2247-2254. Helley, E. J., and Smith, W. (1971). "Development and calibration of a pressure-difference bedload sampler." , US Geol. Survey Open File Report, Menlo Park, CA. Houck, C , Joines, J., and Kay, M. (1995). "A genetic algorithm for function optimization: a Matlab application." , North Carolina State University College of Engineering, Dept. of Industrial Engineering, Technical Report 95-09. Hu, C , and Hui, Y. (1996a). "Bed-load transport. I: Mechanical characteristics."/. Hydraulic Eng. (ASCE), 122(5), 245-254. Hu, C , and Hui, Y. (1996b). "Bed-load transport. II: Stochastic characteristics."/ Hydraulic Eng. (ASCE), 122(5), 255-261. Hubbell, D. W., and Sayre, W. W. (1964). "Sand transport studies with radioactive tracers." / Hydraulics Div., ASCE, 90(HY3), 39-68. Hui, Y., and Hu, C. (1991). "Saltation characteristics of particle motions in water." Shuili Xuebao (J. Hydraulic Eng.) 12, 59-64 (in Chinese). Inggs, M. R., and Lord, R. T. "Interpolating satellite derived wind field data using ordinary kriging, with application to the Nadir Gap." International Geoscience.and Remote Sensing Symposium, 141-143. Isaaks, E. H., and Srivastava, R. M. (1989). An Introduction to Applied Geostatistics, Oxford University Press, New York. Jackson, W. L., and Beschta, R. L. (1982). "A model of two-phase bedload transport in an Oregon Coast Range stream." Earth Surf. Proc. Landf, 7, 517-527. 200 Jorgensen, B. (1997). The Theory of Dispersion Models, Chapman & Hall, London. Joyce, T. M. (1989). "On in situ "calibration" of shipboard ADCPs." J. Atmos. Oceanic Technol., 6, 169-172. Kellerhals, R. (1967). "Stable channels with gravel-paved beds." J. Waterways Harbors Div., Am. Soc. Civil. Eng. ,97,1165-1180. Keshavarzy, A., and Ball, J. E. (1999). "An application of image processing in the study of sediment motion."/. Hydraulic Res., 37(4), 559-576. Klingeman, P. C , and Milhous, R. T. (1971). "Oak Creek vortex bed-load sampler." EOS, Transactions, American Geophysical Union, 52(5), 434. Kohanowich, K. M., Stanton, T. P., and Thornton, E. B. "Acoustic sediment flux measurements from Duck '94." Coastal Dynamics '95, Proceedings of the International Conference on Coastal Research in Terms of Large Scale Experiments, Gdansk, Poland, September 4-8, 1995, 739-748. Kolb, M. (1995). "Experiences with vessel borne ADCPs in shallow waters." Proc. IEEE Fifth Working Conference on Current Measurement, IEEE, 79-82. Kostaschuk, R., and Villard, P. (1996). "Flow and sediment transport over large subaqueous dunes: Fraser River, Canada." Sediment., 43, 849-863. Kostaschuk, R. A., and Ilersich, S. A. (1995). "Dune geometry and sediment transport: Fraser River, British Columbia." River Geomorphology, E. J. Hickin, ed., John Wiley & Sons, Chichester, 19-36. Kostaschuk, R. A., Church, M. A., and Luternauer, J. L. (1989). "Bedforms, bed material, and bedload transport in a salt-wedge estuary: Fraser River, British Columbia." Can. J. Earth Sci., 26, 1440-1452. Kuhnle, R. A., and Southard, J. B. (1988). "Bed load transport fluctuations in a gravel bed laboratory channel." Water Resour. Res., 24(2), 247-260. Lane, S. N., Richards, K. S., and Chandler, J. H. (1995). "Morphological estimation of the time-integrated bed-load transport rate." Water Resour. Res., 31(3), 761-772. Laronne, J. B., and Duncan, M. J. (1992). "Bedload transport paths and gravel bar formation." Dynamics of Gravel-bed Rivers, P. Billi, R. D. Hey, C. R. Thorne, and P. Tacconi, eds., John Wiley & Sons Ltd, Chichester, England, 177-200. Larsen, R. J., and Marx, M. L. (1986). An introduction to mathematical statistics and its applications, Prentice-Hall, Englewood Cliffs, NJ. 201 Lee, H.-Y., Chen, Y.-H., You, J.-Y., and Lin, Y.-T. (2000). "Investigations of continuous bed load saltating process." J. Hydraulic Eng. (ASCE), 126(9), 691-700. Lee, H.-Y., and Hsu, I.-S. (1994). "Investigation of saltating particle motions."/. Hydraulic Eng. (ASCE), 120(7), 831-845. Lemmin, U., and Lhermitte, R. (1999). "Discussion of: ADV measurements of turbulence: can we improve their interpretation? by V.I. Nikora and D.G. Goring." J. Hydraulic Eng. (ASCE), 125, 987-988. Lowe, R. L., Inman, D. L., and Drake, T. G. (1991). "A bedload sensor for wave and current regimes." EOS, Trans. Am. Geophys. Union, 72(44), 230. Matheron, G. (1971). "La theorie des variables regionalisees et ses applications." , Cahiers du Centre de Morphologie Mathematique, Fontainebleau, France. McLean, D. G., and Tassone, B. (1987). "Discussion of'Bed load sampling and analysis'." Sediment Transport in Gravel-bed Rivers, C. R. Thorne, J. Bathurst, C, and R. Hey, D., eds., John Wiley & Sons Ltd., New York, 109-113. McLean, D. G., and Church, M. (1999). "Sediment transport along lower Fraser River: 2. Estimates based on the long-term gravel budget." Water Resour. Res., 35(8), 2549-2559. McLean, D. G., Church, M., and Tassone, B. (1999). "Sediment transport along lower Fraser River: 1. Measurements and hydraulic computations." Water Resour. Res., 35(8), 2533-2548. Medwin, H., and Clay, C. S. (1998). Fundamentals of Acoustical Oceanography, Academic Press, San Diego, CA. Montgomery, D. C , and Peck, E. A. (1982). Introduction to Linear Regression Analysis, John Wiley & Sons, New York. Mosley, M. P. (1981). "The influence of organic debris on channel morphology and bedload transport in a New Zealand forest stream." Earth Surf. Proc. Landf, 6, 571-579. Mulhoffer, L. (1933). "Untersuchungen uber der Schwebstoff und Geschiebefurung des Inn nachst Kirchbichl, Tirol [Investigations into suspended load and bedload of the River Inn, near Kirchbichl, Tirol]". Wasserwirtschaft, Heft 1-6. Neill, C. R. (1987). "Sediment balance considerations linking long-term transport and channel processes." Sediment Transport in Gravel Bed Rivers, C. R. Thorne, J. C. Bathurst, and R. D. Hey, eds., John Wiley, New York, 225-240. 202 Nelson, J. M., Shreve, R. L., McLean, S. R., and Drake, T. G. (1995). "Role of near-bed turbulence structure in bed load transport and bed form mechanics." Water Resour. Res., 31(8), 2071-2086. Nino, Y., Garcia, M. (1994). "Gravel saltation 2. modeling." Water Resour. Res., 30(6), 1915-1924. Nino, Y., Garcia, M., and Ayala, L. (1994). "Gravel saltation 1. experiments." Water Resour. Res., 30(6), 1907-1914. Novak, P. (1957) "Bed load meters - development of a new type and determination of their efficiency with the aid of scale models." Trans. Int. Assoc. Hydraul. Res. 7th Gen. Meeting, Lisbon, A9-1: A9-11. Parker, G., and Klingeman, P. C. (1982). "On why gravel bed streams are paved." Water Resour. Res., 18(5), 1409-1423. Parker, G., Klingeman, P. C , and McLean, D. G. (1982). "Bedload and size distribution in paved gravel-bed streams." J. Hydraulic Engrg., 108(HY4), 544-571. Pinkel, R. (1980). "Acoustic Doppler techniques." Air-Sea Interaction Instruments and Methods, F. Dobson, L. Hasse, and R. Davis, eds., Plenum Press, New York, 171-199. Polonichko, V. (2000). SonTek, San Diego, CA. Powell, D. M., and Ashworth, P. J. (1995). "Spatial pattern of flow competence and bed load transport in a divided gravel bed river." Water Resour. Res., 31(3), 741-752. Reid, I., Layman, J. T., and Frostick, L. E. (1980). "The continuous measurement of bedload discharge."/. Hydraulic Res., 18, 243-249. Rennie, C. D., and Millar, R. G. (2000). "Spatial variability of stream bed scour and fill: a comparison of scour depth in chum salmon (Oncorhynchus keta) redds and adjacent bed." Can. J. Fish. Aquatic Sci., 57(5), 928-938. Rennie, CD., and Millar, R.G. (2001). "Measurements of gravel bedload transport velocity using an acoustic Doppler profiler." Gravel-bed Rivers 2000 CD-ROM, T. Nolan, and C. Thorne, eds., Special Publication of the New Zealand Hydrological Society. Rennie, C. D., Millar, R. G., and Church, M. A. (2002). "Measurement of bed load velocity using an acoustic Doppler current profiler." /. Hydraulic Eng. (ASCE), 128(5), 473-483. Rouhani, S., and Myers, D.E. (1990). "Problems in space-time kriging of geohydrological data." Math. Geol, 22(5), 611-623. 203 Saffman, P.G. (1965). "The lift on a small sphere in a slow shear flow." J. Fluid Meek, 22(2), 385-400. Schmidt, K.-H., and Ergenzinger, P. (1992). "Bedload entrainment, travel lengths, step lengths, rest periods - studied with passive (iron, magnetic) and active (radio) tracer techniques." Earth Surf. Proc. Landf, 17, 147-165. Sear, D. A. (1996). "Sediment transport processes in pool-riffle sequences." Earth Surf. Proc. Landf ,21,241-262. Shen, C , and Lemmin, U. (1999). "Application of an acoustic particle flux profiler in particle-laden open-channel flow." J. Hydraulic Res., 37(3), 407-419. Sidle, R. C. (1988). "Bed load transport regime of a small forest stream." Water Resour. Res., 24(2), 207-218. Simons, D. B., Richardson, E. V., and Nordin, C. F. J. (1965). "Sediment transport in alluvial channels: Bedload equation for ripples and dunes." United States Geological Society Professional Paper, 462H. SonTek. (1998). SonTek ADP™ Acoustic Doppler Profiler Technical Documentation, SonTek, San Diego, CA. Stanton, T. P., and Thornton, E. B. "Sediment fluxes above a mobile sandy bed in the nearshore." Coastal Sediments '99, Volume One, Proceedings of the 4th International Symposium on Coastal Engineering and Science of Coastal Sediment Processes, Hauppauge, Long Island, New York, June 21-23,1999, 241-252. Stelczer, K. (1981). Bed-Load Transport Theory and Practice, Water Resources Publications, Littleton, Colorado. Steveston Marine & Hardware. (2001). Tide Chart Book 2001, R.D. Penhall Ltd., Vancouver. Sutton, D. W., and Jaffe, J., S. (1992). "Acoustic bedload velocity estimates using a broadband pulse-pulse time correlation technique." J. Acoust. Soc. Am., 92(3), 1692-1698. Taylor, J. A., Vincent, C. E., Thorne, P. D., Hardcastle, P. J., Humphrey, V. F., Zhang, J. D., Schaafsma, A., Dohmen-Janseen, C. M., and Perennes, M. (1998). "Three-dimensional sediment transport measurements by acoustics (TRIDISMA)." Oceans '98 Conference Proceedings, Vol. 2, 1108-1114. Ten Brinke, W. B. M., Wilbers, A. W. E., and Wesseling, C. (1999). "Dune growth, decay and migration rates during a large-magnitude flood at a sand and mixed sand-204 gravel bed in the Dutch Rhine river system." Spec. Pubis. Int. Ass. Sediment., 28, 15-32. Theriault, K. B. (1986a). "Incoherent multibeam Doppler current profiler performance: Part I - estimate variance." IEEEJ. Oceanic Eng., OE-11(1), 7-15. Theriault, K. B. (1986b). "Incoherent multibeam Doppler current profiler performance: Part II - spatial response." IEEE J. Oceanic Eng., OE-11(1), 16-25. Thomson, R. E. (1981). Oceanography of the British Columbia Coast, Department of Fisheries and Oceans, Ottawa. Thorne, P. D., Waters, K. R., and Brudner, T. J. (1995). "Acoustic measurements of scattering by objects of irregular shape."/. Acoust. Soc. Am., 97(1), 242-251. van Es, A. J., and Kok, A. R. (1998). "Simple kernel estimators for certain nonparametric deconvolution problems." Statistics and Probability Letters, 39, 151-160. van Rijn, L. C. (1982). "Equivalent roughness of alluvial bed." J. Hydraulic Eng., 108(10), 1431-1456. van Rijn, L. C. (1984a). "Sediment transport. Part I: bed load transport." J. Hydraulic Eng., 110(10), 1431-1456. van Rijn, L. C. (1984b). "Sediment transport. Part III: bed forms and alluvial roughness." J. Hydraulic Eng., 110(12), 1733-1754. Villard, P., and Kostaschuk, R. (1998). "The relation between shear velocity and suspended sediment concentration over dunes: Fraser Estuary, Canada." Marine Geol. 148,71-81. Villard, P., Church, M., and Kostaschuk, R. (in press). "Estimating bed load in sand-bed channels using bottom tracking from an Acoustic Doppler Profiler." Fluvial Sedimentology VII. Wand, M. P. (1998). "Finite sample performance of deconvolving density estimators." Statistics and Probability Letters, 37, 131-139. Weisstein, E. W. (1999). CRC Concise Encyclopedia of Mathematics, Chapman & Hall/CRC, Boca Raton. Wilcock, P. R., and McArdell, B. W. (1997). "Partial transport of a sand/gravel sediment." Water Resour. Res., 33(1), 235-245. 205 Williams, J. J. (1996). "Turbulent flow in rivers." Advances in Fluvial Dynamics and Stratigraphy, P. A. Carling and M. R. Dawson, eds., Jon Wiley & Sons, Chichester, 1-32. Yalin, S. M. (1977). Mechanics of Sediment Transport, Pergamon Press, Oxford. Yalin, S. M. (1972). Mechanics of Sediment Transport, Pergamon Press, Oxford. Yang, C. T., and Sayre, W. W. (1971). "Stochastic model for sand dispersion." J. Hydraulics Division, ASCE, 97(HY2), 265-288. Young, D. S. (1987). "Random vectors and spatial analysis by geostatistics for geotechnical applications." Math. Geol., 19(6), 467-479. Yusuf, F. (2001). "Application of a two dimensional hydrodynamic model to the Fraser River gravel reach" unpublished M.A.Sc. thesis, Dept. Civil Engineering, University of British Columbia, Vancouver. Zar, J. Ff. (1996). Biostatistical Analysis, 3rd Edition, Prentice-Hall, Englewood Cliffs. Zedel, L., and Hay, A. E. (1999). "A coherent Doppler profiler for high-resolution particle velocimetry in the ocean: Laboratory measurements for turbulence and particle flux." J. Atmos. Oceanic Technol, 16, 1102-1117. 206 Appendix A: Rotating particles In Section 1.3.7 it is noted that observations of saltating particles in field conditions do not show significant rotation during the trajectory (Drake et al. 1988). Based on these observations, the contribution of rotation to the surface velocity of particles has been ignored in this thesis. Nino and Garcia (1994) reported that particle rotation rates in their flume experiments were highly dependent on specifics of the particle collision with the bed. Hui and Hu (1991, reported in Lee and Hsu 1994) found a spin rate of about 40 rev/s in their flume tests. The rotation rate varied with particle size and shape. For a 1 mm sand grain this would equal a surface velocity of 0.13 m/s. Clearly, better understanding of bedload particle rotation is required for assessment of acoustic bedload velocity measurements. One would expect rotation of a freely rotating particle about the cross-stream axis due to the mean velocity gradient. Higher velocities at the top of the particle should force it to turn over. Presumably, the lack of rotation observed by Drake et al. (1988) was due to random forces imposed by the turbulent fluid, which prevented consistent rotation about the cross-stream axis of the particle. This appendix offers a derivation of the possible additional surface velocity due to rotation based on the forces imposed by the mean velocity gradient. For a freely rotating particle, the angular rotation (Q) is: (A.1) 207 where k is the velocity gradient (Saffana 1965). The velocity gradient at the location of the saltating particle will be: k=UfD-U(f-l)D D V ' where D is the particle diameter, and fD is the location above the bed of the top of the particle, in terms of particle diameters. A typical particle saltation height is 3D, although observed saltation heights range from 1.2D to 9.2D (see Section 1.3.7). For a semi-log linear velocity profile (Equation 1.20) u = — ln(h) + — ln\ '30^ K \ k s J Substitution of (1.20) into (A.2) yields: — In f f A KD The surface velocity of the particle due to rotation is: 2 Substitution of (A.3) and (A.l) into (A.4) yields: v„ = — In r 4K (A.3) (A.4) (A.5) Assuming u* equal to 0.1 m/s and/equal to 3, vr is 0.025 m/s. For u* equal to 0.2 m/s, vr is 0.050 m/s. The translational velocity (v,) of the particle can also be estimated using (1.20). This assumes that the particle velocity equals the water velocity. In reality, the particle 208 velocity will be less than the water velocity by a slip velocity that is approximately equal to the particle fall velocity (co) (Bagnold 1973). Also, bedload tends to occur during sweeps with high excess instantaneous streamwise (and downward) velocities. Particle velocities are thus not determined by the mean flow velocity. Hu and Hui (1996a) found that the slip velocity (estimated as mean flow velocity at location of particle minus mean particle velocity) was a function of dimensionless shear stress, and varied from about 0.3co to 2.2co. Alternatively, the saltation velocity has been found to equal Cu*. Typical values of C range from 10 for sands (Lee et al. 2000) to 3 for gravels (Nino et al. 1994). Using this estimate of vt, the ratio of particle rotational velocity to translational velocity is: ( f A In U-l. (A.6) v, 4 C K For typical values of/equal to 3 and C equal to 10, the ratio is 0.025, thus the rotational velocity is a small percentage of the translational velocity. However, the largest possible value of this ratio, given the parameter ranges described above (f= 1.2 and C = 3), is 0.36. Further, the rotation rate may depend on the interaction with the bed more than the velocity gradient (Nino and Garcia 1994). It is apparent that assessment of the importance of particle rotation on the measured Doppler velocity merits further consideration. 209 Appendix B: ADP Instrument Settings Please refer to SonTek (1998) for details on instrument operation. This appendix lists the instrument settings used in this thesis. The 1.5 MHz ADP used in this thesis was operated using the firmware CPU Software Version 6.4 and DSP Software Version 4.0, except laboratory trials 1-6, 9-10, and 13-16, which used CPU Software Version 6.1. The DOS data collection system was used (adpsa.exe). When DGPS was integrated with the ADP, the command ADPSA -p was used to start pinging. The DGPS was integrated using the following one-line configuration file (named gpsport.def) which tells the ADP which DGPS script to use and governs the serial connection: $GPGGA 2 9600 N 8 1 All of the instrument settings, except bottom track pulse length, were simply adjusted using the user interface of adpsa.exe. The bottom track pulse length could only be adjusted using the SonTek terminal emulator (sonterm.exe). To enter sonterm type at the DOS prompt: >sonterm The bottom track configuration is seven numbers, one for each of pulselength, blanklength, Nsamp, Ntail, TanW, CorrThresh, MinScore. Only the pulselength was changed from the defaults. To view the current settings, type: >btconf To change the settings, type the following line, with a number substituted for each variable: >btconf set pulselength blanklength Nsamp Ntail TanW CorrThresh MinScore To return to the default settings, type: >btconf init Note that the bottom track configuration settings are not reset to the defaults when the instrument is shut down and powered up again. Table B.l ADP instrument settings Parameter Fraser River Laboratory Norrish Fraser River Fraser River gravel bed (Chapter 2) Creek gravel bed sand bed Sea stationary (Chapters 1 and 3) stationary (Chapter 3) mobile boat (Chapter 4) Reach mobile boat (Chapter 4) Water 0 0 0 0 0 salinity (ppt) Temperature measured measured measured measured measured mode Averaging interval (s) 5 1 1 1 1 Number of 40 a 1 10 40 40 depth cells Depth cell 0.25 4 0.25 0.25 0.25 size (m) Blanking 0.40 20 0.40 0.40 0.40 distance (m) Coordinate ENU XYZ XYZ ENU ENU system Sensor 0.36 0 0.35 0.64 0.46 mounting depth (m) Magnetic 19.6 0 19.6 19.6 19.6 variation (deg.) Bottom track YES YES YES YES YES Bottom track 0.20 0.60 0.20 0.60 0.60 pulse length (m) a Variable, depending on local depth 211 Appendix C: Photographs Figure C. 1: Fraser River at Agassiz-Rosedale Bridge, facing upstream. Former WSC gauge site 08MF035. Width is 510 m. Described as a wandering gravel-bed river (D50,surf = 42 mm), with a slope of 4.8x10-4 (McLean et al. 1999). If bridge spans are numbered from left to right (facing upstream), the thalweg flows through span 5, and a shallow bar occurs at span 2. All 2000 freshet bedload sampling stations were downstream of the bridge at spans 3 and 4 (Chapters 1,3, and 4). Figure C.2: a) ADP, b) 24 size VuV bedload sampler, and c) Helley-Smith bedload sampler. Figure C.3: Boat for Fraser River gravel-bed data collection (Chapters 1,3, and 4). Figure C.4: Laboratory set-up for synthetic bedload. Note: an improved artifical bed was implemented for Trials 7, 8, 11, 12, 17, and 18. (Chapter 2). 213 Figure C.5: Norrish River at Hawkins Pickle Road Bridge, facing downstream from bridge, February 22,2002. WSC gauge site 08MH058 (Chapter 3). River width at bridge is 28 m. Figure C.6: Norrish Creek on Feb 22,2002. (Chapter 3). ADP securely mounted to bridge deck Figure C.7: Boat for Fraser River sand bed Sea Reach data collection (Chapter 4). Figure C.8: Fraser River sand bed Sea Reach, facing upstream (Chapter 4). 215 Appendix D: Error model fits Plots similar to Figure 3.1 are provided for Fraser River and Norrish Creek error model fits for estimated compound Poisson gamma, optimized compound Poisson gamma, beam gamma, and total gamma models. Probability density Probability density CO Probability density 00 I CD CO 0) 3 H —I Q) U & o" 13 Probability density rO 05 CO O O O O O O ' i I I 1 I I 1_ i i DJ 3 OJ cr co o_ c ro Q. c? 3 3 O CD II O OJ II n II p o OJ OJ NJ Oi 3 NJ 3 TJ w OJ co~ P, "co* < o' ro ro Q. w ro co 5' < ro CT o T3 CD —+i OJ OJ T3 p. 3 arti icie< OJ o icie< ro ro «; CD QJ, II 3 3 CD CD OJ 0) Z) 13 £ 8 o o 0 -vl CD O 3 .. 8 w o oo co 3 r* co Q . CD < o —H CT Z3 g CD O CD-CD. TJ 0) —t 0) 3 CD CD 55 Probability density Probability density CO I cu CD 3 0) ZJ 3* 3 o" Probability density 03 CD cu 3 CO CD 3 co" P o >-l P cn CD < a CL a 00 CD CD o_ co T3 CD CD Q. 3 CD P 3 Q. 3 a o 3 3 55* II Probability density 3 _n«? CD ~ - o P P cr co o_ c CD Q. c? —I CD 3 O CD II O b Ol M •a II "ll O OJ II o 3 Sr M o> P " * g 3 -o CD CO P 3 3. < CO 3 3 CD CD P P 3 3 CO < — Q. CD < P % CD « i i . O co • CD O II M O 00 ^ 3 o p •a p 3. cV CD < CD O O 3 s-_ CD £T < CD • P P O a 3 r; §• * » 8 I . — CO o Q_ (D TJ 0) 3 CD I—K CD |3 Probability density Probability density 9*. co" r-t-3. cr C o »—^ o • f l i-l P Co Ct> i-t < CD OP c 1-1 CT> O CO T3 m o M O W O O o cv T3 II o o CD Cv e I-I CD CV P p 3 <-» cv p -7 co 3 2. w El cn DJ < CD —i DJ CQ CD CT CD C L Cn TJ CD CD C L 3 55" o o to -j a " —l . 2.3 i Q J II p isfor load )99 m/s med signal + Probability density <-C O I CO CD 0) 3 0) ZJ 3 CD o" OJ i s OJ =2: CQ 3. CD 3 CT CT CD CD C L OJ to 3 TJ 3 CD V) CD "a CL OJ g P cn" a> o ho o b o b 00 p ^ 2, 2". c 3 OJ 4^ o 3 co" CD OJ O 3 " CT CD OJ 3 cr CD CL O OJ CL CO CQ" 3 DJ Probability density i i OJ 3 OJ CT to O c I— CD Q. 3 O CD II O b cn ro P Q CD II 9: O O Li DJ jv C L O < CD O O 3 co" 3 II o b oo 3 co" "5T CL CD < 3 3 CD CD OJ OJ 3 3 3 II 2. O CO • ro o II NJ o o° 3 co" CL ro < a o CO TJ CT DJ CD OJ *< 3 OJ o < CO' ro ro DJ ^ CQ ro Q o o_ cp_ TJ 0) -1 3 CD CD 3 Probability density co TJ p ' II Probability density CO I CD CD —I —t ZJ CD i—i-O* ZJ Probability density co CD 3 3 O to" CD 3 i l 3 03 tT C/> O c 0 Q . a5 3 8 n o CD -N 00 Q Q D" < cr CD rf 3 CD § II 3 " l l P cr " o CD O 9. .j. oo O O 3 0) < In 3 3 CD CD £0 CO 3 3 cp_ C/> o o O 0 0 3 5r cn Q . CD < CL CD < cr ro 3 3 g O a. cp_ TJ 0) —i 0) 3 CD CD o 3 g-3 ' Probability density Probability density KJ KJ O Probability density Probability density w 9 3 CO" o Probability density <-CO I CD CD 0) 3 0) ZJ W cT CO o" ZJ Probability density ro CD CL co ro CL 3 ro co 3 CD a O 3 CO" Probability density CD cu 3 OJ CT co o c ro a. I 3 O CD II O b NJ 13 II NJ 00 NJ 3 CO" l o CL CD < cr T3 OJ g. o ro < ro o o II II J P NJ , GO OJ < _. tt a> co OJ 5 . ^ P, w 0) < S- W «> 8 R* 3 2. Q-CO 3 ~ » cr -o ro OJ OJ 3-. 3 S" OJ QJ (D 3 5 „ ro 3 CD - ~ - 3 o 2, O CO 3 & 3 3 ro ro OJ OJ 3 3 3 II 2. o ro f° II CJ1 ro < CT CD O Q . cp_ Tl 0) - i 0) 3 (D crT —i c/> Probability density Probability density CO I CD CD 0) 3 — i 0) 3 cT —i 3 CO r-t-o' 3 Probability density <-3 Z! cr CD CD CO 3 CO cr co o c CD g_ 3 II o b NJ cn o. o 0) C l < CD O O Q Q e cr 3 3 CD CD DO 03 3 3 CO & -5T O. ' < | a o co - a cr 00 CD Cd" CO *< 3 CO o < co" CD CD CO w CO CD C l CD f 0 ii cn CO o CD TJ 0) —i 0) 3 CD CD CO to to to Probability density <-CO I CO CD Q) 3 0) in & 3 <—t-O* Z5 Probability density 3 CD OJ 3 OJ CT CO o c f— CD g. j? 3 II O b NJ NJ Q q cr CD § OJ 5> II II o a o o)

3 CD DJ 224 Appendix E: Statistical comparisons of spatially block averaged bedload velocity The spatial block averages for bedload velocity were performed in a 5x5 grid. The blocks are numbered 1 to 25 in ascending columns, starting at the south-west block (ie. the north-west block is number 5 and the north-east block is number 25). The following tables describe the statistical comparisons between blocks for bedload velocity resolved in the downstream direction. Gravel-bed reach of Fraser River For the gravel bed reach, only block #6 was significantly different from the other blocks. blockprimbs Block N Mean Std. Std. Error 95% Confidence Minimum Maximum Deviation Interval for Mean Lower Upper Bound Bound 2 20 .210 .296 .066 .071 .348 -.032 1.368 3 39 .065 .143 .023 .018 .111 -.301 .509 4 17 .169 .215 .052 .058 .279 -.006 .774 6 14 .378 .521 .139 .077 .679 -.060 1.559 7 10 .188 .333 .105 -.050 .426 -.208 .907 8 19 .087 .146 .034 .017 .158 -.125 .318 9 26 .074 .286 .056 -.041 .189 -.229 1.199 12 16 .224 .337 .084 .045 .404 -.282 1.195 13 37 .132 .215 .035 .061 .204 -.150 .940 14 12 .001 .113 .033 -.071 .072 -.127 .197 18 36 .068 .252 .042 -.017 .154 -.384 .759 19 27 .055 .189 .036 -.019 .130 -.556 .493 20 13 .030 .167 .046 -.071 .131 -.128 .469 24 23 .057 .318 .066 -.081 .194 -.488 .862 25 86 .041 .240 .026 -.011 .092 -.477 .700 Total 395 .096 .260 .013 .070 .122 -.556 1.559 225 Sand-bed Sea Reach of Fraser River Descriptives prinibs, m/s 95% Confidence Interval for Mean N Mean Std. Deviation Std. Error Lower Bound Upper Bound Minimum Maximum 1 84 .054 .120 | .013 .028 .080 | -.348 .381 ) 2 86 .061 .100 | .011 .040 .082 | -.206 .368 3 79 .096 .142 | .016 .064: .128 | -.385 .399 4 58 .226 .161 I .021 .184 .269 j -.111 .594 6 84 .466 .203 | .022 .422 .510 | -.110 1.181 7 69 .366 .268 I .032 .302 .431 | -.206 .943 8 81 .141 .247 | .027 .086 .196 | -.298 .937 9 53 .417 .241 I 0 3 3 .350 .483 | -.124 1.001 12 57 .484 .265 I .035 .414 .554 | -.150 1.130 13 124 .333 .293 I .026 .281 .385 | -293 1.178 14 139 .284 .303 I .026 .233 .335 j -.311 1.306 15 91 \ .357 .286 | .030 .297 .416 | -.129 .969 19 101 .290 .316 | .031 .228 .353 | -.336 1.561 20 264 .238 .295 | .018 .202 .273 | -.227 1.304 25 129 .231 .274 I .024 .183 .279 | -.289 1.254 Total 1499 .261 .282 | .007 .247 .275 | -.385 1.561 Test of Homogeneity of Variances primbs, m/s | Levene Statistic j df1 j df2 Sig. j | 12.7241 14 1484! .000! ANOVA primbs, m/s Sum of Squares j df Mean Square F Sig. Between Groups 20.727) 14 1.480 22.316 .000 Within Groups 98.449 1484 6.634E-02 Total 119.176 1498 226 Post Hoc Tests Multiple Comparisons Dependent Variable: primbs, m/s Mean Difference (l-J) Std. Error Sig. 95% Confidence Interval (1) BLOCK (J) BLOCK Lower Bound Upper Bound Tukey u c n 1 2 -6.55E-03 3.95E-02 1.000 -.14 .13. 3 -4.16E-02 4.04E-02 1.000 -.18 9.53E-02 4 -.17f) 4.40E-02 .008 -.32 -2.25E-02 6 -.41 (*) 3.97E-02 .000 -.55 -.28 7 -.310 4.18E-02 .000 -.45 -.17 8 -8.65E-02 4.01 E-02 .697 -.22 4.95E-02 9 -.360 4.52E-02 .000 -.52 -.21 12 -.43(*) 4.42E-02 .000 -.58 -.28! 13 -•280 3.64E-02 .000 -.40 -.16 14 -•23(*) 3.56E-02 .000 -.35 -.11 15 -.30(*) 3.90E-02 .000 -.43 -.17 19 -.24(*) 3.80E-02 .000 -.37 -.11 20 -.18(*) 3.23E-02 .000 -.29 -7.39E-02 25 -.18(*) 3.61 E-02 .000 -.30 -5.43E-02 2 1 6.55E-03 3.95E-02 1.000 -.13 .14 3 -3.50E-02 4.01 E-02 1.000 -.17 .10 4 -.17(*) 4.38E-02 .013 -.31 -1.67E-02 6 -.41(*) 3.95E-02 .000 -.54 -.27 7 -.31 (*) 4.16E-02 .000 -.45 -.16 8 -7.99E-02 3.99E-02 .796 -.22 5.53E-02 9 -.36(*) 4.50E-02 .000 -.51 -.20 12 -.42(*) 4.40E-02 .000 -.57 -.27 13 -.27(*) 3.61 E-02 .000 -.39 -.15 14 -.22(*) 3.53E-02 .000 -.34 -.10 15 -.30(*) 3.87E-02 .000 -.43 -.16 227 19 -23C) 3.78E-02 .000 -.36 -.10 20 -•1&T) 3.20E-02 .000 -.29 -6.83E-02 25 -.17C) 3.59E-02 .000 -.29 -4.86E-02 1 4.16E-02 4.04E-02 1.000 -9.53E-02 .18 2 3.50E-02 4.01 E-02 1.000 -.10 .17 4 -.13 4.45E-02 .186 -.28 2.09E-02 6 -.370 4.04E-02 .000 -.51 -.23 7 -.27(*) 4.24E-02 .000 -.41 -.13 8 -4.49E-02 4.07E-02 .999 -.18 9.32E-02 9 --32H 4.57E-02 .000 -.48 -.17 o 12 --39C) 4.48E-02 .000 -.54 -.24 13 -•24C) 3.71 E-02 .000 -.36 -.11 14 -.19(*) 3.63E-02 .000 -.31 -6.47E-02 15 -.26(*) 3.96E-02 .000 -.40 -.13 19 -.19(*) 3.87E-02 .000 -.33 -6.33E-02 20 -.14(*) 3.30E-02 .002 -.25 -2.97E-02 25 --14D 3.68E-02 .019 -.26 -1.04E-02 1 .17(*) 4.40E-02 .008 2.25E-02 .32 2 •17(*) 4.38E-02 .013 1.67E-02 .31 3 .13 4.45E-02 .186 -2.09E-02 .28 6 -.24(*) 4.40E-02 .000 -.39 -9.10E-02 7 -.14 4.59E-02 .134 -.30 1.56E-02 8 8.52E-02 4.43E-02 .842 -6.51 E-02 .24 A 9 -.19C) 4.89E-02 .009 -.36 -2.46E-02 12 -.26(*) 4.80E-02 .000 -.42 -9.50E-02 13 -.11 4.10E-02 .364 -.25 3.21 E-02 14 -5.77E-02 4.03E-02 .984 -.19 7.88E-02 15 -.13 4.33E-02 .146 -.28 1.61 E-02 19 -6.44E-02 4.24E-02 .974 -.21 7.95E-02 20 -1.16E-02 3.74E-02 1.000 -.14 .12 25 -5.07E-03 4.07E-02 1.000 -.14 .13 6 1 .41 (*) 3.97E-02 .000 .28 .55 228 2 .41 n 3.95E-02 .000 .27 .54 3 •370 4.04E-02 .000 .23 .51 4 •24(*) 4.40E-02 .000 9.10E-02 .39 7 .10 4.18E-02 .521 -4.18E-02 .24 8 33(*) 4.01 E-02 .000 .19 .46 9 4.96E-02 4.52E-02 .999 -.10 .20 12 -1.78E-02 4.42E-02 1.000 -.17 .13 13 •13(*) 3.64E-02 .020 9.83E-03 .26 14 •18(*) 3.56E-02 .000 6.17E-02 .30 15 .11 3.90E-02 .242 -2.27E-02 .24 19 •18C) 3.80E-02 .000 4.68E-02 .30 20 •230 3.23E-02 .000 .12 .34 25 •24(*) 3.61 E-02 .000 .11 .36 1 .31 (*) 4.18E-02 .000 .17 .45 2 .31 n 4.16E-02 .000 .16 .45 3 •27(*) 4.24E-02 .000 .13 .41 4 .14 4.59E-02 .134 -1.56E-02 .30 6 -.10 4.18E-02 .521 -.24 4.18E-02 8 •23(*) 4.22E-02 .000 8.21 E-02 .37 9 -5.05E-02 4.70E-02 .999 -.21 .11 1 12 -.12 4.61 E-02 .400 -.27 3.85E-02 13 3.32 E-02 3.87E-02 1.000 -9.80E-02 .16 14 8.23E-02 3.79E-02 .686 -4.63E-02 .21 15 9.36E-03 4.11 E-02 1.000 -.13 .15 19 7.57E-02 4.02E-02 .863 -6.07E-02 .21 20 •13(*) 3.48 E-02 .018 1.03E-02 .25 25 .13(*) 3.84E-02 .034 4.70E-03 .27 8 1 8.65E-02 4.01 E-02 .697 -4.95E-02 .22 2 7.99E-02 3.99E-02 .796 -5.53E-02 .22 3 4.49E-02 4.07E-02 .999 -9.32E-02 .18 4 -8.52E-02 4.43E-02 .842 -.24 6.51 E-02 6 -•330 4.01 E-02 .000 -.46 -.19 229 7 -.23H 4.22E-02 I 0 0 0 -.37 -8.21 E-02 9 -.28(*) 4.55E-02 | .000 -.43 -.12 12 -•34(*) 4.45E-02 I 0 0 0 -.49 -.19 13 -.19(*) 3.68E-02 | .000 -.32 -6.72E-02 14 -•14(*) 3.60E-02 I .006 -.27 -2.08E-02 15 -•22C) 3.93E-02 | .000 -.35 -8.24E-02 19 -•150 3.84E-02 | .009 -.28 -1.93E-02 20 -9.68E-02 3.27E-02 | .170 -.21 1.42 E-02 25 -9.03E-02 3.65E-02 | .461 -.21 3.36E-02 9 1 .36(*) 4.52E-02 | .000 .21 .52 2 •360 4.50E-02 j .000 .20 .51 3 •320 4.57E-021 .000 .17 .48 4 •190 4.89E-021 .009 2.46E-02 .36 6 -4.96E-02 4.52E-021 .999 -.20 I .10 7 5.05E-02 4.70E-021 .999 -.11 .21 8 •280 4.55E-021 .000 .12 .43 12 -6.73E-02 4.91 E-02 I .990 -.23 9.93E-02 13 8.37E-02 4.23E-02J .811 -5.97E-02 .23 14 .13 4.16E-02J .091 -8.20E-03 .27 15 5.98E-02 4.45E-02 j .992 -9.11 E-02 .21 19 .13 4.37E-021 .202 -2.20E-02 .27 20 •180 3.88E-021 .000 4.74E-02 .31 25 •190 4.20E-021 .001 4.30E-02 | .33 12 1 •43(*) 4.42E-021 .000 .28 •58 2 •420 4.40E-021 .000 .27 .57 3 •390 4.48E-021 .000 .24 .54 4 •260 4.80E-021 .000 9.50E-02 .42 6 1.78E-02 4.42E-02| 1.000 -.13 •17. 7 .12 4.61 E-02 | .400 -3.85E-02 .27 8 •340 4.45E-021 .000 .19 .49 9 6.73E-02 4.91 E-021 .990 -9.93E-02 .23 13 •150 4.12E-02J .020 1.12E-02 .29 230 14 •20C) 4.05E-02 .000 6.28E-02 .34 15 .13 4.35E-02 .185 -2.04E-02 .27 19 •19H 4.27E-02 .001 4.88E-02 .34 20 •25(*) 3.76E-02 .000 .12 .37. 25 •25(*) 4.10E-02 .000 .11 .39 1 .280 3.64E-02 .000 .16 .40 2 •27(*) 3.61 E-02 .000 .15 .39 3 •24(*) 3.71 E-02 .000 .11 .36 4 .11 4.10E-02 .364 -3.21 E-02 .25 6 -•13(*) 3.64E-02 .020 -.26 -9.83E-03 7 -3.32E-02 3.87E-02 1.000 -.16 9.80E-02 13 8 •1?(*) 3.68E-02 .000 6.72E-02 .32 j 9 -8.37E-02 4.23E-02 .811 -.23 5.97E-02 \ 12 -•15(*) 4.12E-02 .020 -.29 -1.12E-02 14 4.91 E-02 3.18E-02 .970 -5.88E-02 .16 15 -2.38E-02 3.56E-02 | 1.000 -.14 9.67E-02 19 4.25E-02 3.45E-02 .997 -7.46E-02 .16 20 9.52E-02O 2.80E-02 .049 1.49E-04 .19 25 .10 3.24E-02 .105 -8.07E-03 .21 1 •23(*) 3.56E-02 .000 .11 .35 2 •22C) '3.53E-02 .000 .10 .34 3 .19(*) 3.63E-02 .000 6.47E-02 .31 j 4 5.77E-02 4.03E-02 .984 -7.88E-02 .19 6 -.18(*) 3.56E-02 .000 -.30 -6.17E-02 7 -8.23E-02 3.79 E-02 .686 -.21 4.63E-02 14 8 .14(*) 3.60E-02 .006 2.08E-02 .27 9 -.13 4.16E-02 .091 -.27 8.20E-03 12 -•20(*) 4.05E-02 .000 -.34 -6.28E-02 13 -4.91 E-02 3.18E-02 .970 -.16 5.88E-02 15 -7.30E-02 3.47E-02 .735 -.19 4.48E-02 19 -6.65E-03 3.37E-02 1.000 -.12 .11 20 4.61 E-02 2.70E-02 .931 -4.54E-02 .14 231 25 5.27E-02 3.15E-02 .942 -5.41 E-02 .16 1 •30H 3.90E-02 .000 .17 .43 2 •30(*) 3.87E-02 .000 .16 •43 3 •260 3.96E-02 .000 .13 .40 4 .13 4.33E-02 .146 -1.61 E-02 .28 6 -.11 3.90E-02 .242 -.24 2.27E-02 7 -9.36E-03 4.11 E-02 1.000 -.15 .13 15 8 •22(*) 3.93E-02 .000 8.24E-02 .35 9 -5.98E-02 4.45E-02 .992 -.21 9.11 E-02 12 -.13 4.35E-02 .185 -.27 2.04E-02 13 2.38E-02 3.56E-02 1.000 -9.67E-02 .14 14 7.30E-02 3.47E-02 j .735 -4.48E-02 .19 19 6.63E-02 3.72E-02 .906 -5.99E-02 .19 20 •120 3.13E-02 .012 1.29E-02 .23 25 •13(*) 3.53E-02 .028 6.04E-03 .25 1 •240 3.80E-02 .000 .11 .37 2 •23C) 3.78E-02 .000 .10 .36 3 •19(*) 3.87E-02 .000 6.33E-02 .33 4 6.44E-02 4.24E-02 .974 -7.95E-02 .21 6 -.180 3.80E-02 .000 -.30 -4.68E-02 7 -7.57E-02 4.02E-02 .863 -.21 6.07E-02 19 8 •150 3.84E-02 .009 1.93E-02 .28 9 -.13 4.37E-02 .202 -.27 2.20E-02 12 -.190 4.27E-02 .001 -.34 -4.88E-02 13 -4.25E-02 3.45E-02 .997 -.16 7.46E-02 14 6.65E-03 3.37E-02 1.000 -.11 .12 15 -6.63E-02 3.72E-02 .906 -.19 5.99E-02 20 5.28E-02 3.01 E-02 t .917 -4.94E-02 .15 25 5.93E-02 3.42E-02 .923 -5.68E-02 .18 20 1 .180 3.23E-02 .000 7.39E-02 .29 2 •ISO 3.20E-02 .000 6.83E-02 .29 3 • 140 3.30E-02 .002 2.97E-02 .25 232 4 1.16E-02 3.74E-02 1.000 -.12 .14 6 -•23(*) 3.23E-02 .000 -.34 -.12: 7 -•130 3.48E-02 .018 -.25 -1.03E-02 8 9.68E-02 3.27E-02 .170 -1.42E-02 .21 9 -.18(*) 3.88E-02 .000 -.31 -4.74E-02 12 -.250 3.76E-02 .000 -.37 -.12 13 -9.52E-02(*) 2.80E-02 .049 -.19 -1.49E-04: 14 -4.61 E-02 2.70E-02 .931 -.14 4.54E-02 15 -.12(*) 3.13E-02 .012 -.23 -1.29E-02 19 -5.28E-02 3.01 E-02 .917 -.15 4.94E-02 25 6.54E-03 2.77E-02 1.000 -8.73E-02 .10 1 •18(*) 3.61 E-02 .000 5.43E-02 .30 2 •17(*) 3.59E-02 .000 4.86E-02 .29 3 .14(*) 3.68E-02 .019 1.04E-02 .26 4 5.07E-03 4.07E-02 1.000 -.13 .14 6 -•24C) 3.61 E-02 .000 -.36 -.11 7 -.13(*) 3.84E-02 .034 -.27 -4.70E-03 25 8 9.03E-02 3.65E-02 .461 -3.36E-02 •21 9 -•19(*) 4.20E-02 .001 -.33 -4.30E-02 12 -.25(*) 4.10E-02 .000 -.39 -.11 13 -.10 3.24E-02 .105 -.21 8.07E-03 14 -5.27E-02 3.15E-02 .942 -.16 5.41 E-02 15 -.13(*) 3.53E-02 .028 -.25 -6.04E-03 19 -5.93E-02 3.42E-02 .923 -.18 5.68E-02 20 -6.54E-03! 2.77E-02 1.000 -.10 8.73E-02 * The mean difference is significant at the .05 level.