Modeling acoustic propagation conditions is an important issue in underwater acoustics and there exist several mathematical/numerical models based on different approaches. Some of the most used approaches are based on ray theory, modal expansion and wave number integration techniques. Ray acoustics and ray tracing techniques are the most intuitive and often the simplest means for modeling sound propagation in the sea. Ray acoustics is based on the assumption that sound propagates along rays that are normal to wave fronts, the surfaces of constant phase of the acoustic waves. When generated from a point source in a medium with constant sound speed, the wave fronts form surfaces that are concentric circles, and the sound follows straight line paths that radiate out from the sound source. If the speed of sound is not constant, the rays follow curved paths rather than straight ones. The computational technique known as ray tracing is a method used to calculate the trajectories of the ray paths of sound from the source.
Ray theory is derived from the wave equation when some simplifying assumptions are introduced and the method is essentially a high-frequency approximation. The method is sufficiently accurate for applications involving echo sounders, sonar, and communications systems for short and medium short distances. These devices normally use frequencies that satisfy the high frequency conditions. This article demonstrates that ray theory also can be successfully applied for much lower frequencies approaching the regime of seismic frequencies.
This article presents classical ray theory and demonstrates that ray theory gives a valuable insight and physical picture of how sound propagates in inhomogeneous media. However, ray theory has limitations and may not be valid for precise predictions of sound levels, especially in situations where refraction effects and focusing of sound are important. There exist corrective measures that can be used to improve classical ray theory, but these are not discussed in detail here. Recommended alternative readings include the books. [1-4] and the articles [5-6].
A number of realistic examples and cases are presented with the objective to describe some of the most important aspects of sound propagation in the oceans. This includes the effects of geographical and oceanographic seasonal changes and how the geoacoustic properties of the sea bottom may limit the propagation ranges, especially at low frequencies. The examples are based on experience from modeling sonar systems, underwater acoustic communication links and propagation of low frequency noise in the oceans. There exist a number of ray trace models, some are tuned to specific applications, and others are more general. In this chapter the applications and use of ray theory are illustrated by using Plane Ray, a ray tracing program developed by the author, for modeling underwater acoustic propagation with moderately range-varying bathymetry over layered bottom with a thin fluid sedimentary layer over a solid half with arbitrary geo-acoustic properties. However, the discussion is quite general and does not depend on the actual implementation of the theory.
2. Theory of ray acoustics
Figure 1 shows a small segment of a ray path and the coordinate system. The segment has horizontal and vertical components
Figure 1 shows that the radius of curvature is
When the sound speed varies with depth the ray angle
The ray parameter
At any point in space, the ray curvature is therefore given by the ray parameter
A ray with horizontal angle
The ray parameter is not constant when the bathymetry varies with range. The change in ray direction is illustrated in Figure 2 showing that after reflection the angle
Consequently, after the ray is reflected, its ray parameter must change from
The coordinates of a ray, starting with the angle
The travel time between the two points is obtained by integrating the quantity 1/
The acoustic intensity of a ray can, according to ray theory, be calculated using the principle that the power within a ray tube remains constant within that ray tube. This is illustrated in Figure 4 showing two rays with a vertical angle separation of
At horizontal distance
Instead of using Eq.(11), it may be more convenient to use the vertical horizontal ray
The last expression in Eq.(13) is obtained by assuming that the ray parameter is constant and by using Snell’s law. The absolute values are introduced to avoid problems with regard to the signs of the derivatives and of sin
With respect to the reference distance
In this treatment the transmission loss includes only the geometric spreading loss. Therefore bottom and surface reflection losses and sea water absorption loss must be included separately.
The geometric transmission loss in Eq.(15) consists of two parts. The first term represents the horizontal spreading of the ray tube and results in a cylindrical spreading loss. The second and third terms represent the vertical spreading of the ray tube and are influenced by the depth gradient of the sound speed.
Eq.(13) predicts infinite intensity under either of two conditions: when
3. A recipe for tracing of rays
A simple receipt for a ray tracing algorithm is to divide the whole water column into a large number of layers, each with the same thickness Δ
These equations give the trajectories and the travel times for any ray’s path to the desired range. By applying Eqs. (13) and (14), the geometrical transmission loss is also determined. The simplicity of this method lies in the approximation of the sound speed profiles with straight-line segments and the ray path’s subsequent decomposition into circular segments. The method’s accuracy is determined by how well the linear fit matches the actual profile. In practice, the sound speed profile is often given as measured sound speeds at relatively few depth points. It is therefore advisable to use an interpolation scheme that is consistent with the usual behavior of the sound speed profile to increase the number of depth points to an acceptable high density.
The examples in this article are generated using the ray trace program PlaneRay that has been developed by the author [7-8]. However, any other ray programs with similar capabilities could have been use and the discussion is therefore valid for ray modeling in general. Other models frequently used and are the Bellhop model , and the models [10-11].
Figure 5 shows an example of ray modeling. The sound speed profile is shown at the left panel and the rays from a source at 50 m depth is shown in the right panel, which also shows the bathymetry and the thickness of the sediment layer over the solid half space.
4. Eigenray determination
To calculate the acoustic field it is necessary to have an efficient and accurate algorithm for determination of eigenrays. An eigenray is defined as a ray that connects a source position with a receiver position. In most case with multipath propagation there are many eigenrays for a given source/receiver configuration, which means that finding all eigenrays is not a trivial task.
The PlaneRay model uses a unique sorting and interpolation routine for efficient determination of a large number of eigenrays in range dependent environments. This approach is described by the two plots in Figure 6, which displays the ray history as function of initial angle at the source. All facts and features of the acoustic fields such as the transmission loss, transfer function and time responses are derived from the ray traces and their history The two plots show the ranges and travel times to where the rays cross the receiver depth line (marked by the red dashed line in Figure 5). A particular ray may intersect the receiver depth line, at several ranges. For instance at the range of 2 km, there are 11 eigenrays and from Figure 6 the initial angles of these rays are approximately found to be 5.9°, 9.6°, 22°, 24° for the positive (down going) rays and−2.0°,−3.6, °− 7-4°− 15.0°− 17.0°− 25.0°,−27.0°, for the negative (up-going waves). However, the values found in this way are often not sufficiently accurate for the determination of the sound field. Further processing may therefore be required to obtain accurate results.
The graphs of Figure 6 are composed of independent points, but it is evident that the points are clustered in independent clusters or groups. This property is used for sorting the points into branches of curves that represents different ray history. These branches are in most case relatively continuous and therefore amenable to interpolation. An additional advantage of this method is that the contribution of the various multipath arrivals can be evaluated separately, thereby enabling the user to study the structure of the field in detail.
In most cases the eigenrays are determined by one simple interpolation yields values that are sufficiently accurate for most application, but the accuracy increases with increasing density of the initial angles at the cost of longer computation times.
Figure 7 shows examples of eigenrays traces with rays a receiver located at 2.5 km from the source for the scenario shown in Figure 5. To this receiver there are a total of 12 eigenrays, spanning the range of initial angles from -30° to 29°.
5. Acoustic absorption in sea water
Sound absorption is important for long range propagation especially at higher frequencies. The absorption increases with frequencies and is dependent on temperature, salinity, depth and the pH value of the water. There exists several expressions for acoustic absorption in sea water; one of the preferred options is the semi-empirical formulae by Francoise and Garrison . Figure 8 shows sound absorption as function of frequency in sea water using this expression for the values given in the figure caption.
6. Boundary conditions at the surface and bottom interfaces
Ray tracing is greatly simplified when no rays are traced into the bottom, but stops at the water-bottom interface. This avoids tracing of multiple reflections in layered bottoms. Instead the boundary conditions at the sea surface and the bottom can be approximately satisfied by the use of plane wave reflection coefficient.
A simple and useful bottom model is assuming a fluid sedimentary layer over a homogeneous solid half space. The reflection coefficient of a bottom with this structure is
In Eq (15) and (16) Z
Figure 9 shows an example of the bottom reflection loss as function of angle and frequency for a bottom with a sediment layer with the thickness
The reflection coefficient of a flat even sea surface is −1 for. For a sea surface with ocean waves there will be diffuse scattering to all other direction than the specular direction, which result in a reflection loss that in the first approximation can be modeled by the coherent rough surface reflection coefficient
In this expression
The reflection loss associated with reflection from a rough sea surface is
The same rough surface reflection coefficient may also be applied to a rough bottom.
Figure 10 shows the rough surface reflection loss as function of grazing angle, calculated for a wave height of 0.5 m and the frequencies of 50 Hz, 100 Hz, 200 Hz and 400 Hz.
7. Synthesizing the frequency domain transfer function and the time responses
The total wave field at any receiving point is calculated in the frequency domain by coherent summation of all the eigenray contributions. The first step in the calculation is to determine the geometrical transmission loss of each of the multipath contributions by applying Eq. (13) and Eq.(14) to the sorted and interpolated range-angle values. The frequency domain transfer function and the transmission loss are obtained by adding the multipath contributions coherently in frequency domain taken into account the phase shifts associated the travel times from the interpolated history of the travel times. The frequency dependent acoustic absorption of sound in water is included at this point in the process. The transfer function
Eq. (26) expresses the transfer function
The synthesis of the received signals is performed in the frequency domain by multiplying the frequency spectrum of the source signal with the transfer function of each of the eigenrays and summing the contributions. The time domain response is obtained after multiplication with the frequency function of a source signal followed an inverse Fourier transform of the product. This requires the choice of a source signal, a sampling frequency (
The total duration of the time window (
It is important to select the values of
Figure 11 shows an example where the transmission loss (in dB) as function of range has been calculated for the frequencies of 100 Hz and 200 Hz. The dashed black line indicates the geometrical spreading loss, which is added for comparison and given by,
This expression yields a transmission loss proportional to 20log(
Figure 12 shows the synthesized time response at receivers spaced at 200 m separations in range up to 6 km. The sound speed and bathymetry is the same as in Figure 5 with the source at 150 m and all receivers at 50 m depth. The time scale is in reduced time to remove the gross transmission delay between the source and receiver. The reduced time is defined as
In Eq. (29),
In the example shown above, the time signal and calculated assuming a narrow band-limited source signal in the form of a Ricker pulse. An example of a Ricker pulse and its frequency spectrum are shown in Figure 13.
The time responses in Figure 12 are sorted according to the history of their eigenrays and color coded to allow for studying the various multipath contributions. This is particularly useful when dealing with transient signal and broad band signal, especially when knowledge of the multipath structure is important. In many such situations only the direct arrival or the refracted arrivals in the water column may carry the useful signals and all the other arrivals represent interference. In this case there are direct arrivals, followed by surface reflected and refracted arrivals at the turning points. Notice the high sound pressure values caused by the caustics at 3 km, 6, km and 7 km, which are apparent in both plots, this issue is discussed in section 8.2.
The red dotted line in Figure 12 represents an estimate of the duration of the cannel impulse response. This time duration is mainly given by the bottom reflection coefficient and the critical angle. Rays that propagate at angles closer to the horizontal plane than the critical angle experience almost no bottom reflection loss and may therefore propagate to long distances. Rays with steeper angles will experience higher reflection losses and die out more rapidly with range. Thus the time duration of the impulse is directly determined by the ratio of sound speeds in the water and the bottom as
This estimate of the time duration of the channel impulse response assumes that the bottom is fluid, homogenous and flat, but the estimate may also be useful in other cases with moderately range dependent depth and with solid or layered bottom.
8. Special considerations
8.1. Frequency of applications
Ray tracing is a high frequency approximation to the solution of the wave equation and in principle more valid for high than for low frequency applications. However, high resolution prediction of higher frequency acoustic fields is difficult both for numerical and physical reasons. Principally most important is the physical limitation caused by the fact that the sound speed and the environment are generally not known in sufficient detail. This can be illustrated by a simple example. Consider coherent communication using a frequency of 10 kHz with wavelength of 10 mm. The required accuracy in order to be correct at a distance of 1 km is that the sound speed is known and stable with a relative error less 10-5, an impossible requirement to satisfy in practice regardless of the numerical accuracy of the computer model.
8.2. Caustics and turning points
As mentioned before, the locations where
Figure 14 shows details of the field at a showing the rays with initial angles in the range of −6° to −1°. The scenario is the same as in of Figure 5, but for clarity the tracing of rays have been stopped after the first bottom reflection and the figure concentrates on the details the field at the caustic at 1760 m range for a ray with initial angle of −5.6°. Figure 15 shows the time responses for ranges in the interval from 1.6 km to 1.9 km. In this case, the source signal is a Ricker pulse with a peak frequency of 200 Hz. There is a first direct arrival (black color) at all ranges. From the range 1760 there is also a refracted arrival a little later than the direct, but with higher amplitude, in particular near the range of 1760 m. Notice the effect of the 90° phase shift for ranges beyond the caustic at 1760 m and that the amplitude at this range is considerable higher than at the other ranges.
8.3. The principle of reciprocity and its validity in ray modeling
The principle of reciprocity is an important and useful property of linear acoustics and systems theory. The principle is very general and valid also in cases where the wave undergoes reflections at boundaries on its path from source to receiver . The reciprocity principle is correctly represented in ray modeling, as can easily be understood from the eigenray plots of Figure 7. The eigenrays from a source position to the receiver position are the same as when source and receiver changes positions. The reflections at the bottom and at the sea surface are also symmetric in angles and consequently the acoustic fields are the same. However, it should be noted that the reciprocity principle applies to a point-to-point situation. This means that, for instance, that the development of the transmission loss as function source-receiver separation is generally not the same for the two directions.
8.4. The validity of using plane wave reflection coefficients
The accuracy of any ray model depends on the validity and limitation of ray theory and the implementation. A fundamental assumption of model is that the interactions with the boundaries are adequately described by plane wave reflection coefficient. In this section the validity of this assumption is investigated.
The general expression for the reflected field is given in text books, for instance in , over horizontal wave numbers
Consider now the situation where is constant and independent of
According to Eq. (33) the reflected wave is the same as the outgoing spherical wave from the image of the source in the mirror position of the real source and modified by the constant reflection coefficient. The situation with a constant reflection coefficient is valid for perfectly flat sea surface where the reflection coefficient is equal to -1 for all angles of incidence. Thus the reflection from a smooth sea surface is accurately described plane wave reflection coefficients.
In the general case, and for reflections from the bottom, the reflection coefficient is not constant and the integral can only be solved approximately or numerically. In order to obtain an approximation of the integral in Eq.(31), the Hankel function is expanded in a power series with the first terms being
Restricting the integral of Eq.(31) to the first term yields
The exponential in the integrand will normally be a rapid varying function and therefore the value of the integral will be small except when the phase term of Eq.(36) is nearly constant. The phase term of Eq.(36) is
The stationary points are defined to the values of the horizontal wave number
The interpretation of this result is quite simple; the reflected wave field is equal to that of the image source multiplied with the reflection coefficient at the specular angle
There are however situations where this approximation is not sufficient in practice. This is discussed in  and in the following their results are cited without proof. The accuracy of the approximation depends on the source or receiver distance from the bottom interface. The result of the analysis is that the distance
With the water parameters of
8.5. Bench marking ray modeling
The wave number integration model OASES  has been used to validate the accuracy and the limitation of the ray trace model using the simple case with constant water depth of 100 m and constant sound speed of 1500 m/s.
Figure 16 show the calculated transmission loss for the frequencies of 25 Hz, 50 Hz, 100 Hz and 200 Hz. The agreements between the results are very god for the higher frequency, but with some discrepancies for the lower frequencies, in particular for 25 Hz. The discrepancy is mainly a phase shift in the interference patterns of the two results, most pronounced for low frequencies and long ranges. This observation agrees with the theory outlined earlier. The seriousness of this discrepancy or errors may not very important in practice since the mean level is nearly the same as shown by the comparison with the OASES model.
9. Case studies
In the following we present two case studies that are relevant application of the modeling techniques descried in this article. The first if these is in connection with acoustic underwater communication and the transmission of digital information. In this case the multipath communication may be a significant problem causing intersymbol interference and significant degradation of reliability and performance. The second case is related to studies on the propagation of low frequency sound and the effect such noise may affect marine life, sea mammals and fish.
9.1. Seasonal variations of communication links
In connection with a study of underwater acoustic communication the propagation over a 6 km track has been modeled for the various seasonal sound speed profiles.
The sound speed some months are shown in Figure 17. The sound speed profiles depend on the sea water temperature, the salinity and the depth. In the present case the sea water temperature variation with depth and the seasons is the main reason for changes in sound speed profile. During winters the surface water is cold and the sound speed is low, in the summer the surface water temperature and the sound speed is higher. The seasonal heating and cooling of the surface water propagates also to deeper depths, but with diminishing temperatures changes. At very large depths the water temperature is nearly the same at all seasons and the sound speed increases linearly and slowly with depth.
Figure 18 shows ray tracing results are for the same profiles as displayed in Figure 17. The purpose of the study was to investigate the possibility of communication to positions beyond a sea mount and to study the multipath arrival structure as function of range and depth.
There is a seamount with a peak at about 3 km from the transmitting station. In order to simplify the interpretation ray tracings in these plots have been terminated after 6 bottom reflections, but all rays are included in the calculation of the acoustic field, but rays with so many bottom reflections, or more, will in most case not be useful for data communication because of the reflection loss and reduced coherence.
Figure 19 shows examples of received time responses at 25 m depth using a Ricker pulse as source signal. The different multipath contributions are color coded for clarity. At distances from the source over 1.5 km the first arrivals is follow paths surface reflected and upward refracted paths
Figure 20 shows the channel responses at a fixed range as function of depth down to 50 m. This figure shows the total response after adding all the individual multi path contributions. The plots demonstrate that the surface channel consists of deep refracted path and a number of paths reflected from the surface and deeper upwards refractions. The stability of these paths may be uncertain and subject to rapid changes in the environmental conditions near the surface due to temperature wind and current.
9.2. Seismic noise propagation
In many areas of the world anthropogenic noise often dominates over the natural ambient noise, especially in the low frequency band from approximately from 10 Hz and upwards to 1000 Hz, or more. This frequency band coincides approximately with the frequencies of perception of sea mammals and fish and may therefore be harmful to their natural activities, or even cause physical damages. An example is the case of the seismic exploration for oil and gas in certain areas where there is important commercial fishing interest. The propagation and distribution of acoustic noise depends the environmental conditions, in particular the oceanographic parameters, the topography of the seafloor and the acoustic properties of the bottom. In this section some of examples are presented to illustrate how the environment may affect the distribution of sound and noise. This study and discussion is also relevant for passive sonar applications to detect and track submerged vehicles and objects base emitted acoustic noise
The effects of bathymetric are illustrated in Figure 21 showing ray traces of upslope and downslope conditions for typical summer conditions at the Halten Bank in the Norwegian Sea. With downslope propagation there is a thinning the ray density with distance and upslope propagation gives a concentration of rays as the water depth diminishes.
Figure 22 and Figure 23 show the calculated sound pressure level as function of range for the downslope and upslope propagation. The sound pulse from an airgun array is modeled as s a Ricker pulse with a peak pressure of 260 dB rel. 1µPa, centered on the frequency of 50 Hz, The horizontal dashed line is the assumed threshold value for fish reaction to sound. The bottom is modeled with a 2 m thick sedimentary layer over solid rock. The sound speed in the sediment layer is 1700 m/s and the density is 1800 kg/m3. The compressional sound speed in the rock is 3000 m/s, and density is 2500 kg/m3. The results in Figure 22 and Figure 23 are obtained under two conditions: (a) with a shear speed of 500 m/s, and (b) with no shear wave in the rock, i.e. the shear speed is zero. The absorptions are assumed to be 0.5 dB per wavelength for all the waves in the sediment layer and the rock. In the first case (a) the bottom reflection loss is as shown in Figure 9 with a significant low frequency reflection loss at angles lower than the critical angle caused by absorptions and conversion to shear wave in the bottom, which draws energy for the reflected wave. In the case of Figure 22 this results in a low-frequency and low-angle reflection loss of about 1 dB. For long ranges and many reflections this adds up to a significant total propagation loss. With no shear conversion the reflection loss is considerably reduced and the sound propagates easier to long ranges. The difference between the sound level at 50 Hz and 100 Hz is partly a result of increase attenuation at the higher frequency and partly that the source level in this case is higher for 50 Hz than for 100 Hz.
Figure 24 and Figure 23 show similar results for downslope and upslope propagation for typical winter conditions represented by a sound speed profile measured in the month of February. For downslope conditions the sound level decrease rapidly with increasing depth and much more rapidly with shear wave conversion (Figure 24a) than without shear (Figure 24b). With upslope propagation (Figure 23) the sound levels are near independent of shear conversion except at the very long rages where the water depth becomes constant. The examples demonstrate that sound propagation in the ocean is strongly influenced by both by the oceanographic conditions and the geophysical properties of the bottom. Reliable prediction of acoustic propagation condition requires modeling tool that can that can handle both bottom and water properties.
The article has outlined the theory of ray modeling and described how the theory can be applied to study acoustic wave propagation in the ocean. The complete acoustic fields are calculated by coherent addition of the contributions of a large number of eigenrays. In this method no rays are traced into the bottom, but the bottom interaction is modeled by plane wave reflection coefficients. Ray tracing is, by definition, frequency independent and therefore the ray trajectories through the water column are valid for all frequencies. Frequency dependency is introduced by reflections from the sea surface and the bottom, including loss associated with absorption and diffuse scattering of a rough ocean and bottom interfaces. Ray tracing is therefore a computational effective method for modeling broad of frequency band wave fields and for calculation of time responses.
Ray tracing is high-frequency approximation to the solution of the wave equation and the accuracy and validity at lower frequencies may be questioned, in particular the use of plane ray reflection coefficient to represent the bottom effects. This problem has been considered both theoretically and by simulations and comparison with more accurate model. The results of this study shows that source and receiver should be at a height above the bottom of at least half a wavelength, but there is no similar requirement to the distance from the sea surface. Less fundamental is the limitation of the numerical accuracy of the determination of the eigenrays, which is most serious in the calculation of the ray amplitude and the transmission loss. These inaccuracies are of more practical nature and can be reduced by refinements in the calculations.
Examples relevant for application in acoustic underwater communication and active sonar have been presented. The propagation of low frequency sound to large distances has been presented showing the effect of the bathymetry and the acoustic properties of the bottom. An important conclusion is the effect of bathymetry and the sound speed structure interacts and that accurate modeling of sound propagation requires information about the oceanography, the bathymetry and the geology of the bottom.