Propagation & channel modeling
API reference
Environment models
Create a generic underwater environment with the given parameters. The following parameters are supported:
bathymetry= bathymetry modelaltimetry= altimetry modeltemperature= temperature modelsalinity= salinity modelpH= pH modelsoundspeed= sound speed profile modeldensity= density modelseabed= seabed sediment modelsurface= surface model
All parameters are optional and have default values.
Propagation models
A fast differentiable mode propagation model that only supports iso-velocity constant depth environments.
ngrid is the number of grid points to use for modal root finding for fluid bottom environments. If ngrid is too small, the mode solver may miss some modes. If ngrid is too large, the mode solver may take a long time to converge. The default value of ngrid of 0 will use a heuristic to automatically determine the number of grid points to use.
nmodes controls the maximum number of modes computed. Setting nmodes to 0 causes all propagating modes to be computed.
A 2D adiabatic mode propagation model based on a range-independent modal propagation model. The adiabatic mode model supports range-dependent bathymetry. Any kwargs passed in are transferred to the underlying model.
dx is the step size in the range direction for integration. dz is the mesh size in the depth direction to compute modes. If dx and/or dz is set to zero, it is automatically determined (~10 points per wavelength).
Adiabatic mode models are usually not reciprocal, i.e., exchanging a source and receiver changes the answer. This is because the bathymetry is assumed to have azimuthal symmetry around the source. In reciprocal mode, the adiabatic computation is modified to ensure reciprocity, but azimuthal symmetry is lost.
Type representing a single acoustic ray arrival.
Properties:
t/time: arrival time (s)ϕ/phasor: complex amplitudens/surface_bounces: number of surface bouncesnb/bottom_bounces: number of bottom bouncesθₛ/launch_angle: launch angle at source (rad)θᵣ/arrival_angle: arrival angle at receiver (rad)path: ray path (optional, vector of 3-tuples or missing)
The properties are accessible with the short names for brevity, and longer more descriptive names where readability is desired or unicode symbols are undesired.
Type representing a single acoustic mode arrival.
Properties:
m/mode: mode numberkᵣ/hwavenumber: horizontal wavenumber (rad/m)ψ(z)/mode_function: mode functionv/group_velocity: group velocity (m/s)vₚ/phase_velocity: phase velocity (m/s)
The properties are accessible with the short names for brevity, and longer more descriptive names where readability is desired or unicode symbols are undesired.
Compute the arrivals at the receiver rx due to the source tx using propagation model pm. Returns an array of arrivals.
For ray models, eigenray paths are typically included in the arrivals. However, if they are not needed, one may set paths=false to allow the propagation model to avoid computing them.
Compute the acoustic field at the receivers rxs due to the source tx using propagation model pm. If rxs denotes a single receiver, the result is a complex scalar. If rxs is an AbstractArray, the result is an array of complex numbers with the same shape as rxs. The amplitude of the field is related to the transmission loss, and the angle is related to the acoustic phase at the source frequency.
Compute the acoustic field at a receiver rxs due to a transmitter tx in the Pekeris waveguide. The field is computed incoherently if mode=:incoherent. Otherwise, the field is computed coherently.
Compute the acoustic field at a receiver rxs due to a transmitter tx in the Pekeris waveguide. The field can be computed incoherently if mode=:incoherent. Otherwise, the field is computed coherently.
Compute the impulse response at the receiver rx due to the source tx using propagation model pm at the given sampling frequency fs. If abstime is true, the result is in absolute time from the start of transmission. Otherwise, the result is relative to the earliest arrival time of the signal at the receiver (possibly with some guard period to accommodate acausal response). ntaps specifies the number of taps in the impulse response. If not specified, the number of taps is chosen automatically based on the arrival times.
Compute the impulse response at the receiver rx due to the source tx using propagation model pm at the given sampling frequency fs.
Several kwargs may be specified:
- If
abstimeistrue(default:false), the result is in absolute time from the start of transmission. Otherwise, the result is relative to the earliest arrival time of the signal at the receiver (with some guard period to accommodate acausal response). ntaps(default:nothingfor automatic) specifies the number of taps in the impulse response.fminandfmaxspecifies the bandwidth of interest (default:0tofs/2). If the impulse response is used to convolve with bandlimited signals, it is recommended that the impulse response bandwidth be reduced to match the signal bandwidth to reduce computational load and improve numerical stability.nmodes(default:nothingfor no limit) is the maximum of modes used in impulse response computation.threshold(default:-60dB) is used to drop attenuated modes to manage computational load.acausal(default:0.02s) controls the computation of acausal impulse response (before the estimated time of arrival of first mode).taper(default:0.1) applies a tukey window to the impulse response to limit it to the estimated delay spread.
The impulse response is computed only for positive frequencies. Such an impulse response is suitable for convolution with passband complex analytic signals. If convolved with real signals, the resulting signal is approximately equivalent to converting the real signal to a complex analytic form and then convolving it with the impulse response.
Channel models and simulation
Compute a channel model from the sources txs to the receivers rxs using propagation model pm. The result is a channel model with the same number of input channels as the number of sources and output channels as the number of receivers. The channel model accepts signals sampled at rate fs and returns signals sampled at the same rate.
An additive noise model may be optionally specified as noise. If specified, it is used to corrupt the received signals.
Propagation model specific keyword arguments kwargs supported by impulse_response() can be passed through when generating a channel.
BasebandReplayChannel(h, θ, fs, fc, step=1; noise=nothing)
BasebandReplayChannel(h, fs, fc, step=1; noise=nothing)Construct a baseband replay channel with impulse responses h and phase estimates θ. The phase estimates are optional. fs is the sampling frequency in Sa/s, fc is the carrier frequency in Hz, and step is the decimation rate for the time axis of h. The effective sampling frequency of the impulse responses is fs ÷ step impulse responses per second.
An additive noise model may be optionally specified as noise. If specified, it is used to corrupt the received signals.
Load a baseband replay channel from a file.
If upsample is true, the impulse responses are upsampled to the delay axis sampling rate. This makes applying the channel faster but requires more memory. rxs controls which receivers to load from the file. By default, all receivers are loaded.
An additive noise model may be optionally specified as noise. If specified, it is used to corrupt the received signals.
Supported formats:
.mat(MATLAB) file in underwater acoustic channel repository (UACR) format. See https://github.com/uwa-channels/ for details.
Simulate the transmission of passband signal x through the channel model ch. If txs is specified, it specifies the indices of the sources active in the simulation. The number of sources must match the number of channels in the input signal. If rxs is specified, it specifies the indices of the receivers active in the simulation. Returns the received signal at the specified (or all) receivers.
fs specifies the sampling rate of the input signal. The output signal is sampled at the same rate. If fs is not specified but x is a SampledSignal, the sampling rate of x is used. Otherwise, the signal is assumed to be sampled at the channel’s sampling rate.
If abstime is true, the returned signals begin at the start of transmission. Otherwise, the result is relative to the earliest arrival time of the signal at any receiver. If noisy is true and the channel has a noise model associated with it, the received signal is corrupted by additive noise.
Simulate the transmission of passband signal x through the channel model ch. If txs is specified, it specifies the indices of the sources active in the simulation. The number of sources must match the number of channels in the input signal. If rxs is specified, it specifies the indices of the receivers active in the simulation. Returns the received signal at the specified (or all) receivers.
fs specifies the sampling rate of the input signal. The output signal is sampled at the same rate. If fs is not specified but x is a SampledSignal, the sampling rate of x is used. Otherwise, the signal is assumed to be sampled at the channel’s sampling rate.
If abstime is true, the returned signals begin at the start of transmission. Otherwise, the result is relative to the earliest arrival time of the signal at any receiver. If noisy is true and the channel has a noise model associated with it, the received signal is corrupted by additive noise.
If start is specified, it specifies the starting time index in the replay channel. If not specified, a random start time is chosen.
Boundary conditions
Compute complex reflection coefficient at a fluid-solid boundary, given:
- angle of incidence
θ(angle to the surface normal) - relative density of the reflecting medium to incidence medium
ρᵣ - relative compressional sound speed of the reflecting medium to incidence medium
cᵣ - relative shear wave speed of the reflecting medium to incidence medium
cᵣₛ - dimensionless compressional absorption coefficient
δ - dimensionless shear absorption coefficient
δₛ
WARNING: Fluid-solid reflection not implemented. Currently ignores shear.
Create a fluid half-space boundary with density ρ, sound speed c, dimensionless absorption coefficient δ, and interfacial roughness σ.
The absorption coefficient δ may also be specified in dB/m/kHz by attaching appropriate units to it (e.g. 13.2u"dB/m/kHz").
Examples
julia> FluidBoundary(1200, 1500)
FluidBoundary(ρ=1200.0, c=1500.0)
julia> FluidBoundary(1200, 1500, 0.1)
FluidBoundary(ρ=1200.0, c=1500.0, δ=0.1)
julia> FluidBoundary(1200, 1500, 0.1, 0.1)
FluidBoundary(ρ=1200.0, c=1500.0, δ=0.1, σ=0.1)
julia> FluidBoundary(ρ=1200, c=1500)
FluidBoundary(ρ=1200.0, c=1500.0)
julia> FluidBoundary(ρ=1200, c=1500, δ=0.1, σ=0.1)
FluidBoundary(ρ=1200.0, c=1500.0, δ=0.1, σ=0.1)
julia> FluidBoundary(ρ=1200.0, c=0.0)
PressureReleaseBoundary
julia> FluidBoundary(ρ=1200.0, c=Inf)
RigidBoundary
julia> FluidBoundary(1.2u"g/cm^3", 1500u"m/s", 13.2u"dB/m/kHz")
FluidBoundary(ρ=1200.0, c=1500.0, δ=0.362803)ElasticBoundary(ρ, cₚ, cₛ)
ElasticBoundary(ρ, cₚ, cₛ, δₚ, δₛ)
ElasticBoundary(ρ, cₚ, cₛ, δₚ, δₛ, σ)
ElasticBoundary(b::FluidBoundary)
ElasticBoundary(b::FluidBoundary, δₛ)
ElasticBoundary(b::FluidBoundary, cₛ, δₛ)Create a solid half-space boundary with density ρ, compressional sound speed cₚ, shear wave speed cₛ, dimensionless compressional absorption coefficient δₚ, dimensionless shear absorption coefficient δₛ, and interfacial roughness σ. If the absorption coefficients or interfacial roughness are unspecified, they are assumed to be 0.
An ElasticBoundary may also be constructed from a b::FluidBoundary by adding shear wave speed cₛ and shear absorption coefficient δₛ. If cₛ is not specified, it is computed using shearspeed(b.cₚ).
Absorption coefficients δₚ and δₛ may also be specified in dB/m/kHz by attaching appropriate units to them (e.g. 13.2u"dB/m/kHz").
Examples
julia> ElasticBoundary(1200, 1500, 500)
ElasticBoundary(ρ=1200.0, cₚ=1500.0, cₛ=500.0)
julia> ElasticBoundary(1200, 1500, 500, 0.1, 0.2)
ElasticBoundary(ρ=1200.0, cₚ=1500.0, cₛ=500.0, δₚ=0.1, δₛ=0.2)
julia> ElasticBoundary(1200, 1500, 500, 0.1, 0.2, 0.1)
ElasticBoundary(ρ=1200.0, cₚ=1500.0, cₛ=500.0, δₚ=0.1, δₛ=0.2, σ=0.1)
julia> ElasticBoundary(ρ=1200, cₚ=1500, cₛ=500)
ElasticBoundary(ρ=1200.0, cₚ=1500.0, cₛ=500.0)
julia> ElasticBoundary(ρ=1200, cₚ=1500, cₛ=500, δₚ=0.1, δₛ=0.2)
ElasticBoundary(ρ=1200.0, cₚ=1500.0, cₛ=500.0, δₚ=0.1, δₛ=0.2)
julia> ElasticBoundary(ρ=1200, cₚ=1500, cₛ=500, δₚ=0.1, δₛ=0.2, σ=0.1)
ElasticBoundary(ρ=1200.0, cₚ=1500.0, cₛ=500.0, δₚ=0.1, δₛ=0.2, σ=0.1)
julia> ElasticBoundary(1.2u"g/cm^3", 1500u"m/s", 500u"m/s")
ElasticBoundary(ρ=1200.0, cₚ=1500.0, cₛ=500.0)
julia> ElasticBoundary(FineSand)
ElasticBoundary(ρ=1484.373, cₚ=1691.954, cₛ=414.413, δₚ=0.01602)
julia> ElasticBoundary(FineSand, 4.3u"dB/m/kHz")
ElasticBoundary(ρ=1484.37, cₚ=1691.95, cₛ=414.4113571750004, δₚ=0.01602, δₛ=0.03265)
julia> ElasticBoundary(FineSand, 500, 0.1)
ElasticBoundary(ρ=1484.373, cₚ=1691.954, cₛ=500.0, δₚ=0.01602, δₛ=0.1)MultilayerElasticBoundary([(h, ρ, cₚ, cₛ), ...])
MultilayerElasticBoundary([(h, ρ, cₚ, cₛ, δₚ, δₛ), ...])
MultilayerElasticBoundary([(h, ρ, cₚ, cₛ, δₚ, δₛ, σ), ...])
MultilayerElasticBoundary([(h, b::FluidBoundary), ...])
MultilayerElasticBoundary([(h, b::ElasticBoundary), ...])Create a multilayer solid boundary with layers defined by a vector of tuples specifying the layer thickness h, density ρ, compressional sound speed cₚ, shear wave speed cₛ, dimensionless compressional absorption coefficient δₚ, dimensionless shear absorption coefficient δₛ, and interfacial roughness σ for each layer. If the absorption coefficients or interfacial roughness is unspecified, they are assumed to be 0. The last entry in the vector is considered the bottom half-space and has an infinite thickness. By convention, we specify a value of Inf for h for that layer.
ρ, cₚ, and cₛ may also be specified with linear variation within a layer by specifying the value as a 2-tuple, with the first entry being the value at the top of the layer and the second entry being the value at the bottom.
The properties of a layer may be specified by giving a FluidBoundary or ElasticBoundary layer instead of ρ, cₚ, cₛ, δₚ, δₛ and σ.
Absorption coefficients δₚ and δₛ may also be specified in dB/m/kHz by attaching appropriate units to them (e.g. 13.2u"dB/m/kHz").
Examples
julia> MultilayerElasticBoundary([
(h = 5.2, FineSand),
(h = Inf, ρ = 2000, cₚ = 2500, cₛ = 500)
])
MultilayerElasticBoundary(2 layers):
(h = 5.2, ρ = 1484.373, cₚ = 1691.954, cₛ = 0.0, δₚ = 0.01602, δₛ = 0.0, σ = 0.0)
(h = Inf, ρ = 2000.0, cₚ = 2500.0, cₛ = 500.0, δₚ = 0.0, δₛ = 0.0, σ = 0.0)
julia> MultilayerElasticBoundary([
(5.2, 1300, 1700, 100),
(Inf, 2000, 2500, 500, 0.1, 0.2)
])
MultilayerElasticBoundary(2 layers):
(h = 5.2, ρ = 1300.0, cₚ = 1700.0, cₛ = 100.0, δₚ = 0.0, δₛ = 0.0, σ = 0.0)
(h = Inf, ρ = 2000.0, cₚ = 2500.0, cₛ = 500.0, δₚ = 0.1, δₛ = 0.2, σ = 0.0)
julia> MultilayerElasticBoundary([
(5.2, 1300, (1700,2000), 100),
(Inf, 2000, 2500, 500, 0.1, 0.2)
])
MultilayerElasticBoundary(2 layers):
(h = 5.2, ρ = 1300.0, cₚ = (1700.0, 2000.0), cₛ = 100.0, δₚ = 0.0, δₛ = 0.0, σ = 0.0)
(h = Inf, ρ = 2000.0, cₚ = 2500.0, cₛ = 500.0, δₚ = 0.1, δₛ = 0.2, σ = 0.0)
julia> MultilayerElasticBoundary([
(5.2, 1300, 1700, 100, 0, 0, 0.1),
(Inf, 2000, 2500, 500, 0.1, 0.2, 0.1)
])
MultilayerElasticBoundary(2 layers):
(h = 5.2, ρ = 1300.0, cₚ = 1700.0, cₛ = 100.0, δₚ = 0.0, δₛ = 0.0, σ = 0.1)
(h = Inf, ρ = 2000.0, cₚ = 2500.0, cₛ = 500.0, δₚ = 0.1, δₛ = 0.2, σ = 0.1)
julia> MultilayerElasticBoundary([
(h = 5.2, ρ = 1300, cₚ = 1700, cₛ = 100, δₚ = 0.1, δₛ = 0.2, σ = 0.1),
(h = Inf, ρ = 2000, cₚ = 2500, cₛ = 500)
])
MultilayerElasticBoundary(2 layers):
(h = 5.2, ρ = 1300.0, cₚ = 1700.0, cₛ = 100.0, δₚ = 0.1, δₛ = 0.2, σ = 0.1)
(h = Inf, ρ = 2000.0, cₚ = 2500.0, cₛ = 500.0, δₚ = 0.0, δₛ = 0.0, σ = 0.0)
julia> MultilayerElasticBoundary([
(5.2u"m", 1.3u"g/cm^3", 1700u"m/s", 100u"m/s", 1.2u"dB/m/kHz", 2.3u"dB/m/kHz"),
(Inf, 2u"g/cm^3", 2500u"m/s", 500u"m/s")
])
MultilayerElasticBoundary(2 layers):
(h = 5.2, ρ = 1300.0, cₚ = 1700.0, cₛ = 100.0, δₚ = 0.03738, δₛ = 0.00421, σ = 0.0)
(h = Inf, ρ = 2000.0, cₚ = 2500.0, cₛ = 500.0, δₚ = 0.0, δₛ = 0.0, σ = 0.0)Pre-defined boundary conditions based on APL-UW Technical Report 9407:
# sea surface boundary conditions
const SeaState0 = WindySurface(0.8)
const SeaState1 = WindySurface(2.6)
const SeaState2 = WindySurface(4.4)
const SeaState3 = WindySurface(6.9)
const SeaState4 = WindySurface(9.8)
const SeaState5 = WindySurface(12.6)
const SeaState6 = WindySurface(19.3)
const SeaState7 = WindySurface(26.5)
const SeaState8 = WindySurface(30.6)
const SeaState9 = WindySurface(32.9)# seabed fluid boundary conditions
const Rock = FluidBoundary(2557.50, 3820.00, 0.01374)
const Pebbles = FluidBoundary(2557.50, 2750.40, 0.01374)
const SandyGravel = FluidBoundary(2549.32, 2073.50, 0.01705)
const VeryCoarseSand = FluidBoundary(2455.20, 1996.64, 0.01667)
const MuddySandyGravel = FluidBoundary(2367.22, 1952.48, 0.0163)
const CoarseSand = FluidBoundary(2282.31, 1910.46, 0.01638)
const GravellyMuddySand = FluidBoundary(2200.47, 1870.42, 0.01645)
const MediumSand = FluidBoundary(1887.44, 1800.29, 0.01624)
const MuddyGravel = FluidBoundary(1652.15, 1741.31, 0.0161)
const FineSand = FluidBoundary(1484.37, 1691.95, 0.01602)
const MuddySand = FluidBoundary(1369.80, 1650.24, 0.01728)
const VeryFineSand = FluidBoundary(1297.16, 1614.79, 0.01875)
const ClayeySand = FluidBoundary(1252.15, 1583.62, 0.02019)
const CoarseSilt = FluidBoundary(1222.49, 1555.35, 0.02158)
const SandySilt = FluidBoundary(1195.89, 1527.85, 0.01261)
const MediumSilt = FluidBoundary(1175.43, 1510.43, 0.00676)
const SandyMud = FluidBoundary(1175.43, 1508.59, 0.00386)
const FineSilt = FluidBoundary(1174.40, 1506.76, 0.00306)
const SandyClay = FluidBoundary(1173.38, 1504.93, 0.00242)
const VeryFineSilt = FluidBoundary(1173.38, 1503.09, 0.00194)
const SiltyClay = FluidBoundary(1172.36, 1501.11, 0.00163)
const Clay = FluidBoundary(1171.34, 1497.44, 0.00148)# seabed elastic boundary conditions
const ElasticRock = ElasticBoundary(Rock, 0.07u"dB/m/kHz")
const ElasticPebbles = ElasticBoundary(Pebbles, 0.07u"dB/m/kHz")
const ElasticSandyGravel = ElasticBoundary(SandyGravel, 13.2u"dB/m/kHz")
const ElasticVeryCoarseSand = ElasticBoundary(VeryCoarseSand, 13.2u"dB/m/kHz")
const ElasticMuddySandyGravel = ElasticBoundary(MuddySandyGravel, 4.8u"dB/m/kHz")
const ElasticCoarseSand = ElasticBoundary(CoarseSand, 13.2u"dB/m/kHz")
const ElasticGravellyMuddySand = ElasticBoundary(GravellyMuddySand, 4.8u"dB/m/kHz")
const ElasticMediumSand = ElasticBoundary(MediumSand, 13.2u"dB/m/kHz")
const ElasticMuddyGravel = ElasticBoundary(MuddyGravel, 4.8u"dB/m/kHz")
const ElasticFineSand = ElasticBoundary(FineSand, 13.2u"dB/m/kHz")
const ElasticMuddySand = ElasticBoundary(MuddySand, 4.8u"dB/m/kHz")
const ElasticVeryFineSand = ElasticBoundary(VeryFineSand, 13.2u"dB/m/kHz")
const ElasticClayeySand = ElasticBoundary(ClayeySand, 4.8u"dB/m/kHz")
const ElasticCoarseSilt = ElasticBoundary(CoarseSilt, 13.4u"dB/m/kHz")
const ElasticSandySilt = ElasticBoundary(SandySilt, 13.4u"dB/m/kHz")Sources and receivers
AcousticSource(pos, frequency; spl=0)
AcousticSource(x, z, frequency; spl=0)
AcousticSource(x, y, z, frequency; spl=0)An source at location pos with nominal frequency and source level spl (dB re 1 µPa @ 1 m). The source is assumed to be omnidirectional and well approximated by a point source. While the source may have some bandwidth, the nominal frequency is used to estimate propagation effects such as absorption, reflection coefficients, etc.
If the location of the source is unknown, it may be specified as nothing. This is useful when the propagation model does not require the source location (e.g., data-driven models).
Noise models
rand([rng::AbstractRNG, ] noise::AbstractNoiseModel, nsamples; fs)
rand([rng::AbstractRNG, ] noise::AbstractNoiseModel, nsamples, nchannels; fs)Generate random noise samples from the noise model noise with the specified size nsamples (can be a 2-tuple for multichannel noise). The noise is returned as a signal sampled at fs. The optional rng argument specifies the random number generator to use.