The 12 December 2017 Baumgarten Gas Hub Explosion: A Case Study on Understanding the Occurrence of a Large Infrasound Azimuth Residual and a Lack of Seismic Observations

The Baumgarten explosion occurred on 12 December 2017 at a gas storage site about 30 km east of Vienna, Austria. Acoustic arrivals from this accidental surface explosion were detected at dozens of stations of the AlpArray seismic network to distances up to 150 km, mainly in easterly directions. Thus it was expected that the Hungarian infrasound array PSZI located about 230 km to the east-southeast of Baumgarten would detect this acoustic wave as well. Standard progressive multichannel correlation processing and frequency-wavenumber analysis identified a signal emerging at 7:57:55 UTC from an azimuth of 296°–300° and with trace-velocity > 400 m/s. The extraordinarily high trace-velocity and excessive backazimuth residual, relative to the explosion site direction of 282°, however cast strong doubts on the arrival’s connection to the Baumgarten event. Accounting for the effect of non-planar geometry of the infrasound array results in a reduction of the azimuth residual by half. Additionally, 2D and 3D raytracing methods are used including the European Centre for Medium-Range Weather Forecasts (ECMWF) atmospheric model to further explain the remaining azimuth residual as well as to elucidate the large trace velocity estimates. The prevailing stratospheric winds in excess of 150 m/s are identified as the underlying cause. Including both factors the initial azimuth residual of up to 18° decreases to ~ 4°, allowing to associate the infrasound signal at PSZI with the Baumgarten event. Finally, the data from a seismic station at 30 km range is re-investigated for magnitude estimation. The local magnitude of ML < 1.0 explains well the scarcity of seismic observations within 50 km range, where three or four stations show signals, mainly consisting of Rg-type surface waves, but no body waves.


Introduction
Accidental explosions are of particular interest in the context of nuclear explosion monitoring of the Comprehensive Nuclear-Test-Ban Treaty (CTBT), as they provide ground truth information. Since the atmosphere is a dynamic and anisotropic medium, propagation effects are highly dependent on location, path, and time (Drob et al. 2003). Ground truth sources make it possible to assess our capabilities in detecting weak infrasound signals, accurately estimating their arrival direction as well as modeling atmospheric sound propagation from source to receiver within the context of a treaty relevant event exhibiting acoustic signatures.
For this reason, a significant number of studies were concerned with such events since the CTBT was opened for signature in 1996. Among these events are the gas pipeline explosion near Ghislenghien in Belgium (Evers et al. 2007), the Buncefield, UK, oil depot explosion (Ottemöller and Evers 2008), the Gerdec, Albania, and Chelopechene, Bulgaria, explosions (Green et al. 2011), and the explosions at a Slovakian ammunition factory (Kristekova et al. 2008). Additional ground truth sources have been rocket engine tests (Pilger et al. 2013;Park et al. 2018), and explosions related to the destruction of old ammunition (Gibbons et al. 2007).
A surface explosion occurred on 12 December 2017 more than 30 km east-northeast of Vienna, Austria, at a gas distribution hub near the village of Baumgarten shortly before 9 am local time (or 08:00 UTC) (see Fig. 1). The geographic location of the gas hub is at 48.32°N and 16.87°E, close to the border between Austria and Slovakia. We also note here that Schneider et al. (2018) cite ground-truth information of this explosion not being related to a typical combustion, but rather an instantaneous release of highly pressurized gas. The closest infrasound arrays are an installation deployed in spring 2017 by the Geodetic and Geophysical Institute (CSFK) at the Piszkes Observatory in Hungary some 230 km east of the explosion site (see Fig. 1) and the IMS station IS26 in Southeast Germany at a distance of about 240 km to the west. Due to the pattern of westerly stratospheric winds during the winter season in Central Europe (Le Pichon et al. 2008;Koch 2010) the expectation of a signal observable at IS26 is very low, but for the station to the east appears rather high.
Seismo-acoustic signals from the Baumgarten explosion were detected at the AlpArray seismic network of permanent as well as temporary stations that have been installed since 2016 (Schneider et al. 2018). The network consists of more than 600 stations with an average station spacing of less than 40 km. Within a distance range of 185 km from the Baumgarten site 68 stations were operational, of which up to 33 stations showed seismo-acoustic arrivals to distances of more than 150 km. These signals, exhibiting acoustic wave speeds with celerities between 270 and 350 m/s, traveled as tropospheric and stratospheric phases to stations to the north, east and south, while for nearly all stations in westerly directions such signals could not be found. Schneider et al. also provided an estimate of the origin time (07:44:16.1 UTC). In this study, we use our origin time estimate of 07:44:20 UTC for analysis determined from acoustic arrivals at four nearby seismic stations and using propagation speeds of 330-350 m/s depending on station azimuth.
We reviewed the data of the infrasound array in Hungary for any signals arriving with acoustic wave speeds and adequate backazimuth. Such a signal could be found shortly before 08:00 UTC, but with a excessively large backazimuth residual and extraordinarily high apparent trace velocity, certainly casting severe doubts on the signal of being associated with the Baumgarten explosion. In addition, this case is further complicated by the fact that a quarry is located along the preliminary backazimuth at a distance of about 40 km. However, telephone conversation with the quarry operator revealed that no blasting had been carried out during the relevant time period. The high trace velocity of the signal was furthermore in disagreement with a source at the quarry; the shape of previous signals from such quarry events were also significantly different from the signal observed on 12 December 2017, both in terms of duration as well as a more gradual signal onset. Map showing the site of the Baumgarten explosion (star) and the locations of the Hungarian infrasound array PSZI and the nearest seismic station MODS in Slovakia (triangles). Also depicted is the city of Vienna and the distance and azimuth to PSZI 4958 K. Koch et al. Pure Appl. Geophys. This study therefore aims to: (1) improve backazimuth estimates from array processing (and provide a reduction of the backazimuth residual) by taking into account local topography (Sect. 2); (2) investigate acoustic wave propagation and its effects on apparent backazimuth and trace velocity (Sect. 3); and (3) derive an estimate of event magnitude to explain the lack of seismic observations throughout the region (Sect. 4). Finally, cross-bearing methods for infrasound source location are applied frequently (e.g. Evers et al. 2007;Landès et al. 2012;Hupe et al. 2019). In addition infrasound arrivals are only used at the International Data Centre of the CTBTO (Comprehensive Test-Ban Treaty Organization) for source location, if a phase is azimuth defining, i.e. shows an azimuth error of less than 10° (Mialle et al. 2019;CTBTO 2012). Hence a major scientific motivation is to present this case of an infrasound arrival with an extraordinary azimuth bias that can indeed be associated to a specific source when including corrections to standard array processing and considering propagation path effects.

Infrasound Observations at PSZI (Hungary) and Results Obtained from Array Processing Techniques
The Piszkes Observatory is located in northern Hungary at latitude 47.92°N and longitude 19.89°E and is therefore 228 km to the east at an azimuth of 100°from the Baumgarten explosion site (Fig. 1). The four infrasound elements were deployed in spring 2017 at distances within about 200 m from the seismic station, with the reference element PSZI1 colocated with seismic station PSZ (Fig. 2a). The individual elements are at elevations between 900 and 950 m due to local topography. This is substantial compared to the overall aperture of the array of about 400 m by 200 m, with the array's baseline dipping towards the west.
The uncertainty of estimates from array processing was calculated using DASE toolkit progressive multichannel correlation (DTK PMCC) (Mialle et al., 2019), as shown in Fig. 2b. This pattern, however, is only valid for a planar setup of the array and does not account for any elevation differences. Due to the array's elongation in the NW-SE direction the azimuth uncertainty is highest in this direction, amounting to up to 10°, and at a minimum in terms of trace velocity error of about 20 m/s. In the perpendicular direction the azimuth uncertainty is estimated at a minimum on the order of 3°-4°, and at a maximum for trace velocity error, reaching up to 60 m/s. For the direction to Baumgarten, i.e. 282°, the azimuth error is larger than 9.5°, while for the tracevelocity the error is below 30 m/s. The data from the PSZI array for the time of the Baumgarten explosion were retrieved from the GEOFON Data Centre (see Acknowledgements) and were analyzed with array processing techniques. The band-pass filtered data (1-5 Hz) for a 2-min window starting at 7:56:30 UTC on 12 December 2017 is shown in Fig. 3a. Using an origin time of the Baumgarten explosion of 07:44:20 UTC this time window corresponds to a celerity window between 320 and 275 m/s. According to Negraru et al. (2010), who assigned different celerity ranges to tropospheric, stratospheric or thermospheric arrivals, this window would most likely contain arrivals from stratospheric propagation paths. Besides a few acoustic disturbances an isolated and clearly impulsive signal at all four acoustic elements is found. The additional seismic channel of the vertical ground motion from station PSZ exhibits a similar arrival. The arrival time of the signal pulse is at 7:57:55 UTC.
Taking an 18 s window around this arrival and performing a frequency wavenumber (FK) analysis (Stammler 1993) we obtain a clear peak in the wavenumber plane at an azimuth of 296°and a corresponding trace (apparent) velocity of 407 m/s (Fig. 3b). This value is substantially higher than expected for a regular tropospheric and/or stratospheric arrival, and would be rather atypical for an acoustic wave generated by a close-by near-surface source such as a blast at the aforementioned quarry. This high trace velocity can only be explained by two scenarios: either (1) unusual atmospheric conditions, generating e.g. an extraordinary strong stratospheric duct for the case of the Baumgarten explosion, or (2) an acoustic source at altitude, e.g. from the sonic boom of an aircraft, with wave incidence from above. Furthermore, for a source in the direction of Vol. 177, (2020) The 12 December 2017 Baumgarten Gas Hub Explosion Baumgarten the determined azimuth implies a residual of 14°, which is very unusual for infrasound arrivals from ground-truth sources (cf. Smets et al. 2015;Pilger et al. 2018;Blixt et al. 2019).
To obtain supporting and potentially improved array processing results for the PSZI data we further apply the PMCC technique (Cansi 1995) to the corresponding time segment, as shown in Fig. 4a. The time segment represents celerity values with respect to the Baumgarten origin time of 275-300 m/s. As a result, we obtain an about 10-15 s long signal with frequencies between 1 and 4-5 Hz. Averaging the azimuth and trace velocity values of all relevant pixels over the signal duration the mean azimuth is determined to 300.5°and the mean trace velocity to 463 m/s. This is more than 50 m/s faster than the estimate from F-K analysis, and the azimuth residual is larger than 18°. Inclined interfaces below seismic array stations (Bondár et al. 1999;Koch and Kradolfer 1997) can cause significant backazimuth bias, and elevation differences were found to produce similar deviations (Edwards and Green 2012) in infrasound beamforming. We further processed the PSZI data with PMCC in a mode, where elevation differences of the array elements are corrected for (Nouvellet et al. 2014). The outcome from this processing is shown in Fig. 4b, where the (mean) azimuth residual is decreased to 10°, and the (mean) trace velocity decreases to 432 m/s. A summary of the results from the different processing steps is given in Table 1.

Atmospheric Propagation Modeling Results
For the atmospheric propagation modeling, a few major considerations from our data analysis have to be taken into account, namely, to match the modeling results with the observations of the trace velocity, the celerity and the observed azimuth deviation from the theoretical backazimuth. The observed trace velocity is related to the sound speed at the ray's turning height, following from the ray parameter being constant along the path and its relationship to horizontal slowness (e.g. Aki and Richards 1980). The celerity is derived from the travel time for the ray between the source and the receiver. For these two observables 1D or 2D ray tracing would be adequate (such as available by e.g. Margrave 2000). However, in order to model the significant azimuth bias from the atmospheric conditions prevailing on 12 December 2017, we used the 3D ray tracing code GeoAc described by Blom and Waxler (2012) and Blom (2013).
The physical parameters characterizing the atmospheric structure are retrieved from the European Centre for Medium-Range Weather Forecasts (ECMWF) six-hourly model analysis data (https:// www.ecmwf.int) for altitudes from the ground to about 60 km, above which they are combined with climatological models for horizontal wind and temperatures [HWM07, Drob et al. (2008); MSISE00, Picone et al. (2002)]. For the range-dependent 2D-ray tracing modeling the profile along the propagation path from the source to the receiver is interpolated with a resolution of 200 m both in range and altitude. For the 3D-modeling with GeoAc we defined a spatial grid of profiles with 0.5°9 0.5°both in latitude and longitude between 47°and 49°N and 16°and 21°E , respectively. A coarser grid was used, compared to Schneider et al. (2018), because it was computationally more efficient and did not affect noticeably stratospheric propagation results. The specific ECMWF model taken was for 06:00 UTC on 12 December 2017 and showed high wind speeds in particular for the zonal direction for both the troposphere as well as the stratosphere. The zonal wind speed significantly exceeded 100 m/s, i.e. is about a third to a half of the standard sound speed of 330-340 m/s near the surface. This wind speed pattern persisted on that day for latitudes between 40°a nd 50°N from the Mid-Atlantic through Eastern Europe and reaching almost Central Asia.
For the region in Austria and Hungary for which atmospheric infrasound propagation is modeled, maximum wind speeds are found in excess of 150 m/ s and in nearly easterly (zonal) direction. The zonal effective sound speed profiles incorporated in the 3D- therefore not allowing any returns from the lower thermosphere.
In the case of arrays, trace velocity is obtained directly from array processing. For other closely spaced stations this quantity can be estimated by the b Figure 3 a Waveforms registered at infrasound array PSZI and the colocated seismic station PSZ (vertical component); b results from F-K analysis of the acoustic data marked by a box in a. The graphical output shows slowness (in s/degree). The shown value translates into an apparent velocity of 407 m/s Figure 4 PMCC results from processing the PSZI data: a analysis was carried out without considering array element altitudes; b analysis results accounting for the differences in altitudes. The arrowed lines mark the average (mean) estimates for azimuth and trace velocity, which are also given. Only pixels are shown that were used for averaging With the 2D ray-tracing along the great-circle path from Baumgarten to PSZI (Fig. 6) both tropospheric as well as stratospheric ducting is expected from the high-velocity wave-guides identified in the effective sound speed profiles (Fig. 5). However, the tropospheric wave undergoes a number of bounces before reaching PSZI and, with the expected attenuation in the lower atmosphere, should therefore not be observable. For the stratospheric wave-guide the ray reaches station PSZI with the second bounce. Due to the favourable atmospheric conditions for stratospheric ducting, the observed signal must be associated with this wave type, as also implied by the corresponding trace velocity found at PSZI (Figs. 3,  4). Since the high stratospheric winds are directed almost to the east, it is expected that they may directly affect the apparent arrival direction, as measured by the backazimuth of the signals. Therefore, we performed additional 3D raytracing with GeoAc to (1) identify the domains of tropospheric and stratospheric ducting and (2) to determine the bearing angle for the eigenray to PSZI. For the latter purpose, standard GeoAc parameters are chosen with take-off angles in the range of 35°-50°and azimuth error up to 9°. The result of this ray-tracing is shown in Fig. 7 for the first two bounce points for the guided waves through the troposphere and stratosphere. As also found and argued by Schneider et al. (2018), atmospheric infrasound propagates mainly to the east stretching from the north to the south, with no arrivals being expected in directions towards the west. For tropospheric waves, the propagation is mainly to the north and east. In south-easterly directions, there is a For the arrival at PSZI, a celerity of 283 m/s is obtained from 810 s of travel time and the distance of 230 km. The 3D ray tracing also resulted in an azimuth deviation of more than 6°between the source bearing and the backazimuth of 288°for the eigenray. Removing this bearing difference from the backazimuth estimate the residual is further reduced from up to 10°down to less than 4°.

Seismological Aspects
Besides the large number of seismo-acoustic arrivals Schneider et al. (2018) studied a few seismic signals that they associate to the Baumgarten explosion as Rg waves with group velocities between 500 and 800 m/s. While two of these recordings are from AlpArray stations, for which the data are not openly available at this time and could not be studied here, we analyze only the data from the open station MODS in Slovakia located at a distance of about 30 km from the explosion site. The seismic signal from this station was shown by Schneider et al. to be the most convincing, as it does not only exhibit a high signal to noise ratio, in the 0.6-1 Hz frequency band, but it also shows almost perfect radial polarisation. Therefore, we use the MODS station to determine a preliminary local magnitude for the Baumgarten explosion to explain the lack or scarcity of seismic arrival observations. For estimating the local magnitude, we carried out Wood-Anderson simulations of the records (Fig. 8), measured the maximum horizontal amplitude according to Richter (1935) and applied standard amplitude-distance attenuation. As a result, we obtained a local magnitude ML of about 0.6-0.9. Of course, this magnitude is affected by shortcomings in the attenuation curve and the source depth. While the use of the attenuation curve for Southern California  Richter (1935). For this simulation, an additional high pass filter was applied to enhance the S wave group 4966 K. Koch et al. Pure Appl. Geophys. may have a small effect at the short range, the implicit source depth of 15 km for earthquakes in California may introduce a stronger bias for an explosion at the surface. However, this bias exists whenever near-surface sources are considered. To verify our magnitude measurement within this framework, we also analysed MODS data from two earthquakes in April 2016 and November 2017, having been reported with magnitudes of 4.0 and 3.1 respectively, and both located about 60 km from the station. For these events, at fixed depths of 8 and 11 km, we determined magnitudes values of 4.0 and 3.0. Thus, the magnitude estimated for the Baumgarten explosion may be reasonable within the 0.3 magnitude unit uncertainty range given above. The magnitude ML of less than 1 can explain well why only a few stations within 50 km range were able to register a seismic signal. According to forensic information cited in Schneider et al. (2018), the source may have been a strong pressure pulse from the sudden release of highly pressurized gas rather than an explosive combustion. Secondly, the sedimentary structure of the Vienna basin (OMV, pers. communication, 2018) is in good agreement with the observed Rg velocities of less than 1 km/s. Slow and soft sediments near the surface may have contributed to inefficient coupling into seismic waves, as reflected by the small seismic magnitude of the Baumgarten event.

Discussion
The Baumgarten explosion has produced a widely observed acoustic signal registered by the seismic stations of the AlpArray network (Schneider et al. 2018). With these seismo-acoustic observations only one parameter can be determined, namely the celerity following from the arrival time of the signal, that represents the mean velocity referred to a path along the ground surface. According to Negraru et al. (2010), propagation within the tropospheric, stratospheric or thermospheric duct can be associated with specific celerity ranges. This was indeed found for the AlpArray observations (Schneider et al. 2018), where observed celerities between 320 and 345 m/s are attributed to tropospheric arrivals and celerities between 280 and 305 m/s to stratospheric returns.
For array measurements additional parameters can be retrieved, i.e. the azimuth of the incoming wave and the trace velocity which is the direct measurement of the effective sound speed at the ray's turning altitude. For the signal detected at PSZI on 12 December 2017, we have found a celerity of 285 m/s which indicated stratospheric ducting. This is confirmed by the propagation modeling and is also supported by the rather large trace velocity of the signal exceeding 400 m/s using two different analysis techniques, thus fitting into the atmospheric specifications of ECMWF for the altitude range between 35 and 40 km.
With respect to the initial azimuth residual of nearly 20°found at PSZI for the acoustic source location at Baumgarten we were able to reduce this error by taking into account the altitude differences between the four array elements. Another large part may be attributed to the atmospheric propagation parameters with the extreme wind conditions in the stratosphere introducing a significant bearing difference of more than 6°to the true backazimuth direction. Similar maximal azimuth differences have also been found by Smets et al. (2015) and Pilger et al. (2018) when comparing azimuth deviations of infrasound signals for ground truth sources in Central and Northern Europe. Besides uncertainties in the atmospheric specifications, leading to errors in the ray tracing calculations, the data processing uncertainty itself, as indicated by the azimuth difference from FK and PMCC analysis (see Figs. 3a, 4a, Table 1) is another factor which may contribute to the remaining residual.
Finally, we discuss the two scenarios mentioned before that may cast doubt on the PSZI signal being indeed related to the Baumgarten explosion. The first case with potential blasting activity at the aforementioned quarry situated along the initially determined backazimuth can be ruled out for a number of reasons: (1) former such blasts created significantly different signal shapes at PSZI than the 12 December signal; (2) the determined trace velocity ([ 400 m/s) for a signal from a distance of about 40 km should not be ducted in the stratosphere, where only such velocities exist from the effective Vol. 177, (2020) The 12 December 2017 Baumgarten Gas Hub Explosion 4967 sound speed profile; and (3) the answer to a corresponding inquiry with the quarry operator was negative. The second case with a source at altitude, encompassing aircraft activity frequently noticed in the PSZI vicinity, is highly unlikely. Supersonic booms have not been observed from nearby activity, but only from remote sources over the North Sea or the Mediterranean, with very different signal shapes and/or arrival directions. This applies to subsonic aircraft signals as well, where time-varying estimates for azimuth and trace velocity are commonly found.

Summary and Conclusions
In this study, we have carried out FK as well as PMCC analyses on an impulsive signal observed at PSZI some 230 km distance east-southeast from the Baumgarten gas hub near Vienna, Austria, where an explosion occurred at about 7:44:20 UTC on 12 December 2017. The resulting celerity of 283 m/s indicates stratospheric ducting of the signal, which is supported by the trace velocity in excess of 400 m/s, corresponding to the velocities of the effective sound speed profile for altitudes of 35 km and higher. The estimated backazimuth of the signal, associated with an excessive azimuth residual, was closely examined and shown to result from the differences in array elements' altitudes as well as from high stratospheric wind speeds impinging on the bearing to the source. The use of a 3D ray tracing code allowed to quantify these effects. In the context of the analysis guidelines applied at the IDC (CTBTO 2012) for infrasound phase, this study is significant, as it demonstrates that infrasound data and their interpretation need special treatment in the presence of adverse factors (local topography, unusual atmospheric winds). Such treatment is currently not applied (Mialle et al. 2019).
We also looked at seismological parameters that may be deduced for this event. Since a wealth of seismo-acoustic observations to distances of some 150 km stand in stark contrast to the seismic observations found for distances of less than 50 km (Schneider et al. 2018), we estimated a local magnitude for the explosion source and found an upper bound for ML of less than 1. This small magnitude explains well the scarcity of seismic observations, with Rg type signals demonstrated in this study, as well as by Schneider et al. (2018).
The Baumgarten explosion is remarkable in the context of CTBT monitoring: it was observed exceptionally well as an acoustic source both by seismic stations as well as by a favorably located infrasound array. It also demonstrated that such a source, even though it caused quite substantial structural damage at the explosion site, did not couple significant energy into seismic waves. So not all explosive sources meet the expectation to be detected by multiple technologies, even when a significant number of regional/national stations are deployed near the source. When considering the IMS network, a sparse teleseismic network designed to detect a 1kT explosion anywhere in the world and often set equal to a magnitude 4 event, the Baumgarten explosion puts a lower bound on the sizes and types of explosions that may be detected at IMS stations.