A depth estimation system for laboratory studies using video imagery.

INTRODUCTION

Knowledge of sea bottom topography is of great importance in coastal engineering (1). Often studies which have great environmental impact require data of the spatial and temporal evolution of the sea bottom in order to quantify the proportion and direction of sediment transport. This information can be used to test or calibrate numerical or empirical models or be used for making decisions concerning, for example, dredging or beach nourishment projects. Traditional survey techniques to obtain bathymetric information are precise, however they are normally expensive and some times impossible to implement given that the wave climate should be calm.

On the other hand, the bathymetric estimation of the nearshore topography traditionally requires surveying equipment, a boat and more than one person, thus making it expensive and often impractical. Although technological advances and the use of Global Positioning Systems (GPS) have led to new techniques being developed (2), it is worth looking for alternative approaches. One of these is the use of depth inversion algorithms combined with remote sensing techniques.

The use of video images, for the quantification of nearshore processes, has become popular in the last decades. This technique is known as photogrammetry defined as the art, science and technology of obtaining reliable information about physical objects and the environment through processes of recording, measuring and interpreting photographic images (3). This method is based on the hypothesis that almost every nearshore phenomena that can be discerned visually can be quantified using video processing techniques (4). Recent methods using photogrammetry for nearshore depth estimation are treated in the present work.

There are several approaches for wave number estimation using sensor arrays (4), (5), but not many have used PCA before (6), (7). Wallace and Dickinson (8) used the Complex Empirical Orthogonal Functions (CEOF) analysis for extraction of wave phase information in atmospheric studies. Horel (9) also used this method for identifying travelling waves in geophysical data sets. Merrifield and Guza (6) defined the wave number as the gradient in the resulting spatial phase from the CEOF analysis.

Water depth estimation is made using spatially extensive arrays in order to determine both wave number components, [k.sub.x] and [k.sub.y], in addition to the frequency information. Depth inversion algorithms were first used during World War II when scientists began to take interest in wave propagation characteristics in order to predict surf conditions and bathymetry for military purposes. The bathymetry and wave characteristics were estimated from the early time-lapse aerial photography of the ocean surface. Holland (10) found that for surface gravity waves, in water depths outside the surf zone region, the linear dispersion relation was highly accurate with a depth estimation of 6% different from the observed depth.

In the last two decades, measurements of shoaling wave characteristics have been made using video techniques following the methodology introduced by Lippmann and Holman (11). Lippmann and Holman (12) found good correlation between time series of pressure sensor and pixel intensity in the surf zone. They also found that the celerity spectra compares favourably with solitary wave theory. They were able to estimate wave directional phase speed in the surf zone using video techniques. Stockdon and Holman (7) estimated the nearshore bathymetry based on video images. The propagation of water waves were measured using pixel intensity time series collected from a cross-shore array. The cross-shore component of the wave number was computed as proposed by Merrifield and Guza (6). Water depth was computed using the linear dispersion relation and the results were compared with bathymetry maps, finding a mean percent error of 13%. More recently, Curtis et al. (13), describe a method for video measurement of wave direction in an idealized coastal inlet physical model, based on the estimated root-mean-square average wave number as proposed by Herbers et al. (5).

The purpose of this study is to develop and validate a simple video-based system for depth estimation in an idealized planar beach. The experiments were carried out in the lab using a numerical wave propagation simulation representing completely controlled conditions. The study first encompasses the image rectification, wave number and wave frequency estimation from pixel intensity time series. Secondly, the estimation of water depth in the numerical model simulation is achieved using the wave characteristics estimated before. Finally, the differences between estimated and real depth are discussed in terms of the possible sources of error.

SYSTEM DESCRIPTION

The video-based system test is comprised of the data processing from a wave propagation numerical simulation for a given deep-water wave height [H.sub.0], period [T.sub.0], incident wave angle 0 and nearshore bathymetry h(x,y). The finite-difference model for the wave propagation is based on the algorithm proposed by Dalrymple (14). This explicit numerical model computes the wave refraction and shoaling of linear and nonlinear water waves over irregular bathymetry. The model uses a central difference scheme in cross-shore direction x, giving an error of order [x.sup.2] (14). The linear solution of this model was chosen to be consistent with the use of the linear dispersion equation for depth estimation. As a result the differences between observed and estimated depths would be directly related to the rectification and wave number estimation method applied. Using the numerical model mentioned before, (1) is then applied to introduce time dependence. Based on earlier observations of correlation between time series of pressure sensor and pixel intensity (12), surface elevation plots (2DH) are created for every time step and further used as a frame in an array (movie).

[eta](x,y,t) = [[H(x,y)]/2]cos([K.sub.x](x,y)x + [K.sub.y](x,y)y-[sigma]t) (1)

where [k.sub.x] and [k.sub.y] represent the cross-shore (x) and long-shore (y) wave number components respectively; H = wave height given by the propagation model described before; = 2 /T where T = wave period; t = time. However, some complication arises as the wave approaches the shore because the wavelength decreases. For this reason, the propagation mesh grid was subdivided into three regions with x = y = 10, 1 and 0.2 m, respectively. The simulation was done in gray-scale, where greater intensities correspond to wave crests and lower intensities to wave troughs. A snapshot from the numerical model simulation previously described is shown in Fig. 1, where it is easily observed that wave length and wave angle decrease as the wave train approaches the shore.

[FIGURE 1 OMITTED]

The resulting intensity pattern from the simulation is similar to the one obtained by Curtis et al. (13) when placing Micro-spheres in the Coastal and Hydraulics Laboratory (CHL) wave basin.

In the field, recorded pixel intensity time series are shifted by a phase lag with respect to wave crest, depending on light source location. However, the extracted phase information from the time series remains the same. Because the video analysis technique for depth estimation relies only upon phase, the use of the wave simulation described above is justified.

Data collection: A wave propagation simulation using an idealized bathymetry was projected onto a wall. A video camera set at an oblique angle filmed the simulated field conditions as shown in Fig. 2. Five ground control points were surveyed over the projection area and the video was recorded for further analysis.

[FIGURE 2 OMITTED]

Data processing

Image extraction: The video was imported from a digital video camera to a dual processor. Each file was opened as uncompressed DV/NTSC footage, rendered and exported as QuickTime files using Sorensen 3 compression technology. This reduced file size with little or no loss of quality. For each movie clip, still frames in GIF format were extracted from the raw footage with a reduction rate from 29.97 frames per second to 8 frames per second, ensuring systematic sampling.

Image rectification: An image is the 2-D representation of a 3-D world, therefore one necessary requirement for quantifying the information recorded is the knowledge of this transformation. The image rectification, defined as the process of making equivalent vertical photographs from tilted photos (3), was performed using a planar projective transformation of the form,

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (2)

or more briefly [X.sub.i] = [H.sub.ij][x.sub.i], where [X.sub.i] and [x.sub.i] represent the rectified and image coordinate vectors respectively and H is the transformation matrix which accounts for rotation, scaling and image distortion. Once the matrix H is obtained, which relates the true coordinates to the image coordinates, the matrix is multiplied by every coordinate point and the corresponding rectified positions are obtained. The intensity found at (x, y, z) is mapped to (X, Y, Z). The use of a planar projective transformation method for image rectification presented some advantages over other methods, in that it does not require any information about the intrinsic and extrinsic parameters of the camera, since they are implicit in the solution of the system of (2). For this reason the method was deemed suitable for small-scale experiments due to the simplicity of surveying a large number of Ground Control Points (GCPs).

The image rectification algorithm, used in the present study, is based on the planar projective transformation method previously described. The inputs are the GCPs coordinates, the desired resolution (m/pixels) and the output image size (m x n). The program requires at least four GCPs to solve (2) for H by a least squares method using Singular Value Decomposition (SVD). Once the transformation matrix H is determined, the program rectifies the sequence of images. The values of each rectified grid coordinate are transformed through (2) to determine its corresponding image location. The pixel intensity value found in that location from the original oblique image is mapped to the rectified grid location in the rectified image matrix. Finally, the rectified image is transformed from RGB to gray-scale image and stored in a mxnxp matrix containing all rectified images.

Wave frequency estimation: For wave frequency and wave number estimation, the video technique first requires the collection of pixel intensity time series at an array of pixel locations. As a wave propagates through a fixed point on the free surface, the intensity of the surface reflection varies from light to dark. As mentioned before, in the field this variation depends on the viewing angle of the camera relative to the wave direction and to the light source location. The forward face of the wave will have a higher intensity than the back face (13) as shown in Fig. 3.

[FIGURE 3 OMITTED]

From the video analysis, only phase information can be obtained and wave height cannot be determined. The wave period estimation is performed using spectral analysis of pixel intensity time series, where the spectral peak frequency is selected.

Wave number estimation: The Principal Component Analysis (PCA) method has been applied in several coastal studies. This method is also known as Empirical Orthogonal Functions (EOF). The objective of PCA is to describe the spatial and/or temporal changes of any process by the minimum number of eigenfunctions. The advantage of this method is that the first mode eigenfunction accounts for the largest portion of the variance of data. However, there is a limitation; PCA detects only standing oscillations, not travelling waves (9). For that reason, the method was extended to be able to detect propagation phenomena. An extension to the PCA method is the Complex Principal Component Analysis (CPCA) or CEOF, from which propagating phenomena can be analyzed.

Stockdon and Holman (7) estimated the wave number components from the analysis of wave phase structure using a frequency domain EOF from a pixel intensity array. The present study used an equivalent CPCA in the time domain instead. The analysis procedure of CPCA is similar to that of PCA, except that the analysis first involves converting the data from a real series of numbers to a complex series.

The Hilbert transformation is applied to the original time series, [u.sub.j](t), in order to obtain the imaginary component. From the resulting complex time series, the phase information between stations is found, [U.sub.j](t) = [u.sub.j](t) + I[^.u](t), where the imaginary part is the Hilbert transformation of the real part. The covariance matrix of [U.sub.j](t) is[C.sub.jk] = <[U.sub.j.sup.*](t)[U.sub.k](t)>[.sub.t].

An important observation at this point is that the CPCA in the time domain is equivalent to that found from a frequency domain PCA with the cross-spectrum integrated over all frequencies (6).

The eigenfunctions, [B.sub.n](x), where n denotes mode number, are obtained from C, [B.sub.n](x) = (x) = PCA[C]and the temporal expansions are given by,

[A.sub.n](t) = [N.summation over (j=1)][U.sub.j](t)[B.sub.n](x)

where both [B.sub.n](x) and [A.sub.n](t) are complex. The complex time series can be represented as:

[U.sub.j](t) = [N.summation over (n=1)][A.sub.n](t)[B.sub.n](x)

Where the first mode explains close to 100% of the variance, the spatial phase can be calculated using,

[[PHI].sub.n] = [-tan.sup.-1][Im[B.sub.n](x)/Re[B.sub.n](x)]

Finally, wave number components are calculated in phase from the cross-shore and long-shore gradients as (6):

[k.sub.x] = d[PHI]/dx (3)

[k.sub.y] = d[PHI]/dy (4)

Where the wave angle is given as:

[theta] = arctan ([k.sub.y]/[k.sub.x]) (5)

Water depth estimation: As a result of the application of the dynamic free surface boundary condition during the formulation of small amplitude wave theory, the so-called dispersion equation is derived,

[[sigma].sup.2] = gk tanh kh

This can also be expressed as,

c = [c.sub.0]tanh kh

This relation describes the way in which a field of propagating waves of many frequencies would separate or disperse due to the different celerities of the various frequency components (2). From a different point of view, the dispersion equation relates wave celerity cand water depth h. Knowing the wave number components and the wave frequency, the water depth h(x,y), can be estimated as:

h(x,y) = [[[tanh.sup.-1]([[sigma].sup.2]/[gk])]/k] (6)

where k = [square root of [k.sub.x.sup.2] + [k.sub.y.sup.2]]; [sigma] = 2[pi]/T; g = gravitational acceleration (9.81 m [sec.sup.-1]) and h = water depth.

All depth prediction methods that use these characteristics of shoaling waves to predict bottom topography are referred to as depth inversion algorithms (DIAs) (1). Once the cross-shore and long-shore components of the wave number and the peak frequency are estimated from the video images, the depth at every position can be calculated using (6).

RESULTS

Performance of the image rectification method: In order to show the performance of the rectification method (planar projective transformation) used here, a small and a large-scale case were examined.

Case 1: Small scale test: A square of known dimensions was drawn on the floor using gray tape. Four GCPs were surveyed and an oblique picture (Fig. 4a) was taken and rectified (Fig. 4b). The outer dimensions of the square were 1x1 m and the tape width was 4.9 cm everywhere. The output image resolution was chosen to be 3.3 mm [pixel.sup.-1] on a 600x600 pixel matrix.

[FIGURE 4 OMITTED]

Known distances from the object were compared with the rectified scaled image and the Root Mean Square error (RMS) was estimated for each measurement.

RMS = [square root of [[[[N.summation over (i = 1)][([[X.sub.i]-[x.sub.i]]/[X.sub.i]).sup.2]]/N]]] (7)

where X is the actual measurement, x is the video estimated and N is the number of observations.

The error analysis consisted of the estimation of the tape thickness at several locations from the rectified image shown in Fig. 4b and subsequent comparison with the surveyed data. Using (7) an RMS error of 5% for length was estimated.

Case 2: Large scale test: Image rectification from a multiple-camera system at El Puntal Spit (Santander, Spain) was used to test the rectification algorithm performance under field conditions. A total of 16 permanent and temporal GCPs were used and its rectification is shown in Fig. 5.

[FIGURE 5 OMITTED]

Although less accurate than commercially available systems, this rectification code allows qualitative analysis from large-scale images. However, the algorithm has a strong dependence on GCPs since the parameter of the camera (extrinsic and intrinsic) are implicit in the solution of (2). Another disadvantage presented while using this method is that it requires at least four GCPs in order to perform the rectification, which is sometimes hard to achieve in the field. Recently, more advanced methods for image rectification have been developed. Holland et al. (4) presented a method that requires two GCPs, or even one GCP, given the intrinsic parameters of the camera through the calibration process.

Due to the good performance of the rectification method in the small-scale case, this was used to rectify the video taken from the numerical wave propagation simulation projected onto a wall in order to determine water depth as sketched in Fig. 2.

Depth estimation Numerical model data: Time series of wave data, generated by the numerical model described in Section 2, were used for depth error estimation by using the wave number estimation method proposed by Merrifield and Guza (6). The case corresponds to water 1 m deep, with a long-crested wave of an 8 sec period, approaching the shore with a 30[degrees] angle over a mild-sloping beach (1:40). The resulting spectrum is shown in Fig. 6 with a peak frequency of 0.125 Hz corresponding with the 8 sec period.

[FIGURE 6 OMITTED]

The wave number components were estimated using CPCA, (3) and (4) where the gradients were computed from deep water to shore. As expected, the magnitude of the cross-shore wave number component [k.sub.x] increases as the wave propagates onshore (Fig. 7a) and the long-shore component [k.sub.y] remains constant at every cross-shore position for the case of straight and parallel contours (Fig. 7b). For this reason the long-shore component only needs to be determined at one cross-shore position in this experiment. Numerical oscillations were observed near the boundary, however they disappear a distance away where estimated values tend towards real values.

[FIGURE 7 OMITTED]

From the wave number components obtained at cross-shore distance x = 1000 m, a wave angle estimation of 27.58[degrees] compares well with the real value of 27.52[degrees] at the same position. Once the cross-shore and long-shore components of the wave number and the peak frequency are known, the water depth can be estimated using (6).

Figure 8 shows water depth estimated from intermediate waters to a position just before the breaker zone along a selected transect. The space average difference between the true and the estimated depths was of [bar.D] = 0.027 m along the profile and its mean relative error [bar.R] = 1.2%. The error in this case results merely from the use of CPCA for wave number estimation, since the data were numerically generated with a linear model.

[FIGURE 8 OMITTED]

Video-based system data: Wave characteristics were the same as used by the numerically generated data case. An oblique video, simulating field conditions, is recorded, digitized and rectified for its analysis. Pixel intensity time series were collected at different cross and long-shore positions from the recorded video of a wave propagation simulation. The sampling frequency was set equal to 2 Hz for a 112-s record.

The wave period was obtained as the average power spectrum from pixel intensity time series collected at different positions (Fig. 9). The resulting period was of T = 8.04 sec which agrees with the real value of 8 sec used for feeding the model.

[FIGURE 9 OMITTED]

In addition, a pixel intensity analysis was done in the space domain. Wave number estimation from the CPCA at two different cross-shore positions was estimated. The proposal being to evaluate the performance of the video technique with increasing distance from the position of the camera Fig. 10. Wave number estimation along a transect close to the camera is expected to be more accurate with other estimates farther away. However, no significant difference was observed for small areas.

[FIGURE 10 OMITTED]

Finally, in Fig. 11 the real and the video-based estimation of the water depth are shown. The mean error along the transects was [bar.D] = 0.116 m, with a mean relative error of [bar.R] 2.3%. In this case, the errors on the video-based water depth estimation for the wave propagation simulation are associated with the image rectification process and with the use of only the first mode of the CPCA for wave number estimation. The difference between error estimation from the numerical model data and the video-based data contributes to the total error introduced from the video technique.

[FIGURE 11 OMITTED]

The use of the linear solution for the wave propagation resulted in an unrealistic decrease in depth estimation error for estimates closer to shore. This situation does not happen in reality, as discussed by Stockdon and Holman (7). In this study non-linear effects, wave current interaction, diffraction, reflection and wave breaking were omitted in the numerical model.

CONCLUSIONS

Although less accurate than using in-situ instruments in the field, the use of photogrammetry for depth estimation from wave shoaling characteristics is an approach that allows a non intrusive and inexpensive monitoring of coastal processes. It is feasible that a system based on a simple rectification algorithm can be developed for both qualitative and small-scale quantitative purposes.

A simple video-based system, using the depth estimation technique developed by Stockdon and Holman (7), was implemented and tested for depth estimation from a numerical wave propagation simulation. The technique is based on complex principal component analysis of pixel intensity time series for wave number estimation. The combined use of numerical modelling and video photogrammetry proposed here is efficient for the qualitative and quantitative evaluation of video-based systems when controlled experiments are desired and/or laboratory facilities such as a wave basin are not available. This approach to test depth inversion algorithms is much less expensive than the use of physical models and provides completely controlled conditions and avoids laboratory problems of light refraction and reflection. The limitations would depend on the complexity of the numerical model chosen to perform the simulation.

Analysis of the performance of the system shows that the depth estimation error related to the wave number estimation technique is of the same order of magnitude as that related with the rectification technique ([approximately equal to]1%).

The technique described in this paper can be easily used for studying spatial and temporal bathymetric changes and the results obtained inexpensively with positive economic and environmental implications.

ACKNOWLEDGMENTS

The first author gratefully acknowledges the financial support provided by the Mexican National Council of Science and Technology (CONACyT) during his master studies at Florida Institute of Technology under the scholarship number 168776. The authors also acknowledge the Coastal and Ocean Research Group of the University of Cantabria for providing video images and ground control point locations from El Puntal Spit (Spain).

REFERENCES

(1.) Grilli, S.T., 1998. Depth Inversion in Shallow Water based on Nonlinear Properties of Shoaling Periodic Waves. Coastal Eng., 35: 185-209.

(2.) Dean, R.G. and R.A. Dalrymple, 1994. Water Wave Mechanics: Hydrodynamics. World Sci.

(3.) Wolf, P.R., 1983. Elements of Photogrammetry, 2nd Edn., New York: McGraw-Hill.

(4.) Holland, K.T., R.A. Holman, T.C. Lippmann, J. Stanley and N. Plant, 1997. Practical use of video imagery in nearshore oceanographic field studies. IEEE J. Oc. Eng., 22: 81-92.

(5.) Herbers, T., S. Elgar and R. Guza, 1995. Generation and Propagation of Infragravity Waves. J. Geophys. Res., 100, C12, 24, 863-24, 872.

(6.) Merrifield, M.A. and R.T. Guza, 1990. Detecting propagating signals with complex empirical orthogonal functions: A Cautionary Note. J. Phys. Oceanogr., 20: 1628-1633.

(7.) Stockdon, H.F. and R.A. Holman, 2000. Estimation of wave phase speed and nearshore bathymetry from video imagery. J. Geophys. Res., 105: 22015-22033.

(8.) Wallace, J.M. and R.E. Dickinson, 1972. Empirical orthogonal representation of time series in the frequency domain, Part I: Theoretical Considerations. J. Applied Meteorol., 11: 887-892.

(9.) Horel, J.D., 1984. Complex Principal Component Analysis: Theory and Examples. J. Climate Applic. Meteorol., 23: 1660-1673.

(10.) Holland, K.T., 1998. Validation of the Linear Dispersion Relation Using Field Observations. Naval Research Laboratory, Mapping, Charting and Geodesy Branch Marine Geosciences Division.

(11.) Lippmann, T.C. and R.A. Holman, 1989. Quantification of sand bar morphology: a video technique based on wave dissipation. J. Geophys. Res., 94: 995-1011.

(12.) Lippmann, T.C. and R.A. Holman, 1991. Phase and angle of breaking waves measured with video techniques, Coastal Sediments, '91, Kraus, N. (Ed.), Am. Soc. Civ. Eng., New York, pp: 542-556.

(13.) Curtis, W.R., K.K. Hathaway, K.T. Holland and W.C. Seabergh, 2002. Video-based wave direction measurements in a scale physical model, Coastal and Hydraulics Engineering Technical Note ERDC/CHL CHETN IV-49, U.S. Army Engineer Research and Development Center, Vicksburg, MS.

(14.) Dalrymple, R.A., 1988. Model for Refraction of Water Waves, J. Wtrwy., Port, Coast. Oc. Eng.., ASCE, 114: 423-435.

(1) Alec Torres and (2) Rodolfo Silva

(1) Departamento de Ciencias y Tecnicas del Agua y del Medio Ambiente, University of Cantabria, ETSI de Caminos, Canales y Puertos, Avda. de los Castros s/n, Santander 39005, Spain

(2) Instituto de Ingenieria, Universidad Nacional Autonoma de Mexico, Cd Universitaria, 04510 DF, Mexico

Corresponding Author: Alec Torres, Departamento de Ciencias y Tecnicas del Agua y del Medio Ambiente, University of Cantabria, ETSI de Caminos, Canales y Puertos, Avda. de los Castros s/n, Santander 39005, Spain