Wave asymmetry and particle orbits in irregular wave models

Abstract Water particle orbits are key elements in the Lagrange wave formulation. The stochastic Miche implementation of the Lagrange model is a linear Gaussian two-dimensional or three-dimensional space–time model which exhibits typical nonlinear wave characteristics when transformed to Eulerian coordinates. This paper investigates the statistical relation between the degree of front–back asymmetry of individual waves and the orbit orientation for the particle located at the wave maximum at the point of observation. It is shown that, in the Lagrangian model with statistical front–back symmetry, for individual waves there is a clear connexion between the degree of individual wave asymmetry and the orientation of the randomly deformed elliptic orbit: a steep front correlates with upward tilt, a steep back is correlated with a downward tilt. This holds both for waves observed in time and in space, and the dependence is stronger for large amplitude waves than for smaller ones. The dependence is strongly dependent on the depth and on the significant steepness and spectral width. Inclusion of the average Stokes drift has a moderate effect on the dependence. For models with forced front–back asymmetry there is both a systematic dependence and a statistical correlation between asymmetry and tilt; for large amplitude waves the systematic relation dominates. The conclusions are based on Fourier simulations of Gauss–Lagrange waves of first and second order with a wind–sea Pierson–Moskowitz spectrum and a narrow swell JONSWAP spectrum.

Water particle orbits are key elements in the Lagrange wave formulation. The stochastic Miche implementation of the Lagrange model is a linear Gaussian two-dimensional or three-dimensional space-time model which exhibits typical nonlinear wave characteristics when transformed to Eulerian coordinates. This paper investigates the statistical relation between the degree of front-back asymmetry of individual waves and the orbit orientation for the particle located at the wave maximum at the point of observation. It is shown that, in the Lagrangian model with statistical front-back symmetry, for individual waves there is a clear connexion between the degree of individual wave asymmetry and the orientation of the randomly deformed elliptic orbit: a steep front correlates with upward tilt, a steep back is correlated with a downward tilt. This holds both for waves observed in time and in space, and the dependence is stronger for large amplitude waves than for smaller ones. The dependence is strongly dependent on the depth and on the significant steepness and spectral width. Inclusion of the average Stokes drift has a moderate effect on the dependence. For models with forced front-back asymmetry there is both a systematic dependence and a statistical correlation between asymmetry and tilt; for large amplitude waves the systematic relation dominates. The conclusions are based on orbits in ocean waves are scarce owing to the difficulty and cost of obtaining accurate data, even with modern techniques (Romero & Melville 2010;Fedele et al. 2013;Benetazzo et al. 2018). Most orbital studies are therefore performed in combination with experiments in a physical or numerical wave tank.
The following studies are of high relevance for the present work. Chen, Hsu & Chen (2010) and Chen et al. (2012) analyse theoretical Lagrangian models for monochromatic waves over uniform and sloping bottoms and compare particle orbits with experimental results. Grue & Jensen (2012) and Grue & Kolaas (2017) report particle velocities for both laboratory waves and directional ocean waves based on Romero & Melville (2010) and reconstruct the orbits with respect to different phases and vertical positions. Starting with the 'Choppy wave' model Nouguier, Guérin & Chapron (2009), Nouguier, Chapron & Guérin (2015 and Guérin et al. (2019) improve on the Lagrange model for semi-regular and irregular waves with consequences for particle orbit studies. Creamer et al. (1989) give an alternative to the second-order technique in the cited works. Experimental studies by van den Bremer et al. (2019) on monochromatic waves focus on mass transport and particle orbits at different depths near wave crests and wave troughs, respectively.
We shall investigate systematic and statistical relations between the degree of front-back asymmetry of individual waves in space and time and the geometry of particle orbits. From a statistical viewpoint one may call the former property an apparent (Euler) characteristic, and the latter an intrinsic (Lagrangian) characteristic. The dominant way to empirically describe a sea state is by its energy spectrum, a typical Eulerian concept. The stochastic connection with a Lagrangian description goes via hydrodynamic equations of varying complexity. In this work we will use Gaussian stochastic models for the Lagrangian vertical and horizontal movements and numerically transform them to an Euler description of the water surface. The link between the two descriptions is described, e.g. in Lindgren & Lindgren (2011) and Guérin et al. (2019). We study unidirectional irregular waves developing in time and space along a straight line over constant depth h.
To finish this introduction we make a comment on the photographs taken by Wallet & Ruellan (1950) on particle trajectories in plane water waves reflected by a partially absorbing barrier, reproduced in Hutter, Wang & Chubarenko (2011, figure 7.8). For a monochromatic wave without reflection the trajectories are very regular in intermediate water, nearly circles at the surface and elliptical near the bottom. With increasing moderate reflection the reflected wave will interfere with the forward travelling wave with increased irregularity and orientation as the result. In the present work the irregularity will come from a continuum of frequencies resulting in a statistical distribution for the trajectories including a statistical relation with the wave asymmetry.

The Gauss-Lagrange models
The theme in this paper is the relation between individual wave front-back asymmetry and the shape and orientation of water particle orbits. The idea is that the relation is an inherent relation, appearing in all theoretical models for irregular waves. The model chosen is the first-order Gauss-Lagrange model, and its second-order extension.
2.1. The free and forced first-order models The first-order Gauss-Lagrange two-dimensional wave model is defined by two correlated Gaussian space-time random fields, describing the vertical and horizontal particle positions as functions of time t and original horizontal position u along the horizontal axis. Together, these fields define the orbital movements of the water particles. We consider here only particles at the free surface but the model extends to general depth. Following Pierson (1961), the energy spectrum S(ω) of the vertical field is called the orbital spectrum. It is not identical to the Euler spectrum, obtained from observations of the ocean surface, but the difference is of no relevance in the present work.
Representing the continuous energy spectrum by a spectrum over discrete frequencies ω j = j ω and wavenumbers k j , we write the models for particles on the free surface as vertical position : Here, h j = cosh(k j h)/ sinh(k j h) is the depth dependent amplitude factor, with dispersion relation ω j = gk j tanh(k j h). The amplitudes A j can be taken as deterministic, A j = 2S(ω j ) ω, but to obtain correct variability in the simulations we must let them be random, A j = a 2 j + b 2 j , with independent normal variables a j , b j , with mean zero and equal variance such that A 2 j has expected value 2S(ω j ) ω. The relative phases φ j shall be independent and uniformly distributed in [0, 2π], and the phase shift between vertical and horizontal movement is π/2 as in (2.2), independent of frequency.
As was shown in Lindgren &Åberg (2009) andLindgren (2010) the Lagrange model with constant phase shift (2.2) will produce waves which are statistically symmetric with respect to wave front and wave back slopes. To get statistically asymmetric waves one can modify (2.2) and introduce a frequency dependent phase shift ρ j , horizontal position : (2.3) A heuristic choice, inspired by the filtering equation and not underpinned by hydrodynamic theory, leads to where the a-value determines the degree of statistical wave asymmetry. The model can be used to reproduce statistical asymmetry that can be observed in real ocean waves. From the pair W(u, t), X(u, t) one can implicitly define a Lagrange wave L(x, t) by that is, its surface height at time t at location X(u, t) is W (u, t). We call the Lagrange model defined by (2.2) and constant phase shift the free Lagrange model and the one defined by (2.3) the model with forced asymmetry, which will be dealt with in § 4.3.

2.2.
The second-order model In the second-order Lagrange model the independent frequency components in (2.1) and (2.2) are combined over pairs of frequencies to give second-order terms, added to the first-order terms. This model exists in many different versions and we present here only the simplest form for the vertical process as an example Here, Z j = a j + ib j are the complex representations of the first-order amplitudes and H + jk , H − jk are complex depth and frequency dependent transfer factors; for examples see Marthinsen & Winterstein (1992), Fouques, Krogstad & Myhaug (2006), Nouguier et al. (2015) and Herbers & Janssen (2016).
When the sum of first-and second-order terms are combined in (2.6) we get the second-order Lagrange model. That model has recently been studied by Guérin et al. (2019), whose theoretical results have bearing also on the particle orbits. To complement their theoretical studies we will simulate the relation between wave asymmetry and orbits in the second-order model and then use the formulation by Biésel (1952) to generate the second-order corrections W 2 (u, t), X 2 (u, t). The routine is implemented in the MATLAB package WAFO and is described in Lindgren & Prevosto (2017). The routine also gives the stochastic local Stokes drift X 2,drift (u, t). In the simulations we will use the average Stokes drift from the integral formula, kωS(ω) dω.

Wave asymmetry and particle orbits
The themes in this paper are the systematic and statistical relations between the degree of front-back asymmetry of irregular random waves, and the shape and orientation of the associated particle orbits. For that reason we make the following definitions of the asymmetry of individual waves, different from those used by, among others, Stansell, Wolfram & Zachary (2003).
For time waves we record the free surface height variation at a fixed location, and identify the sequence of local maxima and minima in the recording. By an individual time wave we mean the part of the signal that lies between two successive local minima. Its front and back slopes are defined as the rate of increase from the initial minimum to the maximum and the rate of decrease from the maximum to the terminal minimum, respectively. As measure of asymmetry we take the logarithm of the ratio between the two. A positive skewness measure means a fast increase and slow decrease.
For space waves we record the surface variation along a straight line at a fixed time. Identifying local minima and maxima we compute front and back slopes, with the front facing the direction towards which waves are travelling. A positive logarithmic skewness measure means a steep front.
For orbit shape and orientation we will study the particle that is located at the maximum of each individual wave in time or space, at the time of observation. This means that we follow that particle in a time interval before and after the time point when the maximum was registered. The orbit is approximated by an ellipse and the direction of its major axis identified. We choose clockwise orientation for both time and space wave orbits, which is consistent with the assumed main wave direction. Figure 1 illustrates the skewness-orbit relation. Figure 1 Figure 2 illustrates the wave shape/orbit relation for space wave in the same sea state as in figure 1. We will compare time and space waves systematically in § 4.1.

Method
It will be shown by simulation that there exist strong relations between wave skewness and orbital orientation, both for models with statistically symmetric first-order waves, (2.2), and for models with different types of wave asymmetry, including second-order waves with and without a Stokes drift, and waves with forced statistical asymmetry, (2.3).
For a pre-defined orbital spectrum S(ω), ω > 0, we shall generate independent realizations of the vertical W(u, t) and horizontal X(u, t) components in the Gauss-Lagrange model along a line 0 ≤ u ≤ U for a time span 0 ≤ t ≤ T, in a discrete space and time grid, u n = n u, t m = m t. In the simulation experiments we will use the routine spec2ldat3DM (Lindgren & Prevosto 2017) in the module lagrange in the MATLAB toolbox WAFO (WAFO-group 2017).
If X(u, t) is a strictly increasing function of u then (2.6) has a unique solution L( fixed we obtain a Lagrange space wave L(u, t 0 ) and keeping x = x 0 fixed we obtain a time wave L(x 0 , t). As noted already by Pierson (1961) it can happen that X(u, t) is not strictly increasing but for the sea states and water depths we consider in the paper these events are very rare.
For each wave realization we identify all local maxima and minima and the corresponding front and back steepnesses. Simultaneously computing the particle orbits from the model, we can find a pair 'before-after asymmetry -orbit shape and orientation' for each local maximum. The procedure is repeated with independent realizations, and it will give an estimate of the statistical relation between wave asymmetry and orbit orientation. Details on the sampling procedure to obtain individual waves and on the different measures will be given in § § 3.1 and 3.2.
As an alternative to the min/max definition one can take the trough/crest definition, where the wave is defined by the maximum and minimum between successive zero crossings. This will lead to more chaotic orbits with extra twists and less characteristic skewness measures. The link between tilt and skewness will be weaker but still significant.

Space waves
For space waves we fix a time point t 0 and observe maxima and minima, M k , m k , in L(x, t 0 ) as a function of x. We denote their locations by x + k , x − k with the wave front facing the positive x-direction and define the front and back steepnesses by (3.1) The front-back asymmetry is recorded in logarithmic scale as A k = log(−s + k /s − k ) and waves with positive A k have a steep front and less steep back; see figure 3 for an illustration. To find the trajectory of the particle that is located at the maximum M k we observe that, when X(u, t 0 ) is strictly increasing, which is the case in all the examples, for each k there is a unique originũ k such that X(ũ k , t 0 ) = x + k , and hence W(ũ k , t 0 ) = M k . The centred trajectory of the top particle therefore has coordinates, as a function of τ , The time interval for the orbit of the top particle will be defined by its movement before and after the observation time point t 0 , as follows. For the maximum at x + k , originating from positionũ k we define d ± k so that W(ũ k , t 0 + τ ) has local minima for t 0 − d − k and t 0 + d + k , these minima being the closest to t 0 . This definition works well if the maximum M k is sufficiently high. It is less suitable for very low maxima for which other methods are needed for an analysis of the orbit-asymmetry relation.

Time waves
For time waves we fix an observation point x 0 ∈ [x min , x max ], and observe the surface variation L(x 0 , t) at that point, identifying all local maxima and minima in that time series, and the corresponding front and back steepnesses.
In each realization, we denote the height of local maxima in L(x 0 , t) by M k , k = 1, . . . , n and their time of occurrence by t + k . We denote the height of the preceding and following local minima by m k , m k+1 and let their time of occurrence be t − k , t − k+1 . The wave front and wave back steepness measures of the time wave are then defined as To find the trajectory of the particle at the maximum M k we identify, when X(u, t) is strictly increasing, the unique originũ k such that X(ũ k , t + k ) = x 0 and hence W(ũ k , t + k ) = L(x 0 , t + k ). The particle trajectory near t + k therefore has coordinates, as a function of τ , In the simulation, the orbit O k (τ ) is sampled at τ j = j t, to give O k (τ j ).
3.3. Orbit characteristics As can be expected, trajectories will be quite irregular for both space and time waves and never exactly elliptical, as is the idealized shape, and they have no unique orientation. The deviation from elliptic shape can be large for small waves, for example with a negative maximum or with a small front or back amplitude, and the meaning of orbit orientation can be questionable. Therefore, the present study will be restricted to waves whose maximum exceeds H s /8, i.e. half the standard deviation of the Gaussian waves, and significant amplitudes, min(M k − m k , M k − m k+1 ) > H s /4. We call these major waves.
As a proxy for orbit orientation, we will fit an ellipse to the sampled trajectory in a finite time interval, −d − k < τ < d + k , around the maximum, and use its orientation as a measure of orbit tilt. A natural choice for the time interval, used in the simulations, is the interval between the two vertical minima closest to the maximum, i.e. d − but that is by no means the only possibility. As can be seen, e.g. in figure 10, our choice results in approximately full round orbits for first-order waves, with more irregularity for second order, and with natural gaps for waves with Stokes drift.
For the fit of the proxy ellipse we will use the MATLAB routine fit_ellipse which makes a least-squares fit of data to a conic representation of an ellipse. This algorithm may fail to find an approximating ellipse, in particular for minor waves. For the major waves analysed in this study it may also in a few cases (less than 0.5 %) report exactly horizontal orbits. These cases are removed in the figures. From the fitted ellipse we extract the tilt of the major axis, counted anticlockwise relative to the horizontal axis, i.e. a positive value means upward tilt. The orbit in figure 3 is an example of the approximation procedure. REMARK 3.1. We have used the orientation of the fitted ellipse as a measure of orbit tilt. It is obviously a statistical measure whose distribution can only be obtained by a simulation experiment. As a possible alternative one could take the orientation of the velocity vector at the point/time of the maximum, a quantity whose distribution can be derived explicitly in the Gauss-Lagrange model. Preliminary studies have indicated that there is indeed a correlation between these two measures of orbit tilt. How it relates to the wave steepness remains to be studied.

Sampling the space and time waves and the influence of Stokes drift
The wave definition (3.1) is unnecessarily complicated for simulation of surface wave asymmetry. When X(u, t 0 ) is strictly increasing there is no need to compute the Lagrange wave L(x, t 0 ) to find its maxima and minima. It suffices to identify maxima and minima in the Gaussian wave W(u, t 0 ) and transform the steepness by the process X(u, t 0 ). Further, an average Stokes drift only affects the particle orbits, not the wave shapes. The time wave simulation sampling of individual waves is more complicated, in particular if Stokes drift is involved. The sampling technique without Stokes drift amounts to finding the local maxima of the Gaussian field W(u, t) along a curveũ(t) defined by X(ũ(t), t) = x 0 , varying Thus, inclusion of Stokes drift means that we change the sampling procedure when we collect time wave maxima. No change is made in the generating Gaussian fields. Clearly, the two sampling procedures pick different maxima with different wave shape, but they come from the same generating fields which are stationary in time and space. The gentle slope of the Stokes drifted curve may have a very small systematic effect on wave asymmetry distribution, and a special simulation study has shown that there is no visible difference in the joint skewness/orbit tilt distribution.
It should be mentioned that the discussion in Guérin et al. (2019) about the drifting particle origin deals with a more complicated model, but the principle is the same.

Examples
We illustrate the technique on two spectra with very different wave characteristics described in table 1, a Pierson-Moskowitz (PM) spectrum and a JONSWAP PM: Pierson-Moskowitz wind-sea spectrum with nominal significant wave height H s = 4.5 m, spectrum peak period T p = 10 s, mean zero crossing period T 02 = 7.5 s and mean zero crossing wavelength 71 m. The Benjamin-Feir index BFI at infinite depth is 0.23 and it is 0.29 at 20 m depth. J20: an extreme and narrow JONSWAP spectrum with H s = 4.5 m and T p = 10 s and a peak enhancement factor γ = 20, mean zero crossing period 9.2 s and mean zero crossing wavelength 112 m. BFI = 0.87 at infinite depth, 0.92 and 1.12 at depths 40 m and 20 m. TABLE 1. Specification of the example spectra. Note that the simulation algorithm truncates the spectra to ω max = 2 which leads to a minor increase in the BFI value. The mean crossing data refer to the truncated spectra.
spectrum (J20). Both have H s = 4.5 m and peak period T p = 10 s. The PM spectrum is representative of a wind-sea spectrum close to saturation. Its Benjamin-Feir index BFI = 0.23 at infinite depth. The J20 spectrum is a narrow-band one, representative of a strong steep swell situation. Its BFI is 0.87 at infinite depth. In the simulation some truncation of the spectra will take place with a minor increase in the BFI values.
We will now illustrate the difference between time and space waves, the effect of varying depth, the influence of wave steepness and spectral peakedness, forced asymmetry and the inclusion of second-order terms and average Stokes drift.

Time waves versus space waves
The waves and orbits shown in figure 1(a-d) are based on a 200 s simulation sampled at 10 Hz. The orbital spectrum is the PM spectrum at infinite depth. From the simulation the major time waves with front and back amplitudes and maximum height exceeding H s /4 and H s /8, respectively, were extracted. That resulted in 9 time waves with positive orbit tilt (figure 1a,b) and 8 with negative tilt (figure 1c,d). The orbits shown in the right column extend over the time period from the preceding minimum to the following minimum with blue circle and red cross indicating the beginning and the end of the orbit. The black dots indicate the maximum of the orbit of the top particle.
The p.d.f. of tilt and skewness in figure 1(e) is based on almost 100 000 pairs of waves and orbits in simulated Lagrange time waves. The corresponding space wave figure 2, is based on 25 000 wave/orbit pairs from the same model and it shows a similar, but not identical, dependence pattern at infinite depth.
The conclusion is clear. There is a statistically significant dependence between time wave front-back asymmetry and the orientation of the involved orbits, as measured by the orientation of the 'best fit' ellipses. Since the orbits are not closed, the term 'ellipse' shall not be taken literally but as a proxy for the irregular shape and the orientation of its greater axis as statistical measure of orbit orientation. Also for space waves there is a clear dependence for this unidirectional wave model. For a model with directional orbital spectrum the relations might be different.

Influence of water depth, wave steepness and spectral peakedness
A question of great interest is how the tilt/skewness dependence depends on water depth. We will study the dependence for time and for space waves for the PM spectrum at depths h = 40, 30, 20 m and for the JONSWAP spectrum at depths h = Inf, 40, 20 m.   We compare the depth dependence for the two spectra, in figures 5 and 6 for PM waves and in figures 7 and 8 for J20 waves. As shown in the leftmost diagrams in all four figures, at large depths both spectra give a bi-modal tilt/skewness p.d.f. for the major waves we study. Both spectra show a gradual change over intermediate depths towards a unimodal diamond shaped distribution. The distribution change for the J20 space waves is more distinct than that for the PM time waves.
One should have in mind that the distributions are based on major waves with large amplitudes and maximum height. The orbit examples in figure 1 show why the proxy ellipse orientation is so strongly coupled to the asymmetry. It is a geometric fact, not related to the stochastic wave model that we use. The change from a bi-modal distribution at deep water to a unimodal one at shallow water is probably an effect of the increase in orbit eccentricity with decreasing depth, which forces the two modes to collapse.

Lagrange models with forced asymmetry
The tilt/skewness dependence illustrated so far is a statistical measure caused by the geometric constraints of the waves. The Lagrange model (2.3) with non-negative a-parameter and frequency dependent phase shift imposes a systematic shift in the particle orbit distribution.
The a-values used in the examples in Lindgren & Åberg (2009), Lindgren (2010 and Lindgren & Lindgren (2011) represent large asymmetries in wave shape and a strong average tilt. We choose a small value, a = 0.1, at depth h = 30 m, and illustrate the tilt/skewness distribution for time waves with a PM spectrum in figure 9. As seen, the orientations of the proxy ellipses have been shifted towards positive, upward tilt, an effect also seen in the joint tilt/skewness density.

The second-order Lagrange model and the influence of Stokes drift
We now show how the second-order correction terms affect the Lagrange trajectories from (2.6) and then also include the Stokes drift from the second-order theory.
For the PM model the average Stokes drift at depth h = 20 m is kωS(ω) dω = 0.11 m s −1 . Figure 10 shows typical tilted space waves and corresponding orbits for first-and second-order models, without and with an added Stokes drift. There is no big difference between first and second order but the effect of the Stokes drift is clearly seen in the rightmost panels, (e,f,k,l). The gap between orbit start and end is compatible with the size of the Stokes drift and the peak period T p = 10 s. Figure 11 illustrates how the second-order orbits are translated and opened by the Stokes drift while they retain a somewhat elliptic shape. One can also note that, due to a smaller frequency bandwidth, the J20 orbits are more regular than the PM orbits. Figure 12 shows the tilt/skewness p.d.f. for the three PM space models. It confirms the similarity between first-and second-order models. It also illustrates the effect on the joint distribution of the introduction of the Stokes drift in space, which makes the orbits look more horizontally elongated, increasing the chance that the proxy ellipse is identified as nearly horizontal with near zero tilt. The remaining dependence between tilt and wave skewness in the right plot is therefore the best illustration of the stochastic covariation between these two characteristics of the wave dynamics. For the time waves in figure 13 the Stokes drift has no effect on the joint distribution, in agreement with the argument in § 3.4. Finally we study the J20 swell spectrum and show, in figures 14 and 15, the results for time and space waves at infinite depth where the average Stokes drift is 5.5 cm s −1 . The first-order model delivers a bi-modal tilt/skewness distribution both in space and in time while the second-order models both give different more uniform tilt distributions. The small Stokes drift has almost no effect for the time waves in figure 15. The strong tilt/skewness correlation is clear here as for the PM model.

Summary and conclusions
In a simulation study we have illustrated and discussed the statistical relation between the front-back asymmetry of irregular ocean waves in the Gauss-Lagrange wave model. We used two different orbital spectra, one Pierson-Moskowitz (PM) spectrum for wind driven sea states and one narrow-band JONSWAP (J20) for a strong steep swell situation. Both spectra have H s = 4.5 m and peak period T p = 10 s. Their Benjamin-Feir values on infinite depth are 0.23 and 0.87, respectively.
We simulated the vertical and horizontal Gaussian components of the Lagrange model and generated observable wave trains in time and in space in which we identified individual waves by the local min-max-min definition. The log ratio between the slopes before and after the local maxima was used to give a measure of wave skewness, front-back for space waves and 'before-after' for time waves. Only major waves with sufficient maximum height and amplitude were considered. For each wave we identified the water particle located at the maximum at the time/location of observation and reconstructed its orbit as a function of time between the preceding and following local minima. To cope with the great variability in the shape of the orbits a standardized procedure was adopted to give a measure of orbit orientation. We did choose a least-squares fit of an ellipse to the orbits and used the inclination of its main axis as a measure of orbit tilt. This measure is sensitive to the geometry of the waves and it also reflects the increasing orbit eccentricity with decreasing depth.
The covariation between the degree of asymmetry of individual waves, as measured by the front-back skewness, and the estimated tilt of the associated particle orbits have been summarized in their bivariate p.d.f., based on between 24 000 and 150 000 pairs of waves/orbits. The purpose of this paper has been to introduce a statistical technique to find intrinsic patterns in irregular ocean waves and test it in one common wave model, thereby complementing the model based and experimental measures mentioned in the introduction. The general conclusion of the study is that, in the Lagrange model, there is a systematic positive relation between the wave skewness and the main orientation of the orbits as estimated in the proposed method. Many factors determine the details of the joint distribution and we discuss the main findings.

The difference between time waves and space waves
Orbit tilt and wave skewness are clearly positively dependent, in both time and space, figures 1 and 2. The dependence can be seen as a summary of the wave and orbit geometry and is not coupled to the specific mathematical model that describes them. One difference between the two types of observation is how one can handle the presence or absence of Stokes drift. We have included examples of its effect on the space waves from the PM spectrum on shallow water, figure 12 and on space waves from the J20 spectrum on infinite depth in figure 14. The technique used in the study allows inclusion of a Stokes drift also in time waves, with some modification of sampling mechanism. The modification is different from that suggested in Guérin et al. (2019). We have not investigated its effect on the distribution but it can be expected to be small.

Depth dependence
For both spectral types there is a typical shift in dependence with decreasing water depth, from a bi-modal p.d.f. in deep water to a unimodal p.d.f. in shallow water. Both spectra show the accompanying increase in BFI values, but further studies are needed to quantify its influence. The change is present in time waves as well as in space waves (figures 6-8) in accordance with increasing eccentricity of the orbits, cf. Wallet & Ruellan (1950). The orientation of more eccentric ellipses will also be more exactly identified by the fitting algorithm than the more irregular orbits at great depth.

Second-order effect and the effect of Stokes drift
Inclusion of a second-order term in the model has a clear but small effect on the shape and extent of tilt/skewness dependence, as shown for PM waves in figures 12(a,b) and 13(a,b). The average Stokes drift, when added to the local model, makes a difference (figure 12c). The effect is explained in figure 10(e, f,k,l) and it is due to the elongation of orbits, that concentrates the orientation towards horizontal. For the JONSWAP spectrum the higher-order models make a bigger change, as illustrated in figures 14 and 15 except for the Stokes drift in time waves.

The dependence of the model
The observed statistical dependence is a consequence of the general wave geometry and it should be rather robust for change of spectrum. Experiments, not illustrated in the paper, with other type of spectra and characteristics confirm the general pattern for the Lagrange model. It might also be insensitive to the specific stochastic wave elevation model used, for example a second-order Eulerian model or that of Creamer et al. (1989). A systematic study of alternative models and of the random variation of the orbits around the typical shape is outside the scope of the present paper. An interesting theme for a future study would be a comparison with a similar experiment in a physical wave tank.

Declaration of interests
The authors report no conflict of interest.