© The Author(s), Published by EDP Sciences, 2025

1 Introduction

The heavy impact of the environment during sound propagation measurements, which creates spurious and unwanted effects on the signals registered, has been at the center of discussion in the acoustician community since the beginning. Embleton [1] was among the earliest researchers to compile a systematic overview of key environmental factors affecting outdoor sound propagation. Subsequent studies have expanded and refined this work, notably through the comprehensive treatment provided in the book by Attenborough [2] and in his more recent one [3]. However, due to the ever-changing environmental conditions, even during the same measurement campaign (wind direction, temperature, turbulence, etc.), maintaining constant environmental conditions and accounting for subtle variations presents significant challenges inpractice. For these reasons, during the years studies were conducted regarding this matter, all aiming to develop frameworks that could correctly model them and/or predict their effects, like [4, 5] where the impact of the atmosphere was investigated. Regarding the influence of the ground, C. Madshus et al. [6], showed the effects of air-ground interaction connected with high-energy blast sound propagation. A more comprehensive approach was proposed by Wilson et al. [7], in which the impact of all these uncertainties was analyzed and quantified using stochastic sampling techniques. A relevant aspect in assessing the uncertainties is that some of the variables (e.g., the noise output of the source, the weather, the background noise) are almost never under control and hard to quantify. In a study from D. Valente et al. [8] it was concluded that long-term average meteorological conditions, seasonality, climate, and time of day are not sufficient to predict variability. Instantaneous meteorological conditions are in fact responsible for significant variability of the sound pressure levels; if the ISO 17201-1 [9] prediction method is used, the potential uncertainties affecting the determination of the source can be up to ±20 [dB]. The study also observed that the variability stabilizes beyond 2 [km].

The subject of uncertainty estimation can become arbitrarily articulated, depending on the use intended for the measurements. An extensive treatment of this topic can be found in the literature [10–12]. In general, the numerical complexity of the problem and the linked experimental workload necessary to characterize environmental conditions and uncertainties resulted in most commercial software and standards offering limited and coefficient-based capabilities for environmental and ground effects prediction. While sufficient for many practical scenarios, such tools may not provide sufficient resolution for detailed academic investigations, where capturing the full variability of environmental effects is essential. On top of it, within the specific topic of firearm noise propagation, not many large datasets are available, rendering stochastic techniques to evaluate uncertainties and environmental effects almost impossible to implement due to lack of data. The research performed by G. Billot et al. [13] contributed to expanding the available firearm noise dataset through multiple measurements conducted to validate the nonlinear propagation model on which the paper builds. While the proposed model accounts for refractive effects caused by wind and temperature gradients, as well as reflections from a perfectly rigid ground, the authors emphasize in their conclusions the need for deeper research into uncertainties related to atmospheric conditions and more refined ground characterizations, further reinforcing the relevance and motivation of the present study.

On this premise, the scope of this research is to better understand and portray the experimental uncertainties regarding the environmental effects of outdoor firearm noise propagation, with the aim to further refine the quantification of the impact of ground and atmosphere. This objective will be pursued with both an experimental and numerical approach. The findings of this study will contribute to the advancement of a physics-based framework for modeling and controlling firearm blast wave propagation, supporting the development of effective mitigation strategies and improved hearing protection measures for operators. The first part (Sects. 2 and 3) will describe in detail the experiment conducted, focusing on the instrumentations and setups adopted, in order to obtain at the end an overall uncertainty budget that can account for both statistical and random errors linked to firearm noise propagation. The second part (Sects. 4 and 5) will focus on improving the 2D Nonlinear Progressive wave Equation (NPE) described in [13] by including a more refined ground model, the porous ground, that will be used to assess the impact of different terrains on the propagation of a blast-like noise. The same model will then be used to evaluate the impact of different wind and atmospheric conditions. Since the primary goal of integrating the updated NPE code is to establish a controlled computational environment that permits the isolation and analysis of the impact of environmental factors rather than to replicate particular measurement results, a direct comparison between the computational simulations and experimental measurements is not used in this study.

2 Experimental setup

This section provides a detailed description of the equipment used during the experimental activity, as well as the setup adopted for outdoor firearm noise measurements. It also outlines the associated uncertainties and limitations relevant to the measurement process.

2.1 Measurement equipment

The developed measurement chain is capable of performing multi-channel continuous sampling. Various transducers, mentioned in Table 1, gather the raw analog signal; their application is contingent upon the kind of firearm being tested and the distance from the noise source. The microphones are always equipped with a windscreen GRAS AM0071 that effectively reduces spurious noise, therefore enhancing outdoor measurement performance, albeit with a slight alteration to the microphone’s frequency response. The gathered signal, after proper amplification, is conditioned and transformed into discretized packets of data. This task is performed thanks to the DAQ (Data AcQuisition system) comprised of PXI-4472/4462 or PXIe-4492 sound modules with a PXI-8106 or PXIe-8861 controller embedded in a PXI-1042Q or PXIe-1082 chasis. The output data, in the form of time-domain pressure histories, undergoes various corrections, including calibration of the microphone’s sensitivity, adjustment to free-field response, background noise compensation, and windscreen effects compensation. Given the absence of a power grid during on-field measurements, the electrical equipment is powered by a portable generator (Honda EU22i, 2.2 kVA).

Since the location of the noise source is known and it can be assumed to be a point source emitting spherically diverging sound waves [14], the microphones can be pointed towards it (free-field type, 0°) or perpendicularly (pressure type, 90°). The low sensitivity of the 1/4^″ pressure microphones and the pencil probe becomes irrelevant when the sensor is located in the vicinity of a firearm, where a high upper end of the dynamic range is the driving requirement. In the present work, free-field microphones are deployed where lower sound pressure levels are expected to avoid the clipping of the signal due to excessive amplitudes and to capitalize on their high sensitivity, which enables them to capture smaller-scale disturbances. The sampling frequency f_s used throughout the present work is 70 [kHz].

2.1.1 Weather measurements

Temperature (T_am [°C]), pressure (p_am [Pa]), humidity (RH %), and wind were evaluated using various instruments. These included real-time operated weather stations, such as the ultrasonic anemometer Vaisala WMT702, which measures the wind speed and direction with a f = 1 [Hz]. According to ISO 20906:2009 [15], the ultrasonic anemometer was positioned at least 1 [m] below the microphone and at least 1.5 [m] horizontally distant from the microphone support mast. Laboratory tests provided evidence that the ultrasonic interference is minimal if the ISO 20906 guidelines are followed. When the anemometer was unavailable, the nearest local weather station served the purpose. In the context of this study, weather data is not used to characterize the instantaneous atmospheric conditions during each individual measurement but rather to assess broader categories of atmospheric conditions typically encountered during outdoor firearm noise measurements. These generalized conditions are then used for the subsequent numerical analyses. For the particular weather conditions, the reader is referred to those published in [13], since the experimental data presented here is either identical or gathered during the same measurement activity.

2.2 Setup

During the experiments, two distinct configurations were used to examine various facets of the muzzle blast sound’s propagation. The particular distances and the type of microphone utilized varied depending on factors such as the terrain at the shooting location and safety issues. The blueprints for all of the settings are based on the two setups that are shown here.

2.2.1 Muzzle blast propagation

In order to study the evolution of the muzzle blast as it propagates, an array of sensors is placed on a straight line having an azimuthal direction θ with respect to the line of fire of the firearm (Fig. 1). The choice of this direction is heavily influenced by the topography of the terrain where the shooting range is located, with care to maximize the flatness of the ground, avoid obstacles or reflective surfaces, enforce the safety of the shooting range, and avoid instrument damage. The microphones are mounted on tripods, typically set at human height (h_r = 1.5 [m] [9, 16]), while the height of the muzzle with respect to the ground h_s varies depending on the type of firearm (e.g., vertical firing angle, shooting position). The source-receiver distance d_i, where i is the index associated with the receiver, is measured along the propagation line between the intersections of the vertical projections of the gun muzzle and of the microphone tip with the ground surface. The distance of the closest receiver from the muzzle d₀ is determined on the basis of two criteria. In practice, d₀ falls anywhere between 2 [m] and 15 [m] depending on the weapon system and dynamic range. The position of the furthest receiver is practically constrained by the considerations on the maximum cable length and is, at best, short of 300 [m]. The setup here described is the same as that adopted in [13].

Figure 1. Outline of the standard experimental configuration for the study of muzzle blast propagation.

2.2.2 Directivity pattern setup

A second setup was adopted (Fig. 2), where the sensors are positioned at a fixed distance d from the receiver along a semicircle, with the muzzle at the center and equal angular spacing between them. Studying the impact of wind on the muzzle blast sound’s propagation is the aim of this setup. Assuming identical measuring conditions, a given microphone should record the same noise level each time; nevertheless, because wind might change from shot to shot, the microphones will record varying levels, with each one experiencing a different wind direction as large as 150°. The radius d was chosen with the same criteria used to obtain d₀ (Sect. 2.2.1). The angular spacing θ₀ typically ranges from 22.5° to 30° and a minimum offset with respect to the line of fire (θ₀ >  15°) is mandatory.

Figure 2. Outline of the standard experimental configuration for the characterization of the directivity pattern.

3 Experimental uncertainty quantification

Any measured quantity has a margin of doubt associated with it, and in the case of environmental noise measurements, the factors affecting the source and thepropagation path usually influence the measurement uncertainty more than instrumentation shortfalls [10]. Traditionally, an error is viewed as having two components: a random component, or Type A uncertainty (Sect. 3.2), and a systematic component, or Type B uncertainty (Sect. 3.1). Regardless of the error type, the most frequently encountered sources of uncertainty in environmental noise measurements can be grouped into three categories (Tab. 2). The data showcased in Figure 3 was captured during a measurement activity where a particularly strong wind was present and demonstrates that, even within a recording period of nearly 3 min, the wind direction sweeps an almost 30° wide span, which can be seen by the growing trend represented by the red line. In this specific case, taking the mean of the recording, represented by the black straight line, is not representative of the growing trend. A smaller window, centered around the actual impulsive event, may prove more effective in capturing the instantaneous effect of the wind.

Figure 3. Example of wind speed and direction (North is 0°) readings at 1 [Hz] performed by a Vaisala WMT702 ultrasonic anemometer. The black straight line is the mean.

3.1 Systematic error in the measurement chain

The Type B uncertainty is systematic, does not change when measurements are repeated (e.g., the measuring equipment is always the same, and therefore the uncertainty tied to it is always the same), and its assessments are generally based on best-guess empirical estimates, published data, calibration certificates, or manufacturers’ data. The quantification of the error relies on an estimate of the upper and lower limits (±x, as shown in Tab. 1), for which the standard uncertainty associated with a rectangular distribution is $u=\frac{x}{\sqrt{3}}$. The expanded uncertainty (U = k u) is finally obtained by associating a coverage factor k accounting for the desired confidence interval. Table 3 includes some of the typical sources of systematic error in the measurement chain. The uncertainty of 3 [dB] associated with the environmental factors is a conservative guess that accounts for the statistical spread associated with long-term changes, otherwise difficult to quantify (season, time of the day, weather conditions). This value is consistent with the expected uncertainty associated with the prediction method NT ACOU 099 [17].

3.2 Random error and statistical relevance of the sample size

Type A uncertainty varies unpredictably from one measurement to another and is linked to short-term changes in the environmental conditions and to poor repeatability associated with human factors. Regarding environmental conditions, they differ from environmental factors cited in the type B errors, as the first ones refer to instantaneous variations, such as local gusts of wind, whereas the second ones are associated with long-term or slowly varying influences, such as seasonal changes. Random error can therefore be referred to as shot-to-shot variation, and its magnitude can be statistically estimated by repeating measurements under the sameconditions. The present section is an attempt to best evaluate the random uncertainty with the available data at hand. Let N be the number of experimentally observed sound levels L _i , …, L _N , repeated under the same conditions (Tabs. 4 and 5). As stated in Section 2.2.1, the data presented in Table 4 are the same exhibited in [13] whereas the ones in Table 5 are new and were obtained with the setup described in Section 2.2.2.

The Type A uncertainty is quantified according to the methodology recommended by G. Taraldsen et al. [18], which consists of a translation of the relative uncertainty on the energy scale into an uncertainty on the decibel scale. The experimentally observed sound pressure level is obtained by superposition of an energy mean value $L_{\overline{W}}$ and of an expanded uncertainty k u _F :

$$\begin{array}{c}\hfill L=L_{\overline{W}}\pm k{u}_F[\mathrm{dB}].\end{array}$$(1)

Let W _i  = W₀10 ^{L _i /10}, where W₀ is a reference value (in this instance being W₀ = 20 [μPa]), and $\overline{W}=\frac{1}{N}\sum _i^N{W}_i$. The energy mean, as defined by G. Taraldsen et al. [18], is the mean sound power level and is calculated as follows:

$$\begin{array}{c}\hfill L_{\overline{W}}=10log(\overline{W}/W_0)[\mathrm{dB}].\end{array}$$(2)

This value differs from the mean $\overline{L}=(L_i,\dots ,L_N)/N$. As a first step to obtain the expanded uncertainty, the standard uncertainty is approximated by the Standard Error of Mean (SEM) of the distribution:

$$\begin{array}{c}\hfill \mathrm{SEM}=\mathrm{SD}/\sqrt{N}[\mathrm{Pa}]\end{array}$$(3)

where ${\mathrm{SD}}^2=\frac{1}{(N-1)}\sum _i^N{(W_i-\overline{W})}^2$ is the square of the standard deviation. While the SD is a descriptive parameter that measures the variability of data by expressing the dispersion of individual observations about the mean, the SEM is an inferential predictive parameter that estimates the variability of the possible values of means of samples. Although opposed by some authors [19], the adoption of SEM is motivated by its use as a way to compute confidence intervals that quantify the precision with which the sample mean $L_{\overline{W}}$ estimates the true population mean [20]. The final measurement uncertainty of $L_{\overline{W}}$ is:

$$\begin{array}{c}\hfill u_F\left(L_{\overline{W}}\right)\equiv 10log\{1+\frac{\mathrm{SEM}}{\overline{W}}\}[\mathrm{dB}].\end{array}$$(4)

The desired expanded uncertainty (U = k u _F ) associated to equation (4) is obtained by assigning the appropriate coverage factor k. Considering the use of samples limited in size in the present work, the student’s t distribution is a safer assumption than a Gaussian distribution, when N< 30. For example, if a 95% confidence interval is sought for a sample of size N = 10, the coverage factor is k = 2.26. This methodology can be indistinctly applied to both L_peak and SEL.

3.2.1 Results of type A uncertainty

Peak sound pressure levelL_peak

Figures 5 and 4 display the random experimental expanded uncertainty calculated for L_peak and allow to examine the effects of the source-to-receiver distance, angle, and the sample size. The variability generally grows higher after the first 50 [m], but more often erratically than monotonically. The reason for this behavior is because environmental effects tend to have a limited influence in the very near field of the firearm due to the dominant role of the direct shockwave and the relatively short propagation distance.

Figure 4. Type A expanded uncertainty associated with the Lpeak, as a function of the source-to-receiver angle θ i . The error band is calculated with a 95% confidence interval student’s t distribution, associated with a sample size N from Table 5.

Figure 5. Type A expanded uncertainty associated with the Lpeak, as a function of the source-to-receiver distance di. The error band is calculated with a 95% confidence interval student’s t distribution, associated with a sample size N from Table 4.

The propagation distance introduces shot-to-shot variability due to atmospheric turbulence and local unsteady variations in the wind profile, which can cause fluctuations in the wavefront and effectively However, the expanded uncertainty never exceeds ±2 [dB] for all the firearms under scrutiny.

Sound exposure level SEL Figure 6 displays the random experimental expanded uncertainty calculated for each separate frequency band and allows to examine the effects of source-to-receiver distance, angle, sample size, and frequency range. Some frequency bands are more affected by statistical variability than others. The range between 250 [Hz] and 800 [Hz] is where the effects of finite-impedance ground absorption are most prominent [22]. This ground-associated drop corresponds to a greater variability in the measurements in that frequency region, as can be seen in the case of the DF30 30 mm in Figure 6b and the LAW M-72 66 mm in Figure 6c. Unsteady refractive effects and turbulence are the main causes of the shot-to-shot variability of the pressure waves’ angle of incidence. The former alters the way the pressure waves interact with the ground, while the latter lessens the coherence between direct and ground-reflected sound components. As a consequence, different destructive and constructive interference patterns yield different sound exposure levels at the receiver. In particular, the unsteady nature of the wind-driven refractive effects can be correlated with volatility in the higher frequency bands. A plausible reason why the high end of the frequency spectrum is subjected to a large uncertainty will be later discussed in Section 5.3. This phenomenon is particularly relevant in the case of the setup described in Section 2.2.2, as can be seen in Figures 6e and 6f, where the changing angle between the mean wind profile and the transmission path that links all the aligned receivers causes uncertainties in the high frequencies. The majority of the observed variability falls within the range indicated by NT ACOU 099 [17] (3 [dB]), even though there are higher uncertainty at bigger distances (d >  100 [m]) and with smaller sample sizes (N <  10). These observations must always be kept in mind when comparing experimental results to numerical predictions in order to put into perspective the accuracy of the numerical solver.

Figure 6. Type A expanded uncertainty associated with one-third octave band SEL, at different source-to-receiver distances d i (a–d) and azimuthal angles θ i (e, f). The error band is calculated with a 95% confidence interval student’s t distribution, associated with a sample of size N.

3.3 Overall uncertainty budget

Since the major possible error sources are now identified and estimated, it is possible to compile an uncertainty budget by considering each separate contribution to the uncertainty chain. All the sources are combined and multiplied by a coverage factor k that assigns a confidence interval. To illustrate the procedure, the case of the measurement of the muzzle blast peak overpressure of the M2 Browning (FN M2HB QCB) .50 caliber machine gun is taken as an example.

The Type B uncertainties are quantified in Table 3, and the Type A uncertainty has already been calculated and shown in Figure 5. If the desired confidence for a sample of size N = 10 is 95%, the coverage factor k for a student’s t distribution is k = 2.26. For example, the overall expanded uncertainty for the measured L_peak at d = 23.8 [m] is calculated as follows:

$$\begin{array}{c}\hfill U_{\mathrm{tot}}=U_{\mathrm{A}}+U_{\mathrm{B}}=k(u_{\mathrm{A}}+u_{\mathrm{B}})=\pm 5.6[\mathrm{dB}].\end{array}$$(5)

If U_tot is computed for all receivers (Fig. 7), the true value of the mean is expected to fall, with a 95% confidence level, within the dashed-delimited band ($U_{\mathrm{tot}}^{\max }=\pm 6.3$ [dB]). In this case, the component that dominates the expanded error is the contribution from the Type B uncertainty. However, if the comparison between the two sources of error is conducted with a different weapon, like the DF30 mm, which shows an overall bigger Type A uncertainty (Fig. 4), the contribution to the overall uncertainty budget would be similar between the two types of errors. Furthermore, the type B uncertainty value is greatly influenced by the 3 [dB] component connected with environmental sources, which is based on empirical and cautious assumptions.

Figure 7. Contribution of Type A and Type B errors to the overall expanded uncertainty associated with the measured decay of the muzzle blast peak overpressure of the M2 Browning .50 caliber.

3.4 Background noise

Among the external elements of perturbation when performing sound measurements, the totality of undesired noise picked up by the microphone falls under the classification of background noise. Regardless of the type of source (e.g., flow-induced noise, generic spurious noise sources, etc.), the approach consists of evaluating the net effect on the measured signal and, whenever possible, compensating for it. In addition, using a windscreen helps reduce the ambient noise. A sample measurement of the ambient noise, typically 30 [s], is recorded by each of the receivers during the shooting session. To make sure that the weapon-borne sound stands out from the external ambient noise, the total sound level in each one-third octave band is compared with that of the sampled background noise. If the L_E of the impulsive event exceeds the L_eq of the background noise by at least 10 [dB] [24], the effect of the latter can be ignored, whereas for less than 3 [dB], the specific band needs to be discarded. For intermediate values, the corrections proposed by A.P. Peterson [23] are applied (Tab. 6).

The spectra measured far from the source (Fig. 8) are most likely to show low sound exposure levels that compete on equal terms with the background noise, especially in the low and high ends of the frequency range.

Figure 8. Example of background noise correction of the measured one-third octave band spectrum of an M2 Browning .50 caliber.

As a result, the tails of the spectrum often need to be corrected or discarded. While low frequencies are of secondary importance for small-caliber weapons, they contain a large quantity of acoustic energy in large-caliber weapons. The disposal of frequencies below 50 [Hz] may therefore complicate the spectral analysis at propagation distances where the muzzle blast sound levels are attenuated to the level of the background noise. The correction or removal of frequencies above 12 [kHz] is not uncommon farther away from the source. Although the investigation on nonlinearity concerns the high-frequency end of the spectrum, the effects can already be observed starting from 2 [kHz]. Understanding these dynamics contributes to more accurate evaluation of experimentaldata.

4 Comparison between experimental uncertainty and discretization errors

This analysis is included to demonstrate that numerical discretization errors are minor compared to experimental uncertainties, thereby validating the suitability of the numerical model for investigating environmental effects.

The numerical model implemented by G. Billot et al. [13] will be used to conduct this comparison. The model is based on the 2D Cartesian nonlinear progressive-wave equation (NPE), originally developed by B.E. McDonald and W. Kuperman [25], which has been modified to a cylindrical version that can properly handle the atmospheric absorption [26], resulting in the following equation:

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill D_{t}R=& -\frac{c_{0}R}{2r}+\frac{1}{2}{\zeta }_{\mathrm{th}}({\partial }_r^2+\frac{1}{r}{\partial }_r)R\hfill \\ \hfill & -{\partial }_r(c_{1}R+c_0\frac{\beta }{2}{R}^2)-\frac{c_0}{2}{\partial }_z^2{\int }_{r_f}^{r}R\mathrm{d}r\hfill \end{array}\end{array}$$(6)

where R = p′/ρ₀ c₀ ² represents the non-dimensional pressure perturbation (p′ [Pa] is the acoustic pressureperturbation with respect to the ambient value, ρ₀ = 1.225 [kg/m³] is the ambient density, and c₀ = 343 [m/s] is the speed of sound), ζ_th [m²/s] the sound diffusivity,c₁ [m/s] the environmental component of the speed of associated to refractive effects of meteorological nature such as temperature and wind driven vertical gradients, and β = γ + 1/2, where γ is the heat capacity ratio. The first term on the right side of the equation is linked to the geometrical spreading, the second one to the atmospheric absorption, the third one to non-linearity and refraction, and the fourth one to the diffraction. The model does not include turbulence effects. As was shown in [13], when analyzing an analytical waveform based on the Howitzer signature blast, and obtained with the method outlined by Rigby et al. [32], increasing grid spacing from Δx₁ = 2.45 [mm] to Δx₂ = 9.8 [mm] results in a loss of peak amplitude and a worse resolution of sharp gradients, with the largest error being 0.4 [dB] in the higher frequencies, and a longer front shock rise time. Hence, it can be concluded that the discretization errors are far less relevant in comparison to the ones caused by experimental uncertainty. This observation can justify the use of this method to further inquire into the impact of the environment on the experimental results, as the errors connected with the method itself are much smaller than the uncertainties that will be investigated. Furthermore, while the model does not offer a fully accurate representation of real-world conditions, it provides a controlled and consistent framework to isolate and analyze specific environmental effects – an approach that would be challenging or unfeasible using experimental data alone due to their inherent variability and complexity.

5 Expected influence of environmental conditions

To better understand the environmental effects of wind, temperature, and finite-impedance ground on the propagation of sound, the numerical model described in Section 4 will be used. However, instead of the rigid ground conditions originally used in [13], a more refined one will be implemented: the Porous Ground. By modifying the certain parameters of the equation, the listed environmental effects will be analyzed one at a time. Compared to the model proposed by Salomons et. al. [27] in a recent study, this one can handle non-linear propagation and treats ground reflection and atmospheric refraction separately with a more refined technique. It is important to note that a direct comparison between the computational simulations and experimental measurements was not performed in this study. The main objective of incorporating the updated NPE code is not to replicate specific measurement results, but to create a controlled computational environment that allows for isolating and analyzing the influence of environmental factors, such as ground impedance and atmospheric conditions, which are difficult to decouple in real-world scenarios. Additionally, the firearms have only been tested on one kind of ground, which severely restricts a thorough comparison between the trials and the model.

5.1 Porous ground for the 2D NPE

With this approach, the ground surface is treated like an isotropic fluid medium capable of transmitting dilatational waves in any direction as the result of an impinging wave on the surface. Therefore, the ground layer will be modeled as a rigid external frame containing an air-filled granular material described by three parameters: Flow Resistivity σ₀ [kPa s/m²], Porosity Ω₀, and Tortuosity Φ. Given the lower energy content of the firearms considered in this study, the assumption of a rigid external ground frame is considered appropriate, despite previous findings on ground elasticity effects for stronger blasts [6]. This approach has already been used in the past in different scenarios [28–31], leading to an altered NPE, the derivation of which is available in [28]:

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill D_{\mathrm{t}}R=& -\frac{c_0}{\sqrt{\mathrm{\Phi }}}{\partial }_x[(1-\sqrt{\mathrm{\Phi }})R+\frac{\beta }{2}{R}^2]\hfill \\ \hfill & -\frac{c_0}{2\sqrt{\mathrm{\Phi }}}\int {\partial }_z^{2}R\mathrm{d}x-\frac{{\sigma }_0{\mathrm{\Omega }}_0}{2\mathrm{\Phi }{\rho }_0}(1+\frac{\xi c_0}{\sqrt{\mathrm{\Phi }}}|R|)R.\hfill \end{array}\end{array}$$(7)

The ground layer is hence treated as a separate domain where sound propagation is described by equation (7), whereas the interaction between this domain and the external one is accounted for by the boundary conditions. The correction term $(1-\sqrt{\mathrm{\Phi }})$ replaces the refractive perturbations c₁, and the sound speed within the porous ground is reduced by the tortuosity ($c_{\mathrm{ground}}=c_0/\sqrt{\mathrm{\Phi }}$ [m/s]). Consequently, if Φ = 1 and losses through the medium are neglected (σ = 0 [kPa s/m²]), the equation reduces to the atmospheric NPE. The Forchheimer’s nonlinearity is represented by ξ [s/m⁻¹] and is relevant only when dealing with significant gradients of the particle velocity in the ground; therefore, it can be discarded since the correct modelling of the wavelength within the ground is outside the scope of this research. Details regarding the implementation of the Porous Ground boundary condition within the 2D NPE model are provided in Appendix A, where the specific treatment and numerical approach are discussed in full.

5.2 Ground effects

To simulate the effects of the ground, an analytical source based on the Howitzer signature blast was developed using the methodology outlined by Rigby et al. [32], resulting in a blast-like impulsive source of 4 [kPa] peak amplitude and 1.5 [ms] signal duration located at z = 3 [m]. An array of virtual receivers is positioned at a fixed height z = 1.5 [m] and at x = [5, 15, 30, 50, 95] [m]. The size of the atmospheric moving window is (3 × 18) [m] (the first value refers to the x direction, the second to the z), while the ground layer is (3 × 2) [m]. Each cell of the rectangular structured mesh has an aspect ratio of 3 (Δx = Δz/3 = 3.4 [mm]) and the time step is Δt = Δx c₀ = 0.01[ms]. Four different ground conditions are considered: asphalt (σ₀ = 200000[kPa s/m²], Φ = 100, Ω₀ = 0.01), grass (σ₀ = 200[kPa s/m²], Φ = 2.22, Ω₀ = 0.45), snow (σ₀ = 10[kPa s/m²], Φ = 1.43, Ω₀ = 0.07), and grass up to 10 [m] and asphalt thereafter. The values considered to define these ground types can be found in [2, 33]. The numbers used specifically for asphalt pertain to hot-rolled asphalt, which is effectively a non-porous surface similar to a rigid condition. It is chosen to better understand the difference with more porous terrains, despite the fact that it is an unlikely case. Within the parameters of this paper, the results would not be substantially impacted by the adoption of a more realistic asphalt condition. No wind or temperature refraction (c₁ = 0 [m/s]) effects are included, and only the linear flow resistivity is retained.

5.2.1 Wavefront

In all cases (Fig. 9), the pressure waves are convected towards the right of the computational window, owing to the linear dependency of the particle velocity with amplitude (u ≈ p′/ρ₀c₀) that underpins the nonlinear steepening effects. The reflected wave initially forms an angle ≈45° to then achieve an almost complete superposition with the main wave at x = 95 [m]. The presence of grass (Fig. 9b) and snow (Fig. 9c) modifies the pressure distribution compared to the classic blast signature retained by asphalt (Fig. 9a) and the combination of grass and asphalt (Fig. 9d).

Figure 9. Color map of the 2D pressure distribution for different surfaces: (a) Asphalt, (b) Grass, (c) Snow, and (d) Grass up to 10 [m], then asphalt.

5.2.2 Spectrum

The ground effects are also obvious in the frequency spectrum (Fig. 10), particularly in the mid-frequency (250 [Hz]–1.2 [kHz]), and in the shift of the primary attenuation from 630 [Hz] at x = 30 [m] to 2 [kHz] at x = 95 [m]. To isolate the frequency-dependent effect, a free-field solution (σ₀ = 0, Φ = 0, Ω₀ = 0) is imposed, and in Figure 11 the excess attenuation for each one-third octave band f_i is shown. The findings are comparable to those in [27], however, in this case, the extra attenuation shifts to higher frequencies with increasing distance since non-linearity is incorporated into the propagation model. The asphalt behaves similarly to a perfectly rigid surface, with the reflection that doubles the pressure of the incident wave, resulting in an increase of 6 [dB] at lower frequencies, whereas the case of grass followed by asphalt diverges from this trend after x = 15 [m], where the delay anticipated by the time domain signal corresponds to a spectrum that is shifted out of phase with respect to the asphalt scenario. Finally, the excess attenuation manifests in the form of an alternating sequence of constructive and destructive interferences, sometimes referred to as the comb filter effect.

Figure 10. One-third octave band sound exposure levels for an analytical waveform interacting with different surfaces.

Figure 11. Excess attenuation with respect to the free-field solution for an analytical waveform interacting with different surfaces.

5.2.3 Amplitude

Starting at the peak amplitude when the waves are separate (x = 5 [m], Fig. 12), the ground reflection, at about half the amplitude of the main peak for asphalt, approaches the main wave. As they merge, constructive interference amplifies the reflection, making it the largest peak between 15 [m] and 30 [m]. By x = 95 [m], the reflection integrates fully into a single lobe.

Figure 12. Time-domain pressure histories for an analytical waveform interacting with different surfaces.

The non-linear effects, which are more prominent

in higher intensity signals, are the reason why the main peak attenuates more

(−10.6 [dB]) than the one of the reflected wave (−6.6 [dB]) between x = 5 [m] and x = 30 [m]. During the assimilation, the

duration of the positive lobe increases to a maximum to then be restored to a value similar to the initial one. In the case of the grass layer, the merging of the waves never fully happens, the zero-crossing is delayed, and the recovery to the static pressure takes longer. Finally, the snow case is out of phase compared to the others, and the negative phase is dampened to the point that, after 30 [m], the waveform shape hardly resembles that of a blast wave.

The impact of ground type on blast-like noise propagation, particularly for sound-dampening terrains (snow, grass), revealed significant signal attenuations that shift from the mid (≈630 [Hz]) to the high (≈2 [kHz]) frequencies with increasing distance. Attenuations of this magnitude can approach 10 [dB] at x >  50 [m]. Except in a few cases, the ground impact is typically higher when compared to those associated with experimental uncertainty (for instance, the overall uncertainty budget for the M2 Browning U_tot = ±5.6 [dB]).

5.3 Wind and temperature

Wind and temperature influence the transmission path by curving the otherwise straight trajectory of the sound waves. If vertical profiles of wind speed, direction, and temperature are measured with a meteorological mast, the effective sound speed profiles needed for the NPE model can be calculated directly with curve fitting methods. In this work, the Harmonoise engineering technique [34] will be employed, which is based on the Monin-Obukhov similarity theory [35]. The atmosphere is therefore assumed to be stratified but horizontally homogeneous and steady, which makes it possible to express the sound speed c(z) as the sum of the adiabatic sound speed c_ad(z) and the wind speed component u(z):

$$\begin{array}{c}\hfill \begin{array}{cc}& c_{\mathrm{ad}}(x)\approx c_0+a_c\ln (1+z/z_0)+b_{c}z\hfill \\ \hfill & u(z)\approx a_u\ln (1+z/z_0)+b_{u}z.\hfill \end{array}\end{array}$$(8)

In equation (8), c₀ denotes the adiabatic sound speed at ground level, and a _c , a _u , b _c , b_u are empirical coefficients that describe how sound speed and wind speed vary with height. The variable z₀ represents the aerodynamic roughness length, and z-pagination is the vertical height above the ground. Once a range of values of wind speed at z = 10 [m] and a stability class (selected from Tab. B.3 in the Appendix B) are chosen, it is possible to determine the surface-layer scaling parameters (friction velocity u_*, Monin-Obukhov length L, and temperature scale T_*) and use them to approximate the coefficients (a_c, a_u, b_c, b_u) with the Businger-Paulson profiles [34].

In Table 7 is reported the set of representative cases that will be modeled, where the value corresponding to each of the classes can be seen from the tables in Appendix B (Tabs. B1–B.3). The cases include various combinations of wind direction, wind speed, and atmospheric stability, selected to reflect typical outdoor conditions under which measurements may occur. These seven cases are plotted for a 25 [m] tall domain (Fig. 13). Only the weak downwind (WD) condition qualifies as a favorable measuring condition [9, 10, 17]. In contrast, under a homogenous atmosphere (HA), the temperature gradient bends acoustic waves upward, causing overall lower sound levels near the ground. These refractive profiles are incorporated into the NPE with c₁(z), using the same analytical source and numerical parameters described in Section 5.2. However, since a shift of the acoustic waves within the window is expected, the domain size will vary (3–5 [m] x-direction, and 25–30 [m] z-direction). In order to isolate the influence of wind and temperature, the ground is considered as perfectly rigid. The atmospheric conditions analyzed in this section are best interpreted as instantaneous or short-term variations, and therefore more closely associated with the Type A uncertainties introduced earlier, rather than with long-term or averaged conditions typically linked to Type B uncertainties.

Figure 13. Vertical profiles of refractive sound speed for the conditions described in Table 7.

5.3.1 Peak and minimum sound pressure level

Figure 14 shows the effects of different atmospheric conditions on L_peak and L_p, min, the latter defined as

Figure 14. Evolution with source-receiver distance of Lpeak and Lp, min metrics for an analytical waveform subjected to different atmospheric conditions.

$$\begin{array}{c}\hfill L_{p,min}=20log\left(\frac{|p_{min}|}{2\times {10}^{-5}}\right).\end{array}$$(9)

where p_min is the minimum acoustic pressure. The influence of refraction appears only after 65 [m], with the upwind condition lowering it, whereas the downwind one does the opposite. Also, the discrepancy between weak and moderate downwind cases is little, but it is more noticeable for the two upwind conditions.

5.3.2 Amplitude

By looking at Figure 15, it can be seen that the downwind refractive condition bends the acoustic rays downward towards the ground surface, promoting the formation of multiple peaks coalescing together, and the more intense the downwind condition is, the earlier the wavefront reaches the receiver. In the case of the upwind refraction, the front shock is thicker, and the traveling wave is slowed down.

Figure 15. Time-domain waveforms numerically under the effect of the vertical refractive profiles defined in Table 7.

5.3.3 Spectrum

From Figure 16 it can be noticed that the downwind conditions enhance the sound levels in the low end of the frequency domain, creating a gap of as much as 10 [dB] between the MD and the WU cases at 200 [Hz] and after 200 [m] of propagation. However, the largest contribution by refractive effects to the frequency spectrum can be seen above 500 [Hz]. Since the downward refracting profile counteracts the decrease in the front shock rise time, the refraction serves as a modulator for the waveforms’ steepening. This intensifies the phenomenon of nonlinear steepening, which was shown by G. Billot et al. [13] to persist even at longer distances, and slows the high-frequency spectrum’s decay. For the initial and boundary conditions chosen at the start of the analysis, the refractive effects become tangible earlier than 100 [m] and for frequencies lower than 2 [kHz].

Figure 16. Comparison of the one-third octave band sound exposure levels for an analytical waveform subjected to different atmospheric conditions.

6 Conclusions

The objective of this research was to examine the uncertainty linked with the acoustic measures of outdoor firearm noise propagation. This was achieved combining experimental and numerical approaches, tackling both the experimental limitations given by the equipment and the uncertainties associated with the environmental conditions, and evaluating numerically the impact of different ground surfaces and atmospheric conditions. The uncertainties were classified into systematic ones (Type B) and random errors (Type A). Regarding the L_peak uncertainty, Type B resulted in an expanded uncertainty of U_B = ±4.1 [dB] (Tab. 3 with k = 2.26) for all receivers and guns, while type A (Figs. 5 and 4) ranged from U_A <  ±0.5 [dB] (LAW M-72 66 mm) to U_A = ±2 [dB] (DF30 30 mm) depending on the different firearm as well as the distance of the receiver from the source and his angle. The impact of random errors associated with one-third octave bands SEL (Fig. 6) increases with distance and is most relevant between 250 and 800 [Hz], with the Howitzer 105 mm (Fig. 6d) reaching a max value of U_A = ±10 [dB] at x = 276 [m] and 1 [kHz] while the LAW M-72 66 mm (Fig. 6c) never exceeds U_A = ±2 [dB]. Even though the behavior of the overall uncertainty budget observed for the M2 Browning .50 caliber (Fig. 7) suggests that the systematic errors are more relevant, the value is significantly affected by the 3 [dB] component associated with environmental sources, which is based on empirical and precautionary assumptions. In addition, by considering a different gun, the result would differ. Moreover, at long distances from the source and in contained frequency ranges, the random error is the most significant one. This implies that minimizing shot-to-shot variability can substantially reduce data fluctuations and improve measurement consistency. Automated firing mechanisms, enhanced environmental monitoring, and optimized measurement setups could help minimize variability.

An assessment of the experimental uncertainties alongside the discretization errors associated with the 2D Nonlinear Progressive Wave Equation (NPE) model, as used by G. Billot et al. [13], suggests that the former are the primary contributors to the overall error budget. While the discretization errors are minimal (around 0.4 [dB]), the random environmental factors introduce larger variability, particularly at longer propagation distances and higher frequencies. This highlights the critical need to address environmental randomness to improve the reliability of measurements. The subsequent implementation of a detailed ground model showed the heavy impact of the ground type on the blast-like noise propagation, particularly for sound-dampening terrains (snow, grass) in both the amplitude and frequency content of the signal, resulting in strong signal attenuations that shift from the mid (≈630 [Hz]) to the high (≈2 [kHz]) frequencies with increasing distance, and that can be up to 10 [dB] at x >  50 [m]. This underlines the need to always take into account the surface on which the measurements are conducted, as it can strongly affect the results, even more so than the uncertainties underlined in the previous chapters. Regarding the impact of wind and temperature, it becomes more relevant the farther away from the noise source, and can lead to completely opposite results in the case of upwind or downwind conditions. In contrast, the difference between weak and moderate wind conditions for the same direction has a smaller effect. Accurate identification of wind direction is therefore a primary consideration when assessing the atmospheric conditions during outdoor noise measurements. Furthermore, the variability observed in the numerical simulations is consistent with the shot-to-shot variability identified in the experimental uncertainty analysis. Ultimately, this study highlights the value of addressing uncertainty in firearm noise propagation from two complementary angles: experimental and numerical. While the two approaches were not directly coupled, their combined use provides a more comprehensive understanding of the challenges involved in accurately assessing firearm noise in realistic conditions. This dual methodology supports the development of more reliable prediction tools and better-informed regulatory and mitigation strategies for both military and civilian applications.

Acknowledgments

The authors wish to thank and acknowledge the contribution of the Artillery Battalion of the Belgian Defence, which made the firing ranges accessible to us; the Material Evaluation Center of the Directorate General Material Resources; W. Deweerdt, M. Van Cautter, J. Grossen, J. Weckx and A. Vanhove of the Royal Military Academy for their support; the master interns who contributed to the measurement campaign and the post-processing: S. Deschutter from the Vrije Universiteit Brussel (Belgium); and P. Aubret and S. Benet from the Ecole speciale militaire de Saint-Cyr (France).

Funding

The present work was supported by the scientific research funding of the Royal Higher Institute for Defence under Grant No. MSP/21-02 (MilSOund).

Conflicts of interest

The authors hereby declare no conflict of interest to disclose.

Data availability statement

Data are available on request from the authors.

References

Appendix A Porous ground implementation on the 2D NPE

The implementation will be in analogy with the methodology used for the atmospheric NPE [13], but with a smaller number of grid points in the z-direction. To prevent spurious numerical reflection at the bottom interface of the ground, the same absorbing buffer layer operator prescribed in [13] for the upper boundary will be used. The boundary conditions that need to be enforced at the interface are the continuity of pressure (A.1a) and of normal particle speed (A.1b).

$$\begin{array}{cc}\left[R^a\right]& =\left[R^g\right]\end{array}$$(A.1a)

$$\begin{array}{cc}[\sqrt{\Phi }{\partial }_z{R}^a-\frac{{\sigma }_0{\Omega }_0}{{\rho }_0{c}_0}\int {\partial }_z{R}^{a}dx]& =\left[{\Omega }_0{\partial }_z{R}^g\right]\end{array}$$(A.1b)

Superscripts a and g respectively refer to quantities belonging to the atmospheric and the ground domains, whereas the square brackets denote the field quantity evaluated on the interface.

In order to incorporate the continuous boundary interface conditions into the algorithm, equations (A.1) need to be discretized. For this purpose, the interface is placed midway between the grid rows with indices j = 1, where the bottom of the atmospheric computational domain is located, and j = 0 (Fig. A.1).

Figure A.1. Grid illustration at the air-ground interface. Auxiliary non-physical values (empty circles) are generated by the boundary interface conditions (Eq. (A.3)).

For the generic layer l, the following approximations define the pressure and its spatial derivative at the interface:

$$\begin{array}{cc}\left[R^l\right]& =\frac{R_{i,1}^l+R_{i,0}^l}{2}\cdot .\end{array}$$(A.2a)

$$\begin{array}{cc}\left[{\partial }_z{R}^l\right]& =(R_{i,1}^l-R_{i,0}^l)\Delta z^{-1}.\end{array}$$(A.2b)

The two unknown non-physical auxiliary points R_i, 0^a and R_i, 1^g generated by these approximations are obtained by replacing equation (A.2) into equation (A.1):

$$\begin{array}{cc}{R}_{i,0}^a& =\frac{(A-G)}{(A+G)}{R}_{i,1}^a+\frac{2G}{(A+G)}{R}_{i,0}^g+\frac{S}{(A+G)}\\ & +\sum _{m=N_x}^{i+1}(R_{m,1}^a-R_{m,0}^a)\end{array}$$(A.3a)

$$\begin{array}{cc}{R}_{i,1}^g& =\frac{(G-A)}{(A+G)}{R}_{i,0}^g+\frac{2A}{(A+G)}{R}_{i,1}^a+\frac{S}{(A+G)}\\ & +\sum _{m=N_x}^{i+1}(R_{m,1}^a-R_{m,0}^a)\end{array}$$(A.3b)

where

$$\begin{array}{c}A=\sqrt{\Phi }+\frac{1}{2}SG={\Omega }_{0}S=\frac{{\sigma }_0{\Omega }_0\Delta x}{c_0{\rho }_0}\cdot \end{array}$$(A.4)

R_i, 1^a and R_i, 0^g are the field values obtained from the solution of the NPE in its respective media, for the bottom and top rows.

Appendix B Wind and stability classes

All Tables

All Figures

	Figure 1. Outline of the standard experimental configuration for the study of muzzle blast propagation.
In the text

	Figure 2. Outline of the standard experimental configuration for the characterization of the directivity pattern.
In the text

	Figure 3. Example of wind speed and direction (North is 0°) readings at 1 [Hz] performed by a Vaisala WMT702 ultrasonic anemometer. The black straight line is the mean.
In the text

	Figure 4. Type A expanded uncertainty associated with the Lpeak, as a function of the source-to-receiver angle θ i . The error band is calculated with a 95% confidence interval student’s t distribution, associated with a sample size N from Table 5.
In the text

	Figure 5. Type A expanded uncertainty associated with the Lpeak, as a function of the source-to-receiver distance di. The error band is calculated with a 95% confidence interval student’s t distribution, associated with a sample size N from Table 4.
In the text

	Figure 6. Type A expanded uncertainty associated with one-third octave band SEL, at different source-to-receiver distances d i (a–d) and azimuthal angles θ i (e, f). The error band is calculated with a 95% confidence interval student’s t distribution, associated with a sample of size N.
In the text

	Figure 7. Contribution of Type A and Type B errors to the overall expanded uncertainty associated with the measured decay of the muzzle blast peak overpressure of the M2 Browning .50 caliber.
In the text

	Figure 8. Example of background noise correction of the measured one-third octave band spectrum of an M2 Browning .50 caliber.
In the text

	Figure 9. Color map of the 2D pressure distribution for different surfaces: (a) Asphalt, (b) Grass, (c) Snow, and (d) Grass up to 10 [m], then asphalt.
In the text

	Figure 10. One-third octave band sound exposure levels for an analytical waveform interacting with different surfaces.
In the text

	Figure 11. Excess attenuation with respect to the free-field solution for an analytical waveform interacting with different surfaces.
In the text

	Figure 12. Time-domain pressure histories for an analytical waveform interacting with different surfaces.
In the text

	Figure 13. Vertical profiles of refractive sound speed for the conditions described in Table 7.
In the text

	Figure 14. Evolution with source-receiver distance of Lpeak and Lp, min metrics for an analytical waveform subjected to different atmospheric conditions.
In the text

	Figure 15. Time-domain waveforms numerically under the effect of the vertical refractive profiles defined in Table 7.
In the text

	Figure 16. Comparison of the one-third octave band sound exposure levels for an analytical waveform subjected to different atmospheric conditions.
In the text

	Figure A.1. Grid illustration at the air-ground interface. Auxiliary non-physical values (empty circles) are generated by the boundary interface conditions (Eq. (A.3)).
In the text