A method for estimating the yield of explosions from shock-wave and acoustic-wave measurements is presented. The method exploits full waveforms by comparing pressure measurements against an empirical stack of prior observations using scaling laws. The approach can be applied to measurements across a wide-range of source-to-receiver distances. The method is applied to data from two explosion experiments in different regions, leading to mean relative errors in yield estimates of 0.13 using prior data from the same region, and 0.2 when applied to a new region.

Accurate characterization of the explosive yield of human-related or natural blasts provides important information on the size, mechanism, and time evolution of the event. For example, seismoacoustic methods have been used to determine the yield of a gas pipeline explosion (Evers et al., 2007) and to discriminate between events during an attack on a munitions dump in Baghdad (Aleqabi et al., 2016). The explosive yield of volcanic explosions reflects both the fundamental physics of eruptions as well as the type and range of the resulting hazards (e.g., Johnson and Miller, 2014; Spina et al., 2015). In addition, estimating the yield of explosions from shock wave or acoustic wave measurements has important applications for explosion forensics. Current methods for yield estimation that use shock wave or acoustic measurements are based on using one of two parameters: overpressure (e.g., Kinney and Graham, 1985; Koper et al., 2002) or impulse (e.g., Koper et al., 2002; Ford et al., 2014; Kim and Rodgers, 2016). A major challenge of these methods is that these parameters are also sensitive to propagation effects. Static propagation effects include the transition from the nonlinear shock wave to the acoustic wave regime and ground properties (e.g., terrain and acoustic absorption). In addition, acoustic waves are sensitive to meteorological effects over short temporal and spatial scales. These localized variations in atmospheric properties are difficult to capture even with state-of-the-art meteorological models. Thus, even pressure wave propagation models that incorporate terrain effects can be difficult to implement under real-world conditions. In this paper, we describe a yield estimation technique that is applicable over a large range of source-to-receiver distances and works for both shock waves and acoustic waves.

Most practical methods for estimating the explosive yield using blast waves or acoustic waves have been developed for observations over limited source-to-receiver ranges (or more precisely, over limited scaled ranges, where the scaled range adjusts the absolute range to account for differences in yield). For example, the scaling laws defined in Kinney and Graham (1985) are only applicable to shock waves. Semi-empirical models such as the ANSI model (ANSI, 1983) and the blast operational overpressure model (BOOM) (Douglas, 1987) are fit to observational data at very specific scaled ranges in specific locations. Thus, they do not work well for observations outside those ranges and/or locations. Methods that utilize full waveform modeling (Kim and Rodgers, 2016) suffer from problems related to under-parameterization of local atmospheric effects. The high-resolution atmospheric data necessary to predict full waveforms is rarely available in the field, particularly for unexpected events.

The approach introduced in this paper utilizes the full pressure waveform, avoiding the reduction of information content that parametric methods make. The method relies upon recent high-density deployments of acoustic sensors for ground-truth explosion experiments, which are described below. The method is only now possible because of such data-rich experiments, which significantly extend the observational record of pressure data from explosions.

The method described in this paper exploits a rich dataset of shock wave and acoustic observations from explosions of known yield. The method is based on the fact that observations from multiple surface explosions, with known but different yields, can be plotted on a single scaled-time/scaled-range plot (e.g., Fig. 1). Prior to estimating the yield of a new explosion, observations from previous ground-truth experiments are stacked by averaging sampled pressure observations in scaled time and scaled range bins. In this study, we generated stacks by averaging pressure observations from 0 to 700 m/kg1/3 and −0.005 s/kg1/3 to 0.03 s/kg1/3 in bins of 20 m/kg1/3 and 0.0001 s/kg1/3 (Fig. 1). The minimum and maximum scaled times and ranges were chosen to capture the full set of pressure-time measurements from the observational dataset used in this study. The bin dimensions used are a compromise between scaled time and scaled range resolution and noise reduction through averaging. This was achieved through visual inspection by the analyst.

Fig. 1.

Picture of the stack showing raw waveforms for HRR-1, HRR-4, HRR-5, and HRR-6 (left panel), averaged waveforms (center panel), and peak overpressure amplitudes of the averaged waveforms (right panel). Scaled times are plotted relative to analyst-derived arrival time picks. The stack in the middle panel was used to estimate the yield for HRR-3.

Fig. 1.

Picture of the stack showing raw waveforms for HRR-1, HRR-4, HRR-5, and HRR-6 (left panel), averaged waveforms (center panel), and peak overpressure amplitudes of the averaged waveforms (right panel). Scaled times are plotted relative to analyst-derived arrival time picks. The stack in the middle panel was used to estimate the yield for HRR-3.

Close modal

Following Kinney and Graham (1985), scaled time and distance (tsc,rsc) are defined relative to the observed time and distance (t,r), by

(1)
(2)

where W is the explosion yield in kilograms (TNT equivalent), and

(3)
(4)

scale for the observed ambient pressure (Pobs) and temperature (Tobs) to reference pressure Pref=1013.25mbar and temperature Tref=15°C. The use of temperature and pressure in these scaling equations relates to the source term (the efficiency of the acoustic source) and not to path terms (the propagation of the signal).

To estimate the yield of a new explosion, we perform a grid search whereby observed data at a given range are converted to scaled range by assuming a trial yield. This waveform is compared against the corresponding stack record, with data transformed into units of scaled time relative to an arrival time pick. The estimated yield is thus the trial yield that results in the minimum normalized root-mean-square residual over all observations. The residual from M sensors, where each sensor is at a different geographic location, can be defined as

(5)

where p=[p1,p2,...,pn] are scaled-time binned pressure observations at a given station, p¯(W)=[p¯1,p¯2,...,p¯n] is the empirical stack for the scaled range closest to the true scaled range for the trial yield, and pimax and pimin are maximum and minimum values of the pressure at the ith sensor. The factor of (pimaxpimin) normalizes each root-mean-square (RMS) residual to ensure that near-field and far-field observations have the same weight. The residual defined by (1/N)j=1N(pjp¯(W)j)2 in Eq. (5) effectively provides the RMS difference between an observed waveform, converted into scaled-time for a trial yield and averaged in each scaled-time bin, and the average binned reference waveform for the corresponding scaled-range. The total residual in Eq. (5) is thus the mean of these RMS differences, with each RMS difference normalized by absolute amplitude.

The Humming Roadrunner (HRR) experiment was a series of chemical explosions conducted at the White Sands Missile Range (WSMR) in New Mexico (these events were also studied by Kim and Rodgers, 2016, and Bonner et al., 2013 and a more complete description of the dataset is provided in these references). Five shots (Table 1 and Fig. 2) were conducted at or above the ground-surface and provide the primary testing dataset for this study. For each shot, a dense network of overpressure and acoustic sensors were deployed, and digitized at 1000 Hz (Fig. 2). The five shots were conducted in different locations and at different times in a region of mountainous topography and thus include complex meteorological and terrain effects. Notably, HRR-3 occurred at a time when a strong storm moved through the study region (Table 2). The waveforms associated with HRR-3, discussed in detail below, are distinct from the waveforms associated with the other shots.

Table 1.

Comparison between estimated and true yields for each HRR and HTA shot. HTA estimated yields are evaluated using a stack formed from HRR data.

Shot numberTNT equivalent yield (kg)Estimated yield (kg)Estimated yield (Kim and Rodgers)Relative error(this study)Relative error(Kim and Rodgers)
HRR-1 18 144 21 210 17 300 0.17 0.05 
HRR-3 9072 6870 9300 0.24 0.03 
HRR-4 9072 9700 12 100 0.07 0.25 
HRR-5 45 359 49 400 61 000 0.09 0.34 
HRR-6 45 359 49 400 58 800 0.09 0.30 
HTA-1 227 355 — 0.36 — 
HTA-2 227 233 — 0.03 — 
HTA-3 227 175 — 0.23 — 
Shot numberTNT equivalent yield (kg)Estimated yield (kg)Estimated yield (Kim and Rodgers)Relative error(this study)Relative error(Kim and Rodgers)
HRR-1 18 144 21 210 17 300 0.17 0.05 
HRR-3 9072 6870 9300 0.24 0.03 
HRR-4 9072 9700 12 100 0.07 0.25 
HRR-5 45 359 49 400 61 000 0.09 0.34 
HRR-6 45 359 49 400 58 800 0.09 0.30 
HTA-1 227 355 — 0.36 — 
HTA-2 227 233 — 0.03 — 
HTA-3 227 175 — 0.23 — 
Fig. 2.

Maps of the source and sensor locations for the HRR (top) and HTA (bottom) experiments. In the top panel, colored symbols represent sources (cyan circle = HRR-1, red triangle = HRR-3 and HRR-4, blue square = HRR-5 and HRR-6) and open symbols represent sensor locations (circles = HRR-1, triangles = HRR-3 and HRR-4, and squares = HRR-5 and HRR-6). Stations to the NW and SE of HRR-3 used in the analysis presented in Fig. 3 are filled red and yellow, respectively. In the bottom panel, the circle is the source location for each shot, and open circles show sensor locations. In each panel, topography is shaded with a light source.

Fig. 2.

Maps of the source and sensor locations for the HRR (top) and HTA (bottom) experiments. In the top panel, colored symbols represent sources (cyan circle = HRR-1, red triangle = HRR-3 and HRR-4, blue square = HRR-5 and HRR-6) and open symbols represent sensor locations (circles = HRR-1, triangles = HRR-3 and HRR-4, and squares = HRR-5 and HRR-6). Stations to the NW and SE of HRR-3 used in the analysis presented in Fig. 3 are filled red and yellow, respectively. In the bottom panel, the circle is the source location for each shot, and open circles show sensor locations. In each panel, topography is shaded with a light source.

Close modal
Table 2.

Meteorological data acquired at the time of each shot. Wind speeds for HRR shots were taken at 10 m above the ground.

Shot numberWind speed (m/s)Wind direction (deg.)Pressure (mbar)Temperature (C)
HRR-1 2.5 29.7 848.5 31.6 
HRR-3 11.9 121.2 856.8 23.9 
HRR-4 0.3 194.2 856.9 25.2 
HRR-5 3.1 54.2 855.8 32.5 
HRR-6 3.3 188.3 856.5 29.7 
HTA-1 4.5 — 817 3.3 
HTA-2 2.7 — 818 2.6 
HTA-3 0.0 — 813 0.7 
Shot numberWind speed (m/s)Wind direction (deg.)Pressure (mbar)Temperature (C)
HRR-1 2.5 29.7 848.5 31.6 
HRR-3 11.9 121.2 856.8 23.9 
HRR-4 0.3 194.2 856.9 25.2 
HRR-5 3.1 54.2 855.8 32.5 
HRR-6 3.3 188.3 856.5 29.7 
HTA-1 4.5 — 817 3.3 
HTA-2 2.7 — 818 2.6 
HTA-3 0.0 — 813 0.7 

The Humming Tarantula (HTA) experiment was a series of much smaller explosions detonated at the Energetic Materials Research and Testing Center near Socorro, New Mexico, in January and February of 2015 (Fig. 2). The HTA experiment comprised several explosions of different height-of-burst and depth-of-burial and was designed to understand the effects of these parameters on the relative coupling of energy into seismic and acoustic waves. Here, we only use the above-ground shots, HTA1–HTA3, which were conducted at height-of-bursts ranging from 0.6 to 4.8 m. For each shot, a network of overpressure instruments and infrasound sensors were deployed (overpressure were deployed in the near-field to measure shock waves, while infrasound sensors were deployed in the far-field).

The method described here assumes that an explosive event has been detected and located, and that the onset times of shock wave or acoustic wave signals at each station have been determined, since waveforms are aligned by arrival time pick for calculating the residual in Eq. (5). In this study, arrival times have been picked by an analyst, although in principle they could be determined automatically. For the purpose of this short letter, we do not consider the effect of errors associated with the arrival time picks, or of errors associated with event locations, since we assume that arrival time picks are accurately determined by an analyst and that the event location is known from other means.

To assess the method described above on the HRR data, we performed a leave-one-out test whereby we estimated the yield for each of the five HRR explosions in Table 1 using the other four HRR explosions to generate the empirical stack. The stack was then used to estimate the normalized RMS residual for a set of 50 trial yields, logarithmically spaced from 100 to 100 000 kg, with the estimated yield taken as the yield resulting in the minimum residual. The resultant yield estimates, shown in Table 1, result in similar relative errors, |(WestWtrue)/Wtrue|, to the results published by Kim and Rodgers (2016) which was based on finite difference modeling—our yield estimates were closer for three shots, while the results of Kim and Rodgers were closer for two shots.

To further explore the method, we apply the technique separately to observations taken along two profiles to the NW and SE of HRR-3 (the stations used are colored red and yellow in Fig. 2). The observations of this shot exhibited very strong azimuthal effects that are caused by a storm front that passed through the study region at the time of the event. Radiosonde measurements, taken 30 min before the shot, show quite different propagation environments to the NW and SE (Fig. 3). These different propagation environments are manifest in the overpressure observations, which decrease with range to the NW but remain steady to the SE due to a strong directional wind-driven waveguide (Fig. 3). Despite such strong azimuthal differences, our method results in yield estimates of 8697 kg when using only the stations to the NW (a relative error of only 0.04) and 6135 kg when using only the stations to the SE (a relative error of 0.32). For reference, using all the data to estimate yield for this shot results in a relative error of 0.24 (Table 1). As shown in Fig. 3, in addition to having very different overpressures in opposite directions, the stack waveforms and waveforms measured from HRR-3 exhibit strong differences. However, despite there being strong differences in waveforms, the RMS residual between the observed and stack waveforms is minimized when using a trial yield close to the true yield. This suggests that our method is sensitive to the long-wavelength properties of the signal, with high-frequency effects averaged out. Because the method is based on a differential measurement, rather than an absolute measurement (like overpressure), the template matching method continues to give reasonable results for HRR-3, despite the strong propagation effects.

Fig. 3.

Application of the template matching method to HRR-3 observations to the NW (red triangles in Fig. 2) and the SE (yellow triangles in Fig. 2) of the shot. (a) Effective sound speed profiles to the NW (dashed line) and the SE (solid line) from a radiosonde that was launched adjacent to the shot site 30 min prior to launch. (b) Overpressure measurements taken along the NW (circles) and SE (triangles) profiles. (c) A comparison between the observed data (solid lines) and stack data (dashed lines) along the NW profile, with observations scaled by the true yield. (d) As (c) for the SE profile.

Fig. 3.

Application of the template matching method to HRR-3 observations to the NW (red triangles in Fig. 2) and the SE (yellow triangles in Fig. 2) of the shot. (a) Effective sound speed profiles to the NW (dashed line) and the SE (solid line) from a radiosonde that was launched adjacent to the shot site 30 min prior to launch. (b) Overpressure measurements taken along the NW (circles) and SE (triangles) profiles. (c) A comparison between the observed data (solid lines) and stack data (dashed lines) along the NW profile, with observations scaled by the true yield. (d) As (c) for the SE profile.

Close modal

To assess the sensitivity of the method to the number of observations, we performed a Monte Carlo analysis, focusing on HRR-4. By sampling Mi stations at random (where Mi ranges from 3 to 20) from the 68 observations of HRR-4, and repeating the experiment 20 times for each value of Mi, we estimated the mean and standard deviation of the resultant yield estimate from the 20 trials. The results (see supplemental figure1) indicated that the mean yield estimate is insensitive to the number of observing stations because the mean yield estimate shows no systematic trend with the number of observing stations and is consistently higher than the true yield with a relative error ranging from 0.07 to 0.12. However, the standard deviation of the estimates decreases from about 2800 kg with only three stations, to about 1400 kg with 20 stations. These results suggest that the method works well for even small numbers of observations; increasing the number of observations improves the consistency of the result obtained through multiple Monte Carlo trials.

To test the transportability of the method, where a stack formed using experimental data in one region is applied to estimate yields using data from a different region and time of year, we used a stack of HRR shot data to estimate the yields of the HTA shots. The results, shown in Table 1, indicate that the stack generated for the HRR experiment can be successfully transported to estimate yields for the HTA experiment (which has a different set of source-receiver paths with different terrain effects, and occurred at a different time of year with different meteorological effects). This result suggests that stacks formed from datasets in one region, or across multiple regions, may be transportable to new regions. However, it is important to note that the datasets used in this study do not include shots performed under unusual conditions such as temperature inversions, although the analysis of HRR-3 shows strong wind focusing to the SE.

The method outlined in this paper uses the full pressure time-series recorded from explosions to estimate yield. The method can be applied to pressure measurements over a large range of source-to-receiver distances that include shock wave and acoustic wave measurements, thus mitigating a major limitation of existing empirical methods. Because we were able to use an empirical stack developed for set of source-receiver paths during the summertime, and apply it to successfully estimate yields given a completely different set of source-receiver paths during the wintertime, we suggest that the method may be less sensitive to specific source-receiver effects than conventional parameter-based methods. Empirical evidence has found that impulse, a more averaged property of the blast wave, may be less sensitive to atmospheric conditions than instantaneous measurements such as overpressure (Ford et al., 2014). Through a detailed analysis of HRR-3 data, we suggest that using the full waveform, as we do in this paper, further mitigates the effects of meteorology on yield estimates.

We are grateful to Weston Geophysical for providing the measurement data that was taken by multiple people and organizations during the Humming Roadrunner experiment. We thank Bob Reinke at the Defense Threat Reduction Agency and Chris Young at Sandia for helpful discussions. We also thank the editor, Vladimir Ostashev, and anonymous reviewers who provided valuable feedback on an earlier version of this manuscript. Sandia National Laboratories is a multimission laboratory managed and operated by National Technology and Engineering Solutions of Sandia, LLC., a wholly owned subsidiary of Honeywell International, Inc., for the U.S. Department of Energy's National Nuclear Security Administration under Contract No. DE-NA0003525.

1

See supplementary material at http://dx.doi.org/10.1121/1.4984121E-JASMAN-141-506706 for an illustration of the sensitivity of the results to the number of observations.

1.
Aleqabi
,
G.
,
Wysession
,
M.
, and
Ghalib
,
H.
(
2016
). “
Characterization of seismic sources from military operations in urban terrain (MOUT): Examples from Baghdad
,”
Bull. Seism. Soc. Am.
106
,
23
41
.
2.
ANSI
(
1983
). S2.20-1983,
Estimating Airblast Characteristics for Single Point Explosions in Air, With a Guide to Evaluation of Atmospheric Propagation and Effects
(
American National Standards Institute
,
New York
).
3.
Bonner
,
J.
,
Russell
,
D.
, and
Reinke
,
R.
(
2013
). “
Modeling surface waves from aboveground and underground explosions in alluvium and limestone
,”
Bull. Seism. Soc. Am.
103
,
2953
2970
.
4.
Douglas
,
D.
(
1987
). “
Blast operational overpressure model (BOOM): An airblast prediction method
,” Air Force Weapons Laboratory,
Report No. AFWL-TR-85-150
.
5.
Evers
,
L. G.
,
Ceranna
,
L.
,
Haak
,
H. W.
,
Le Pichon
,
A.
, and
Whitaker
,
R.
(
2007
). “
A seismoacoustic analysis of the gas-pipeline explosion near Ghislenghien in Belgium
,”
Bull. Seism. Soc. Am.
97
,
417
425
.
6.
Ford
,
S. R.
,
Rodgers
,
A.
,
Xu
,
H.
,
Templeton
,
D.
,
Harben
,
P.
,
Foxall
,
W.
, and
Reinke
,
R.
(
2014
). “
Partitioning of seismoacoustic energy and estimation of yield and height-of-burst/depth-of-burial for near-surface explosions
,”
Bull. Seism. Soc. Am.
104
,
608
623
.
7.
Johnson
,
J.
, and
Miller
,
A.
(
2014
). “
Application of the monopole source to quantify explosive flux during vulcanian explosions at Sakurajima Volcano (Japan)
,”
Seism. Res. Lett.
85
,
1163
1176
.
8.
Kim
,
K.
, and
Rodgers
,
A.
(
2016
). “
Waveform inversion of acoustic waves for explosion yield estimation
,”
Geophys. Res. Lett.
43
,
6883
6890
, doi:.
9.
Kinney
,
G.
, and
Graham
,
K.
(
1985
).
Explosive Shocks in Air
(
Springer-Verlag
,
New York
).
10.
Koper
,
K.
,
Wallace
,
T.
,
Reinke
,
R.
, and
Leverette
,
J.
(
2002
). “
Empirical scaling laws for truck bomb explosions based on seismic and acoustic data
,”
Bull. Seism. Soc. Am.
92
,
527
542
.
11.
Spina
,
L.
,
Taddeucci
,
J.
,
Cannata
,
A.
,
Gresta
,
S.
,
Lodato
,
L.
,
Privitera
,
E.
,
Scarlato
,
P.
,
Gaeta
,
M.
,
Gaudin
,
D.
, and
Palladino
,
D.
(
2015
). “
Explosive volcanic activity at Mt. Yasur: A characterization of the acoustic events (9–12th July 2011)
,”
J. Volc. Geo. Res.
302
,
24
32
.

Supplementary Material