Vibration Data Processing for Bedload Monitoring in Underwater Environments

: This paper proposes a novel data processing framework dedicated to bedload monitoring in underwater environments. After calibration, by integration the total energy in the nominal bandwidth, the proposed experimental setup is able to accurately measure the mass of individual sediments hitting the aluminum plate. This require a priori knowledge of the vibration trasients in order to match a predeﬁned dictionary. Based on unsupervised hierarchical agglomeration of complex vibration spectra, the proposed algorithms allow to accurately localize the transients corresponding to the shocks created by sediment impacts on a steel plate.


Introduction
Underwater bedload transport surveys is important for assessing stability issues such as reservoir silting or channel self-cleaning. To this purpose, sediment traps are currently used to derive the sediment balances.
Monitoring the bedload transport in underwater environments is crucial for understanding stability issues such as reservoir silting or channel self-cleaning. In this context, sediment traps are widely used to derive the sediment balances. The Birkbeck sampler has become one of the preferred method for in situ bedload measurements [1]. Being an integral measure of the transported sediment volume in a certain time interval, this direct method suffers from one main drawback: it has a finite capacity. Without the use of devices such as a sludge pump to empty the accumulating sediment, the sampler fills and eventually ceases to yield data.
An alternative nonintrusive bedload monitoring instrument is the buried geophone station [2], [3], [4]. Under protection of steel plates, several geophones can provide continuous and automatic measurements, even during large floods. Every stone passing the steel geophone equipped plate generates an impulse recorded by the signal acquisition board. Using a calibrated voltage thresholding scheme the grain impacts are recorded and counted. One can notice that this instrument is highly sensitive to the environmental noise, since it is operating at low frequencies (< 1 KHz).
In [5], Bogen and Moen used piezoelectric acoustic transducers for bedload monitoring stations in Norway. The reported results present a clear correspondence with the sediment flow, but no estimation scheme is proposed.
In this paper, the potential of ultrasound (US) transducers coupled with piezoelectric accelerometers for bedload monitoring stations is further explored in a controlled laboratory environment. The obtained results provid interesting sensitivity capabilities, in accordance with the ground truth. The main objectives are to validate the operating frequency of the transducers, to provide consistent calibration and processing and to evaluate the uncertainty of the derived measurements. To pursue these objectives, an experimental platform, illustrated in Fig. 1-(a) has been developed and tested in the GIPSA-Lab controlled tank system. The following sensors have been employed: • four ultrasound contact transducers PAC R50 (EM, R2, R3, R4), • 2 calibrated piezoelectric accelerometers Endevco 233E (A1, A1), • 1 geophone R.T. Clark (G).
The EM US transducer has been excited using the ENI 2100L power amplifier and the Picotest G5100A arbitrary signal generator. The received signals from A1, A2, R3, R4 have been conditioned using the Nexus low noise amplifier and recorded using the PXIe-1082 NI system. After calibration, by integrationg the total energy in the nominal bandwidth, the proposed experimental setup is able to accurately mesure the mass of individual sediments hitting the aluminium plate. This require a priori knowledge of the vibration trasients in order to match a predefined dictionary. We propose a novel framework for detecting and localizing the transients (in the passive configuration) corresponding to the shocks created by sediment impacts on the steel plate. The proposed algorithm is based on hierarchical agglomeration of complex vibration spectra. The multivariate segmentation algorithm proposed in [6] is selected: the multimodal signals are analyzed by exploiting the asymptotic distribution of the covariance matrix of the complex spectra.
Within this context, this manuscript is synthesizing the results presented in [7,8] in order to propose an unified acoustic data processing and analysis framework for the monitoring the bedload transport in underwater environments. The paper is structured as follows. Section 2 introduces the mesurement methodology, Section 3 describes the proposed data processing framework, while Section 4 presents some qualitative and quantitative performance assessment. Section 5 concludes the paper.

Calibration and measurement
For the proposed bedload monitoring system, the static calibration is a direct relationship between the output and the input of the system determined experimentally.
Static calibration requires a reference standard, to which the output is reported. For our system of sensors, the calibration process determines a relationship between the input voltage controlling the EM US transducer and the output characteristic, called static characteristic curve or calibration. Fig.  2 shows the static characteristic in water, measured by transmitting a sinusoid signal and looking at the reception of the peak value of the Fourier transform. The value obtained supports the idea that our system is linear in all cases. In circumstances where the assumption of stationary is unsustainable, linear dynamic models need to be replaced by more sophisticated models. The dynamic characteristic of each sensor is related to the performance when the input EM voltage signal is a function of time. The speed of response of a system refers to the ability to respond to sudden changes in the amplitude of the input signal. Thus, the response rate is calculated based on the parameters describing the system performance under the transition: constant measuring of the time delay (latency measurement), setup time (settling time), the dead time and range dynamics. The obtained dynamic characteristic are illustrated in Fig. 3 and match the previously derived static characteristic. The transfer function of an acoustic transducer can be computed by positioning two transducer elements in transmission configuration: • by continuous excitation of one transducer with wideband white Gaussian noise, the other responds in its nominal frequency band; • by considering the performance parameters from the dynamic characteristic, one can define an impulsive excitation for estimating the system transfer function. Accordingly, by taking the impulse duration (short sine waveform) to be smaller than the system stabilization time, the transducer is sensing this calibration signal as an impulse and not as a step function.
In this manuscript, the second method for estimating the transfer function is employed and the results are illustrated on Fig. 4. After using controlled steel balls for calibration, two types of sediments were investigated and confidence intervals for the mass measurements have been obtained. Notice that the maximum speed in water have been reached before hitting the multi-sensor bedload monitoring platform. Fig. 5 shows obtained bedload mass measurements in a passive configuration in the nominal bandwidth of the impact (< 14 KHz). One can observe that the bedload monitoring platform can provide  competitive results, both in terms of measurement error and confidence intervals in controlled noise free environment.

Hierarchical agglomeration of complex vibration spectra
In order to correctly employ the mesurement system presented in Section 2, reliabale a priori knowledge of the vibration trasients in order to match the predefined calibration dictionary. Thus a preliminary trasient segmentation step is required for the recorded acoustic data.
The first stage of the proposed segmentation scheme is pre-processing. As each acquisition is triggered on the amplitude level of the sensor exhibiting the fastest response time (the accelerometers in our case), the obtained transient signal is not symmetric. In order to improve the spectral representation, each received signal is circularly shifted to the right and filtered in the common bandwidth. The complex spectrum is obtained by applying the Fast Fourier Transform (FFT) on each signal. One advantage of this representation resides in the fact that it is independent of the time of arrival of each transient.

SIRV spectral estimation
In the next step, a multivariate random vector is obtained by concatenating all the complex spectra of the available signals. For each sensor, let k be the m × 1 complex target vector corresponding to the same frequency range. One way to model the statistical properties of this multivariate random vector leads to Spherically Invariant Random Vectors (SIRV) [9]. It is defined as the product of a square root of a positive random variable τ (variations in power) with an independent circular complex Gaussian vector z with zero mean and covariance matrix [M] = E{zz H } (Gaussian kernel): where the superscript H denotes the complex conjugate transposition and E{·} the mathematical expectation. SIRV representation is not unique, so a normalization condition is necessary. The ML estimator of the normalized covariance matrix under the deterministic texture case is the solution of the following recursive equation : with Tr([M] FP ) = p.
Its asymptotic distribution can be assimilated to the Wishart Probability Density Function (PDF) [9]. In this way, each random vector is described by its normalized covariance matrix, which is independent on the total power at the reception. This will form the feature space for the transient segmentation algorithm.

Hierarchical segmentation
The hierarchical segmentation algorithm from [6] is adapted to the vector of complex spectra. The segmentation strategy is a classical iterative merge algorithm. At each iteration, the two segments which minimize the Stepwise Criterion (SC) are merged. The basic principle of the proposed hierarchical segmentation can be divided into three steps : 1. Definition of an initial partition (which is formed by the acquired buffers, in our case). 2. For each segments pair, SC is computed. Then, the two segments which minimize the criterion are found and merged. 3. Stop if the maximum number of merges is reached, otherwise go to step 2.

Similarity measure
At each iteration, merging two segments yields a decrease in the log-likelihood function. The stepwise criterion is based on this consideration. The hierarchical segmentation algorithm merges the two segments S i and S j which minimizes the loss of likelihood of the partition (which is defined as the sum of likelihoods of partition's segments). The stepwise criterion (SC i,j ) can be expressed as: where MLL(·) denotes the segment maximum log-likelihood function. It is the log-likelihood of the segment (samples in each segment are considered independent realizations) with respect to the assumed probability density function (the Wishart distribution in our case) whose parameters are estimated in the maximum likelihood (hence, the name) sense. Its expression is given by: θ S represents the set of distribution parameters (normalized covariance matrix in our case).

Generalized Maximum Log-Likelihood (GMLL)
In general, the normalized covariance matrix is unknown. One solution consists in replacing the SIRV parameters by their estimates. After replacing the covariance matrix [M] by its respective ML estimators, the stepwise criterion becomes: where GMLL(S) is the generalized maximum log-likelihood function for segment S. For the Wishart PDF, the generalized maximum log-likelihood function for segment S is:

Stop criterium
The use of the L-method [10] has been considered in this paper. This method employs the very error (quality) function that is used to perform cluster merging during the hierarchical segmentation algorithm, specifically the Log-Likelihood Function (LLF) of the partition (i.e. the sum of the MLL values for all the segments of the partition). As this is readily computed during the proposed method, no further computational effort is required. The knee of this error function is identified and the optimal number of clusters is chosen at that point. The knee of a curve is somewhat similar to the point of maximum curvature.

Results and Discussion
The proposed experimentations were carried in the GIPSA-Lab controlled tank system. Figs. 6 and 7 illustrate the transient signals obtained from the four sensors (accelerometers and ultrasound) in passive configuration, after right circular shift and common band pass filtering. In order to build the testing database, five metallic balls with different diameters have been selected for the experiments. For each ball, 25 independent impacts have been recorded using the bedload transport monitoring platform from Fig. 1. The balls have been manually released form the surface of the water at approximately the same position on the metallic plate. An independent buffer has been recorded for each impact with the four sensors simultaneously and coherently.
After concatenating the 25 × 5 buffers, Fig. 8 shows the derived signals used to test the proposed hierarchical agglomeration algorithm.
The obtained results are presented in the confusion matrix from Tab. 1. There are at least two facts to be noticed. Firstly, the diagonal structure of the confusion matrix reveals a rather good classification accuracy. Secondly, the L-method stoped the segment merging at 8 classes. This is explained by the fact that the impacts of the ball B4 have been splitted in two classes, while some impacts from B1 and B3 generated two additional classes. This is in agreement with the subjective visual assessment of the signals from Fig. 8. Table 1. Quantitative performance assessment: confusion matrix. SC (from 1 to 8) are the classes obtained by the proposed algorithm, while TC (from 1 to 5) are the ground truth based classes. Finally, an overall quantitative indicator of the performance of the proposed algorithm is the classification accuracy A, defined as:  with TP Bi the true positive for the metallic ball Bi, TN Bi the true negative and T = 125 the total population. After fusing the classes SC 5 and SC 6, SC 2 and SC 8, SC 4 and SC 7, we can compute the total classification accuracy according to Eq. 7. The obtained value is A = 90.9%.

Conclusion
This paper proposed a new framework for detecting and localizing the transients corresponding to the shocks created by sediment impacts on the steel plate. The proposed classification strategy is based on hierarchical agglomeration of complex vibration spectra: the multimodal signals were analyzed by exploiting the asymptotic distribution (SIRV) of the normalized covariance matrix. Qualitative and quantitative performance assessment has been carried out using vibration signals recorded by the multi-sensor bedload transport monitoring platform form the GIPSA-Lab.
Future work will enroll in two main directions. Firstly, we will try to explore as much as possible all the benefits of the proposed classification algorithm in real life scenarios. Secondly, we will continue with improving the description of transient vibration signals by using non-stationary time-frequency representations instead of the Fourier transform.