A new accurate retrieval of bromine monoxide inside minor volcanic plumes from Sentinel-5 Precursor/TROPOMI

<jats:p>Abstract. Bromine monoxide (BrO) is a key radical in the atmosphere, influencing the chemical state of the atmosphere, most notably the abundance of ozone. The main effect of BrO onto tropospheric ozone concentrations occurs in bromine release events in polar regions, salt pans and volcanic plumes. Ozone depletion caused by halogen release has been observed and modeled for such conditions, in particular inside volcanic plumes. Furthermore, the molar bromine to sulphur ratio in volcanic plumes is a proxy for the magmatic composition of a volcano and potentially an eruption forecast parameter. The integrated column of BrO in the atmosphere, which in turn serves as an estimate for the bromine content, can be detected simultaneously with SO2 via spectroscopic measurements using the Differential Optical Absorption Spectroscopy (DOAS). Thus, a direct derivation of the BrO/SO2 ratio can be performed from a single measurement. Satellite spectroscopic observations offer the potential to observe and monitor volcanic bromine release globally. The detection of BrO in volcanic plumes is limited by the precision and sensitivity of the retrieval, which so far only allowed for the detection of BrO during major eruptions, leading to a potential sampling bias when looking at the BrO/SO2 ratio. The The TROPospheric Monitoring Instrument (TROPOMI) onboard Sentinel-5 Precursor (S-5P) however, with its unprecedented spatial resolution of up to 3.5 × 5.5 km2 and high signal-to-noise ratio, enables the detection of BrO in minor eruptions or even quiescent degassing. In this study, we investigate, how far the BrO retrieval can be improved using TROPOMI data and how well BrO can be detected, even in small eruptions and during quiescent volcanic degassing. There are two steps, for which improvements in accuracy are investigated and applied: the improvement and quantitative determination of (1) the detection limit of the DOAS BrO column retrieval and (2) the correction of non-volcanic background BrO signal. First, the DOAS retrieval settings are varied and their influence on accuracy and precision is investigated with respect to the detection limit and potential systematic influences. Based on these results, we propose a dedicated DOAS evaluation scheme optimized for the detection of BrO in volcanic plumes. For the DOAS retrieval, we propose the use of a large fit window from 323–360 nm, yielding a factor of 1.8 lower statistical uncertainty compared to previous BrO DOAS algorithms, while not enhancing systematic influences. Second, the effect of the background BrO is reduced by a latitude dependent empirical correction scheme correlated to cloud information as well as information on the ozone column. Via these improvements, the combined statistical and systematic uncertainties of the resulting BrO vertical column density is in the order of 7 × 1012 molecules cm−2, which allows for the detection of even slightly enhanced BrO amounts inside minor eruptive plumes of bromine-rich volcanoes.                         </jats:p>

section_content: Bromine monoxide (BrO) has been found to be a catalyst for ozone (O 3 ) destruction (e. g. in the polar regions, see Barrie et al., 1988; Simpson et al., 2007). Volcanic bromine release is thought to also affect atmospheric ozone chemistry (von Glasow, 2010; Surl et al., 2021), and strong eruptions ejecting bromine into the tropical stratosphere can potentially impact global stratospheric ozone abundance over a ten-year period (Brenna et al., 2019). BrO in volcanic plumes was first detected at Soufriere Hills (Bobrowski et al., 2003) using ground-based Differential Optical Absorption Spectroscopic (DOAS, Perner and Platt, 1979) measurements. BrO and SO 2, which can be detected by DOAS simultaneously, have been measured over the last 20 years at many volcanoes using ground-based instruments and the BrO/SO 2 ratio has proven to be a useful tool both to determine the composition of a volcanic plume (Bobrowski and Platt, 2007; Kern et al., 2009; Hörmann et al., 2013) as well as a proxy for changes in the volcanic system (Bobrowski and Giuffrida, 2012; Lübcke et al., 2014; Bobrowski et al., 2015; Dinger et al., 2018; Warnach et al., 2019). 

After the unsuccessful attempt to detect BrO in volcanic plumes using satellite data from GOME and SCIAMACHY (Afe et al., 2004), volcanic BrO was first detected from GOME-2 in the plume of Kasatochi in 2008 by Theys et al. (2009). A global detection algorithm of BrO and the BrO/SO 2 ratios inside volcanic plumes from GOME-2 was developed by Hörmann et al. (2013). 

DOAS retrievals of BrO from satellites have varied over time, both w. r. t. the wavelength fit range as well as other DOAS fit settings. The first detection of tropospheric BrO from the GOME instrument by Wagner and Platt (1998), employs a wavelength range of 345 -359. 5 nm. Other studies use almost identical fit ranges for spectra from GOME (Richter et al., 1998; Hegels et al., 1998; Chance, 1998; Richter et al., 2002) and for spectra from OMI (Chance, 2002). For SCIAMACHY and GOME-2 measurements, the short wavelength limit of the fit range was lowered to 336 nm, and the upper limit varied between 347 nm (Afe et al., 2004; Begoin et al., 2010), 351 -352 nm (Valks et al., 2009; Theys et al., 2009), and 360 nm (Heue et al., 2011; Sihler et al., 2012; Hörmann et al., 2013). The short wavelength limit of 336 nm was also found best in a long-term study combining GOME, SCIAMACHY, and GOME-2 measurments (Bougoudis et al., 2020). For recent global studies of BrO from GOME-2 (Theys et al., 2011) and TROPOMI (Seo et al., 2019) the lower boundary was further lowered and a fit range of 332 -359 nm and 334. 6 -358 nm was chosen respectively. 

The DOAS BrO fit used in this study is based on the BrO fit developed by Sihler et al. (2012) for the retrieval of BrO in arctic regions, which was also used for the detection of BrO in volcanic plumes by Hörmann et al. (2013). There, a fit range from 336 -360 nm was chosen. The omission of shorter wavelengths was mainly attributed to the strong ozone absorption in the UV at high latitudes, caused by the large ozone columns and high solar zenith angles (SZAs). 

In this study, we report on the further improvement of the BrO DOAS retrieval with a specific focus on the detection of BrO in volcanic plumes and propose a new accurate BrO retrieval for the detection within minor volcanic plumes. Retrieval improvements are investigated by varying the DOAS fit settings, by both considering lower wavelengths and the inclusion/exclusion of the formaldehyde (HCHO) cross section in the DOAS fit. The aim of this task was to achieve both a low statistical uncertainty and low systematic errors, especially low interferences caused by clouds, ozone and HCHO. 

There are several reasons why -in contrast to previous retrievals -lower wavelengths can be included in current BrO DOAS retrievals for the detection of volcanic plumes. First, the introduction of the so-called Pukite terms (Puk, īte et al., 2010) can largely compensate the spectral effects of the strong ozone absorption. Second, volcanic plumes mainly occur at lower latitudes, where ozone columns and the SZAs are considerably lower compared to high latitudes. Third, the use of an earthshine reference spectrum recorded at the equator reduces the optical density of ozone to zero at the equator, further reducing the potential ozone interference at low latitudes. 

The question whether it is advantageous to include HCHO in the BrO DOAS fit is long debated, as the interference between BrO and HCHO is a well known difficulty both in BrO and HCHO DOAS retrievals (Theys et al., 2011; Vogel et al., 2013; De Smedt et al., 2018). The structure of the absorption cross-section of BrO and HCHO are very similar (cf. Fig. 1 ). 

Especially the fact, that the largest absorption peak of BrO overlaps with a major absorption peak of HCHO at 339 nm, can cause a spectral interference between both trace gases. For arctic applications, such as Wagner and Platt (1998), Richter et al. (1998), but also for the more recent study by Sihler et al. (2012), HCHO was not included as an absorber in the DOAS retrieval, because HCHO is almost not abundant in arctic regions. However, in equatorial regions, where HCHO can reach high column densities due to biogenic emissions or biomass burning, spectral interference with HCHO can become important. HCHO is therefore included within the DOAS analysis of the global BrO detection algorithms of GOME-2 by Theys et al. (2011) and the TROPOMI algorithm of polar BrO by Seo et al. (2019). Interferences between BrO and HCHO are observed for equatorial regions by Theys et al. (2011). Using 332 nm compared to 336 nm (used in Sihler et al., 2012) as the lower fit boundary reduces these interferences (Theys et al., 2011). Hence, both Theys et al. (2011) and Seo et al. (2019) employ fit ranges with a lower wavelength limit of 332 nm and 334. 6 nm respectively. 

The stratospheric background of BrO together with the potential BrO background of the free troposphere is up to one order of magnitude larger than volcanic BrO columns of small eruptions. Therefore, an imperfect subtraction of the background column can yield uncertainties in the order of small volcanic signals. Hence, an accurate background correction is essential. 

Fortunately, horizontal gradients of the stratospheric BrO VCD occur on a much larger scale than the extension of volcanic plumes. Moreover, while stratospheric gradients mainly originate from changes in the tropopause height (Sihler et al., 2012), volcanic plumes are dispersions of point sources. Thus, the patterns are spatially independent from each other and a separation can be achieved by the application of a spatial high-pass filter. Moreover, volcanic plumes can be independently detected due to their enhanced SO 2 column density. Background correction on this basis was for example implemented by Hörmann et al. (2013) by first masking the volcanic plume, which is well detectable via the SO 2 signal, and subtracting a 2D polynomial fitted to the area around the volcanic plume. While our BrO DOAS retrieval is in principle applicable to the complete globe, similar to the retrievals presented in Theys et al. (2011); Sihler et al. (2012) ; Seo et al. (2019), it is specifically optimized for the most accurate detection of BrO originating from volcanic plumes, which are local/regional events whose spatial extent is typically of the order of 10 -100 km, only rarely exceeding 1000 km. 

This paper is structured as follows: in Sect. 2, we describe the TROPOMI instrument and the basics of our BrO DOAS retrieval from satellite. The investigation of different DOAS fit settings is described in Sect. 3. For these fit settings the systematic influence of clouds, ozone and HCHO are discussed, corrected, and the remaining systematic influences are quantified in Sect. 4. In combination with an investigation of the statistical uncertainties (Sect. 5. 1), the most suitable fit settings for an accurate retrieval within minor volcanic plumes is proposed in Sect. 5. The proposed retrieval algorithm is tested on multiple volcanic plumes of different emission strengths and at different latitudes in Sect. 6 to demonstrate its accuracy and efficacy. 

Lastly, conclusions are drawn in Sect. 7.
section_name: Introduction

section_name: Methods

section_content: The TROPospheric Monitoring Instrument (TROPOMI, Veefkind et al., 2012) is onboard the satellite Sentinel-5 precursor (S-5P) of the European Space Agency (ESA). S-5P is a sun-synchronous, polar-orbiting satellite with an orbital inclination of 98. 7 • on a low-Earth orbit (altitude: 824 km) launched on 13 October 2017. Its equator crossing time is roughly at 13 : 30 h mean local solar time and it has a repeat cycle of 17 days. S-5P and its single payload TROPOMI are designed to determine and monitor the atmospheric composition. TROPOMI uses a push-broom method to scan the Earth's surface. The opening angle of 108 • leads to a ground swath of approx. 2600 km (separated to 450 pixels) and allows for an almost complete daily coverage of the whole globe. 

TROPOMI is a hyperspectral imaging spectrometer measuring radiance spectra via four different detectors in the ultraviolett 270 -320 nm (UV), UV-visible 310 -500 nm (UVIS), near-infrared 675 -775 nm (NIR), and short-wave-infrared 2305 -2385 nm (SWIR). A detailed overview of the spectral characteristica for each detector is given in Veefkind et al. (2012).
section_name: The TROPOMI instrument

section_content: For the DOAS retrieval of BrO the L1B radiance spectra from the UVIS channel are used. The UVIS channel has a along-track pixel size of 7 km, which was reduced to 5. 5 km on 6 August 2019. There are 450 across-track pixels, whose size varies from 3. 5 km at nadir up to 14. 4 km at large viewing angles. Thus, the highest spatial resolution is 5. 5 × 3. 5 km 2 = 19. 25 km 2.
section_name: TROPOMI UVIS spectra used

section_content: For the investigation of systematic errors and their correction several additional input data-sets are used in this study. 

First, for the investigation of the influence of clouds onto the retrieval and for its correction, there are two cloud datasets used: the cloud height (CH) and the cloud fraction (CF). Both data-sets are calculated using the FRESCO algorithm (Compernolle et al., 2021) and provided within the NO 2 operational product (van Geffen et al., 2021; NO2, 2021). While the cloud fraction is directly inferred within the FRESCO algorithm, the cloud height is calculated from the FRESCO cloud pressure (CP), surface pressure (SP) and surface altitude (SA) all provided within the NO 2 product (van Geffen et al., 2021; NO2, 2021) via the barometric pressure formula. To ensure the use of the latest FRESCO product, we used the reprocessed version of the NO 2 product, provided via the sentinel-5 precursor expert hub (https://s5pexp. copernicus. eu/). For this study, we choose the FRESCO cloud products, as they are easily accessible via the operational NO 2 product and are therefore available for the complete TROPOMI sensing period. However, also other cloud products are available (Latsch et al., 2022) and could in principle be used. For example the use of the cloud fractions provided by the Mainz Iterative Cloud Retrieval Utilities (MICRU, Sihler et al., 2020) was tested in Warnach (2022) and found suitable for application within the retrieval presented here. However, unfortunately it is not available for the complete TROPOMI measurement period. 

Second, for the study of the influence of HCHO onto the BrO retrieval, the operational TROPOMI L2 HCHO product is used (De Smedt et al., 2018). 

Third, for the estimation of the stratospheric column strength the O 3 VCD is used. The O 3 VCD is derived directly from the BrO DOAS fit (cf. Sect. 4. 2 for a detailed description).
section_name: Additional input data

section_content: The DOAS method was introduced by Perner and Platt (1979) and is based on the attenuation of light in the atmosphere described by the Beer-Lambert-Bouguer law. The DOAS method makes use of the fact that molecules are distinguishable by their unique absorption structure originating from electronic, vibrational, and rotational excitations. These structures usually contain high frequency structures, which can be distinguished from the broad band extinction features of atmospheric scattering processes such as Mie and Rayleigh scattering. Essentially, the optical density τ of the absorption of the trace gas in question is retrieved by comparing the measurement spectrum's intensity I(λ) with a so called «reference spectrum» I 0 (λ):

The result of the DOAS analysis will be the difference of the columns of the trace gases between reference and measurement spectrum, the slant column density (SCD). DOAS thus requires only knowledge about the absorption cross-sections σ of the relevant absorbing molecules. Traditional DOAS algorithms employ a combination of a linear (for the contribution of the ref. 

spectra) and non-linear (for spectral shift and distortion compensation) fit in order to determine the column of the absorber. 

However, for satellite applications, where huge data-sets are recorded daily, a linearisation of the non-linear DOAS fit was developed by Beirle et al. (2013) and is used here. The fitting algorithm used here was implemented and described in detail by Borger et al. (2020). The slant column density (SCD) retrieved in the DOAS fit is the integrated concentration along the light path. As the light path depends on the measurement geometry, i. e. viewing angle, solar angle, and the radiative transfer, the One simple approximation for the AMF is the geometrical AMF, depending only on the solar zenith angle (SZA) and viewing zenith angle (VZA):

The geometric AMF is a valid approximation for stratospheric light paths in the UV if the SZA does not exceed 70 • (Burrows et al., 2011) and accounts for the satellite and solar viewing geometry. For tropospheric columns, usually a more complex derivation using radiative transfer calculations is performed, including the trace gas profile (i. e. the plume height in the case of a volcanic plume), surface albedo and cloud cover. For trace gases close to the surface the AMF can be up to one order of magnitude lower compared to the geometric approximation (Theys et al., 2017). 

In this study we employ only geometric AMFs, since (1) the focus of the study is to look at systematic influences and background corrections and (2) as the main focus of volcanic BrO is to derive the BrO/SO 2 ratio, for which the AMF divergence will cancel out in first approximation.
section_name: The DOAS fitting routine

section_content: Traditionally, DOAS retrievals from satellites use the direct solar spectrum recorded on each day as a reference spectrum. This ensures that the DOAS retrieval can yield absolute slant columns of the respective trace gas. However, the direct solar spectra are obtained through a different optical input channel and traverses different optical parts (see Veefkind et al., 2012). Thus, these differences will cause stripes in the across-track direction (Richter and Wagner, 2001; Chance, 2007). While this can be corrected via so called »destriping« algorithm (e. g. Chance, 2007; Boersma et al., 2011; Hörmann et al., 2016), in recent times it has been found advantageous to use an earth-shine spectrum as reference spectrum for the retrieval (Theys et al., 2017; Seo et al., 2019). Ideally, such an earthshine spectrum is calculated from radiances obtained at regions where the respective trace gas is not abundant. Thus, regions over the equatorial Pacific are used (Sihler et al., 2012; Hörmann et al., 2013; Seo et al., 2019), far off potential anthropogenic pollution sources and remote from volcanoes. In this study, we expand the region an earth-shine spectrum is calculated using the complete equatorial latitude band as a reference region (20 • S to 20 • N and

For SO 2, the assumption that remaining SO 2 is negligible holds true for this reference area. For BrO however, there is still a significant stratospheric background signal over the equatorial region with VCDs of roughly 3. 5 × 10 13 molecules cm -2, as stated by Richter et al. (2002) and also found within this study. Thus, using the earthshine reference will still yield differential SCDs instead of absolute SCDs. However, the stratospheric signal will introduce an offset in the resulting slant column densities, which is only slowly changing in spatial direction. In this study, the focus is on local volcanic plumes. are singular events, whose (1) spatial patterns usually have a high contrast compared to the background, and (2) spatial extent is almost always much smaller compared to stratospheric patterns. Thus, the stratospheric signal can be removed effectively by performing a local background correction for SO 2 and BrO (see Sect. 6. 2). 

The earthshine reference spectra are calculated in the following way: Firstly, all radiances within the global latitude band [180 • W to 180 • E] around the equator [20 • S to 20 • N] are selected, as mentioned above. In a second step, the radiances are normalized with respect to their maximum intensity in order to prevent that scenes or geometries which result in a brighter backscatter radiation, e. g. cloudy scenes, to be systematically weighted stronger. In the last step, the mean of the normalized radiances is taken as the daily reference spectrum in the DOAS analysis. This process is performed separately for each spectra of the 450 across-track detector pixels of the TROPOMI instrument and for each day. A comparison of results using either an earthshine or a direct solar spectrum for measurements over the pacific region [±20 • N, 135 -105 • W] on 1 October 2018 is shown in Fig. 2. It can be seen that the retrieved VCDs show no difference (here the stratospheric influence in the irradiance data is eliminated for comparison by subtracting the mean BrO VCD). The fit root-mean-square (RMS), however, is about 25 % lower for the earthshine fit and roughly 6 × 10 -3. This RMS distribution is in very good agreement with RMS reported over a similar pacific region by Seo et al. (2019) (cf. Fig. 11b ) who employed a DOAS earthshine fit independently from the fit presented in this study.
section_name: Selection of reference region

section_content: In the following sections, different wavelength ranges of the BrO DOAS fit are investigated. The performance with respect to the statistical variation, systematic spatial patterns of various origin and interference with HCHO is checked by comparing the results of four fit wavelength ranges, where the upper fit range boundary is unchanged at 360 nm, as varying the upper fit range boundary was found to be less critical for the fit performance (Vogel et al., 2013; Seo et al., 2019). This upper wavelength fit range boundary is thus chosen at 360 nm as chosen in Sihler et al. (2012), which is used as a baseline BrO DOAS fit. The lower wavelength boundary is varied, yielding the following wavelength ranges:

-336 -360 nm: This smallest fit wavelength range includes four BrO absorption bands (see Fig. 1 ) and is the baseline wavelength range used by Sihler et al. (2012); Hörmann et al. (2013 Hörmann et al. (, 2016)). 

-332 -360 nm: This fit wavelength range is chosen, because it adds one more BrO absorption band compared to the baseline, while adding no additional HCHO absorption band (see Fig. 1 ). This lower wavelength boundary is also used -323 -360 nm: This largest fit wavelength range was chosen, such that it includes seven BrO absorption bands while avoiding interference with very strong O 3 or SO 2 bands at wavelength below 320 nm. It was not used in DOAS satellite retrievals so far due to the strong O 3 absorption for high latitudes. 

- 323 -328. 5 & 332 -360 nm: This encompasses also the largest fit wavelength range, but the HCHO absorption band around 329 nm, which overlaps with a BrO absorption band, is excluded in order to avoid HCHO interference (see Fig. 1 ). 

All four fit wavelength ranges are implemented with two different fit settings, which differ only in the inclusion and exclusion of the HCHO absorption cross section in the DOAS fit. Thus, there is a total set of eight different fits compared. The fits including the HCHO absorption cross-section are labeled »HCHO«, the ones excluding HCHO are labeled »SR«(standard retrieval), i. e. the fit HCHO 332 -360 nm corresponds to a DOAS fit evaluated in the wavelength range 332 nm to 360 nm with the HCHO absorption cross-section included, while the fit labeled SR 332 -360 nm corresponds to the same DOAS fit just without the HCHO absorption cross-sections included. These labels will be used throughout this paper. The overview over the DOAS fit settings, and the included absorption cross-sections for both SR and HCHO fits are listed in Table 1. Since the sensitivity study described in Sect. 4 suggests that the fit SR 323 -360 nm is the most suitable fit (cf. Sect. 5), most of the results in the main part of the paper will be shown only for this fit. The corresponding results for the other fits can be found in the appendix Sect. A1.
section_name: Fit test setup

section_content: The global BrO VCD map for 1 October 2018 for the fit SR 323 -360 nm is shown in Fig. 3. Enhanced BrO occurs at high latitudes and is caused by the stratospheric signal, related to variations in the tropopause height (Sihler et al., 2012; Schöne, 2023). Hence, the variations of the BrO dSCD is directly linked to variations in the ozone column, which in turn can be used as a proxy for the strength of the stratospheric BrO column (Sihler et al., 2012). The enhanced BrO columns in southern high latitudes can be additionally caused by tropospheric BrO events occurring during polar spring (Sihler et al., 2012; Herrmann et al., 2022). In the equatorial region there are also distinct spatial patterns in the BrO columns not originating from volcanic plumes. They might be related to clouds and will be investigated in more detail in Sect. 4. 1.
section_name: Overview of the fit performance

section_content: In this section, the systematic influences of following three effects onto the BrO retrieval are investigated for the eight different BrO fit settings. Potential systematic influences of clouds and ozone in the polar as well as equatorial region onto the volcanic BrO data are described (Sect. 4. 1) and corrected (Sect. 4. 2). The question of spectral interference of HCHO is investigated and quantified (Sect. 4. 3). Lastly, the overall effect of systematic influences is quantified (Sect. 4. 4). Based on the latter as well as the statistical uncertainties (quantified in Sect. 5. 1), an optimal fit setting is chosen in Sect. 5 and the combined fit uncertainty is quantified for different volcanic scenarios globally.
section_name: Investigation of systematic effects

section_content: The aforementioned structures in the global BrO map (cf. Fig. 3 ) are probably correlated to clouds (equatorial region) and the ozone column (high-latitude regions). To showcase this, we look at two regions: First, a pacific equatorial region ([±20 • N,

, where the cloud structures are assumed to dominate the systematic structures, and second a high-latitude region ([40

, where changes in the stratospheric column height are assumed to dominate the systematic structures in the BrO column. The former is depicted in Fig. 4 including both FRESCO CF (Fig. 4a ) and FRESCO CH (Fig. 4c ) as well as the BrO VCD for fit SR 323 -360 nm (Fig. 4b ). Furthermore, also the BrO VCD after applying the correction scheme for clouds and ozone (described in Sect. 4. 2) is included in Fig. 4d. BrO partial column below the clouds or a purely spectral response or has other origins. We therefore test a correction of these structures via an empirical scheme based on the two cloud parameters. 

It is noteworthy that this decrease in the BrO column for cloudy scenes is an indication for cloud shielding of tropospheric

BrO columns and can be used to derive and estimate for the tropospheric column, as done by Theys et al. (2011). However, as it is unclear if the cloud relation does not have a different origin and in this study we simply want to remove non-volcanic influences on the BrO column, this question is not indulged further here. The high-latitude example is depicted in Fig. 5, including the FRESCO CF (Fig. 5a ) as well as the O 3 VCD (Fig. 5c ) and the BrO VCD (Fig. 5b ). There is a clear gradient in BrO VCD visible along 50 • N, where BrO VCDs increase by more than 1 × 10 13 molecules cm -2. This follows the gradient line of O 3 along the same latitude, indicating that this increase is most probably due to a higher stratospheric column. This gradient is overlayed by structures in the BrO VCD map south of 50 • N, e. g. between 40 -45 • N and 85 -95 • W. These coincide with cloud cover (indicated by the cloud fraction). In contrast to the equatorial region, the sign of the relation is inverted and high cloud fraction results in elevated BrO VCDs. Both these gradients and structures in the BrO VCD are reduced after applying the correction scheme (see Sect. 4. 2), as can be seen in the corresponding map of the corrected BrO VCD (Fig. 5d ). 

In the following, we will investigate the cloud and ozone dependency for different latitudes. In order to investigate this The mean BrO VCD is systematically decreasing, both for increasing CH and increasing CF. The corresponding plots for all eight fits can be found in the appendix Fig. A1. The relation of the mean BrO VCD both to CH and CF appears to be almost linear for all eight fit ranges, but with varying strength. Therefore this linearity seems independent of the fit settings chosen 275 here allowing to derive a correction term, which will be described in Sect. 4. 2. 

While the cloud influence is very strong and rather linear for the equatorial region, there are also systematic structures in mid-and high latitudes (compare Fig. 3 ). Since these structures are not all related to the presence of clouds, but also to the O 3 columns (cf. equator case is also found for neighbouring latitude bands from 30 • S to 10 • S and from 10 • N to 30 • N. For mid-to high latitudes north and south, this dependency of BrO on CH and CF is changing and/or is overlayed by other features. 

While the general CH dependency of the BrO VCD (higher clouds, lower BrO VCD) is prevailing for high latitudes, the CF dependency is changing its sign, so that it is inverted for 70 • N -90 • N, i. e. a lower cloud fraction leads to a lower BrO VCD. This is already observed in the case depicted in Fig. 5, indicating that this is a general finding applicable to high-latitudes and not specific to the case illustrated in Fig. 5. In addition, there is a positive BrO signal overlayed for cloud heights between 1000 and 5000 m, independent of cloud fraction. This dependency is clearly visible for northern latitudes (30

and appears to shift to higher cloud heights for higher latitudes, i. e. this shifts from 0 -3000 m, to 1000 -4000 m and to In order to correct for these systematic effects on the BrO VCDs, different correction approaches are investigated in the next section (Sect. 4. 2). Since the dependency on cloud parameters is latitudinal dependent, the correction approach is applied to different latitude bands independently. Furthermore, since the systematic structures in high latitudes might be correlated to the stratospheric column, the inclusion of the O 3 VCD in the correction scheme as a first order indicator for the extent of the stratospheric column is tested. The three correction schemes are referred to as »ozone«, »latitude«, and »ozone latitude« in the following (see Table 2 ). All schemes include the two cloud parameters cloud fraction and cloud height derived from FRESCO. In addition, the schemes »ozone« and »ozone latitude« also include the O 3 VCD, taken from the BrO fit. Since there are two different O 3 absorption cross-sections (see below) as well as the two »Pukite«-pseudo absorbers included in the fit (see Table 1 ) to account for the spectral absorption signal of O 3, the contributions of the individual fitted O 3 terms are combined to the complete slant column following Puk, īte and Wagner (2016) :

where the slant columns S O3,223K, and S O3,243K correspond to the cross sections for two temperatures while the two last terms S O3, λ, S O3, σ correspond to the two Pukite-term pseudo-absorbers (Puk, īte et al., 2010), λ to the wavelength, and σ O3

to the ozone absorption cross-sections. For simplicity, the for the wavelength λ and the cross-section σ O3 was taken only at the wavelength of the first maximum of the absorption cross section closest to the lower wavelength boundary (e. g. for the fit wavelength range 323 --360 nm: λ = 325. 0 nm and σ O3 = 1. 5 × 10 -20 cm 2 molecules -1 ). Each of the three BrO correction schemes is employed as a multi-dimensional polynomial P including the respective fitted parameters (see Table 2 ). The polynomial P comprising of p fitted parameters x 1. . . x p and their corresponding polynomial degrees n 1. . . n p is defined by in the following way:
section_name: Influence of clouds and ozone

section_content: where c i1,. . .,ip are the resulting fitted coefficients describing the dependency of the fitted quantity (in our case the BrO VCD)

to the fitting parameters (in our case cloud parameters). 

This means for example, that for the correction scheme »latitude« (see Table 2 ), where the fitted parameters are the cloud height and the cloud fraction with polynomial degree of two for both parameter, the polynomial P looks like:

The principle output of this fit are the coefficients c i1,. . .,ip. Afterwards the polynomial P is subtracted from the uncorrected BrO VCDs forming the corrected BrO VCDs. 

In this section the polynomial P is tested for the complete global data-set. However, for the the schemes »latitude« and »ozone latitude« there are separate sets of coefficients c i1,. . .,ip for each 20 • latitude band. Thus, for the retrieval of BrO inside a volcanic plume (Sect. 6) the correction coefficients can be calculated specifically for the latitude band around the volcano's latitude, ensuring a most accurate correction specific for each volcano's latitude. 

In the following, the correction schemes will be applied and discussed for the fit SR 323 -360 nm only. The correction schemes yield similar results for the other fits, but the fit SR 323 -360 nm shows the best performance with respect to HCHO influence (see Sect. 4. 3) and statistical uncertainties (see Sect. 5. 1) and is therefore chosen as the best fit in the end. 

The corrected global maps -after applying the correction schemes »ozone«, »latitude« and »ozone latitude« -are schown in Fig. 7a -c respectively for the settings SR 323 -360 nm (excluding HCHO). Note that the scale of the BrO VCD color-code differs from the uncorrected map (cf. Fig. 3 ). 

The systematic structures in the equatorial region are effectively removed both by »latitude« and »ozone latitude« schemes (Fig. 7b, c ). For the »ozone« correction scheme (Fig. 7a ), however, the systematic structures prevail reduced in magnitude. 

The systematic structures in mid-to high-latitude are reduced for all three correction schemes. However, for the »latitude« scheme, which is the only one where the O 3 VCD is not included as a parameter, there are still significant structures visible for These are probably tropospheric signals in the polar spring (compare Wagner and Platt, 1998; Sihler et al., 2012). If volcanic plumes are analysed for such conditions, the correction schemes developed in this study will not be able to separate the BrO in the volcanic plume form such events of enhanced tropospheric BrO. 

For both, the »latitude« and the »ozone latitude« schemes, the correction is done separately for 20 • latitude bands. It can be seen that on the edges of the latitude bands jumps in the BrO VCD are occurring (see Fig. 7b, c ). However, for the »ozone latitude« scheme (Fig. 7c ) these only occur at high latitudes. As the latitude band can be adjusted for each volcanic plume, they can be chosen so that the volcano is at the center of the latitude band so that the volcanic plume is not located on the edges. 

Furthermore, only a 20 • latitude band around the volcano should be considered. Only very large plumes -recurring once every several years -might exceed this limit. For these cases the correction coefficient can also be applied For the development of a global BrO product including these corrections these jumps might be avoided by an interpolation of the correction parameter with latitude. 

In order to quantify the efficacy of the different correction schemes, the mean of the corrected BrO VCD is calculated for each CH-CF bin for each latitude band independently. The corresponding plots are depicted in Fig. 8 for all three correction schemes as well as the bins without any correction. Ideally, there should be no or little dependency on the cloud parameters remaining after correction and no structures independent of the cloud parameter appearing. 

All three correction schemes reduce the cloud dependency. It is noteworthy that the color-scale of the uncorrected plots is twice as large as the corrected ones. 

It can be seen that the correction schemes working on each latitude band separately (»latitude« and »ozone latitude«) are much better at reducing the dependencies, both on cloud fraction as well as on cloud height compared to the global »ozone« scheme, which produces new CH-CF gradients in equatorial regions. However, these gradients are weaker than the original gradients. The »latitude« correction reduces most cloud dependencies for latitude bins between 90 • N and 50 • S. However, cloud independent signals remain in mid-to high-latitudes. The inclusion of ozone in the correction scheme (»ozone latitude«) is further reducing the cloud dependencies, but most importantly drastically reduces the cloud independent structures. For this scheme maximal differences in the CH-CF plots are typically 3 × 10 12 molecules cm -2 for latitudes between 90 • N to 50 • S and 5 × 10 12 molecules cm -2 for latitudes south of 50 • S. 

Thus, including ozone in the latitude dependent correction scheme (»ozone latitude«) appears to be the best correction both with respect to consistency in the global map (see Fig. 7c ) and in the CH-CF bin plots (see Fig. 8 ). Therefore, this correction scheme is chosen for this study. In this approach the correction scheme is performed on separate 20 • latitudinal bands. For the calculation of the correction factors, we thus propose to use the data on the 20 • latitudinal band around the latitude of the volcanic plume. 

From now on, all BrO VCDs used in this study will be by default corrected using the »ozone latitude« correction scheme. The main questions is whether the inclusion of the HCHO absorption cross-section in the fit is necessary for the analysis of BrO in volcanic plumes, i. e. whether the exclusion of HCHO in the fit introduces a systematic offset in the BrO results in presence of enhanced HCHO absorption. One reason why it would be advantageous not to include the HCHO absorption cross-section is that the analysis yields a larger noise in the retrieved BrO SCDs if HCHO is included (as shown in Sect. 5. 1). 

Equally important is the question, which fit wavelength range results in the lowest systematic influence of HCHO onto the BrO results and to quantify the strength of a potential false BrO signal in the presence of enhanced HCHO absorption. 

As information on the atmospheric HCHO column, the HCHO SCDs from the operational S-5P/TROPOMI HCHO L2 product (De Smedt et al., 2018) are used. The cloud information (cloud fraction and cloud height) is taken from FRESCO as in the previous sections. 

The HCHO SCD derived from the operational S-5P/TROPOMI HCHO L2 product for the selected region is plotted in Fig. 9 (noted »HCHO L2 SCD«, topmost left). There is a region of high HCHO columns exceeding 3 × 10 16 molecules cm -2 on the African continent around 5 • S and 29 • E. Most importantly, this region is mostly cloud free (see Fig. 9, noted »FRESCO cloud fraction«, topmost right), ensuring minimal influence of clouds. 

The BrO VCDs are corrected using the »ozone latitude« correction scheme described in Sect. 4. 2. The corresponding maps for the eight fits are plotted in Fig. 9. The maps of the fits including HCHO (HCHO, Fig. 9 left column) do not reveal a distinct positive BrO VCD signal in the region of the HCHO event. To the contrary, there rather seems to be a systematic negative offset in this region, most visible for the fits HCHO 323 -360 nm and HCHO 323 -328. 5 & 332 -360 nm. The fits excluding HCHO (SR, Fig. 9 right column), show a positive offset in the region of high HCHO SCDs. This is most prominent for the fit range SR 332 -360 nm (fourth from top, right column), where the BrO VCDs exceed 1 × 10 13 over a large area. 

In order to investigate the potential systematic response of the BrO retrieval onto HCHO, the pixels are separated into bins of HCHO SCD and cloud fraction. In order to separate the HCHO effect from potential cloud influence (which is also an response, while the fits starting at 323 nm show a slightly weaker response. The fits not including the HCHO cross-section (SR, dashed lines in Fig. 10 ) show a positive response to HCHO, i. e. BrO SCDs increase for increasing HCHO SCDs. However, this only occurs for high HCHO SCDs and varies strongly for the different fit ranges. The effect is strongest for the fit range SR 332 -360 nm, where BrO VCDs exceed 1 × 10 13 molecules cm -2, whereas it is much less pronounced for the other fit ranges,

where maximal BrO VCDs are in the range of 5 × 10 12 molecules cm -2. In addition, for low HCHO SCDs the BrO VCD changes only slightly or remains constant. This is most prominent for the fit range SR 323 -360 nm (green dashed line). Here, the BrO VCDs are zero for HCHO SCDs below 1. 5 × 10 16 molecules cm -2. For higher HCHO SCDs, the BrO VCDs slowly increase, reaching about 5 × 10 12 molecules cm -2 for maximal HCHO SCDs of about 4 × 10 16 molecules cm -2. Interestingly the fit range SR 323 -328. 5 & 332 -360 nm, in which an overlapping BrO/HCHO absorption band is left out compared to SR 323 -360 nm and thus a weaker response would be expected, shows a higher false positive signal. 

In order to include also cloudy scenes, the mean BrO SCD is calculated for the complete HCHO-CF bins. The corresponding plots for the eight fits are shown in Fig. A3, where the four different fit wavelength ranges are plotted from left to right. The fits excluding the HCHO cross-section are in the upper row, while the fits including the HCHO cross-section are plotted in the lower row. 

Both effects, the decreasing of BrO SCDs with increasing HCHO SCD for the fits including the HCHO cross-section (HCHO) as well as the opposite behavior of the fits not including the HCHO cross-section (SR), prevail also for cloudy scenes, but are sometimes weaker. Nevertheless, for volcanoes in regions where HCHO columns are high, it should be taken care for spatial HCHO patterns overlapping with the volcanic plumes (identified by its SO 2 signal). 

It can be thus concluded, that the inclusion of the HCHO cross-section in the BrO DOAS retrieval is not necessary. Even though there can be a systematic false positive BrO response produced by very large HCHO columns, this effect is very low for the fit range SR 323 -360 nm, with no response for HCHO SCDs below 1. 5 × 10 16 molecules cm -2 and a response not exceeding 5 × 10 13 molecules cm -2 for very large columns. HCHO columns above 1. 5 × 10 16 molecules cm -2 are only observed for tropical biomass burning events, mostly ocurring around the equator, i. e. the systematic impact of HCHO on the BrO retrieval using the SR 323 -360 nm fit is negligible for most cases.
section_name: The resulting O

section_content: In this section, the remaining potential systematic influence of clouds and other systematic effects onto the retrieval are investigated, i. e. the potential »false« BrO signal which could be created by clouds, ozone or other sytematic patterns. For this purpose the BrO VCD (after the correction with the »ozone latitude« correction scheme) is investigated for the fit SR 323 -360 nm. 

The uncertainties of any data-set is comprised of both statistical and systematic uncertainties. By averaging over neighbouring pixels the statistical variation can be gradually reduced and often the remaining uncertainties can be attributed to systematic effects. The statistical variation can be quantified by a Gaussian fit over the BrO distribution of the binned pixels. This Gaussian fit is performed separately for each latitude band as well as for each binning factor and the resulting standard deviations are plotted in Fig. 11. It can be seen in Fig. 11 that the standard deviation is decreasing slower for increasing increasing binning factors (square root of the number of pixels binned) compared to the purely statistical model (black line in Fig. 11a ) for all latitudes. For binning factors exceeding 20 (i. e. more than 400 pixels are binned together) the statistical model predicts a standard deviation of less than 5 × 10 11 molecules cm -2, while the measured data shows a standard deviation between 1 and 6 × 10 12 molecules cm -2 depending on the latitude. Thus, this remaining variance can be attributed to systematic effects. Around the equator (10 • S and 10 • N, Fig. 11b, violet square) the minimum standard deviation is almost 1 × 10 12 molecules cm -2, and this is chosen to be representative for the equator region. For the mid-latitudes below 50 • latitude (triangles in Fig. 11 ), the bands on 10 -30 • N (green triangle), 30 -50 • N (light blue triangle), and 10 -30 • S (yellow triangle), the standard deviation is in the range of 1-1. 5×10 12 molecules cm -2, while it is 2. 5×10 12 molecules cm -2 for the band 30-50 • S (light red triangle). We therefore choose the uncertainty to be 1. 5 × 10 12 molecules cm -2 for mid-latitudes. For the high-latitudes north of 50 • N (dark red diamond and dark blue circle) and south of 50 • S (blue diamond), there is a large difference. The northern high latitudes (dark red diamond and dark blue circle) show a standard deviation of 2. 5 × 10 12 molecules cm -2 and 2. 5 × 10 12 molecules cm -2, similar to the latitude band 30 -50 • S (light red triangle). The same latitude band in the south (50 -70 • S, blue triangle),

where the polar spring induces tropospheric BrO enhancements, which are not accounted for by the correction, yields a much higher standard deviation of 4. 0 × 10 12 molecules cm -2. For high-latitudes, we therefore choose an uncertainty of 2. 5 × 10 12 molecules cm -2 under normal conditions and an uncertainty of 4. 0 × 10 12 molecules cm -2 for polar spring conditions. 

The latitude band 70 -90 • S did not yield a meaningful gaussian fit due to the strong interference of enhanced tropospheric BrO event in polar spring. These uncertainties are then used as general systematic uncertainties in Table 4. 

It is noteworthy that for high latitude bands of 70 • S -50 • S, 50 • S -30 • S and 50 • S -70 • N, there is even an increase in the standard deviation for binning factors above 20. This is due to the lower number of total pixel in this region. Due to the limited light availability in these regions (solar zenith angle > 90 • ), which will lead to an even smaller number of data-points after binning. This effect can be seen for all latitudes, where a slight increase in the standard deviation compared to lower binning factors is found.
section_name: Quantifying the systematic effects

section_content: In addition to the systematic influences onto the BrO retrieval, the statistical uncertainties of the eight BrO DOAS fits also has to be taken into account in deciding which fit settings is optimal.
section_name: Proposed fit settings for the retrieval of BrO in volcanic plumes

section_content: In order to investigate the statistical performance of the BrO DOAS fit, the distribution of BrO VCDs in the equatorial region

is examined for all eight fit settings (see Fig. 12 ). The distributions are separated into the fits excluding the HCHO cross-section (SR, left), and the ones including the HCHO cross-section (HCHO, right), so that the same fit wavelength ranges have the same colors in both subplots. Increasing the fit wavelength range towards lower wavelength leads to a significant reduction in the statistical variation for both fit settings. For the SR fits the smallest fit range (SR 336 -360 nm, which is the baseline fit range from Sihler et al. ( 2012 )) shows a standard deviation σ of 1. 03 × 10 13 molecules cm -2, which subsequently decreases to 7. 49 × 10 12 molecules cm -2 (SR 332 -360 nm) and finally 6. 34 × 10 12 molecules cm -2 (SR 323 -360 nm) when the wavelength range is expanded. This is to be expected as a larger wavelength range increases the information content. 

In comparison to the SR fits, the HCHO fits have a larger standard deviation for all wavelength ranges, leading to a minimal σ of 6. 92 × 10 12 molecules cm -2 (compared to 6. 34 × 10 12 molecules cm -2 for the corresponding SR fit) and a maximum σ of 1. 71 × 10 13 molecules cm -2 (compared to 1. 03 × 10 13 molecules cm -2 for the corresponding SR fit). 

Thus, in order to achieve the lowest statistical uncertainties, the largest fit wavelength range excluding the HCHO crosssection should be chosen (SR 323 -360 nm). 

Looking at the whole globe for the fit range SR 323 -360 nm shows that the standard deviation is around 6. 5×10 12 molecules cm -2

for the mid-low latitude band between 70 • S and 50 • N (see Fig. 13a ). In the arctic region, the standard deviation is gradually decreasing towards high northern latitudes, which can be due to the higher AMFs in high latitudes, i. e. the variation of the SCD is divided by a larger number. In the Antarctic region, the standard deviation also decreases for a latitude of 80 • S. However, the distribution is not well represented by a gaussian distribution (see Fig. 13b, blue lines). There is a positive tail of elevated BrO VCDs due to enhanced tropospheric BrO columns in polar spring events. Furthermore, the distribution is not centered around zero, i. e. the correction term is offset. Concluding, the statistical variation of the fit SR 323 -360 nm can be best estimated by the standard deviation of 6. 5 × 10 12 molecules cm -2 or less for all latitudes.
section_name: Quantifying the statistical uncertainties

section_content: The investigation of the eight different fits can be summarized by the following:

-The statistical uncertainty decreases if (1) the wavelength fit range of the fit is increased and (2) the HCHO cross section is not included in the DOAS fit. 

-The systematic influence of clouds and the stratospheric background can be well corrected by the correction scheme described in Sect. 4. 1 for all eight fits. Thus, there is no preference for a fit range or fit setting. 

-The influence of HCHO is lowest for the fit SR 323 -360 nm. The BrO column shows no cross-sensitivity for moderate HCHO columns and the lowest increase for very high HCHO columns. 

Thus, the fit range SR 323 -360 nm is chosen for the BrO DOAS evaluation, since it shows the lowest statistical variation and is least susceptible to HCHO interference. The complete overview of the DOAS fit settings for this fit is listed in Table 3. 

The chosen best fit settings are recommened also for further studies on BrO in volcanic plumes.
section_name: Summary of the findings

section_content: In this section, the error estimations from the previous sections are summarised for all different cases identified throughout these sections and the total error is calculated. Three latitude regions are selected: Equator is the latitude band between 20 • S and 20 • N, mid-latitudes stretch from the equatorial region to 50 • N and 50 • S, and high latitudes are further north/south. In addition, there are two special cases distinquished: The polar spring in high latitudes with increased tropospheric BrO columns 510 and biomass burning events which are typically located close to the equator. Since biomass burning can seldomly also occur in mid-latitude or boreal summer, these cases are included as well. 

Finally, for each of these cases the statistical uncertainty and the systematic uncertainties are added quadratically: 10 ). The statistical uncertainties are derived for each latitude region in Sect. 5. 1 (depicted in Fig. 13 ). 

The added systematic, statistical, and the combined overall uncertainties are shown for all cases in Table 4. Under normal conditions (that means in the absence of high HCHO columns and polar tropospheric BrO enhancements), the combined uncertainty is between 6. 6 and 7. 0 × 10 12 molecules cm -2. Such conditions are representative for most of the volcanoes, e. g. 

the island volcanoes in the equatorial region (Galapagos, Vanuatu, Mayon) but also all mid-as well as high-latitude volcanoes such as Etna, Kilauea, and La Palma as well as Kamchatka, Iceland and the Aleutian volcanic regions. 

The volcanoes with a risk of a general influence of biomass burning are located in the equator region in the vicinity of rain forest. These are most notably the African rift valley volcanoes Nyiragongo/Nyamuragira, the northern Andes volcanoes -e. g. Nevado del Ruiz -and the Central American volcanic Arc volcanoes -e. g. Masaya and Fuego. Here, an uncertainty of 7. 4 × 10 12 molecules cm -2 has to be assumed. While mid-to high-latitude volcanic regions can also be affected by high HCHO columns which would increase the combined uncertainty by roughly 0. 8 × 10 12 molecules cm -2 to 7. 6 to 8. 2 × 10 12 molecules cm -2, this is expected only under very rare circumstances.
section_name: Estimation of the combined statistical and systematic uncertainties

section_content: This section shows the final selected retrieval of BrO inside volcanic plumes. The SO 2 DOAS retrieval is quickly explained (cf. Sect. 6. 1) as well as the additional local background correction (cf. Sect. 6. 2). The retrieval steps are: The proposed BrO DOAS fit in the chosen wavelength range of 323 -360 nm (cf. Sect. 5) is performed and the latitudinal background correction (cf. 4. 2) is applied. The plume is masked based on the SO 2 signal, whose retrieval is explained in Sect. 6. 1. This mask is used for a secondary local background correction (cf. Sect. 6. 2). The results of all the different steps is discussed in detail in Sect. 6. 3. 

6. 1 SO 2 retrieval SO 2 is one major constituent of volcanic gas plumes, accounting for 1-25 % of the volcanic gas emissions (Textor et al., 2004). 

In contrast to BrO, there is no significant stratospheric column of SO 2. As SO 2 has also differential absorption structures in the UV, it is widely used as a indicator for the detection of volcanic plumes (e. g. Hörmann et al., 2013; Fioletov et al., 2016; Carn et al., 2017). 

DOAS SO 2 retrievals from satellites were performed starting with the launch of the Global Ozone Monitoring Experiment (GOME) in 1995 (Eisinger and Burrows, 1998; Thomas et al., 2005; Khokhar et al., 2005) ; later DOAS retrievals were implemented also for the SCanning Imaging Absorption SpectroMeter for Atmospheric CHartographY (SCIAMACHY) (Afe et al., 2004), GOME-2 (Rix et al., 2012; Hörmann et al., 2013), the Ozone Monitoring Instrument (OMI) (Theys et al., 2015) and most recently for TROPOMI (Theys et al., 2017). As SO 2 columns in volcanic plumes vary over multiple orders of magnitude, recent algorithms introduced a scheme using multiple fits in different wavelength ranges with varying SO 2 absorption strength (Hörmann et al., 2013; Theys et al., 2015; Warnach, 2022). In this study, an SO 2 retrieval is employed purely to identify the volcanic plume for an accurate local background correction of BrO (see next section). Thus, only the fit range of 312 -324 nm, already sensitive to small SO 2 columns, is employed. 

This fit range is used in the SO 2 verification algorithm for TROPOMI (developed by C. Hörmann, MPIC, Mainz), and refined by S. Warnach, MPIC, within the verification for the upcoming Sentinel-5 mission (van Roozendael and the Sentinel 5 Verification team) as well as on synthetic spectra within the SO 2 verification for the upcoming Sentinel-4 mission (Wagner and the Sentinel 4 Verification team). In all these verification exercises very good agreement was found to the operational algorithm of S-5P/TROPOMI (Theys et al., 2017). 

The DOAS fit settings for the SO 2 retrieval includes the SO 2 absorption cross-section (203 K, Bogumil et al., 2003), and two Serdyuchenko et al., 2014). The other fit parameters (pseudo-absorbers), and the Fraunhofer reference are identical to those of the BrO fit. A detailed description of the SO 2 retrieval can be found in Warnach (2022).
section_name: Retrieval of BrO inside volcanic plumes

section_content: In order to correct for a potential remaining spatial pattern in the background distribution of SO 2 as well as BrO, a spatial polynomial correction is applied. This spatial polynomial was already successfully implemented as a local background correction in Hörmann (2013). Even though the BrO and SO 2 retrieval and correction schemes presented in this study are much less affected by background gradients due to the inclusion of Pukite-Terms for O 3 and the aforementioned cloud-ozone correction scheme, there are still potential small spatial patterns remaining, which should be further reduced. 

The retrieved VCD V m can be decomposed into the background VCD V b and the volcanic VCD V v :

The spatial polynomial correction scheme makes use of the fact that volcanic plumes can usually be considered as a localized phenomenon of several 100 km extent within a large scale background distribution pattern on the scale of 1000 km for both BrO and SO 2. Furthermore, both volcanic and background distributions are independent from each other and it can be assumed that the background is smooth with respect to the scale of the volcanic plume and background information gathered spatially around the plume can be interpolated over the plume region. 

The background VCD V b for each pixel k is modeled by two-dimensional polynomial of degree three:

where x is the pixel latitude and y the longitude corrected by a latitudinal dependent squeeze: In order to calculate the polynomial, first the volcanic plume is masked. For this a SO 2 VCD threshold is chosen, above which 575 a pixel is masked as influenced by the volcanic plume. This is chosen in this study as four times the standard deviation of SO 2 (for TROPOMI typically in the range of 2 × 10 16 molecules cm -2, compare Warnach, 2022). The polynomial coefficients c i,j

are then determined for this subset of data for both BrO and SO 2 independently and applied to each pixel of the complete data-set, leading to the background VCD. The resulting (volcanic) VCD is obtained by subtracting the background VCD V b from the retrieved VCD V m. 580
section_name: Local background correction

section_content: The complete sequence of processing steps of the retrieval of BrO for a volcanic plume is shown for the example of the An example of a small scale eruption at high latitudes is the plume of Sheveluch, Kamtchatka on 18 April 2019. While only low SO 2 VCDs are observed -not exceeding 8 × 10 16 molecules cm -2 (see Fig. 16a ) -there is still a clearly enhanced BrO

VCDs visible coinciding with the SO 2 pattern (see Fig. 16b ) in the order of several ×10 13 molecules cm -2. However, there is also elevated background BrO of up to 3 × 10 13 molecules cm -2 visible north of the plume over the island. These background structures are clearly identified by the correction term (see Fig. 16d ) and the resulting BrO VCD map (see Fig. 16c ) shows a smooth background of roughly zero. 

It is noteworthy that this is an example where the correction term actually shows a large gradient at the location of the plume (at 57 • N and 161 • -162 • E) where the correction term varies by more than 5 × 10 12 molecules cm -2. This gradient is caused by a cloud edge in this area (see Fig. B1c, d ).
section_name: Results

section_content: The plume of Ambrym volcano, Vanuatu on 30 September 2018

Another minor eruption where a BrO signal is detected is the eruption of Ambrym volcano on 30 September 2018 (see Fig. 17 ). 

While the SO 2 signal can be detected over several 100 km, maximum SO 2 columns do not exceed 1. 5 × 10 17 molecules cm -2 (see Fig. 17a ). Nevertheless, BrO VCDs exceed 3 × 10 13 molecules cm -2 indicating that the volcanic emissions are relatively BrO-rich (see Fig. 17b ). In contrast to the high-latitude Sheveluch eruption (see 6. 3. 1), the BrO background in this case is much weaker (see Fig. 17b ), since the plume is located close to the equator where the reference spectrum is taken. Furthermore, there are only low clouds in the area (see Fig. B2d ), hence cloud influence is low. Nevertheless, the correction term applies minor corrections for cloudy regions, e. g. for the high cloud at 14 • S 168 • E Fig. B2d, where the correction factor accounts for the negative BrO VCD (see Fig. 17d ), and slightly homogenizes the background (see Fig. 17c ). An example for a medium eruption, with a very localized volcanic plume is the Mt. Etna eruption on 29 January 2019 (see Fig. 18 ). This eruption is part of a phase of eruptive activity which started with weak degassing in December 2018 culminating in an explosive VEI-2 eruption on 24 December 2018 (Calvari et al., 2020) accompanied by lava flows (GVP, 2019). Activity subsided to strombolian activity throughout January and February 2019 (GVP, 2019). The eruptive event on 29 January 2019 falls within a week of increased activity during the last week of January (GVP, 2019). The event is a medium eruption, with 620 maximum SO 2 VCD of 5 × 10 17 molecules cm -2, but with a very localized plume (see Fig. 18a ). SO 2 VCDs exceeding
section_name: 6.3.2

section_content: We present a new retrieval scheme designed specifically for the accurate detection of BrO inside minor volcanic plumes from satellite spectra. Expanding the DOAS fit wavelength range to lower wavelength of 323 nm and thus a fit range of 323 -360 nm reduces the statistical uncertainty, without increasing any systematic influences. This is made possible by the introduction of 630 so called »Pukite« terms (Puk, īte et al., 2010), which account for second order effects of the strong absorption of O 3 at lower wavelengths. In addition, this wavelength range was shown to be least susceptible for a systematic influence caused by inter-
section_name: Conclusions

section_content: Acknowledgements. We thank ESA and the S-5P/TROPOMI level 1 and level 2 team for providing the L1B spectra and the NO2 L2 data. 

Moreover, we would like to thank them for the achieved high data quality of TROPOMI and its products. We thank DLR for financial support under funding number 50EE1811B. We thank the Max Planck Computing and Data Facility (MPCDF) for providing access to the high-performing computing cluster used for the DOAS evaluation and post-processing calculations. We thank Jonas Kuhn and Moritz Schöne, both at the IUP Heidelberg, and Steffen Ziegler, MPIC, for valuable discussions about the BrO DOAS retrieval. We thank Natural Earth (naturalearthdata. com) for making their 10m coastline and borders publicly available, which was used to create the coastline and borders in the Earth's maps included in this paper.

section_content: Data availability. The TROPOMI BrO VCDs (uncorrected, corrected and the correction terms) as well as the SO2 VCDs are available upon request.

section_content: Appendix A: Detailed data for systematic influences A1 Cloud and ozone influences BrO VCD as a function of CF and CH at the equatorial region on 1 October 2018 for all eight fit settings (see Fig. A1 ). 

The relation is stronger for larger fit wavelength ranges (left side) and seems also slightly increased for fits including HCHO (bottom row). The response is weakest for the fits with wavelength range 332 -360 nm both for SR and HCHO fits. Competing interests. The authors declare that they have no conflict of interest.