Image-Based Bed Material Mapping of a Large River

The composition or bed material plays a crucial role in the physical hydromorphological processes of fluvial systems. However, conventional bed material sampling methods provide only pointwise information, which can be inadequate when investigating large rivers of inhomogeneous bed material characteristics. In this study, novel, image-based approaches are implemented to gain areal information of the bed surface composition using two different techniques: monocular and stereo computer vision. Using underwater videos, captured in shorter reaches of the Hungarian Danube River, a comparison of the bed material grain size distributions from conventional physical samplings and the ones reconstructed from the images is carried out. Moreover, an attempt is made to quantify bed surface roughness, using the so-called Structure from Motion image analysis method. Practical aspects of the applicability of image-based bed material mapping are discussed and future improvements towards an automatized mapping methodology are outlined.


Introduction
There is a permanent interaction mechanism between the water flow and the riverbed, which on the one hand determines the sediment transport capacity of the river yielding to different morphodynamic processes, and on the other hand, the riverbed itself, through the hydraulic resistance, influences the flow features. This back-and-forth effect eventually leads to a complex and continuously changing system of the rivers, possessing great practical importance. For instance, the conditions of fluvial navigation strongly depend on the riverbed morphology in the low water regime, when the ship draught can easily exceed the available flow depth. Bank infiltrated drinking water production depends on, e.g., porosity, stratification, permeability, armoring, colmation, and characteristic sediment grain sizes [1,2]. Moreover, besides the human-related importance, the flow-bed interactions also affect the conditions of aquatic habitats. For example, grain sizes, hiding-factor, inhomogeneity, bedforms, and porosity are all essential for the spawning fish or macroinvertebrates. When analyzing reach scale hydromorphological processes, there is a strong need to gain detailed information on the flow features, the riverbed morphology, the sediment transport, and as a crucial boundary condition: the composition of the riverbed.
Despite its importance, the conventional bed material sampling methods, in a lot of cases, are time-consuming, expensive, and not representative enough, especially when dealing with large rivers of inhomogeneous bed composition. Such problems can easily arise when investigating the transition zones of rivers, i.e., where gravel dominated bed is transformed into sand bed, showing a highly mixed composition, which can vary both in time and space. In these cases, the point-like, local data

Study Sites
The field tests were performed along the Hungarian section of the Danube River. This reach of the Danube is a transition zone in terms of bed slope and bed material, where the bed material composition shows strongly varying fractions from silt to gravel even within short distances [20]. For practical reasons, two study sites were chosen, located approximately 100 km apart. The first study site is located close to the settlement Gönyű between rkm (river kilometer) 1791.2 and 1790.6 ( Figure 1). The second one is located at Göd between rkm 1669.4 and 1667.2 ( Figure 2). Both study areas can be characterized by an approximately 300 m channel width, a bed slope of 15 cm/km, whereas the mean discharges are around 2200 m 3 /s and 1500 m 3 /s (representing the larger branch along a 33 km long island), respectively. At both sites, 3-3 cross-sections were studied in more detail.  Cross-sectional ADCP (Acoustic Doppler Current Profiler) and video measurements (moving boat measurements) were carried out in the sections marked with red.

Materials and Methods
The field tests can be separated into two groups. The over-all flowchart of the locations, aims, applied methods, and expected results can be seen in Figure 3. Its elements are further discussed in Sections 3.1 and 3.2.

Applied Monocular Computer Vision Approach for Local Bed Composition Investigation
First, the applicability of the monocular computer vision approach was investigated at Study site I. Here, videos were captured of the riverbed from which individual images were assessed ( Figure 4). A GoPro Hero 4 camera in a waterproof case was attached to a~50 kg heavy isokinetic weight, together with a dive flashlight (see the setup in Figure 5). The measurement setup contained a known-sized object fixed right under the camera, aligned with the bottom of the weight so that it touches the riverbed when the setup lands, providing a reference length scale for the image-based analysis (see the screw in Figure 4). The camera was set to 1920 × 1080 HD resolution and 48 fps.  Besides the imaging, physical bed material samples were also taken, using the conventional bucket sampler in the same verticals, where the videos were recorded. The physical samples were needed for validation purposes of the monocular computer vision approach. The physical bed material samples were brought to the laboratory for drying, and sieving to provide their grain size distributions (GSDs). For the image analysis, the so-called transferable wavelet-method was chosen, using the Matlab code of Buscombe [7], which does not need further calibrations, but the adequate measurement instrumentation. This method converts the input image into a greyscale picture and then analyzes the grey-value (intensity) of each pixel in given rows and columns. For each row and columns, the variation of intensity is then treated as signals. Using continuous Morlet-type wavelet transforms (instead of Fourier transforms), the signals are decomposed, and the power spectra of the transforms are calculated, finally resulting in the grain size distribution in pixels. The goal of this study was to test the method in the field performing underwater imaging in a river reach, where strongly inhomogeneous bed material composition presents. Furthermore, the potential in extending the pointwise information into areal, map-like information of the bed composition was investigated.
To get the GSD in the real length scale, the size of the reference object was used. Eventually, sieving and an image-based GSD could be provided for each sample, enabling the comparison of the physical and indirect methods. However, it is relevant to consider that while sieving gives a volume distribution, the image processing results in an area (or grid-by-number) distribution, so a transformation of one of the GSDs was necessary to provide comparable datasets [7,21]. In order to transform the image-based areal information to volume distribution, the characteristic grain shape should be known, i.e., the thickness of the sediment particles. Results from a recent shape analysis were therefore used ( Figure 6) based on the categorization method of Zingg, which suggested that the study site can be characterized with disk shape sediments. The average b/a and c/b values of the dataset were 0.76 and 0.62, respectively ( Figure 6, black X). Grains more or less tend to lie with their shortest (c) axes perpendicular to the riverbed, meaning that the intermediate (b) and longest (a) axes are presented in the photos [22]. Hence, we used the intermediate axes retrieved from the image processing, to calculate their other axes by using the average ratio-pairs for disk shape particles ( Figure 6, black X). Assuming a constant density of 2.65 g/cm 3 , the weight was calculated, providing a by-weight distribution. Performing the previous steps on all the image-based results, the direct comparison of the physical and image-based GSD could be carried out.

Applied Stereo Computer Vision Techniques for Estimation the Bed Surface Roughness
In the second part of the study, the stereo computer vision techniques were studied. For this purpose, cross-sectional video-recordings were carried out in study site II by fixing the imaging equipment on a moving vessel. The camera settings were similar to the case of Section 3.1. During the surveys, a real-time kinematic GPS recorded the actual vessel coordinates. In this way, the grain composition information constructed from the images could be located in the absolute coordinate system. The imaging instrumentation used in the first case study ( Figure 5) had to be slightly modified due to the fact that the object for scale referencing needed to touch the riverbed. This restraint needed to be avoided, therefore the system was improved using laser lights instead (four laser pointers fixed on the frame), with a known distance between each other. The laser lights could be well-identified on the images, enabling the transformation of the image series from the distorted image coordinate system to a two-dimensional (2D) system.
The image processing method called Structure-from-Motion (SfM) was then applied to the image series. The basic idea of the method is to look for and link the common features, points between the input images. If pictures are taken from different but overlapping angles and directions, it will result in a 3D surface model of the given object. The method is using the Scale Invariant Feature Transform (SIFT) algorithm, known from computer vision, to first retrieve and store the strong features and key points from the images, defining their position, magnitude, orientation, and scale. After this, their neighborhood and other features are examined and paired. As a result, a sparse point cloud is generated [24,25]. The next step includes using one of the Multi-View Stereo (MVS) methods to create a dense point cloud [26]. Finally, the points will act as nodes to generate the digital model of the object. Instead of using several cameras for reconstructing the bed surface, we exploited the moving character of the camera, which enabled us to capture images of the same locations from different angles, and sometimes even different distances. For this, the Agisoft Metashape software [27] was used, which is capable of the above-mentioned steps.
Based on the reconstructed surface models, an attempt was made to calculate bed surface roughness. Flow resistance is influenced by skin or grain friction and form drag from developing bedforms. Here, the investigation was carried out for evaluating the former one, by calculating the average roughness height. For this purpose, we followed the theories laid down by mechanical engineering and processing technology, where the surfaces are represented as a composition of different sinusoidal waves. The high frequency waves are defining the roughness, medium frequencies are the waviness, and low frequencies consist of the form of the given surface ( Figure 7) [28][29][30][31]. The suitable decomposition of the surface profile was performed along transects, using the International Organization for Standardization (ISO) 16610-21 closed profile Gaussian filter. With this approach, the so-called Gaussian mean-line (consisting of the low and medium frequencies) could be produced. Calculating the difference between the raw profile and the Gaussian mean-line gives the roughness profile of the section [28]. . The approach to estimate bed surface roughness from the reconstructed bed surface models: along a normal section, surfaces can be decomposed to high frequency (roughness) and medium frequency (waviness) sinusoidal waves (or even lower frequency ones to get the surface form, not shown in the figure). [30].
For this purpose, one must define the value of the so-called cut-off wavelength (λ c ). This value is, in fact, a filter value and is used for the Gaussian-filter function as an input. This parameter distinguishes between roughness (shorter wavelengths than λ c ) and waviness (longer wavelengths than λ c ). For the suitable determination of this parameter, recommendations are given in the mentioned standard (e.g., it should be at least 2.5 times the average spacing between adjacent peaks of the measured profile on the sampling length). For the calculation, a Matlab code was written, based on the ISO 16610-21 closed profile Gaussian filter [31], and a case-specific cut-off wavelength was chosen, after taking the standard recommendations into account and visualizing the data. After the separation of the roughness profile, the average roughness height (Ra) could be calculated with Equation (1) to describe the roughness of the given section: where L [m] is the sampling length (length of the investigated section) and Z(x) [m] is the roughness value of the given x point.

Comparative Assessment of Grain Size Distributions
In the following, a comparison between the GSDs from physical samplings together with the conventional sieving methods and the GSDs constructed from the monocular image-based technique is introduced. As for the latter, two types of composition curves are assessed, i.e., the raw ones resulted directly by the imaging method and the ones after the area-volume transformation, as described above. The median grain sizes (D 50 ) and the standard deviations are also assessed.
Depending on the bed composition, the characteristic grain sizes, and shapes, the agreements between the methods show different behavior. There were samples, where the grain size distributions constructed from the images showed a good match with the sieving based GSDs. Such an example can be seen in Figure 8, taken from the sampling point 1/1 (see the sampling location in Figure 1). After transforming the image-based result, the GSD curves showed an even better match on the lower tail of the curves (<30%). However, they presented coarser fractions on the higher part. A potential explanation can be that grain-shapes of the coarser fractions differ from the assumed disk-shapes. The standard deviation also showed a better match with the sieving after the transformation, meaning the shapes of the curve are better captured than before the transformation. Thus, in these cases, the fractions above 10 mm are better estimated without the area-volume transformation. Furthermore, there were several samples where the area-volume transformation improved the agreement, such as in the case of point 2/2 ( Figure 9). These samples are generally taken from the coarsest zones of the riverbed, where no sand and fine gravel fractions present. The coarsest fractions are still not captured accurately, but the method works best when no fine fractions can be found on the bed surface. The most problematic zones are those locations, where fine, sand fractions are present. It is known that in image-based methods, the minimum resolvable grain size depends on the image resolution [4], and the validity of the hereby used digital image analysis method starts from approximately 0.7 mm (it is 2-3 px normally, [6,32]) due to the applied camera settings and parameters. To see whether the amount of the coarser fractions are well reconstructed by the image-based method, we removed the fractions that were finer than 1 mm from the sieving based GSD. An example from point 3/4 can be seen in Figure 10. Seeing how the curve moved closer to the image-based curve, we considered it as a proof for the expectations and earlier experience regarding the limitations of the image-based methods (originating from simply resolution limits). The reason for the original image-based GSD showing better match with the sieving (after removing sand), than the transformed result, might be because, in this case, the shortest c-axes (non-visible for the camera) of the grains is not that significant in regards to volume (and mass), so instead of the applied ones, other ratios, or even shape, would be more realistic. Involving all the samples, a comparison of the median grain diameter (Figure 11a) as well as the standard deviation values (Figure 11b) were performed between the image-based and sieving based methods. Overall, an adequate agreement was found in the gravel fractions and weak estimation for sand fractions.

Bed Surface Model and Estimation of Bed Roughness
As described above, videos captured at the second study site along the river cross-sections were used for performing a stereo computer vision-based, SfM technique. In the following, the results prepared from survey Section 2 (see location in Figure 2) are shown. First, the series of images extracted from the video were merged to provide a large, continuous image of the riverbed for further analysis ( Figure 12). For the sake of illustration, only an approximately 21 m long section (evaluation length) of the whole cross-section is analyzed and showed, but the method, considering the computational demand of the image processing, can certainly be extended for longer profiles. It is worth noting that the computational analysis of this section took approximately 10 h on a normal PC. After the creation of the bed surface model, a polyline was drawn along the modeled section to extract the total profile and calculate the roughness heights. The next step was to use closed profile Gaussian-filter to get the mean-line (low-and medium-frequency elements of the total profile). The difference between the filtered and original profile gives the roughness profile ( Figure 13). During the video recording, it was seen that in the middle of this section (from 10 to 12 m in Figure 13), the riverbed was coarser, gravel dominated than in other parts (e.g., between 2 and 4 m), where rather the sand dominated. This phenomenon was also visible after creating the 3D model ( Figure 14).  The difference in riverbed composition was also visible in the roughness profile ( Figure 13). In the coarser part, the total height (vertical distance between a maximum valley and maximum peak) was larger, than in the sand dominated part. For the coarser part Ra gravel = 0.0050 m, for the sand part Ra sand = 0.0017 m, while for the whole 21 m long section Ra total = 0.0029 m was calculated, as average roughness height values. Indeed, the riverbed section, where gravel was dominating, the characteristic grain diameter fell in the range of 0.5-1 cm, whereas the typical sand fractions have an order of magnitude lower diameter. This difference in the roughness values could well be captured, underlining that the SfM based fine-scale bed surface model can be used for separating different bed material classes.

Discussion
In this study, attempts were made to implement image-processing methods for under-water morphological field measurements of large rivers by emphasizing the grain size distribution and roughness of the riverbed as parameters of interest. The field tests were done in the transition zone of the Hungarian section of the Danube River, where the conventional sampling methods are generally inadequate for the representative characterization of bed composition. In contrast with earlier studies applying image-based methods for bed material analysis (e.g., [2][3][4][5][6][7][8][9][10][11][12][13][14][15][16][17]), in this study, the investigated river sections show strong inhomogeneity in terms of bed material, ranging from silt to gravel dominated zones. The key question of this research was whether applying a suitable underwater imaging equipment, combining with different image processing methods, could provide areal kind information of the riverbed composition with acceptable quality for further analysis.
As for the developed field imaging equipment, it can be stated that the instrumentation, together with the applied methods, are easy and quick to use. The moving vessel measurements would be the best option for further studies in the future, with the improved instrumentation setup (i.e., lasers), because it can provide adequate input for spatial referencing. It is crucial, however, to note that the introduced underwater imaging is feasible only in low water situations when the water transparency enables clear images. Moreover, there should be no bedload transport during imaging, because the moving particles can bias the resulted grain composition. However, in such circumstances, a suitable image post-processing method could exclude moving particles from the images, as was already introduced by the authors [33].
The applied transferable wavelet-method, as a monocular approach, showed promising results when comparing the constructed GSDs with the ones from sieving analysis. It is, however, important to note that in the case of sand and finer fractions, the method is not able to identify them, resulting in biased GSDs. This insufficiency can be explained with the resolution limitations of image-based methods. To overcome this problem, a combination of the monocular and stereo vision approaches may be necessary, where the latter can account for the separation of fine and coarse particle dominated zones, as discussed later on. The issue of the theoretical differences between sieving and image-based analysis, i.e., volume and areal distribution, was also addressed. A transformation procedure was tested to account for the conversion of areal to volume distributions considering a typical grain shape at the study site. In some cases, the applied estimation gave better results, suggesting that further and more sophisticated involvement of the grain shapes can improve the results.
As for the stereo vision method, a great advantage is the possibility of merging continuously collected images of the riverbed and analyzing grain size distribution for areas with arbitrary sizes, which can eventually yield map-like GSD information. The stereo vision method enables the generation of the 3D model of the riverbed along the camera-path. The herein applied SfM method is reasonably robust. However, it is a key feature to find suitable measurement circumstances in terms of e.g., relative speed of the vessel, camera resolution, framerate, or the distance of the camera from the riverbed. If these parameters are well chosen, this method can indeed provide quantitative information on the bed surface roughness.
The grain-scale 3D model of the recorded bed surface contains relevant information for a better understanding of flow riverbed interactions but can be of great importance for parameterizing computational flow models. Moreover, as mentioned before, flow resistance is dependent not only on the skin friction, but on the form drag as well. The herein applied filtering method is capable of separating the roughness (skin friction), waviness (forms, form drag), and the channel shape from each other, based on the so-called cut-off wavelengths. Hence, it could be further used in the future for calculating other parameters.
As for future improvements, the exploitation of increasing camera resolution and framerate will lead to more and more accurate assessments with e.g., further decreasing the lower-limit for sensible particle sizes for both the wavelet (GSD) and the SfM (grain-scale 3D model and surface roughness) methods. A logical continuation of the research can be the combination of the wavelet-method with Deep Learning algorithms [18,19] that are capable of recognizing overall textures and patterns. With the combination of the two approaches, it could also contribute to a more accurate separation of sand and gravel patches. Regarding the distribution comparison and grain shape assessment, further studies are planned to carry out more sophisticated transformations of the curves (e.g., not just one average shape and axes-ratios, but different ones for each size class). In case of surface roughness, the average roughness height (root mean square deviation) will be calculated from the presented roughness profiles and attempt will be made to correlate these values with Nikuradse roughness heights and compare the this-way-calculated D 50 values to ones that are retrieved from sieving or the wavelet method itself. Moreover, besides the presented line roughness (2D), area roughness (3D) [17] will be calculated as well. The final goal would be the automatization of the wavelet-method to retrieve the grain size distributions from any part of the riverbed, mapped by the moving camera(s).

Conclusions
In this paper, the first results of the assessment of the riverbed composition were introduced based on video imaging of the riverbed. Two different approaches were tested in shorter sections of a large river: monocular and stereo vision methods. The former provided information about the characteristic grain size of the bed surface, whereas the second was used to create fine-scale, 3D bed surface models based on, which, information on the local roughness could be extracted. The results clearly demonstrated the capabilities of both methods, eventually complementing each other. The monocular method showed adequate accuracy for gravel fractions, whereas the grain composition of finer fractions could not be reconstructed well, due to the resolution of the images. On the other hand, the tested stereo vision method could well distinguish between sand and gravel fractions, and so once this separation can be done using the second technique, the first one can be used for more detailed assessment.
The greatest advantage of the introduced image analysis methods is that compared to the conventional, physical bed sampling procedures where local, pointwise information can be extracted from the time-consuming field campaigns, and expensive laboratory analysis, here, after the automatization of the introduced data processing, grain composition information (even along longer streamwise and transversal river transects) can be gained. When the physical characteristics of the bed material strongly vary even in shorter reaches, the quantitative description of these gradients can be crucial.
The results of this study can be well exploited in river engineering problems, when fine scale computational hydro-morphological modelling is performed for shorter river reaches, and map sort of information is needed of the sediment characteristics for model parameterization. Such 3D numerical models are more and more frequently used, e.g., for supporting restoration measures. Furthermore, when, for instance, habitat assessment of rivers is carried out, these image-based methods can also support the description of the varying nature of the benthic zones.
The permanent development of computational capacity, video resolution as well as the new developments in artificial intelligence methods, it is foreseen that these techniques become more and more suitable even for large scale analysis of rivers, moreover, the assessment of static images can be further extended to dynamic images, i.e., the quantification of the moving particles and so the investigation of sediment transport on the grain scale can be a realistic goal.

Funding:
The study was partly funded by European Union-European Regional Development Fund via the Interreg Austria-Hungary SEDDON II project. The second author acknowledges the support of the ÚNKP-19-4 New National Excellence Program of the Ministry for Innovation and Technology and the Bolyai János fellowship of the Hungarian Academy of Sciences. This research was partly supported by MTA TKI of the Hungarian Academy of Sciences. Support of grant BME FIKP-VÍZ by EMMI is also kindly acknowledged. The authors acknowledge the funding of the OTKA FK 128429 grant.