Next Article in Journal
Comments on “Wind Gust Detection and Impact Prediction for Wind Turbines”
Next Article in Special Issue
Multi-Frequency, Multi-Sonar Mapping of Shallow Habitats—Efficacy and Management Implications in the National Marine Park of Zakynthos, Greece
Previous Article in Journal
A National Assessment of Wetland Status and Trends for Canada’s Forested Ecosystems Using 33 Years of Earth Observation Satellite Data
Previous Article in Special Issue
Side-Scan Sonar Image Mosaic Using Couple Feature Points with Constraint of Track Line Positions
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Modeling Microseism Generation by Inhomogeneous Ocean Surface Waves in Hurricane Bonnie Using the Non-Linear Wave Equation

by
Joshua D. Wilson
General Dynamics—Applied Physical Sciences, Groton, CT 06340, USA
Work Conducted at the Massachusetts Institute of Technology, Department of Mechanical Engineering, Cambridge, MA 02139, USA.
Remote Sens. 2018, 10(10), 1624; https://doi.org/10.3390/rs10101624
Submission received: 2 July 2018 / Revised: 30 September 2018 / Accepted: 4 October 2018 / Published: 12 October 2018
(This article belongs to the Special Issue Advances in Undersea Remote Sensing)

Abstract

:
It has been shown that hurricanes generate seismic noise, called microseisms, through the creation and non-linear interaction of ocean surface waves. Here we model microseisms generated by the spatially inhomogeneous waves of a hurricane using the non-linear wave equation where a second-order acoustic field is created by first-order ocean surface wave motion. We treat range-dependent waveguide environments to account for microseisms that propagate from the deep ocean to a receiver on land. We compare estimates based on the ocean surface wave field measured in hurricane Bonnie in 1998 with seismic measurements made roughly 1000 km away in Florida.

1. Introduction

Hurricanes generate seismic noise, commonly referred to as microseisms, in the 0.1 to 0.6 Hz frequency range. On a few occasions measurements of these microseisms have been used to track hurricanes [1,2,3,4,5,6,7,8]. More generally these microseisms are a significant cause of noise in seismic measurements [7,9,10] and raise the detection threshold for monitoring earthquakes [11].
Physically the process for generating microseisms can be broken down into four steps: (1) The spatially and temporally varying hurricane winds create complicated wave patterns on the ocean surface; (2) In locations where these surface waves are in opposition, their non-linear interactions generate acoustic waves; (3) These acoustic waves propagate through the water and couple into the ocean bottom; (4) The resulting seismic waves propagate along the sea bottom where they may be detected at great distances by seismic sensors. In this paper we develop an approach for estimating the microseisms generated by a hurricane and their seismo-acoustic propagation.
Our modeling approach is based on the non-linear wave equation where a second-order acoustic field is generated by a source distribution that depends on the first-order ocean surface wave motion. The acoustic (or seismic) field at a receiver can then be expressed as the integral over the source distribution multiplied by the waveguide Green function. This approach is ideal for hurricane-generated microseisms since it can be applied to spatially inhomogeneous surface wave fields. Also, this approach may be used in range-dependent waveguide environments as is the case when a hurricane at sea generates microseisms that propagate up the continental margin to a receiver on land.
As an example, we use the ocean surface directional wave spectrum from hurricane Bonnie [12,13] to estimate the microseismic source levels. We than calculate the Green’s function to describe the seismo-acoustic propagation from Bonnie’s location in the North Atlantic to a seismometer in Florida. Our calculated levels are then compared against actual measurements.
Our approach builds on a long history of microseism modeling efforts. Initially, microseism generation was modeled for idealized situations [14] where the ocean surface wave field was taken to be spatially homogeneous and the ocean waveguide to be range-independent. Previous microseismic models, however, cannot be applied in typical hurricane scenarios because the surface wave fields are inhomogeneous in that the wave height spectra in different parts of the hurricane can vary both in magnitude and directionality. In some cases this is because spatial homogeneity is assumed over infinite [15,16,17,18,19,20] or very large [14,21,22] surface wave areas. One model that does account for finite microseismic source generation regions is limited by the assumption that the receiver is at the center of the source area [23,24]. Also, these previous models are not applicable to the range-dependent environments typical in hurricane measurements where microseisms generated by a hurricane over the deep ocean are measured by sensor on land. Our non-linear wave equation method is shown to agree with earlier models [18,20] if we make the same simplifying assumptions that the source generation region is spatially homogeneous and that the waveguide can be modeled as an infinite half-space.
While the treatment of inhomogeneous surface wave fields and range-dependent waveguides presented here had not been previously considered, the medium non-linear wave equation used here has a parallel derivation to the wave equations or hydrodynamic equations used in earlier works. Longuet-Higgins [14], Hasselmann [17] and Brekhovskikh [15] base their derivations on perturbing Bernoulli’s equation while Hughes [18] and Kibblewhite and Wu [19] begin by separately perturbing the equations of momentum, continuity and state. We show that the non-linear wave equation used here is equivalent to perturbation expansions [18] for the physical parameters typical in hurricane microseisms. Lloyd [20] and Cato [23] base their derivations on Lighthill’s equation; however, Lloyd [20] shows that both perturbation and Lighthill approaches yield the same end result. We also note that the non-linear wave equation can be derived as a second-order approximation of Lighthill’s equation [25].
These second-order non-linear theories should not be confused with the linear theories proposed by Banerji [26] and Bowen et al. [27]. These linear theories claim that the first-order motion of a surface gravity wave creates a first order pressure fluctuation on the sea floor regardless of how deep the ocean is. This is contrary to classic surface wave theory which shows that first-order wave motion decays exponentially with depth such that, in deep water, the first-order pressure fluctuation on the sea floor goes to zero [14,28,29,30].
Our derivation is in agreement with earlier works [14,31], which show that microseisms are generated by the non-linear interaction of ocean surface waves with roughly the same wavelength but opposing propagation directions. Measurements [13,32] and models [12] of surface directional wave height spectra in hurricane Bonnie show complex patterns with the opposing surface waves necessary to generate microseisms. Based on the wave height spectra in hurricane Bonnie, we calculate the microseismic source levels generated by the non-linear interaction of the ocean surface waves.
It has been shown both theoretically [10,14] and experimentally [4,11,33,34,35,36,37] that microseisms propagate as Rayleigh waves along the sea floor. Here we develop a normal-mode model for the propagation of Rayleigh waves along the seafloor assuming a water depth that changes with range. The Green’s function from this model is used to estimate the propagation of the microseism from deep water to a seismometer on land.
Given the wave height spectra in hurricane Bonnie, the calculation of the microseismic source level generated by the non-linear wave interaction, and the propagation of those waves along the sea floor, we estimate the noise levels at a seismometer in Florida. These estimates are then compared against measured data.
It should be noted that there are other sources of microseisms besides hurricanes. Microseisms have been associated with other non-hurricane storms, coastal regions where complex wave patterns are prevalent, and the typical ambient ocean wave spectra [11,38,39,40,41,42,43,44]. These other sources of microseisms are not considered here.

2. Materials and Methods

2.1. Ocean Surface Gravity Waves

Hurricanes are characterized by high winds that can vary quickly with position, both in direction and speed, as shown in Figure 1A. In addition a hurricane may move at speeds up to 15 m/s [45] so that the winds at any location can also change with time. These spatially and temporally varying winds generate complex ocean surface wave directional spectra Figure 2 with wave heights that can exceed 10 m [12,13,32].
For surface gravity waves that are homogeneous over an infinite ocean surface area we can express the complex surface wave height as the linear superposition of plane waves [46] where κ is the sea surface wavenumber and σ is the corresponding frequency where in deep water σ 2 = g κ . Throughout this paper the surface wavenumber κ is expressed either in cartesian coordinates as ( κ x , κ y ) or in polar coordinates as ( κ , θ ) . Also, a cartesian coordinate system is used for position where r = ( ρ , z ) = ( x , y , z ) where z is defined downward from the sea surface. The complex surface wave height can then be written as
ζ ^ ( ρ , t ) = - A ( κ ) e - i σ t + i κ · ρ d 2 κ
where A ( κ ) is the surface wave height spectrum.
Since the complex wave height ζ ^ ( ρ , t ) is produced by the contributions of many independent random physical processes, we can assume by Central Limit Theorem that the statistical distribution of ζ ^ ( ρ , t ) is Gaussian [47]. In addition we may define the reference position from which wave height is measured such that ζ ^ ( ρ , t ) has zero mean. Given this assumption, in Appendix A.1 we show that A ( κ ) must also be a zero-mean complex Gaussian random variable [46,47].
Since the spectrum A ( κ ) is zero-mean and homogeneous over an infinite ocean surface area, the second moment of the spectrum can be written as (see Appendix A.1) [16,47,48]
< A ( κ ) A * ( κ ) > = S ( κ ) δ ( κ - κ ) ,
where < > denotes the expectation. This expression shows that, for stationary processes, the wavenumber components κ and κ decorrelate. Also, since the spectrum A ( κ ) is Gaussian, the fourth moment can then be written as [16]
< A ( κ 1 ) A ( κ 2 ) A * ( κ 3 ) A * ( κ 4 ) > = δ ( κ 1 - κ 3 ) S ( κ 1 ) δ ( κ 2 - κ 4 ) S ( κ 2 ) + δ ( κ 1 - κ 4 ) S ( κ 1 ) δ ( κ 2 - κ 3 ) S ( κ 2 ) .
In a hurricane the surface waves are not homogeneous over an infinite area, but rather the surface wave spectrum A ( κ ) changes gradually with position in the storm. To characterize the spatially varying surface wave field in the hurricane, we divide the sea surface into a grid made up of finite regions m with surface areas S m over which the wave spectrum A m ( κ ) can be taken to be homogeneous. The complex wave height in a particular finite region m then can be expressed as
ζ ^ m ( ρ , t ) = - A m ( κ ) e - i σ t + i κ · ρ d 2 κ
In Appendix A.2 we show that the finite size of the regions m introduces a ‘windowing’ effect when calculating the moments of the spectrum A m ( κ ) . However, we also show that this windowing effect can be neglected if the dimensions L x and L y of the region are much greater than the wavelength λ g of the surface wave. For the hurricane waves of interest in this paper, the surface wavelengths λ g range from roughly 100 to 300 m requiring that the dimensions L x and L y be in the order of 1 km or greater. Provided this condition is met we can make the approximation.
< A m ( κ ) A n * ( κ ) > S m n ( κ ) δ ( κ - κ )
where S m n ( κ ) is the cross-spectral density of the surface wave fields at locations m and n.
As before we can assume by Central Limit Theorem that ζ ^ m ( ρ , t ) is a zero-mean complex Gaussian random variable so that A m ( κ ) is also a zero-mean complex Gaussian random variable. Because of this we can write the fourth moment as
< A m ( κ 1 ) A m ( κ 2 ) A n * ( κ 3 ) A n * ( κ 4 ) > δ ( κ 1 - κ 3 ) S m n ( κ 1 ) δ ( κ 2 - κ 4 ) S m n ( κ 2 ) + δ ( κ 1 - κ 4 ) S m n ( κ 1 ) δ ( κ 2 - κ 3 ) S m n ( κ 2 ) .
This expression for the fourth moment of the wave height spectra will be used later in Section 2.3.
From Equation (4) we can write the real part of the surface wave height as
ζ m ( ρ , t ) = 1 2 - ( A m ( κ ) e - i σ t + i κ · ρ + A m * ( κ ) e i σ t - i κ · ρ ) d 2 κ
We can also express the real surface wave particle velocity v as a linear superposition of plane waves such that
v m ( r , t ) = - - σ 2 ( κ x κ i ^ x + κ y κ i ^ y ) ( A m ( κ ) e - i σ t + i κ · ρ + A m * ( κ ) e i σ t - i κ · ρ ) e - κ z - i σ 2 i ^ z ( A m ( κ ) e - i σ t + i κ · ρ - A m * ( κ ) e i σ t - i κ · ρ ) e - κ z d 2 κ .
where i ^ x , i ^ y and i ^ z represent unit vectors in the x, y and z directions respectively.
Experiments by Forristall et al. [49] show that the linear wave theory of Equations (1), (4), (7) and (8) is adequate to describe surface waves even in high sea states with significant wave breaking as in a hurricane. The use of linear surface wave theory can also be justified since the typical surface wave heights in a hurricane are an order of magnitude less than the wavelength [12,13,32,49].
Wright et al. [13] and Walsh et al. [32] measured the spatial variation of the directional wave spectra A m ( κ ) in hurricane Bonnie using an aircraft-mounted Scanning Radar Altimeter (SRA). These measurements, however, are limited to the locations and times of the aircraft flights and do not give a complete picture of the wave field in the hurricane.
To quantify the surface wave field at any location and time in hurricane Bonnie, Moon et al. [12] used the Wave Watch III (WW3) ocean surface wave model [50]. These model results, shown in Figure 2, are in close agreement with measurements from [13,32] and with data from buoys and oceanographic stations. In this paper we use the surface directional wave spectra A m ( κ ) in hurricane Bonnie calculated using WW3 [12].

2.2. Non-Linear Wave Equation

The non-linear wave equation [51,52] describes the second-order acoustic wave field generated by first order fluid motion. For the long surface waves in a hurricane viscosity may be neglected so that the non-linear wave equation can be written as [53]
2 p 2 ( r , t ) - 1 c 2 2 p 2 ( r , t ) t 2 = Γ ρ c 4 2 p 1 2 ( r , t ) t 2 - ρ 2 c 2 2 t 2 ( v 1 ( r , t ) · v 1 ( r , t ) ) - ρ 2 2 ( v 1 ( r , t ) · v 1 ( r , t ) ) ,
where the first order velocity v 1 and pressure p 1 terms on the right-hand side generate a second-order acoustic pressure field p 2 on the left-hand side, c is the sound speed, ρ is density and Γ is the coefficient of nonlinearity [51,53].
In this work we will consider the case where the first order velocity v 1 , from Equation (8), is due to surface wave motion. The first order pressure p 1 can be found using conservation of momentum where - p 1 = ρ v 1 / t . Using the relationship between first order velocity and pressure we can compare the relative magnitudes of the source terms in Equation (9). Given typical parameters κ = 0 . 02 to 0 . 7 rad/m, σ = 0 . 4 to 0 . 8 rad/s, c 1500 m/s, ρ 1000 kg/m 3 , and Γ - 2 . 6 [54] we find that the ρ 0 2 2 ( v 1 ( r , t ) · v 1 ( r , t ) ) term of Equation (9) exceeds the other source terms by several orders of magnitude. Dropping the lesser terms yields
2 p 2 ( r , t ) - 1 c 2 2 p 2 ( r , t ) t 2 = - ρ 2 2 ( v 1 ( r , t ) · v 1 ( r , t ) ) = q ( r , t ) .
This equation is equivalent to the perturbation expression used in the previous microseism derivation of Hughes [18]. This can be seen by taking Equation (10) of [18] and assuming irrotational flow ( ( v 1 ( r , t ) · v 1 ( r , t ) ) = 2 ( v 1 ( r , t ) · ) v 1 ( r , t ) ) and continuity ( · ( ρ 0 v 1 ( r , t ) ) = 0 ). Taking the Fourier Transform of Equation (10) we find the frequency domain Helmholtz equation for the second-order field
2 P 2 ( r , ω ) + k 2 P 2 ( r , ω ) = q ( r , ω ) .
We will begin by solving for the pressure field generated by a finite volume V m with a source distribution q m which depends on the fluid velocity v 1 , m from Equation (8). The volume V m extends vertically from the ocean surface at z = 0 to a depth well below the surface wave region ( z > λ g = 2 π / κ ). The volume extends over the horizontal area S m with dimensions L x and L y centered at ρ m = ( x m , y m ) as defined in Section 2.1. We can solve the Helmholtz (Equation (11)), using the frequency domain Green function g ( r r , r s , ω ) such that,
P 2 , m ( r r , ω ) = V m g ( r r , r s , ω ) q m ( r s , ω ) d 3 r s ,
where r s = ( x s , y s , z s ) and r r = ( x r , y r , z r ) are the source and receiver locations respectively.
It should be noted that, although the ocean surface is perturbed by the presence of surface waves, we choose to define the upper boundary of our volume by the quiescent surface z = 0 following the example of Lloyd [20]. He argues that the Green function g “has the form g = g 0 + g 1 + g 2 + where g 0 is the Green function for the quiescent ocean and where g j , j 1 , is of degree j in a wave height perturbation parameter < h 2 > 1 / 2 . The first-order term g 1 represents Bragg reflected by the moving spatially periodic components”. He shows that the terms g j for j 1 may be omitted in the excitation of p 2 . He goes on to apply the Green function for a quiescent ocean half-space and, even with this approximation, matches the result of earlier work by Hughes [18]. In Appendix C we will demonstrate agreement between our end result and that of both [20] and [18].
Equation (12) can be rewritten as
P 2 , m ( r r , w ) = 0 d z s - w ( ρ s - ρ m ) g ( r r , r s , w ) q m ( r s , w ) d 2 ρ s ,
where the window function w ( ρ s - ρ m ) is unity for - L x / 2 < x s - x m < L x / 2 and - L y / 2 < y s - y m < L y / 2 and zero elsewhere. The double integral over ρ s represents the field at a receiver r r due to a distribution of sources q over a finite area defined by the window function w. If the receiver r r is in the far field of the finite area ( | r r - r m | > L x 2 / λ a and | r r - r m | > L y 2 / λ a where λ a is the acoustic wavelength and where r m = ( ρ m , z s ) = ( x m , y m , z s ) ), we can make the plane wave approximation,
g ( r r , r s , w ) g ( r r , r m , w ) e - i k r · ( ρ s - ρ m )
where the acoustic wavenumber is defined in cartesian coordinates as k = ( k r , k z ) = ( k x , k y , k z ) where ω / c = k x 2 + k y 2 + k z 2 = k . In free space this would simplify to
g ( r r , r s , w ) = 1 4 π | ρ r - ρ s | e i k | ρ r - ρ s | 1 4 π | ρ r - ρ m | e i k | ρ r - ρ m | e - i k r · ( ρ s - ρ m )
however, we will use the more general form given in Equation (14). Inserting Equation (14) into Equation (13) yields
P 2 , m ( r r , w ) 0 d z s g ( r r , r m , w ) - w ( ρ s - ρ m ) e - i k r · ( ρ s - ρ m ) q m ( r s , w ) d 2 ρ s ,
Again the double integral over ρ s analogous to a horizontal planar array of acoustic sources with broadside directed down toward the ocean bottom. In the following analogy we can define this double integral as the array output [55]
B m ( k , z s ) = - w ( ρ s - ρ m ) e - i k r · ( ρ s - ρ m ) q m ( r s , w ) d 2 ρ s ,
The array output in Equation (17) can be evaluated for any horizontal wavenumber vector k r ; however, only those wavenumber vectors where | k r | ω / c lead to a propagating acoustic field. This can be illustrated if we consider the microseismic field generated by a horizontal plane of sources near the sea surface as described in Equation (17). In previous works [14] it has been shown that microseisms generated at the ocean surface are first transmitted downward to the sea floor where they then propagate along the bottom as Rayleigh waves. This vertical propagation from the sea surface to the sea floor can be modeled to first order using the free-space Green function
g ( w ) = 1 4 π d e i k z d
where d represents the water column depth. If | k r | exceeds ω / c then k z becomes imaginary leading to an exponential decay with depth d. For example if | k r | = 2 ω / c then k z = i 3 ω / c so that, for ω = 1 . 3 rad/s, c = 1500 m/s and an ocean depth d = 1000 m, the field is attenuated by e - k z z r 0 . 2 or 13 dB. Greater depths will lead to even more attenuation. Because of this evanescent attenuation we will only consider the acoustic field where | k r | ω / c .
Assuming the acoustic pressure field given in Equation (16) is temporally stationary, the power spectral density may be written as [56]
S P m , n ( r r , ω ) δ ( ω - ω ) = < P 2 , m ( r r , w ) P 2 , n * ( r r , w ) > 0 0 d z s d z s g ( r r , r m , w ) g * ( r r , r n , w ) × < B m ( k , z s , m ) B n * ( k , z s , n ) > ,
where < > represents the expectation and r n = ( ρ n , z s ) = ( x n , y n , z s ) . In the next section we derive expressions for the second moment < B m ( k , z s , m ) B n * ( k , z s , n ) > of the array output and the power spectral density S P m , n due to ocean surface waves.
The power spectral density S P m , n in Equation (19) represents the field at a receiver from a single pair of finite volumes V m and V n . In the next section we determine the microseismic field generated by a hurricane where the surface wave field is inhomogeneous and extends over a large region hundreds of kilometers across. To do this we divide the hurricane region into finite volumes V m , as shown in Figure 3, and then sum their contributions to find the total power spectral density of the received field written as
S P ( r r , ω ) = m , n S P m , n ( r r , ω ) .

2.3. Power Spectral Density Due to Ocean Surface Gravity Waves

Taking the ocean surface wave velocity as the first order field we can determine the source term from Equation (11) where q ( r s , ω ) is the Fourier transform
q m ( r s , ω ) = 1 2 π - q m ( r s , t ) e i w t d t
of
q m ( r s , t ) = - ρ 2 2 ( v 1 , m ( r s , t ) · v 1 , m ( r s , t ) ) .
Given Equation (8) for the first order velocity we find
2 ( v 1 , m ( r s , t ) · v 1 , m ( r s , t ) ) = - - σ 1 σ 2 4 × [ A m ( κ 1 ) A m ( κ 2 ) e - i σ + t e i κ + · ρ s F 1 ( κ 1 , κ 2 ) + A m * ( κ 1 ) A m * ( κ 2 ) e i σ + t e - i κ + · ρ s F 1 ( κ 1 , κ 2 ) + A m ( κ 1 ) A m * ( κ 2 ) e - i σ - t e i κ - · ρ s F 2 ( κ 1 , κ 2 ) + A m * ( κ 1 ) A m ( κ 2 ) e i σ - t e - i κ - · ρ s F 2 ( κ 1 , κ 2 ) ] × e - ( κ 1 + κ 2 ) z s d 2 κ 1 d 2 κ 2 ,
where σ + = ( σ 1 + σ 2 ) , σ - = ( σ 1 - σ 2 ) , κ + = ( κ 1 + κ 2 ) , and κ - = ( κ 1 - κ 2 ) . Note that the Laplacian operator ( 2 ) in Equation (23) leads to the direction cosine functions
F 1 ( κ 1 , κ 2 ) = ( κ + ) 2 - ( κ x + ) 2 - ( κ y + ) 2 - 1 + κ 1 x κ 2 x κ 1 κ 2 + κ 1 y κ 2 y κ 1 κ 2 = 2 κ 1 κ 2 [ 1 - cos ( θ 1 - θ 2 ) ] [ - 1 + cos ( θ 1 - θ 2 ) ] ,
and
F 2 ( κ 1 , κ 2 ) = ( κ - ) 2 - ( κ x - ) 2 - ( κ y - ) 2 1 + κ 1 x κ 2 x κ 1 κ 2 + κ 1 y κ 2 y κ 1 κ 2 = 2 κ 1 κ 2 [ - 1 + cos ( θ 1 - θ 2 ) ] [ 1 + cos ( θ 1 - θ 2 ) ] ,
where κ x + = ( κ 1 x + κ 2 x ) , κ y + = ( κ 1 y + κ 2 y ) , κ x - = ( κ 1 x - κ 2 x ) , and κ y - = ( κ 1 y - κ 2 y ) .
Substituting Equation (23) into the definition of q m ( r s , t ) from Equation (22) and taking the Fourier transform yields,
q m ( r s , ω ) = - ρ 8 - - σ 1 σ 2 × [ A m ( κ 1 ) A m ( κ 2 ) F 1 ( κ 1 , κ 2 ) δ ( ω - σ + ) e i κ + · ρ s + A m * ( κ 1 ) A m * ( κ 2 ) F 1 ( κ 1 , κ 2 ) δ ( ω + σ + ) e - i κ + · ρ s + A m ( κ 1 ) A m * ( κ 2 ) F 2 ( κ 1 , κ 2 ) δ ( ω - σ - ) e i κ - · ρ s + A m * ( κ 1 ) A m ( κ 2 ) F 2 ( κ 1 , κ 2 ) δ ( ω + σ 1 - ) e - i κ - · ρ s ] × e - ( κ 1 + κ 2 ) z s d 2 κ 1 d 2 κ 2 .
By inserting Equation (26) into the definition of the array output B m ( k , z s ) from Equation (17) we obtain
B m ( k , z s ) = - ρ 8 - - σ 1 σ 2 × [ A m ( κ 1 ) A m ( κ 2 ) F 1 ( κ 1 , κ 2 ) δ ( ω - σ + ) - - w ( ρ s - ρ m ) e i κ + · ρ s e - i k r · ( ρ s - ρ m ) d 2 ρ s + A m * ( κ 1 ) A m * ( κ 2 ) F 1 ( κ 1 , κ 2 ) δ ( ω + ( σ 1 + σ 2 ) ) - - w ( ρ s - ρ m ) e - i ( κ 1 + κ 2 ) · ρ s e - i k H · ( ρ s - ρ m ) d 2 ρ s + A m * ( κ 1 ) A m ( κ 2 ) F 2 ( κ 1 , κ 2 ) δ ( ω - ( σ 1 - σ 2 ) ) - - w ( ρ s - ρ m ) e i ( κ 1 - κ 2 ) · ρ s e - i k H · ( ρ s - ρ m ) d 2 ρ s + A m * ( κ 1 ) A m ( κ 2 ) F 2 ( κ 1 , κ 2 ) δ ( ω + ( σ 1 - σ 2 ) ) - - w ( ρ s - ρ m ) e - i ( κ 1 - κ 2 ) · ρ s e - i k H · ( ρ s - ρ m ) d 2 ρ s ] × e - ( κ 1 + κ 2 ) z s d 2 κ 1 d 2 κ 2 .
Note that the wave height spectra A m and A n in Equation (27) can be brought out of the spatial integral since they are constant over the ocean surface area S m defined by the window function w. Integrating over ρ s then leads to
B m ( k , z s ) = - ρ 8 - - σ 1 σ 2 × [ A m ( κ 1 ) A m ( κ 2 ) F 1 ( κ 1 , κ 2 ) δ ( ω - σ + ) W ( k r ) δ ( k r - κ + ) + A m * ( κ 1 ) A m * ( κ 2 ) F 1 ( κ 1 , κ 2 ) δ ( ω + σ + ) W ( k r ) δ ( k r + κ + ) + A m ( κ 1 ) A m * ( κ 2 ) F 2 ( κ 1 , κ 2 ) δ ( ω - σ - ) W ( k r ) δ ( k r - κ - ) + A m * ( κ 1 ) A m ( κ 2 ) F 2 ( κ 1 , κ 2 ) δ ( ω + σ - ) W ( k r ) δ ( k r + κ - ) ] × e - ( κ 1 + κ 2 ) z s d 2 κ 1 d 2 κ 2 ,
where ∗ represents the two-dimensional convolution written as,
F ( k r ) G ( k r ) = - F ( k r ) G ( k r - k r ) d 2 k r
and where
W ( k r ) = 2 sin ( k x L x 2 ) k x 2 sin ( k y L y 2 ) k y .
Note that, when integrating over ρ s , the ρ m term in the window function w ( ρ s - ρ m ) of Equation (27) introduces a phase factor e - i k r · ρ m . This cancels the other phase term e i k r · ρ m in Equation (27). From this we see that the source array output given in Equation (28) has no phase dependence on location ρ m .
The delta functions in Equation (28) show that the horizontal component k r of the acoustic wavenumber vector is generated by and equal to either the sum κ + or difference κ - of the surface wavenumber. In Section 2.2 we discussed how the non-propagating field where | k r | > ω c is not significant and may be ignored. Therefore, when integrating over κ 2 we need only consider a range of integration Ω corresponding to propagating acoustic waves. For the sum terms in Equation (28) Ω corresponds to the integration region where | κ + | ω c and for the difference terms Ω corresponds to | κ - | ω c .
Here we will define the source area dimensions L x and L y to be ‘acoustically small’ ( L x , L y < < λ a ). The acoustic wavelengths of the microseismic field in water range from roughly 6 to 11 km requiring that the dimensions L x and L y be on the order of 1 km or less. Note that this requirement on L x and L y , combined with the requirement that L x , L y > > λ g from Section 2.1, means that λ g < < L x , L y < < λ a . If L x , L y < < λ a , the function W ( k r ) L x L y in Equation (30) within the range of acoustically propagating wavenumbers k r ω / c . With this approximation Equation (28) can be written as
B m ( k , z s ) = - ρ L x L y 8 - Ω σ 1 σ 2 × [ A m ( κ 1 ) A m ( κ 2 ) F 1 ( κ 1 , κ 2 ) δ ( ω - σ + ) + A m * ( κ 1 ) A m * ( κ 2 ) F 1 ( κ 1 , κ 2 ) δ ( ω + σ + ) + A m ( κ 1 ) A m * ( κ 2 ) F 2 ( κ 1 , κ 2 ) δ ( ω - σ - ) + A m * ( κ 1 ) A m ( κ 2 ) F 2 ( κ 1 , κ 2 ) δ ( ω + σ - ) ] × e - ( κ 1 + κ 2 ) z s d 2 κ 1 d 2 κ 2 .
The power spectral density of the pressure in Equation (19) contains the second moment or variance of the array output B m ( k , z s ) which from Equation (31) can be expressed as,
< B m ( k , z s ) B n * ( k , z s ) > = ρ 2 L x 2 L y 2 64 - Ω - Ω σ 1 σ 2 σ 3 σ 4 × < [ A m ( κ 1 ) A m ( κ 2 ) F 1 ( κ 1 , κ 2 ) δ ( ω - ( σ 1 + σ 2 ) ) + A m * ( κ 1 ) A m * ( κ 2 ) F 1 ( κ 1 , κ 2 ) δ ( ω + ( σ 1 + σ 2 ) ) + A m ( κ 1 ) A m * ( κ 2 ) F 2 ( κ 1 , κ 2 ) δ ( ω - ( σ 1 - σ 2 ) ) + A m * ( κ 1 ) A m ( κ 2 ) F 2 ( κ 1 , κ 2 ) δ ( ω + ( σ 1 - σ 2 ) ) ] × [ A n * ( κ 3 ) A n * ( κ 4 ) F 1 ( κ 3 , κ 4 ) δ ( ω - ( σ 3 + σ 4 ) ) + A n ( κ 3 ) A n ( κ 4 ) F 1 ( κ 3 , κ 4 ) δ ( ω + ( σ 3 + σ 4 ) ) + A n * ( κ 3 ) A n ( κ 4 ) F 2 ( κ 3 , κ 4 ) δ ( ω - ( σ 3 - σ 4 ) ) + A n ( κ 3 ) A n * ( κ 4 ) F 2 ( κ 3 , κ 4 ) δ ( ω + ( σ 3 - σ 4 ) ) ] > × e - ( κ 1 + κ 2 ) z s e - ( κ 3 + κ 4 ) z s d 2 κ 1 d 2 κ 2 d 2 κ 3 d 2 κ 4 .
Since the wave height spectra A m and A n are uncorrelated zero-mean Gaussian random variables as discussed in Section 2.1, we can substitute Equation (6) into Equation (32) and integrate over κ 3 and κ 4 to obtain
< B m ( k , z s ) B n * ( k , z s ) > = ρ 2 L x 2 L y 2 32 - Ω σ 1 2 σ 2 2 S m n ( κ 1 ) S m n ( κ 2 ) × [ F 1 2 ( κ 1 , κ 2 ) δ ( ω - σ + ) δ ( ω - σ 1 + ) + F 1 2 ( κ 1 , κ 2 ) δ ( ω + σ + ) δ ( ω + σ + ) + F 2 2 ( κ 1 , κ 2 ) δ ( ω - σ - ) δ ( ω - σ - ) + F 2 2 ( κ 1 , κ 2 ) δ ( ω - σ - ) δ ( ω + σ - ) + F 2 2 ( κ 1 , κ 2 ) δ ( ω + σ - ) δ ( ω - σ - ) + F 2 2 ( κ 1 , κ 2 ) δ ( ω + σ - ) δ ( ω + σ - ) ] × e - ( κ 1 + κ 2 ) z s e - ( κ 1 + κ 2 ) z s d 2 κ 1 d 2 κ 2 .
Then by integrating κ 2 in Equation (33) over the integration region Ω , corresponding to | κ + | ω c and | κ - | ω c , we find
< B m ( k , z s ) B n * ( k , z s ) > = 2 ρ 2 π L x 2 L y 2 k 2 - σ 1 4 κ 1 4 × [ S m n ( κ 1 ) S m n ( - κ 1 ) δ ( ω - 2 σ 1 ) δ ( ω - 2 σ 1 ) + S m n ( κ 1 ) S m n ( - κ 1 ) δ ( ω + 2 σ 1 ) δ ( ω + 2 σ 1 ) + 4 S m n ( κ 1 ) S m n ( κ 1 ) δ ( ω ) δ ( ω ) ] × e - 2 κ 1 ( z s + z s ) d 2 κ 1 .
Since ω c < < κ this integration leads to the approximations κ 1 - κ 2 in the sum terms and κ 1 κ 2 in the difference terms.
The term 4 S m n ( κ 1 ) S m n ( κ 1 ) δ ( ω ) δ ( ω ) of Equation (34) leads to a constant ( ω = ω = 0 ) value in the power spectral density that is irrelevant. Also, for the remainder of the paper, we suppress the subscripts of κ 1 and σ 1 to simplify notation.
The wave height cross-spectral density S m n ( κ ) can be converted to the frequency-angle ( σ - θ ) domain S m n ( σ , θ ) by evaluating the Jacobian κ σ κ where σ 2 = g κ = g | κ | [57]. This yields the relation S m n ( σ , θ ) = κ σ κ S m n ( κ ) = 2 σ 3 g 2 S m n ( κ ) . Making this substitution in Equation (34) and integrating over σ yields,
< B m ( k , z s ) B n * ( k , z s ) > ρ 2 π L x 2 L y 2 k 2 ω 9 512 g 2 δ ( ω - ω ) × e - ω 2 z s 2 g ( z s + z s ) 0 2 π S m n ω 2 , θ S m n ω 2 , θ + π d θ .
Note that the delta functions in Equation (34) introduce a frequency doubling effect in Equation (35) where the frequency ω of the acoustic wave is twice the frequency σ of the surface gravity wave.
Substituting Equation (35) into Equation (19) and integrating over ω leads to the power spectral density of the pressure field
S P m , n ( r r , ω ) = ρ 2 π L x 2 L y 2 k 2 ω 9 512 g 2 0 0 d z s d z s × g ( r r , r m , w ) g * ( r r , r n , w ) e - ω 2 z s 2 g e - ω 2 z s 2 g × 0 2 π S m n ω 2 , θ S m n ω 2 , θ + π d θ = 0 0 d z s d z s g ( r r , r m , w ) g * ( r r , r n , w ) × S B m n ( ω , z s ) ,
where we define the ‘microseismic source cross-spectral density’ as
S B m n ( ω , z s ) = - < B m ( k , z s ) B n * ( k , z s ) > d ω = ρ 2 π L x 2 L y 2 k 2 ω 9 512 g 2 e - ω 2 z s 2 g e - ω 2 z s 2 g × 0 2 π S m n ω 2 , θ S m n ω 2 , θ + π d θ .
Summing the contribution from all finite volumes V m gives us the total received power spectral density
S P ( r r , ω ) ρ 2 π L x 2 L y 2 k 2 ω 9 512 g 2 m , n 0 0 d z s d z s × g ( r r , r m , w ) g * ( r r , r n , w ) e - ω 2 z s 2 g e - ω 2 z s 2 g × 0 2 π S m n ω 2 , θ S m n ω 2 , θ + π d θ m , n 0 0 d z s d z s g ( r r , r m , w ) g * ( r r , r n , w ) × S B m n ( ω , z s ) .
Equation (38) provides an analytic expression for microseisms generated by inhomogeneous ocean surface wave fields, as opposed to previous formulations where spatial homogeneity is assumed over infinite [15,17,18,19,20] or very large [14,21,22] surface wave areas. Also the Green functions g ( r r , r m , w ) and g * ( r r , r n , w ) in Equation (38) may be calculated for any arbitrary range-dependent or range-independent ocean waveguide using standard propagation models. Later we will use a formulation for range-dependent Rayleigh wave propagation to calculate the Green functions and model microseismic propagation in a typical North Atlantic waveguide environment.
While Equation (38) is applicable to inhomogeneous surface wave fields and arbitrary ocean waveguides, it can also be applied to simpler homogeneous surface wave fields and Green functions. In fact, in Appendix C we find that, given an infinite half-space and assuming the surface wave spectrum is range-independent ( S m n ( ω / 2 , θ ) = S ( ω / 2 , θ ) ), the power spectral density of the pressure field in Equation (38) simplifies to the results derived by Hughes [18] and by Lloyd [20].

2.4. Microseismic Source Levels in Hurricane Bonnie

In this section we show how Equation (37) can be used to find the microseismic source cross-spectral density generated by the surface wave field in a realistic hurricane. In these examples we use the wave spectra from hurricane Bonnie. Through this analysis we demonstrate the relationship between the wind speeds, surface wave spectra, and microseismic source spectra in a hurricane.
Beginning on 22 August 1998, hurricane Bonnie traveled along the east coast of the Bahamas, Florida, and South Carolina [58]. During this time microseisms were recorded by an Incorporated Research Institutions for Seismology (IRIS) seismometer at the Disney Wilderness Preserve in Florida ( 28 . 1 N, 81 . 4 W). Here we model the microseisms generated by hurricane Bonnie that would be measured by this seismometer.
To evaluate the microseismic source cross-spectral density S B m n of Equation (37) requires knowledge of the surface wave height spectrum A m . Moon et al. [12] calculate the surface wave height spectrum for the wind speeds measured in hurricane Bonnie using the Wave Watch III (WW3) wind-wave modeling program [50]. They also show how their results are in close agreement with aircraft-based measurements of surface wave height spectrum by Wright et al. [13] and Walsh et al. [32].
Figure 1A shows the surface wind speed in hurricane Bonnie on 24 August 1998. Illustrated is the typical hurricane structure with the high wind speeds of the ‘eye wall’ surrounding the ‘eye’ at the hurricane’s center. From these wind speeds, WW3 is used to determine the wave height spectra A m ( κ ) as defined in Equation (5) (Figure 2). In Figure 2 we see that hurricanes generate waves that propagate in many different and often opposing directions. For example in Figure 2c we see a peak in the surface wave spectrum (in red) with waves traveling North ( ( κ x , κ y ) ( 0 . 0 , 0 . 04 ) ), West ( ( κ x , κ y ) ( - 0 . 04 , 0 . 0 ) ), and South ( ( κ x , κ y ) ( 0 . 0 , - 0 . 04 ) ).
From the wave height spectra calculated using WW3, the microseismic source cross-spectral density S B m n may be calculated from Equation (37). We define the ‘microseismic source level’ as
L B m n ( ω , z s ) = 10 log ( S B m n ( ω , z s ) L x 2 L y 2 1 m 4 1 P a 2 1 s )
expressed in dB re Pa 2 s / m 4 where the somewhat arbitrary length scales L x and L y have been factored out. Figure 1B shows L B m m ( ω , 0 ) at ω = 1 . 26 rad/s ( f = 0 . 2 Hz) on 24 August 1998, the day where we found the source level to be highest. We see that the peak of L B m m at location (c) is not at the same location as the maximum wind speed (a). This is because, while location (a) has large wave heights (≈37 dB re m 2 / ( rad 2 / m 2 ) ), the waves are all propagating primarily to the South without any opposing waves. At location (c), however, there are opposing waves (propagating both North and South) which generate microseisms as expressed in Equation (37). Also the source level is higher in the low-wind-speed eye (b) than at the maximum wind speed location (a) due to the opposition in the surface waves (Figure 2b).
The arrows in Figure 1 represent the direction hurricane Bonnie was moving and we see that there is a small peak in the source level well ‘behind’ the hurricane at location (d) even though the wind speed there is relatively low. This is again due to the opposing waves at this location (Figure 2d). This illustrates the complex relationship between wind speed and wave spectra, where the wave spectra is a function not only of wind speed but also of hurricane geometry and translation speed [12], and the complex relationship between wave spectra and microseismic source level, where the level depends on the non-linear interaction of opposing waves.

3. Results

In this section we model the propagation of the microseisms generated by hurricane Bonnie through the North Atlantic waveguide. The model results are then compared with seismic data gathered in Florida.
To determine the microseisms received at the sensor in Florida the Green functions g ( r r , r m , w ) and g * ( r r , r n , w ) of Equation (38) are calculated using the adiabatic Rayleigh wave propagation model derived in Appendix B. The bathymetry between the hurricane, located over the roughly 5 km deep Hatteras Abyssal plain on 24 August, and the seismometer, located in Florida is characterized by a gentle upslope. We will model this up-sloping environment with the simplified geometry shown in Figure 4. The sound speeds and densities for this model waveguide are based on typical values for the deep North Atlantic Ocean [59] and bottom [60].
With this environment the Green functions g ( r r , r m , w ) and g * ( r r , r n , w ) are calculated using Equation (A47). The sea-floor depth at the center of hurricane Bonnie is given in Table 1, as well as the range between the seismometer and the hurricane center. The seismometer in Florida is buried in the ground to a depth of 162 m as shown in Figure 3.
Given the source cross-spectral density from Equation (38) and Green functions from Equation (A47), the power spectral density of the field, from Equation (38), received by the seismometer in Florida is calculated. Figure 5 shows the power spectral density of the horizontal velocity (Figure 5a–d) and vertical velocity (Figure 5e–h) of the earth’s crust at the Florida seismometer based on our model (the red dashed curves) and on measured data (in blue).
The dotted blue lines in Figure 5 represent portions of the measured data that are dominated by non-hurricane related ambient noise. The level of this ambient noise is determined by measuring the noise levels during the week before and the week after the hurricane passed. The un-corrupted microseism signals are taken to be those that exceed the ambient by at least 3 dB while anything below that is considered to be corrupted.

4. Discussion

The theoretical results and measured data show reasonable agreement with peaks in the spectrum at roughly f = 0 . 2 Hz and a peak level between 50 and 60 dB re ( nm / s ) 2 / Hz . Note that the peak of the wave height spectral in Figure 2 is at roughly κ = 0 . 04 rad/m which equates to a surface wave frequency σ = 0 . 63 rad/s. From the frequency doubling effect discussed in Section 2.3, this gives a peak in the acoustic field at ω = 1 . 26 rad/s or f = 0 . 2 Hz which we see in the received field in Figure 5.
There is a lesser peak in the model result at around 0.1 Hz that is not obvious in the measured data. This peak in the model comes from the model for the Rayleigh wave propagation and suggests perhaps a limitation of the single layer (fluid water layer over a solid half-space) approximation adopted here.
It should be noted that the horizontal- and vertical-velocity measurements are consistent with Rayleigh waves. There is no evidence that other seismic wave types play a significant role.

5. Conclusions

Here we present an analytic model, based on the non-linear wave equation, to describe the microseisms generated by a hurricane. This model is ideal for hurricane-generated microseisms since it can be used to calculate the acoustic field due to spatially inhomogeneous surface waves. Also, this model may be used in range-dependent waveguide environments as is the case when a hurricane at sea generates microseisms that propagate up the continental margin to a receiver on land. This modeling is useful because microseisms are a primary cause of noise in seismic measurements [7,9,10] that raise the detection threshold for monitoring earthquakes [11] and tsunamis.
Based on the ocean surface directional wave spectrum in hurricane Bonnie [12,13], we predict the microseismic source levels generated by the non-linear interaction of the ocean surface waves. We then model propagation of the microseismic field from hurricane Bonnie in the North Atlantic to a seismometer in Florida to hindcast the measured signal. We find that these results compare reasonably well with seismic measurements.

Funding

Funding for open access publication provided by General Dynamics—Applied Physical Sciences.

Acknowledgments

Portions of this work were funded by the Office of Naval Research Graduate Traineeship Award in Ocean Acoustics. Wave Watch 3 (WW3) surface wave models for hurricane Bonnie were provided by Il-Ju Moon of the University of Rhode Island. The IRIS seismometers are part of the Global Seismic Network (GSN) and are installed, maintained and operated by the United States Geological Survey (USGS) Albuquerque Seismological Laboratory (http://aslwww.cr.usgs.gov and http://www.liss.org). The GSN is a cooperative scientific facility operated jointly by the Incorporated Research Institutions for Seismology (IRIS), the USGS, and the National Science Foundation (NSF).

Conflicts of Interest

The authors declare no conflict of interest.

Appendix A. Spectral Properties Surface Waves

Appendix A.1. Homogeneous Surface Wave Fields

In Equation (1) we define the complex surface wave height in the form of an inverse Fourier transform such that
ζ ^ ( ρ , 0 ) = - A ( κ ) e i κ · ρ d 2 κ
The surface wave height spectra A ( κ ) can then be written as the Fourier transform
A ( κ ) = 1 4 π 2 - ζ ^ ( ρ , 0 ) e - i κ · ρ d 2 ρ
We assume that surface wave height has zero mean ( < ζ ^ ( ρ , 0 ) > = 0 ). Taking the expectation of Equation (A2)
< A ( κ ) > = 1 4 π 2 - < ζ ^ ( ρ , 0 ) > e - i κ · ρ d 2 ρ
we can see that the wave height spectra must also have zero mean ( < A ( κ ) > = 0 ).
The second moment < A ( κ ) A * ( κ ) > can also be derived using Equation (A2) such that
< A ( κ ) A * ( κ ) > = < 1 4 π 2 - ζ ^ ( ρ , 0 ) e - i κ · ρ d 2 ρ 1 4 π 2 - ζ ^ * ( ρ , 0 ) e i κ · ρ d 2 ρ > = 1 16 π 4 - - < ζ ^ ( ρ , 0 ) ζ ^ * ( ρ , 0 ) > e - i κ · ρ e i κ · ρ d 2 ρ d 2 ρ
If we define ρ = ρ - ρ we can write this as
< A ( κ ) A * ( κ ) > = 1 16 π 4 - - R ( ρ ) e - i κ · ρ e i κ · ( ρ - ρ ) d 2 ρ d 2 ρ = 1 16 π 4 - R ( ρ ) e - i κ · ρ d 2 ρ - e - i κ · ρ e i κ · ρ d 2 ρ = S ( κ ) δ ( κ - κ )
where R ( ρ ) = < ζ ^ ( ρ , 0 ) ζ ^ * ( ρ , 0 ) > is the wave height correlation and the power spectral density S ( κ ) is its Fourier transform. This expression shows that, for homogeneous surface wave fields, the different wavenumber components are uncorrelated.
In Section 2.1 we also assume, based on the Central Limit Theorem, that the complex wave height ζ ^ ( ρ , t ) is a Gaussian random variable. Since the spectrum A ( κ ) is the Fourier transform of the complex wave height (Equation (A2)), it must also be a Gaussian random variable. The fourth moment of a Gaussian random variable can be written in terms of its second moments such that
< A ( κ 1 ) A ( κ 2 ) A * ( κ 3 ) A * ( κ 4 ) > = δ ( κ 1 - κ 3 ) S ( κ 1 ) δ ( κ 2 - κ 4 ) S ( κ 2 ) + δ ( κ 1 - κ 4 ) S ( κ 1 ) δ ( κ 2 - κ 3 ) S ( κ 2 ) .
We can also examine the cross correlation between two surface wave processes ζ ^ m and ζ ^ n where
ζ ^ m ( ρ , 0 ) = - A m ( κ ) e i κ · ρ d 2 κ
end
ζ ^ n ( ρ , 0 ) = - A n ( κ ) e i κ · ρ d 2 κ
The surface wave height spectra can then be written as the Fourier transform
A m ( κ ) = 1 4 π 2 - ζ ^ m ( ρ , 0 ) e - i κ · ρ d 2 ρ
and
A n ( κ ) = 1 4 π 2 - ζ ^ n ( ρ , 0 ) e - i κ · ρ d 2 ρ
Again we assume that surface wave heights and the wave height spectra have zero mean so that the second moment < A m ( κ ) A n * ( κ ) > can be written as
< A m ( κ ) A n * ( κ ) > = < 1 4 π 2 - ζ ^ m ( ρ , 0 ) e - i κ · ρ d 2 ρ 1 4 π 2 - ζ ^ n * ( ρ , 0 ) e i κ · ρ d 2 ρ > = 1 16 π 4 - - < ζ ^ m ( ρ , 0 ) ζ ^ n * ( ρ , 0 ) > e - i κ · ρ e i κ · ρ d 2 ρ d 2 ρ = 1 16 π 4 - - R m n ( ρ ) e - i κ · ρ e i κ · ( ρ - ρ ) d 2 ρ d 2 ρ = 1 16 π 4 - R m n ( ρ ) e - i κ · ρ d 2 ρ - e - i κ · ρ e i κ · ρ d 2 ρ = S m n ( κ ) δ ( κ - κ )
where R m n ( ρ ) = < ζ ^ m ( ρ , 0 ) ζ ^ n * ( ρ , 0 ) > is the wave height cross correlation and the cross-spectral density S m n ( κ ) is its Fourier transform.
Again we assume that the spectra A m ( κ ) and A n ( κ ) are Gaussian random variables so that the fourth moment can be written as
< A m ( κ 1 ) A m ( κ 2 ) A n * ( κ 3 ) A n * ( κ 4 ) > = δ ( κ 1 - κ 3 ) S m n ( κ 1 ) δ ( κ 2 - κ 4 ) S m n ( κ 2 ) + δ ( κ 1 - κ 4 ) S m n ( κ 1 ) δ ( κ 2 - κ 3 ) S m n ( κ 2 ) .

Appendix A.2. Inhomogeneous Surface Wave Fields

Now we consider the case where the surface wave height is defined over a finite area. The complex surface wave height in a finite region centered at the origin can be written as
ζ ^ w ( ρ , 0 ) = ζ ^ ( ρ , t ) w ( ρ )
where the window function w is unity for - L x / 2 < x < L x / 2 and - L y / 2 < y < L y / 2 and zero elsewhere. Taking the Fourier transform (Equation (A2)) of both sides of Equation (A13) yields
A w ( κ ) = A ( κ ) W ( κ ) = - A ( κ ) W 0 ( κ - κ ) d 2 κ
The second moment of A w can then be written as
< A w ( κ ) A w * ( κ ) > = - W ( κ - κ ) d 2 κ - W * ( κ - κ ) d 2 κ < A ( κ ) A * ( κ ) >
which, from Equation (A5) becomes
< A w ( κ ) A w * ( κ ) > = - W ( κ - κ ) d 2 κ - W * ( κ - κ ) d 2 κ S ( κ ) δ ( κ - κ ) = - W ( κ - κ ) W * ( κ - κ ) S ( κ ) d 2 κ
As the dimensions L x and L y become much larger than a wavelength λ g the functions W begin to approximate delta functions such that
< A w ( κ ) A w * ( κ ) > δ ( κ - κ ) δ ( κ - κ ) S ( κ ) d 2 κ δ ( κ - κ ) S ( κ ) = < A ( κ ) A * ( κ ) >
We can now consider the more complicated case of the cross power spectral density between the surface wave heights in two finite regions m and n.
ζ ^ w , m ( ρ , 0 ) = ζ ^ ( ρ , 0 ) w ( ρ - ρ m )
and
ζ ^ w , n ( ρ , 0 ) = ζ ^ ( ρ , 0 ) w ( ρ - ρ n )
which, after taking the Fourier transform, become
A w , m ( κ ) = A m ( κ ) * ( W ( κ ) e - i κ · ρ m )
and
A w , n ( κ ) = A n ( κ ) * ( W ( κ ) e - i κ · ρ n )
As before, we can take the second moment which becomes
< A w , m ( κ ) A w , n * ( κ ) > = - W ( κ - κ ) e - i ( κ - κ ) · ρ m d 2 κ - W * ( κ - κ ) e i ( κ - κ ) · ρ n d 2 κ < A m ( κ ) A n * ( κ ) > = - W ( κ - κ ) e - i ( κ - κ ) · ρ m d 2 κ - W * ( κ - κ ) e i ( κ - κ ) · ρ n d 2 κ S m n ( κ ) δ ( κ - κ ) = - W ( κ - κ ) e - i ( κ - κ ) · ρ m W * ( κ - κ ) e i ( κ - κ ) · ρ n S m n ( κ ) d 2 κ
Again, as the dimensions L x and L y become much larger than a wavelength λ g the functions W begin to approximate delta functions such that
< A w , m ( κ ) A w , n * ( κ ) > e - i ( κ - κ ) · ( ρ m - ρ n ) δ ( κ - κ ) S m n ( κ )
This expression is similar to the previous Equations (A5) and (A17) with the addition of an exponential phase term. Since the delta function is zero for all κ - κ 0 and since the phase term is unity when κ - κ = 0 this expression can be reduced to
< A w , m ( κ ) A w , n * ( κ ) > δ ( κ - κ ) S m n ( κ ) = < A m ( κ ) A n * ( κ ) >

Appendix B. Adiabatic Propagation of Generalized Rayleigh Waves in a Range-Dependent Ocean Environment

To determine the microseisms received by a sensor either in the ocean or in the earth’s crust, the Green function of Equation (38) must be calculated. It has been shown both theoretically [10,14] and experimentally [4,33,34,35,36,37] that microseisms propagate as Rayleigh waves along the sea floor.
Stoneley [61] shows how the depth or thickness of the ocean layer affects Rayleigh wave phase speed and how, for large depths or frequencies, there may be multiple propagating Rayleigh modes. Press and Ewing [10], Press and Ewing [62] and Ellis and Chapman [63] later show that the Rayleigh wave field can be expressed as a sum of normal-mode contributions. Unlike ‘classic’ Rayleigh waves which are defined to propagate along a vacuum/elastic boundary, Ewing et al. [64] coin the term ‘generalized Rayleigh waves’ to describe Rayleigh waves that propagate along an elastic boundary, such as at the sea floor, under a finite thickness fluid layer like the ocean. The modal expressions of [10,62,63] are only applicable, however, to environments where the ocean depth is constant. Varying water depth at the continental margin can affect the Rayleigh wave propagation.
To solve the range-dependent propagation problem, Arvello and Überall [65] divide the waveguide into range-independent segments, each with its own numerically calculated eigenvalues and eigenfunctions, and apply the adiabatic normal-mode approximation to evaluate the seismic and acoustic propagation between segments. The computational nature of their approach, however, does not provide the physical insight given by the analytic results of the earlier range-independent studies [10,61,62,63,64]. Here we provide an analytic solution for the range-dependent case where a fluid ocean of varying depth overlays an elastic ocean floor, as opposed to the computational solution method given by [65].
We consider the simple case where a monopole source is located in a homogeneous fluid medium overlaying an elastic half space. Unlike previous range-independent studies, here the depth or thickness of the fluid layer may change with range from the source. For the example in this paper we model an up-sloping environment where depth decreases with range, however, other depth profiles may also be considered.
It should be noted that, in addition to modeling microseisms, this range-dependent Rayleigh wave model also has application in other seismic research where a sensor on land measures signals generated by sources at sea. For example undersea earthquakes [66] and underwater explosions, like those monitored under the Comprehensive Nuclear Test Ban Treaty [67,68], also produce Rayleigh waves on the sea floor. The propagation of these Rayleigh waves is also often measured to infer the geologic characteristics of the ocean floor using surface-wave tomography [69,70,71].
We begin by considering the range-independent form of the normal-mode solution for a source at a depth z 0 in a constant-depth fluid layer overlaying an elastic halfspace [65]
ϕ 1 ( ρ , z ) = n u 1 , n ( z s ) u 1 , n ( z r ) H 0 ( 1 ) ( k n ρ ) ,
ϕ 2 ( ρ , z ) = n u 1 , n ( z s ) u 2 , n ( z r ) H 0 ( 1 ) ( k n ρ ) ,
and
ψ 2 ( ρ , z ) = n u 1 , n ( z s ) v 2 , n ( z r ) H 0 ( 1 ) ( k n ρ ) ,
where ρ = | ρ | is the horizontal distance from source to receiver and z r and z s are the receiver and source depth respectively. The variable ϕ 1 represents the displacement potential due to compressional waves in the fluid medium, and ϕ 2 and ψ 2 represent the compression and shear potentials in the solid medium. The variable k n is the horizontal component of the wavenumber and u 1 , n , u 2 , n and v 2 , n are the modal eigenfunctions for each mode n. The potentials are related to the horizontal displacement d r , vertical displacement d z , vertical normal stress p z z and vertical shear stress p z r by [62,64]
d r = ϕ r + 2 ψ r z d z = ϕ z + 2 ϕ z 2 + ω 2 β 2 ψ p z z = λ 2 ϕ + 2 μ w z p z r = μ ( q z + w r ) .
The Lame’s constants λ and μ are related to the compressional wave speed in the fluid v 1 and the compression and shear wave speeds in the solid α 2 and β 2 by
v 1 2 = λ 1 / ρ 1 α 2 2 = ( λ 2 + 2 μ 2 ) / ρ 2 β 2 2 = μ 2 / ρ 2 ,
where ρ 1 and ρ 2 are the densities in the fluid and elastic media respectively [62,64].
Pierce [72] has shown that for slowly varying range-dependent waveguides the modal sums of Equations (A25)–(A27) can be written as
ϕ 1 ( ρ , z ) = e - i π 4 2 π ρ × n u 1 , n ( z s ) u 1 , n ( ρ , z r ) e - i 0 ρ k n ( ρ ) d ρ k n ( ρ )
ϕ 2 ( ρ , z ) = e - i π 4 2 π ρ × n u 1 , n ( z s ) u 2 , n ( ρ , z r ) e - i 0 ρ k n ( ρ ) d ρ k n ( ρ )
and
ψ 2 ( ρ , z ) = e - i π 4 2 π ρ × n u 1 , n ( z s ) v 2 , n ( ρ , z r ) e - i 0 ρ k n ( ρ ) d ρ k n ( ρ )
based on the adiabatic mode approximation. This approximation requires that the change in the waveguide environment as a function of range ρ is negligible over a wavelength scale and that there is no coupling between modes. In these expressions both the modal eigenfunctions u 1 , n , u 2 , n and v 2 , n and the horizontal wavenumber k n are allowed to change with range ρ from the source.
For an iso-speed fluid layer of depth H ( ρ ) overlaying an elastic half-space the general solutions for the mode functions u 1 , n , u 2 , n and v 2 , n can be written in terms of sinusoids and exponentials [62,63] such that Equations (A30)–(A32) become
ϕ 1 ( ρ , z ) = e - i π 4 2 π ρ × n Φ 1 sin ( ξ n ( 0 ) z s ) sin ( ξ n ( ρ ) z r ) e - i 0 ρ k n ( ρ ) d ρ k n ( ρ )
ϕ 2 ( ρ , z ) = e - i π 4 2 π ρ × n Φ 2 sin ( ξ n ( 0 ) z s ) e η n ( ρ ) z r e - i 0 ρ k n ( ρ ) d ρ k n ( ρ )
and
ψ 2 ( ρ , z ) = e - i π 4 2 π ρ × n Ψ 2 sin ( ξ n ( 0 ) z s ) e ζ n ( ρ ) z r e - i 0 ρ k n ( ρ ) d ρ k n ( ρ )
where
ξ n ( ρ ) = ω 2 v 1 2 - k n 2 ( ρ ) = k n ( ρ ) c n 2 ( ρ ) v 1 2 - 1 η n ( ρ ) = k n 2 ( ρ ) - ω 2 α 2 2 = k n ( ρ ) 1 - c n 2 ( ρ ) α 2 2 ζ n ( ρ ) = k n 2 ( ρ ) - ω 2 β 2 2 = k n ( ρ ) 1 - c n 2 ( ρ ) β 2 2
represent the vertical components of the wavenumber vector for each mode n. The phase speed for each mode is defined as c n ( ρ ) = ω / k n ( ρ ) where ω is the frequency in radians/Section.
The amplitudes Φ 1 , Φ 2 and Ψ 2 of the mode functions in Equations (A33)–(A35) must be chosen such that the boundary conditions
( p z z ) f = ( p z z ) s ( p z r ) s = 0 ( d z ) f = ( d z ) s at z = H ( ρ )
are satisfied at the fluid solid interface z = H ( ρ ) . The subscripts f and s indicate that the variable is to be evaluated in the fluid or solid layer respectively. These boundary conditions yield a system of three equations for the unknown amplitudes Φ 1 , Φ 2 and Ψ 2 whereby any two of the unknowns may be solved in terms of the third. The third amplitude is then normalized to some convenient value. Later we will choose a normalization such that our solution is consistent with that of [62] for constant bathymetry H ( ρ ) = H . The solution to this set of equations only exists, however, if the determinant
- λ 1 ( ξ n 2 + k n 2 ) sin ( ξ n H ) ( - λ 2 ( η n 2 + k n 2 ) + 2 μ 2 η n ) e - η n H - 2 μ 2 k n 2 e - ζ n H 0 - 2 η n e - η n H ( 2 ζ n 2 + ω 2 β 2 2 ) e - ζ n H ξ n cos ( ξ n H ) η n e - η n H - k n 2 e - ζ n H = 0
equals zero [62,64] where the dependence of k n , ξ n , η n , ζ n and H on range ρ has been suppressed to simplify notation. Equation (A38) reduces to the equation
ρ 1 ρ 2 ω 4 β 2 4 η n ( ρ ) ξ n ( ρ ) tan ( ξ n ( ρ ) H ( ρ ) ) - 4 k n 2 ( ρ ) η n ( ρ ) ζ n ( ρ ) - 2 k n 2 ( ρ ) - ω 2 β 2 2 2 = 0
where the roots k n ( ρ ) correspond to the horizontal wavenumbers for each mode n. For the case of constant bathymetry H ( ρ ) = H , Equation (A39) is identical to those given in earlier works [10,62,64]; however, for this more general range-dependent case the roots k n ( ρ ) can now vary with range ρ . Unlike a ‘classic’ Rayleigh wave which has only one root or propagating value, ‘generalized’ Rayleigh waves may have multiple roots k n ( ρ ) depending on frequency and fluid layer thickness as shown in Equation (A39) for the range-dependent case and as shown by [10] for the range-independent case. Note that Equation (A39) has no roots for k n ( ρ ) > ω / v 1 .
Given the roots k n ( r ) from Equation (A39), we can express the potentials as
ϕ 1 ( ρ , z ) = 2 H ( ρ ) 2 π ρ e - i π 4 n Φ 1 ( k n ( ρ ) , H ( ρ ) ) × sin ( ξ n ( 0 ) z s ) sin ( ξ n ( ρ ) z ) e - i 0 ρ k n ( ρ ) d ρ k n ( ρ )
ϕ 2 ( ρ , z ) = 2 H ( ρ ) 2 π ρ e - i π 4 n Φ 2 ( k n ( ρ ) , H ( ρ ) ) × sin ( ξ n ( 0 ) z s ) e η n ( ρ ) z e - i 0 ρ k n ( ρ ) d ρ k n ( ρ )
and
ψ 2 ( ρ , z ) = 2 2 π ρ e - i π 4 n Ψ 2 ( k n ( ρ ) , H ( ρ ) ) × sin ( ξ n ( 0 ) z s ) e ζ n ( ρ ) z e - i 0 ρ k n ( ρ ) d r k n ( ρ )
where the solutions for the modal amplitudes are
Φ 1 ( k n ( r h o ) , H ( ρ ) ) = - ρ 1 ρ 2 c n 4 ( ρ ) β 2 4 η n ( ρ ) ξ n ( ρ ) k n ( ρ ) H ( ρ ) c n 2 ( ρ ) / v 1 2 - 1 cos ( ξ n ( ρ ) H ( ρ ) ) D ( k n ( r h o ) , H ( ρ ) )
Φ 2 ( k n ( ρ ) , H ( ρ ) ) = - ρ 1 ρ 2 c n 2 ( ρ ) β 2 2 ( 2 - c n 2 ( ρ ) / β 2 2 ) k n ( ρ ) H ( ρ ) c n 2 ( ρ ) / v 1 2 - 1 D ( k n ( r h o ) , H ( ρ ) )
and
Ψ 1 ( k n ( ρ ) , H ( ρ ) ) = - 2 ρ 1 ρ 2 c n 2 ( ρ ) β 2 2 η n ( ρ ) ξ n ( ρ ) D ( k n ( r h o ) , H ( ρ ) )
and where
D ( k n ( r h o ) , H ( ρ ) ) = ρ 1 ρ 2 c n 4 ( ρ ) β 2 4 [ sin ( ξ n ( ρ ) H ( ρ ) ) c n 2 ( ρ ) / v 1 2 - 1 1 - c n 2 ( ρ ) / α 2 2 1 + 1 - c n 2 ( ρ ) / α 2 2 c n 2 ( ρ ) / v 1 2 - 1 - k n ( ρ ) H ( ρ ) 1 - c n 2 ( ρ ) / α 2 2 c n 2 ( ρ ) / v 1 2 - 1 sec ( ξ n ( ρ ) H ( ρ ) ) ] - 4 [ 1 - c n 2 ( ρ ) / β 2 2 1 - c n 2 ( ρ ) / α 2 2 + 1 - c n 2 ( ρ ) / α 2 2 1 - c n 2 ( ρ ) / β 2 2 + 2 1 - c n 2 ( ρ ) / α 2 2 1 - c n 2 ( ρ ) / β 2 2 - 2 ( 2 - c n 2 ( ρ ) / β 2 2 ) ] cos ( ξ n ( ρ ) H ( ρ ) )
Equations (A40)–(A46) provide an analytic model for the propagation of ‘generalized’ Rayleigh waves in a range-dependent environment. Note, however, that the mode amplitudes have been normalized such that Φ 1 ( k n ( ρ ) , H ( ρ ) ) , Φ 2 ( k n ( ρ ) , H ( ρ ) ) and Ψ 1 ( k n ( ρ ) , H ( ρ ) ) are the same as Φ 1 ( k n , H ) , Φ 2 ( k n , H ) and Ψ 1 ( k n , H ) of [62] for the case of constant bathymetry H ( ρ ) = H (Equation (38) of [62] is missing a factor of 2 which was later corrected in Equations (4)–(184) of [64]). The Green function used in this paper can be written as
g ( r r , r s , w ) = 4 π f p z z ( | ρ r , ρ s | , z r , z s )
where the 4 π f accounts for the normalization adopted by [62].

Appendix C. Range-Independent Half-Space

To compare the expression in Equation (38) with derivations by others we can simplify our solution for the case of a range-independent surface wave field over an infinite ocean half-space. The Green function for a source near the free surface of an infinite ocean half-space can be written as a dipole
g ( r r , r m , ω ) = e i k R m 2 π R m k z s c o s ( β m ) ,
where R m = | r r - r m | is the distance between source and receiver positions r m and r r respectively and β m = Atan ( | ρ r - ρ m | / | z r - z m | ) is the angle from vertical. From Equation (38) we write the power spectral density of the pressure field, substituting Equation (A48) as
S P ( r r , ω ) ρ 2 π L x 2 L y 2 k 2 ω 9 512 g 2 S ω 2 , θ S ω 2 , θ + π d θ × 0 0 d z s d z s e - ω 2 z s 2 g e - ω 2 z s 2 g × m , n e i k R m 2 π R m k z s c o s ( β m ) e - i k R n 2 π R n k z s c o s ( β n )
For densely spaced r m and r n the sums can be approximated as integrals on ρ s and ρ s
S P ( r r , ω ) ρ 2 π L x 2 L y 2 k 2 ω 9 512 g 2 S ω 2 , θ S ω 2 , θ + π d θ × 0 0 d z s d z s e - ω 2 z s 2 g e - ω 2 z s 2 g × - - d 2 ρ s L x L y e i k R 2 π R k z s c o s ( β ) × - - d 2 ρ s L x L y e - i k R 2 π R k z s c o s ( β )
For k z r > > 1 the integral over ρ s can be approximated as z s L x L y e i k ( z r - z s ) + i π / 2 and the integral over ρ s as z s L x L y e - i k ( z r - z s ) - i π / 2 so that
S P ( r r , ω ) ρ 2 π k 2 ω 9 512 g 2 S ω 2 , θ S ω 2 , θ + π d θ × 0 0 d z s d z s z s z s e - ω 2 z s 2 g e - ω 2 z s 2 g e i k z s e - i k z s .
Integrating over z s then yields
S P ( r r , ω ) ρ 2 π k 2 ω g 2 32 S ω 2 , θ S ω 2 , θ + π d θ
since ω 2 2 g > > k at the frequencies considered here. This result is identical to Equation (33) of [18] and Equation (35) of [20] where S ( σ , θ ) equates to 8 ω f 3 g 2 X ( k 1 ) G ( θ ) from Hughes and to 8 σ 3 g 2 f ( κ ) from Lloyd (including the 1 2 correction of Hughes derivation discussed by Lloyd).

References

  1. Chi, W.-C.; Chen, W.-J.; Kuo, B.-Y.; Dolenc, D. Seismic Monitoring of Western Pacific Typhoons. Mar. Geophys. Res. 2010, 31, 239–251. [Google Scholar] [CrossRef]
  2. Davy, C.; Barruol, G.; Fontaine, F.R.; Sigloch, K.; Stutzmann, E. Tracking Major Storms from Microseismic and Hydroacoustic Observations on the Seafloor. Geophys. Res. Lett. 2014, 10, 8825–8831. [Google Scholar] [CrossRef]
  3. Gerstoft, P.; Shearer, P.M.; Harmon, N.; Zhang, J. Global P, PP and PKP Wave Microseisms Observed from Distant Storms. Geophys. Res. Lett. 2008, 35, L23306. [Google Scholar] [CrossRef]
  4. Gilmore, M.H. Microseisms and ocean storms. Seismol. Soc. Bull. 1946, 36, 89–119. [Google Scholar]
  5. Lin, J.; Wang, Y.; Wang, W.; Li, X.; Fang, S.; Chen, C.; Zheng, H. Seismic Remote Sensing of Super Typhoon Lupit (2009) with Seismological Array Observations in NE China. Remote Sens. 2018, 10, 235. [Google Scholar] [CrossRef]
  6. Lin, J.; Lin, J.; Xu, M. Microseisms Generated by Super Typhoon Megi in the Western Pacific Ocean. J. Geophys. Res. Oceans 2017, 10, 9518–9529. [Google Scholar] [CrossRef]
  7. Ramirez, J.E. An experimental investigation of the nature and origin of microseisms at St. Louis, Missouri. Bull. Seismol. Soc. Am. 1940, 30, 35–84. [Google Scholar]
  8. Zhang, J.; Gerstoft, P.; Bromirski, P.D. Pelagic and Coastal Sources of P-Wave Microseisms. Geophys. Res. Lett. 2010, 37, L15301. [Google Scholar] [CrossRef]
  9. Lee, A.W. The effect of geologic structure upon microseismic disturbance. R. Astron. Soc. Mon. Not. Geophys. Suppl. 1932, 3, 83–105. [Google Scholar] [CrossRef]
  10. Press, F.; Ewing, M. A Theory of Microseisms with Geologic Applications. Trans. Am. Geophys. Union 1948, 29, 163–174. [Google Scholar] [CrossRef]
  11. Webb, S.C. Broadband seismology and noise under the ocean. Rev. Geophys. 1998, 36, 105–142. [Google Scholar] [CrossRef] [Green Version]
  12. Moon, I.-J.; Ginis, I.; Hara, T.; Tolman, H.L.; Wright, C.W.; Walsh, E.J. Numerical simulation of sea surface directional wave spectra under hurricane wind forcing. J. Phys. Oceanogr. 2003, 33, 1680–1706. [Google Scholar] [CrossRef]
  13. Wright, C.W.; Walsh, E.J.; Vandemark, D.; Krabill, W.B.; Garcia, A.W.; Houston, S.H.; Powell, M.D.; Black, P.G.; Marks, F.D., Jr. Hurricane directional wave spectrum spatial variation in the open ocean. J. Phys. Oceanogr. 2001, 31, 2472–2488. [Google Scholar] [CrossRef]
  14. Longuet-Higgins, M.S. A theory of the origin of microseisms. Philos. Trans. R. Soc. Lond. Ser. A 1950, 243, 1–35. [Google Scholar] [CrossRef]
  15. Brekhovskikh, L.M. Underwater sound waves generated by surface waves in the ocean. Izv. Atmos. Ocean. Phys. 1966, 2, 582–587. [Google Scholar]
  16. Guralnik, Z.; Bourdelais, J.; Zabalgogeazcoa, X.; Farrell, W.E. Wave-Wave Interactions and Deep Ocean Acoustics. J. Acoust. Soc. Am. 2013, 134, 3161–3173. [Google Scholar] [CrossRef] [PubMed]
  17. Hasselmann, K. A statistical analysis of the generation of microseisms. Rev. Geophys. 1963, 1, 177–210. [Google Scholar] [CrossRef]
  18. Hughes, B. Estimates of underwater sound (and infrasound) produced by nonlinearly interacting ocean waves. J. Acoust. Soc. Am. 1976, 60, 1032–1039. [Google Scholar] [CrossRef]
  19. Kibblewhite, A.C.; Wu, C.Y. The generation of infrasonic ambient noise in the ocean by nonlinear interactions of ocean surface waves. J. Acoust. Soc. Am. 1989, 85, 1935–1945. [Google Scholar] [CrossRef]
  20. Lloyd, S.P. Underwater sound from surface waves according to the Lighthill-Ribner theory. J. Acoust. Soc. Am. 1981, 69, 425–435. [Google Scholar] [CrossRef]
  21. Kibblewhite, A.C.; Wu, C.Y. The theoretical description of wave-wave interactions as a noise source in the ocean. J. Acoust. Soc. Am. 1991, 89, 2241–2252. [Google Scholar] [CrossRef]
  22. Kibblewhite, A.C.; Wu, C.Y. Acoustic source levels associated with the nonlinear interactions of ocean waves. J. Acoust. Soc. Am. 1993, 94, 3358–3378. [Google Scholar] [CrossRef]
  23. Cato, D.H. Sound generation in the vicinity of the sea surface: Source mechanisms and the coupling to the received sound field. J. Acoust. Soc. Am. 1991, 89, 1076–1095. [Google Scholar] [CrossRef]
  24. Cato, D.H. Theoretical and measured underwater noise from surface wave orbital motion. J. Acoust. Soc. Am. 1991, 89, 1096–1112. [Google Scholar] [CrossRef]
  25. Westervelt, P.J. Parametric Acoustic Array. J. Acoust. Soc. Am. 1963, 35, 535–537. [Google Scholar] [CrossRef]
  26. Banerji, S.K. Microseisms associated with disturbed weather in the Indian seas. Philos. Trans. R. Soc. Lond. Ser. A 1930, 229, 287–328. [Google Scholar] [CrossRef]
  27. Bowen, S.P.; Richard, J.C.; Mancini, J.D.; Fessatidis, V.; Crooker, B. Microseism and infrasound generation by cyclones. J. Acoust. Soc. Am. 2003, 113, 2562–2573. [Google Scholar] [CrossRef] [PubMed]
  28. Banerji, S.K. Theory of microseisms. Proc. Indian Acad. Sci. 1935, 1, 727–753. [Google Scholar]
  29. Darbyshire, J. Microseisms, in the Sea, Ideas and Observations on Progress in the Study of the Seas; Hill, M.N., Ed.; Interscience: New York, NY, USA, 1966; p. 701. [Google Scholar]
  30. Wilson, J.D. Classifying Hurricanes with Undersea Sound. Ph.D. Thesis, Massachusetts Institute of Technology, Cambridge, UK, 2006. [Google Scholar]
  31. Miche, M. Mouvements ondulatoires de la mer en profonder constante ou decroissante. Ann. Ponts Chaussées 1944, 114, 131–164. [Google Scholar]
  32. Walsh, E.J.; Wright, C.W.; Vandemark, D.; Krabill, W.B.; Garcia, A.W.; Houston, S.H.; Murillo, S.T.; Powell, M.D.; Black, P.G.; Marks, F.D., Jr. Hurricane directional wave spectrum spatial variation at landfall. J. Phys. Oceanogr. 2002, 32, 1667–1684. [Google Scholar] [CrossRef]
  33. Bernard, P. Short Review and Recent Results on Microseisms. Geophys. Surv. 1983, 5, 395–407. [Google Scholar] [CrossRef]
  34. Bradner, H.; De Jerphanion, L.G.; Langlois, R. Ocean Microseism Measurements with a Neutral Buoyancy Free-Floating Midwater Seismometer. Bull. Seismol. Soc. Am. 1970, 60, 1134–1150. [Google Scholar]
  35. Goodman, D.; Yamamoto, T.; Trevorrow, M.; Abbott, C.; Turgut, A.; Badiey, M.; Ando, K. Directional Spectra Observations of Seafloor Microseisms from an Ocean-Bottom Seismometer Array. J. Acoust. Soc. Am. 1989, 86, 2309–2317. [Google Scholar] [CrossRef]
  36. Gutenberg, B. Microseisms in North America. Bull. Seismol. Soc. Am. 1931, 21, 1–24. [Google Scholar]
  37. Lee, A.W. On the direction of Approach of Microseismic Waves. Proc. R. Soc. Lond. Ser. A 1935, 886, 183. [Google Scholar] [CrossRef]
  38. Haubrich, R.A. Earth noise, 5 to 500 millicycles per second 1. Spectral stationarity, normality, and nonlinearity. J. Geophys. Res. 1965, 70, 1415–1427. [Google Scholar] [CrossRef]
  39. Haubrich, R.A.; MacKenzie, G.S. Earth noise, 5 to 500 millicycles per second 2. Reaction of the earth to oceans and atmosphere. J. Geophys. Res. 1965, 70, 1429–1440. [Google Scholar] [CrossRef]
  40. Haubrich, R.A.; McCamy, K. Microseisms: coastal and pelagic sources. Rev. Geophys. 1969, 7, 539–571. [Google Scholar] [CrossRef]
  41. Kibblewhite, A.C.; Ewans, K.C. Wave-wave interactions, microseisms, and infrasonic ambient noise in the ocean. J. Acoust. Soc. Am. 1985, 78, 981–994. [Google Scholar] [CrossRef]
  42. Nichols, R.H. Infrasonic ambient ocean noise measurements: Eleuthera. J. Acoust. Soc. Am. 1981, 69, 974–981. [Google Scholar] [CrossRef]
  43. Schimmel, M.; Stutzmann, E.; Ardhuin, F.; Gallart, J. Polarized Earth’s ambient microseismic noise. Geochem. Geophys. Geosyst. 2005, 12, 1–14. [Google Scholar] [CrossRef]
  44. Webb, S.C. The equilibrium oceanic microseism spectrum. J. Acoust. Soc. Am. 1992, 92, 2141–2158. [Google Scholar] [CrossRef]
  45. Holland, G.J. (Ed.) Global Guide to Tropical Cyclone Forecasting; World Meteorological Organization: Geneva, Switzerland, 1993. [Google Scholar]
  46. Komen, G.J.; Cavaleri, L.; Donelan, M.; Hasselmann, K.; Hasselmann, S.; Janssen, P.A.E.M. Dynamics and Modeling of Ocean Waves; Cambridge University Press: Cambridge, UK, 1994; pp. 16–27. [Google Scholar]
  47. Kinsman, B. Wind Waves, Their Generation and Propagation on the Sea Surface; Prentice-Hall: Englewood Cliffs, NJ, USA, 1965; pp. 336–352. [Google Scholar]
  48. Mandel, L.; Wolf, E. Optical Coherence and Quantum Optics; Cambridge University Press: Cambridge, UK, 1995; pp. 58–59. [Google Scholar]
  49. Forristall, G.Z.; Ward, E.G.; Cardone, V.J.; Borgmann, L.E. The directional spectra and kinematics of surface gravity waves in tropical storm Delia. J. Phys. Oceanogr. 1978, 8, 888–909. [Google Scholar] [CrossRef]
  50. Toleman, H.L. User Manual and System Documentation of WAVEWATCH-III; Version 1.18; Techreport Note. 166; Ocean Modeling Branch, NCEP National Weather Service, NOAA, US Department of Commerce: College Park, MD, USA, 1999; p. 110.
  51. Dean, L.W., III. Interactions between sound waves. J. Acoust. Soc. Am. 1962, 34, 1039–1044. [Google Scholar] [CrossRef]
  52. Westervelt, P.J. Scattering of sound by sound. J. Acoust. Soc. Am. 1957, 29, 199–203. [Google Scholar] [CrossRef]
  53. Morse, P.M.; Ingard, K.U. Theoretical Acoustics; Mc-Graw-Hill: New York, NY, USA, 1968; pp. 245–866. [Google Scholar]
  54. Beyer, R.T. Nonlinear Acoustics; Acoustical Society of America, Woodbury: New York, NY, USA, 1997; p. 101. [Google Scholar]
  55. Van Trees, H.L. Optimum Array Processing. Part IV of Detection, Estimation and Modulation Theory; Wiley-Interscience: New York, NY, USA, 2002; pp. 95–118. [Google Scholar]
  56. Papoulis, A.; Pillai, S.U. Probability, Random Variables and Stochastic Processes; McGraw-Hill: New York, NY, USA, 1965; p. 515. [Google Scholar]
  57. Tucker, M.J.; Pitt, E.G. Waves in Ocean Engineering; Elsevier: New York, NY, USA, 2001; pp. 33–34. [Google Scholar]
  58. Avila, L.A. Preliminary Report, Hurricane Bonnie 19–30 August 1998; National Hurricane Center: Miami, FL, USA, 1998.
  59. Makris, N.C.; Avelino, L.Z.; Menis, R. Deterministic Reverberation from Ocean Ridges. J. Acoust. Soc. Am. 1995, 97, 3547–3574. [Google Scholar] [CrossRef]
  60. Lay, T.; Wallace, T.C. Modern Global Seismology; Academic Press: New York, NY, USA, 1995; pp. 252–254. [Google Scholar]
  61. Stoneley, R. The Effect of the Ocean on Rayleigh Waves. Geophys. Suppl. Mon. Not. R. Astron. Soc. 1926, 1, 349–356. [Google Scholar] [CrossRef]
  62. Press, F.; Ewing, M. Propagation of Explosive Sound in a Liquid Layer Overlying a Semi-Infinite Elastic Solid. Geophysics 1950, 15, 426–446. [Google Scholar] [CrossRef]
  63. Ellis, D.D.; Chapman, D.M.F. A Simple Shallow Water Propagation Model Including Shear Wave Effects. J. Acoust. Soc. Am. 1985, 78, 2087–2095. [Google Scholar] [CrossRef]
  64. Ewing, W.M.; Jardetzky, W.S.; Press, F. Elastic Waves in Layered Media; McGraw-Hill: New York, NY, USA, 1957; pp. 157–189. [Google Scholar]
  65. Arvello, J.I.; Überall, H. Adiabatic Normal-Mode Theory of Sound Propagating Including Shear Waves in a Range-Dependent Ocean Floor. J. Acoust. Soc. Am. 1990, 88, 2316–2325. [Google Scholar] [CrossRef]
  66. Sutherland, F.H.; Vernon, F.L.; Orcutt, J.A.; Collins, J.A.; Stephen, R.A. Results from OSNPE; improved teleseismic earthquake detection at the seafloor. Bull. Seism. Soc. Am. 2001, 94, 1868–1878. [Google Scholar] [CrossRef]
  67. Mitchell, B.J. Monitoring the Comprehensive Nuclear-Test-Ban Treaty; Preface. Pure Appl. Geophys. 2001, 158, 1339–1340. [Google Scholar] [CrossRef]
  68. Ritzwoller, M.H.; Levshin, A.L. Monitoring the Comprehensive Nuclear-Test-Ban Treaty; Introduction. Pure Appl. Geophys. 2001, 158, 1341–1348. [Google Scholar]
  69. Barmin, M.P.; Ritzwoller, M.H.; Levshin, A.L. A fast and reliable method for surface wave tomography. Pure Appl. Geophys. 2001, 158, 1351–1375. [Google Scholar] [CrossRef]
  70. Pilidou, S.; Priestley, K.; Debayle, E.; Gudmundsson, Ó. Rayleigh wave tomography in the North Atlantic: High resolution images of the Iceland, Azores and Eifel mantle plumes. Lithos 2005, 79, 453–474. [Google Scholar] [CrossRef]
  71. Singh, D.D. Rayleigh wave group-velocity studies beneath the Indian Ocean. Bull. Seismol. Soc. Am. 2005, 95, 502–511. [Google Scholar] [CrossRef]
  72. Pierce, A.D. Extension of the Method of Normal Modes to Sound Propagation in an Almost-Stratified Medium. J. Acoust. Soc. Am. 1965, 37, 19–27. [Google Scholar] [CrossRef]
Figure 1. (A) Wind speed in m/s and (B) microseismic source level L B m m ( ω , z s = 0 ) at ω = 1 . 26 rad/s ( f = 0 . 2 Hz) in dB re Pa 2 s/m 4 from Equation (39) at 1200 on 24 August 1998 as a function of latitude and longitude. The arrow indicates the direction hurricane Bonnie was moving. The letters (a–d) represent features of interest; (a) indicates the location of maximum wind speed, (b) indicates the eye where wind speed is zero, and (c,d) indicate peaks in the microseismic source level. This figure shows that, while a hurricane can produce significant microseismic source levels (B), these level do not directly follow the wind speeds (A) in the hurricane.
Figure 1. (A) Wind speed in m/s and (B) microseismic source level L B m m ( ω , z s = 0 ) at ω = 1 . 26 rad/s ( f = 0 . 2 Hz) in dB re Pa 2 s/m 4 from Equation (39) at 1200 on 24 August 1998 as a function of latitude and longitude. The arrow indicates the direction hurricane Bonnie was moving. The letters (a–d) represent features of interest; (a) indicates the location of maximum wind speed, (b) indicates the eye where wind speed is zero, and (c,d) indicate peaks in the microseismic source level. This figure shows that, while a hurricane can produce significant microseismic source levels (B), these level do not directly follow the wind speeds (A) in the hurricane.
Remotesensing 10 01624 g001
Figure 2. The wave height power spectral level ( 10 log S m m ( κ ) /(m 2 /(rad 2 /m 2 ))) in dB re m 2 / ( rad 2 / m 2 ) at the locations of interest (ad) given in Figure 1 at 1200 on 24 August 1998. The peak in the spectra is at roughly κ = 0 . 04 rad/m or σ = 0 . 63 rad/s which corresponds to an acoustic frequency of ω = 1 . 26 rad/s or f = 0 . 2 Hz. At some locations (bd) there are waves propagating with opposing wavenumber vectors κ , while at other locations (a), most of the waves propagate in the same direction. From Equation (37) we expect these locations with opposing waves to produce the greatest microseismic source levels and shown in Figure 1B.
Figure 2. The wave height power spectral level ( 10 log S m m ( κ ) /(m 2 /(rad 2 /m 2 ))) in dB re m 2 / ( rad 2 / m 2 ) at the locations of interest (ad) given in Figure 1 at 1200 on 24 August 1998. The peak in the spectra is at roughly κ = 0 . 04 rad/m or σ = 0 . 63 rad/s which corresponds to an acoustic frequency of ω = 1 . 26 rad/s or f = 0 . 2 Hz. At some locations (bd) there are waves propagating with opposing wavenumber vectors κ , while at other locations (a), most of the waves propagate in the same direction. From Equation (37) we expect these locations with opposing waves to produce the greatest microseismic source levels and shown in Figure 1B.
Remotesensing 10 01624 g002
Figure 3. Geometry of the hurricane wave field and ocean waveguide (not to scale). In our model the waveguide environment may be range-dependent and in this paper we consider the example of upslope propagation from the deep North Atlantic to Florida. The range and ocean depth parameters R and d are given in Table 1. The depth of the receiver below the ocean bottom of 162 m corresponds to the depth below the earth surface of the actual seismometer in Florida. The compression wave speeds c p are 1500 and 5200 m/s in the water and bottom respectively. The shear wave speed c s in the bottom is 3000 m/s. The densities ρ in the water and bottom are 1.0 and 2.5 g/cm respectively.
Figure 3. Geometry of the hurricane wave field and ocean waveguide (not to scale). In our model the waveguide environment may be range-dependent and in this paper we consider the example of upslope propagation from the deep North Atlantic to Florida. The range and ocean depth parameters R and d are given in Table 1. The depth of the receiver below the ocean bottom of 162 m corresponds to the depth below the earth surface of the actual seismometer in Florida. The compression wave speeds c p are 1500 and 5200 m/s in the water and bottom respectively. The shear wave speed c s in the bottom is 3000 m/s. The densities ρ in the water and bottom are 1.0 and 2.5 g/cm respectively.
Remotesensing 10 01624 g003
Figure 4. The ocean depth (solid line) between hurricane Bonnie and the seismometer in Florida at noon on 24 August. Also shown is the ocean depth for the idealized up-sloping environment used to calculate the Green functions in Section 3 (dashed line). The scale of the figure makes the actual slope appear to change rapidly; however, the maximum slope of the ocean floor is roughly 3 .
Figure 4. The ocean depth (solid line) between hurricane Bonnie and the seismometer in Florida at noon on 24 August. Also shown is the ocean depth for the idealized up-sloping environment used to calculate the Green functions in Section 3 (dashed line). The scale of the figure makes the actual slope appear to change rapidly; however, the maximum slope of the ocean floor is roughly 3 .
Remotesensing 10 01624 g004
Figure 5. Horizontal (ad) and vertical (eh) velocity spectra modeled (red dashed) and measured (blue) at the seismometer in Florida at four times on 24 August 1998. The dotted blue lines represent portions of the data that where dominated by non-hurricane related ambient noise. Note that the peak in the spectra is at roughly f = 0 . 2 Hz. This frequency corresponds to the peak in the wave height spectrum at κ = 0 . 04 rad/m seen in Figure 2.
Figure 5. Horizontal (ad) and vertical (eh) velocity spectra modeled (red dashed) and measured (blue) at the seismometer in Florida at four times on 24 August 1998. The dotted blue lines represent portions of the data that where dominated by non-hurricane related ambient noise. Note that the peak in the spectra is at roughly f = 0 . 2 Hz. This frequency corresponds to the peak in the wave height spectrum at κ = 0 . 04 rad/m seen in Figure 2.
Remotesensing 10 01624 g005
Table 1. Parameters for Hurricane Bonnie on 24 August 1998.
Table 1. Parameters for Hurricane Bonnie on 24 August 1998.
TimeHurricane Center PositionRange fromOcean Depth at
Lat ( N)Lon ( W)Sensor (km)Hurricane Center (km)
000024.871.810285.1
060025.272.19835.5
120025.672.49395.5
180026.172.88815.2

Share and Cite

MDPI and ACS Style

Wilson, J.D. Modeling Microseism Generation by Inhomogeneous Ocean Surface Waves in Hurricane Bonnie Using the Non-Linear Wave Equation. Remote Sens. 2018, 10, 1624. https://doi.org/10.3390/rs10101624

AMA Style

Wilson JD. Modeling Microseism Generation by Inhomogeneous Ocean Surface Waves in Hurricane Bonnie Using the Non-Linear Wave Equation. Remote Sensing. 2018; 10(10):1624. https://doi.org/10.3390/rs10101624

Chicago/Turabian Style

Wilson, Joshua D. 2018. "Modeling Microseism Generation by Inhomogeneous Ocean Surface Waves in Hurricane Bonnie Using the Non-Linear Wave Equation" Remote Sensing 10, no. 10: 1624. https://doi.org/10.3390/rs10101624

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop