Acoustic Sensors for Air and Surface Navigation Applications

This paper presents the state-of-the-art and reviews the state-of-research of acoustic sensors used for a variety of navigation and guidance applications on air and surface vehicles. In particular, this paper focuses on echolocation, which is widely utilized in nature by certain mammals (e.g., cetaceans and bats). Although acoustic sensors have been extensively adopted in various engineering applications, their use in navigation and guidance systems is yet to be fully exploited. This technology has clear potential for applications in air and surface navigation/guidance for intelligent transport systems (ITS), especially considering air and surface operations indoors and in other environments where satellite positioning is not available. Propagation of sound in the atmosphere is discussed in detail, with all potential attenuation sources taken into account. The errors introduced in echolocation measurements due to Doppler, multipath and atmospheric effects are discussed, and an uncertainty analysis method is presented for ranging error budget prediction in acoustic navigation applications. Considering the design challenges associated with monostatic and multi-static sensor implementations and looking at the performance predictions for different possible configurations, acoustic sensors show clear promises in navigation, proximity sensing, as well as obstacle detection and tracking. The integration of acoustic sensors in multi-sensor navigation systems is also considered towards the end of the paper and a low Size, Weight and Power, and Cost (SWaP-C) sensor integration architecture is presented for possible introduction in air and surface navigation systems.


Introduction
The directionality of acoustic waves has been long used for localization by human beings. The term 'echolocation' was coined by Donald R. Griffin [1], where he discusses ship captains exploiting sound to ascertain the ship's surroundings and avoid obstacles in low visibility environments. Acoustic sensors provide a low Size, Weight, and Power (SWaP) solution, which is low cost, scalable and robust. Moreover, acoustic sensors have the capability to provide high-resolution spatial information at short distance range. Radio-based localization techniques like Global Navigation Satellite Systems (GNSS) are prone to data degradations in dense urban environments and indoors [2]. On the other hand, electromagnetic techniques suffer from interference from other sources as well as metal structures. Optical navigation sensors are also still relatively expensive and their performance deteriorates

•
Bats can lower their call intensity as they approach strong reflective objects to prevent the echo sound pressure level from becoming too large. Bats can exhibit very high-resolution of target detection with time difference discriminations of 10-12 nanoseconds [3,6,7]. The duration of echolocation can vary considerably, with individual clicks being approximately~50-100 µs long to constant frequency signals which are longer than 30 ms. Table 1 lists various bat species and their call type, based on their diet. • As seen in the Table 1, the echolocation call can consist of a single frequency or multiple frequencies comprising a harmonic series. The pulse interval of the call also varies with proximity to the target. As bats approach their target, the repetition rate of their calls increases to get faster localisation updates. Also, the pulse interval of the call gives an indication of the maximum range from which bats can detect objects. Gleaning bats can passively listen to prey-generated sounds to localize their prey by interrupting echolocation or drastically reducing call intensity shortly before capturing prey [4]. • Big brown bats (Eptesicus fuscus) might roughly localise the position of its prey by listening to conspecific-generated echoes [5]. This is true for bottlenose dolphins (Tursiops truncates) as well [8].
Certain bats like Mexican free-tailed bat tend to increase their emission rate when flying in pairs. However, when flying in bigger groups, bats tend to decrease their emission rates, thereby reducing mutual interference [9]. This temporal modulation of emission in bats is similar to timing algorithms employed in electronic communication systems. These algorithms, also referred to as back-off algorithms, introduce probabilistic delays in resending packets lost due to interference [10].  [11,12]. Bats exhibit a variety of behaviour while coping with interference from environmental noise as well as both calls and echoes from nearby bats. While some species treat the presence of nearby conspecifics as any other source of noise or object in their field of view [13], certain species of bats like the free-tailed bats (Molossidae) compensate for interference by calling louder or by varying the frequency or duration of echolocation pulses [14].
Echolocation in bats refers to sensing of the environment based on Time of Arrival (TOA) of sound waves emitted by them. This helps bats navigate as well track their prey during the night. The strength of the received signal is indicative of the size of the target. Also, analysis of frequency spectrum of the echo gives an idea of the surface texture of the target. Most bat echolocation calls are ultrasonic, ranging from 20-200 kHz and the sound intensity can reach up to 130 dB. There have been attempts to investigate echolocation abilities in human beings as well, especially in visually impaired. In [15], the echolocation abilities of blind and sighted humans are reviewed, suggesting enhanced auditory abilities in visually impaired than normally sighted humans. The effect of prior visual experience on sound localisation in late blind individuals has been studied in detail in [16].
There has been some research focus on analyzing the flight dynamics of bats and attempts made to emulate the same [17]. Bats use either their tongue or vocal chords to produce sonar signals [8]. Bats can vary the frequency, signal duration, signal intensity, harmonic composition and pulse interval according to their surroundings. Bats use narrowband signals for ranging of distant targets and broadband for localization. Some species of bats also account for Doppler shift by varying their call frequency [3,18]. Attempts have been made to develop biomimetic sonars inspired by bats' external ears or pinnae (Figure 1), for localization and mapping, referred to as BatSLAM [19].

Sound Propagation
Acoustic waves are longitudinal waves that require a material medium to propagate. Fundamentally, sound can be defined as mechanical energy transmitted by pressure waves in a material medium. Acoustic waves are mechanical waves, i.e., they involve rapid to and fro displacements or vibrations of molecules in the medium. The velocity of sound in a medium (c m ) varies with the bulk modulus (B) and density (ρ) of the medium as shown in Equation (1). Sound travels faster in a medium with high bulk modulus or stiffness, like solids as compared to a medium with lower bulk modulus, like fluids: The attenuation rate of sound waves varies with frequency, with higher frequencies attenuating at a faster rate. Attenuation can occur either due to reflection/scattering at interfaces or absorption [21]. However, higher frequencies, having short wavelengths, reflect strongly from small objects. Reflection from surfaces causes interference with the incident sound wave, which could be constructive or destructive. Interference depends upon the frequency of sound as well as the difference between the path length of direct and reflection paths [22]. Furthermore, the speed of sound in air varies with temperature, pressure, humidity, and wind, thereby affecting the propagation of sound. The generic equation for sound propagation can be given by: where L p (r): The sound pressure level at distance r from the source (dB); L w : The sound power level of the source (dB); A i : The combination of modifying factors that either attenuate or enhance the transmission of the sound energy as it propagates from source to receiver. Acoustic sources have both far-field and near-field regions. Wavefronts produced by the sound source in near-field are not parallel and the intensity of the wave oscillates with range and the angle between source and receiver. However, in the far-field, wavefronts are nearly parallel, with intensity varying only with range to a centroid between sound sources, in accordance with the inverse squared rule. The near-field distance r n f is given by: where D is the equivalent aperture of the transmitter given by: where k is the wave number and θ 3dB is the half power beam angle. The wavefront for a sound source radiating equally in all directions is a sphere of radius r, whose intensity I from the source of power W is given by: Assuming a point source of sound in an unbounded homogenous atmosphere, the propagation of sound is affected by just two attenuating effects. While the first attenuation effect is geometric, which is solely dependent on the distance from the sound source, the second attenuating effect is the atmospheric absorption. Sound propagates due to the oscillation of air molecules about their mean position; with a higher frequency of sound leading to a higher rate of oscillation. This vibration of the air molecules leads to loss of energy through two dissipative mechanisms. While one of the mechanisms comprises of frictional losses, which includes both viscous action and heat conduction, the other mechanism involves the interaction of water vapour with the resonance of oxygen and nitrogen molecules. Hence, there are heat conduction losses, shear viscosity losses, and molecular relaxation losses [23].

Sound Attenuation in Atmosphere
However, in practical situations, the propagation of sound in the atmosphere is affected by additional factors like ground effects, attenuation due to finite barriers and buildings, reflections, wind, and temperature gradient effects, and atmospheric turbulence. The atmospheric sound attenuation factors are discussed in detail in Sections 3.1.1-3.1.5.
Geometrical divergence refers to the spherical spreading in the free field from a point sound source. The attenuating effect is the geometric attenuation, which results from the spreading of the radiated sound energy over a sphere of increasing diameter as the wavefront propagates away from the source. Equation (6) shows the relationship between the sound power level of the source, L w , and sound pressure level, L p (r), at a distance r from that source. The variation of sound pressure level with distance from the source is shown in Equation (7). Unlike atmospheric attenuation, geometric attenuation is independent of the frequency of the propagating sound wave. From Equation (7), it can be inferred that sound intensity or sound pressure level, L p decreases by 6 dB per doubling of distance away from the source: The geometrical divergence, in dB, is given by: where d is the distance from the sound source to receiver (m) and d 0 is the reference distance which is 1 m from an omnidirectional point sound source.

Atmospheric Absorption (A atm )
Air absorption becomes significant at higher frequencies and at long ranges, thereby acting as a low-pass filter at long range. The pressure of a planar sound wave at a distance x from a point of pressure P 0 is given by: The attenuation coefficient, α, for air absorption depends on frequency, humidity, temperature, and pressure, with its value being calculated using Equations (10)-(12) [24].
where f is the frequency, T is the absolute temperature of the atmosphere in kelvins, T 0 is the reference value of T (293.15 K) and f r,N and f r,O are relaxation frequencies associated with the vibration of nitrogen and oxygen molecules, respectively, and are given by: where p s is local atmospheric pressure, p 0 is the reference atmospheric pressure (101,325 Pa) and H is the percentage molar concentration of water vapour in the atmosphere which is given by: where r h is the relative humidity and ρ sat is given by: where C sat is given by: Similarly, ρ sat can also be written as [25,26]: (16) Figure 2 shows the variation of absorption coefficient α with the frequency of sound at 293.15 K, one atmospheric pressure, 20% relative humidity and H being 4.7 × 10 −3 . As there are two relaxation frequencies associated with oxygen and nitrogen, the frequency dependence of the attenuation coefficient for sound in the air has three distinct regions. At very low frequencies, where the sound frequency is much lower than that associated with nitrogen molecules, the attenuation is dominated by vibrational relaxation of nitrogen molecules (α 1 ). The frequency dependence is quadratic with an apparent bulk viscosity associated with the nitrogen relaxation. In the intermediate region, the frequency is substantially larger than that associated with nitrogen relaxation, but still substantially less than that associated with oxygen relaxation, with quadratic frequency dependence, smaller coefficient and apparent bulk viscosity that is associated with oxygen relaxation (α 2 ). In the higher frequency region, there is quadratic dependence again, although with an even smaller coefficient, and with the intrinsic bulk viscosity associated with molecular rotation (α 3 ) [24].
where is the relative humidity and is given by: where is given by: Similarly, can also be written as [25,26]: . .
(16) Figure 2 shows the variation of absorption coefficient with the frequency of sound at 293.15 K, one atmospheric pressure, 20% relative humidity and being 4.7 × 10 −3 . As there are two relaxation frequencies associated with oxygen and nitrogen, the frequency dependence of the attenuation coefficient for sound in the air has three distinct regions. At very low frequencies, where the sound frequency is much lower than that associated with nitrogen molecules, the attenuation is dominated by vibrational relaxation of nitrogen molecules ( ). The frequency dependence is quadratic with an apparent bulk viscosity associated with the nitrogen relaxation. In the intermediate region, the frequency is substantially larger than that associated with nitrogen relaxation, but still substantially less than that associated with oxygen relaxation, with quadratic frequency dependence, smaller coefficient and apparent bulk viscosity that is associated with oxygen relaxation ( ). In the higher frequency region, there is quadratic dependence again, although with an even smaller coefficient, and with the intrinsic bulk viscosity associated with molecular rotation ( ) [24]. Having calculated the absorption coefficient for a given temperature, pressure, relative humidity and percentage molar concentration of water vapor in the atmosphere, the attenuation of sound due to atmospheric absorption, during propagation through a distance (m) is given by:  Having calculated the absorption coefficient for a given temperature, pressure, relative humidity and percentage molar concentration of water vapor in the atmosphere, the attenuation of sound due to atmospheric absorption, during propagation through a distance d (m) is given by: Adding the effect of a bounding ground plane to the sound propagation model allows for sound to propagate directly from source to receiver as well as through secondary propagation path resulting from a reflection off the ground plane as shown in Figure 3. This secondary propagation path can result in interference effects between the direct and reflected waves at the receiver. The interference effect can be constructive or destructive, depending on the relative amplitudes and phase of the direct and reflected waves. The relationship between the direct and reflected waves depends on a variety of factors including the difference between the direct and reflected path lengths, which is a function of source and receiver separation distance (dp) as well as their height above the ground (h s and h r ), the wavelength of sound and reflective properties of the ground which can cause variations in the phase and amplitude of the reflected sound wave. Adding the effect of a bounding ground plane to the sound propagation model allows for sound to propagate directly from source to receiver as well as through secondary propagation path resulting from a reflection off the ground plane as shown in Figure 3. This secondary propagation path can result in interference effects between the direct and reflected waves at the receiver. The interference effect can be constructive or destructive, depending on the relative amplitudes and phase of the direct and reflected waves. The relationship between the direct and reflected waves depends on a variety of factors including the difference between the direct and reflected path lengths, which is a function of source and receiver separation distance (dp) as well as their height above the ground (ℎ and ℎ ), the wavelength of sound and reflective properties of the ground which can cause variations in the phase and amplitude of the reflected sound wave. Based on the acoustical properties of ground surfaces, they are classified into three types based on the ground factor ( ). Hard ground, which has a of 0, includes concrete, paving, water, ice and all other low porosity ground surfaces. On the other hand, porous ground, which has a of 1, consists of grass, trees, foliage and all other ground surfaces which are suitable for the growth of vegetation. Surfaces which are a combination of both hard and porous grounds, i.e., having a value of ranging from 0 to 1, where represents the fraction of the ground surface that is porous, are known as mixed ground. The total ground attenuation is obtained by summing up the attenuation for the source region specified by the ground factor , for the middle region specified by the ground factor , and for the receiver region specified by the ground factor , as shown in Equation (18):

Screening ( )
An object shall be considered as a screening obstacle ( Figure 4) if it meets the following requirements:


The object has a surface density of at least 10 kg/m 2 ;  The surface of the object is closed without cracks or gaps;  The horizontal dimension of the object normal to the source-receiver line ( + ) is larger than the acoustic wavelength λ at the nominal midband frequency for the octave band of interest, i.e., + > λ ( Figure 5).
The barrier diffraction could be either single diffraction in case of thin barriers ( Figure 5) or double diffraction in thick barriers. In case of more than two barriers, the barrier attenuation can be approximated to be a case of double diffraction, by choosing the two most effective barriers, Based on the acoustical properties of ground surfaces, they are classified into three types based on the ground factor (G). Hard ground, which has a G of 0, includes concrete, paving, water, ice and all other low porosity ground surfaces. On the other hand, porous ground, which has a G of 1, consists of grass, trees, foliage and all other ground surfaces which are suitable for the growth of vegetation. Surfaces which are a combination of both hard and porous grounds, i.e., having a value of G ranging from 0 to 1, where G represents the fraction of the ground surface that is porous, are known as mixed ground. The total ground attenuation is obtained by summing up the attenuation A s for the source region specified by the ground factor G s , A m for the middle region specified by the ground factor G m , and A r for the receiver region specified by the ground factor G r , as shown in Equation (18):

Screening (A bar )
An object shall be considered as a screening obstacle ( Figure 4) if it meets the following requirements: • The object has a surface density of at least 10 kg/m 2 ; • The surface of the object is closed without cracks or gaps; • The horizontal dimension of the object normal to the source-receiver line (l i + l r ) is larger than the acoustic wavelength λ at the nominal midband frequency for the octave band of interest, i.e., l i + l r > λ ( Figure 5).
The barrier diffraction could be either single diffraction in case of thin barriers ( Figure 5) or double diffraction in thick barriers. In case of more than two barriers, the barrier attenuation can be approximated to be a case of double diffraction, by choosing the two most effective barriers, neglecting the effect of the others. The effect diffraction (in dB) for downward sound propagation over the top edge and the vertical edge, respectively, is given by Equations (19) and (20) [27]. D Z is the barrier attenuation for each octave band and A gr is the ground attenuation in the absence of the barrier, as described in Section 3.1.3.

Wind and Temperature Gradient Effects
As a result of uneven heating of the Earth's surface, the atmosphere is constantly in motion. The turbulent flow of air across the rough solid surface of the Earth generates a boundary layer. The lower part of the meteorological boundary layer, called the surface layer, extends over 50-100 m in typical daytime conditions [24]. Turbulent fluxes vary by less than 10% of their magnitude in the surface layer, but the wind speed and temperature gradients are the largest. Turbulence can be modelled as a series of moving eddies with a distribution of sizes. Various turbulence models like Gaussian, Von Kármán, and Kolmogorov are used in atmospheric acoustics [28][29][30]. It has been shown that turbulence effects decrease with increase in elevation of sound sources from the ground [31].
As the temperature decreases with height, in the absence of wind, this causes sound waves to bend, or refract, upwards. Wind velocity either adds or subtracts from the velocity of sound,

Wind and Temperature Gradient Effects
As a result of uneven heating of the Earth's surface, the atmosphere is constantly in motion. The turbulent flow of air across the rough solid surface of the Earth generates a boundary layer. The lower part of the meteorological boundary layer, called the surface layer, extends over 50-100 m in typical daytime conditions [24]. Turbulent fluxes vary by less than 10% of their magnitude in the surface layer, but the wind speed and temperature gradients are the largest. Turbulence can be modelled as a series of moving eddies with a distribution of sizes. Various turbulence models like Gaussian, Von Kármán, and Kolmogorov are used in atmospheric acoustics [28][29][30]. It has been shown that turbulence effects decrease with increase in elevation of sound sources from the ground [31].
As the temperature decreases with height, in the absence of wind, this causes sound waves to bend, or refract, upwards. Wind velocity either adds or subtracts from the velocity of sound,

Wind and Temperature Gradient Effects
As a result of uneven heating of the Earth's surface, the atmosphere is constantly in motion. The turbulent flow of air across the rough solid surface of the Earth generates a boundary layer. The lower part of the meteorological boundary layer, called the surface layer, extends over 50-100 m in typical daytime conditions [24]. Turbulent fluxes vary by less than 10% of their magnitude in the surface layer, but the wind speed and temperature gradients are the largest. Turbulence can be modelled as a series of moving eddies with a distribution of sizes. Various turbulence models like Gaussian, Von Kármán, and Kolmogorov are used in atmospheric acoustics [28][29][30]. It has been shown that turbulence effects decrease with increase in elevation of sound sources from the ground [31].
As the temperature decreases with height, in the absence of wind, this causes sound waves to bend, or refract, upwards. Wind velocity either adds or subtracts from the velocity of sound, depending upon whether the source is upwind or downwind of the receiver, height above ground and temperature inversions. Wind effects tend to dominate over temperature effects when both are present. The general relationship between the speed of sound profile c(z), the temperature profile T(z) and wind speed profile u(z) in the direction of sound propagation, for a height z, is given by [32]:

Other Sound Attenuation Factors
Various standardised techniques for measuring the attenuation of sound outdoors due to atmospheric absorption effects have been developed in ISO 9613-1 and ISO 9613-2 [27,33]. This includes attenuation of sound due to miscellaneous effects like foliage, housing or industrial sites. The analytical models presented rely on the values of various atmospheric parameters like temperature, pressure, relative humidity, wind speed and time of the day. Besides, fog and precipitation can also affect the attenuation of sound [34]. Experiments show that precipitation affects the temperature variation, hence indirectly affecting sound attenuation outdoors [31]. Hence, to summarize, the attenuation of sound in atmosphere can be given by: where A misc is the sound attenuation due to other miscellaneous effects like wind and temperature gradient effects, precipitation, foliage, and housing or industrial sites. Table 2 gives the different ranging parameters involved in the design of the acoustic sensors. The ranging equation is given by Equation (23), where R m is the measured range, R a is the actual range from the transmitter (x T , y T , z T ) to the receiver (x R , y R , z R ) and ε is the error in the measured range. The error term ε comprises mainly of error due to multipath (ε Mp ), Doppler shift (ε Ds ) and atmospheric effects (ε Atm ):

Echolocation Errors
where:

Doppler Effect
The Doppler effect is caused by the perceived change in sound frequency due to relative motion between sound source and receiver. Figure 6 shows the elevation angle (E n ) of the nth transmitter (Tr n ) to the receiver (R) as well as the relative bearing (χ n ), the tangential velocity of the transmitter ( → v T ) and the azimuth of the Line of Sight (LOS) projection (χ n ).
→ v 0 is the velocity of the receiver and the case of results in a null Doppler shift as no component of receiver velocity vector is in the direction of LOS to the transmitter. As is evident from Equation (26), the Doppler shift is inversely proportional to both elevation and azimuth angle. The pitch increases as the relative distance between the sound source and receiver decreases. The change in observed frequency of sound is given as: where → v n = nth transmitter velocity component along the LOS; → v a = receiver velocity projection along the LOS; c = speed of sound (ms −1 ); f = sound frequency (Hz); E n = elevation angle of the nth transmitter; χ n = azimuth of the LOS projection.   The moving sound source emits crests every ′ units of time. The Doppler-shifted frequency from a moving sound source emitting frequency is given by: Considering a simplified case where transmitter motion and LOS from the transmitter to the receiver are coplanar, the sound field at times t and (t + t ) for a moving sound source (Tr) which has moved to a new point (Tr ) in time t , is shown in Figure 7.  Considering a simplified case where transmitter motion and LOS from the transmitter to the receiver are coplanar, the sound field at times and ( + ′) for a moving sound source ( ) which has moved to a new point ( ') in time ′, is shown in Figure 7. The moving sound source emits crests every ′ units of time. The Doppler-shifted frequency from a moving sound source emitting frequency is given by: The moving sound source emits crests every t units of time. The Doppler-shifted frequency f from a moving sound source emitting frequency f s is given by: where M is the Mach number for the sound source and θ is the direction of the receiver to the sound source at the time (t + t ).
Assuming no relative motion between transmitter-receiver and the speed of the transmitter being Mc, if it takes time t for sound to reach the receiver from the transmitter, the error in range due to Doppler shift is given by: where θ is the direction of receiver motion relative to the LOS between the transmitter and the receiver.

Multipath
An important property of a medium that influences the strength or amplitude of reflected waves is acoustic impedance. Acoustic impedance can be defined as the product of the density of the medium (ρ) and speed of sound (c m ) in the medium: Acoustic impedance gives a measure of the sound transmitted and reflected back at the interface of two mediums. The ratio of the reflected pressure amplitude, P r , to the incident pressure amplitude, P i , called amplitude reflection coefficient, is given by: Assuming a homogenous medium, sound propagation, especially at high frequencies, can be assumed to be a straight line from the source to the receiver [35]. Assuming a range independent geometry for a homogenous medium, the sound waves are subject to multiple reflections, as shown in Figure 8. As most of the reflecting surfaces are irregular, sound waves experience a significant amount of scattering or diffraction on reflection. where is the direction of receiver motion relative to the LOS between the transmitter and the receiver.

Multipath
An important property of a medium that influences the strength or amplitude of reflected waves is acoustic impedance. Acoustic impedance can be defined as the product of the density of the medium ( ) and speed of sound ( ) in the medium: Acoustic impedance gives a measure of the sound transmitted and reflected back at the interface of two mediums. The ratio of the reflected pressure amplitude, , to the incident pressure amplitude, , called amplitude reflection coefficient, is given by: Assuming a homogenous medium, sound propagation, especially at high frequencies, can be assumed to be a straight line from the source to the receiver [35]. Assuming a range independent geometry for a homogenous medium, the sound waves are subject to multiple reflections, as shown in Figure 8. As most of the reflecting surfaces are irregular, sound waves experience a significant amount of scattering or diffraction on reflection. Using ray-tracing [36,37], the reflection point and the defined point , as shown in Figure 9, should satisfy the equation: where sound waves are emitted from point to the receiver location after reflection at point . is defined as a point on the reflecting surface and n is a unit vector normal to that surface. The line equation connecting and is given by: where is a parameter between 0 and 1. Combining Equations (31) and (32): Using ray-tracing [36,37], the reflection point S and the defined point V, as shown in Figure 9, should satisfy the equation: where sound waves are emitted from point Tr to the receiver location R after reflection at point S. V is defined as a point on the reflecting surface and n is a unit vector normal to that surface. The line equation connecting Tr and R image is given by: where m is a parameter between 0 and 1. Combining Equations (31) and (32): Assuming specular reflection, the extra path length is given by: The starting point for performing ray-tracing of acoustic waves is Helmholtz equation, which can be written in Cartesian coordinates = ( , , ) as: where is the total pressure, is the angular frequency of the source located at . A solution of the Helmholtz equation, which is called the ray series, is used to obtain the ray equations [38]. The ranging error due to multipath is given by: where: = location of ith acoustic transmitter; = ith reflection point; = location of ith acoustic receiver; = number of reflecting surfaces.

Atmospheric Effects
In the "International Standard Atmosphere" (ISA), the troposphere extends up to 11 km. The temperature gradient in the troposphere can be assumed to be constant. Following this assumption, the ranging error due to atmospheric effects is given by: where is the variation of the speed of sound due to the wind. The variation of the speed of sound due to temperature [39] is given by: where: = speed of sound at sea-level; = height above sea-level; = is the gradient of the speed of sound, where is the sea-level temperature ( K ) and = is the variation of Assuming specular reflection, the extra path length L mS is given by: The starting point for performing ray-tracing of acoustic waves is Helmholtz equation, which can be written in Cartesian coordinates x = (x, y, z) as: where p is the total pressure, ω is the angular frequency of the source located at x 0 . A solution of the Helmholtz equation, which is called the ray series, is used to obtain the ray equations [38]. The ranging error due to multipath is given by: where: Tr i = location of ith acoustic transmitter; S i = ith reflection point; R i = location of ith acoustic receiver; N = number of reflecting surfaces.

Atmospheric Effects
In the "International Standard Atmosphere" (ISA), the troposphere extends up to 11 km. The temperature gradient in the troposphere can be assumed to be constant. Following this assumption, the ranging error due to atmospheric effects is given by: where c w is the variation of the speed of sound due to the wind. The variation of the speed of sound due to temperature [39] is given by: where: c 0 = speed of sound at sea-level; H = height above sea-level; dc dH = λc 0 2T 0 is the gradient of the speed of sound, where T 0 is the sea-level temperature (K) and λ = dT dH is the variation of temperature with height. The variation of the speed of sound due to wind is given by: where: c r = speed of sound relative to air; δ = angle of wavefront normal with the horizontal; v w = horizontal wind velocity. The magnitude of the horizontal wind velocity near the Earth's surface is predominantly determined by the prevailing horizontal pressure gradient in the atmosphere and the surface friction [39]. The surface friction arises from the relative motion between air and the ground surface and has to be accounted for heights up to 1000 m.

Ranging Error Analysis
The range measured by acoustic sensors can be written as: where t is the time of flight. R can be rewritten as: The uncertainty in range measurement can be obtained by calculating the deviation in range measurement error from all error sources. The cumulative deviation of ranging error is given by: where: The error budgeting has been numerically validated in a case study, taking realistic values of variables, as shown in Table 3. This analysis gives an error of 0.24 m for a range of 10 m, with more than 95% of the error being due to multipath. However, in real life situations, the hardware limitations associated with the practical realization of acoustic sensors can introduce additional errors, which have to be considered as well in the overall error budgeting.

Monostatic Approach
A major limitation of the multistatic system is that it works in predetermined environments only. To overcome this limitation, a monostatic approach can be applied for obstacle detection and tracking. The generalized SONAR (Sound Navigation and Ranging) equation is given by: where: SNL = signal-to-noise ratio of the returning echo; SL = source sound level; 2TL = two-way transmission losses; TS = target strength; NL = noise level; DI = directivity index. This approach utilizes a collocated transceiver. The transmission losses are a function of spherical spreading and atmospheric attenuation. The basic range equation for the monostatic approach is: where ρ k t is the actual measurement, t k is the nominal time of reception, t t is the nominal time of emission, and c is the speed of sound.

Multistatic Approach
In a familiar environment, for example in a building, a multistatic sensor approach can be applied that utilises base stations (BS) at known fixed locations in the test environment. Similar work has been done so far [40], which utilises transmitters at known locations to calculate the mobile station (MS) or receiver coordinates. A platform fitted with an acoustic receiver (R), as shown in Figure 10, can utilise distance measurements from multiple transmitters to calculate its position. An optimised arrangement of transmitters ensures that the Position Dilution of Precision (PDOP), as shown in Equation (58), is kept within an acceptable threshold [41][42][43]. The optimised geometry also results in keeping the cost of the system down: where σ 2 x , σ 2 y , and σ 2 z are the variances of the estimation errors along each axis and h 1 , h 2 , h 3 , and h 4 are altitudes of the tetrahedron formed by joining unit vectors from four BS and one MS [44].
where ρ k t is the actual measurement, t k is the nominal time of the receiver clock k at reception, t t is the nominal time of the transmitter clock s at emission, and c is the speed of sound. After taking into consideration the clock biases, propagation delay in the air, multipath error and random measurement noise, the complete expression for range measurement becomes: where: ρ k t (t r,k ) = geometric distance (m); dt k = receiver clock error (s); dt t = transmitter clock error (s); P k,t t (t k ) = propagation delay in air (standard conditions) (m); d k,t (t k ) = error due to hardware code delay at the receiver (m); d t t (t k ) = error due to hardware code delay at the transmitter (m); d k,T t (t k ) = multipath error (m); ε t = random measurement noise (m). The multipath error depends on the relative geometry of the transmitter and the receiver with respect to the surrounding reflective surface as well as its reflective properties. The coordinates of the receiver as well as the timing information are derived from the simultaneous/sequential observation of four (or more) transmitters. Assuming a constant clock error for measurements to all transmitters and neglecting all other error terms, the following system of equations is obtained: Considering this system of equations, a minimum of four BSs are required to yield four equations, which are solved to provide the positioning solution. One approach to solving the time-independent non-linear system of equations is to linearize them by using a recursive least squares algorithm for positioning. With more than four BS being used for calculating the position of the MS, the problem is over-determined and can be solved in the least squares sense to yield an optimal estimate of the MS location. In practice, the solution is obtained iteratively starting from an initial guess of the receiver position (x 0 , y 0 , z 0 ) as shown in Figure 11, where: where ∆x, ∆y and ∆z are the differences between the true solution and the initial guesses (corrections to the nominal values). Equation (62) can also be written as: Considering this system of equations, a minimum of four BSs are required to yield four equations, which are solved to provide the positioning solution. One approach to solving the time-independent non-linear system of equations is to linearize them by using a recursive least squares algorithm for positioning. With more than four BS being used for calculating the position of the MS, the problem is over-determined and can be solved in the least squares sense to yield an optimal estimate of the MS location. In practice, the solution is obtained iteratively starting from an initial guess of the receiver position ( , , ) as shown in Figure 11 where ∆ , ∆ and ∆ are the differences between the true solution and the initial guesses (corrections to the nominal values). Equation (62) can also be written as: After simplification and disregarding the higher order error terms: Equation (68) can be written as: where: After simplification and disregarding the higher order error terms: Equation (68) can be written as: where: Denoting the ignored error term by e, Equation (69) can be written as: or: This is an unconstrained minimization problem. J needs to be minimized with respect to n unknown values of ∆r. Thus: Pre-multiplying Equation (74) by −1/2 A T A −1 yields the least squares solution for linear algebraic equations, which can be written as: where: A * * is called the pseudo inverse of A for overdetermined systems. Note that the order of matrix ∆r remains 3 × 1, irrespective of the number of rows in A, as long as there are a minimum of three rows. The number of rows in A denotes the number of BS's that are in range. To correct the initial guess, a linear (first-order Taylor) expansion of ρ k is adopted. Let n (x n , y n , z n ) denote a preliminary best estimate value. The linearised equation is: where R n s is the nominal range measurement to the sth BS and ∆R s is the difference between the actual and nominal range measurements. In addition to the ranging errors that affect the sensor position accuracy, the relative geometry of BS and MS affects the estimation. The linearised equations to s BSs are given by: where B is the matrix of coefficients of the linear set of equations, R ρ is the covariance matrix of the pseudorange errors, and σ 0 2 is a scale factor known as the a priori variance of unit weight.
Various other multilateration algorithms have been proposed in the literature which present numerical algorithms based on iteratively solving simultaneous position equations based on an initial estimate [45][46][47]. A comparison of a linear least squares estimator, a non-linear least squares estimator, and an iteratively reweighted least squares technique was also made [46], where it was shown that the performance of non-linear least squares method was the best. In [48], a probabilistic model of the error in distance measurement for trilateration using Extended Kalman Filtering (EKF) was developed. However, most of the iterative algorithms converge locally, hence are sensitive to the initial estimate. An inappropriate initial estimate can lead to convergence to the wrong local optimum in the vicinity of the initial estimate, described as mirroring in [40]. Besides, iterative algorithms are computationally intensive. Various closed-form trilateration algorithms have also been discussed in the literature [49][50][51][52]. In [53], a closed-form algorithm for solving least squares trilateration problem is discussed, which works for an overdetermined system as well as in cases where there is no real solution. Additionally, global optimization techniques can also be explored, but their computational complexity can render them unsuitable for real-time applications.

Combination of Multistatic and Monostatic Approaches
A third approach combines the monostatic and multistatic approaches. This hybrid approach enables relative navigation among platforms equipped with monostatic transceivers. Swarms of platforms communicate with each other via on-board acoustic sensors to calculate their position as well as augment the knowledge of their environment. The transmitters of the multistatic system are fixed at known locations. The transmitters relay their unique identification information and transmission time along with the acoustic signal to enable positioning based on multilateration. Figure 12 shows the schematic of relative navigation in case of three independent transceivers. T i is the ith transmitter fixed at a known location while TR j is the jth transceiver on-board the platform, with r ij being the distance between the ith transmitter and the jth transceiver at a given instant of time. transmitter fixed at a known location while is the jth transceiver on-board the platform, with being the distance between the ith transmitter and the jth transceiver at a given instant of time.

Overview of State-of-the-Art Acoustic Sensors
Acoustic sensors, mostly ultrasonic, have been widely used in navigation, especially indoor navigation and personal mobility. In [54], an ultrasonic indoor positioning system based on Time Difference of Arrival (TDOA) is discussed. The proposed system consists of fixed active ultrasonic transmitters and passive receivers arranged in an omnidirectional hexagonal pattern on the mobile platform. The indoor navigation system discussed in [40] is based on Time of Arrival (TOA) to calculate the coordinates of the receiver in real time. Synchronization between the transmitters at fixed known positions on the ceiling and the mobile receiver enables Time Division Multiple Access (TDMA) implementation for positioning of the receiver using multilateration. In [55], a combination of globally and locally-referenced Local Positioning Systems (LPS) is introduced to cover an extensive indoor environment. In addition, certain acoustic LPS employ Direct Sequence Code Division Multiple Access (DS-CDMA) [56] techniques, where the receiver on-board the robot determines its position using TDOA between a reference beacon and other beacons, with each beacon transmitting a unique 255-bit Kasami code. More advanced Code Division Multiple Access (CDMA) based acoustic positioning systems determine TOA at the receiver using signal correlation [57], used in all modern Radio Frequency (RF)-based communication systems, including the GNSS.
Certain navigation systems like [58][59][60][61] use both RF as well as ultrasonic signals to calculate the position of the mobile platform. Cricket indoor location system [58] uses RF signals for synchronizing the ultrasonic transmitters and receiver on the mobile platform. In Dolphin [59], the

Overview of State-of-the-Art Acoustic Sensors
Acoustic sensors, mostly ultrasonic, have been widely used in navigation, especially indoor navigation and personal mobility. In [54], an ultrasonic indoor positioning system based on Time Difference of Arrival (TDOA) is discussed. The proposed system consists of fixed active ultrasonic transmitters and passive receivers arranged in an omnidirectional hexagonal pattern on the mobile platform. The indoor navigation system discussed in [40] is based on Time of Arrival (TOA) to calculate the coordinates of the receiver in real time. Synchronization between the transmitters at fixed known positions on the ceiling and the mobile receiver enables Time Division Multiple Access (TDMA) implementation for positioning of the receiver using multilateration. In [55], a combination of globally and locally-referenced Local Positioning Systems (LPS) is introduced to cover an extensive indoor environment. In addition, certain acoustic LPS employ Direct Sequence Code Division Multiple Access (DS-CDMA) [56] techniques, where the receiver on-board the robot determines its position using TDOA between a reference beacon and other beacons, with each beacon transmitting a unique 255-bit Kasami code. More advanced Code Division Multiple Access (CDMA) based acoustic positioning systems determine TOA at the receiver using signal correlation [57], used in all modern Radio Frequency (RF)-based communication systems, including the GNSS.
Certain navigation systems like [58][59][60][61] use both RF as well as ultrasonic signals to calculate the position of the mobile platform. Cricket indoor location system [58] uses RF signals for synchronizing the ultrasonic transmitters and receiver on the mobile platform. In Dolphin [59], the relative position of a node placed on the ceiling is calculated by trilateration using TOA of ultrasonic signals from three nodes. The RF signal from the nodes contains node location information as well as time synchronization. A commercial-off-the-shelf (COTS) hardware components based architecture for a self-driving miniature vehicle is developed using real-time operating system (RTOS) in [62]. There 1:10 scale self-driving vehicle has three ultrasonic sensors, three infrared sensors and a camera. In [63], an ultrasonic relative positioning system is developed for a group of ground-based robots. Moving in formation, the robots can determine the distance and orientation of nearby robots based on TOF evaluation of ultrasonic pulses as well as a RF communication link. In [60], personal tracking is achieved with the help of a portable unit called Bat, which consists of a radio transceiver, controlling logic and an ultrasonic transducer. Ultrasound receiver units are placed on the ceiling at known points and are interconnected. The Bats are triggered by a radio message from the base station, which synchronizes the trigger with the receivers. In [61], a multi-block navigation system is developed which utilizes RFID transmitters to trigger relevant beacons to send an ultrasonic signal for localization of the mobile object. Table 4 lists some representative COTS ultrasonic ranging sensors operating at frequencies ranging from 40 kHz to 500 kHz. It can be observed that as the operating frequency of the ultrasonic ranging sensor increases, the detection range decreases mainly due to higher attenuation. As described in Section 2, bats can vary the frequency of their echolocation calls based on the distance from their prey or an obstruction, which allows them to achieve optimal range and angular resolution performance in various conditions. Various acoustic navigation aids for visually impaired have been investigated. An ultrasonic FM mobility aid is proposed in [64]. While, in [65], a bio-inspired mobility aid for visually impaired based on echolocation of bats is discussed. This device uses downswept FM ultrasound emissions to detect obstacles, which are perceived as localized sound images corresponding to the direction and size of the obstacles. In [66], a digital signal processor based ultrasonic navigation aid for the visually impaired is presented. The nearest obstacle in front of the user is detected using 2-dimensional echolocation and a binaural audio feedback is provided. GuideCane [67] and UltraCane [68] use ultrasonic sensors in the cane to help visually impaired navigate quickly and safely in the presence of obstacles and hazards.

Integration of Acoustic Sensors in Multi-Sensor Navigation Systems
Currently, air and surface vehicles rely on a combination of GNSS, Inertial Measurement Unit (IMU) and in some cases, Vision-Based Navigation (VBN) and Aircraft Dynamics Model (ADM) virtual sensor for air vehicle navigation and guidance [2]. Integration of ANS (Acoustic Navigation System) with the existing NGS (Navigation and Guidance System) enables accurate and reliable positioning, even in low visibility indoor environments, using low Size, Weight and Power, and Cost (SWaP-C) sensors. GNSSS signals are prone to data degradations or complete loss of signal in dense urban environments as well as indoors due to multipath effects, interference, or antenna obscuration [36,69]. High precision IMUs employ relatively high cost, weight and volume inertial components (e.g., ring laser and fibre optic gyroscopes), rendering them unsuitable for many ground and air vehicles like cars, Remotely Piloted Aircraft Systems (RPAS), etc. Although commercially available Micro Electromechanical Systems (MEMS) devices can support the design of low cost IMUs, their accuracy decreases steeply with operating time. The performance of VBN sensors is affected by low visibility conditions which could be due to night/low light conditions, smoke, fog, precipitation as well as presence of transparent or opaque objects [70]. An ADM virtual sensor is dependent on the aircraft's physical sensors and is based on the assumption of the aircraft being a rigid body with a constant and static mass distribution [71]. An integration of sensor data achieved through Multi-Sensor Data Fusion algorithms can lead to an improved positioning solution for air and ground platforms. The performance analysis of different multi-sensor integrated architectures, including the Extended Kalman Filter (EKF)-based VBN-IMU-GNSS (EVIG), the EKF-based VBN-IMU-GNSS-ADM (EVIGA) and the Unscented Kalman Filter (UKF)-based VIGA (UVIGA) have shown that these integration schemes provide improvements in position, velocity, and attitude (PVA) data in all flight phases in air vehicles, when compared to individual sensor measurements. In particular, the EVIGA and UVIGA systems achieve horizontal/vertical position accuracies in line with International Civil Aviation Organization (ICAO) CAT-I and CAT-II requirements [71]. Figure 13 presents the NGS architecture including the ANS and thereby resulting in an ANS-VIGA (AVIGA) integration scheme.  [71]. Figure 13 presents the NGS architecture including the ANS and thereby resulting in an ANS-VIGA (AVIGA) integration scheme. In this integrated system, the data output rate for ANS is 2 Hz, 20 Hz for VBN and GNSS at 2 Hz to augment the MEMS-IMU running at 100 Hz. The IMU position and velocity information are compared to the Global Positioning System (GPS) position and velocity and form the measurement input of the data fusion block. A similar process is also applied to the attitude data, the differences of which are used in the data-fusion algorithms. The data-fusion algorithm provides estimates of the PVA errors, which are then removed from the sensor measurements to obtain the corrected navigation states. The corrected PVA as well as the estimates of the accelerometer and the gyroscope biases are used to update the IMU raw measurements. The attitude data provided by ADM augmentation and by the IMU are compared to feed the data-fusion block at 100 Hz. The attitude data provided by the VBN and IMU sensors are compared at 20 Hz. The best estimate of attitude is compared with the IMU attitude to obtain the corrected attitude. By employing a UKF, the AVIGA system performance in terms of attitude data accuracy can be increased in addition to a significant extension of the ADM validity time. Additionally, a pre-filter can be also used to pre-process the Six Degrees-of-Freedom (6-DoF) ADM navigation solution. Another pre-filter can be employed to process the ANS data to remove any outliers. In order to select the navigation sensors based on its performance (accuracy, availability, continuity, and integrity), a sensor selection and prioritisation approach is employed using Adaptive Boolean Decision Logic (ABDL), as shown in Figure 13. Multiple operating modes can be implemented using the ABDL, wherein each mode is a unique combination of sensors. The selection of sensors is dictated by the flight phase, integrity alerts (both In this integrated system, the data output rate for ANS is 2 Hz, 20 Hz for VBN and GNSS at 2 Hz to augment the MEMS-IMU running at 100 Hz. The IMU position and velocity information are compared to the Global Positioning System (GPS) position and velocity and form the measurement input of the data fusion block. A similar process is also applied to the attitude data, the differences of which are used in the data-fusion algorithms. The data-fusion algorithm provides estimates of the PVA errors, which are then removed from the sensor measurements to obtain the corrected navigation states. The corrected PVA as well as the estimates of the accelerometer and the gyroscope biases are used to update the IMU raw measurements. The attitude data provided by ADM augmentation and by the IMU are compared to feed the data-fusion block at 100 Hz. The attitude data provided by the VBN and IMU sensors are compared at 20 Hz. The best estimate of attitude is compared with the IMU attitude to obtain the corrected attitude. By employing a UKF, the AVIGA system performance in terms of attitude data accuracy can be increased in addition to a significant extension of the ADM validity time. Additionally, a pre-filter can be also used to pre-process the Six Degrees-of-Freedom (6-DoF) ADM navigation solution. Another pre-filter can be employed to process the ANS data to remove any outliers. In order to select the navigation sensors based on its performance (accuracy, availability, continuity, and integrity), a sensor selection and prioritisation approach is employed using Adaptive Boolean Decision Logic (ABDL), as shown in Figure 13.
Multiple operating modes can be implemented using the ABDL, wherein each mode is a unique combination of sensors. The selection of sensors is dictated by the flight phase, integrity alerts (both preventive and reactive) and other factors. Table 5 lists some of the typical characteristics of COTS sensors used for navigation in small sized platforms.

Conclusions and Recommendations for Future Research
This paper presents acoustic wave propagation and its applications in air and surface vehicles navigation. Taking inspiration from echolocating mammals, especially bats, novel acoustic navigation techniques are introduced. Various sound wave attenuation factors like geometric divergence, atmospheric absorption, ground effect, screening, and wind and temperature gradient effects are discussed in detail. Mathematical error modeling for acoustic range measurements, taking into consideration Doppler shift, multipath and atmospheric effects, is presented. Different acoustic sensor arrangements for navigation, monostatic, multistatic and a combination of the two, are discussed. Also, the state-of-the-art in acoustic sensors and their applications in navigation are presented. Finally, the integration of acoustic sensors in multi-sensor navigation systems is discussed.
The use of acoustic sensors for air and surface vehicle navigation holds very high potential, especially when addressing low SWAP-C requirements. Current research trends indicate that future acoustic sensors will implement sophisticated bio-inspired features, which may significantly enhance their range and resolution performance [75,76]. In future applications, acoustic sources might also be used as Signals of Opportunity (SoOP) for navigation in GNSS challenged/denied environments [77]. Although acoustic sensor networks have been utilized in the design of smart city networks [78] and for intelligent transport systems [79], their potential has not yet been fully exploited in air and surface navigation systems. So, acoustic sensors are expected to attract much research and industry efforts in the near future, due to a growing interest in alternative sources of navigation data in urban/indoor environments for a variety of civil and military applications. In order to undertake the development of either monostatic or multistatic acoustic sensors, a careful analysis of aero-acoustic effects is required, especially for aviation applications. In particular, it is essential to assess the magnitude of potential in-band and out-of-band interferences (e.g., engine, propellers, aerodynamic noise) and implement appropriate engineering solutions that mitigate/prevent the negative effects of such interferences on sensor performance. For instance, the use of larger propellers in rotary-wing aircraft might significantly reduce in-band interferences for acoustic sensors operating at higher frequencies. Similar considerations apply to the acoustic emissions produced by gas-turbine engines (i.e., size and number of rotor blades). In this case, however, noise reduction and frequency tailoring can be obtained by adopting other engineering solutions, such as chevrons, liners and noise absorptive materials [80].