5. Velocimetric method
In contrast with the relatively direct character of the
stereoscopic method, the velocimetric technique is indirect in
nature. The technique seeks to reconstruct the surface topography on
the basis of measurements of the horizontal velocity field. For that
purpose, the method exploits the Bernoulli relation applicable to free
surface streamlines. This relation takes a particularly simple form if
the flow can be approximated as a superposition of small perturbations
upon a rapid mean current. To carry out the procedure, images from the
top camera are first used to track particles on the water free surface
and estimate their horizontal velocities. The elevation field is then
obtained after suitable longitudinal and transverse averaging of the
individual velocity vectors. These steps are detailed in the next
paragraphs.
5.1 Elevation-from-velocity principle
The following two assumptions are taken to apply locally (on the
scale of a few wavelengths) to antidune flows: 1) the time evolution
of the loose sediment bed is sufficiently slow that it appears
stationary to the rapidly flowing fluid, hence the flow can be taken
as quasi-steady; 2) the flow can be considered inviscid (but it
does not have to be assumed irrotational). The dynamics of the free
surface can then be described by requiring conservation of mechanical
energy along streamlines (see e.g. Batchelor 1967). Each surface
streamline is subject to the following set of ordinary differential
equations:
, | (16) |
where ζ is a
curvilinear coordinate measured along the streamline projection in the
horizontal plane (x, y), u = (u,
v, w) is the water velocity, and η is the free surface
elevation. Equations (15a,b,c) describe the streamline as a parametric
curve , along which the Bernoulli
equation (16) is written. Vertical velocity w can be eliminated
from the Bernoulli equation using (15c). It follows that if the
horizontal components of the surface velocity field
(u, v) of a steady flow are known (and if
dissipation can be neglected), one can solve for the free surface
elevation η
along each streamline, up to a constant of integration. The constant
can be determined by invoking a suitable boundary condition (point of
known elevation along the streamline) or stationarity condition (known
average elevation). This constitutes the basis of the velocimetric
technique proposed here for the measurement of free surface
topography. On a streamline-by-streamline basis, the above principle
can be applied to a wide class of flows beyond the present antidune
case: examples include flows over oscillatory bottoms, uniform
currents past isolated accidents of topography, or shooting flows
controlled upstream by a vane of known opening.
Before proceeding, it is useful to reduce the system to a
dimensionless form. Given a characteristic scale [u] for the
velocity, a characteristic wavenumber can be defined as (Reynolds
1965)
. | (17) |
Conversion to dimensionless variables can then be achieved by the
transformation:
resulting in the following set of dimensionless equations:
. | (21) |
5.2. Approximation for small amplitude oscillations over a rapid mean current
An
especially simple situation arises if the free surface can be
approximated as small amplitude oscillatory perturbations superposed
upon a rapid mean flow. This is the case for the present antidune
flows, for which one can further assume that: 1) the mean flow is
oriented in the longitudinal x direction; 2) the mean flow is
approximately uniform, with only a weak mean velocity gradient along
the transverse y direction. In that case, it is possible to
decompose velocities into mean and perturbation components and
introduce a small parameter ε
such that
,
, | (22a,b) |
where is the mean current velocity, and
measures the strength of the perturbation
velocities u', v' (having characteristic scale
[u']) relative to the mean current velocity (having characteristic scale [u]). One can
then substitute expansions
 | (23) |
 | (24) |
in eqs (20)-(21) to obtain, at order 0,
, | (25) |
and at order 1,
 | (26) |
where it was assumed that transverse variations in mean velocity
are at most
O(ε ).
One can finally translate back to dimensional variables and obtain the
key result
, | (27) |
which is particularly simple both to use and to interpret:
superimposed upon a mean surface level , the free surface topography perturbations are proportional to the perturbations in
longitudinal velocity u', but opposite in phase. For torrential
flows, the water velocity is slowest at the crests and fastest at the
troughs, an observation pointed out already by Leonardo Da Vinci in
his notebooks several centuries ago (Da Vinci ca. 1500/1975).
In practice, the ratio of fluctuation velocities to mean velocities
obtained in the present
experiments is about 0.1, with mean velocities of about 1 m/s. This
leads to errors on the perturbed elevation of the order of
(1 m/s)2(0.1)2/(10 m/s2) ≈ 1 mm. More precise estimates
obtained by calculating second order corrections are listed in Table
1. These second order corrections are comparable in amplitude to
measurement errors. Because they are also highly sensitive to
measurement noise, only the first order estimate is retained in what
follows. Whatever the order of approximation, an obvious limitation of
the velocimetric method is that only the topography perturbations are
accessible, not the mean elevation. If needed, this has to be obtained
by other means.
5.3. Horizontal velocities in world coordinates
Horizontal velocities are obtained as
follows from the image sequences acquired with the top camera
(viewpoint V = T). First, particle positions on the digital
images are pinpointed using the algorithm
of section 3.1. Particles are then tracked in the image plane using
the Voronoï match algorithm of section 3.2. Index pairs
j(i) are obtained, identifying matches , likely to
correspond to the same physical particles on successive images.
To convert image coordinates into world coordinates , the image positions are first projected onto a
plane corresponding to the mean water surface . This is carried out as follows. To each particle
image viewed under camera viewpoint T is
associated a ray
 | (28) |
where and are given by equations (3a,b). On the other hand,
plane Π0 can
be specified by its Cartesian equation , where is the unit
normal to the plane and its signed
distance from the origin (with an arbitrary
point belonging to the plane). Upon substitution, the parameter value
corresponding to the point of
intersection of the ray with the plane is given by
. | (29) |
In the experimental configuration, plane Π0 is close to horizontal
with Cartesian equation . Once world
coordinates have been approximated by their
projections in plane Π0, the horizontal
components of velocity vector
ui,m are estimated from:
,
, | (30) |
where Δt is the time interval separating successive
frames.
The above procedure would be sufficient to obtain accurate
measurements of the horizontal velocities if either the water free surface was perfectly
planar, or if we had been able to place the camera far above the water
surface with a precisely vertical optical axis. In practice however,
the free surface is perturbed by the antidune oscillations and the top
camera is oriented slightly obliquely (as illustrated on Fig. 2). As a
result, a non-negligible bias results when the projection is performed
onto a hypothetical non-perturbed plane. To correct for this bias, an
iterative procedure is used: horizontal velocities derived from
particle motions projected onto the mean water surface are used to
construct a first estimate of the perturbed free surface
topography. Particle motions are projected again onto the estimated
perturbed topography, then a new horizontal velocity field is
obtained, and so on. Since free surface slopes depart only moderately
from horizontal, this procedure converges very quickly and two
iterations were found sufficient.
5.4. Averaging of the velocity field
An averaged horizontal velocity field remains to be extracted from
the full set of horizontal particle velocities defined at horizontal positions . Rather than combining right away the vectors
corresponding to different times tm
into one single set to be indiscriminately binned and averaged,
the regularisation procedure adopted here seeks to preserve the
streamline structure which underlies the various individual velocity
vectors. This structure presents important implication for the error
formation process, and it is sought to account for it by splitting the
averaging procedure into two steps: 1) along-trajectory averaging; 2)
transverse averaging. For each of the two steps, a balance between
attenuation and noise error can again be sought in order to minimise
the combined error.
Introducing hats to distinguish measured estimates from hatless
"true" quantities, consider measurements acquired at successive time instants
tm along a single particle trajectory. We have the
standard particle tracking velocimetry (PTV) estimate (e.g. Adrian
1991)
 | (31) |
where measurement approximates the
true particle position at time
tm, and measurement approximates the true instantaneous velocity at intermediate time . A simple stochastic model for error formation
can be derived by assuming (Wernet and Pline 1993) that measured
positions correspond to the actual particle
position to which is added Gaussian white
noise, i.e.
 | (32) |
where σx is the root mean square (rms)
error on particle position, and where the
γ m are independent Gaussian random
variables of zero mean and unit standard deviation. Invoking (31) and
(32) as well as Lagrangian relation
x(t) = ∫ u(t)dt,
velocity estimate can be expressed
in terms of the actual velocity history u(t) as
, | (33) |
where
. | (34) |
This constitutes an observation equation (Honerkamp 1998) for the
PTV measurement process, i.e. an idealised description of the
connection between the true signal and the measured one. From
expression (33), two discrepancies are seen to arise between
measurements and the actual instantaneous
velocities u(t). The first discrepancy is a corrupting
high frequency noise due to random position error. The second
discrepancy results from averaging over interval Δt, which attenuates high frequency
velocity variations.
Minimisation of the combined error can again be sought by balancing
the effects of noise and attenuation. A filtered signal can be obtained by convoluting signal with a rectangular window of width
(2n+1)Δt:
. | (35) |
(The same operation can also be applied to the transverse velocity
measurements to yield along-trajectory filtered
values .) Substitution of (35) into (33)
yields an observation equation for the filtered velocities
 | (36) |
where
. | (37) |
Relation (36) is equivalent to the earlier observation equation
(33), except that it is applied to the longer time interval
(2n+1)Δt. From (36),
the root mean square (rms) error on the filtered velocity due to
position error is given by
. | (38) |
Because measurements are obtained along a trajectory, position errors
cancel each other along the way and only the endpoint errors
matter. As a result, the rms error on velocity (38) due to position
errors decreases in proportion to 1/(2n+1), a much faster
decrease with filter width than the factor which would be obtained if errors on velocity were
uncorrelated Gaussian random variables. It is to account for this
feature of the error formation process that a two step procedure
starting with along-trajectory filtering is selected in preference to
indiscriminate binning of the whole set of velocity measurements.
Because the present experiments feature variations of relatively long
wavelength, one can approximate , yielding for the
random position error the rough estimate
. | (39) |
The resulting values for the different runs are listed in Table 1.
Converted back to pixel coordinates, the particle positioning error is
around 0.1 pixel. Due to the limited resolution of the camera, such
small errors in pixel units still translate into significant errors on
velocities.
The corresponding attenuation error can be estimated by idealising
the velocity history of a particle along its streamline as a harmonic
signal, i.e.
 | (40) |
where angular frequency ω ~ 2π[u]/[λ], amplitude , and characteristic scales [u], [λ], and [u'] can be estimated
for a given run (see Table 1). Alternatively, equation (17) can be
invoked to approximate . Convoluting
signal (40) with a rectangular window of width (2n+1)Δt results in filtered signal
 | (41) |
where the damped amplitude
is given by (e.g. Jähne 1995)
 | (42) |
in which again sinc(x) = sin(x)/x. Developing
(42) by Taylor around the origin and setting , the rms attenuation error on the filtered signal
can be approximated by
. | (43) |
Balancing position noise error (38) and attenuation error (43), one
obtains for the conditions of the present experiments an optimal
filter width (2n + 1) = 13. The estimated error on
trajectory-filtered velocities due to either noise or attenuation is
of the order of 0.01 m/s.
The last remaining operation is transverse averaging. This is
performed in two steps. First, filtered along-trajectory velocities
are sampled at the intersections of the trajectories with transverse
planes placed at equal intervals xi
along the x-axis. This yields for each
xi a set of irregularly distributed data
points . The second step of the procedure
consists of converting these measurements into regularly spaced data
at positions yj of a uniform
grid. The averaging is performed by way of overlapping bins of equal
width Δy, where the bin width Δy is again
chosen to minimise the combined error resulting from noise and
attenuation. Noise at this stage comes both from the imperfect PTV
measurements and from fluctuations in time of the physical
trajectories, and is assumed to be Gaussian. Using the same error
balancing approach as before (see section 4.3), the optimal bin width
is found to be Δy ≈ 35 mm. The
estimated error due to either noise or attenuation at the transverse
averaging stage is again of the order of 0.01 m/s. Noise and
attenuation errors incurred at both the streamline and transverse
averaging stages yield an overall rms velocity error of the order of
0.02 m/s. This amounts to errors of the order of 2 % relative to
the mean velocity (see Table 1). Upon invoking the approximate
Bernoulli equation (27), elevation errors due to velocity measurement
errors are of the order of 2 mm.
To
these must be added the neglected contribution from second order
perturbations. Assuming the errors to be independent, an overall
error is again estimated by adding the squares of the rms
contributions and taking the square root. The resulting values for
the various runs are given in Table 1. Reconstruction errors arising
from interpolation onto a regular grid are not significant for the
velocimetric method because interpolation is conducted to first order
accuracy. The interpolants used are piecewise linear rather than
piecewise constant. The latter approach, less accurate but more
robust, is used for binning in the stereo method due to the presence
of stereo mismatches.
| |