Stereoscopic measurement of free surface flows

The present paper adresses the problem of combined three-dimensionnal measurements of shape and velocity of moving free surfaces. A measurement method based on the cross-correlation of image pairs obtained from a calibrated stereoscopic vision system is presented. The underlying concept of the method consists in the generation of parametric shape and displacement forms which are directly projected on the camera models. This procedure is then integrated in an iterative optimization process so that elevation, orientation, curvature and displacement of each surface subset are accurately estimated. Application is made on inclined plane flows of complex non-Newtonian fluids as an alternative to conventionnal rheometric devices.


Introduction
Advanced two-and three-dimensional Digital Image Correlation (DIC) methods are widely used in the material and solid mechanics community with applications to shape and deformation measurements (see [8] for a review).In the fluid mechanics community, the Particle Image Velocimetry (PIV) technique and its derivatives have become unavoidable tools for the two-dimensionnal analysis of a wide spectra of flows, and three-dimensionnal applications are becoming more and more common.Despite the many practical differences between DIC and PIV methods, their fundamentals are identical, and a common application can be found in the study of free surface flows.
Most common 3D surfacic DIC and PIV methods extend the concept of optimized 2D correlation of image subsets [1] in the context of stereovision so that 3D coordinates as well as the three displacement components can be obtained after triangulation [5].
For free surface or interfacial fluid flows, the three-dimensionnal Interfacial PIV technique (3D-IPIV) [18] derived from the Stereo PIV (SPIV) [14] and refined SPIV algorithms [11] have been succesfully applied to the simultaneous measurement of both free surface shape and velocities.Other PIV-based techniques have been adapted to the measurement of free surface using an extension of the synthetic Schlieren method [12].Alternative, non correlation-based techniques have been applied to the measurement of both shape and velocity of free surface flows [2], or to the 3D tracking of particles near a free surface [10], while a number of developments mainly focus on the measurement of the sole free surface [17,19,6].
In an effort to provide measurement solutions adapted to the analysis of free surface flows, a DIC/PIV method is proposed that identifies parametric free surface shapes and tracks their displacements and deformations using the direct projection of the physical scene onto a calibrated stereovision system.For validation purposes, the technique is applied to an inclined plane flow of a polyoxyde ethylene solution.
14th International Conference on Experimental Mechanics

Measurement method
The proposed measurement technique belongs to the family of DIC/PIV methods and has been developped for hydrodynamic and rheological applications as an alternative to the traditional SPIV-based algorithms used in the fluid mechanics community.The initial formulation used to represent and identify the measured surfaces is based on an Eulerian representation in Cartesian co-ordinates of the height of a moving fluid column and its successive derivatives.Formally, such an approach can be compared to the developments proposed by Helm et.al [7].

Procedure
The measurement method is organized as follows: -The fluid domain is represented by Cartesian coordinates (x, y, z), where (x, y) correspond to a fixed reference plane, typically the flat bed of a flume.-Each measurement point consists initially in a set of physical 2D coordinates (x 0 , y 0 ) representing a position of interest in the reference plane, and for which the fluid height, shape and flow velocity are unknown.-A square-base free surface subset surrounding (x 0 , y 0 ) is identified and parametrized by its height h(x 0 , y 0 ), orientation ∂h/∂α) α∈(x,y) and curvature ∂ 2 h/∂α∂β -The displacement (D x , D x , D z ) and linear deformation ∂D α /∂β) (α,β)∈(x,y,z)×(x,y) of the surface subset are measured.
For the measurement of free surfaces, the method is based on the reconstruction of image patterns generated from a parametric shape predictor and projected onto the camera models.The projected image coordinates are then textured according to their position in the experimental images using bilinear interpolation.The candidate parametric shape is therefore represented by patterns of gray levels corresponding to the candidate 3D scene as viewed from each camera.Cross-correlation is applied between these patterns so that a cross-correlation coefficient, or a cross-correlation map, is obtained.Thus, the measurement can be validated if the shape predictor leads to a cross-correlation coefficient approaching one.The full procedure is integrated into an iterative optimization algorithm which generates the shape parameters and whose objective function is the cross-correlation coefficient itself or a combination of the cross-correlation coefficients obained from more than one pair of cameras.
The measurement of three-component surface velocities is performed once surface subsets have been parametrized.Then, a similar procedure including both the shape and displacement parameters is used, in which cross-correlation is computed for each pair of sequential images.Thus, the objective function of the optimization algorithm is the sum of the cross-correlation coefficients returned from each camera sequence.

Image interpolation
Classicaly, each grayscale image recorded by camera #i is represented by a discrete distribution of gray levels gi (ũ i , ṽi ) recorded at pixel (ũ i , ṽi ), and which will be interpolated as a continuous bivariate function g i (u i , v i ) of real-valued sensor coordinates (u i , v j ) using bilinear terms.

Free-surface identification
A shape parameter p s comprised of the height, orientation and curvature of the unknown free surface subset located above a measurement point (x 0 , y 0 ) is defined as: so that, for a point (x, y) = (x 0 + ∆x, y 0 + ∆y) located in a neighbourhood of (x 0 , y 0 ), the candidate free surface elevation is expressed as the second order bivariate Taylor development of h in terms of (∆x, ∆y): Given the calibration of camera #i with a pinhole camera model of matrix M i , the sensor coordinates are obtained from the physical coordinates (x, y, z) using the usual projection formula: Thus, for a stereoscopic system comprised of two cameras (#1 and #2), the normalized cross-correlation coefficient corresponding to the free-surface candidate at point (x 0 , y 0 ) is defined as follows: where Ω(x 0 , y 0 ) is a discrete square neighbourhood of the point (x 0 , y 0 ) defining a list of physical coordinates (x, y) over which the grayscale functions values g i (u i , v i ) are evaluated, averaged as ḡi , and cross-correlated.

Displacement measurement
In a similar manner, a displacement parameter p d is defined as follows: The candidate displaced positions of the initial surface subset are computed from the initial positions defined in Equation 2 using an additional first order development of (x, y, z) in terms of (∆x, ∆y): The real valued sensor coordinates for the final images and the corresponding grayscale function values are then computed from these displaced positions using the calibation matrices.
Once an initial free-surface subset has been identified, the correlation coefficient between two successive images obtained by camera #i is defined for the validated shape parameter p s and the candidate displacement parameter p d as: 14th International Conference on Experimental Mechanics where the superscript 0 refers to constant values calculated using the fixed parameter p s in the initial image.The variables without superscript refer to the final image and are obtained using both p s and the candidate displacement parameter p d according to Equations 2 and 6.
For a binocular stereoscopic system, the sum of the correlation coefficients for the displacement parameter is then: which is therefore used as the objective function of the optimization algorithm.Once the displacement parameter has been validated, velocities (V x , V y , V z ) are classically estimated using the displacements (D x , D y , D z ) measured between two set of images divided by the time lag separating the two exposures.

Optimization algorithm
The Nelder and Mead simplex algorithm [13] is used for maximization of the correlation functions 2 and 6 against parameters p s and p d , respectively.The optimization algorithm is initiated using a first approximation of p s or p d and stops when the simplex reaches a limit size.
Such a measurement method presents several advantages: first, it allows the search for any type of parametric shape, and not only the localization of matching areas between two or more images, so that, for example, surface orientation, curvature and velocity gradients can be directly obtained.Second, the projection of parametric shapes to generate image patterns allows measurements at user-chosen locations, whereas classic stereo-matching techniques need an interpolation of the triangulated result to measure the surface at specific points.Third, the projection of parametric shapes is not affected by the respective position and orientations of the cameras, which suppresses some issues associated to image deformation and the need to rectificate images [4,11].
Concerning the reconstruction of image subsets according to a shape parameter, two aspects can be discussed.For surface subsets that are well aligned with their projection line on the camera sensor, the interpolated grayscale information will show poor contrast, as the corresponding pixels will be confined in a thin strip of the sensor.On the other hand, an image reconstruction based on surface subsets allows to take account of the physical size of the surface markers for the definition of the measurement neighbourhood, so that their apparent size can be controlled and will be identical on each camera.This point may prove useful in pratice for particle markers, such as those used in PIV, as their number, size and appearence in the correlation window are crucial parameters [20].
Intrinsequely, this principle of measurement is not limited to stereoscopic system, as the optimization step can handle any combination of cross-correlated images patterns which can be generated from more than two cameras, whithout extending the number of optimization parameters.

Inclined plane channel
The experimental device consists in an inclined rectangular channel realized in PMMA and designed to generate impulsively started free surface flows, such as encountered in the case of dam-break situations (Figure 1).The channel is 350mm in width, 690mm in length and 100mm in height, with an open-end that ensures a free flow condition along its lower edge and a compartment closed by a sluice gate in its higher part, which is used as a fluid tank.An adjustable resting frame allows a precise inclination α of the channel, which was set to 10 o in the experiment.

Material
The polyoxyde ethylene solution is prepared from 300g powder of polyoxyde ethylene powder poured into 2700g of deionized water, leading to a mass concentration of 10%.Then the mixture is let at rest 12002-p.4 14th International Conference on Experimental Mechanics for two weeks to ensure the complete diffusion of polymer chains into the water.Homogeneity of the fluid during the experiment is ensured by a manual mixing of the solution prior to the measurement process.The density of the fluid is measured using a picnometer and a precision weighing machine and is equal to 1200 kg/m 3 .

Acquisition devices
This flow is recorded using two Jai CV-M2 CCD Cameras equipped with 50mm f/22 Nikon lenses and linked to the R&D VISION HIRIS real-time image acquisition software.The cameras are disposed in a classic stereoscopic arrangement, symetrically oriented at ±30 o from the normal to the channel floor around the main flow direction, and elevated of approximatively 0.8m, as shown in Figure 1.Sheimpflug mounts are used to align the focal plane of each camera with the channel floor.

Camera calibration
Each camera of the binocular stereoscopic system is independently calibrated using the method of Faugeras and Toscani [3].The calibration markers consist in 260 circular dots regularly distributed on a two-level Cartesian grid Lavision Typ#20 (x − y steps of 15mm, z-levels of 10.5mm and 12 mm) which are automatically detected using specific in-house pattern recognition algorithms [15].The calibration pattern conforms to the coordinate system depicted in Figure 1, so that camera calibration can be directly included in the surface identification algorithm.An additional calibration has been carried-out using a more refined 0.5mm-step Cartesian grid positioned at height z locations, in order to assess the linearity of the optical system over the whole depth of measurement.Optical distortion was found to be negligible, so that the sole Typ#20 calibration pattern associated to the Faugeras-Toscani calibration method could be used satisfyingly.

Seeding and illumination
Similarly to most of the measurement methods using cross-correlation of images, speckle image patterns are used in order to enhance the correlation signal.For that purpose, light sub-millimetric dark organic markers are randomly spread upstream of the measurement region during the displacement, while the flow is illuminated by a centered 600W light normally aligned with the channel floor.An average of 10 4 particles is recorded in each image, and contrast between the lightly colored fluid and the darker markers was found sufficient to efficiently apply the cross-correlation procedures (Fig. 2).Despite the need to spread the markers over the flow during the experiment, this surface seeding technique has the advantage of preserving the rheological properties of the fluid, which would not necessarily be the case if volumic seeding was used.

Results
A series of 500 double images of the polyoxyde ethylene flow is acquired at a rate of 1Hz, ensuring time-resolved measurements of the surface motion.The interrogation region is a 100mm×120mm cartesian grid situated at the free end of the channel and centered along the direction transverse to the flow.In order to assess the feasability of rheometric analyses using this experimental device, the interrogation region has been subdivided into the five 20mm square zones depicted in Figure 1.
Synthetic particle images complying with the experimental devices have been analysed in order to assess the accuracy of the measurement.These preliminary results indicate that the free surface as well as millimetric displacements can be measured with a 50µm accuracy for slopes not exceeding 20%, which is far beyond the actual measured values.
The identification of the surface shape and motion has been carried out on overlapping 12.5mm square subsets over a 21×25 cartesian grid of 5mm cell size.The experimental results show that the local free surface shape and local displacements are actually obtained with sub-millimetric accuracy, as it is observable in Figure 3, in which the disjoint surface subset and velocity distribution exhibit a good spatial continuity.According to the theory of inclined plane flows, in which only the streamwise component should be non-vanishing [9], a discrepancy criterion between this component and the norm of the velocity can be estimated.Here, the discrepancy is below 0.3% for the central zones and below 0.6% for the side zones, showing good agreement with theory.This series of data allows the extraction of height-velocity characteristics of the flow, leading to the identification of the rheometric law of the solution in terms of shear rate γ and shear stress τ at mid-height positions: Figure 4 shows the rheograms resulting from the complete series of measurements, in which the free surface height decreases from 15mm to 5mm and the velocity from 1mm/s to 0. These result can be compared to rheograms obtained over a wider range of shear rates with a classical rheometer using a coaxial cylinder geometry (Gemini, Malvern). 12002-p.6 14th International Conference on Experimental Mechanics The subdivision of the interrogation area into separate zones indicate that edge effects significantly affect the measurements in the left and right zones for high shear rates, showing that the uniform plane flow assumption can only be verified near the central line of the flume.The slight discrepancy of the rheograms obtained in the front, central and rear zones also indicate a slight increase of the shear rate near the free end of the flume.

Conclusion
A stereoscopic, cross-correlation-based, measurement technique has been developed to identify simultaneously the three-dimensional free surface and the three-component surface velocities of materials flowing down an inclined channel.The measurement method is based on the optimisation of parametric shape and displacement terms using an iterative algorithm.Tests based on synthetic images and experimental flows have led to statisfying results, showing a sub-millimetric accuracy of both shape and displacement measurements.Rheograms from the inclined plane flow and a conventional rheometer have been succesfully compared and validate the principle of a rheometric application of this technique.However, the theoretical assumptions needed for the obtention of rheometric laws are hardly met over the whole channel, and slight discrepancies have been found to significantly alter the rheograms.An effort to enhance the experimental procedure and further assessment of the measurement accuracy is then needed.This procedure has also been applied to clay suspensions [9] and seems to be an promising solution for the study of coarse suspensions, such as cohesive sediments, whose rheological characteristics are difficult to obtained using conventional rheometers.

Fig. 3 .
Fig. 3. Sample result of surface and velocity measurements corresponding to images of Figure 2. The average norm of the velocity vectors is 1mm/s.

Fig. 4 .
Fig. 4. Rheograms obtained from inclined plane measurements of the polyoxyde ethylene flow (left) and comparison with data obtained from a conventional rheometer (right).