Articles | Volume 8, issue 4
Research article
21 Dec 2020
Research article |  | 21 Dec 2020

Inertial drag and lift forces for coarse grains on rough alluvial beds measured using in-grain accelerometers

Georgios Maniatis, Trevor Hoey, Rebecca Hodge, Dieter Rickenmann, and Alexandre Badoux

Quantifying the force regime that controls the movement of a single grain during fluvial transport has historically proven to be difficult. Inertial micro-electromechanical system (MEMS) sensors (sensor assemblies that mainly comprise micro-accelerometers and gyroscopes) can used to address this problem using a “smart pebble”: a mobile inertial measurement unit (IMU) enclosed in a stone-like assembly that can measure directly the forces on a particle during sediment transport. Previous research has demonstrated that measurements using MEMS sensors can be used to calculate the dynamics of single grains over short time periods, despite limitations in the accuracy of the MEMS sensors that have been used to date. This paper develops a theoretical framework for calculating drag and lift forces on grains based on IMU measurements. IMUs were embedded a spherical and an ellipsoidal grain and used in flume experiments in which flow was increased until the grain moved. Acceleration measurements along three orthogonal directions were then processed to calculate the threshold force for entrainment, resulting in a statistical approximation of inertial impulse thresholds for both the lift and drag components of grain inertial dynamics. The ellipsoid IMU was also deployed in a series of experiments in a steep stream (Erlenbach, Switzerland). The inertial dynamics from both sets of experiments provide direct measurement of the resultant forces on sediment particles during transport, which quantifies (a) the effect of grain shape and (b) the effect of varied-intensity hydraulic forcing on the motion of coarse sediment grains during bedload transport. Lift impulses exert a significant control on the motion of the ellipsoid across hydraulic regimes, despite the occurrence of higher-magnitude and longer-duration drag impulses. The first-order statistical generalisation of the results suggests that the kinetics of the ellipsoid are characterised by low- or no-mobility states and that the majority of mobility states are controlled by lift impulses.

1 Introduction

River sediment transport is a critical process in landscape evolution (Tucker and Hancock2010), controls river morphology and ecology (Recking et al.2015) and affects river engineering (Van Rijn1984). The study of the two-way relationship between transport processes and the corresponding morphology has a long history (e.g. Gilbert and Murphy1914), using approaches and mathematical conceptualisations that range from deterministic (e.g. Gilbert and Murphy1914; Shields1936; Ali and Dey2016) to probabilistic (e.g. Einstein1937; Grass1970; Ancey et al.2008).

Fluvial sediment transport is a complex two-phase flow defined by (a) hydraulics (Kline et al.1967; Nelson et al.1995; Papanicolaou et al.2002), (b) sediment properties and arrangement (Ashida and Michiue1971; Komar and Li1988; Kirchner et al.1990; Buffington et al.1992; Hodge et al.2013; Prancevic and Lamb2015), (c) flow history across timescales (Shvidchenko and Pender2000; Diplas et al.2008; Valyrakis et al.2010; Phillips et al.2018; Masteller et al.2019) and (d) biological and chemical processes that can rearrange or stabilise sediment (Johnson et al.2011; Vignaga et al.2013; Johnson2016).

To analyse the motion of a grain resting on a riverbed that is sheared by a turbulent flow (Dey and Ali2018), a large group of laboratory and theoretical studies use an implicit (fixed) reference frame. Historically, such analyses have been deterministic (implementing a single threshold shear stress or force at which grains are entrained (Gilbert and Murphy1914; Shields1936; Yalin1963; Iwagaki1956; Ikeda1982; Dey1999). However, stochastic descriptions arguably capture better the complex particle–fluid interplay since near-bed turbulence which drives gain motion is inherently stochastic (Einstein1937; Grass1970; Papanicolaou et al.2002; Marion and Tregnaghi2013). Coupled with advances in monitoring techniques (e.g. Papanicolaou et al.2002; Fathel et al.2016), stochastic treatments have led to Lagrangian, primarily numerical, formulations being applied to the full range of motion (McEwan et al.2004; Bialik et al.2015). The term Lagrangian means that sediment flow is observed from the perspective of individual mobile sediment grains and not a domain fixed in time and space (as in traditional Eulerian approaches using fixed xy and z coordinates; Ballio et al.2018). The turbulence impulse approach (Diplas et al.2008; Valyrakis et al.2010; Celik et al.2010) accounts for both the magnitude and the duration of the hydraulic forcing and is often categorised as a stochastic approach (Dey and Ali2018) since the stochastic nature of local turbulence is accounted for by the integration of turbulent forces acting on the grain over time. However, the grain forces are treated deterministically through a detailed treatment of the force balance during incipient motion. Finally, the spatio-temporal approach (Coleman and Nikora2008) is different as the equations of motion are applied separately for the fluid (in a spatially averaged domain) and the sediment particles, linking the mode of transport with the scales of turbulence (Bialik et al.2015).

Field experiments tracing individual sediment grains are in principle Lagrangian (Hassan and Roy2016). Lagrangian analytical models have been developed following advances in monitoring techniques that allow tracking of individual grains, including magnetic (e.g. Schmidt and Ergenzinger1992; Hassan et al.2009) and RFID tracers (e.g. Schneider et al.2014; Tsakiris et al.2015). An important milestone in the development of Lagrangian approaches for sediment transport was the introduction of discrete particle modelling techniques in simulations (McEwan et al.2001, 2004; Schmeeckle and Nelson2003), which opened up the prospect for upscaling the Lagrangian metrics.

Lagrangian measurements find direct application in coarse-grain gravel bed and bedrock river environments (e.g. Hassan et al.1992, 2009; Ferguson et al.2002; Hodge et al.2011; Liedermann et al.2012), and the morphological impact of Lagrangian dynamics in those environments is pronounced (Hodge et al.2011). For example, the inertia of the typically larger particles transported in these streams has been identified as one of the factors contributing to inaccurate predictions of transport rates (Buffington and Montgomery1997; Bunte et al.2004; Singh et al.2009). Equally important is the lack of information on the energy transfer between these large particles and the riverbed, particularly during impact. Recent experiments (Gimbert et al.2019) show how this energy transfer can be inferred from seismic measurements, opening the way for testing hypotheses that relate to river reach scale processes (e.g. Burtin et al.2014). Finally, for a complete understanding of these interactions, the rotational component of grain movement cannot be ignored (Niño and García1998).

A particular advance in monitoring technology has been the development of sediment grain scale inertial sensors which record at high frequency the accelerations and angular velocities experienced by grains during entrainment and motion (Kularatna et al.2006; Akeila et al.2010; Frank et al.2015; Maniatis et al.2013; Gronz et al.2016; Maniatis et al.2017). These applications became possible after the development of compact MEMS (micro-electromechanical mechanical system) inertial measurement units (IMUs), assemblies of 3D MEMS accelerometers and 3D MEMS gyroscopes, which overcome many technical difficulties posed by older instrumentation (mainly caused by limited storage capacity; e.g. Ergenzinger and Jupner1992; Spazzapan et al.2004). The focus of those works is the development of an IMU-based sensor assembly (IMU enclosed in a grain or a purpose-specific grain-shaped artificial enclosure) that can successfully measure grain dynamics.

MEMS-IMU sensors measure the acceleration and angular velocity of the grain, which can be used to calculate the net force acting on the grain. For the most accurate measurement of this force, sensors should be located at the grain's centre of mass. Data collected from within grains undergoing transport has potential to describe the timing of motion, forces acting on the grain and grain location. As a grain moves, its centre of mass moves, and so the reference point for the force measurements is mobile. The latter means that the IMU measurements need to be transformed to a frame of reference that can be understood by an observer. Generally, an IMU accelerometer is a non-inertial frame fixed within the mobile body frame of the sensor assembly.

In theory, the accelerations recorded by the IMU could be integrated to calculate grain velocity and integrated again to reveal location, a process referred to as dead reckoning. One long-term goal for this approach would be to use IMUs to track a large number of grains through fluvial systems. However, real fixed (“strap-down”) IMUs based on MEMS are not suitable for these integrations since the data contain several sources of uncertainty, including signal noise and nanoscale misalignment of sensor axes. With sensors that are cheap enough to be deployed in large numbers, the accumulation of errors means that they cannot be used for 3D tracking of long-term unconstrained motions (Woodman2007; Kok et al.2017; VectorNav2016). This problem is well known in the fields of navigation and electrical engineering, and the modelling of IMU errors is a significant research area (Zekavat and Buehrer2011) since the applications of this technology are numerous (Gebre-Egziabher et al.1998; Grewal et al.2007).

Despite the limitations, IMU sensors have been considered to be a suitable technique for measuring grain motion (e.g. Gronz et al.2016). Recently Gimbert et al. (2019) used accelerometers to measure particle bed impacts in order to complement seismic measurements. However, for field deployments, relevant sensors have only been used as start-and-stop motion sensors (Olinde and Johnson2015).

The first goal of this paper is to introduce a simple rigid-body model that connects measurements derived from an idealised IMU with existing models for grain motion. For this model to be successful, it is necessary to map the IMU body frame dynamics to the reference frame of motion (flume or riverbed). This resolution allows the inertial measurements to be related to the forcing on the particle and defines an explicit and unambiguous threshold for particle motion. The second goal is to introduce the calculation of inertial impulses above the drag and lift thresholds of motion, following the example of Diplas et al. (2008), for grain entrainments and short-transport events. The calculations are performed for observations made in a set of flume entrainment experiments using two sensor assemblies: one spherical and one ellipsoidal. We apply the same analysis to a time series of successive transport events measured with the ellipsoid sensor in a steep Alpine river (Erlenbach, Switzerland), calculating the force regime and the generated impulses during grain motion. Finally, we discuss how the combined dataset (flume and field experiments) for the ellipsoid sensor can be used for bootstrap calculations leading to the generalisation of the derived measurements.

2 Force measurements and a Newton–Euler regime for sediment motion

An IMU accelerometer records accelerations of the grain, which can be converted to forces acting on the grain by multiplying by grain mass. However, those accelerations are both fixed (e.g. due to gravity) and variable (e.g. due to fluid forces applied to the particle). The measurements are difficult to interpret because fixed and variable accelerations cannot be decoupled after the accelerometer begins to move, and these accelerations are perceived differently from different reference frames. If the sensor is static and gravity is recorded from the reference frame of the sensor, it is then possible to derive a reliable orientation measurement. If the sensor moves and the sensor frame accelerations are compensated for gravity (by removing gravity from the raw accelerometer measurement), then the sensor will record the 3D components of the resultant or net force that mobilises the particle. This resultant is the force that can be observed from an observer who is static in relation to the particle. In Appendix A, we provide a technical explanation of the differences between non-inertial (mobile sensor frame) and inertial (as observed from a fixed frame external to the sensor) accelerations and demonstrate the necessary transformations.

The physics of IMU sensors also define the main difference between IMUs and other force sensors that have been used to monitor grain motion, such as load cells (e.g. Schmeeckle et al.2007; Lamb et al.2017). Load cells measure the total vertical and horizontal forces applied to the grain prior to the onset of motion, but their application is also limited because they prevent the grain from being fully entrained. In the case of a grain sitting on a horizontal bed and being subjected to an increasing lift force, a load cell records an increasing force up until the point at which the lift force is equal to the weight and the grain is entrained. In contrast, the gravity-compensated IMU will record zero acceleration (and therefore zero net force) until the point of entrainment, after which the grain starts to accelerate. In flowing water, grains do not transition instantly from being completely stationary to being fully entrained (Garcia et al.2007), and the IMU will also record the movement of grain vibrations prior to entrainment. Consequently, load cells and IMUs measure different parts of the transport process; load cells record the forces that are applied to the grain prior to entrainment, i.e. to the point where the forces balance the weight of the grain, whereas IMUs measure the forces that accelerate the grain from this balanced position.

The following derivations rely on the transformation of the IMU accelerations from the mobile reference frame of the monitored particle (frame b) to a local static frame where the applied hydraulic, particle and gravitational forces are analysed (frame r hereafter). The x and y directions of the r frame (rx and ry hereafter) are, respectively, parallel and perpendicular to the mean direction of the flow and the river or flume bed. The z direction (rz hereafter) is normal to the flow and the riverbed. In addition, an independent static frame i is used, the origin of which coincides with the centre of the Earth, as a stable reference for the initial alignment of frames b and r (Appendix A). br transformations rely on the successive tracking of the relative orientation between the two frames, which are represented here using quaternions (Appendices A and G1).

Figure 1Forces in the Newton–Euler model. The diagram shows the linear forces applied on the centre of the mass of the target particle. FDnet is the net drag force, FLnet is the net lift force and FG is the force of gravity (Eq. 1), all analysed in the local static frame r. The rotational component (torques generated by the surface force frot around the centre of the mass of the target particle) is analysed on the body frame of the particle (b, here depicted as aligned with the static frame r; the general case is discussed in Appendix A). ϕ is the downstream slope. For the results presented in this work rx is the downstream direction. Also, for the flume experiments there is no cross-stream slope (θ=0).


Combined measurements of grain acceleration and angular velocity allow direct calculation of the forces and turning moments acting at the grain's centre of mass. This type of model formulation is the Newton–Euler model in the rigid-body-dynamics literature (O'Reilly2008). For a spherical particle resting on an inclined bed, irrespective of the degree of exposure to the flow, the Newton–Euler regime is defined as follows:


Fnet and Tnet are net forces and torques applied on and around the centre of mass of the particle. For the right-hand elements of Eqs. (1) and (2), m is the particle mass and Icm is its moment of inertia (Appendix C). arx, ary and arz are linear accelerations resolved in the static frame r (Fig. A1) measured using accelerometer measurements and Eq. (A4). αbx, αby and αbz are rotational accelerations resolved in the mobile frame but fixed on the particle frame b such that

(3) α b x α b y α b z T = d ω b x d t d ω b y d t d ω b z d t T ,

where ωbx, ωby and ωbz are angular velocities within the body frame reference frame as recorded by gyroscope measurements.

The left-hand elements of Eqs. (1) and (2) describe the forces applied to the particle. Hereafter, the lower-case vectors and scalars (e.g. terms f and frot) will refer to the interactions between hydraulic forces (turbulence) and particle forces (support forces and friction) that are not measured directly by an IMU and which cannot be decoupled using IMU measurements. The vector f is the linear component of those interactions, applied on the centre of mass of the particle. The vector frot is the component of those interactions applied on the surface of the particle, and it generates the torques around the centre of its mass (tangential to the particle's radius). The vector FGr defines the force of gravity rotated in the r frame and compensated for hydrodynamic effects (e.g. buoyancy) as

(4) F Gr = R i r q W s i f v = F Gr x F Gr y F Gr z .

The matrix R(irq) denotes the orientation of the slope in 3D (the relative orientations between the frames i and r, Eq. (A5). Wsi is the immersed weight of the spherical particle equal to mbg, where mb is its immersed mass (Papanicolaou et al.2002) and g is the acceleration of gravity, both acting at its centre of mass. Wsi has a constant direction in frame i and is rotated in r using R(irq). fv=1+[0.5ρ/(ρs-ρ)] accounts for the hydrodynamic mass effect (Papanicolaou et al.2002; Celik et al.2010), and ρ and ρs are the densities of the water and the particle, respectively. For our ellipsoid, we calculate FGr assuming a sphere with the same volume as the ellipsoid, since resolving the hydrodynamic mass coefficient for an ellipsoid in 3D is beyond the scope of this work. The magnitude of FGr can be measured by deploying an IMU as an orientation sensor (Appendix A1).

The magnitude of FGr components along the drag direction (rx) is given by the vector

(5) F G D = F Gr x .

For the direction normal to the bed (rz, direction of lift hydraulic force) the FGr component is the vector

(6) F G L = F Gr z .

Similarly, we can use the linear accelerations of Eq. (1) to separate the components of the net force (Fnet) along the drag and lift directions as

(7) F Dnet = m a r x


(8) F Lnet = m a r z .

2.1 Inertial measurements and the threshold of motion

The linear force threshold of motion is defined for ar=[0 0 0]T, representing the explicit state where the forces are balanced (a resultant force Fnet of Eq. (1) equal to 0) and the particle is not moving. Since all the linear forces are transformed in the r frame, exceeding this threshold relates to a motion along the rx drag direction (threshold defined at FDnet=0) and/or a motion along the rz lift direction (threshold defined at FLnet=0). Given the force balance of Eq. (1), the threshold of motion is also the point where the resultant of the combined forces represented by f (which includes the hydraulic lift and drag forces) balance the rotated gravity forces. For the drag direction, the component of f is given by the vector

(9) f D = F Dnet - F G D ,

and the threshold of motion is at the point where fD=-FGD.

If FDnet exceeds 0 then fD exceeds the component of gravity FGD towards the positive (+) drag direction (downstream). Consequently, the condition FDnet>0 is the exact equivalent of the condition fD>-FGD. We use the second description (for both the drag and the lift components) in order to highlight the fact that the f interaction primarily balances the components of gravity defined by the mean orientation of the bed.

For the lift direction, the component of f is given by

(10) f L = F Lnet - F G L .

Similarly to the drag components, the condition FLnet>0 is the exact equivalent of the condition fL>-FGL (the particle moving towards the rz positive direction, upwards).

The rotational threshold of motion is defined by the state where the balance of torques around the centre of mass of the particle are balanced (Tnet=0). In the Newton–Euler model introduced here, the sum of torques captures a spinning rotational component defined by the product of the moment of inertia of the rotating particle (Icm) and the body frame angular accelerations (Eqs. 2 and 3). It is important to note that this rotation differs from the orbital rotation defined around the centre of mass of a supporting particle, which is a common description in the hydraulic literature (e.g. Papanicolaou et al.2002).

The torques are analysed in the body frame of the particle (b) and a non-directional description of the rotation threshold is given by the norm of angular accelerations (αb) exceeding 0. After Eqs. (2) and (3) the critical condition for particle rotation is given by

(11) T net = I cm α b = I cm d ω b d t = f rot R 0 .

It is useful to re-state here that all the derivations so far are for a spherical particle (the equivalent of Eq. (11) for an ellipsoid is given in Appendix C). Using Eq. (11), it is possible to estimate the magnitude of the tangential component frot (Appendix C). The calculation reveals that, for the scale of the particles discussed here, the effect of the tangential force, in terms of force magnitude applied to the particle, is negligible in comparison to the linear forces. This is explained by the dependency of the sum of torques (Tnet) on the moment of inertia Icm of the particle, which is generally a very small number even for relatively coarse particles (e.g. for the spherical particle introduced in Sect. 3 it is equal to 0.00085 kg m−2). For this reason, we will only focus on the linear net force and the interaction f, which we calculate as the difference between the net force and the gravitational components. However, we demonstrate the scale difference between the rotational and the linear components of particle motion in Fig. 3 and Appendix C.

2.2 General kinematics and impulse

For completeness we note that, if ar=[arxaryarz]T is the 3D vector of linear acceleration of Eq. (1) in the r reference frame, grain linear velocities can be calculated by single integration. The velocity in the r frame is υr(t)=υr(0)+0tar(t)dt) and the total kinetic energy can by calculated from these velocities as K=12mvr2+12Icmωb2, with 12mvr2 being the translational and 12Icmωb2 the rotational component. The condition K=0 represents the threshold of rolling for any rigid body. However, this is not a direct equivalent to the typical rolling mode of entrainment because of the differences in the definition of rotation (spinning vs. orbital) and the fact that here there is no assumption about the slipping condition as K can be used to describe a rolling-with-slipping motion.

To account for both the duration and the magnitude of a force, the impulse I for duration δt starting from the time ti is defined as

(12) I = t i t i + δ t F ( t ) d t .

The subsequent analysis focuses on the calculation of impulses for specific time durations δt (Diplas et al.2008; Celik et al.2010; Valyrakis et al.2010). We calculate the impulse of the resultant interaction f above the threshold of motion and towards the positive drag and lift direction. This happens when the balancing gravity components FGr (Eqs. 5 and 6) are exceeded, thus for the conditions

(13) f D > - F G D


(14) f L > - F G L ,

which, as explained above, can only be satisfied when the magnitude of the net force exceeds 0 (Fnet>0, net force measured using gravity-compensated accelerometer data, Eq. A4). In practice, we only account for the net forces that exceed the noise threshold of the accelerometer sensor (Fig. 3, YEI2014; Maniatis2016). Finally, it is important to note that the calculated impulses are transferred to the particle from fluid turbulence and coherent flow structures; however this transfer is not described in this work. Here, the impulses capture directly the flow–particle interaction. The terms impulse, net impulse and inertial impulse are used interchangeably hereafter, and they relate to impulses of the non-zero net force.

3 Laboratory and field experiments

Two sensor assemblies were deployed, one sphere and one ellipsoid (described in Maniatis2016). The 90 mm diameter, 1.019 kg, sphere is solid aluminium with a symmetrical cavity for the IMU centred at the origin of the sphere. The ellipsoid (axes 100, 70, 30 mm), made of the same material, weighs 0.942 kg. The cavity in the ellipsoid was designed to ensure that the IMU axes align with the principal axes of the whole device. The density of both devices after the cavity cut is 2670±3 kg m−3, approximating the density of quartz (2650 kg m−3). The measuring unit is the TSS-DL-HH-S sensor from YEI-TechnologiesTM (YEI2014), equipped with a gyroscope (±2000 s−1 sensitivity) and an accelerometer with a maximum range of ±400 g. The acceleration range is one of the main reasons for the selection of this IMU as lower-range accelerometers, particularly those in the very common ±20 g range, are not suitable for capturing forcing in natural environments (Maniatis et al.2013). The factory maximum sampling frequency is 250 Hz. The nominal sampling frequency of the sensor (used for all the flume and field experiments presented in this work) is 50 Hz, which permits constant use for approximately 5 h (LiPo rechargeable battery). This frequency was chosen after considering displacements of 15 particle diameters per second reported by Drake et al. (1988). This is adequate for capturing the dynamics of the particles if collisions and strong interactions with the bed are excluded (different types of higher-frequency piezoelectric sensors are more suitable for the full measurement of impacts). The measuring unit was calibrated through a series of shaking table and rolling drop experiments, which are described in Maniatis (2016, Sect. 6), along with the corresponding filtering workflow. During this calibration the noise threshold for the accelerometer was defined at 1.1±0.23 N, and the value of 1.2 N is used for the following presentation. The gyroscope has a much lower noise threshold, which was calculated to be <0.0001 s−1.

3.1 Laboratory entrainment experiments and separation of impulsive events

All flume experiments took place in an ArmfieldTM flume of 7 m effective length and 0.9 m internal width and at a fixed bed slope of S=0.02. This slope corresponds to an orientation by the quaternion irq=[0.92, 0.14, −20, 0.28] (Appendix A). irq was measured by aligning the X-IMU axis with the centreline of the flume, the Y-IMU axis with the transverse direction and the Z-IMU axis with the direction normal to the bed. The positive X-IMU axis coincides with cross-sectionally averaged flow direction. A bed of plastic hemispheres of the same diameter (90 mm) as the spherical device was constructed. The hemispheres were glued to form a 0.5 m (L) × 0.9 m (W) section, and the whole section was placed at the point which allowed for the test particle to be 4.5 m from the upstream boundary of the flume. A thin layer of 1.5 mm uniform sand was glued to both the hemisphere section and the upstream flume surface. The section upstream of the hemispheres was filled with very densely packed non-uniform, rounded gravel (D50=0.015 m), enabling the development of turbulent flow. No sediment transport occurred from the upstream gravel section during the entrainment experiments. The hemispheres were glued in positions that produced a 0.045 m protrusion of the sphere above the top of the hemispheres, although the sand layer made the protrusion higher, by partially filling the gaps between the hemispheres, and non-uniform around the test particle (≈0.050 m). The ellipsoid was only supported by the hemispheres and was fully exposed to flow, thus having a protrusion equal to its c axis (≈0.03, Fig. 2).

Figure 2Laboratory setting and initial alignment. The diagram shows the arrangement of the bed of hemispheres and the test position (4.5 m downstream of the entrance of the flume). Upstream of the bed of hemispheres, the flume was filled with densely packed gravel (D50=0.015 m, black layer upstream of the hemisphere section in the graph) to enable the development of the flow. The sphere diameter is 0.09 m, and ellipsoid b axis is 0.07 m. The height of the centroid above the crest of the glued hemispheres is approximately 0.02 m for the sphere and 0.015 m for the ellipsoid. The depths correspond to the hydraulics at the point of entrainment, summarised in Table D1. r stands for riverbed (flume in this case) reference frame, b for body frame and i for the inertial reference frame (irq=[0.92, 0.14, −20, 028]), Appendix A). S is the slope of the channel (0.02). The sketch also depicts the initial alignment for each device (bi for the sphere and br for the ellipsoid). While the spherical particle was in full contact with the bed (hemispheres and coating/filling gravel), the setting could result to a film of water flowing underneath the ellipsoid.


For the experiments, the spherical device was placed on the flume centreline in a saddle position between four bed hemispheres, and the three sensor axes were aligned with the inertial frame i (ibq=[1, 0, 0, 0], Appendix B). After positioning the sensor, the discharge was increased at a constant rate of 0.028 L s−2 until the particle was entrained. Acceleration and rotation were measured for the duration of the flow increase, throughout the sensor movement, and for a further 10 s after it stopped moving. All the experiments were videoed at 60 frames per second using a standard GoPro Hero 7. The same experimental protocol was followed for the ellipsoid device, differing only in that the particle was initially aligned to the frame of reference of the flume bed (ibq=irq=[0.92, 0.14, −20, 0.28]). This resulted in the X-sensor long axis coinciding with the flume's x direction (the direction of the flow, Fig. 2).

Ten entrainment experiments were conducted with each device. We define entrainment as when the particle moves by one particle diameter or b axis length for the ellipsoid. Having identified the time when the grain has moved by this distance from the video, the timing of the vibrations which directly and continuously preceded entrainment was also determined from the video. Many periods of vibration which do not lead to entrainment were recorded by the sensor and are visible on the videos. Drag and lift forces as well as the duration and the inertial impulses for cases where the grain started to move (fD>-FGD and/or fL>-FGL) were calculated using the derivations of Sect. 2.1. For the spherical sensor, the gravitational components are FGD=3.99 N and FGL=-7.25 N (Fig. 3). The equivalent values for the ellipsoid, using the geometry of a sphere of equal volume, were FGD=5.11 N and FGL=-9.28 N. The critical discharge for the sphere was 24.8±1.8 L s−1, which corresponds to a measured depth of 0.095±0.015 m, measured from the top of the supporting hemispheres. The critical discharge for the ellipsoid was 45.2±2.2 L s−1, which corresponds to a measured depth of 0.12±0.02 m (τ*=0.01 and 0.02, respectively; τ*=HS(ρp-ρf)D; Appendix D). The hydraulic parameters at mean critical discharges are calculated in Table D1. The flume experiments have a particle-diameter-to-flow-depth ratio close to 1 (d/H1) despite the fact that the particles were fully submerged at the critical discharge for all the experiments (Fig. 2). The tested conditions are relevant to the entrainment of coarse particles in steep mountain streams, but they should not be directly generalised to other bedload transport regimes despite the generality of the Newton–Euler model presented in Sect. 2.

Figure 3Example flume and field experiments (a, b). Calculated drag (a) and lift forces (b) for the same flume experiment using the spherical sensor. The vectors FDnet and fD record the magnitudes of the net force and the hydraulics–resistance interaction on the 2D plane parallel to the flume bed, respectively. The vectors FLnet and fL record the net force and the hydraulic resistance interaction, respectively, along the normal to the flume bed direction. The gravitational components (FGD=3.7 N and FGL=-7.25 N) were calculated using Eq. (4). The vertical dashed line (t=10.1 s) shows the exact point of the entrainment of the particle as determined from the video recording. The NT lines indicates the noise threshold of the accelerometer sensor; the calculated impulses concern sequences of point forces that exceed NT. (c) Net drag forces recorded during one field experiment (Erlenbach). The vertical line indicates the time of 1 s after the release of the sensor (see Sect. 3.3). The calculation of the tangential force magnitude (frots and frote for the sphere and the ellipsoid, respectively) is based on Eqs. (C1) and (C5).


3.2 Probabilistic impulse threshold for motion

Entrainment was observed independently from video recordings which were synchronised with the experiments from the start of the flow increase (Sect. 3.1). These observations are used to calculate statistically the probability of entrainment as a function of the impulse of the interaction f exceeding gravity. Following the framework presented in Maniatis et al. (2017), the exact time of entrainment was noted in the video recording, and the derived inertial impulses were separated into a binary pre- and post-entrainment dataset. Logistic regression was used to describe the probability of entrainment, with Pr>0.5 defining the threshold of motion. Following the conceptualisation of Grass (1970), exceeding that threshold relates to impulses that are able to fully dislodge the particle, in contrast to the conditions that relate to pre-entrainment vibrations. The difference here is that this conceptualisation is applied to events in which impulses exceed the thresholds defined in Eqs. (5) and (6). Video recording was not possible in the field setting, so this calculation is only presented for the laboratory experiments.

3.3 Field testing

Field experiments took place within a 5 m long straight and confined reach of the Erlenbach mountain stream in Switzerland, approximately 15 m upstream of the concrete channel section and 55 m upstream of the sediment retention basin in which continuous bedload transport measurements have been made during the past 30 years (Turowski et al.2011; Rickenmann et al.2012). The stream has a step-pool morphology allowing the sensor to be retrieved from pools, so the ellipsoid sensor was submerged on a bare bedrock section close to the edge of a step (average slope S=0.1, cross-averaged flow depth H=0.1 m, τ*=0.095), aligned to the same orientation as the riverbed (ibq=irq=[0.50 −0.39, 0.34, −0.68], assumed parallel to the banks and the cross-averaged flow direction, at the approximate centre line of the stream) and allowed to be transported until it stopped moving and remained immobile for at least 10 s (Fig. 3c). In all the experiments the sensor was entrained fully and immediately as there was no vibration in situ. The first 1 s of each transport event was removed from the data as the effect of holding and releasing the sensor was still present. Ten transport events were recorded and processed similarly to the flume experiments. The corresponding fD>-FGD and fL>-FGL (Fig. 3c) were calculated similarly to the flume experiments. The average travel distance for each transport event was 2±0.43 m (from the point of release to the point of deposition, tape measurements), and the average event duration, after the first second of release, was 3±0.6 s (Fig. 3). To establish a representative orientation for the reach in relation to the orientation of gravity, the IMU was aligned parallel to the approximate centreline (X-IMU axis, X+= cross-sectionally averaged flow direction), transverse (Y-IMU axis) and normal to the bed (Z-IMU axis) directions within the stream. The hydraulic parameters are summarised in Table D1.

4 Results

The flume experiments demonstrate the differences between the spherical and the ellipsoid particle during incipient motion (Figs. 4 and 5). For the sphere, drag and lift impulses over the gravity forces (fD>-FGD and fL>-FGL, Eqs. 13 and 14) occur for similar durations and generate impulses of similar magnitude (fD>-FGD impulse median = 0.45 N s−1, fL>-FGL impulse median = 0.46 N s−1). The relationship between the duration of exceedance events and the generated impulse follows an approximately linear trend, although variability is marginally higher for the relationship between drag impulses I and corresponding durations (t). For the relationship I vs. t, R2=0.78 (p value <6.52×10-16) for the drag events and 0.89 (p value <4.2×10-9) for the lift events (Fig. 4a).

Figure 4Inertial impulses and duration of threshold exceedance events for laboratory experiments. Impulses of all the inertial forces that exceeded the gravity forces. For the spherical particle (a) the drag impulse median (median of fD>-FGD events) is 0.01 N s−1 higher than the lift impulse median (median of fL>-FGL events). For the ellipsoid (b) the equivalent difference is 0.013 N s−1. The relationship between the duration (t) and the impulse (I) during the exceedance events is linear for both the sphere and the ellipsoid (and for both fD>-FGD and fL>-FGL) events. The entrainment of the ellipsoid is more dependent on short and low lift impulses than the sphere, demonstrating the effect of shape on the inertial dynamics.


Figure 5Probabilistic inertial impulse threshold for laboratory experiments. Logistic regression of the probability of entrainment for the spherical (a) and ellipsoid (b) particles. The calculation is based on the combination of video recordings and inertial impulse measurements during drag and lift threshold exceedance (fD>-FGD and fL>-FGL events). For the sphere there is little statistical difference between the calculated inertial impulses as over 95 % of the values relate to an entrainment event (the probability threshold 0.5 is always exceeded). However, there is significant variability in this calculation (wide 95 % confidence intervals), which indicates the wide range of impulses that can lead to entrainment of the sphere. For the ellipsoid, the probabilistic lift inertial impulse threshold relaxes to 0.27 N s−1 (blue vertical line, b) and the drag threshold relaxes to 0.74 N s−1 (red vertical line, b).


The results from the ellipsoid sensor demonstrate a strong influence of the lift forces. Exceedance impulses occur for similar durations and magnitudes; however there is a strong bias of the lift distribution towards the shorter and low-impulse events. The drag duration and impulse distributions include more and higher magnitude outliers than the lift distributions (fD>-FGD impulse median = 0.08 N s−1, fL>-FGL impulse median = 0.022 N s−1). For the ellipsoid, the relationship I vs. t has R2=0.95 (p value <2.2×10-16) for the drag events and 0.67 (p value <2.2×10-16) for the lift events (Fig. 4b). For all these threshold-exceeding events the sensor was vibrating until entrainment as observed from both video and IMU data.

Using the video recording observations, the impulse thresholds for entrainment were approximated with logistic regression. The probability of 0.5 corresponds to the threshold impulse for which the probability changes from the particle being more likely to be at rest to being more likely to be entrained. In this context, with this approximation we calculate a gradational threshold of entrainment (Begin and Schumm1979) and not an absolute one. The probability of entrainment as a function of impulse (Fig. 5a and b) highlights the control of short lift events on the entrainment of the ellipsoid. The impulse threshold for the sphere is close to 0, as all the approximated probabilities exceed 0.5. However, there is significant variability in this calculation (wide 95 % confidence intervals), which indicates the wide range of impulses that can lead to entrainment of the sphere. In contrast, the entrainment of the ellipsoid demonstrates a dependency on lift impulses as the lift threshold is lower (ellipsoid lift impulse threshold =0.27±0.03 N s−1) and the drag threshold is approximated with less confidence (ellipsoid drag impulse threshold =0.74±0.27 N s−1).

Finally, the results from the field experiments (Fig. 6) indicate similar statistical behaviour to the laboratory experiments but higher variability. Drag forces are of higher magnitude and duration than the lift forces (fD>-FGD impulse median = 0.13 N s−1, fL>-FGL impulse median = 0.08 N s−1), but there is an abundance of low-magnitude lift impulses that affect strongly the motion of the ellipsoid. In the Erlenbach, the duration of the exceedance events is also a proxy for the generated impulse, with the relationship I vs. t being more linear for the lift events compared to the drag events (I vs. t R2=0.66 and p value <2.2×10-16 for the drag events, and R2=0.88 and p value <2.2×10-16 for the lift events).

Figure 6Inertial Impulses and duration of threshold exceedance events for field experiments. Impulses of all the inertial forces that exceeded the gravity forces. During short-transport events (average travel distance =2±0.43 m) the drag inertial impulse median (median of fD>-FGD events) is 0.05 N s−1 higher than the lift inertial impulse median (median of fL>-FGL events). The relationship between the duration (t) and the inertial impulse (I) is linear for both fD>-FGD and fL>-FGL events. During in situ transport the drag forces are of comparable magnitude and duration; however, short and low-magnitude fL>-FGL impulses have a strong influence on the motion of the ellipsoid.


5 Discussion

5.1 IMU sensors and geomorphological applications

The advantage of using an IMU sensor for capturing grain motion is that the sensor solves a complex force and torque balance and removes any ambiguity in whether or not a test particle is in motion, as motion leads to the explicit thresholds Fnet and/or Tnet exceeding 0. Entrainment is captured directly and, assuming correct sensor calibration, robustly. IMUs can be a useful tool for geomorphologists since they offer a realistic prospect for monitoring particle motion during transport without invasive apparatus, which is not possible with standard equipment, especially in field applications (e.g. PIT tracers). At the same time, it is important to recognise that exceedance of the explicit thresholds above does not always produce complete dislodgment of the particle and also does not directly describe the modes of transport in the context that is commonly assumed for sediment hydraulics (e.g. differences in spinning and orbital rotations, Sect. 2.1). For a complete understanding and effective prediction of grain motion both the hydraulic and the particle forces need to be measured, analysed and decoupled from the inertial forces we measure in this study.

Further, there has been a recent rapid increase in use of IMU sensors, but most off-the-shelf IMU sensors are not suitable for the range of forces characterising natural sediment transport, especially if the focus is on particle interaction or impacts (Maniatis et al.2013). In addition, the physics of IMU sensors are complex, and a number of common assumptions about their use do not always hold. For example, while dead reckoning appears to allow positions to be recovered by double integration of linear accelerations, uncertainties introduced during the production of IMUs (mostly nanometre-scale imperfections on the alignment of the MEMS) lead to extreme uncertainty in positional estimates. A second issue involves calibrating IMU accelerometers, which has often been done using free-fall drop experiments. An accelerometer in free fall will measure zero acceleration despite being subjected to the acceleration of gravity, as gravity in the context of the body frame of the accelerometer is a so-called fictitious force (Appendices A and B). Consequently, the force or impact results of a free-fall drop experiment which relies solely on IMU measurements are highly dependent on how quickly the sensor is programmed to enter and wake up from the free-fall detection state (Clifford2006). It is possible to approximate the height of the free fall using the approximate time of the free-fall state. However, the measurement of the impact force needs a very detailed description of both the impact surface and the low-level code that controls all the basic operations of the sensors (on–off routines, logging, storage handling etc.). This low-level programming is a black box for proprietary off-the-shelf sensors and for users without suitable programming skills. Finally, it is not possible to derive directional information, even for forces, for long mobile periods without complementary corrections (Kok et al.2017) or without a detailed presentation of the reference frames involved and their initial alignment.

Here, we calibrated and deployed a commercial IMU sensor following standard procedures (Maniatis2016), but the precise corrections used are sensor specific, and similar procedures should be followed again for any other IMU sensor. The calibration of force measurements is likely to be standardised and simplified in the near future as the use of IMU sensors develops further. Similar standardisation for the direction of forces is potentially further away as it requires using IMU sensors that rely on optical technology and which are currently not manufactured with physical dimensions or within a price range that is accessible for sediment transport studies (De Agostino et al.2010).

5.2 Relationship with previous work and first-order statistical generalisation

Two aspects of this study are particularly important to address before we make comparisons with previous studies. The first is that we made inertial measurements from within the sediment particles, which are fundamentally different from measurements of fluid turbulence that are often used for predicting sediment motion. The second is that the flow regimes under which we made measurements, with varying shallow flows, differ from those in many studies of sediment motion. Both of these aspects provide new insights into sediment movement, but they require care in making direct comparisons with studies that have used different approaches and/or hydraulic conditions. In addition, it is useful to note that grain protrusion is not discussed in this work, despite being an important control on grain motion and particularly entrainment (e.g. Dey and Ali2018), since the presented laboratory and field experiments only correspond to particles that are highly exposed to the turbulent flow.

This work uses a theoretical framework which has the potential to enhance the mathematical modelling of sediment transport. The Newton–Euler model of Sect. 2, in conjunction with the quaternion transformations of Appendix A, can be read as a 3D and unrestricted Lagrangian–Eulerian model for sediment transport. In our analysis, particle dynamics are transformed from a Lagrangian domain (and the mobile body frame of the particle b) to a static Eulerian domain (frame r), which is most commonly used for the analysis of turbulent flow. Ballio et al. (2018) analyse the topic in detail and provide a comprehensive 1D Lagrangian–Eulerian model which also accounts for the intermittency of sediment transport using a binary classification of mobile and non-mobile states. Our presentation can be used to define 3D Lagrangian dynamics, including rotation, in full and then to transform the corresponding kinematic properties to the Eulerian domain for direct comparison with the turbulent forces. We acknowledge that the verification of the 3D Lagrangian–Eulerian model is heavily dependent on the inertial measurements, and particularly the constant tracking of relative orientation between the frames. However, it is possible to predict that future calibration experiments deploying IMUs will be used this way to parametrise simulations.

Previous laboratory studies using fixed force meters attached to grains (Schmeeckle et al.2007; Cameron et al.2019) report nominal drag and lift forces of the order of 10−1 N for 0.008–0.025 m diameter particles. The measured force normalised by the submerged weight of the average-diameter particle is of the order of 0.08. In the flume results presented here, the inertial drag and lift forces during entrainment are of the order of 101 N, with normalised values for the sphere and ellipsoid being 15.6 and 17.1, respectively. Differences in force magnitude are expected since (a) the inertial sensor is freely mobile, enabling the inertia of the moving particle to be recorded, and (b) the hydraulic regimes in both Schmeeckle et al. (2007) and Cameron et al. (2019) correspond to flows deeper than the particle diameters, with fully developed boundary layers.

Static vibration sensors were also deployed by Lamb et al. (2017), who attached them to a wide range of test particles (D= 0.075–0.218 m). They reported drag forces of the order of 101–102 N and lift forces of the order of 101 N (drag forces from 5.4 to 40.7 when normalised by the submerged weight of the average-diameter particle). Also, the hydraulic regime used by Lamb et al. (2017) is comparable to both the laboratory and the field experiments presented here. The two studies yield very similar forces at entrainment, using different measurement methods, thereby validating the results from our inertial sensor. However, in general it is important to consider the type of sensor (static or restricted vs. mobile), the data processing model and the experimental protocol (varied or steady flow) when different force measurements are compared.

Lamb et al. (2017) also observe a predominance of negative lift forces, especially for partially submerged particles, which may have significant morphological implications, such as potentially explaining, along with turbulence intensities, lower-than-predicted sediment fluxes in steep mountain streams. In our work, the inertial negative lift forces are measured (Appendix E), but the exceedance events (fL>-FGL) are only calculated for the positive lift forces (FLnet>0). The negative lift forces are components of the resultant force, which can have a strong hydraulic component, as argued in Lamb et al. (2017), but they can also be a reaction to positive lift forces during the motion of the particle (and especially the motion of the ellipsoid, Figs. 3 and E2), which requires further investigation.

The laboratory inertial impulse calculations demonstrate that, for unrestricted entrainments, there are observable differences between spherical and ellipsoid particles, with the latter being more sensitive to the lift forces at entrainment threshold conditions. Those differences support previous work on the effect of shape on the response of particles in various hydraulic regimes (e.g. Komar and Li1986; Demir2000) and on the mode of near-bed transport. Measurements of the effect of particle shape can now be made directly using inertial sensors.

The corresponding inertial calculations from the field also demonstrate that the ellipsoid is highly sensitive to lift forces and impulses. We observe higher mean magnitude lift forces on the ellipsoid in the field (2.57 N) than in the laboratory (2.13 N) (Fig. E1). The greater negative lift forces during motion suggest that the particle has a reaction to the positive lift forces, specifically those that exceed the threshold of motion and lead to transport. A similar increase in magnitude was observed for the drag components, with the instantaneous forces being up to 10 times the mean (Fig. E1). The magnitude and duration of the exceedance events in the field are comparable to the observations in the laboratory experiments (Figs. 4 and 6; the relationship between impulse and duration of events is more variable for the field experiments). However, the extreme forces are short lived and so generate very small impulses close to 0 N s−1. Differences in force magnitude and duration can relate to transitions, as described by Shih and Diplas (2018), from hydraulic “impulse controlled” transport, as in our flume experiments, to “force-magnitude controlled” transport, corresponding to the dynamics recorded in Erlenbach. However, an important difference between laboratory and field experiments lies in the scales of turbulence (Coleman and Nikora2008; Singh et al.2009), which requires further investigation since detailed flow measurements were not made during the presented experiments (e.g. PIV measurements).

Overall, differences of particle inertial dynamics during grain entrainment and translation are important because they can potentially enhance predictions for grain particle travel distances with measurements from the field and particularly for large distances (Hassan et al.1991, 2013). Measuring those differences is the most direct insight we can have for studying the effect of several morphological controls (e.g. grain arrangement, burial depths, sediment sorting) until the high-frequency 3D measurement of tracer positions during transport becomes possible.

Considerable effort has been applied to define distributions for hydraulic impulses during the entrainment of spherical particles and relate them to critical thresholds (Diplas et al.2008; Valyrakis et al.2010; Celik et al.2010; Valyrakis et al.2010, 2011). In addition, there have been recent efforts towards upscaling the effect of hydraulic impulses to fully developed bedload equations (Shih and Diplas2018) and results pointing towards evaluating the morphological impact of different hydraulic impulse regimes (Phillips et al.2018). These works highlight the importance of deriving general statistical descriptions for grain inertial impulses.

Our data provide new insights into the roles of drag and lift impulses in entrainment. To begin to generalise these findings and to assess the interactions between lift and drag forces, a bootstrapping method is used here. We approximate the distributions of inertial lift and drag impulses for an ellipsoid particle and from a combination of laboratory and field measurements. This analysis is also the first step towards calculating the combined behaviour of the drag and lift distributions, which can lead to the definition of joint distributions that have stronger explanatory and perhaps predictive value.

To combine the results from the laboratory and field ellipsoid experiments, the impulsive exceedance events were normalised since the conditions are different for the laboratory and the natural conditions. Normalisation used the mean impulses for all drag and all lift forces. Also, the mean impulses were calculated separately for the laboratory and field experiments (I^=I/Iexp). After the normalisation, the laboratory and field results are combined into one dataset of normalised drag (I^Drag) and normalised lift (I^Lift) impulses, which are assumed to be uncorrelated. The latter is justified by the fact that the point lift and drag forces are statistically uncorrected as shown in Appendix E (Fig. E2).

Figure 7Bootstrap normalised impulse sampling (lift and drag). I^Drag and I^Lift are fitted with a gamma and a lognormal distribution, respectively (Appendix F2). The normalised drag and lift thresholds (red dashed lines) are calculated using the probabilistic drag and lift inertial impulse thresholds presented in Fig. 5, which were divided by the recorded mean inertial drag and lift impulse from the laboratory experiments.


The fitting of the representative distributions for I^Drag and I^Lift permits bootstrap sampling from those distributions. Figure 7 shows 50 000 random I^Drag and I^Lift combinations, sampled from the selected gamma and lognormal distributions (for I^Drag and I^Lift, respectively; Appendix F). After taking into account the normalised drag and lift impulse thresholds as defined in the Results (Fig. 5), we conclude that the probability for the exceedance of the lift threshold is approximately 0.02, the probability for the exceedance of the drag threshold is approximately 0.005 and the probability of both thresholds being exceeded simultaneously is approximately 0.0001. The calculation confirms the observation that the transport of the ellipsoid particle is defined by states of zero or very low mobility (97.8 % for the calculated combinations corresponds to dynamics that are below the normalised probabilistic impulse threshold for entrainment, Fig. 5). In total 80 % of the 1371 threshold exceedance events correspond to lift threshold exceedances, 19 % to drag threshold exceedances and 0.4 % to exceedances of both thresholds. The calculation suggests that the majority of the mobility states of the ellipsoid will relate to the action of lift forces. The very small probability for the simultaneous exceedance of both thresholds is another possible effect of the particle's shape as spherical particles will protrude more and are more likely to be equally affected by both drag and lift components (Fig. 5a). This type of calculations requires sample sizes that only advanced instrumentation, such as that presented in this work, can deliver. Similar frameworks can be used for meta-analysis of existing results and to inform the design of future experiments.

6 Conclusions

This work introduces a framework that can be used to derive and interpret IMU measurements in sediment transport studies. The derivation of inertial measurements from mobile sediment grains requires a physical model that links the inertial dynamics with existing force (or moment) balance equations for sediment transport. The types of sensors and associated smart-pebble assemblies currently deployed for the measurements of grain inertial dynamics are not suitable for 2D or 3D tracking of grain position. However, it is possible to measure net forces and impulses if the necessary transformations are applied consistently.

Field and laboratory measurements of inertial lift and drag impulses highlight the different entrainment behaviours of a spherical and an ellipsoidal particle. The lift net force is dominant during the unrestricted entrainment of the ellipsoid, while there is no statistical difference between the effects of lift and drag impulses on the entrainment of the sphere. The drag component can be stronger during transport; however short impulses influence the motion of the ellipsoid significantly.

The continuous improvement of the sensor technology along with the better understanding of the physics described by inertial measurements can lead to a unified treatment of the resultant grain dynamics during bedload transport. These are the dynamics that represent exactly the interaction of hydraulic and sediment forces in different regimes and can enhance the parametrisation of important hydro-morphological controls.

Appendix A: Frames of reference, rotations and IMU measurements

To discuss the measurements recorded by an IMU, and particularly the measurements from an accelerometer and a gyroscope, it is necessary to introduce three basic frames of reference and select one of the many representations for arbitrary rotations in 3D.

The following assumptions and simplifications are used throughout this study:

  • Due to the small-scale (10−1 to 101 m, typically) motion of sediment grains, an Earth frame (one that coincides with the inertial frame as defined below but rotates with the Earth) is not defined. Also, the angular velocity of the Earth (approximately 7.29×10-5 rad s−1) is ignored.

  • For the same reason, the non-gravitational fictitious forces (such as the Coriolis effect) are ignored.

  • For the mathematical derivations, ideal IMUs (no error accumulation is considered) and perfectly aligned sensor assemblies are assumed. The errors associated with IMUs and especially with the magnitude of the integration errors are presented in relevant electrical engineering sources (e.g. Kok et al.2017).

We define the body frame b as the coordinate frame of the moving IMU. For an ideal IMU the origin of this frame is located exactly at the centre of both the accelerometer and the gyroscope, and this centre falls precisely on the centre of mass of the complete sensor assembly (Maniatis et al.2013, 2017). Frame b is mobile but fixed in the sensor assembly (fixed axes representation).

The local geographical frame r is the stationary frame within which hydrodynamics are analysed. This is the reference frame used implicitly for all the single-grain motion studies (Dey and Ali2018). For laboratory experiments, the rxry plane is exactly parallel to the flume bed, and the rz direction is normal to the bed. For the field experiments this alignment will be an approximation due to variations of the local topography. The inertial frame i is a stationary frame. Strap-down IMUs measure acceleration and angular velocity changes in response to this frame, and its origin lies at the centre of the Earth.

Figure A1Frames of reference. r stands for riverbed (flume in this case) reference frame, b for body frame and i for the inertial reference frame.


Transforming information between these three reference frames is non-trivial. A widely used method to represent the change between frame is the application of quaternions (Hamilton1844; Diebel2006). Quaternions are an extension of complex numbers used in the description of 3D mechanics, particularly 3D rotations. They are considered the most efficient description of unrestricted 3D rotations, as they are free from numerical errors that occur when other representations are used (such as the Gimbal lock error associated with rotations expressed by Euler angles, Appendix G2). A typical introduction to quaternions can be found in Valenti et al. (2015), and we follow that primer for a brief introduction to quaternion algebra in Appendix G1.

A unit quaternion ABq defines a rotation from frame A to frame B, and successive rotations are represented by quaternion multiplication. For each ABq, a direction cosine matrix (DCM) R(ABq) is defined as a function of ABq components (Eq. G11), which also represents a rotation from frame A to B. If Bv and Av are observations of the vector v in frames B and A, respectively, they are related through the following typical matrix operation:

(A1) B v = R A B q A v .

If the frames A and B are relatively static (such as the inertial frame i in relation to the local geographic frame r), then both ABq and R(ABq) are explicit. If B is rotating in relation to A (such as the body frame b in relation to the inertial frame i), ABq and the corresponding R(ABq) need to be recursively updated. The transition quaternion q̃ between two successive poses is defined by the applied angular velocity as

(A2) q ̃ = cos ω δ t 2 sin ω δ t 2 ω b x ω sin ω δ t 2 ω b y ω sin ω δ t 2 ω b z ω T ,

where ωbx, ωby and ωbz are angular velocities observed along the bx, by and bz body frame axes, respectively, by the 3D gyroscope; ω=ωb=ωbx2+ωby2+ωbz2 is the norm of angular velocities; and δt is the time of rotation, set here equal to the frequency of the IMU measurements.

Equation (A2) is part of the direct multiplication method (Whitmore2000; Zhao and van Wachem2013), and the updated quaternion ABq is derived as

(A3) A B q = A B q q ̃ ,

with the operation  denoting quaternion multiplication (Eq. G4). After each update ABq is set as ABq.

Inertial accelerometers measure the proper acceleration ab applied within the body frame b. These accelerations will include gravitational acceleration, a uniform force in the inertial frame i. To derive the linear acceleration in the inertial frame i, it is necessary to rotate the body frame measurement to the inertial frame and then to subtract gravitational acceleration. For this rotation the recursively updated R(biq) DCM is used after calculating the biq using Eq. (A3). The linear acceleration in the local geographical frame r is then given by

(A4) a r = R i r q R b i q a b - g ,

where ab=[abxabyabz]T is the vector of the b frame accelerometer measurements, and R(irq) the explicit DCM that rotates the accelerations from i to r derived from the quaternion irq. ar defines the magnitude and direction of the resultant force or net force.

A1 Accelerometers as orientation sensors

In Eq. (A4), the raw accelerometer measurement is “compensated” for the action of the gravitational field in order to extract a measurement for the resultant force in the r frame. To explain the need for this calculation, it is important to note that accelerometers measure proper acceleration. Proper acceleration is different from coordinate acceleration, which is generally defined as the rate of velocity (in a fixed inertial frame). In practice, if an accelerometer is placed on a flat surface, it will measure an upwards (positive) acceleration equal to 1 g (9.81 m s−2). In free fall, an accelerometer will measure zero acceleration because in the non-inertial frame of the sensor there is neither a static (e.g. gravity) nor a dynamic (motion or vibrational) force applied and the accelerometer will only “feel” the terminal velocity of the impact if it lands on a surface. When the accelerometer is subjected to an external force and it begins to move (in relation to a static inertial frame), there is no way to separate the components of the static and the dynamic forces, unless the relative orientation of the frames is constantly monitored (in the manner we describe in Eq. A4). At the same time, when the sensor is static, the non-compensated signal (the raw acceleration from the sensor) can provide an estimate for the relative orientation of the sensor's frame (bi or br quaternions in relation to the frame to which static forces are applied (gravity frame, i)). Valenti et al. (2015) provide the solution for calculating the ibq quaternion from raw accelerometer measurements, which, adapted to the notation used here, is as follows:

(A5) i b q = a b z + 1 2 - a b y 2 a b z + 1 a b x 2 a b z + 1 0 T , a b z 0 - a b y 2 1 - a b z 1 - a b z 2 0 a b x 2 1 - a b z T , a b z < 0 .

In this work we used Eq. (A5) to estimate the initial alignment of the frames for both the flume and the field experiments. In addition, the YEI sensor implements an onboard calibration routine which can verify the initial alignment using sensor fusion and a series of nonlinear filters (YEI2014).

Appendix B: 3D IMU measurements

Figure B1Example incipient motion IMU data. Measurements from the incipient motion experiments using the spherical sensor. (a) Unfiltered and uncompensated inertial acceleration measurements. The sensor is initially aligned to gravity, which results in the z axis of the accelerometer measuring a mean value of 9.81 m s−2. (b) Angular velocity (rad s−1) measurements derived from the gyroscope. (c) Linear acceleration along the three-body-frame axis. This is the result of removing gravity from the inertial measurements shown in (a) and applying a FFT high-pass filter as described in Maniatis (2016) (Sect. 6). Panel (d) shows the kinetic energy calculation after integrating once the signal presented in panel (c) and applying the formula K=12mvr2+12Icmωb2 as described in Sect. 2.


Appendix C: The magnitude of the tangential forces and resultant torques

The rotational component captured by the gyroscope relates to the moments of the forces applied on the surface of the measured particle via Eq. (2). For a sphere, frot is defined by rearranging Eq. (11) as

(C1) f rots = I cm d ω b d t R ,

where R is the radius of the particle. The relationship holds because the moment of inertia (Icm) of the sphere is uniform. This implies that the shape of the particle does not affect the direction of the resultant force. For the spherical sensor presented in this work the moment of inertia is calculated using the formula for a sphere (Icm=52mR2) and is equal to 0.00085 kg m−2.

For the ellipsoid, the same calculation is significantly more complicated. Firstly, the moment of inertia is not uniform. In this work, we implemented a numerical calculation during the design phase of the enclosure (using Solidworks, Système2016) for the principle axes of the ellipsoid. The principle axes coincide with the bx, by and bz body frame axes. Those are the equivalent of the axes a, b and c of the ellipsoid (where a=0.1 m, b=0.07 m and c=0.03 m), respectively; they are fixed in the body frame and aligned with the IMU axes (see Appendix A). The principle components of inertia for the ellipsoid sensor were calculated as Ixx=0.00057 kg m−2, Iyy=0.00060 kg m−2 and Izz=0.00094 kg m−2. The non-principle components of inertia were calculated as being of the order of 10−8 kg m−2, and they were ignored. The balance of moments for the principle axes is given by the following system of equations (often discussed as Euler equations, O'Reilly2008):


For the calculation of the tangential force, a good approximation (ignoring the secondary moments of inertia) is given by dividing the principle moments by the half length of the corresponding principle axis of the ellipsoid, giving

(C5) f rote = 2 M x a 2 + 2 M y b 2 + 2 M z c 2 .

Equations (C1) and (C5) were used to calculate the magnitude of the rotational component for the example experiments presented in Fig. 3. Figure C1 shows the same calculations at a scale that reveals their fluctuation.

Figure C1Force magnitude of the rotational component. (a) Flume experiment (incipient motion) using the spherical sensor (Fig. 3a). (b) Field experiment (Erlenbach) using the ellipsoid sensor (Fig. 3e). The differences between the linear and the rotational components are between 1 and 2 orders of magnitude. frot is of the order of 10−3 to 10−2 N (sphere in the flume and ellipsoid in Erlenbach, respectively), while the linear forces (applied directly on the centre of the mass of the particle) are of the order of 100 to 101 N (Appendix E).


Appendix D: Hydraulic parameters

Table D1The parameters are estimated as follows: ρpρf is the ratio of an experimental particle density to fluid density (ρf=1000 kg m−3), Q is the flow rate at the point of entrainment for the flume experiments and during the transport events in Erlenbach, H is the flow depth (measured from the bottom of the bed to the water surface for Erlenbach and from the top of the hemisphere bed for the flume experiments), Ub=Q/A is the bulk mean velocity (A is the cross-sectional area of the flow), W is the channel width and Sb is bed slope. The value 0.105 (or 0.1) is also the average bed slope of the lowermost natural reach in Erlenbach of about 30 m length upstream of the stream gauging station. F=Ub/(gd)0.5 is the Froude number; Rb=(Ubd)/ν is the bulk Reynolds number (where ν is the fluid kinematic viscosity at 0 C for Erlenbach and at 25 C for the flume experiments). For the calculations referring to the ellipsoid, we assume D equal to the particle's b axis.

Download Print Version | Download XLSX

Appendix E: Summary of statistics for inertial lift and drag forces

Figure E1Histogram of inertial forces from all experiments. The inertial dynamics show that net lift (FLnet) and drag (FDnet) forces consistently fluctuate around 0. The vertical lines indicate the corresponding medians. The mean force magnitude for the sphere (in the lab) is FDnet=0.57 N and FLnet=0.62 N for the drag and lift directions, respectively. For the ellipsoid in the laboratory experiment the mean drag force magnitude is FDnet=0.74 N and the mean lift force magnitude is FLnet=2.13 N. Finally for the ellipsoid in Erlenbach mean FDnet=1.22 N and FLnet=2.57 N.


Figure E2Lift vs. drag force magnitude correlation (flume experiments). Regression analysis applied to the magnitude of calculated forces (drag and lift) shows a moderate correlation for the spherical particle (statistically significant Pearson's R=0.3) and a weak correlation for the ellipsoid in both the laboratory and the field experiments – statistically significant Pearson's R=0.046 for the lab measurements, ellipsoid (Lab), and insignificant Pearson's R=-0.00055 for the field ones, ellipsoid (Erlenbach). The latter supports the assumption of statistical independence between the two components for the ellipsoid, justifying the randomisation presented in Sect. 5 (Fig. 7). Pearson's R is an unbiased metric for this sample size.


Appendix F: Normalised impulse: selection of representative drag and lift distributions

Three types of right-tail distributions were considered (Cullen et al.1999; Appendix F2) as well-fitting candidates for both I^Drag and I^Lift: the Weibull, the gamma and the lognormal distributions. Goodness-of-fit analysis (Table F1) shows that I^Drag is approximated marginally better by a gamma distribution (median shape = 0.529, median scale = 0.5) and the I^Lift is approximated better by a lognormal distribution (median meanlog =−0.66, median sdlog = 1.13).

Figure F1a shows the Cullen and Frey diagram for the identification of candidate distributions for the normalised drag impulses (I^Drag). Being in the “beta” region, the diagram indicates a right-tail distribution (beta distributions are restricted between 0 and 1, and there is no physical relationship with the possible values for impulses). The skewness–kurtosis relationship (blue dots) indicates a right-tail distribution as a candidate as well. Figure F1b–e show the graphical comparison between three candidate distributions (Weibull, gamma and lognormal). Weibull and gamma distributions outperform the lognormal distribution on the tails of the histogram (Q–Q plot). The median values are also captured better by the Weibull and gamma distributions (P–P plot). Finally, the histogram and CDF diagrams confirm that the lognormal distribution is the least representative of I^Drag. The gamma distribution marginally outperforms the Weibull distribution for I^Drag both in graphical and goodness-of-fit comparisons, and it is selected for the bootstrap calculation of Fig. 7.

Figure F2a shows the Cullen and Frey diagram for the identification of candidate distributions for the normalised lift impulses (I^Lift). The skewness–kurtosis relationship (blue dots) indicates a right-tail distribution as a candidate. Figure F2b–e show the graphical comparison between three candidate distributions (Weibull, gamma and lognormal). The lognormal distribution outperforms the other candidates at the tail (Q–Q plot) and the median regions (P–P plot). Finally, the histogram and CDF diagrams confirm that the lognormal distribution is the best representative of I^Lift, and it is selected for the bootstrap calculation of Fig. 7.

Figure F1Choice of distribution for drag impulses.


Figure F2Choice of distribution for lift impulses.


Table F1Fitted distribution statistics.

Download Print Version | Download XLSX

Table F2Statistics for selected distributions.

Download Print Version | Download XLSX

Figure F3Bootstrap parameters for selected distributions. The graphs quantify the stability of the selected distributions for drag and lift impulses. For the gamma distribution (drag impulses) 1000 bootstrapped parameters were cross-compared, revealing a range of 0.1 for the shape parameter and 0.2 for the rate parameter (a). For the lognormal distribution (lift impulses) 1000 bootstrapped parameters were cross-compared, revealing a range of 0.15 for the log-mean parameter and 0.14 for log standard deviation parameter (b). The differences are marginal, indicating good stability of the selected distributions for the scaling of the data.


Appendix G: Quaternions and rotations

G1 Summary of quaternion algebra

Quaternions can be written in the form

(G1) q = q 1 + q 2 i + q 3 j + q 4 k ,

where q1q4 are the components of quaternion q (and ik and j are unit imaginary numbers).

The quaternion conjugate is given by

(G2) q = q 1 - q 2 i - q 3 j - q 4 k .

The sum of two quaternions is then

(G3) q + w = q 1 + w 1 + q 2 + w 2 i + q 3 + w 3 j + q 4 + w 4 k ,

and quaternion multiplication is defined as

(G4) q w = q 1 w 1 - q 2 w 2 - q 3 w 3 - q 4 w 4 + q 1 w 2 + q 2 w 1 + q 3 w 4 - q 4 w 3 i + q 1 w 3 - q 2 w 4 + q 3 w 1 + q 4 w 2 j + q 1 w 4 + q 2 w 3 - q 3 w 2 + q 4 w 1 k .

The quaternion norm is therefore defined by

(G5) n ( q ) = q q = q 2 2 + q 1 2 + q 3 2 + q 4 2 .

With little manipulation, the quaternions can be directly related to four-element vectors.

Quaternions can be interpreted as a scalar plus a vector by writing

(G6) q = q 1 + q 2 i + q 3 j + q 4 k = ( s , v ^ ) ,

where s=q1 and v^=q2i+q3j+q4k. In this notation, quaternion multiplication has the form

(G7) q 1 q 2 = s 1 , v ^ 1 s 2 , v ^ 2 = s 1 s 2 - v ^ 1 v ^ 2 , s 1 v ^ 2 + s 2 v ^ 1 + v ^ 1 v ^ 2 .

Finally, the rotation about the unit vector n^ by an angle θ can be computed using the quaternion

(G8) q = ( s , v ) = cos 1 2 θ , n ^ sin 1 2 θ .

The components of this quaternion are called Euler parameters. After rotation, a point p=(0p) is then given by

(G9) p = q p q - 1 = q p q

since n(q)=1.

A concatenation of two rotations, first q1 and then q2, can be computed using the identity

(G10) q 2 q 1 p q 1 q 2 = q 2 q 1 p q 1 q 2 = q 2 q 1 p q 2 q 1 .

Finally, the transformation that gives the equivalent DCM for a quaternion q=q1+q2i+q3j+q4k is given by

(G11) R ( q ) = q 1 2 + q 2 2 - q 3 2 - q 4 2 , 2 q 2 q 3 - q 4 q 1 , 2 q 2 q 4 + q 3 q 1 2 q 2 q 3 + q 4 q 1 , q 1 2 - q 2 2 + q 3 2 - q 4 2 , 2 q 3 q 4 - q 2 q 1 2 q 2 q 4 - q 3 q 1 , 2 q 3 q 4 + q 2 q 1 , q 1 2 - q 2 2 - q 3 2 + q 4 2 .

G2 The Gimbal lock (adapted from Maniatis, 2016)

To demonstrate the advantage of quaternions, we randomly rotate the static vector of gravity. In an orthogonal Cartesian frame where the origin of the z axis is the centre of the Earth, gravity is measured as [Gx, Gy, Gz]=[0, 0, 9.81] m s−2. If we assume a rigid body rotating freely and randomly in this frame, we can do the rotation calculations. Avoiding further mathematisation, the series of the calculations is the following:

  • Randomisation of the body frame angular velocities of the rigid body ωx, ωy and ωz in a [−2π2π] range. A frequency of 100 Hz is used.

  • Calculation of successive quaternions using direct multiplication for random angular velocities.

  • Calculation of the direction cosine matrix from Euler angles.

  • Rotating the vector of gravity in the body frame of the rigid body using both the direction cosine matrix and common matrix vector multiplication.

The vector expressed in the body frame is shown in Fig. G1. The results are different, and Gimbal lock (an inconsistent axis change when the second rotation approaches ±π/2) occurs after the 450th iteration, which corresponds to 8 s in simulation time.

Figure G1Random rotation of the static vector of gravity. [gx,gy,gz]=[0,0,9.81] m s−2 in the gravity frame of reference as expressed in the body frame of a randomly rotating rigid body (dt=0.01 s). Panel (a) demonstrates the rotation calculations with the usage of Euler angles. Gimbal lock occurs after the 400 iterations. Panel (b) shows the same rotation series calculated via quaternions. No Gimbal lock occurs and the result is easy to interpret as it is based on the use of the measured body frame angular velocities.


Data availability

The dataset used in this paper is curated by Georgios Maniatis in the following Github repository: (Maniatis2020).

Code availability

All the calculations were performed using the R statistical software and specifically the libraries orientlib (Murdoch2003) for the rotation calculations and fitdistrplus (Delignette-Muller and Dutang2015) for fitting statistical distributions. Both libraries are open and freely available.

Author contributions

GM designed and calibrated the sensor, designed and implemented the experiments, performed all the physical and statistical calculations, and produced the first draft of this paper. TH supervised the laboratory experiments, reviewed several versions of the manuscript, and contributed significantly to the interpretation and contextualisation of the results. RH contributed to the design of the laboratory experiments, reviewed several versions of the manuscript, and contributed to the interpretation and contextualisation of the results. DR contributed to the design, supervised and assisted with the field experiments, reviewed several versions of the manuscript, and contributed to the interpretation and contextualisation of the results. AB contributed to the design of the field experiments and reviewed several versions of the manuscript.

Competing interests

The authors declare that they have no conflict of interest.


The flume experiments were conducted in the School of Engineering of the University of Glasgow. The authors thank Tim Montgomery (University of Glasgow) and Tobias Nicollier (WSL) for their assistance with the laboratory and field experiments, respectively; three anonymous reviewers, who significantly improved the manuscript with their contribution; and Jens Turowski, who provided feedback on previous versions of the manuscript. Finally, Georgios Maniatis thanks Katerina Georgiou for her assistance with the design and production of the figures and the typewriting of the manuscript.

Financial support

Georgios Maniatis was supported by a University of Glasgow Lord Kelvin Adam Smith Scholarship (charity number SC004401) at the time of the flume experiments. The field experiments were supported by an Early Career Research award from the British Society for Geomorphology (Quantification of erosional processes using inertial sensors across environments, 2017).

Review statement

This paper was edited by Claire Masteller and reviewed by three anonymous referees.


Akeila, E., Salcic, Z., and Swain, A.: Smart pebble for monitoring riverbed sediment transport, IEEE Sensors J., 10, 1705–1717, 2010. a

Ali, S. Z. and Dey, S.: Hydrodynamics of sediment threshold, Phys. Fluids, 28, 075103,, 2016. a

Ancey, C., Davison, A., Böhm, T., Jodeau, M., and Frey, P.: Entrainment and motion of coarse particles in a shallow water stream down a steep slope, J. Fluid Mech., 595, 83–114, 2008. a

Ashida, K. and Michiue, M.: An investigation of river bed degradation downstream of a dam, in: in Proceedings of 14th Int. Association of Hydraulic Research Congress, vol. 3, Wallingford, UK, 247–255, 1971. a

Ballio, F., Pokrajac, D., Radice, A., and Hosseini Sadabadi, S. A.: Lagrangian and Eulerian description of bed load transport, J. Geophys. Res.-Earth, 123, 384–408, 2018. a, b

Begin, Z. and Schumm, S.: Instability of alluvial valley floors: a method for its assessment, T. ASAE, 22, 347–350, 1979. a

Bialik, R. J., Nikora, V. I., Karpiński, M., and Rowiński, P. M.: Diffusion of bedload particles in open-channel flows: distribution of travel times and second-order statistics of particle trajectories, Environ. Fluid Mech., 15, 1281–1292, 2015. a, b

Buffington, J. M. and Montgomery, D. R.: A systematic analysis of eight decades of incipient motion studies, with special reference to gravel-bedded rivers, Water Resour. Res., 33, 1993–2029, 1997. a

Buffington, J. M., Dietrich, W. E., and Kirchner, J. W.: Friction angle measurements on a naturally formed gravel streambed: implications for critical boundary shear stress, Water Resour. Res., 28, 411–425, 1992. a

Bunte, K., Abt, S. R., Potyondy, J. P., and Ryan, S. E.: Measurement of coarse gravel and cobble transport using portable bedload traps, J. Hydraul. Eng., 130, 879–893, 2004. a

Burtin, A., Hovius, N., McArdell, B. W., Turowski, J. M., and Vergne, J.: Seismic constraints on dynamic links between geomorphic processes and routing of sediment in a steep mountain catchment, Earth Surf. Dynam., 2, 21–33,, 2014. a

Cameron, S., Nikora, V., and Marusic, I.: Drag forces on a bed particle in open-channel flow: Effects of pressure spatial fluctuations and very-large-scale motions, J. Fluid Mech., 863, 494–512, 2019. a, b

Celik, A. O., Diplas, P., Dancey, C. L., and Valyrakis, M.: Impulse and particle dislodgement under turbulent flow conditions, Phys. Fluids, 22, 046601,, 2010. a, b, c, d

Clifford, M.: Detecting Freefall with Low-G Accelerometers, Sensor and Analog Products Division, Tempe, AZ, 2006. a

Coleman, S. and Nikora, V.: A unifying framework for particle entrainment, Water Resour. Res., 44, W04415,, 2008. a, b

Cullen, A. C., Frey, H. C., and Frey, C. H.: Probabilistic techniques in exposure assessment: a handbook for dealing with variability and uncertainty in models and inputs, Springer Science & Business Media, USA, 1999. a

De Agostino, M., Manzino, A. M., and Piras, M.: Performances comparison of different MEMS-based IMUs, in: IEEE/ION Position, Location and Navigation Symposium, 4–6 May 2010, Palm Springs, California, 187–201, 2010.  a

Delignette-Muller, M. L. and Dutang, C.: fitdistrplus: An R package for fitting distributions, J. Stat. Softw., 64, 1–34, 2015. a

Demir, T.: The influence of particle shape on bedload transport in coarse-bed river channels, PhD thesis, Durham University, Durham, 2000. a

Dey, S.: Sediment threshold, Appl. Math. Model., 23, 399–417, 1999. a

Dey, S. and Ali, S. Z.: Advances in modeling of bed particle entrainment sheared by turbulent flow, Phys. Fluids, 30, 061301,, 2018. a, b, c, d

Diebel, J.: Representing attitude: Euler angles, unit quaternions, and rotation vectors, Matrix, 58, 1–35, 2006. a

Diplas, P., Dancey, C. L., Celik, A. O., Valyrakis, M., Greer, K., and Akar, T.: The role of impulse on the initiation of particle movement under turbulent flow conditions, Science, 322, 717–720, 2008. a, b, c, d, e

Drake, T. G., Shreve, R. L., Dietrich, W. E., Whiting, P. J., and Leopold, L. B.: Bedload transport of fine gravel observed by motion-picture photography, J. Fluid Mech., 192, 193–217, 1988. a

Einstein, H. A.: Bedload Transport as a Probability Problem, Mitteilung der Versuchsanstalt für Wasserbau an der Eidgenössischen Technischen Hochschule, Tech. Hochsch., Zurich, Switzerland, p. 110, 1937. a, b

Ergenzinger, P. and Jupner, R.: Using COSSY (CObble Satellite SYstem) for measuring the effects of lift and drag forces, Erosion and Sediment Transport Monitoring Programmes in river Basins, IAHS Publications, Oslo, 41–50, 1992. a

Fathel, S., Furbish, D., and Schmeeckle, M.: Parsing anomalous versus normal diffusive behavior of bed load sediment particles, Earth Surf. Proc. Land., 41, 1797–1803,, 2016. a

Ferguson, R. I., Bloomer, D. J., Hoey, T. B., and Werritty, A.: Mobility of river tracer pebbles over different timescales, Water Resour. Res., 38, 1045,, 2002. a

Frank, D., Foster, D., Sou, I. M., Calantoni, J., and Chou, P.: Lagrangian measurements of incipient motion in oscillatory flows, J. Geophys. Res.-Oceans, 120, 244–256, 2015. a

Garcia, C., Cohen, H., Reid, I., Rovira, A., Úbeda, X., and Laronne, J. B.: Processes of initiation of motion leading to bedload transport in gravel-bed rivers, Geophys. Res. Lett., 34, L06403,, 2007. a

Gebre-Egziabher, D., Hayward, R. C., and Powell, J. D.: A low-cost GPS/inertial attitude heading reference system (AHRS) for general aviation applications, in: Position Location and Navigation Symposium, IEEE 1998, 20–23 April 1996, Palm Springs, CA, USA, 518–525, 1998. a

Gilbert, G. K. and Murphy, E. C.: The transportation of debris by running water, US Government Printing Office, Washington, D.C., USA, 86 pp., 1914. a, b, c

Gimbert, F., Fuller, B. M., Lamb, M. P., Tsai, V. C., and Johnson, J. P.: Particle transport mechanics and induced seismic noise in steep flume experiments with accelerometer-embedded tracers, Earth Surf. Proc. Land., 44, 219–241,, 2019. a, b

Grass, A. J.: Initial instability of fine bed sand, J. Hydraul. Div., 96, 619–632, 1970. a, b, c

Grewal, M. S., Weill, L. R., and Andrews, A. P.: Global positioning systems, inertial navigation, and integration, John Wiley & Sons, USA, 2007. a

Gronz, O., Hiller, P. H., Wirtz, S., Becker, K., Iserloh, T., Seeger, M., Brings, C., Aberle, J., Casper, M. C., and Ries, J. B.: Smartstones: A small 9-axis sensor implanted in stones to track their movements, Catena, 142, 245–251, 2016. a, b

Hamilton, W. R.: II. On quaternions; or on a new system of imaginaries in algebra, London Edinburgh Dublin Philos. Mag. J. Sci., 25, 10–13, 1844. a

Hassan, M. A. and Roy, A. G.: Coarse particle tracing in fluvial geomorphology, in: Tools in Fluvial Geomorphology, edited by: Kondolf, G. M. and Piégay, H., John Wiley & Sons, Ltd, UK,, 2016. a

Hassan, M. A., Church, M., and Schick, A. P.: Distance of movement of coarse particles in gravel bed streams, Water Resour. Res., 27, 503–511, 1991. a

Hassan, M. A., Church, M., and Ashworth, P. J.: Virtual rate and mean distance of travel of individual clasts in gravel-bed channels, Earth Surf. Proc. Land., 17, 617–627, 1992. a

Hassan, M. A., Church, M., Rempel, J., and Enkin, R. J.: Promise, performance and current limitations of a magnetic Bedload Movement Detector, Earth Surf. Proc. Land., 34, 1022–1032, 2009. a, b

Hassan, M. A., Voepel, H., Schumer, R., Parker, G., and Fraccarollo, L.: Displacement characteristics of coarse fluvial bed sediment, J. Geophys. Res.-Earth, 118, 155–165, 2013. a

Hodge, R. A., Hoey, T. B., and Sklar, L. S.: Bed load transport in bedrock rivers: the role of sediment cover in grain entrainment, translation, and deposition, J. Geophys. Res., 116, 1–19, 2011. a, b

Hodge, R. A., Sear, D. A., and Leyland, J.: Spatial variations in surface sediment structure in riffle–pool sequences: a preliminary test of the Differential Sediment Entrainment Hypothesis (DSEH), Earth Surf. Proc. Land., 38, 449–465, 2013. a

Ikeda, S.: Incipient motion of sand particles on side slopes, J. Hydraul. Div., 108, 95–114, 1982. a

Iwagaki, Y.: Basic studies on the critical tractive force (1), Trans. JSCE, 31, 1–20, 1956. a

Johnson, J. P. L.: Gravel threshold of motion: a state function of sediment transport disequilibrium?, Earth Surf. Dynam., 4, 685–703,, 2016. a

Johnson, M. F., Rice, S. P., and Reid, I.: Increase in coarse sediment transport associated with disturbance of gravel river beds by signal crayfish (Pacifastacus Leniusculus), Earth Surf. Proc. Land., 36, 1680–1692, 2011. a

Kirchner, J. W., Dietrich, W. E., Iseya, F., and Ikeda, H.: The variability of critical shear stress, friction angle, and grain protrusion in water-worked sediments, Sedimentology, 37, 647–672,, 1990. a

Kline, S., Reynolds, W., Schraub, F., and Runstadler, P.: The structure of turbulent boundary layers, J. Fluid Mech., 30, 741–773,, 1967. a

Kok, M., Hol, J. D., and Schön, T. B.: Using inertial sensors for position and orientation estimation, arXiv preprint: arXiv:1704.06053, 2017. a, b, c

Komar, P. D. and Li, Z.: Pivoting analyses of the selective entrainment of sediments by shape and size with application to gravel threshold, Sedimentology, 33, 425–436, 1986. a

Komar, P. D. and Li, Z.: Applications of grain-pivoting and sliding analyses to selective entrapment of gravel and to flow-competence evaluations, Sedimentology, 35, 681–695, 1988. a

Kularatna, N., Melville, B., Akeila, E., and Kularatna, D.: Implementation aspects and offline digital signal processing of a smart pebble for river bed sediment transport monitoring, in: 5th IEEE Conference on Sensors, Nashville, Tenesse, USA, 1093–1098, 2006. a

Lamb, M. P., Brun, F., and Fuller, B. M.: Direct measurements of lift and drag on shallowly submerged cobbles in steep streams: Implications for flow resistance and sediment transport, Water Resour. Res., 53, 7607–7629, 2017. a, b, c, d, e

Liedermann, M., Tritthart, M., and Habersack, H.: Particle path characteristics at the large gravel-bed river Danube: results from a tracer study and numerical modelling, Earth Surf. Proc. Land., 38, 512–522, 2012. a

Maniatis, G.: Eulerian-Lagrangian definition of coarse bed-load transport: theory and verification with low-cost inertial measurement units, PhD thesis, University of Glasgow, Glasgow, 2016. a, b, c, d

Maniatis, G.: ESD, Inertial drag and lift forces for coarse grains measured using in-grain accelerometer, Zenodo,, 2020. a

Maniatis, G., Hoey, T., and Sventek, J.: Sensor Enclosures: example Application and Implications for Data Coherence, J. Sensor Actuat. Netw., 2, 761,, 2013. a, b, c, d

Maniatis, G., Hoey, T. B., Hassan, M. A., Sventek, J., Hodge, R., Drysdale, T., and Valyrakis, M.: Calculating the explicit probability of entrainment based on inertial acceleration measurements, J. Hydraul. Eng., 143, 04016097,, 2017. a, b, c

Marion, A. and Tregnaghi, M.: A new theoretical framework to model incipient motion of sediment grains and implications for the use of modern experimental techniques, in: Experimental and Computational Solutions of Hydraulic Problems, Springer, Łochów, Poland , 85–100, 2013. a

Masteller, C. C., Finnegan, N. J., Turowski, J. M., Yager, E. M., and Rickenmann, D.: History-Dependent Threshold for Motion Revealed by Continuous Bedload Transport Measurements in a Steep Mountain Stream, Geophys. Res. Lett., 46, 2583–2591, 2019. a

McEwan, I., Habersack, H., and Heald, J.: Discrete particle modelling and active tracers: new techniques for studying sediment transport as a Lagrangian phenomenon, in: Gravel bed rivers V, edited by: Mosley, M. P., Hydrological Society, Wellington, New Zealand, 339–360, 2001. a

McEwan, I., Sørensen, M., Heald, J., Tait, S., Cunningham, G., Goring, D., and Willetts, B.: Probabilistic modeling of bed-load composition, J. Hydraul. Eng., 130, 129–139, 2004. a, b

Murdoch, D.: Orientlib: An R package for orientation data, J. Stat. Softw., 8, 1–11, 2003. a

Nelson, J. M., Shreve, R. L., McLean, S. R., and Drake, T. G.: Role of near-bed turbulence structure in bed load transport and bed form mechanics, Water Resour. Res., 31, 2071–2086, 1995. a

Niño, Y. and García, M.: Using Lagrangian particle saltation observations for bedload sediment transport modelling, Hydrol. Process., 12, 1197–1218, 1998. a

Olinde, L. and Johnson, J. P.: Using RFID and accelerometer-embedded tracers to measure probabilities of bed load transport, step lengths, and rest times in a mountain stream, Water Resour. Res., 51, 7572–7589, 2015. a

O'Reilly, O.: Intermediate Dynamics for Engineers: A Unified Treatment of Newton–Euler and Lagrangian Mechanics, Cambridge University Press, Cambridge,, 2008. a, b

Papanicolaou, A., Diplas, P., Evaggelopoulos, N., and Fotopoulos, S.: Stochastic incipient motion criterion for spheres under various bed packing conditions, J. Hydraul. Eng., 128, 369–380, 2002. a, b, c, d, e, f

Phillips, C., Hill, K. M., Paola, C., Singer, M., and Jerolmack, D.: Effect of flood hydrograph duration, magnitude, and shape on bed load transport dynamics, Geophys. Res. Lett., 45, 8264–8271, 2018. a, b

Prancevic, J. P. and Lamb, M. P.: Particle friction angles in steep mountain channels, J. Geophys. Res.-Earth, 120, 242–259, 2015. a

Recking, A., Piton, G., Vazquez-Tarrio, D., and Parker, G.: Quantifying the morphological print of bedload transport, Earth Surf. Proc. Land., 41, 809–822,, 2015. a

Rickenmann, D., Turowski, J. M., Fritschi, B., Klaiber, A., and Ludwig, A.: Bedload transport measurements at the Erlenbach stream with geophones and automated basket samplers, Earth Surf. Proc. Land., 37, 1000–1011, 2012. a

Schmeeckle, M. W. and Nelson, J. M.: Direct numerical simulation of bedload transport using a local, dynamic boundary condition, Sedimentology, 50, 279–301, 2003. a

Schmeeckle, M. W., Nelson, J. M., and Shreve, R. L.: Forces on stationary particles in near-bed turbulent flows, J. Geophys. Res.-Earth, 112, F02003,, 2007. a, b, c

Schmidt, K.-H. and Ergenzinger, P.: Bedload entrainment, travel lengths, step lengths, rest periods – studied with passive (iron, magnetic) and active (radio) tracer techniques, Earth Surf. Proc. Land., 17, 147–165, 1992. a

Schneider, J. M., Turowski, J. M., Rickenmann, D., Hegglin, R., Arrigo, S., Mao, L., and Kirchner, J. W.: Scaling relationships between bed load volumes, transport distances, and stream power in steep mountain channels, J. Geophys. Res.-Earth, 119, 533–549, 2014. a

Shields, A.: Application of similarity principles and turbulence research to bed-load movement, Technical Report, Soil Conservation Service, Pasadena, California, 1936. a, b

Shih, W. and Diplas, P.: A unified approach to bed load transport description over a wide range of flow conditions via the use of conditional data treatment, Water Resour. Res., 54, 3490–3509, 2018. a, b

Shvidchenko, A. B. and Pender, G.: Flume study of the effect of relative depth on the incipient motion of coarse uniform sediments, Water Resour. Res., 36, 619–628, 2000. a

Singh, A., Fienberg, K., Jerolmack, D. J., Marr, J., and Foufoula-Georgiou, E.: Experimental evidence for statistical scaling and intermittency in sediment transport rates, J. Geophys. Res.-Earth, 114, F01025,, 2009. a, b

Spazzapan, M., Petrovčič, J., and Mikoš, M.: New tracer for monitoring dynamics of sediment transport in turbulent flows, Acta Hydrotech., 22, 135–148, 2004. a

Système: Dassault Systèmes, SolidWorks Software webpage, available at: (last access: 1 July 2020), 2016. a

Tsakiris, A. G., Papanicolaou, A., Moustakidis, I., and Abban, B. K.: Identification of the Burial Depth of Radio Frequency Identification Transponders in Riverine Applications, J. Hydraul. Eng., 141, 04015007,, 2015. a

Tucker, G. E. and Hancock, G. R.: Modelling landscape evolution, Earth Surf. Proc. Land., 35, 28–50, 2010. a

Turowski, J. M., Badoux, A., and Rickenmann, D.: Start and end of bedload transport in gravel-bed streams, Geophys. Res. Lett., 38, L04401,, 2011. a

Valenti, R. G., Dryanovski, I., and Xiao, J.: Keeping a good attitude: A quaternion-based orientation filter for IMUs and MARGs, Sensors, 15, 19302–19330, 2015. a, b

Valyrakis, M., Diplas, P., Dancey, C. L., Greer, K., and Celik, A. O.: Role of instantaneous force magnitude and duration on particle entrainment, J. Geophys. Res.-Earth, 115, F02006,, 2010. a, b, c, d

Valyrakis, M., Diplas, P., and Dancey, C. L.: Entrainment of coarse grains in turbulent flows: An extreme value theory approach, Water Resour. Res., 47, W09512,, 2011.  a

Van Rijn, L. C.: Sediment transport, part I: bed load transport, J. Hydraul. Eng., 110, 1431–1456, 1984. a

VectorNav: Inertial Measurement Units and Inertial Navigation, VecotrNav webpage, available at:, last access: 9 May 2016. a

Vignaga, E., Sloan, D. M., Luo, X., Haynes, H., Phoenix, V. R., and Sloan, W. T.: Erosion of biofilm-bound fluvial sediments, Nat. Geosci., 6, 770–774, 2013. a

Whitmore, S. A.: Closed-form integrator for the quaternion (Euler angle) kinematics equations, US Patent 6,061,611, 2000. a

Woodman, O. J.: An introduction to inertial navigation, Technical Report UCAMCL-TR-696, 14, University of Cambridge, Computer Laboratory, Cambridge, 2007. a

Yalin, M. S.: An expression for bed-load transportation, J. Hydraul. Div., 89, 221–250, 1963. a

YEI: 3-Space Sensor, User's Manual, YEI Trechnology, available at: (last access: 18 December 2020), 2014. a, b, c

Zekavat, R. and Buehrer, R. M.: Handbook of Position Location: Theory, Practice and Advances, in: vol. 27, John Wiley & Sons, USA, 2011. a

Zhao, F. and van Wachem, B.: A novel Quaternion integration approach for describing the behaviour of non-spherical particles, Acta Mechanica, 224, 3091–3109, 2013. a

Short summary
One of the most interesting problems in geomorphology concerns the conditions that mobilise sediments grains in rivers. Newly developed smart pebbles allow for the measurement of those conditions directly if a suitable framework for analysis is followed. This paper connects such a framework with the physics used to described sediment motion and presents a series of laboratory and field smart-pebble deployments. Those quantify how grain shape affects the motion of coarse sediments in rivers.