On Evolution of Young Wind Waves in Time and Space

The mechanisms governing the evolution of the wind-wave field in time and in space are not yet fully understood. Various theoretical approaches have been offered to model wind-wave generation. To examine their validity, detailed and accurate experiments under controlled conditions have to be carried out. Since it is next to impossible to get the required control of the governing parameters and to accumulate detailed data in field experiments, laboratory studies are needed. Extensive previously unavailable results on the spatial and temporal variation of wind waves accumulated in our laboratory under a variety of wind-forcing conditions and using diverse measuring techniques are reviewed. The spatial characteristics of the wind-wave field were determined using stereo video imaging. The turbulent airflow above wind waves was investigated using an X-hot film. The wave field under steady wind forcing as well as evolving from rest under impulsive loading was studied. An extensive discussion of the various aspects of wind waves is presented from a single consistent viewpoint. The advantages of the stochastic approach suggested by Phillips over the deterministic theory of wind-wave generation introduced by Miles are demonstrated. Essential differences between the spatial and the temporal analyses of wind waves’ evolution are discussed, leading to examination of the applicability of possible approaches to wind-wave modeling.


Introduction
Although the process of the generation of sea waves by wind has fascinated the human mind since ancient times, attempts to understand the processes leading to the excitation of waves by wind started only about 100 years ago. Our understanding of physics governing the interaction of atmosphere and ocean in general, and, more specifically, of mechanisms leading to the excitation of waves on the water surface and their evolution and growth due to the action of wind remains incomplete. Extensive reviews of the current state of knowledge on the subject were presented in [1][2][3][4], and the effects related to its limited fetch and duration are emphasized in [5,6].
Many important questions related to the extremely complex process of energy and momentum exchange between wind and waves still do not have fully satisfactory answers [6,7]. Accurate and reliable measurements in the airflow over the waves and in the water below remain a very demanding task. Field measurements in oceans and lakes are advantageous, as the collected data provide information on the measured parameters at natural scales. However, unsteady and unpredictable environmental and wind-forcing conditions are just two of the many formidable difficulties that accompany field experiments and complicate obtaining reliable information. Even when the requirement of steadiness of wind forcing can be seen as effectively satisfied, field measurements are often performed at a single location, thus preventing studies of wave evolution with fetch and limiting the possibilities for direct estimates of wind-wave growth rates. To mitigate some of the difficulties, Donelan et al. [8,9] carried out wind-wave measurements in relatively small lakes, while the applied in this work, as well as in the study of the effect of wind gusts on the developed wind-wave field [44]. Veron and Melville [45] applied a number of sensors to study waves driven by accelerating wind. Uz et al. [46,47] investigated the wave field in a tank under steady, abruptly starting and modulated wind speeds. A range of wavelengths was investigated using a scanning laser slope gauge. Hwang and Wang [12] and Hwang et al. [48] reported on field experiments under unsteady forcing and suggested appropriate space-time conversion to relate fetch and duration-limited results.
Wind velocity over the ocean varies both in magnitude and in direction. The variations occur in time and in space, and span scales that differ by orders of magnitude. The unsteadiness and inhomogeneity of wind waves in nature contributes to the significant scatter of results gained in field experiments. The theoretical analysis of the wind-wave generation poses numerous problems, even for steady forcing, due to a lack of adequate understanding of the mechanisms governing the energy and momentum transfer from wind to waves and of the relative importance of the details of those mechanisms as compared to the nonlinear interactions among wind waves [37,49]. It is often assumed in theoretical and numerical studies that the wave field is spatially homogeneous, so that under steady forcing, waves evolve in time only. The accumulation of quantitative data on the spatial as well as on the temporal variability of wind waves, in particular under unsteady forcing, is thus crucial for the improvement of wave forecasting models under realistic conditions.
To obtain such experimental data, a relatively small wind-wave facility has been constructed in our laboratory that enables measurements of waves under controlled conditions for steady or unsteady wind forcing. The short time and length scales associated with wind waves in a moderately scaled facility enable the accumulation of statistically reliable parameters characterizing a random wind-wave field in experiments of reasonable duration. However, it should be stressed that while experiments performed under controlled conditions offer clear advantages over field measurements, waves excited in any laboratory tank apparently are very different from those in the ocean. Wavelengths even in the largest available facilities are shorter by orders of magnitude than the lengths of typical sea waves. Ocean waves exhibit considerable angular spreading, whereas in wind-wave tanks, the deviation of the waves from the downwind direction, as well as the possible effects of sidewalls, often do not attract sufficient attention. Young wind waves in the laboratory are usually steeper than those in the ocean, resulting in stronger nonlinearity. The frequency of relatively short waves in the laboratory may be affected by wind-induced shear current; thus, nonlinear effects such as the Stokes drift current and frequency shift are more essential than for longer waves in the ocean. Wind waves in the laboratory may contain considerable energy at frequencies and wavelengths where the effects of capillarity become essential. Due to these and additional differences between water waves in the laboratory and in nature, the relevance of conclusions based on small-scale experiments to full-scale phenomena has to be addressed.
To this end, the laboratory results should be compared with the data obtained in field experiments using appropriate scaling. Kitaigorodskii [50] was among the first to suggest the application of dimensional analysis to study the spatially evolving wind-wave field. Following his approach, the principal parameters, such as the peak wave frequency f p and the characteristic wave energy m 0 = (η 2 )-in which η(x, t) is the instantaneous surface elevation at the measurement location x along the wind direction (the fetch), and the bar denotes time average-depend on the dimensionless fetchx = gx/U 2 . Here, g is the acceleration due to gravity and U is the characteristic wind velocity. For historical reasons, wind velocity at a certain elevation z above the water surface, usually U 10 = U (z = 10 m), was taken as the reference value. However, the momentum exchange between the wind and water surface is best represented by the friction velocity u * = (τ s /ρ a ) 1/2 , in which τ s is the mean shear stress at air-water interface, and ρ a is the air density. The dependencies of the peak frequency f p and wave energy m 0 on fetch x can be written in the following dimensionless form: Atmosphere 2019, 10, 562 4 of 54 g 2 m 0 /u 4 * = f 2 (gx/u 2 * ), (2) Similarly, the dependence of the power frequency spectrum F(f ) on fetch can be presented as: Demonstration of the applicability of conclusions based on results obtained in laboratory experiments to much larger scales require the application of suitable scaling and understanding the limitations imposed by the facility size.
This paper summarizes and discusses the results of an extensive research project carried out in our laboratory. The project dealt with detailed measurements of wind-wave evolution in a laboratory tank under diverse operational conditions. Both steady and impulsive wind-forcing cases were considered. The results of those measurements appeared in numerous publications. In this review, they are systematically presented in a consistent and organized form from a coherent viewpoint that emerged from the understandings gained over the course of this project. An extensive discussion of those understandings follows the presentation of the accumulated experimental results. In particular, the essential differences between the spatial and the temporal evolution of wind waves are specified, with consequences associated with the applicability of Miles' and Phillips approaches to the determination of initial wind wave growth in time and in space. The significance of those understandings to the modeling of wind waves is also discussed.
The paper is structured as follows. The experimental facility and the data-processing procedures are described in Section 2. The results of measurements of waves excited by steady wind are presented in Section 3. Measurements of the parameters characterizing steady turbulent airflow above water waves and the momentum exchange between wind and waves are described, including coupling between the turbulent airflow above the waves and wind waves. Then, the main findings on the spatial evolution of the random stationary wave field and its statistical parameters are presented. Particular attention is given to the three-dimensional and random nature of wind waves. Section 4 is devoted to experimental results on the spatial and temporal variation of waves under effectively impulsive or nearly impulsive wind forcing. Then, the total body of the accumulated experimental results is discussed in Section 5 in view of understandings gained over the course of the whole project. The aspects considered here include the relevance of the obtained laboratory results to larger scales in the ocean, the applicability of various theoretical approaches to wind-wave modeling, and the limitations of the prevailing theoretical approaches. The main conclusions are briefly summarized in Section 6.

Experimental Facility and Techniques
The wind-wave flume in Water Waves Laboratory, Tel-Aviv University (TAU) described in detail in Liberzon and Shemer [51] consists of a closed-loop wind tunnel installed atop of a test section that is 5 m long, 0.4 m wide, and 0.5 m high (see Figure 1). It was designed and constructed with the purpose to enable measurements at any desired location along the test section by a variety of sensors and experimental methods. A computer-controlled blower enables the maximum wind speed in the test section, which may exceed 15 m/s. The rectangular air inlet and outlet openings in the tank are 0.4 m wide and 0.25 m high. The wind tunnel is equipped with 1 m 3 settling chambers at the inlet and the exit of the test section; a 5-cm thick honeycomb with 5-mm hexagon cells and a converging nozzle with an area reduction ratio of about 4 at the inlet assure uniform airflow in the test section. Sidewalls and the bottom of the wave tank are made of clear glass to enable flow visualization of the wave field from all directions. Transparent removable Perspex plates with a partially sealed slot along the centerline of the test section to facilitate inserting the sensors serve as the tank's roof. An instrument carriage built of aluminum profiles that is movable along the test section can be positioned at any fetch. When necessary for imaging purposes, the Perspex plates at the imaged fetch are replaced by a clear glass sheet to minimize optical disturbances. The water depth in the test section is maintained at about separate vertical precision stages with the positioning accuracy of 0.05 mm. The stages are driven by a pair of identical step motors and controlled by a PC. The shorter 15-cm stage is used for the positioning and calibration of the wave gauges, while airflow sensors are mounted on the 40-cm long one. A Pitot tube with an outer diameter of 1 mm connected to a sensitive pressure transducer with a resolution of 2.5 . 10 -5 Pa was used to measure the mean airflow velocity profile, while a different 3mm Pitot tube was used to monitor the airflow during the wind-wave experiments. Static pressure fluctuations were measured using a special probe, as described in the research of Liberzon and Shemer [52]. Capacitance-type wave gauges made of 0.3-mm or 0.5-mm anodized tantalum wires were used to measure the instantaneous surface elevation η. To measure the vertical profiles of the turbulent air velocity fluctuations in the mean flow, u', and the vertical, w', directions, as well as the mean and/or ensemble-averaged Reynolds stress , an X-hot film probe was used; for details, see [53]. To eliminate wetting of the air sensors, a maximum wave height detector consisting of two adjacent bare wires located at fixed vertical displacement below the sensitive air sensors was used to estimate the maximum possible crest height at each location and wind velocity by an iterative sequence of measurements, as described in [51]. These data, in conjunction with the wave gauge records, enable the determination of the vertical coordinates of the air sensors relative to the mean water surface level at each fetch and wind velocity. Maintaining constant temperature in the test section is essential for velocity measurements in turbulent airflow using thermal anemometry. Therefore, a heat exchanger connected to an external water chiller and controlled by a temperature controller was installed and set to keep constant air temperature independent of wind velocity. Wave gauges provide information on the temporal local variation of the surface elevation, but not on the shape of the water surface. One approach to gain insight into the 2D structure of the windwave field is based on measurements of the temporal variation of the instantaneous surface elevation slope in the direction along the tank and the spanwise direction, ∂η/∂x (t) and ∂η/∂y (t). This is done by applying the laser slope gauge (LSG) [54][55][56][57]. The LSG is a non-intrusive instrument with a very fast response; it was used in numerous laboratory studies of water waves [58][59][60][61].
The instrument is based on the deflection by the inclined water surface of the laser beam directed vertically from below the transparent tank bottom. The laser beam passing through the water is refracted by the surface waves and then collimated by the Fresnel lens, forming an image on the 25 by 25 cm 2 diffusive screen in its focal plane. In our facility, the LSG is placed on a separate frame that can be mounted at any location along the test section. The instrument consists of four major components: a laser diode, a Fresnel lens, a diffusive glass screen and a position sensing detector (PSD). The estimated smallest wavelength that can be resolved by the instrument did not exceed about 2 mm. For the adopted LSG configuration, the coordinates of the beam image on the screen are the function of the instantaneous water slope only [62]. More details on the LSG assembly in our experimental facility, its calibration, and its application for wind-waves studies are presented in [63].
Stereo video imaging is a different approach that is capable of providing the temporal variation of the instantaneous water surface elevation over the imaged area. The images were recorded by two synchronized 4-Mpixel video cameras capable of recording at 180 frames per second (fps). The data Wave gauges provide information on the temporal local variation of the surface elevation, but not on the shape of the water surface. One approach to gain insight into the 2D structure of the wind-wave field is based on measurements of the temporal variation of the instantaneous surface elevation slope in the direction along the tank and the spanwise direction, ∂η/∂x (t) and ∂η/∂y (t). This is done by applying the laser slope gauge (LSG) [54][55][56][57]. The LSG is a non-intrusive instrument with a very fast response; it was used in numerous laboratory studies of water waves [58][59][60][61].
The instrument is based on the deflection by the inclined water surface of the laser beam directed vertically from below the transparent tank bottom. The laser beam passing through the water is refracted by the surface waves and then collimated by the Fresnel lens, forming an image on the 25 by 25 cm 2 diffusive screen in its focal plane. In our facility, the LSG is placed on a separate frame that can be mounted at any location along the test section. The instrument consists of four major components: a laser diode, a Fresnel lens, a diffusive glass screen and a position sensing detector (PSD). The estimated smallest wavelength that can be resolved by the instrument did not exceed about 2 mm. For the adopted LSG configuration, the coordinates of the beam image on the screen are the function of the instantaneous water slope only [62]. More details on the LSG assembly in our experimental facility, its calibration, and its application for wind-waves studies are presented in [63].
Stereo video imaging is a different approach that is capable of providing the temporal variation of the instantaneous water surface elevation over the imaged area. The images were recorded by two synchronized 4-Mpixel video cameras capable of recording at 180 frames per second (fps). The data processing and the reconstruction of the instantaneous surface elevation were carried out using methods described in detail in [16,64,65]. In those studies, stereo video imaging was used to investigate the sea surface; application of this technique to wind waves in a laboratory facility has certain advantages, as it allows controlled forcing and extensive (in terms of the dominant wavelengths) spatial coverage. The application of stereo video imaging in a closed laboratory requires adequate illumination. A uniform pulsed illumination of the wind-wave field was provided by two 14" by 14" white light-emitting diode (LED) panels. Two LED panels were used. The images span about 15 cm across the tank around its centerline; they cover about 40 cm, or a few dominant wavelengths, along the tank, with the exact number depending on the fetch and wind speed. This coverage allows studying the spatial inhomogeneity of the wind-wave field. However, reconstruction of the instantaneous water surface shape from stereo imaging is extremely computer-time consuming, somewhat limiting the number of imaged pairs of frames that can be processed for each wind-forcing condition and fetch in a reasonable time. More details on the stereo imaging technique applied in our experiments are given in [66].
The current induced by wind and waves at the water surface was measured using the particle tracking velocimetry (PTV) technique. Black 6-mm diameter thin paper disks were used as tracers and spread over the water surface; a video camera operating at 60 fps with resolution of 1280 by 720 pixels located above the water surface recorded the position of the floaters that move with the water surface. The camera imaged a water surface area of about 15 cm in the wind direction and about 10 cm across the test section. The recorded video clips allowed monitoring each particle's trajectory and calculating each particle's velocity components along and across the imaged section. The particles' positions were collected from each frame, and consequently, two velocity components-along the tank, U s and across the tank, V s -were tracked across all frames of appearance. Numerous measurements were performed to obtain statistically reliable mean values for different operational conditions [51,67].
All experiments were controlled by a LabView program that enables running the experiment, including setting the wind speed in the tunnel, the vertical positioning of the air sensors, the calibration of the sensors before and after data acquisition sessions, and the data acquisition, in a fully automatic mode without human intervention. More information about the experimental facility, calibration of the sensors, and of the detailed experimental and data-processing procedures is provided in [51,63,66,67].

Characterization of the Air Flow above Wind Waves: Mean Velocity Profiles and Vertical Distributions of Reynolds Stresses
The interfacial shear stress τ s and thus the friction velocity u * are usually determined in field as well as in laboratory experiments by performing measurements in the airflow. An accurate measurement of u * is of crucial importance for the comparison of data obtained under various experimental conditions and for the evaluation of predictions of the theoretical models. The structure of the mean flow above the air-water interface is analogous to that of the turbulent flow above a stationary rough surface so that the logarithmic turbulent velocity profile U(z) outside the viscous sublayer is usually represented as: where κ = 0.41 is the von Karman constant and z 0 is the characteristic roughness. The drift velocity induced by the wind U s at the local mean water surface z = 0 should be accounted for and subtracted Atmosphere 2019, 10, 562 7 of 54 from the measured mean wind velocity U(z) in Equation (4) [33,51]. Both the friction velocity u * and the roughness z 0 can be determined by fitting the measured velocity profile to Equation (4). Alternatively, the friction velocity can be estimated by application of the eddy correlation method [4]. To this end, instantaneous turbulent velocity fluctuations u and w in the horizontal (with wind) and vertical directions are measured by an X-hot-film probe. A profile of the turbulent Reynolds shear stress −ρu w (z) extrapolated to the mean water surface elevation defines the friction velocity . Values of the friction velocity in air u * determined by thermal anemometry were reported by [33,68,69], among others. In [53], measurements of instantaneous surface elevation and mean air velocity turbulent velocity fluctuation profiles were performed at seven fetches (1.0 m ≤ x ≤ 3.4 m) and at six different blower settings corresponding to the maximum wind velocities U max in the test section ranging from 5.5 to 11.2 m/s. No significant variation of U max was detected; deviations from the mean value along the test section did not exceed a few percentage points. Therefore, the average for all the fetch values of U max were used as the reference velocity. At each fetch and wind speed, the data acquisition was initiated only after the stationary wind and wave conditions, as well as the constant temperature, were attained.
To get accurate values of the friction velocity, the density of the measuring locations was higher closer to the air-water interface. At each measuring location and airflow rate, the sampling of X-hot-film data at the rate of 120 Hz continued for at least 5 minutes (min) to accumulate a sufficient amount of data for statistical analysis. Sensitivity to the sampling frequency was checked by performing experiments where the hot-film data were sampled at the rate of 1200 Hz/channel. No significant effect of the sampling frequency on the statistical turbulent flow parameters was detected. At each blower setting, the hot-film calibration was repeated after the completion of measurements at all elevations. In the Pitot-tube measurements of the vertical mean air velocity profiles, the vertical positions ranged from 1 mm to 130 mm relative to the highest wave crest, with the sampling duration at each point being 1 min.
The mean velocity distributions for two representative fetches (x = 260 cm and x = 240 cm) are presented in Figure 2. Additional experimental data on the vertical air velocity profiles can be found in [51,53]. Coles [70] showed that for solid stationary surfaces, starting from z + = zu * /ν ≥ 50 (approximately about 20% of the turbulent boundary layer thickness), the profile is no longer of logarithmic form, and is affected by the flow in the outer ('wake') layer. Thus, the boundary layer thickness over the mean water surface δ (taken here as the elevation where the mean velocity attains 99% of its maximum value) was chosen in Figure 2 as a parameter to normalize the elevation over the mean water surface. The values of δ range from about 50 mm at shorter fetches to about 80-90 mm at x = 340 cm; the variation of δ with wind velocity at each fetch is insignificant. The advantage offered by the automatic experimental procedure made possible air velocity measurements at about 100 vertical positions. For lower wind velocities and thus lower wave crests, measurements were performed quite close to the water surface. Since waves grow with the fetch x and with the wind velocity U max , the lowest measuring location had to be shifted upward (relative to the mean water level) accordingly. It is possible to distinguish in the profiles plotted in Figure 2 between two regions. The lower part of each profile, up to about 0.4δ, i.e., for maximum elevations not exceeding about 15 mm for lower fetches to about 30 mm at the far end of the test section, is linearly fitted and exhibits logarithmic velocity distribution, while above this region, deviation from the fit can be noticed [70,71]. When the log fit is performed on the whole set of data that extends to the outer part of the boundary layer, the resulting slope of the line is somewhat different, and the values of u * and z 0 resulting from such a fit are less accurate. The logarithmic fit curves for a given U max at panels in Figure 2 have similar slopes. The same observation holds for additional fetches, as shown in [53]. In view of Equation (4), this indicates that the friction velocity u * does not vary notably with fetch. The vertical shift of the profiles with x and with U max implies that the effective surface roughness z 0 in the present experiments generally increases with fetch and with the wind velocity. The vertical profiles of the root mean square (RMS) values of the deviations of the instantaneous velocity components from the local mean values in both longitudinal (u ) and vertical (w ) directions measured with an X-type hot-film anemometer are presented in Figure 3. The effect of the wind velocity in the test section represented by U max is studied in panels (a) and (c) for a fixed fetch x = 300 cm. In Figure 3b,d, the profiles of the normalized by U max RMS values of two components of the turbulent fluctuations are plotted for various fetches, and a constant value of U max = 7.7 m/s are plotted. Those results are related to the classical measurements of turbulent flow over a solid surface in Figure 3a,c. Here, the corresponding vertical profiles of turbulent velocity fluctuations measured in turbulent boundary layer over a smooth plate by Klebanoff [72] and over a rough plate by Corrsin and Kistler [73] (see also Hinze [74]) are presented as well; those measurements were performed at comparable free stream air velocities.
blower settings corresponding to the maximum wind velocities Umax in the test section ranging from 5.5 to 11.2 m/s. No significant variation of Umax was detected; deviations from the mean value along the test section did not exceed a few percentage points. Therefore, the average for all the fetch values of Umax were used as the reference velocity. At each fetch and wind speed, the data acquisition was initiated only after the stationary wind and wave conditions, as well as the constant temperature, were attained.
To get accurate values of the friction velocity, the density of the measuring locations was higher closer to the air-water interface. At each measuring location and airflow rate, the sampling of X-hotfilm data at the rate of 120 Hz continued for at least 5 minutes (min) to accumulate a sufficient amount of data for statistical analysis. Sensitivity to the sampling frequency was checked by performing experiments where the hot-film data were sampled at the rate of 1200 Hz/channel. No significant effect of the sampling frequency on the statistical turbulent flow parameters was detected. At each blower setting, the hot-film calibration was repeated after the completion of measurements at all elevations. In the Pitot-tube measurements of the vertical mean air velocity profiles, the vertical positions ranged from 1 mm to 130 mm relative to the highest wave crest, with the sampling duration at each point being 1 min. The mean velocity distributions for two representative fetches (x = 260 cm and x = 240 cm) are presented in Figure 2. Additional experimental data on the vertical air velocity profiles can be found in [51,53]. Coles [70] showed that for solid stationary surfaces, starting from (approximately about 20% of the turbulent boundary layer thickness), the profile is no longer of logarithmic form, and is affected by the flow in the outer ('wake') layer. Thus, the boundary layer thickness over the mean water surface δ (taken here as the elevation where the mean velocity attains 99% of its maximum value) was chosen in Figure 2 as a parameter to normalize the elevation over the mean water surface. The values of δ range from about 50 mm at shorter fetches to about 80-90 mm at x = 340 cm; the variation of δ with wind velocity at each fetch is insignificant. The advantage offered by the automatic experimental procedure made possible air velocity measurements at about 100 vertical positions. For lower wind velocities and thus lower wave crests, measurements were performed quite close to the water surface. Since waves grow with the fetch x and with the wind velocity Umax, the lowest measuring location had to be shifted upward (relative to the mean water As in Figure 2, the distance z from the mean surface level is normalized by the boundary layer thickness δ. The values of the normalized turbulent velocity fluctuations u and w' in Figure 3 exhibit qualitative and to some extent quantitative similarity, with those in the turbulent boundary layer over a solid surface, although they are generally somewhat higher even than those in the flow over a rough rigid plate. The difference can be mainly attributed to the higher turbulent level in our wind tunnel, as can be seen from comparison of the turbulent intensities outside the boundary layer. Due to the limitations imposed by water waves, measurements were mostly made at somewhat higher elevations. Note also that in closer vicinity to the interface, the RMS values of the velocity fluctuations as z tends to zero over a solid surface vanish, while over the waves, they remain finite due to the orbital wave motion at the air-water interface.
Atmosphere 2019, 10, x FOR PEER REVIEW 9 of 52  [72] in a turbulent boundary layer over a flat plate, and of Corssin and Kistler [73] for longitudinal velocity fluctuations over rough and smooth plates are plotted as well.
As in Figure 2, the distance z from the mean surface level is normalized by the boundary layer thickness δ. The values of the normalized turbulent velocity fluctuations u' and w' in Figure 3 exhibit qualitative and to some extent quantitative similarity, with those in the turbulent boundary layer over a solid surface, although they are generally somewhat higher even than those in the flow over a rough rigid plate. The difference can be mainly attributed to the higher turbulent level in our wind tunnel, as can be seen from comparison of the turbulent intensities outside the boundary layer. Due to the limitations imposed by water waves, measurements were mostly made at somewhat higher elevations. Note also that in closer vicinity to the interface, the RMS values of the velocity fluctuations as z tends to zero over a solid surface vanish, while over the waves, they remain finite due to the orbital wave motion at the air-water interface.  [72] in a turbulent boundary layer over a flat plate, and of Corssin and Kistler [73] for longitudinal velocity fluctuations over rough and smooth plates are plotted as well.
Measurements of the Reynolds shear stresses u w allow the direct evaluation of the friction velocity u * (the so-called eddy correlation method). The measured variation of the Reynolds shear stress u w with height is presented in Figure 4 for four wind velocities U max ; additional data can be found in [53]. The vertical profiles of the Reynolds stress are affected by the mean pressure gradient [52]. The pressure drop over the test section is quite modest, ranging from a few Pascals for low flow rates to about 20 Pa for U max ≈ 9 m/s. The Reynolds shear stresses within the boundary layer vary linearly with elevation for all wind flow rates and fetches, vanishing in the central part of the test section, and changing sign with approach to the roof of the channel. The data points within the boundary layer were linearly fitted, with fits extrapolated to the mean water level z = 0. The limiting value of the Reynolds shear stress at z = 0 can be seen as the squared value of the friction velocity u * (the so-called eddy correlation method). The very limited scatter of the limiting values of u w extrapolated to z = 0 in each frame of Figure 4 indicates that the values of the friction velocity measured by the eddy correlation method in the present experiments are practically independent of fetch. However, the interfacial shear stress increases notably with the wind velocity.

Measurements of the Reynolds shear stresses
uw  allow the direct evaluation of the friction velocity * u (the so-called eddy correlation method). The measured variation of the Reynolds shear stress uw  with height is presented in Figure 4 for four wind velocities Umax; additional data can be found in [53]. The vertical profiles of the Reynolds stress are affected by the mean pressure gradient [52]. The pressure drop over the test section is quite modest, ranging from a few Pascals for low flow rates to about 20 Pa for Umax ≈ 9 m/s. The Reynolds shear stresses within the boundary layer vary linearly with elevation for all wind flow rates and fetches, vanishing in the central part of the test section, and changing sign with approach to the roof of the channel. The data points within the boundary layer were linearly fitted, with fits extrapolated to the mean water level z = 0. The limiting value of the Reynolds shear stress at z = 0 can be seen as the squared value of the friction velocity * u (the so-called eddy correlation method). The very limited scatter of the limiting values of uw  extrapolated to z = 0 in each frame of Figure 4 indicates that the values of the friction velocity measured by the eddy correlation method in the present experiments are practically independent of fetch. However, the interfacial shear stress increases notably with the wind velocity.
The eddy correlation method for the determination of the shear stress is considered more reliable The eddy correlation method for the determination of the shear stress is considered more reliable as compared to determination of u * from the velocity profiles, as presented in Figure 2 and in the application of Equation (4), since the Reynolds stress is measured directly. However, thermal anemometry is more complicated and time consuming; thus, the measurements in [53] were carried out for a limited set of parameters. To validate the accuracy of both approaches, the friction velocities estimated from the logarithmic fit of mean velocity profiles and the values derived by application of the eddy correlation method are compared in Figure 5a for all fetches and wind flow rates.
Atmosphere 2019, 10, x FOR PEER REVIEW 11 of 52 estimated from the logarithmic fit of mean velocity profiles and the values derived by application of the eddy correlation method are compared in Figure 5a for all fetches and wind flow rates. It is obvious from Figure 5a that these two independent methods yield very similar results, with the scatter within about 5-10% and data points spread at both sides of the 45° line. Thus, it is demonstrated that estimating the values of * u from the logarithmic fit of the mean wind velocity profile that contains a sufficiently large number of data points close to the air-waterinterface can serve as a reliable method for the determination of the friction velocity. However, it should be stressed that measurements should be restricted to elevation satisfying the condition z/δ < 0.4.
As can be seen in Figure 5b, the values of * u do not depend notably on the location along the test section (with a possible exception of very short fetches), and are roughly proportional to Umax. The fit in this figure corresponds to the relation * u = 0.73Umax. These results are in agreement with [51]. An attempt has been made to examine the effect of the surface drift velocity, Us, on the values of * u and z0, by substituting U(z) -Us instead of the mean velocity U(z) in Figure 1(a). The measurements of Us were performed in [51] using particle tracking velocimetry and demonstrated that the drift velocity constitutes about 3% of Umax. However, this modification had only a minor effect on * u and z0 that did not cause deviation exceeding 2-3% for * u , and practically did not affect the effective roughness. The values of the effective water surface roughness, z0, that were also retrieved from those fitted velocity profiles as shown in Figure 2, as well as more details on the measurements of friction velocity are presented in [53].

Statistical Parameters of Waves under Steady Wind Forcing
Measurements of the principal statistical parameters of wind waves were carried out in Zavadsky et al. [75] at eight fetches (100 cm ≤ x ≤ 380 cm). The blower operated at nine airflow rates, resulting in the maximum wind velocities in the test section in the range 3.3 m/s ≤ Umax ≤ 12.3 m/s and friction velocities 0.17 m/s ≤ * u ≤ 1.28 m/s. At each fetch and wind velocity Umax, the instantaneous variation of the surface elevation was continuously recorded for 90 min at the sampling rate 120 Hz, thus enabling the accumulation of very large ensembles of experimental data.
The power spectra of the surface elevation variation with time, η(t), were calculated by dividing each record into windows that contained 8192 points (duration of each window about 70 s), with 50% overlap, resulting in 94 segments for each fetch and airflow rate. The variation with fetch of the surface elevation power spectra averaged over all windows is presented in Figure 6a,b for two representative airflow rates. It is obvious from Figure 5a that these two independent methods yield very similar results, with the scatter within about 5-10% and data points spread at both sides of the 45 • line. Thus, it is demonstrated that estimating the values of u * from the logarithmic fit of the mean wind velocity profile that contains a sufficiently large number of data points close to the air-waterinterface can serve as a reliable method for the determination of the friction velocity. However, it should be stressed that measurements should be restricted to elevation satisfying the condition z/δ < 0.4.
As can be seen in Figure 5b, the values of u * do not depend notably on the location along the test section (with a possible exception of very short fetches), and are roughly proportional to U max . The fit in this figure corresponds to the relation u * = 0.73U max . These results are in agreement with [51]. An attempt has been made to examine the effect of the surface drift velocity, U s , on the values of u * and z 0 , by substituting U(z) -U s instead of the mean velocity U(z) in Figure 1a. The measurements of U s were performed in [51] using particle tracking velocimetry and demonstrated that the drift velocity constitutes about 3% of U max . However, this modification had only a minor effect on u * and z 0 that did not cause deviation exceeding 2-3% for u * , and practically did not affect the effective roughness. The values of the effective water surface roughness, z 0 , that were also retrieved from those fitted velocity profiles as shown in Figure 2, as well as more details on the measurements of friction velocity are presented in [53].

Statistical Parameters of Waves under Steady Wind Forcing
Measurements of the principal statistical parameters of wind waves were carried out in Zavadsky et al. [75] at eight fetches (100 cm ≤ x ≤ 380 cm). The blower operated at nine airflow rates, resulting in the maximum wind velocities in the test section in the range 3.3 m/s ≤ U max ≤ 12.3 m/s and friction velocities 0.17 m/s ≤ u * ≤ 1.28 m/s. At each fetch and wind velocity U max , the instantaneous variation of the surface elevation was continuously recorded for 90 min at the sampling rate 120 Hz, thus enabling the accumulation of very large ensembles of experimental data.
The power spectra of the surface elevation variation with time, η(t), were calculated by dividing each record into windows that contained 8192 points (duration of each window about 70 s), with 50% overlap, resulting in 94 segments for each fetch and airflow rate. The variation with fetch of the surface elevation power spectra averaged over all windows is presented in Figure 6a,b for two representative airflow rates. It is obvious from these panels in Figure 6 that at each flow rate, the total wave energy increases along the channel, while the peak frequency decreases with fetch; therefore, both the horizontal and vertical axes of the two panels of Figure 6 are scaled differently. A comparison of the panels in Figures 6a and 6b also demonstrates downshifting of the peak frequency at a given fetch with increasing wind speed. These results are in agreement with [51]. The characteristic frequency of the wind-wave field can be defined either as the frequency at the spectral peak, fp, or as the dominant frequency, fdom, based on the spectral moments: ( ) .
As mentioned in the Introduction, the zero-th moment m0 defines the characteristic wave amplitude Note that in Equation (5), the integration is carried out over the free waves' domain only; this domain for each spectrum, ωmin < ω < ωmax, is taken within ±60% of the peak frequency.
Since fdom is an integral value and thus less prone to variations in the experimentally estimated wave power spectra; it is a more robust quantity at each fetch and airflow rate. However, the frequency of the spectral peak fp is a physical and visual parameter that is more intuitively straightforward than fdom, and therefore is usually chosen as the characteristic of the wave field. To mitigate the scatter in the power spectra as seen in Figure 6, the value of fp is defined as the frequency of the peak of the parabolic fit performed within about ±0.3 Hz around the maximum values in the It is obvious from these panels in Figure 6 that at each flow rate, the total wave energy increases along the channel, while the peak frequency decreases with fetch; therefore, both the horizontal and vertical axes of the two panels of Figure 6 are scaled differently. A comparison of the panels in Figure 6a,b also demonstrates downshifting of the peak frequency at a given fetch with increasing wind speed. These results are in agreement with [51]. The characteristic frequency of the wind-wave field can be defined either as the frequency at the spectral peak, f p , or as the dominant frequency, f dom , based on the spectral moments: As mentioned in the Introduction, the zero-th moment m 0 defines the characteristic wave amplitude Note that in Equation (5), the integration is carried out over the free waves' domain only; this domain for each spectrum, ω min < ω < ω max , is taken within ±60% of the peak frequency.
Since f dom is an integral value and thus less prone to variations in the experimentally estimated wave power spectra; it is a more robust quantity at each fetch and airflow rate. However, the frequency of the spectral peak f p is a physical and visual parameter that is more intuitively straightforward than f dom , and therefore is usually chosen as the characteristic of the wave field. To mitigate the scatter in the power spectra as seen in Figure 6, the value of f p is defined as the frequency of the peak of the parabolic fit performed within about ±0.3 Hz around the maximum values in the experimentally determined spectra. The comparison of the values of f dom and f p carried out in [75] demonstrated that although the difference between these two parameters is small, f dom is consistently larger than f p due to the greater contribution of higher frequencies in the computation of the spectral moments m 1 , as shown in Equation (5), with the slope of the linear fit of 0.96.
The surface elevation power spectra are now examined in greater detail in Figure 6c,d. In Figure 6c, the spectra obtained at several selected fetches at a constant air flow rate are plotted in a logarithmic scale as a function of the dimensional frequency, whereas in Figure 6d, the power spectra are normalized by their peak values and plotted for several wind velocities as a function of f/f p . The lines corresponding to the power law: with n = 4 are also plotted in Figure 6c,d.
The spectra in the panels c and d of Figure 6 are characterized by peaks at the second and the third harmonics of the peak frequency f p ; at low wind velocities, even weaker peaks corresponding to the fourth harmonic can be identified. Those peaks signify the contribution of nonlinear bound waves. At frequencies exceeding the third or the fourth harmonics of f p , the peaks associated with bound waves cannot be distinguished; those frequencies correspond to the so-called quasi-saturated high-frequency tail. The normalized power spectra of the surface elevation for a range of wind flow rates in Figure 6d seem to collapse on a single curve for f < 3f p .
The tail behavior is only in qualitative agreement with that corresponding to n = 4; a notable spread of the tail slope is clearly visible in Figure 6d. Therefore, a closer look at the power-law dependence on the frequency of the high-frequency 'tail' of the spectra is presented in Figure 7a. The linear fit range 3 < f/f p < 6 was chosen in this figure, as the spectra in this part are free of noticeable higher harmonic peaks, which can affect the slope. The power n that defines the exponent in the power law shown in Equation (6) is defined by the slope of the fit. The values of n are usually somewhat above n = 3 and remain mostly below n = 4; they decrease with wind velocity and have a tendency to increase with fetch. At the largest wind velocity, the measured value of n is the lowest, and is close to the theoretical prediction of n = 17/6, which equals approximately 2.83, for purely capillary waves [76], it increases with fetch from n = 3.1 for x = 100 cm to n = 3.74 for the longest fetch, x = 380 cm.  [75] demonstrated that although the difference between these two parameters is small, fdom is consistently larger than fp due to the greater contribution of higher frequencies in the computation of the spectral moments m1, as shown in Equation (5), with the slope of the linear fit of 0.96. The surface elevation power spectra are now examined in greater detail in Figure 6c,d. In Figure  6c, the spectra obtained at several selected fetches at a constant air flow rate are plotted in a logarithmic scale as a function of the dimensional frequency, whereas in Figure 6d, the power spectra are normalized by their peak values and plotted for several wind velocities as a function of f/fp. The lines corresponding to the power law: with n = 4 are also plotted in Figure 6c,d.
The spectra in the panels c and d of Figure 6 are characterized by peaks at the second and the third harmonics of the peak frequency fp; at low wind velocities, even weaker peaks corresponding to the fourth harmonic can be identified. Those peaks signify the contribution of nonlinear bound waves. At frequencies exceeding the third or the fourth harmonics of fp, the peaks associated with bound waves cannot be distinguished; those frequencies correspond to the so-called quasi-saturated high-frequency tail. The normalized power spectra of the surface elevation for a range of wind flow rates in Figure 6d seem to collapse on a single curve for f < 3fp.  The tail behavior is only in qualitative agreement with that corresponding to n = 4; a notable spread of the tail slope is clearly visible in Figure 6d. Therefore, a closer look at the power-law dependence on the frequency of the high-frequency 'tail' of the spectra is presented in Figure 7a. The linear fit range 3 < f/fp < 6 was chosen in this figure, as the spectra in this part are free of noticeable higher harmonic peaks, which can affect the slope. The power n that defines the exponent in the power law shown in Equation (6) is defined by the slope of the fit. The values of n are usually somewhat above n = 3 and remain mostly below n = 4; they decrease with wind velocity and have a tendency to increase with fetch. At the largest wind velocity, the measured value of n is the lowest, and is close to the theoretical prediction of n = 17/6, which equals approximately 2.83, for purely capillary waves [76], it increases with fetch from n = 3.1 for x = 100 cm to n = 3.74 for the longest fetch, x = 380 cm.
The measured characteristic wave parameters can be presented as relations among dimensionless groups. For steady wind forcing, those parameters depend on fetch only. The friction velocity * u is the quantity directly related to momentum transfer from air to water, and as stated in  The measured characteristic wave parameters can be presented as relations among dimensionless groups. For steady wind forcing, those parameters depend on fetch only. The friction velocity u * is the quantity directly related to momentum transfer from air to water, and as stated in Section 3.1, for a given airflow rate in the test section, it is essentially independent of fetch. The characteristic wave amplitude can be represented either by the RMS value of the surface elevation (η 2 ) 1/2 , or based on the zero spectrum momentum, m 1/2 0 . The definition based on the spectral momentum seems to be advantageous as it contains free wave components only, which is in agreement with most nonlinear wind-wave theories; see e.g., [49]. The comparison of (η 2 ) 1/2 and m 1/2 0 for all experimental conditions demonstrates that wave energy is mainly contained in the spectral range around the peak frequency; the values of m 1/2 0 are smaller than the corresponding RMS values of the surface elevation by less than about 3%.
In view of the relations shown in Equations (1) and (2), the characteristic amplitude and the peak frequency can be presented in dimensionless form asη = gm 1/2 0 /u 2 , where the dimensionless fetchx = gx/u 2 * . Note that for a given dimensional distance from the inlet x, the dimensionless fetchx decreases with u * , and thus with the increased airflow rate in the tank. The dependencies of the dimensionless amplitudesη and dimensionless peak frequenciesf p on the dimensionless fetchx are plotted in Figure 8a  m . The definition based on the spectral momentum seems to be advantageous as it contains free wave components only, which is in agreement with most nonlinear wind-wave theories; see e.g., [49]. The comparison of / gx u . Note that for a given dimensional distance from the inlet x, the dimensionless fetch x decreases with * u , and thus with the increased airflow rate in the tank.
The dependencies of the dimensionless amplitudes η and dimensionless peak frequencies ˆp f on the dimensionless fetch x are plotted in Figure 8a Figure 8a,b are in agreement with the wind-wave channel measurements at a comparable range of dimensionless fetches by Mitsuyasu [77] and Mitsuyasu and Honda [78], who suggested the following relations: Similar results were obtained in field measurements by the JONSWAP group [79] and by Kahma [80]. The points in Figure 6a,b that do not fit the general trend and are marked by filled symbols correspond to lower wind velocities and thus higher dimensionless fetches x . The peak frequencies at lower wind velocities mostly exceed 5 Hz at all fetches; the frequencies of the corresponding short waves are strongly affected by Doppler shift due to the surface drift current [51]. The effect of the drift current that is mostly pronounced at low airflow rates seems to cause the deviation from the power law of the energy growth and peak frequency variation along the test section; therefore, these values were not included in the data fits in Figure 6. The relation between the dimensionless wave amplitude η and peak frequency ˆp f can be obtained from the linear fits in Figure 6: The majority of points (denoted by open symbols) in panels a and b in Figure 8 fit a linear dependence in log-log coordinates, and thus exhibit power-law behavior. The dimensionless frequencŷ f p decreases with the dimensionless fetch asx −0.27 , while the dimensionless wave amplitudesη grow asx 0.506 . The scatter of data represented by open symbols around the fit line can be attributed in part to local variations in the friction velocities. The power dependencies presented in Figure 8a,b are in agreement with the wind-wave channel measurements at a comparable range of dimensionless fetches by Mitsuyasu [77] and Mitsuyasu and Honda [78], who suggested the following relations: Similar results were obtained in field measurements by the JONSWAP group [79] and by Kahma [80]. The points in Figure 6a,b that do not fit the general trend and are marked by filled symbols correspond to lower wind velocities and thus higher dimensionless fetchesx. The peak frequencies at lower wind velocities mostly exceed 5 Hz at all fetches; the frequencies of the corresponding short waves are strongly affected by Doppler shift due to the surface drift current [51]. The effect of the drift current that is mostly pronounced at low airflow rates seems to cause the deviation from the power law of the energy growth and peak frequency variation along the test section; therefore, these values were not included in the data fits in Figure 6. The relation between the dimensionless wave amplitudê η and peak frequencyf p can be obtained from the linear fits in Figure 6: This result differs somewhat from the generally accepted −3/2 power law H s = BT s −3/2 proposed by Toba [81] that is based on extensive empirical data; here, H s is the significant wave height, and T s is defined by the peak frequency.
The dimensionless spectral width ν defined by the spectral moments in Equation (5) as: is an integral characteristic of the shape of the free-wave part of the power spectrum at lower frequencies around f p . The values of ν calculated for all wind velocities in the present experiment are plotted in Figure 8c as a function of the dimensionless fetchx. For each U max , the spectral width decreases with x, the rate of decrease is higher at lower fetches and wind velocities. The values of the spectral width vary from about ν = 0.25 to ν = 0.145; they become nearly constant around ν = 0.17 for higher fetches and stronger winds. Thus, the process of the evolution of wind-wave spectra with fetch is characterized by a decrease in both the peak frequency f p and the spectral width ν.
The deviation of the surface elevation behavior in a random wind-wave field from Gaussianity is examined in Figure 9. The lack of exact Gaussianity can be characterized by coefficients of skewness, λ 3 , and kurtosis, λ 4 , which is defined as: for normally distributed random variables, λ 3 = 0 and λ 4 = 3. The values of λ 3 and λ 4 calculated from the measured data for all cases are plotted in Figure 9a,b.
Atmosphere 2019, 10, x FOR PEER REVIEW 15 of 52 This result differs somewhat from the generally accepted −3/2 power law Hs = BTs −3/2 proposed by Toba [81] that is based on extensive empirical data; here, Hs is the significant wave height, and Ts is defined by the peak frequency.
The dimensionless spectral width ν defined by the spectral moments in Equation (5) as: is an integral characteristic of the shape of the free-wave part of the power spectrum at lower frequencies around fp. The values of ν calculated for all wind velocities in the present experiment are plotted in Figure 8c as a function of the dimensionless fetchx . For each Umax, the spectral width decreases with x, the rate of decrease is higher at lower fetches and wind velocities. The values of the spectral width vary from about ν = 0.25 to ν = 0.145; they become nearly constant around ν = 0.17 for higher fetches and stronger winds. Thus, the process of the evolution of wind-wave spectra with fetch is characterized by a decrease in both the peak frequency fp and the spectral width ν. The deviation of the surface elevation behavior in a random wind-wave field from Gaussianity is examined in Figure 9. The lack of exact Gaussianity can be characterized by coefficients of skewness, λ3, and kurtosis, λ4, which is defined as:  Figure 9a. It characterizes the vertical wave asymmetry and is positive for all records except for very short fetches and the lowest wind velocities, where the wave heights are extremely small (usually less than 1 mm). The measured values of skewness tend to increase with fetch and with wind velocity, attaining values close to 0.5 at larger fetches and wind velocities. Positive values of λ 3 characterize an essentially nonlinear wave field where the crest heights exceed those of troughs due to the contribution of the second-order bound waves that cause an effective increase in wave crests and decrease in troughs [82]. The kurtosis coefficient characterizes the deviation of wave height probability distribution from Gaussian. The values of λ 4 plotted in Figure 9b are below the Gaussian limit of three for all wind velocities and at all fetches, except for the low wind velocity at moderate fetches, where it is somewhat higher. For wind velocities exceeding about 7 m/s, the kurtosis coefficient seems to be nearly constant at λ 4 ≈ 2.4-2.6 and essentially independent of both the fetch and the wind velocity. These values of λ 4 indicate a relative deficit of large-amplitude waves in the distribution. The values of kurtosis presented in Figure 9(b) are in qualitative agreement with the experimental results of Shemer et al. [82] obtained for mechanically generated unidirectional random waves at comparable values of the spectral width ν. However, they disagree with the results of wind-wave simulations by Annenkov and Shrira [83] based on the Zakharov equation, which obtained λ 4 > 3.
The relatively low probability of very steep waves in the present experiments, as suggested by the results of Figure 9b, is further validated by plotting the probability distribution of the wave height along the flume, as shown in Figure 9c. The results presented for the extreme values of the airflow rate in these experiments are also representative for all other wind velocities. For linear narrow-band Gaussian waves, the wave heights H follow the Rayleigh exceedance distribution [84]: Longuet-Higgins [85] also has shown that nonlinearities can cause the deviation of surface elevation distributions F(H) from the Gaussian statistics. Numerous laboratory experiments, in situ measurements, and numerical simulations indeed support this conclusion (see, e.g., [82,86,87] and the references therein).
The exceedance distributions obtained in the present experiments are compared in Figure 9c,d with the Rayleigh distribution. Each individual wave height was calculated as the difference between the maximum and minimum between two consecutive positive zero crossings. Note that the statistics in those panels is based on ensembles that contain more than 10 4 individual waves. For all experimental conditions, notable deviations from the Rayleigh distribution are indeed clearly visible in those panels.
For the Rayleigh distribution, the significant wave height H s ≈ 4 η 2 1/2 . The probability of small and moderate waves with heights H < H s , is higher in Figure 9c,d than that corresponding to the Rayleigh distribution. However, the probability of waves higher than H s is overestimated by Equation (10). These deviations from the Rayleigh distribution are stronger for the higher wind velocity in Figure 9d as compared to the weaker wind case in Figure 9c. No dependence of the measured F(H) on fetch can be clearly identified in those panels. Few extreme waves recorded at the low wind flow rate in Figure 9c appear exclusively at short fetches characterized by vanishing wave heights. Extremely steep (the so-called freak) waves are usually defined as waves with heights exceeding 2H s or 8 η 2 1/2 . Thus, it follows from the present results that the probability of such waves is considerably overestimated by the Rayleigh distribution; those waves are totally absent in all ensembles accumulated at the high wind velocity, as shown in Figure 9d.
The asymmetry of wind waves was further investigated by Zavadski and Shemer [88]. The positive values of the skewness coefficient λ 3 in Figure 9a indicate asymmetry relative to the horizontal axis. The asymmetry of individual waves relative to the vertical axis can be characterized in terms of the distances between a successive wave crest, trough, and the next crest [89]. However, this approach is applicable mainly for the characterization of individual two-dimensional waves; it is ill-suited to describe the statistically representative asymmetry of random wind waves. Elgar [90] and Elgar and Guza [91] demonstrated that the asymmetry A relative to the vertical axis is related to the imaginary part of the frequency bispectrum of the surface elevation variation in time, η(t). They showed that the imaginary part of the bispectrum can be presented by the Hilbert transform H(η(t)), and that the asymmetry A can be seen as the normalized skewness of H(η(t)): where the overbar denotes time averaging.
The asymmetry of the surface elevation calculated using Equation (11) is plotted in Figure 10a for several fetches and wind velocities. The positive values of A(η) correspond to the forward face of the wave being steeper, on average, than the rear slope. For a given wind velocity, the asymmetry decreases with fetch. The asymmetry of wind waves was further investigated by Zavadski and Shemer [88]. The positive values of the skewness coefficient λ3 in Figure 9a indicate asymmetry relative to the horizontal axis. The asymmetry of individual waves relative to the vertical axis can be characterized in terms of the distances between a successive wave crest, trough, and the next crest [89]. However, this approach is applicable mainly for the characterization of individual two-dimensional waves; it is ill-suited to describe the statistically representative asymmetry of random wind waves. Elgar [90] and Elgar and Guza [91] demonstrated that the asymmetry A relative to the vertical axis is related to the imaginary part of the frequency bispectrum of the surface elevation variation in time, η(t). They showed that the imaginary part of the bispectrum can be presented by the Hilbert transform H(η(t)), and that the asymmetry A can be seen as the normalized skewness of H(η(t)): where the overbar denotes time averaging.
(a) (b) (c) The asymmetry of the surface elevation calculated using Equation (11) is plotted in Figure 10a for several fetches and wind velocities. The positive values of A(η) correspond to the forward face of the wave being steeper, on average, than the rear slope. For a given wind velocity, the asymmetry decreases with fetch.
In addition to wave heights, the exceedance of probability distributions for crest heights and trough depths were also computed in [88]. The crest and trough parameters were calculated as local maxima and local minima between two consecutive positive zero crossings. In those experiments, the total sampling duration was about 3 hours at each fetch and wind velocity, allowing to record ensembles containing the number of individual waves (0.3-0.8· 10 5 ) that are even larger than those used for computing statistics in Figure 9. Therefore, probability distribution tails could be accurately estimated.
For linear narrow-band Gaussian waves, the version of Equation (11) for the wave crests ηc and troughs ηt for the Rayleigh exceedance distributions is: The measured exceedance distributions for crests, troughs, and wave heights are plotted in Figure 10b,c for two representative cases. All the measured exceedance distributions of ( ) and ( , ) deviate notably from the Rayleigh shape. As in Figure 9c,d, the probability of extremely steep (rogue) waves, with heights exceeding 2 ̅̅̅ 1/2 by a factor of eight, is below the Rayleigh prediction by orders of magnitude in all cases. The probability of extremely high crests exceeding 4 2 ̅̅̅ 1/2 is vanishingly small compared to that corresponding to the Rayleigh distribution. However, the probability of very high crests is by orders of magnitude larger than that of very deep troughs. In addition to wave heights, the exceedance of probability distributions for crest heights and trough depths were also computed in [88]. The crest and trough parameters were calculated as local maxima and local minima between two consecutive positive zero crossings. In those experiments, the total sampling duration was about 3 h at each fetch and wind velocity, allowing to record ensembles containing the number of individual waves (0.3-0.8 × 10 5 ) that are even larger than those used for computing statistics in Figure 9. Therefore, probability distribution tails could be accurately estimated.
For linear narrow-band Gaussian waves, the version of Equation (11) for the wave crests η c and troughs η t for the Rayleigh exceedance distributions is: The measured exceedance distributions for crests, troughs, and wave heights are plotted in Figure 10b,c for two representative cases. All the measured exceedance distributions of F(η) and F(η c,t ) deviate notably from the Rayleigh shape. As in Figure 9c,d, the probability of extremely steep (rogue) waves, with heights exceeding η 2 1/2 by a factor of eight, is below the Rayleigh prediction by orders of magnitude in all cases. The probability of extremely high crests exceeding 4η 2 1/2 is vanishingly small compared to that corresponding to the Rayleigh distribution. However, the probability of very high crests is by orders of magnitude larger than that of very deep troughs. Therefore, the vertical asymmetry due to nonlinear bound waves is well pronounced, resulting in significant skewness, λ 3 . Additional data and a comparison of results with some theoretical distributions that account for nonlinear effects are presented in [88]. It is shown there that the nonlinear exceedance distributions still overestimate the probability of very steep waves.

Coupling between Waves and Airflow
Surface current is created due to the shear stress induced by the wind blowing over the water surface, as well as due to waves' nonlinearity; the water surface velocity provides the boundary condition for the forcing wind. To assess the effect of the drift current on the wind waves' parameters, measurements of the surface water velocity were performed in [51] using particle tracking velocimetry (PTV). Small light-weighted paper disks dispensed a few millimeters above the water surface through a long dry pipe served as tracers and were imaged by a video camera. Water velocity was measured at several fetches up to x = 410 cm from the entrance to the test section; at each fetch and wind velocity, about 350-400 were used. Special care was taken to provide adequate illumination; for details, see [66]. The 640 × 480-pixel camera imaged a 15-cm long area at the rate of 60 fps. The recorded video clips allowed monitoring each particle's trajectory and calculating each particle's velocity components along and across the imaged section. Two surface velocity components, U s and V s , were tracked across all the frames of appearance using the PTV processing software. While the PTV-derived mean values of the crosswind surface velocity V s were found to be close to zero, as expected ( Figure 2), the along-wind component U s can be quite significant, and was found to constitute about 3% of the maximum wind velocity in the test section, U max . Note that the values of the water surface velocity have been subtracted from the measured air velocities U(z) while performing the fit according to Equation (4). The logarithmic velocity profiles plotted in Figure 2 account for this adjustment; however, the importance of this correction proved to be mostly negligible.
The presence of water surface drift may be more substantial for the sufficiently short slow-propagating gravity-capillary waves. In an absence of current, the radian frequency ω = 2πf is related to the wave number k by the gravity-capillary dispersion relation: where σ is the water-air surface tension. The shear current close to the water surface can affect the dispersion relation due to the Doppler shift. To evaluate this effect, the wave phase velocity (celerity) of short waves, c = ω/k, was measured directly using two consecutive probes, thus allowing determining the wave number from independently measured wave celerities and frequencies. Each frequency component of the spectrum was examined separately. The celerity for each frequency for which the cross-correlation between the two gauge outputs was sufficiently high was determined by calculating the complex cross-spectrum of the two wave gauges' output. The ensemble of measured values of c determined at the local peak frequency f p at numerous fetches and wind flow rates was compared to the celerities of the corresponding gravity-capillary waves calculated based on Equation (14), c gc . The obtained results were fitted to an empirical dispersion relation of the form: Subsequently, an empirical dispersion relation in Equation (15) was suggested in [88] with coefficients a = 0.006 m and b = −2.2·10 −5 m, which are somewhat different from those in [51]. The relation in [90] was obtained using a different approach; it fits the data of [51] well, but covers a wider range of wave numbers: 20 m −1 < k < 150 m −1 . This relation is in agreement with results reported elsewhere regarding the presence of the wind-induced shear flow [28,92]. More details on experiments that resulted in the empirical dispersion relation in Equation (15) are given in [51,88].
Contrary to studies in which wave-related properties in the airflow were investigated over deterministic mechanically generated waves, the turbulent flow over a naturally evolving random wind-wave field is considered here. While the records of the surface elevation show quasi-periodic behavior with the local dominant typical frequency, air velocity components u(t) and w(t) exhibit high-frequency fluctuations that are typical for turbulent flow. The frequency spectra of u(t) and w(t) exhibit peaks at the same frequency f p as in the spectra of η(t); however, the velocity fluctuations spectra are notably wider than those of η; their amplitudes remain significant for the whole range of frequencies [53]. While the turbulent fluctuations in u and w are modulated by wind waves at frequencies close to f p , even near the air-water interface and in the presence of steep and relatively narrow-band wind waves, the fluctuations in both velocity components and the surface elevation variation are not closely related. The relation between water waves and air velocity fluctuations weakens further at higher elevations, shorter fetches, or lower wind flow rates.
To extract the quantitative relations between various random parameters in the absence of fixed frequency and phase references, cross-spectral analysis was applied. The extent of similarity between two time-dependent signals f (t) and g(t) at various radian frequencies ω can be characterized by the cross-correlation coefficient r fg (ω) calculated as the absolute value of the normalized cross-spectrum Γ fg (ω) [93]. The phase shift θ fg (ω) between the two signals f (t) and g(t) at frequency ω is defined as the argument of the corresponding component of the complex cross-spectrum value at this frequency, S fg (ω), so that θ fg is positive when f leads g, and negative otherwise. Apparently, it only makes sense to consider the phase shift between the signals at frequencies where they are sufficiently well correlated, as represented by the value of r fg (ω).
The records accumulated in [53] were divided into windows containing 1024 data points for each signal, corresponding to the window length of about 8 s (frequency resolution of about 0.12 Hz), with 50% overlap. About 70 independent estimates of the cross-spectra were obtained from the available records for each sensor location and airflow rate; the averaged over all windows cross-spectra were computed at all flow conditions and for different measuring locations.
The cross-correlation coefficients between the fluctuating parts of the two velocity components, u and w, and the surface elevation at wind velocities below U max = 7.7 m/s and fetches shorter than 220 cm were too low to enable extracting meaningful information. The cross-correlation coefficients between the horizontal velocity fluctuations u at a number of vertical locations and the surface elevation, r ηu are presented at fetch x = 260 cm and wind velocity U max = 11.2 m/s in Figure 11a. Figure 11b depicts the corresponding cross-correlation coefficient for the vertical velocity component, r ηw , at x = 340 cm. In both panels, a well-defined peak around the fundamental wave frequency f p is observed; the values of f p decrease with fetch. The maximum of the cross-correlation coefficient with the surface elevation is notably higher for the vertical velocity fluctuations as compared to the horizontal ones. This difference in the correlation level stems from the fact that the wave-induced fluctuations of velocity fluctuations in airflow above water have similar amplitudes in the longitudinal and vertical directions, whereas the velocity fluctuations in the horizontal direction due to the background turbulence not related to the underlying wave field are significantly higher than the vertical ones. The values of the correlation coefficients decay with height. Secondary peaks in the frequency dependence of the cross-correlation coefficient can be noticed at the second harmonic of the peak wave frequency, and to some extent even at the third harmonic. These secondary peaks indicate that velocity fluctuations in air are also affected by higher-order bound wave harmonics.
The resulting phase shifts θ ηu and θ ηw for U max = 11.2 m/s and two fetches, x = 260 cm and x = 340 cm, are presented in Figure 11c. For each velocity component, the phase shifts relative to the surface elevation are clustered in tight "clouds". The larger scatter in the θ ηu "cloud" as compared to that of θ ηw can be attributed to the lower correlation level of the longitudinal velocity component with the surface elevation variation, as compared to r ηw , as seen in Figure 11a,b. The measured phase shifts seem to be independent of elevation. For both fetches in Figure 11c, the average phase shift between η and the longitudinal velocity fluctuations is close to −30 • , and between η and w, it is close to 60 • . Figure 11d demonstrates that similar results are obtained for other wind velocities and fetches.
The results of Figure 11c,d can be seen as phase shifts between the wave-coherent velocity fluctuations in air and the water surface elevation η. The phase shift that is close to 90 • between the two velocity components results in a vanishing phase-coherent shear Reynolds stress. To validate this conclusion, an attempt was made in [53] to examine directly the correlation between the surface elevation η(t) and the fluctuating instantaneous contributions to the Reynolds shear stress, −u w (t). Even for the longest fetch and highest wind velocity employed in that study, and thus the longest and highest wind waves, no significant correlation between those two time series has been found; the correlation coefficient was small at all frequencies, including the vicinity of f p . Thus, it seems that not only do the the wave-induced velocity fluctuations in air not contribute to the mean Reynolds shear stress, there is also no detectable wave-induced time-dependent variation of −u w with the phase of the dominant wave. Even for the longest fetch and highest wind velocity employed in that study, and thus the longest and highest wind waves, no significant correlation between those two time series has been found; the correlation coefficient was small at all frequencies, including the vicinity of fp. Thus, it seems that not only do the the wave-induced velocity fluctuations in air not contribute to the mean Reynolds shear stress, there is also no detectable wave-induced time-dependent variation of uw   with the phase of the dominant wave. Measurements of the fluctuations of the static pressure in the airflow above the water surface were made in [51] using the static pressure probe described in [52]. Some representative results showing the dependence on frequency of the cross-correlation coefficients between the static pressure fluctuations in air and the water surface elevation, rpη, are given in [51] for several elevations above the mean water surface z. The values of rpη are generally quite small for most frequencies, the correlation between surface elevation and pressure variations becoming significant around the peak Figure 11. Characterization of wave-coherent air velocity fluctuations. (a) Cross-correlation coefficients r ηu at various elevations z above the highest crest at x = 260 cm and U max = 11.2 m/s; (b) as in (a) for r ηw, x = 340 cm and U max = 11.2 m/s; (c) Phase shifts θ ηu and θ ηw in the vicinity of the peak frequency at x = 260 cm and U max = 11.2 m/s; (d) Phase shifts between the two velocity components and the surface elevation at the peak wind-wave frequency for various fetches and wind velocities.
Measurements of the fluctuations of the static pressure in the airflow above the water surface were made in [51] using the static pressure probe described in [52]. Some representative results showing the dependence on frequency of the cross-correlation coefficients between the static pressure fluctuations in air and the water surface elevation, r pη , are given in [51] for several elevations above the mean water surface z. The values of r pη are generally quite small for most frequencies, the correlation between surface elevation and pressure variations becoming significant around the peak frequency f p and only sufficiently close to the air-water interface. Note that for lower wind velocities in the test section and thus smaller wind-wave amplitudes, the wave-coherent pressure fluctuations were correspondingly weaker, rendering the less reliable results on cross-coherence between air pressure fluctuations and the surface elevation variations. The dependence on frequency of the cross-correlation coefficient between air pressure and η strongly resembles the shape of wave frequency spectra, cf. Figure 6. Close to the air-water interface, significant correlation seems to exist also at frequencies corresponding to the second, and to some extent, even to the third harmonic of f p , thus reflecting the contribution of bound waves to wave-coherent air pressure fluctuations.
The phase difference values between pressure fluctuations and surface elevation were obtained from the complex pressure-surface elevation cross-spectra for those frequencies only at which significant correlation was detected (defined as r pη > 0.6). The dependence of θ pη on frequency (around the peak frequency at about 3 Hz) is presented in Figure 12b,c for two fetches and the highest wind flow rate in experiments. For all experimental conditions at which significant correlation between pressure and surface elevation fluctuations was detected, the values of θ pη were smaller than 180 • and showed no notable variations with z or with fetch. Thee phase angles of pressure fluctuations presented in Figure 12a,b are in good agreement with the results reported in [9,31,94]. The phase difference values between pressure fluctuations and surface elevation were obtained from the complex pressure-surface elevation cross-spectra for those frequencies only at which significant correlation was detected (defined as rpη > 0.6). The dependence of θpη on frequency (around the peak frequency at about 3 Hz) is presented in Figure 12b,c for two fetches and the highest wind flow rate in experiments. For all experimental conditions at which significant correlation between pressure and surface elevation fluctuations was detected, the values of θpη were smaller than 180° and showed no notable variations with z or with fetch. Thee phase angles of pressure fluctuations presented in Figure 12a,b are in good agreement with the results reported in [9,31,94].
The RMS values of the pressure fluctuations in each record were extrapolated to the undisturbed water surface level, and used in [51] as characterizing the magnitude of wave-coherent pressure fluctuations at the local dominant wave frequency for each airflow rate in the test section, while the values of θpη represented their phases. Those results allowed estimating the form drag τform [4, 5,13], the amount of work being done by pressure on the water surface, and the growth rate of the wave energy at the dominant frequency. Detailed results are given in [51].

The Two-Dimensional Features of the Wind-Wave Field
Measurements of wind-generated waves, such as described above, are mostly based on wave gauges that provide information on the temporal variation of the surface elevation at the sensor The RMS values of the pressure fluctuations in each record were extrapolated to the undisturbed water surface level, and used in [51] as characterizing the magnitude of wave-coherent pressure fluctuations at the local dominant wave frequency for each airflow rate in the test section, while the values of θ pη represented their phases. Those results allowed estimating the form drag τ form [4, 5,13], the amount of work being done by pressure on the water surface, and the growth rate of the wave energy at the dominant frequency. Detailed results are given in [51].

The Two-Dimensional Features of the Wind-Wave Field
Measurements of wind-generated waves, such as described above, are mostly based on wave gauges that provide information on the temporal variation of the surface elevation at the sensor location. However, the wind-wave field is never unidirectional. The diverse approaches that have been applied in attempts to extract information on the two-dimensional structure of the wave field were briefly reviewed in [66]. To this end, two essentially different methods described in Section 2 were used in our studies to gain information on the 2D spatial structure of wind waves: stereo video imaging and the laser slope gauge (LSG). Due to improvements in computational power and computational algorithms, the application of stereo video imaging in field experiments with cameras either installed on offshore platforms or on a moving vessel expanded significantly in recent years. The comparison of video imaging results with those obtained using a standard instrumentation that was carried out in those experiments provided evidence that at open sea scales, stereo wave imaging is effective to retrieve data on medium-to-short wavelengths; see [17] and additional references therein. Stereo video imaging, similarly to alternative optical methods of reconstruction of the instantaneous water surface shape, requires considerable computational resources that currently limit their effective spatial resolution and the overall duration of continuous imaging. The question of accuracy of the surface shape obtained using stereo imaging not yet been fully resolved.
In [66], stereo video imaging was applied in a small-scale laboratory facility to study the two-dimensional structure of the wind-wave field under controlled conditions. Fast video recording enables extracting both the temporal and spatial information. The application of various instruments under identical forcing conditions enabled the verification of reliability and determination of the limitations of the applied methods. Details of the verification procedure, various sources of errors, and estimates of accuracy of the stereo video imaging reconstruction process are given in [66]. For stereo imaging of the wind-wave field, all other sensors were removed, and the video images were acquired in a separate set of experiments at three fetches, x = 110 cm, x = 220 cm, and x = 340 cm, and three wind speeds, U max = 6.5 m/s, 8.5 m/s, and 10.5 m/s. A typical reconstructed "snapshot" representing the instantaneous surface elevation in the surface area that is common to both cameras is presented in Figure 13. In [66], stereo video imaging was applied in a small-scale laboratory facility to study the twodimensional structure of the wind-wave field under controlled conditions. Fast video recording enables extracting both the temporal and spatial information. The application of various instruments under identical forcing conditions enabled the verification of reliability and determination of the limitations of the applied methods. Details of the verification procedure, various sources of errors, and estimates of accuracy of the stereo video imaging reconstruction process are given in [66]. For stereo imaging of the wind-wave field, all other sensors were removed, and the video images were acquired in a separate set of experiments at three fetches, x = 110 cm, x = 220 cm, and x = 340 cm, and three wind speeds, Umax = 6.5 m/s, 8.5 m/s, and 10.5 m/s. A typical reconstructed "snapshot" representing the instantaneous surface elevation in the surface area that is common to both cameras is presented in Figure 13. Surface elevation variation with time at a fixed location, η(x,y,t), can be retrieved from the sequence of reconstructed 3D surfaces, as shown in Figure 13. The root mean square (RMS) values of the surface elevation / were computed at all wind velocities and fetches from a time series taken from the virtual wave gauges located across the central part of the reconstructed images. These records were compared in [66] with the wave-gauge results at an identical wind velocity and fetch (but not acquired simultaneously with the video imaging). Comparison of the two signals indicated that they have comparable wave heights and dominant frequencies; the records derived from the stereo imaging were usually somewhat noisier. Taking the wave gauge results as a reference, the stereo imaging method overestimates / by up to 30% at low wind velocities and thus small wave amplitudes; however, for stronger winds and longer fetches resulting in the surface elevation fluctuations exceeding 5 mm, the relative error was reduced to under 15%. In view of the inaccuracy of stereo imaging of low-amplitude waves, frequency spectra of the surface elevation variation obtained by the stereo imaging and by the wave gauge were compared in [66] only for the moderate and high wind forcing and at longer fetches. Although not identical, the spectra obtained by different techniques exhibit similarity in all examined cases, in particular within about 30% of the peak frequency. Away from the spectral peak, and thus for much smaller wave amplitudes, the frequency spectra derived from sequences of reconstructed stereo video records become less reliable. Thus, it transpires that at the currently available resolution, stereo imaging is capable of providing information on the gross features of the wave field, while more fine details may be inaccurate. Thus, the major advantage offered by the stereo imaging over alternative measuring methods is were computed at all wind velocities and fetches from a time series taken from the virtual wave gauges located across the central part of the reconstructed images. These records were compared in [66] with the wave-gauge results at an identical wind velocity and fetch (but not acquired simultaneously with the video imaging). Comparison of the two signals indicated that they have comparable wave heights and dominant frequencies; the records derived from the stereo imaging were usually somewhat noisier. Taking the wave gauge results as a reference, the stereo imaging method overestimates η 2 1/2 by up to 30% at low wind velocities and thus small wave amplitudes; however, for stronger winds and longer fetches resulting in the surface elevation fluctuations exceeding 5 mm, the relative error was reduced to under 15%. In view of the inaccuracy of stereo imaging of low-amplitude waves, frequency spectra of the surface elevation variation obtained by the stereo imaging and by the wave gauge were compared in [66] only for the moderate and high wind forcing and at longer fetches. Although not identical, the spectra obtained by different techniques exhibit similarity in all examined cases, in particular within about 30% of the peak frequency. Away from the spectral peak, and thus for much smaller wave amplitudes, the frequency spectra derived from sequences of reconstructed stereo video records become less reliable. Thus, it transpires that at the currently available resolution, stereo imaging is capable of providing information on the gross features of the wave field, while more fine details may be inaccurate.
Thus, the major advantage offered by the stereo imaging over alternative measuring methods is its ability to characterize the global structure of the wind-wave field. This advantage was exploited in [66] to assess the spatial coherence of the wave field applying the spatial autocovariance function in the along-wind and crosswind directions, which is defined as: respectively, as a quantitative parameter. The autocovariance function f (∆x) was computed along the imaged area for constant values of the spanwise coordinate y; the reference value of x was taken at the right (upwind) end of the image. The resulting values of f (∆x) were averaged over various y coordinates in the central part of the image, and then over all the reconstructed surfaces. Similarly, in computations of the autocovariance function in the crosswind direction, g(∆y), the averaging was first performed over different constant values of x within the reconstructed surface elevation map; subsequent averaging was carried out over the whole set of the instantaneous reconstructed surfaces.
The computed autocovariance functions f (∆x) and g(∆y) for different fetches x and wind velocities U max are plotted in Figure 14. In all examined cases, the autocovariance functions decay fast and effectively attain zero after the second local minimum. In analogy with the autocovariance function for a deterministic monochromatic wave, the values of ∆x at the first local minimum in all the panels of Figure 14 may be seen as corresponding to one-half of the local dominant wavelength λ d . Since the wave field is spatially inhomogeneous and the characteristic wavelengths λ d increase with the fetch for each wind velocity, the autocovariance functions f (∆x) were computed by direct integration without invoking the Fourier transform that implies periodicity. The estimated wavelengths range from 7 cm to 27 cm, which according to the dispersion relation in Equation (15) correspond to dominant frequencies ranging from 2.4 Hz to 4.7 Hz.
The locations of minima and maxima of the autocovariance functions in the crosswind direction, g(∆y), plotted in Figure 14d-f do not differ notably from those in panels a to c of that Figure, although the functions g(∆y) decay faster (as emphasized by the different horizontal scales adopted for f (∆x) and g(∆y)). The functions f (∆x) and g(∆y) effectively vanish at distances comparable to the dominant wavelength, indicating that the spatial variations of the surface elevation become statistically independent at relatively short distances. Those short correlations lengths in both along-wind and crosswind directions are in accordance with the snapshot in Figure 13 that shows short-crested and three-dimensional wind waves in the present experiments. the imaged area for constant values of the spanwise coordinate y; the reference value of x was taken at the right (upwind) end of the image. The resulting values of (∆ ) were averaged over various y coordinates in the central part of the image, and then over all the reconstructed surfaces. Similarly, in computations of the autocovariance function in the crosswind direction, (∆ ), the averaging was first performed over different constant values of x within the reconstructed surface elevation map; subsequent averaging was carried out over the whole set of the instantaneous reconstructed surfaces.  Figure 14. In all examined cases, the autocovariance functions decay fast and effectively attain zero after the second local minimum. In analogy with the autocovariance function for a deterministic monochromatic wave, the values of ∆x at the first local minimum in all the panels of Figure 14 may be seen as corresponding to one-half of the local dominant wavelength λd. Since the wave field is spatially inhomogeneous and the characteristic wavelengths λd increase with the fetch for each wind velocity, the autocovariance functions f(∆x) were computed by direct integration without invoking the Fourier transform that implies periodicity. The estimated wavelengths range from 7 cm to 27 cm, which according to the dispersion relation in Equation (15) correspond to dominant frequencies ranging from 2.4 Hz to 4.7 Hz.
The locations of minima and maxima of the autocovariance functions in the crosswind direction, g(∆y), plotted in Figure 14d-f do not differ notably from those in panels a to c of that Figure, although the functions g(∆y) decay faster (as emphasized by the different horizontal scales adopted for f(∆x) and g(∆y)). The functions f(∆x) and g(Δy) effectively vanish at distances comparable to the dominant wavelength, indicating that the spatial variations of the surface elevation become statistically independent at relatively short distances. Those short correlations lengths in both along-wind and crosswind directions are in accordance with the snapshot in Figure 13 that shows short-crested and three-dimensional wind waves in the present experiments. Integral length scales that characterize the spatial extent of coherence of the wind-wave field can in downwind, Λ x , and crosswind, Λ y , directions can be defined as: respectively. Both integral length scales Λ x and Λ y normalized by the local characteristic scale of wind waves, λ d , are plotted in Figure 15. The dominant wind-wave lengths λ d were determined for each fetch and wind velocity from the dominant frequency f d using the empirical dispersion relation shown in Equation (15). The results of Figure 15 demonstrate that the normalized integral scales Λ x and Λ y are notably smaller than unity. This implies that although the wind waves propagate in a relatively narrow channel with mean wind velocity directed along the test section, the waves cease to be coherent at distances much shorter than the dominant wavelength λ d . Since for all wind velocities and fetches, the dominant wavelengths remain smaller than the channel width, the waves in the central part of the test section are thus essentially unaffected by the presence of sidewalls. In view of results presented in Figure 15b, the data taken along the lines y = const. that are separated by more than about 5 cm can be considered as statistically independent.
velocities and fetches, the dominant wavelengths remain smaller than the channel width, the waves in the central part of the test section are thus essentially unaffected by the presence of sidewalls. In view of results presented in Figure 15b, the data taken along the lines y = const. that are separated by more than about 5 cm can be considered as statistically independent.
(a) (b) To complete the discussion on wind waves' coherence, the temporal autocovariance function was defined in analogy with Equation (16) as: and computed from wave gauge measurements, laser slope gauge, and stereo video imaging [66]. The wave gauge-derived absolute values of h(∆t) at the first minimum are significantly larger than the corresponding values of f(∆x). The fact that the wavenumber spectrum is wider than the frequency spectrum (approximately by a factor of two for pure gravity waves) contributes to this difference. Moreover, a fetch-limited wave field is statistically stationary, but spatially inhomogeneous. The spatial variation of the wind-wave field over the expanse of the image causes an additional decrease in the autocovariance function f(∆x). The temporal autocovariance functions retrieved from the virtual wave gauge at a single point of the reconstructed surfaces from stereo imaging and the corresponding autocovariance function calculated from the measurements carried out by a real wave gauge yielded close locations of the first local minima. The temporal correlation estimated by all measurement methods vanishes for time delays ∆t of the order of the corresponding wave period. An attempt was made in [66] to compute the directional wave spectra S(k,θ), in which k = (kx 2 +ky 2 ) 1/2 is the wavenumber and θ = tan −1 (ky/kx) defines the wave vector angle using sequences of the threedimensional instantaneous water surfaces reconstructed from stereo video imaging. Then, the To complete the discussion on wind waves' coherence, the temporal autocovariance function was defined in analogy with Equation (16) as: and computed from wave gauge measurements, laser slope gauge, and stereo video imaging [66]. The wave gauge-derived absolute values of h(∆t) at the first minimum are significantly larger than the corresponding values of f (∆x). The fact that the wavenumber spectrum is wider than the frequency spectrum (approximately by a factor of two for pure gravity waves) contributes to this difference. Moreover, a fetch-limited wave field is statistically stationary, but spatially inhomogeneous. The spatial variation of the wind-wave field over the expanse of the image causes an additional decrease in the autocovariance function f (∆x). The temporal autocovariance functions retrieved from the virtual wave gauge at a single point of the reconstructed surfaces from stereo imaging and the corresponding autocovariance function calculated from the measurements carried out by a real wave gauge yielded close locations of the first local minima. The temporal correlation estimated by all measurement methods vanishes for time delays ∆t of the order of the corresponding wave period. An attempt was made in [66] to compute the directional wave spectra S(k,θ), in which k = (k x 2 +k y 2 ) 1/2 is the wavenumber and θ = tan −1 (k y /k x ) defines the wave vector angle using sequences of the three-dimensional instantaneous water surfaces reconstructed from stereo video imaging. Then, the resulting power spectra were averaged over the whole set of reconstructed surfaces, as shown in Figure 13. The spectra at all fetches and for all wind velocities are presented in Figure 16. The directional ambiguity that is common to all instantaneous spatial images can be resolved in most cases, since the direction of wind is from right to left.
Atmosphere 2019, 10, x FOR PEER REVIEW 25 of 52 resulting power spectra were averaged over the whole set of reconstructed surfaces, as shown in Figure 13. The spectra at all fetches and for all wind velocities are presented in Figure 16. The directional ambiguity that is common to all instantaneous spatial images can be resolved in most cases, since the direction of wind is from right to left. imaging of the wave slopes in two perpendicular directions, also report on asymmetry in the directional spectra at early stages of wind waves' evolution. Deviations of the wave propagation direction from that of the wind were observed in field experiments [96]. As waves grow and become longer at high wind speeds and larger fetches, the asymmetry of the spectra disappears; the contours of the directional spectra are at 180° (the wind direction). The wave numbers associated with the dominant waves in those cases are in general agreement with the results obtained from the At low wind velocities and shorter fetches, the wave amplitudes are small; the corresponding spectra are characterized by the presence of waves in the crosswind direction. Thus, the question of directional ambiguity remains open for those cases; most probably, waves propagating in both positive and negative y directions are present. The direction of the dominant wave for moderate winds deviates slightly to the left of the wind direction; this may suggest that the reconstruction procedure for low elevations introduces directional inaccuracy. However, Caulliez and Collard [95], who performed reconstruction of the decimeter-scale water-wave surface based on the simultaneous imaging of the wave slopes in two perpendicular directions, also report on asymmetry in the directional spectra at early stages of wind waves' evolution. Deviations of the wave propagation direction from that of the wind were observed in field experiments [96]. As waves grow and become longer at high wind speeds and larger fetches, the asymmetry of the spectra disappears; the contours of the directional spectra are at 180 • (the wind direction). The wave numbers associated with the dominant waves in those cases are in general agreement with the results obtained from the autocorrelation function, the relative error being within about 15%. At highest wind velocity and the longest fetch, the wavelength is overestimated in the spectral representation; this is apparently due to insufficient image size, which in this case becomes comparable with the wavelength. The spectra are characterized by directional spreading of up to 60 • around the wind direction and thus represent an essentially three-dimensional and short-crested wave field.

Umax
A different experimental approach that sheds light on the three-dimensional structure of wind waves is based on application of the laser slope gauge (LSG). Temporal records of the variation of the instantaneous surface slope components provided by LSG allow characterizing the direction of the vector normal to the instantaneous surface shape at the measuring location. The projection of the vector normal to the surface on the horizontal plane allows determining the azimuthal angle of the instantaneous slope, θ = tan −1 (η y / η x ). The probability distribution functions of the instantaneous azimuthal angles are presented in Figure 17a for three fetches and a single wind velocity U max = 8.5 m/s, and in Figure 17b for the fetch x = 220 cm and different wind velocities. All the probability distribution functions plotted in Figure 17 are nearly symmetric around the channel axis. The directional spreading is essential; the probability of the azimuthal angle in the crosswind direction (θ = ± 90 • ) is significant. All the curves in Figure 17 have two distinctive peaks at θ = 0 • and θ = 180 • ; the prevailing direction of the waves is aligned with that of the wind. The probability that the azimuthal angle θ = 0 • that is obtained for a forward-leaning slope is notably lower than that of the rear inclination. This front-back asymmetry of the wave shape based on the LSG measurements observed in Figure 17 is in agreement with the asymmetry A(η) of the surface elevation measured by a wave gauge, as plotted in Figure 10a. The probability distributions of the slope inclination direction in Figure 17a do not vary notably with fetch; similar behavior was obtained for all other wind velocities. The effect of the wind velocity on the variability of the slope direction is much more noticeable, as shown in Figure 17b. The probability distribution functions of the azimuthal angle θ become flatter, and the probability of the instantaneous slope in the crosswind direction increases for stronger winds.
to insufficient image size, which in this case becomes comparable with the wavelength. The spectra are characterized by directional spreading of up to 60° around the wind direction and thus represent an essentially three-dimensional and short-crested wave field.
A different experimental approach that sheds light on the three-dimensional structure of wind waves is based on application of the laser slope gauge (LSG). Temporal records of the variation of the instantaneous surface slope components provided by LSG allow characterizing the direction of the vector normal to the instantaneous surface shape at the measuring location. The projection of the vector normal to the surface on the horizontal plane allows determining the azimuthal angle of the instantaneous slope, θ = tan −1 (ηy/ ηx). The probability distribution functions of the instantaneous azimuthal angles are presented in Figure 17a for three fetches and a single wind velocity Umax = 8.5 m/s, and in Figure 17b for the fetch x = 220 cm and different wind velocities. All the probability distribution functions plotted in Figure 17 are nearly symmetric around the channel axis. The directional spreading is essential; the probability of the azimuthal angle in the crosswind direction (θ = ± 90°) is significant. All the curves in Figure 17 have two distinctive peaks at θ = 0° and θ = 180°; the prevailing direction of the waves is aligned with that of the wind. The probability that the azimuthal angle θ = 0° that is obtained for a forward-leaning slope is notably lower than that of the rear inclination. This front-back asymmetry of the wave shape based on the LSG measurements observed in Figure 17 is in agreement with the asymmetry A(η) of the surface elevation measured by a wave gauge, as plotted in Figure 10a. The probability distributions of the slope inclination direction in Figure 17a do not vary notably with fetch; similar behavior was obtained for all other wind velocities.
The effect of the wind velocity on the variability of the slope direction is much more noticeable, as shown in Figure 17b. The probability distribution functions of the azimuthal angle θ become flatter, and the probability of the instantaneous slope in the crosswind direction increases for stronger winds. The measured by LSG temporal dependencies of downwind ∂η/∂x(t) and crosswind ∂η/∂y(t) slope components allow comparing their frequency spectra with those of the surface elevation η(t). Such a comparison is carried out in Figure 18 for the fetch x = 220 cm and different wind velocities Umax. The power spectra of η(t) are also plotted in this figure for comparison. As in the spectra of η(t) The measured by LSG temporal dependencies of downwind ∂η/ ∂x(t) and crosswind ∂η/ ∂y(t) slope components allow comparing their frequency spectra with those of the surface elevation η(t). Such a comparison is carried out in Figure 18 for the fetch x = 220 cm and different wind velocities U max . The power spectra of η(t) are also plotted in this figure for comparison. As in the spectra of η(t) (Figure 6), the total power of the slope components spectra increases with U max . Apart from the principal peak at f p , the slope spectra are characterized by secondary peaks representing bound waves; those secondary peaks seem to be less prominent in the spectra of ∂η/ ∂y. The peak frequencies decrease with the wind velocity.
Atmosphere 2019, 10, x FOR PEER REVIEW 27 of 52 ( Figure 6), the total power of the slope components spectra increases with Umax. Apart from the principal peak at fp, the slope spectra are characterized by secondary peaks representing bound waves; those secondary peaks seem to be less prominent in the spectra of ∂η/∂y. The peak frequencies decrease with the wind velocity. For all wind velocities, the peak frequencies fp of the slope components do not differ much from those of the surface elevation spectra. Nevertheless, it can be noticed that the values of fp of the downwind surface slope spectra in Figure 18b exceed somewhat those of the surface elevation spectra in Figure 18a. This can be attributed to the fact that the power of the downwind surface slope, ∂η/∂x(t), is proportional to (akx) 2 , in which kx is the longitudinal component of the wave vector, rather than of the wave amplitude squared, a 2 . The peak frequencies of the crosswind surface slope ∂η/∂y(t) do not differ significantly from the peak frequencies of ∂η/∂x(t). Contrary to the surface elevation spectra, the spectral intensities of the slope components in the peak region seem to be only weakly dependent on U, in particular for longer fetches. At the same fetch and wind velocity, the power of the downwind slope component around fp exceeds significantly that of the corresponding crosswind component; however, at frequencies exceeding 10 Hz, they become comparable. Thus, it appears that short wind waves in the capillary range may not have a preferred propagation direction.
The high frequency part of the power spectra of the surface slope components ∂η/∂x and ∂η/∂y can be presented by the power law shown in Equation (6), similarly to the surface elevation η; see Figure 7b. The slope exponents corresponding to the power law fit of all the tails of the spectra of ∂η/∂x and ∂η/∂y are summarized in Figure 19 for all the fetches and wind velocities employed in [88]. The exponents obtained in those experiments for the surface elevation spectra are also plotted for comparison. At lower wind velocities and shorter fetches, wave amplitudes are quite small; therefore, the results exhibit considerable scatter. For all the cases considered, the exponents of the spectral tail of η are close to n = 3, which is in agreement with Figure 7b for this frequency range. The spectral tail slope exponents for both downwind ∂η/∂x and crosswind ∂η/∂y surface slope components do not differ significantly, and with increasing wind velocity, seem to converge to a value close to unity.
As noticed above, at each frequency f, the ratio of the spectral component powers of the surface slope ∂η/∂x and of the surface elevation η is k 2 (f). It follows from the linear dispersion relation for gravity-capillary waves shown in Equation (14) that at frequencies exceeding about 10 Hz, the contribution of the surface tension term in Equation (14) becomes dominant; thus, k ~ ω 2/3 , and one can expect -nη = 4/3. The difference between exponents in Figure 19 apparently exceeds this value. This discrepancy can be understood in view of the fact that the frequencies of short wind waves are affected by Doppler shift, which is induced by the shear stress and Stokes drift current, as discussed in Section 3.3 in relation to the empirical dispersion relation shown in Equation (15). For all wind velocities, the peak frequencies f p of the slope components do not differ much from those of the surface elevation spectra. Nevertheless, it can be noticed that the values of f p of the downwind surface slope spectra in Figure 18b exceed somewhat those of the surface elevation spectra in Figure 18a. This can be attributed to the fact that the power of the downwind surface slope, ∂η/ ∂x(t), is proportional to (ak x ) 2 , in which k x is the longitudinal component of the wave vector, rather than of the wave amplitude squared, a 2 . The peak frequencies of the crosswind surface slope ∂η/ ∂y(t) do not differ significantly from the peak frequencies of ∂η/ ∂x(t). Contrary to the surface elevation spectra, the spectral intensities of the slope components in the peak region seem to be only weakly dependent on U, in particular for longer fetches. At the same fetch and wind velocity, the power of the downwind slope component around f p exceeds significantly that of the corresponding crosswind component; however, at frequencies exceeding 10 Hz, they become comparable. Thus, it appears that short wind waves in the capillary range may not have a preferred propagation direction.
The high frequency part of the power spectra of the surface slope components ∂η/ ∂x and ∂η/ ∂y can be presented by the power law shown in Equation (6), similarly to the surface elevation η; see Figure 7b. The slope exponents corresponding to the power law fit of all the tails of the spectra of ∂η/ ∂x and ∂η/ ∂y are summarized in Figure 19 for all the fetches and wind velocities employed in [88]. The exponents obtained in those experiments for the surface elevation spectra are also plotted for comparison. At lower wind velocities and shorter fetches, wave amplitudes are quite small; therefore, the results exhibit considerable scatter. For all the cases considered, the exponents of the spectral tail of η are close to n = 3, which is in agreement with Figure 7b for this frequency range. The spectral tail slope exponents for both downwind ∂η/ ∂x and crosswind ∂η/ ∂y surface slope components do not differ significantly, and with increasing wind velocity, seem to converge to a value close to unity.
As noticed above, at each frequency f, the ratio of the spectral component powers of the surface slope ∂η/ ∂x and of the surface elevation η is k 2 (f ). It follows from the linear dispersion relation for gravity-capillary waves shown in Equation (14) that at frequencies exceeding about 10 Hz, the contribution of the surface tension term in Equation (14) becomes dominant; thus, k~ω 2/3 , and one can expect n η x -n η = 4/3. The difference between exponents in Figure 19 apparently exceeds this value. This discrepancy can be understood in view of the fact that the frequencies of short wind waves are affected by Doppler shift, which is induced by the shear stress and Stokes drift current, as discussed in Section 3.3 in relation to the empirical dispersion relation shown in Equation (15). Similarly to the results presented in Figure 9a The front-rear wave asymmetry characterized by the skewness of the surface slope, A(∂η/∂x), decreases as the waves grow larger with fetch, which is in agreement with the behavior of A(η) in Similarly to the results presented in Figure 9a,b, the deviation from the Gaussianity of the surface slope components can be characterized by their higher statistical moments. The skewness coefficient of both LSG-measured surface slope components, λ 3x = η 3 x /η 2 provides quantitative estimates of wave slope asymmetry about the along-wind and crosswind axes. The normalized values of the skewness of the downwind component of the surface slope, η x , characterize the vertical wave asymmetry, and the vertical asymmetry associated with the slope A(η x ) = −λ 3x is obtained from measurements that are independent of those used to compute asymmetry A(η) according to Equation (12) and as plotted in Figure 10a. The variation of −λ 3x with wind velocity and fetch is presented in Figure 20 together with the skewness coefficients, −λ 3y . Negative values of slope skewness are associated with shorter and thus steeper forward wave faces. The absolute values of λ 3y are small (|λ 3y | < 0.1), thus indicating that waves are on average vertically symmetric in the crosswind direction.
Atmosphere 2019, 10, x FOR PEER REVIEW 28 of 52 Figure 19. Comparison of exponents of the spectral tails for surface elevation and two components of the surface slope.
Similarly to the results presented in Figure 9a,b, the deviation from the Gaussianity of the surface slope components can be characterized by their higher statistical moments. The skewness coefficient of both LSG-measured surface slope components, 3 = 3 ̅̅̅ / 2 ̅̅̅  The front-rear wave asymmetry characterized by the skewness of the surface slope, A(∂η/∂x), decreases as the waves grow larger with fetch, which is in agreement with the behavior of A(η) in The front-rear wave asymmetry characterized by the skewness of the surface slope, A(∂η/ ∂x), decreases as the waves grow larger with fetch, which is in agreement with the behavior of A(η) in Figure 10a. The vertical wave asymmetry determined from the skewness of the surface slope decreases with the wind velocity U. The asymmetry magnitudes A(η) are of the same order as those in A(∂η/ ∂x), although they are somewhat different. The dependence of A(∂η/ ∂x) on wind velocity is more pronounced than that of A(η). Additional results on the statistical parameters of wind waves and the discussion of those findings can be found in [75,88].

Spatial and Temporal Variation of Waves under Effectively Impulsive Wind Forcing
Results presented in Section 3 deal with the fetch-limited spatial evolution of a wave field under wind forcing that is statistically steady. While the interfacial shear stress under steadily blowing wind was shown to be approximately homogeneous, wave parameters vary significantly along the test section; thus, the wave field, while statistically stationary, is essentially inhomogeneous. For unsteady wind forcing, the situation is more complicated, since the statistical characteristics of waves in this case depend on time as well as on space. However, only limited experimental data on waves under time-dependent wind under controlled conditions are currently available; these works were briefly reviewed in the Introduction. This present section is based on the experiments carried out by Zavadsky and Shemer [67]; it presents time-resolved statistically reliable results on waves excited from rest by a (nearly) impulsive wind forcing. These findings were further discussed in [97].

Experimental Procedure
To accumulate statistically significant results on a wave field that varies spatially as well as temporally, special experimental and data-processing procedures were developed. Data from multiple realizations of the wave field under identical forcing conditions were recorded, thus allowing the computation of ensemble-averaged values. Each experimental run was initiated when the water surface was calm. To assess the response of the system, both the output voltage of the blower that represents the airflow rate and the wind velocity in the test section measured by a Pitot tube were recorded. The output voltage of the blower varies linearly at its maximum slew rate of 1 V/s until the prescribed steady state is attained. The experiments were carried out at three fetches: x = 120 cm, x = 220 cm, and x = 340 cm. At each fetch, measurements were performed for the following wind velocities in the test section: U = 6.5 m/s, 7.5 m/s, 8.5 m/s, 9.5 m/s, and 10.5 m/s. The results of simultaneous measurements of the blower output voltage and of the wind velocity U(t) are plotted in Figure 21. The instant of the blower initiation by the computer serves as the time reference. The duration of the ramp in the blower-driving signal varies from 3 s to 5 s. At each instant, the wind velocity lags slightly behind the blower output; however, the delay does not exceed 1 s. Note that velocities below about 1.5 m/s are not measured adequately by a Pitot tube. Once the target velocity is attained, the blower maintains a constant airflow rate in the test section for 120 s; it is then shut down.
The temporal variations of the instantaneous surface elevation η, of the surface slopes components ∂η/∂x and ∂η/∂y, of the mean wind velocity U(t), and of the voltage output were recorded at sampling frequency f s = 300 Hz; data sampling started prior to activation of the blower. The time interval between successive runs was set at about 6 min; that proved to be sufficient for the decay of waves and disturbances in the test section. The total duration of a single run, including the waiting period, thus exceeded 8 min. At least 100 independent realizations for each set of operational parameters (the target wind velocity and the fetch) were recorded. Thus, the total duration of the experiments at a single fetch and wind target velocity exceeded 13 h; for all wind velocities at a single fetch, measurements lasted for more than two days. Such long continuous experimental sessions and synchronization between the blower operation and data sampling were possible, since the whole procedure was fully automated and controlled by a computer using a single LabView program, virtually without human intervention. Ensemble averaging of the accumulated data was carried out for each fetch x and target wind velocity U over all realizations as a function of time elapsed since the initiation of the blower in each realization. The characteristic amplitudes of the surface elevation η and of the slope components ∂η/∂x and ∂η/∂y at each instant can be represented by their corresponding ensemble-averaged RMS values calculated over the whole set of realizations. A different procedure was applied to determine the dominant wave frequency as a function of time elapsed from the blower activation. Since frequency Fourier analysis is inapplicable for unsteady forcing, continuous real Morlet wavelet transform was used, which decomposes a time-varying function into wavelets and offers good time and frequency localization. The wavelet 'spectrum', or map, was calculated for each realization; then, the resulting maps were averaged over the whole set of realizations. The scale corresponding to the maximum intensity of the ensemble-averaged map defines the dominant pseudo-frequency fdom at each instant. More details on the determination of fdom(t) are given in [67].

Quantitative Characterization of Wind-Wave Evolution with Time
Under steady wind forcing, the wave field in the present experimental facility is essentially threedimensional, as discussed in Section 3.4. It is instructive to examine the three-dimensional structure of the evolving in time wind-wave field by comparing the temporal variation of the simultaneously acquired surface elevation η(t) and of both components of surface slope, ∂η/∂x(t) and ∂η/∂y(t), as shown in Figure 22a. The records in this figure were scaled and shifted vertically for convenience.
For the initial t < 3.5 s (after the activation of the blower), there are no recognizable fluctuations in η and in both surface slope components. Thus, a comparison of Figures 21 and 22a indicates that the wave field for the conditions of Figure 22a starts to develop only after the constant wind velocity in the test section has been attained. For wind velocities lower than 8.5 m/s, forcing can effectively be presented by a step function. First visible disturbances at the water surface appear only about 4 s following the activation of the blower; additional lag of few seconds is observed until the steady state level of the surface elevation fluctuations is attained. However, the quasi-steady level of fluctuations of the surface slope component, ∂η/∂x, is attained much faster, already at t ≈ 4 s, just about 1 s after the inception of the first visible disturbances. The fluctuations of ∂η/∂y lag somewhat behind the development of slopes in the wind direction, so that for just less than 1 s, the wave field remains approximately unidirectional. Ensemble averaging of the accumulated data was carried out for each fetch x and target wind velocity U over all realizations as a function of time elapsed since the initiation of the blower in each realization. The characteristic amplitudes of the surface elevation η and of the slope components ∂η/∂x and ∂η/∂y at each instant can be represented by their corresponding ensemble-averaged RMS values calculated over the whole set of realizations. A different procedure was applied to determine the dominant wave frequency as a function of time elapsed from the blower activation. Since frequency Fourier analysis is inapplicable for unsteady forcing, continuous real Morlet wavelet transform was used, which decomposes a time-varying function into wavelets and offers good time and frequency localization. The wavelet 'spectrum', or map, was calculated for each realization; then, the resulting maps were averaged over the whole set of realizations. The scale corresponding to the maximum intensity of the ensemble-averaged map defines the dominant pseudo-frequency f dom at each instant. More details on the determination of f dom (t) are given in [67].

Quantitative Characterization of Wind-Wave Evolution with Time
Under steady wind forcing, the wave field in the present experimental facility is essentially three-dimensional, as discussed in Section 3.4. It is instructive to examine the three-dimensional structure of the evolving in time wind-wave field by comparing the temporal variation of the simultaneously acquired surface elevation η(t) and of both components of surface slope, ∂η/∂x(t) and ∂η/∂y(t), as shown in Figure 22a. The records in this figure were scaled and shifted vertically for convenience.
For the initial t < 3.5 s (after the activation of the blower), there are no recognizable fluctuations in η and in both surface slope components. Thus, a comparison of Figures 21 and 22a indicates that the wave field for the conditions of Figure 22a starts to develop only after the constant wind velocity in the test section has been attained. For wind velocities lower than 8.5 m/s, forcing can effectively be presented by a step function. First visible disturbances at the water surface appear only about 4 s following the activation of the blower; additional lag of few seconds is observed until the steady state level of the surface elevation fluctuations is attained. However, the quasi-steady level of fluctuations of the surface slope component, ∂η/∂x, is attained much faster, already at t ≈ 4 s, just about 1 s after the inception of the first visible disturbances. The fluctuations of ∂η/∂y lag somewhat behind the development of slopes in the wind direction, so that for just less than 1 s, the wave field remains approximately unidirectional. To support this observation further, the corresponding parameters averaged over an ensemble of 100 independent realizations are plotted in Figure 22b for the maximum wind velocity U = 10.5 m/s and fetch x = 340 cm. Ensemble averaging is performed relative to the instant of the blower initiation and denoted by angle brackets. It can be seen that the values of the characteristic surface slope <(∂η/∂x) 2 > 1/2 attain their maximum value practically simultaneously with the wind velocity in the test section. The slope oscillations in the crosswind direction lag by less than 1 s, while the surface elevation oscillations develop much slower. The delay in the emergence of crosswind undulations can plausibly be associated with the Langmuir circulation excited in water by accelerated wind [45,98,99], although the attempts to detect this circulation in [67] were indecisive. Quasi-steady levels of slope fluctuations in both directions are of comparable values and were attained at similar times, which were notably earlier than that of η.
The temporal evolution of the main wave parameters ensemble-averaged as shown in Figure  22b is studied in more detail in Figure 23. The results are plotted for three fetches and three wind velocities; additional data can be found in [67]. The instantaneous wave amplitudes are characterized by the RMS values of the surface elevation, <η 2 > 1/2 ; the values of <(∂η/∂x) 2 > 1/2 serve as an indicator of the instantaneous wave steepness (as seen in Figure 22, the behavior of <(∂η/∂y) 2 > 1/2 is similar to that of <(∂η/∂x) 2 > 1/2 , and is therefore not plotted). The dominant frequency is based on the ensemble-averaged wavelet maps obtained for the surface elevation records η(t); the dominant wavelength λ(t) is calculated from fdom(t) using the empirical dispersion relation shown in Equation (15).
The temporal resolution of 3.3 ms in Figure 23 enables distinguishing between the stages in the wave field evolution. At all fetches and wind velocities, the water surface remains effectively undisturbed for about 4 s following the initiation of the blower; the initial ripples during this period do not yet appear. This lag between the response of the water surface and initiation of the blower is sufficient for air in the test section to accelerate to its final prescribed velocity. The only possible exception is the highest flow rate employed in this study, as seen in Figure 21. Once wind is switched on in an effectively impulsive manner, the rapid growth of very short waves is observed along the whole length of the test section. The RMS values of the surface elevation of those ripples do not exceed approximately 1 mm. As already noticed in Figure 22b, this fast initial wave growth is accompanied by a much faster nearly impulsive increase in the RMS values of the characteristic surface slope <(∂η/∂x) 2 > 1/2 , which rapidly approach the steady-state values. Those values of <(∂η/∂x) 2 > 1/2 vary in the range of 0.17 to about 0.22; at a given wind velocity, they are nearly independent of the fetch and increase somewhat with U. To support this observation further, the corresponding parameters averaged over an ensemble of 100 independent realizations are plotted in Figure 22b for the maximum wind velocity U = 10.5 m/s and fetch x = 340 cm. Ensemble averaging is performed relative to the instant of the blower initiation and denoted by angle brackets. It can be seen that the values of the characteristic surface slope <(∂η/∂x) 2 > 1/2 attain their maximum value practically simultaneously with the wind velocity in the test section. The slope oscillations in the crosswind direction lag by less than 1 s, while the surface elevation oscillations develop much slower. The delay in the emergence of crosswind undulations can plausibly be associated with the Langmuir circulation excited in water by accelerated wind [45,98,99], although the attempts to detect this circulation in [67] were indecisive. Quasi-steady levels of slope fluctuations in both directions are of comparable values and were attained at similar times, which were notably earlier than that of η.
The temporal evolution of the main wave parameters ensemble-averaged as shown in Figure 22b is studied in more detail in Figure 23. The results are plotted for three fetches and three wind velocities; additional data can be found in [67]. The instantaneous wave amplitudes are characterized by the RMS values of the surface elevation, <η 2 > 1/2 ; the values of <(∂η/∂x) 2 > 1/2 serve as an indicator of the instantaneous wave steepness (as seen in Figure 22, the behavior of <(∂η/∂y) 2 > 1/2 is similar to that of <(∂η/∂x) 2 > 1/2 , and is therefore not plotted). The dominant frequency is based on the ensemble-averaged wavelet maps obtained for the surface elevation records η(t); the dominant wavelength λ(t) is calculated from f dom (t) using the empirical dispersion relation shown in Equation (15).
The temporal resolution of 3.3 ms in Figure 23 enables distinguishing between the stages in the wave field evolution. At all fetches and wind velocities, the water surface remains effectively undisturbed for about 4 s following the initiation of the blower; the initial ripples during this period do not yet appear. This lag between the response of the water surface and initiation of the blower is sufficient for air in the test section to accelerate to its final prescribed velocity. The only possible exception is the highest flow rate employed in this study, as seen in Figure 21. Once wind is switched on in an effectively impulsive manner, the rapid growth of very short waves is observed along the whole length of the test section. The RMS values of the surface elevation of those ripples do not exceed approximately 1 mm. As already noticed in Figure 22b, this fast initial wave growth is accompanied by a much faster nearly impulsive increase in the RMS values of the characteristic surface slope <(∂η/∂x) 2 > 1/2 , which rapidly approach the steady-state values. Those values of <(∂η/∂x) 2 > 1/2 vary in the range of 0.17 to about 0.22; at a given wind velocity, they are nearly independent of the fetch and increase somewhat with U. The wavelet-derived characteristic frequency of the initial ripples exceeds about fdom = 15 Hz at all fetches and wind velocities. During the fast growth stage of the ripples, their frequencies decrease. The eventual quasi-steady characteristic frequencies of the wind waves decrease with fetches; for a given fetch, they decrease with wind velocity. The variation of the dominant wavelength, λdom, is associated with that of fdom, according to the dispersion relation shown in Equation (15). Contrary to the fast variation of the surface slope components up to their steady-state values, the evolution of <η 2 > 1/2 , fdom, and λdom is much slower, and at all fetches and wind velocities, it occurs at time scales that are identical to those that characterize the temporal variation of <η 2 > 1/2 . The time scales of the evolution increase with fetch and wind velocity.

A Simple Scenario for Intial Wind-Wave Evolution
The evolution in time of the characteristic wave amplitude, <η 2 > 1/2 , and of the instantaneous dominant frequency, fd, are compared in Figure 24 at two fetches x and for two values of wind velocity U. The comparison reveals that for a given U, at the instant when the instantaneous values of <η 2 > 1/2 attain the value corresponding to the steady state at the shorter fetch, the corresponding dominant frequencies fd also attain the steady-state value at the shorter fetch. Similar results were obtained for the other pairs of fetches and wind velocities considered in this study. The wavelet-derived characteristic frequency of the initial ripples exceeds about f dom = 15 Hz at all fetches and wind velocities. During the fast growth stage of the ripples, their frequencies decrease. The eventual quasi-steady characteristic frequencies of the wind waves decrease with fetches; for a given fetch, they decrease with wind velocity. The variation of the dominant wavelength, λ dom , is associated with that of f dom , according to the dispersion relation shown in Equation (15). Contrary to the fast variation of the surface slope components up to their steady-state values, the evolution of <η 2 > 1/2 , f dom , and λ dom is much slower, and at all fetches and wind velocities, it occurs at time scales that are identical to those that characterize the temporal variation of <η 2 > 1/2 . The time scales of the evolution increase with fetch and wind velocity.

A Simple Scenario for Intial Wind-Wave Evolution
The evolution in time of the characteristic wave amplitude, <η 2 > 1/2 , and of the instantaneous dominant frequency, f d , are compared in Figure 24 at two fetches x and for two values of wind velocity U. The comparison reveals that for a given U, at the instant when the instantaneous values of <η 2 > 1/2 attain the value corresponding to the steady state at the shorter fetch, the corresponding dominant frequencies f d also attain the steady-state value at the shorter fetch. Similar results were obtained for the other pairs of fetches and wind velocities considered in this study. These observations can be explained adopting the following simplified pattern of the spatial and temporal evolution of waves under impulsive wind forcing. As seen in Figure 23, it can be assumed that the initial disturbances at the water surface upon the application of wind are approximately homogeneous along the whole test section. This initial wave field contains a broad spectrum of harmonics. Under the action of wind, different harmonics propagate with their respective group velocities cg, and grow in time and in space, until they attain equilibrium values of amplitude and steepness for the given steady wind forcing. Longer waves apparently propagate faster and attain larger amplitudes for a given characteristic steepness; however, attaining equilibrium for longer wave harmonics requires longer times and occurs at larger fetches. Therefore, once a harmonic with the frequency fi attains its equilibrium amplitude at a certain fetch x(fi), the growth stage terminates; its amplitude remains constant, and does not vary fast at longer fetches, xi > x(fi). Adopting this conjecture, the duration of the growth stage at each fetch and wind velocity, tgr(x), depends on the dominant frequency fdom, and can thus be estimated as the time required for the harmonic at fdom to reach the fetch x: where the group velocity is determined using the dispersion relations shown in Equations (14) and (15). The values of the total duration of the growth stage, tgr, calculated from Equation (19) using the empirical values of group velocity for each dominant wavenumber, are compared in Figure 25 with the actual duration of the growth stage ttot estimated from the experimental data presented in Figure  23 and in [67]. Note that the initial time reference point in Figure 25 is taken at the instant when the blower capacity attains its maximum, i.e., 3.5-5.5 s after the activation of the blower, depending on the U, as shown in Figure 21. Figure 25 demonstrates that the relation shown in Equation (19) yields reasonable estimates of the duration of the initial growth stage for all fetches and wind velocities. It should be stressed that this scenario apparently neglects the nonlinear effects. The results of Figures  24 and 25 imply that while young wind waves are strongly nonlinear, the initial growth of short harmonics up to attaining the quasi-equilibrium state is not yet affected strongly by nonlinearity, and can be approximated as a linear process. These observations can be explained adopting the following simplified pattern of the spatial and temporal evolution of waves under impulsive wind forcing. As seen in Figure 23, it can be assumed that the initial disturbances at the water surface upon the application of wind are approximately homogeneous along the whole test section. This initial wave field contains a broad spectrum of harmonics. Under the action of wind, different harmonics propagate with their respective group velocities c g , and grow in time and in space, until they attain equilibrium values of amplitude and steepness for the given steady wind forcing. Longer waves apparently propagate faster and attain larger amplitudes for a given characteristic steepness; however, attaining equilibrium for longer wave harmonics requires longer times and occurs at larger fetches. Therefore, once a harmonic with the frequency f i attains its equilibrium amplitude at a certain fetch x(f i ), the growth stage terminates; its amplitude remains constant, and does not vary fast at longer fetches, x i > x(f i ). Adopting this conjecture, the duration of the growth stage at each fetch and wind velocity, t gr (x), depends on the dominant frequency f dom , and can thus be estimated as the time required for the harmonic at f dom to reach the fetch x: where the group velocity is determined using the dispersion relations shown in Equations (14) and (15). The values of the total duration of the growth stage, t gr , calculated from Equation (19) using the empirical values of group velocity for each dominant wavenumber, are compared in Figure 25 with the actual duration of the growth stage t tot estimated from the experimental data presented in Figure 23 and in [67]. Note that the initial time reference point in Figure 25 is taken at the instant when the blower capacity attains its maximum, i.e., 3.5-5.5 s after the activation of the blower, depending on the U, as shown in Figure 21. Figure 25 demonstrates that the relation shown in Equation (19) yields reasonable estimates of the duration of the initial growth stage for all fetches and wind velocities. It should be stressed that this scenario apparently neglects the nonlinear effects. The results of Figures 24 and 25 imply that while young wind waves are strongly nonlinear, the initial growth of short harmonics up to attaining the quasi-equilibrium state is not yet affected strongly by nonlinearity, and can be approximated as a linear process. Atmosphere 2019, 10, x FOR PEER REVIEW 34 of 52 Figure 25. The calculated according to the model duration of the initial wave growth process, tgr, compared to the measured growth duration, ttot.

The Relation of the Recent Experimental Results to the Existing Theoretical Approaches
The extensive experimental data obtained in [67] are now analyzed in view of previous experimental and theoretical results. Notable changes of the slopes of the curves representing the temporal variation of the different parameters can be identified in Figure 23. Those sharp slope changes occur essentially simultaneously for all parameters plotted in that figure, thus suggesting division of the temporal evolution of the wave field into distinct stages. To define those stages, the variation with the elapsed time of a single representative parameter is plotted in Figure 26. For that purpose, the variation with time elapsed from the initiation of the blower of the characteristic instantaneous wave amplitude, <η 2 > 1/2 , was selected. By the time the wavelets become detectable at the elapsed time t1, the wind velocity in the test section has already attained its prescribed value; see Figure 21. For t > t1, the short ripples initially grow fast, but the growth rate decreases abruptly at t = t2. The characteristic wave amplitudes at that instant are still very small-less than 0.5 mm; however, the values of their representative steepness at t = t2 already approach their maximum limit (cf. Figure 23). The wave growth continues at a significantly slower rate for t2 < t < t3. Then, at t3, the rate of wave amplitude growth increases again, until the temporal wave evolution process ceases at t = t4. This instant represents the total duration of the wave growth at the prescribed location and wind velocity; it is denoted in Figure 25 as ttot. A Figure 25. The calculated according to the model duration of the initial wave growth process, t gr , compared to the measured growth duration, t tot .

The Relation of the Recent Experimental Results to the Existing Theoretical Approaches
The extensive experimental data obtained in [67] are now analyzed in view of previous experimental and theoretical results. Notable changes of the slopes of the curves representing the temporal variation of the different parameters can be identified in Figure 23. Those sharp slope changes occur essentially simultaneously for all parameters plotted in that figure, thus suggesting division of the temporal evolution of the wave field into distinct stages. To define those stages, the variation with the elapsed time of a single representative parameter is plotted in Figure 26. For that purpose, the variation with time elapsed from the initiation of the blower of the characteristic instantaneous wave amplitude, <η 2 > 1/2 , was selected.

The Relation of the Recent Experimental Results to the Existing Theoretical Approaches
The extensive experimental data obtained in [67] are now analyzed in view of previous experimental and theoretical results. Notable changes of the slopes of the curves representing the temporal variation of the different parameters can be identified in Figure 23. Those sharp slope changes occur essentially simultaneously for all parameters plotted in that figure, thus suggesting division of the temporal evolution of the wave field into distinct stages. To define those stages, the variation with the elapsed time of a single representative parameter is plotted in Figure 26. For that purpose, the variation with time elapsed from the initiation of the blower of the characteristic instantaneous wave amplitude, <η 2 > 1/2 , was selected. By the time the wavelets become detectable at the elapsed time t1, the wind velocity in the test section has already attained its prescribed value; see Figure 21. For t > t1, the short ripples initially grow fast, but the growth rate decreases abruptly at t = t2. The characteristic wave amplitudes at that instant are still very small-less than 0.5 mm; however, the values of their representative steepness at t = t2 already approach their maximum limit (cf. Figure 23). The wave growth continues at a significantly slower rate for t2 < t < t3. Then, at t3, the rate of wave amplitude growth increases again, until the temporal wave evolution process ceases at t = t4. This instant represents the total duration of the wave growth at the prescribed location and wind velocity; it is denoted in Figure 25 as ttot. A By the time the wavelets become detectable at the elapsed time t 1 , the wind velocity in the test section has already attained its prescribed value; see Figure 21. For t > t 1 , the short ripples initially grow fast, but the growth rate decreases abruptly at t = t 2 . The characteristic wave amplitudes at that instant are still very small-less than 0.5 mm; however, the values of their representative steepness at t = t 2 already approach their maximum limit (cf. Figure 23). The wave growth continues at a significantly slower rate for t 2 < t < t 3 . Then, at t 3 , the rate of wave amplitude growth increases again, until the temporal wave evolution process ceases at t = t 4 . This instant represents the total duration of the wave growth at the prescribed location and wind velocity; it is denoted in Figure 25 as t tot . A summary of the estimated from the experimental results transition times t i , i = 1, . . . , 4 for three fetches and different wind velocities can be found in [67]. The distinct stages in the wind-wave evolution with the transition at the characteristic times t i identified in Figure 26 can plausibly be attributed to the different mechanisms that govern the wave growth at each stage. An attempt is made to relate the experimental results to the existing theoretical models.
The initial growth of ripples excited by such an impulsive wind forcing is plotted in semi-logarithmic coordinates in Figure 27a for lower wind velocities. Results for stronger winds are presented in [67]. Figure 27a demonstrates that the growth with time of the energy of the initial ripples <η 2 > with time is exponential; thus, it can be presented as: where η 2 0 is the initial disturbance energy, and β is the growth parameter. The slopes in Figure 27a do not vary much with fetch and wind velocity U; similar results were presented in [67] for additional values of x and U. The values of β in Equation (20) for all cases considered do not differ significantly from β ≈ 10 s −1 . The duration of the exponential growth in Figure 27a does not exceed about 200-300 ms, requiring a sufficiently high sampling rate of 300 Hz, which was used in this study for the unambiguous determination of the growth parameter β. The exponential growth of initial short waves in Figure 27a may be attributed to the linear instability mechanism that becomes essential at instant t 1 ; see Figure 26. The growth rate parameters β obtained in [99] indeed compare favorably with the computations and the experimental results reported by Kawai [28]. Since the forcing conditions in those two studies were different, extrapolation of the results of [28] was needed to enable quantitative comparison. summary of the estimated from the experimental results transition times ti, i = 1,…,4 for three fetches and different wind velocities can be found in [67]. The distinct stages in the wind-wave evolution with the transition at the characteristic times ti identified in Figure 26 can plausibly be attributed to the different mechanisms that govern the wave growth at each stage. An attempt is made to relate the experimental results to the existing theoretical models.
(a) (b) Figure 27. Appearance of initial ripples: (a) Exponential growth of ripples' energy; (b) the growth of ripples' slope components. Vertical broken lines denote transition instants t1 and t2, as defined in Figure 26.
The initial growth of ripples excited by such an impulsive wind forcing is plotted in semilogarithmic coordinates in Figure 27a for lower wind velocities. Results for stronger winds are presented in [67]. Figure 27a demonstrates that the growth with time of the energy of the initial ripples <η 2 > with time is exponential; thus, it can be presented as: where 〈 0 2 〉 is the initial disturbance energy, and β is the growth parameter. The slopes in Figure 27a do not vary much with fetch and wind velocity U; similar results were presented in [67] for additional values of x and U. The values of β in Equation (20) for all cases considered do not differ significantly from β ≈ 10 s −1 . The duration of the exponential growth in Figure 27a does not exceed about 200-300 ms, requiring a sufficiently high sampling rate of 300 Hz, which was used in this study for the unambiguous determination of the growth parameter β. The exponential growth of initial short waves in Figure 27a may be attributed to the linear instability mechanism that becomes essential at instant t1; see Figure 26. The growth rate parameters β obtained in [99] indeed compare favorably with the computations and the experimental results reported by Kawai [28]. Since the forcing Figure 27. Appearance of initial ripples: (a) Exponential growth of ripples' energy; (b) the growth of ripples' slope components. Vertical broken lines denote transition instants t 1 and t 2 , as defined in Figure 26.
The exponential growth stage of the initial ripples that terminates at the instant t 2 can be attributed to the significant wave nonlinearity manifested by the instantaneous values of <(∂η/∂x) 2 > 1/2 that attain their maximum quasi-steady values at times close to t 2 . An additional factor may be instrumental in rendering the application of the linear viscous instability theory of unidirectional waves inapplicable at t > t 2 . As seen in Figure 22, the growth of <(∂η/∂x) 2 > 1/2 somewhat precedes that of <(∂η/∂y) 2 > 1/2 . A closer look at the variation with the elapsed time of these two parameters during the initial evolution stages is plotted in Figure 27b for the conditions corresponding to those in Figure 26. The slopes in the crosswind direction indeed grow visibly slower than those in the along-wind direction, and attain their quasi-steady-state values later. The initial disturbances during the exponential growth stage at t < t 2 can be seen as largely unidirectional, as assumed in all linear stability studies, including [28]. However, this assumption ceases to be valid at longer times. The apparent failure of the linear instability approach to describe the evolution of waves under impulsive wind forcing beyond the very initial phase suggests invoking the essentially nonlinear wave generation theory offered by Phillips [21] for the times following the appearance of the initial ripples. The theory suggests a two-stage temporal evolution of wind waves due to random pressure fluctuations. In our experiments, the 'initial growth stage' according to the approach by Phillips may be attributed to elapsed times t 2 < t < t 3 . Phillips does not provide an expression for the description of the wind-wave growth at this stage; rather, his theory states that the mean values of <η 2 > grow linearly with time. To examine this conjecture, the measured dependence of the mean wave energy <η 2 > on time is plotted in Figure 28a as a function of (t-t 2 ). In spite of a considerable scatter in the data that is more prominent for stronger winds, the wave energy indeed seems to increases linearly with time. The duration of this growth stage decreases with the increase in the wind velocity U; see also Figure 23. theory offered by Phillips [21] for the times following the appearance of the initial ripples. The theory suggests a two-stage temporal evolution of wind waves due to random pressure fluctuations. In our experiments, the 'initial growth stage' according to the approach by Phillips may be attributed to elapsed times t2 < t < t3. Phillips does not provide an expression for the description of the wind-wave growth at this stage; rather, his theory states that the mean values of <η 2 > grow linearly with time. To examine this conjecture, the measured dependence of the mean wave energy <η 2 > on time is plotted in Figure 28a as a function of (t-t2). In spite of a considerable scatter in the data that is more prominent for stronger winds, the wave energy indeed seems to increases linearly with time. The duration of this growth stage decreases with the increase in the wind velocity U; see also Figure 23. At t = t3, the wave energy growth rate increases notably. Adopting the two-stage theory of Phillips, the wind-wave field evolution at t > t3 may be associated with 'the principal stage'. This stage terminates when the quasi-steady state is attained at t ≈ t4 for each fetch and wind velocity. Waves at fetches x = 220 cm and x = 340 cm have lengths exceeding 7 cm, as shown in Figure 23, and thus can be considered as purely gravity waves practically unaffected by capillarity. Thus, the application of the following relation obtained by Phillips for the temporal growth of gravity waves at the principal stage of development: is justified for our experimental conditions and remote fetches. In Equation (21), 2 ̅̅̅ represents the mean square pressure fluctuations at the surface, g is the acceleration due to gravity, and ρw is the water density. For the resonance conditions that govern wave field evolution at this stage, the convection velocity of the pressure fluctuations, Uc, is assumed to correspond to the phase velocity of water waves. It is further assumed that the turbulent pressure fluctuations are related to interfacial shear stress, = * 2 , where * is the friction velocity at the air-water interface, so that 2 ̅̅̅~ * 4 . The convection velocity of the pressure fluctuations, Uc, is related in [21] to wind velocity, Uc ≈ U. The following relation is adopted between the friction velocity and the wind velocity U=18 * . Thus, the relation shown in Equation (21) is rewritten as: Thus, the theory of Phillips predicts the linear time growth of the gravity wind-wave energy, while the growth rate is proportional to U 3 . However, the characteristic wave velocities for shortgravity waves in the wavelength domain pertinent to the present experiments, 10 cm < λ < 30 cm, were measured in the present experimental facility [51], and were found to vary only weakly; they At t = t 3, the wave energy growth rate increases notably. Adopting the two-stage theory of Phillips, the wind-wave field evolution at t > t 3 may be associated with 'the principal stage'. This stage terminates when the quasi-steady state is attained at t ≈ t 4 for each fetch and wind velocity. Waves at fetches x = 220 cm and x = 340 cm have lengths exceeding 7 cm, as shown in Figure 23, and thus can be considered as purely gravity waves practically unaffected by capillarity. Thus, the application of the following relation obtained by Phillips for the temporal growth of gravity waves at the principal stage of development: is justified for our experimental conditions and remote fetches. In Equation (21), p 2 represents the mean square pressure fluctuations at the surface, g is the acceleration due to gravity, and ρ w is the water density. For the resonance conditions that govern wave field evolution at this stage, the convection velocity of the pressure fluctuations, U c , is assumed to correspond to the phase velocity of water waves. It is further assumed that the turbulent pressure fluctuations are related to interfacial shear stress, τ = ρ a u 2 * , where u * is the friction velocity at the air-water interface, so that p 2 ∼ u 4 * . The convection velocity of the pressure fluctuations, U c , is related in [21] to wind velocity, U c ≈ U. The following relation is adopted between the friction velocity and the wind velocity U = 18u * . Thus, the relation shown in Equation (21) is rewritten as: Thus, the theory of Phillips predicts the linear time growth of the gravity wind-wave energy, while the growth rate is proportional to U 3 . However, the characteristic wave velocities for short-gravity waves in the wavelength domain pertinent to the present experiments, 10 cm < λ < 30 cm, were measured in the present experimental facility [51], and were found to vary only weakly; they remain substantially lower than the wind velocity, and are essentially independent of U. Thus, the assumption of constant convection velocity U c is more appropriate for the conditions of the present experiments. When this assumption is incorporated into Equation (21), the relation shown in Equation (20) is replaced by: In the framework used in Phillips' approach, the dimensional coefficient C in Equation (23) is supposed to be constant. The validity of Equation (23) is examined in Figure 28b, where the measured variation of the normalized by U 4 wave energy, η 2 /U 4 is plotted as a function of the time elapsed since the initiation of the principal stage of wind-wave development, t-t 3 . The results for two fetches (x = 220 cm and x = 340 cm) and three wind velocities (U = 8.5 m/s, U = 9.5 m/s and U = 10.5 m/s) are well approximated by straight lines. Thus, the linear with time growth of the wind-wave energy during the 'principal stage' of evolution is confirmed experimentally as well. Moreover, the results of Figure 28b strongly support the U 4 dependence of the growth rate on wind velocity. For the longer fetch, x = 340 cm, practically identical slopes are obtained for different wind velocities. For the shorter fetch, x = 220 cm, the scatter in the slopes is somewhat larger; there is no identifiable trend in the normalized wave-energy growth rate on the wind velocity U. For both fetches, the slopes in Figure 28b do not differ significantly.
A different attempt to gain insight into a general pattern of wave energy growth, from the initial calm water surface to the eventual steady state at each fetch, is made in Figure 28c. The plotted ensemble-averaged values of <η 2 (t)> at fetches x = 220 cm and x = 340 cm in this panel are normalized by the corresponding steady-state mean values for each one of the three higher wind velocities in the experiments (U = 8.5 m/s, 9.5 m/s, and 10.5 m/s). The similarity of all curves at each fetch for different wind velocities U is obvious in Figure 28c. Therefore, the results for all wind velocities at each fetch are averaged; the results are plotted in Figure 28c by solid lines. The shapes of both solid curves are qualitatively similar; the three distinct stages in the wave evolution identified in Figure 26 are visible in this panel. The exponentially growing initial short wavelets appear simultaneously at both fetches and last up to t 2 ≈ 5 s. Then, a relatively slow and close to linear growth of the energy of gravity-capillary waves follows; this stage terminates at t 3 ≈ 8 s for the shorter fetch x = 220 cm, and at t 3 ≈ 13 for x = 340 cm. The growth rates during the time interval t 2 < t < t 3 are apparently fetch-dependent. Starting from t 3 , the growth rate of <η 2 (t)> increase sharply; the rates of wave energy change ∂<η 2 (t)>/∂t become essentially constant and apparently similar for both fetches. Finally, the fetch-dependent quasi-steady state is attained at t = t 4 . Thus, the transition times between the development stages t 2 and t 3 as well as the total duration of the evolution are independent of wind velocity U and depend on fetch only.

The Three-Dimensional Aspects of the Evolving Wind-Wave Field
As discussed in relation to Figure 27b, the evolving under impulsive wind forcing wave field becomes essentially three-dimensional quite fast: approximately at the end of the initial exponential growth stage. The temporal variation of the 3D surface geometry during the whole wave growth can be characterized in terms of the probability density function (PDF) of the instantaneous azimuthal angle determined from the instantaneous components of the surface slope measured by LSG, ∂η/∂x, and ∂η/∂x, as performed in Figure 17 for wind waves under steady forcing. The PDF obtained for ensembles of instantaneous azimuthal angles accumulated at several instants after activation of the blower are presented in Figure 29a for a short fetch of x = 120 cm and a low wind velocity of U = 6.5 m/s, and in Figure 29b  In both panels of Figure 29, the probability distribution of the slope inclination angle θ is symmetric around the two peaks at θ = 0° and θ = 180°. The difference in peak heights indicates the upwind-downwind asymmetry of the wave shape; a flatter distribution corresponds to wider directional spreading. As in Figure 17, all the inclination directions of the instantaneous surface have a non-negligible probability at all instants and for both panels of Figure 29. The probability distributions of the angles θ vary strongly with time. The initial recognizable ripples are strongly asymmetric, with the leeward part notably longer than the upwind part. As the waves grow with time, the probability distributions become more symmetric; near perfect front-back symmetry is observed at t = 4.4 s in Figure 29a, and at t = 4.2 s in Figure 29b. Once the steady state is attained, the directional spreading becomes wider, and the downwind part of the wave is typically shorter than its upwind part, which is in agreement with Figure 17.

Relation to Wind Waves at Natural Scales
A systematic study of the evolution of surface water waves under steady or impulsively applied wind forcing in our experimental facility resulted in the extensive body of experimental results summarized in Sections 3 and 4. The relevance of laboratory-scale experiments to wind waves in natural environment, as questioned in the Introduction, is now discussed in view of those results.
First, the effect of the test-section dimensions has to be considered. The length of the dominant waves λ = 2π/k encountered in the experiments for all the wind-forcing conditions employed and at practically at all fetches satisfies deep-water condition, kd > π, in which d is the water depth [100]. For shorter fetches and relatively weak wind forcing, the dominant waves are affected by surface tension as well as by the induced surface current. These effects are taken into account in applying the empirical dispersion relation for gravity-capillary waves, as shown in Equations (14) and (15). However, for longer waves observed at larger fetches and stronger steady wind forcing, those effects become negligible, and the deep-water dispersion relation for gravity waves is satisfied around the peak spectral frequency. In that sense, the waves in the ocean and in the laboratory differ only in scale. The size of the facility apparently limits the maximum length of wind waves. The lateral In both panels of Figure 29, the probability distribution of the slope inclination angle θ is symmetric around the two peaks at θ = 0 • and θ = 180 • . The difference in peak heights indicates the upwind-downwind asymmetry of the wave shape; a flatter distribution corresponds to wider directional spreading. As in Figure 17, all the inclination directions of the instantaneous surface have a non-negligible probability at all instants and for both panels of Figure 29. The probability distributions of the angles θ vary strongly with time. The initial recognizable ripples are strongly asymmetric, with the leeward part notably longer than the upwind part. As the waves grow with time, the probability distributions become more symmetric; near perfect front-back symmetry is observed at t = 4.4 s in Figure 29a, and at t = 4.2 s in Figure 29b. Once the steady state is attained, the directional spreading becomes wider, and the downwind part of the wave is typically shorter than its upwind part, which is in agreement with Figure 17.

Relation to Wind Waves at Natural Scales
A systematic study of the evolution of surface water waves under steady or impulsively applied wind forcing in our experimental facility resulted in the extensive body of experimental results summarized in Sections 3 and 4. The relevance of laboratory-scale experiments to wind waves in natural environment, as questioned in the Introduction, is now discussed in view of those results.
First, the effect of the test-section dimensions has to be considered. The length of the dominant waves λ = 2π/k encountered in the experiments for all the wind-forcing conditions employed and at practically at all fetches satisfies deep-water condition, kd > π, in which d is the water depth [100]. For shorter fetches and relatively weak wind forcing, the dominant waves are affected by surface tension as well as by the induced surface current. These effects are taken into account in applying the empirical dispersion relation for gravity-capillary waves, as shown in Equations (14) and (15). However, for longer waves observed at larger fetches and stronger steady wind forcing, those effects become negligible, and the deep-water dispersion relation for gravity waves is satisfied around the peak spectral frequency. In that sense, the waves in the ocean and in the laboratory differ only in scale. The size of the facility apparently limits the maximum length of wind waves. The lateral dimensions of the test section may plausibly impose limitations on the angular distribution of wave propagation directions. In a deterministic three-dimensional wave field, the conditions imposed by the sidewall boundary may indeed affect the waves strongly; a discrete spectrum in the spanwise direction is often observed under certain forcing conditions [101,102]. However, the waves in our experiments are essentially three-dimensional and random; see Figure 13. The random nature of the waves results in short correlation lengths, as demonstrated in Figures 14 and 15, thus effectively preventing the waves in the central part of the test section from "feeling" the presence of the walls. It is often assumed that the so-called wave age, c/ u * , largely determines the behavior of wind waves. The friction velocities u * in laboratory experiments and in the ocean typically are of the same order of magnitude, while wave phase velocities c are much smaller, even in the largest available facilities as compared to nature. Thus, there is no way to simulate the wave age of long ocean waves in any reasonably scaled facility. Thus, the applicability of the laboratory results to wave scales typical for natural reservoirs should be judged based on the quantitative comparison of the measured statistical parameters after the application of appropriate scaling.
The short initial ripples that appear during the early evolution stage of waves generated by impulsive wind forcing in the laboratory hardly differ from those observed in nature. The presence of sidewalls is not essential during the initial evolution stages due to the randomness and three-dimensionality of the waves. Therefore, the results of experiments obtained in laboratory facilities of limited size can be directly applied to gain understanding of the early stages of a wave field emerging in nature under gusty wind forcing.
The statistical and spectral parameters of the wave field in our facility at later stages of evolution along the test section under steady wind forcing were compared with the available theoretical and experimental data accumulated in both laboratory and field studies. The friction velocity u * was used as a characteristic wind velocity, being a more natural characteristic than the U 10 velocity that is often used in field measurements. Therefore, special efforts were invested in the accurate determination of the friction velocity in our facility under different steady forcing conditions, as described in Section 3.1. Dependencies on the dimensionless fetch of the integral dimensionless wave parameters, such as the characteristic surface elevation amplitude, the dimensionless peak frequency, and the spectral width, are plotted in Figure 8. Except for lower wind velocities and very short fetches characterized by the presence of gravity-capillary waves of small amplitudes, the resulting empirical relations are in agreement with laboratory results in much larger facilities [77,78] and with field observations [79,80]. These relations are also in general agreement with the 3/2 power law suggested by Toba [81] based on field data analysis. Moreover, as specified in Section 3.2, deviations of the surface elevation variation from Gaussianity characterized by higher distribution moments of skewness, λ 3 , and kurtosis, λ 4 -see Equation (9) and Figure 9a,b-as well as by the wave exceedance distributions shown in Figure 9c,d, are also in general agreement with data available elsewhere.
The comparison of spectral shapes rather than integral parameters such as the dimensionless spectral width ν as defined by Equation (8) provides a more sensitive measure to assess the similarity between the waves in the laboratory and in nature. The power spectra in the present experiments that covered up to five decades were obtained by averaging a large number of statistically independent estimates, taking advantage of steady conditions in the experimental facility. It is impossible to accumulate data ensembles of comparable size in field experiments.
The observed decrease in both the peak frequency and the spectral width with the increasing steady wind speed, or with fetch, is characteristic for random wind-generated waves. Higher frequency components grow initially more rapidly with their amplitudes attaining saturation, which is followed It should be stressed that the similarity of wind-wave spectra measured in the laboratory to the JONSWAP shape was noticed already by Plant [104]. The spectral shape defined by Equation (27) was fitted to the data plotted in Figure 30. The coefficients γ, σR, and σL were obtained from the fit of the spectral shape shown in Equation (27) for each fetch using the corresponding dimensionless spectral widths ν averaged over wind velocities; see Figure 8c. Good agreement of the fitted shape with the data at both sides of the peak value f  = 0 is obtained. The values of the peak enhancement coefficient γ remain nearly constant at about γ = 3.58 for all fetches. The coefficients σR ≈ 0.136 and σL ≈ 0.148 were also essentially independent of fetch. In view of insensitivity of the coefficients of Equation (27) to fetch and wind velocity, a single dimensionless normalized JONSWAP-type spectrum averaged over all fetches, and using an averaged value of ν = 0.173, is plotted in Figure 30. It should be stressed that the similarity of wind-wave spectra measured in the laboratory to the JONSWAP shape was noticed already by Plant [104]. The spectral shape defined by Equation (27) was fitted to the data plotted in Figure 30. The coefficients γ, σ R , and σ L were obtained from the fit of the spectral shape shown in Equation (27) for each fetch using the corresponding dimensionless spectral widths ν averaged over wind velocities; see Figure 8c. Good agreement of the fitted shape with the data at both sides of the peak value ∆f = 0 is obtained. The values of the peak enhancement coefficient γ remain nearly constant at about γ = 3.58 for all fetches. The coefficients σ R ≈ 0.136 and σ L ≈ 0.148 were also essentially independent of fetch. In view of insensitivity of the coefficients of Equation (27) to fetch and wind velocity, a single dimensionless normalized JONSWAP-type spectrum averaged over all fetches, and using an averaged value of ν = 0.173, is plotted in Figure 30.
The agreement obtained for the longer waves at the left side of the spectrum is quite good; the fit for the higher-frequency part of the spectrum overestimates somewhat the experimental data for ∆f > 1.4. The sharpness of the spectral peak γ in the fit is only slightly higher than the value γ = 3.3 reported for mature waves in JONSWAP experiments [79]. However, the coefficients σ in the fit of Figure 30 that characterize spectral shape steepness are somewhat higher than the generally accepted σ R ≈ 0.09 and σ L ≈ 0.07. The peak power of the JONSWAP spectrum is determined by the coefficient β, which for the adopted spectral shape and the measured significant wave heights presented in [75] has a value of β ≈ 0.2. This value is very close to β ≈ 0.21, which was suggested by Goda [103]. As discussed in relation to Figure 7, the behavior of the high-frequency spectral tail also does not differ significantly from that obtained in previous theoretical and field studies.
It transpires from this discussion that while certain limitations imposed by size inevitably exist, experiments on wind-generated water waves in facilities of modest dimensions may yield results that are applicable at considerably larger scales. Experiments in small-size and mid-size laboratory facilities, besides being relatively inexpensive, allow a detailed investigation with high temporal and spatial resolution of the diverse processes that are associated with the generation of water waves by wind under repeatable steady or unsteady forcing conditions in a controlled environment. These benefits of laboratory experiments cannot be met in field measurements.

Evolution of Wind Waves in Space and Time
Wind blowing over a water surface in nature is always turbulent, as both the magnitude and the direction of its mean velocity vary randomly in time as well as in space. The resulting wind-generated wave field is extremely complicated; this complexity hampers its quantitative analysis. Wind-generated waves are essentially random not only in a large-scale natural environment, but also in much smaller laboratory facilities. To extract reliable statistical data from measurements, advantage was taken in this study of the ability to control wind forcing. Two generally adopted simplified forcing settings were considered. In Section 3, extensive results on wind waves under steady wind were presented; while Section 4 dealt with waves excited by wind that could be seen as impulsively imposed (the relatively minor effects of the finite response time of the system were discussed in detail in [67]).
In laboratory experiments on wind waves, steady wind forcing is most often applied. The assumption of steady mean wind velocity is also routinely adopted in field measurements, with various degrees of validity of this approximation. Since realistic wind-forcing conditions baffle effective analysis, theoretical attempts to study the development of the wave field are also often limited to steady wind forcing. A simplification frequently applied in theoretical investigations is the so-called fetch-limited case, in which the wave field is assumed to be excited by a steady wind starting from a certain location-say, the shoreline for wind blowing into the sea. The experiments discussed in Section 3 can be seen as corresponding to the fetch-limited case. It is important to stress that under steady wind forcing, the wave field is statistically stationary; thus, the application of temporal Fourier analysis and the computation of frequency spectra of various parameters characterizing the waves are justified, and are indeed commonly applied in laboratory and field experiments. Waves under steady forcing evolve with fetch, as is obvious from the data presented in Section 3. Those spatial variations along the test section in our experimental facility are clearly visible in a video record presented in the Supplementary Material (Video S1). The essential inhomogeneity and lack of spatial periodicity in all the mean statistical parameters of the wave field under steady wind forcing clearly prevent the application of the spatial Fourier analysis.
Another generally adopted simplified theoretical approach is the so-called duration-limited case, where the temporal evolution of waves under impulsive wind forcing is studied. This approach seems to match the experiments under impulsive wind forcing that are described in Section 4. The video record of the evolution of wind waves under such forcing is also presented in the Supplementary Material (Video S2). However, there is an important difference between the experiments and the assumptions adopted in the duration-limited approach. As made obvious in Figure 23, the waves in the experiments evolve in time as well as in space; thus, the wave field is neither stationary nor homogenous, rendering the application of both spatial and temporal Fourier analysis formally unjustified. For that reason, the dominant instantaneous local 'frequency' in the present experiments was determined using wavelet analysis and averaging over the whole ensemble of numerous independent realizations of the wave field, rather than applying the Fourier transform. However, the theoretical studies of the duration-limited case usually assume that the whole wave field is a spatially homogeneous wave field, and its wave-vector spectrum evolves in time [49].
One can argue that any random signal does not entirely satisfy conditions of periodicity, stationarity, and homogeneity; nevertheless, the frequency and wave-vector spectra of such signals of diverse nature are routinely calculated; see examples in Figures 6, 16 and 18 of this paper. However, the curves and contours plotted in those and similar figures represent estimates of the spectra at a fixed location ( Figures 6 and 18), or over a certain area ( Figure 16). Under approximately steady forcing, the measured parameters represent statistically stationary quantities within reasonable accuracy, often allowing the averaging of results based on numerous records taken at different times. One of many examples of variation of the wave-vector spectra measured by airborne remote sensing radar in the near-shore ocean area and demonstrating wave refraction in the process of shoaling is given in [105]. However, if the rate of change is fast, the assumption of stationarity of short segments may lead to inaccurate conclusions, as discussed in [67] regarding the experimental results on the initial growth of waves under impulsive forcing by Kawai [28]. If the records are sufficiently long, or the area covered is large enough, relative to the appropriate time/length scales, they can be subdivided in domains where stationarity/homogeneity can be assumed. Approximate local homogeneity in a statistically steady wave field allows invoking the kinetic equation [106] to analyze the corresponding frequency and wave-vector spectra and additional statistical properties [107].
The situation in the case of unsteady wind forcing in general, and impulsive forcing in particular, is more complicated than the steady forcing case. It is clearly seen in Figure 23 that with a possible exception of the very short initial stage where the first ripples appear more or less uniformly over the whole length of the test section, the ensemble-averaged wave parameters vary with time as well as with fetch. The broad spectrum of wind waves, wave dispersion, and the downshifting of the spectral peak in the process of wave evolution along the tank preclude using a single group velocity c g for the straightforward application of the Galilean transformation to enable studying temporal variations while eliminating spatial ones. It seems that the only conceivable wave field that can be considered spatially homogeneous while evolving in time is that excited in a closed reservoir by omnidirectional wind or mechanical forcing, so that waves have no preferred propagation direction. Such systems indeed were studied extensively theoretically, numerically, and experimentally in the framework of the statistical theory on weakly nonlinear dispersive waves-the so-called wave, or weak, turbulence; see, e.g., Nazarenko [108,109] and the references therein. However, there is no experimental evidence that the preferred wind direction ceases to affect the propagation direction, even for waves much shorter than the dominant one. The simultaneous measurements of two slope components in [67] indicate that slopes associated with the wind direction prevail even at higher frequencies and for stronger wind. Unlike in hydrodynamic turbulence, it has not been demonstrated so far convincingly that free waves in the energy-containing part of the spectrum can be seen as approximately isotropic in the vast majority of the water-wave systems in nature and even in laboratory experiments. Thus, the duration-limited case is a useful abstract concept; however, the applicability of findings based on wave turbulence for waves generated by wind with a preferred propagation direction needs careful experimental verification.

Wind-Wave Growth: The Approach of Phillips vs. Miles-Inspired Theories
All the groundbreaking theoretical studies [18][19][20][23][24][25][26][27][28] investigated wave evolution in time in the framework of the duration-limited case, i.e., without dependence on fetch. Linear monochromatic water waves and unidirectional mean airflow were usually considered; the outcome of those theories is the exponential growth rate in time of waves of a given length for a prescribed vertical mean air velocity profile. However, the essentially three-dimensional structure of the wind-wave field was emphasized in numerous experimental studies, e.g., [67,96,99]. Contrary to the linear and deterministic approach adopted by Miles and his followers, the theory of Phillips [21] is fundamentally random and nonlinear; it also takes into account the directional spreading of wind and of waves. Instead of considering the growth of a single wave harmonic of small amplitude, the nonlinear theory of Phillips predicts a linear with time increase of wave energy in two distinct phases.
It should be noted that the growth of the mean squared values that characterize the disturbance under random forcing that is linear in time, as observed in the experiments on the temporal growth of wind waves summarized in Section 4, is quite common in diverse physical systems under random forcing, such as Brownian motion, random walks, diffusion processes, etc. Closer to the subject of the present study, it is noteworthy to mention the paper by Milewski et al. [110]. They considered a simple four-wave 'toy' model where dispersive waves can exchange energy via resonant nonlinear interactions; some modes in their model were forced by a white noise, while others may be damped. The motivation for [110] came from the fact that the energy exchange among different harmonics in water surface waves occurs through resonant quartets, such as those investigated in wave turbulence studies. However, they stress that the complexity of the systems investigated in the framework of the kinetic equation makes many of the hypotheses underlying these results necessarily heuristic. It was demonstrated theoretically in [110] that the energy of such nonlinear wave systems with negligible dissipation grows linearly with time in a qualitative agreement with the predictions of [21].
In spite of the enormous effort invested studies of wind-wave excitation over the last century, there is surprisingly little direct reliable experimental evidence regarding the dominance of any specific physical mechanism that actually governs wind-wave growth. The theory of Phillips [21], while extensively cited, was actually disregarded over the years, which was probably due to its complexity; the study of Zavadsky and Shemer [67] is the first known to us serious attempt to examine his approach experimentally. The linear deterministic theory of Miles [20] and its subsequent developments are widely accepted as the major guidelines for describing momentum and energy exchange between the ocean and the atmosphere. However, the experimental evidence is scarce. Direct measurements of the exponential growth in time of specific wave harmonics were probably only performed by Larson and Write [41] and Plant and Wright [42]. Kawai [28] reported on the exponential in time growth of the total wave energy observed in his experiments rather than on the evolution of distinct harmonics, as considered in his theoretical analysis. To those few examples, one can add the studies of Hristov et al. [14,15], who provided circumstantial evidence for the existence of a critical layer in wave-induced air velocity profiles for sufficiently fast waves moving with phase velocities c that were high enough to satisfy the condition for wave age c/u * > 15; there is no evidence of the critical layer for slower and younger waves. The existence of the critical layer by itself still does not prove definitely that the Miles mechanism indeed dominates the wave growth, even at large ocean scales. Perrard et al. [111] recently suggested that the exponential growth of the energy of initial ripples as demonstrated in Figure 27a may be triggered by the nearly isotropic pressure fluctuations generated in the airflow boundary layer. Their analysis can be seen as a support for the Phillips [21] approach.
In the absence of convincing experimental results demonstrating an exponential in time growth of a single wave harmonic, the generally accepted notion that Miles-type instability is the main mechanism responsible for wind waves' growth is based largely on the results obtained in analyses of the spatial wave evolution under steady, or approximately steady, wind forcing. For example, Caulliez et al. [98] reported on measurements of the spatial exponential growth of the total energy of initial ripples; however, the stage of the exponential growth apparently terminated when the characteristic wave height remained very small.
Plant [104] compiled results of numerous field and laboratory experiments, and offered the following relation for the temporal growth rate β defined in Equation (18) of waves propagating in the mean wind direction: where the radian wave frequency ω and the phase velocity c characterize the dominant wave. This expression is in qualitative, albeit not quantitative, agreement with Miles [23].
The data accumulated in [51] were used to determine the variation of various spectral frequency components in the vicinity of the dominant frequency for each prescribed flow rate, which was assumed to be exponential with the fetch growth rates. The analysis of the behavior of separate harmonics demonstrated that the significant growth is mainly restricted to components corresponding to lower frequencies and thus longer waves. To estimate the growth of energetic components along the test section, the power of each harmonic, η 2 i was plotted and fitted to the assumed exponential in space growth as exp(γ i x). Only those frequency components at which clear exponential growth with a constant rate γ was observed along at least three examined fetches were taken into account in the fit. The spatial growth was often followed by saturation; fetches exceeding that at which the inception of saturation was detected were excluded from the fit. The directly measured values of the dimensional spatial growth rates γ obtained by this procedure are plotted in Figure 31a as a function of the inverse wave age u * /c. Note that larger values of u * /c correspond to slower-moving shorter waves. Curves corresponding to the empirical relation shown in Equation (26), as well as the prediction of Miles [23], are also plotted in this panel. To this end, the spatial growth rate was estimated as γ = β/ c gr invoking the group velocity calculated from the empirical dispersion relation shown in Equation (15). The data accumulated in [51] were used to determine the variation of various spectral frequency components in the vicinity of the dominant frequency for each prescribed flow rate, which was assumed to be exponential with the fetch growth rates. The analysis of the behavior of separate harmonics demonstrated that the significant growth is mainly restricted to components corresponding to lower frequencies and thus longer waves. To estimate the growth of energetic components along the test section, the power of each harmonic, 2 ̅̅̅ was plotted and fitted to the assumed exponential in space growth as exp(γix). Only those frequency components at which clear exponential growth with a constant rate γ was observed along at least three examined fetches were taken into account in the fit. The spatial growth was often followed by saturation; fetches exceeding that at which the inception of saturation was detected were excluded from the fit. The directly measured values of the dimensional spatial growth rates γ obtained by this procedure are plotted in Figure 31a as a function of the inverse wave age * ⁄ . Note that larger values of * ⁄ correspond to slower-moving shorter waves. Curves corresponding to the empirical relation shown in Equation (26), as well as the prediction of Miles [23], are also plotted in this panel. To this end, the spatial Figure 31. (a) Symbols-dimensionless spatial exponential growth rates under steady wind forcing estimated from experiments for different frequency spectral harmonics, lines-empirical fit by Plant [104] and Miles [23] theory; (b) Measured spectra at U = 7.5 m/s and three fetches.
The experimental results show that the dimensionless spatial growth rate increases with wind velocities as they correspond to higher friction velocities u * which is in qualitative but not quantitative agreement with the predictions based on Equation (28). However, for a constant airflow velocity rate, higher values of γ/k correspond to lower frequencies, which is in a sharp qualitative contrast with Equation (28). The measured growth rates also differ significantly from the predictions. An insight into the reasons for this disagreement can be obtained by examining the spatial variation of the high-frequency part of the predominantly free-wave part of the spectrum that is plotted for selected fetches in Figure 31b for wind velocity U = 7.5 m/s. The decay with x of wave amplitudes of some frequencies in this figure cannot result from dissipation only; moreover, waves at all fetches are subject to nearly identical wind forcing. Thus, the spectra plotted in Figure 31b suggest that essential coupling exists between nonlinearity, wind input, and wave damping during wind waves' evolution with fetch x. Therefore, estimates of the growth or decay rates of individual frequency harmonics that disregard nonlinear interactions as attempted in Figure 31a are probably futile.
The apparent disagreement in Figure 31a between the directly measured growth rates and the estimates provided by the empirical fit in Equation (28), which is based on the Miles approach, is not surprising. Plant and Wright [42] stated that a nonlinear energy transfer from shorter to longer waves is essential, even in the initial stage of wind-wave evolution. Zakharov et al. [49] argued that nonlinearity in fact dominates wave field evolution, while this process is less sensitive to details of energy input and damping. The inapplicability of Miles' critical layer mechanism for short waves where the critical layer may be within the viscous sublayer has been discussed in [36,49]. The direct relation between the spatial, γ, and the temporal, β, exponential growth rates is only applicable for spatially homogeneous waves, and does not hold when spatial wave growth is considered. The exponential in space growth of the total wave energy, rather than that of a single frequency harmonic, as sometimes assumed, does not follow from the linear instability theory. In fact, the well-known relation shown in Equation (6) predicts that the total wave energy m 0 defined by Equation (5) increases nearly linearly with the fetch. In that sense, the spatial wave growth under steady wind forcing that has been experimentally observed in numerous studies at least partially corroborates the Phillips theory for wave field evolution under impulsive forcing.

A Challenge for Wind-Wave Modeling
The existing large-scale wind-wave models such as WAVEWATCH, WAM, and similar are being continuously upgraded over decades, and provide an ever-improving forecast of waves in the ocean. Those models are based on a linear water wave theory that indeed performs sufficiently well; wind input and wave energy dissipation are usually introduced using empirical expressions. The nonlinear effects are incorporated in more recent versions of those models; however, for the relatively weak nonlinearity that characterizes ocean waves, the importance of nonlinearity emerges only at large distances or long times. The relatively coarse resolution of those models does not permit an assessment of the relative contribution of nonlinear terms, as well as an accurate description of the wind input and dissipation adopted in the models.
The detailed experimental data presented in Sections 3 and 4 provide a basis for the quantitative comparison of various approaches to wind-wave modeling with measurements. Beyond the accuracy and detail of data that cannot be attained in field measurements, the wind waves studied are characterized by strong nonlinearity; see, e.g., Figure 22b and additional evidence in [51,66,67,75]. Thus, the appropriately scaled evolution of the wave field due to nonlinear interactions occurs notably faster in those laboratory experiments than in the ocean [112][113][114], facilitating the assessment of the contribution of nonlinearity.
However, the appropriate modeling of developing wind waves as those studied here is far from being straightforward. The random wave field is usually characterized in the Fourier space. The frequency spectra, sometimes complemented by directional parametrization, are usually applied in analyzing the experimental results, while wave-vector spectra are most often considered in theoretical studies. The linear dispersion relation shown in Equation (11) is employed to relate these two formulations. The evolution of wind waves is habitually described theoretically using the kinetic wave equation [4,106] that describes the slow temporal variation of the ensemble-averaged wave-vector spectra of a random wave field due to nonlinear interactions among resonant wave-vector quartets, with additional terms often included to account for wind input and wave dissipation. Annenkov and Shrira [115] performed a comparison of the temporal evolution of a wind-wave field computed using different versions of the kinetic equation with what they called the direct numerical simulation (DNS) of a random wind-wave field based on Monte Carlo simulation using the Zakharov [112,114] equation. In view of the significant differences found between the results obtained by the application of different models, they called for a revision of the fundamentals of the kinetic equation. However, all the models employed in [115] employed wave-vector spectra, and thus were in principle incapable of describing the spatial variation of the wave field excited by wind with a preferred direction. A revision of the kinetic equation that is more profound than that discussed in [115] is required to reproduce the spatial variation of wind waves as observed in the experiments.
Various additional advanced numerical methods that do not rely on the kinetic equation were applied to simulate the evolution of wind waves. For example, direct numerical simulations of the Navier-Stokes equation were employed for that purpose in [116][117][118][119]. Spectral methods invoking periodic boundary conditions were used in all those studies, therefore effectively preventing the spatial development of the wave field as observed in experiments. In fact, there is no variation in time without variation in space due to the nonlinear interactions for any water-wave system with a final spectral width and a preferred propagation direction, such as swell. The main wave characteristics of such an evolving nonlinear dispersive system at the exit from the computational domain are different from those at the inlet; the wave field is necessarily inhomogeneous and non-periodic in space.
For steady forcing with a preferred wind direction, the situation is somewhat simpler, since although the wind-wave field is still inhomogeneous, it is stationary, and thus allows the application of temporal Fourier analysis, yielding frequency spectra. This situation is analogous to the study of deterministic wavemaker-generated nonlinear waves in laboratory flumes. The appropriate nonlinear model for a broad-banded wave spectrum is the well-known Zakharov equation [112][113][114] that describes the evolution in time of the wave-vector amplitudes of a spatially homogenous nonlinear wave field, and is closely related to the kinetic equation. As is the case with the kinetic equation, the solutions of the Zakharov equation cannot be directly verified by experiments that study spatial evolution under steady or periodic mechanical forcing. Therefore, the spatial version of the Zakharov equation [120][121][122] was suggested, in which the evolution in space of the frequency spectrum in a periodic in time wave field is considered. The results of those model computations were found to be in good agreement with numerous experiments carried out in various facilities. It seems that considering frequency rather than wave-vector harmonics is advantageous for modeling wind waves. Although Fourier-based models are limited to steady forcing, they have some promise of success in modeling wind-wave evolution that is consistent with the physics of the process.

Conclusions
The generation of water waves by wind is the most common phenomenon, yet due to its complexity, it is not yet well understood. Numerous physical models aimed at describing the various mechanisms governing the generation of waves by wind have been offered over the years, and different models that attempt to predict the evolution of the wind-wave field were developed. The hypotheses and assumptions underlying those theories are necessarily heuristic and often not sufficiently well substantiated; due to considerable experimental difficulties, they largely have not been verified in quantitative experiments.
The present paper summarizes and discusses the results of a systematic extensive experimental project carried out in our laboratory in recent years in which waves excited by wind in a modestly sized facility were investigated under controlled conditions. Exploiting the advantages of the relatively short time and length scales of the waves in the tank, as well as of a fully controlled and computer-automated experimental procedure, a detailed set of previously unavailable experimental results has been accumulated.
The general case in which waves are affected by wind with velocity magnitude and direction that may vary in time and in space, as usually occurs in nature, is far too complicated to enable extracting new understandings from the experimental data. Two relatively simple cases were therefore considered: (i) waves generated by steady wind forcing, and (ii) waves that appear after an effectively impulsive wind is applied over an initially calm and smooth water surface. These two cases may be related to the fetch-limited and the duration-limited approximations, respectively, that are routinely considered in theoretical studies.
Wave generation by a steady wind has been studied first. The turbulent airflow in the test section above the water surface was well documented. The friction velocities that characterize the interaction of wind and waves at the air-water interface were estimated independently from the logarithmic velocity profiles and from the vertical distribution of the Reynolds stresses. Good agreement between results obtained by the two methods was attained. Measurements of the airflow and the wave field were performed at numerous fetches x and for a range of constant wind velocities U. Numerous statistical parameters that characterize the wave field were determined, including the spatial variation of the characteristic wave amplitude and the frequency spectra. The measured dependencies of the spectral peak frequency and the wave energy on the dimensionless fetch are very close to those obtained elsewhere, with the wave energy growing nearly linearly with fetch. Both the vertical and the horizontal wave asymmetry was estimated; higher statistical moments that characterize the deviation of the wave field from Gaussianity were also presented and supported by calculation of the probability exceedance distributions. Particular attention has been given to the investigation of the three-dimensional aspects of the wave field using a two-component laser slope gauge and stereo video imaging. The essentially random nature of wind waves manifests itself in the correlation length and time-scales that were estimated from the experiments and proved to be short compared to the dominant wind-wave length. The strong three-dimensionality of the wind-wave field and the virtual absence of correlation attained already at short distances indicate that the evolving wave field along the whole test section is effectively insensitive to the existence of sidewalls in the tank. While the characteristic wavelengths and frequencies vary with fetch, the wave steepness η 2 x + η 2 y 1/2 remains nearly constant. Since the wave steepness is the measure of wave nonlinearity, its relatively high values of~0.2 obtained in the present experiments indicate that nonlinearity is significant. The frequency spectra of the temporal variation of both slope components η x and η y were computed and compared with the spectra of the surface elevation η. The appropriately normalized spectra of η(t) were fitted to the standard JONSWAP shape; the fitting coefficients were close to the usually adopted values. The exponential decay of the spectral power with frequency was demonstrated. The frequency spectra of η x and η y become similar at frequencies well above the dominant frequency, thus suggesting that wind waves at high frequencies do not have a preferred propagation direction. The difference of the spectral slope exponents between the surface elevation spectra and those calculated from the LSG data varies with the wind velocity and to a lesser extent with the fetch, and reflects the gradual transition from gravity-capillary to purely capillary waves.
Then, the temporal and spatial evolution of wind waves emerging on a water surface that was initially at rest due to effectively impulsive wind forcing was studied. It should be stressed that temporal rather than spatial variation was considered in the theoretical studies of wind-wave growth surveyed above. Waves generated by impulsive wind forcing in our facility do not differ significantly from those that appear at the initially calm water surface upon the application of wind gusts in larger facilities and in nature. Measurements were carried out for a range of wind velocities and at multiple fetches. The fully automatic experimental procedure controlled by a single computer and the short time scales characterizing the evolution of very young wind waves in a moderate-sized experimental facility enable performing measurements during numerous independent realizations of a random wave field under identical wind-forcing conditions. This experimental approach allowed the computation of ensemble-averaged time-resolved representative wave parameters, as a function of time elapsed since the activation of wind forcing; the major stages in the wind-wave evolution were identified. Since Fourier analysis is inapplicable for the emerging nonstationary and inhomogeneous wave field, the local instantaneous dominant 'frequencies' were determined by ensemble-averaging Morlet wavelet maps.
The first observable surface waves are mostly unidirectional, but very fast become essentially three-dimensional. These initial gravity-capillary ripples grow exponentially with time in agreement with the linear viscous instability mechanism analyzed by Kawai [28]. Their instantaneous frequencies decrease fast; the ripples rapidly attain wave steepness close to that observed under steady wind forcing, being only weakly dependent on fetch and wind velocity. The exponential growth lasts less than 300 ms and then terminates, while the ripples still have submillimeter amplitudes. The growth of wind waves following the appearance of the initial wavelets is notably different. The dominant wave frequency continues to decrease with time at all fetches and wind-forcing conditions, while the wave energy growth with time becomes close to linear, as predicted by Phillips [21] during the initial stage of the waves' development under random forcing. Then, the rate of change of wave parameters increases sharply again, during what can be identified as the principal stage of wind-wave development in the notation of Phillips, until the quasi-steady state is attained. Some of the assumptions regarding the structure of the turbulent airflow above the water surface and on the resonant waves' propagation velocities invoked by Phillips [21] were modified, taking advantage of the information accumulated in our previous studies. These modifications led to the expected wave growth during the principal evolution stage as <η 2 >~U 4 t, which is in agreement with the experimental results.
These results demonstrate that at the scale of the present experiments, the deterministic linear theory that originated from Miles [20] can only explain the growth of initial ripples, while their subsequent evolution is better described in the framework of the random, three-dimensional, and nonlinear approach suggested by Phillips [21]. Therefore, the Phillips theory may serve as an appropriate starting point to understand the physical mechanisms that govern the evolution of wind waves from initial ripples to longer short gravity waves.
Finally, various existing methodologies to model the development of the wind waves generated by wind with a preferred direction were discussed in view of the insights gained in the course of this project. The statistical parameters o5f developing waves always vary in space; thus, the wave field is inhomogeneous, and its evolution cannot be described in terms of wave-vector spectra. The evolution in time of waves that propagate in a certain general direction cannot be considered separately and decoupled from the spatial variation. There is no temporal evolution without spatial evolution in wind waves. Therefore, computations of the temporal evolution of wind waves applying numerical methods that invoke periodic boundary conditions also necessarily yield results that are limited to omnidirectional wind forcing. No such model is capable of predicting the evolution patterns observed in the experiments described in this paper. The assumptions of homogeneity and the spatial periodicity of a wave field cannot probably be fully avoided in any model analysis. However, one has to keep in mind that the resulting model predictions are based on an inevitably flawed approach; thus, the limits of applicability of those simulations to wind waves excited by a realistic directional wind need verification.