Smoke-charged vortices in the stratosphere generated by wildfires and their behaviour in both hemispheres : comparing Australia 2020 to Canada 2017

The two most intense wildfires of the last decade that took place in Canada in 2017 and Australia in 2019-2020 were followed by large injections of smoke in the stratosphere due to pyroconvection. It was discovered by Khaykin et al. (2020, doi:10.1038/s43247-020-00022-5) and Kablick et al. (2020, doi:10.1029/2020GL088101) that, after the Australian event, part of this smoke self-organized as anticyclonic confined vortices that rose in the mid-latitude stratosphere up to 35 km. Based on CALIOP observations and the ERA5 reanalysis, this new study analyzes the Canadian case and find, similarly, that a large plume penetrated the stratosphere by 12 August 2017 and got trapped within a meso-scale anticyclonic structure which travelled across the Atlantic. It then broke into three offsprings that could be followed until mid-October performing three round the world journeys and rising up to 23 km. We analyze the dynamical structure of the vortices produced by these two wildfires and demonstrate how they are maintained by the assimilation of data from instruments measuring the signature of the vortices in the temperature and ozone field. We propose that these vortices can be seen as bubbles of low absolute potential vorticity and smoke carried vertically across the stratification from the troposphere inside the middle stratosphere by their internal heating, against the descending flux of the Brewer-Dobson circulation.


Introduction
A spectacular consequence of large summer wildfires in midlatitude forests is the generation of pyro-cumulus clouds that can reach the lower stratosphere during extreme events (Fromm et al., 2010). The combustion products and accompanying tropospheric compounds (organic and black carbon, smoke aerosols, condensed water, carbon monoxyde, low ozone, ...) that are lifted to the stratosphere can survive several months and be transported over considerable distances (Fromm et al., 2010;Yu et al., 2019;Kloss et al., 2019;Bourassa et al., 2019) filling the mid-latitude band and sometimes reaching the tropics (Kloss et al., 2019). As black carbon is highly absorptive of the incoming solar radiation, the resulting heating produces buoyancy (Ditas et al., 2018) and an additional lift of the non diluted parts of the plume by several kilometers in the stratosphere (Khaykin et al., 2018;Yu et al., 2019). This effect enhances the dispersion and, by increasing the altitude, ensures a longer life time in the stratosphere (Yu et al., 2019), enhancing the radiative effect on climate that has been estimated of comparable amplitude to moderate volcanic eruptions (Peterson et al., 2018;. The Australian wildfire of the summer season 2019-2020 have exceeded in magnitude all previously known events . A striking discovery has been the observation that a part of the stratospheric smoke plumes self-organized as anticyclonic vortices that persisted between one to three months (Kablick et al., 2020; and, for the most intense one (nicknamed below as Koobor, following an aboriginal legend), rose up to 35 km , an altitude never reached by tropospheric aerosols since the Pinatubo eruption. It was conjectured by  that aerosol heating was essential in maintaining the structure and providing the lift. In turn, the vortex created a confinement that preserved the embedded smoke cloud from being rapidly diluted within the environment.
It is a natural extension to investigate the occurrence of such vortices after previous wildfires, in particular over the arXiv:2011.13239v1 [physics.ao-ph] 26 Nov 2020 2Lestrelin, Legras, Podglajen & Salihoglu: Smoke-charged vortices in the stratosphere: Canada 2017 vs Australia 2020 last 15 years during which the required satellite instruments are available. As the strongest previously recorded wildfire of the last decade is the 2017 Canadian event that took place in British Columbia (Hanes et al., 2019), this work is revisiting this case which has already been documented in several studies (Khaykin et al., 2018;Ansmann et al., 2018;Peterson et al., 2018;Yu et al., 2019;Kloss et al., 2019). In particular, a stratospheric rise of up to 30 K day −1 was diagnosed based on satellite observation by Khaykin et al. (2018) who also reports a compact smoke cloud at 19 km over the Haute-Provence observatory, Southern France, on 29 August 2017. An other goal of this study is to complement  by expanding their diagnostics and interpretations on the 2020 case.
Section 2 describes the data and methods used in this work. Section 3 describes the new vortices found after the 2017 Canadian fire and their evolution. It includes a detailed discussion of previous results. Section 4 describes the structure of the vortices based on the 2017 and 2020 cases. Section 5 offers conclusions.
2 Data and methods

Satellite data from CALIOP
Launched in April 2006, the CALIPSO mission (Winker et al., 2010) took on board the Cloud-Aerosol Lidar with Orthogonal Polarization (CALIOP) which is a two-wavelength polarization lidar that performs global profiling of aerosols and clouds in the troposphere and lower stratosphere. We used the total attenuated 532 nm backscatter level 1 (L1) product in the latest available version V4.10 (Powell et al., 2009). The nominal along track horizontal and vertical resolutions are, respectively, 1 km / 60 m between 8.5 and 20.1 km and 1.667 km / 180 m between 20.1 and 30.1 km. The L1 product oversamples the layers above 8.5 km with an uniform horizontal resolution of 333 m. We computed the scattering ratio by dividing the total attenuated backscatter by the calculated molecular backscattering following Vernier et al. (2009) and Hostetler et al. (2006) using the meteorological metadata of the L1 product. To reduce the noise, an horizontal median filter with 40 km width (121 pixels) was applied to the data. In order to separate clouds from aerosols, we also used the level 2 (L2) total scattering aerosol coefficient at 532 nm which is available at 5 km resolution. Both daytime and nighttime measurements have been used. Since the daytime measurements are noisier, they can only be used when the aerosol signal is strong enough, mainly during the first weeks after the release of the plume.
As in , CALIOP inspection has been the first step in order to identify the potential vortices. The L1 sections have been systematically screened from 12 August to 15 October 2017, selecting those which contained isolated compact patterns above 11 km, identified as aerosols by the L2 product. The location and size of the retained patches were then determined by visually matching a rectangular box to the observation, as illustrated in Fig. A1. It is usually very easy to see the boundaries. At the next stage, the patches have been associated to the vortices detected as described below and a further inspection has rejected cases (less than 20% of those retained at the first stage) that correspond to tails left behind by vortices.
Due to an increased solar activity, CALIOP operations have been suspended between 5 and 14 September. This has been an obstacle to establish a continuity based on CALIOP observations alone but the tracking of the vortices filled that gap as described below.

Reanalysis
To track the stratospheric wildfire vortices and diagnose their dynamical structure, we used the ERA5 reanalysis (Hersbach et al., 2020) which is the last generation global atmospheric reanalysis of the European Center for Medium Range Weather Forecast (ECMWF).  used instead the ECMWF operational analysis and forecast. Both are based on the ECMWF Integrated Forecast System (IFS). In the ERA5, the native horizontal resolution of the IFS model is about 31 km and it has 137 levels in the vertical with spacing that varies from less than 400 m at 15 km to about 900 m at 35 km, within the relevant altitude range of this study. We used an extracted version of the ERA5 with 1 • resolution in latitude and longitude, at full vertical resolution and 3-hourly resolution. This choice was dictated by practical considerations and by the fact that 1 • is sufficient to describe synoptic-scale features and that vortices do not travel across more than a few grid points in 3 hours.  used temperature, ozone and vorticity to investigate the structure of the vortices. In this study, we also used potential vorticity (PV). PV is a Lagrangian invariant for inviscid and adiabatic flows (Ertel, 1942); furthermore, it provides a compact and complete picture of the balanced part of the flow (Hoskins et al., 1985). An interesting property of the ERA5 is that its dynamical core preserves PV much better than previous reanalyses (Hoffmann et al., 2019). Although ERA5 potential vorticity can be directly retrieved from the ECMWF archive on a given set of potential temperature levels, we instead recalculated it from the retrieved vorticity, temperature and total horizontal wind fields on model levels, in order to benefit from the full vertical resolution offered by ERA5. For that purpose, we used the definition of PV for the primitive equations in spherical coordinates and model hybrid vertical coordinate η, which is given by : where a is the radius of the Earth, λ is the longitude, φ the latitude, g is the free-fall acceleration, ζ η is the model level vertical component of the vorticity, f the Coriolis parameter, (u, v) the horizontal velocity, p the pressure and θ the potential temperature. The gradients are estimated using centered differences on the retrieved longitude-latitude-η grid. See Eq. (3.1.4) of Andrews et al. (1987) for an analogous formula in barometric altitude coordinate. While the formulation of potential vorticity in Eq. (1) is the most commonly used, it bears the disadvantage of a large background vertical gradient, which proves inconvenient to track and characterize structures along their ascent. To overcome this issue, we used the alternative formulation of Lait (1994), discussed by Müller and Günther (2003): where P is the Ertel PV, = 4 in the Australian case and = 9 2 in the Canadian case, and θ 0 = 420 K. Compared to P , Π is still an adiabatic invariant and exhibits a reduced vertical gradient, so that the vortices characterized by anticyclonic Π anomalies can be unambiguously distinguished during their ascent.

Assimilation increment
The ERA5 is constrained by observations over repeated 12hour assimilation cycles. The assimilation increment is defined over each cycle as the difference between the new analysis and the first guess provided as a final stage of a free forecast run of the model initialized from the previous analysis 12 hours before. In the ERA5, the assimilation increments can be calculated on each day at 6:00 and 18:00 UTC. In order to diagnose how the observations are forcing the vortices, we calculated the assimilation increments of temperature, vorticity, potential vorticity and ozone. Temperature and ozone determine the radiances that are measured by spaceborne instruments and are also directly accessible from in situ instruments. On the contrary, potential vorticity cannot be directly retrieved from any instrument. These three parameters are updated by the assimilation system in order to reduce the difference between observed quantities (typically radiances but also deviations of the GPS signal path) and simulated quantities (radiances that a satellite flying "above the model" would see). It is tempting to see the temperature assimilation increment as an additional heating but this is incorrect as the increment is calculated from an adjusted state, resulting from the iterations of the assimilation, in which both temperature and motion respond to the forcing by the observations that depend only on the temperature.
Neither the ECMWF operational analysis nor the ERA5 assimilate aerosol observations. The smoke plumes are hence totally missing in the IFS where stratospheric aerosols were only accounted by mean climatological distribution during the periods of investigation.

Vortex tracking approach
Once we caught a first occurrence of an isolated bubble of aerosols with CALIOP, we searched the ERA5 data for the occurrence of a corresponding vortex that showed up as an isolated pattern in the Π and ozone maps. While  used only relative vorticity to track the vortices for commodity, we used both Π and the ozone anomaly defined as the deviation with respect to the mean at the same latitude and altitude. The tracking was made every six hour by following a local extremum within a box of usually 12 degrees in longitude, 5 degrees in latitude and a range of at least 30 K of potential temperature in the vertical. Once a vortex was caught, the box was moved forward to the step n + 1 according to the vortex motion between steps n − 1 and n. In very few instances, the tracking was guided by reducing the size of the box. This was needed in particular at the birth of a vortex, or when it split in two parts and in the final stage when one of the methods lost the track before the other.

General description
Although the 2017 Canadian fire started in June (Hanes et al., 2019), it is only by 12 August that the fire reached an intensity such that a large PyroCb developed and reached the lower stratosphere leaving a smoke plume that could be followed by satellite sensors (Khaykin et al., 2018;Peterson et al., 2018). From the inspection of scattering ratio along CALIPSO orbits, we found our first distinctly isolated cloud on 16 August at 9:38 UTC and 12.6 km altitude, over the north of Canada at 63 • N and 102 • W, where the tropopause altitude was 11 km. In the following days and until the CALIOP interruption on 4 September, we could track this cloud, labelled as bubble O, and its offsprings almost every day. Figure 1 shows several typical patterns during that period.
As indicated by the longitudes on Fig. 1, the bubble O moved eastward across the Atlantic. During the early days, it was captured each day by at most two CALIPSO orbits (one day orbit and one night orbit), the adjacent orbits only showing filamentary non compact structures that could be easily distinguished. The bubble rose rapidly reaching 18 km by 25 August (i.e. an ascent rate larger than 0.5 km day −1 ). One surprising feature is that once it reached the eastern Atlantic coast by 27 August, several bubbles could be tracked on multiple close orbits, or even exceptionally, on the same orbit. As we shall see in the following, this corresponds to the dissociation of bubble O into several offsprings that we label as A, B1 and B2 for which several views are shown in Fig. 1. Another issue came soon after as CALIOP operations were suspended for the period 5 -14 September due to increased solar activity. When CALIOP observations resumed after 14 September, the A, B1 and B2 bubbles that were now located at 20 km or above could be found again, over central Asia for A and over the Pacific for B1 and B2, and they were followed during their subsequent journey until mid October. Figure 2 shows a selection of views during that period.

Early evolution
The smoke bubbles have been attributed to vortices based on the ERA5 reanalysis which were available during their whole life cycles. Starting on 14 August, a kernel of almost zero PV and low ozone could be found by 72 • N and 115 • W at an altitude of 12 km just above the tropopause. In the following days, this anomaly developed vertically and connected with the bubble O location as identified from CALIOP (see animation in Sec. S2 of the Supplement). By 17 August, as it crossed the Hudson Bay, it exhibited a well developed intrusion that reached 14 km in the PV longitudinal and latitudinal section. As seen in the animation, this intrusion, while still rising, was subsequently stretched by the vertical shear and split by 19 August into an upper part isolated in the stratosphere at about 15 km while the lower part near the tropopause was further stretched and disappeared. The upper part which was associated with bubble O was tracked in ERA5 ozone and PV fields and from CALIOP until the end of August (Fig. S1 of the Supplement). As it crossed the Atlantic, it got trapped inside a thalweg by 20 August and travelled with it until it reached the European coast. Due to the wind shear prevailing in the associated jet, the vortex O elongated in latitude across the isentropes (see the animation) until it got dissociated over western Europe into the three parts A, B1, B2, as described in the next section. Figure 3 displays the series of event that led to the dissociation of the vortex O over the period 22 August -1st September into its three offsprings that were subsequently followed during one month and half. The CALIPSO orbits which are displayed here are more numerous than those retained for the tracking of the vortices as we are also interested in the tails of PV patches. The presence of a smoke bubble detected by CALIOP along an orbit is shown as a red segment. We see that on all orbits, there are bubbles that match each of the intersections with the low PV patches. The sequence begins, on 22 August ( Fig. 3a) with the elongated structure of the vortex O that emerged at 420 K on the eastern flank of the thalweg within which it crossed the Atlantic. The elongation is due to the intense vertical shear in the jet stream. The initiation of vortex A is already visible as the wrapping up pattern on the north east end of vortex O. On 23 August (Fig. 3a), the pattern of vortex A is lost at level 420 K but is now visible on level 435 K. Two days later ( Fig. 3c), the vortex A was a developed structure reaching 60 • N at level 455 K which was overpassed by CALIPSO ( Fig. 4a). At the same time the vortex O maintained a core near 45 • N (Fig. 1d). In the following days, on 26 and 27 August ( Fig. 3d-e), the vortex A fully separated from vortex O as it moved eastward and rose (Fig. 1e). On 28 and 29 August ( Fig. 3f-g), the vortex A moved away while the vortex O was again elongated and started to split into a western part (Fig. 4b) and an eastern part (Fig. 4c). On 30 and 31 August ( Fig. 3h-i), the western part separated while the eastern part folded itself and separated in two more parts that provided vortex B1 as the north component ( Fig. 1) and vortex B2 as the south component (essentially the continuation of the vortex O) which were both seen on the same CALIOP overpass of 1st September (Figs. 3j & 4d). The vortices B1 and B2 fully separated in the following days while rising and starting to move slowly eastward. The western component could also be followed until September 4, accompanied by patches seen by CALIOP. It remained however below 460 K and could not be linked with certainty to any structure seen from CALIOP after 15 September.

Late evolution
During the observations period that followed the recovery of CALIOP after 15 September and until mid-October, almost all observations of compact smoke bubbles at 20 km and above could be attributed to one of the vortices A, B1 and B2. Those which were discarded are under the shape of filaments or tails at altitudes corresponding to the bottom of the bubbles. As in Fig A perfect match of the horizontal location is not expected as there is no reason for CALIPSO orbit to cross every time the bubble by its center. Nevertheless, we see a very good agreement between the ERA5 trajectory and the location of the 27 smoke bubbles seen by CALIOP that are attributed to vortex A, and the same holds for the two other vortices B1 and B2 (Figs. S2 & S3 of the Supplement) with, both, 21 attributions. The evolution of vortices A, B1 and B2 is also available as animations in Sec. S3 of the Supplement. Vortex A was the first to separate from the mother O by 22 August and moved first eastward until it reached central Asia by 1st September Figure 3. Sequence of PV charts showing the dissociation of vortex O into vortices A, B1 and B2. The map is plotted on the potential temperature surface corresponding to the core of vortex O or its continuation B2 in the ERA5 tracking. The orange lines are plotted at the isentropic level of vortex A, specified within parenthesis, and show the contour of -2 PVU for vortex A. The black, green and purple lines show the intersecting orbits of CALIPSO and the red segments show the parts of the orbits occupied by the bubble (only the bubble core extent is displayed here). The maps are plotted at the hour matching best the selected CALIOP occurrences. The purple orbits are those corresponding to the sections shown in Fig. 1. The green orbits are those corresponding to the sections shown in Fig. 4. Day-time and night-time orbits are used, the day-time orbit are running from south-east to north-west while the night-time orbits are running from north-east to south-west. Starting from 1st September, the remaining of the vortex O is relabelled as B2. . at a longitude of 60 • E and a latitude of 55 • N. It then got trapped into the region of slow motion that extends between the two centers of the Asian monsoon anticyclone (AMA), and started to drift slowly southward while staying at about the same longitude and maintaining its ascent rate. A week later it reached 35 • N where is was caught into the easterly circulation and started to move westward crossing Africa and continuing its path which could be tracked until the west Pacific. Figure 7 shows a composite image of successive images of the localized ozone hole from the ERA5. It was easier to track the vortex using the ozone field than the PV. In particular the PV signal almost vanished as it passed over Africa (see the vortex A animation in the Supplement) while the ozone signature WaS always very clear. We attribute this vanishing to the strong infrared emissivity of the Sahara that limits the sensitivity of the IASI sounders which are important to sense the thermal signature of the vortex . As ozone leaves also a signature in the ultraviolet seen by other assimilated instruments such as GOME-2 , it is less affected. The fact that the PV signature reintensifies as soon as the vortex is over the Atlantic supports this hypothesis. During the last stages of the vortex, by mid-October, we also see a premature loss of the PV signal while the ozone signal is still detectable. This pattern is shared by the two other vortices B1 and B2 and it differs from the 2020 case where such effect is not observed for any of the three vortices. Besides the increase of the IASI fleet from two to three instruments, we do not see any drastic change in the observation system between the two events. While the vortex A was completing its transition to the tropics, the two other vortices B1 and B2 travelled eastward within the westerly flow on the north side of the AMA reaching together the Pacific on 18 September. The vortex B1 crossed the Pacific at mid latitude and got lost near Hudson Bay after crossing most of North America by 14 October. The vortex B2 travelling at higher latitude completed a full round the globe travel during the same period and got lost over central Asia by 11 October. Figure 6 summarized the trajectories of the vortices O, A, B1 and B2 from their birth to their loss. The total paths of the four vortices are 13000, 42400, 28500 and 33400 km, respectively. They rose up to 23 km for vortex A and to 21 km for vortices B1 and B2.

Comparison with previous studies
Several previous studies have discussed various stages of the smoke cloud evolution described above although none made the link with a PV structure. A number of comments are here in order: -The ascent rate was the strongest during the initial stage of vortex O when it rose from its origin just above the tropopause. Previous studies based on CALIOP or limb sounding instruments reported an ascent rate using the upper envelop of the bubble. Yu et al. (2019) reported an ascent from 14 to 20 km from 15 to 24 October and  Khaykin et al. (2018) reported an ascent of 30 K of potential temperature per day, that is 3 km per day, between 16 and 18 August. We used the centroid of the PV and O3 anomalies to define the ascent and found (see Fig. S2  We investigate here the structure of the vortices by the mean of a composite analysis. The dynamical fields surrounding the vortex centroids were regridded regularly in Cartesian geometry within a moving frame which follows the centroids ascent and meridional displacement. As in , the composited model fields were then averaged in time to generate a composite analysis that filters out noise and vacillations unrelated to the vortices and emphasizes their mean structure. Figure 8 depicts the composite of Lait PV Π (defined in Eq. (2)) and temperature anomaly following vortex A (bottom panels) and the main vortex Koobor (top panels) generated by the 2020 Australian fires and tracked by . In the natural coordinates, without expanding the vertical direction, both vortices appear as isolated pancakes of anomalous anticyclonic Lait PV (i.e., negative in the northern hemisphere and positive in the southern one). The analysed temperature anomaly consists in a vertical dipole surrounding the PV with negative temperature anomaly above and positive temperature below it, where the centers are located on the upper and lower edges of the Lait PV distribution. As noted by , this relationship between temperature stratification and vor- Figure 7. Sequence of ozone anomalies along the trajectory of vortex A for selected dates from 31 August to 11 October. The background is the mean ozone anomaly over this period at 460 K. For each selected day, the box shows the ozone anomaly on that day at at 0:00 UTC at the level of the vortex centroid according the the ozone anomaly tracking, reported along the blue line. ) and e) and c) and e) are respectively horizontal, zonal and meridional sections through the centroid of the vortex. The green lines are contours of anticyclonic vertical vorticity (corresponding to 3 · 10 −5 and 5 · 10 −5 s −1 for the top panels, and −1.5 · 10 −5 and −2.5 · 10 −5 s −1 for the bottom panels). The black contours are temperature anomaly with respect to the zonal mean. Note that the displayed horizontal range is reduced by a factor of 2 for vortex A.
ticity is qualitatively consistent with the thermal wind balance, and characteristic of anticyclonic vortices in the quasigeostrophic (QG) equations (Dritschel et al., 2004), but also of balanced anticyclones in the full primitive equations (e.g. Hoskins et al., 1985). In the case of the smoke-charged vortices investigated here, it turns out that their intensity is slightly beyond that of the typical geostrophic flow. The typ-ical Rossby number of the structure can be expressed as follows: with U the maximum horizontal wind speed, L h the horizontal length scale (defined as the diameter of the ring of wind speed maximum) andζ = U/L h the average relativite vorticity. Ro is about 0.06 in the 2017 case and 0.35 for 2020. While in 2020 the Rossby number is beyond the typical QG regime (Ro < 0.1), the aspect ratio α, defined as (L z being the vertical extent of the vorticity contour at maximum wind speed) obeys the stratified QG scaling α f N in both cases, despite the 2017 vortex A being about 8 times smaller in volume than its gigantic 2020 counterpart.
To investigate further the dynamical regime in which the identified vortices evolve (or their representation in the IFS), Table 1 presents their typical sizes, amplitude characterized by the Rossby number Ro and the Froude number: as well as the absolute vorticity amplitude at the vortex centroid normalized by f : which measures the inertial (in)stability of the flow. It should be noted that despite the similarity, ζ a f is not redundant with Ro defined in Eq. 3 and characterizes the extremum rather than the structure average of the vorticity.
Keeping in mind the limited vertical and horizontal resolution of the IFS, it can be noticed that the Froude and Rossby numbers are always of the same order, leading to a Burger number Bu ∼ 1 as typically encountered in geophysical flows. Furthermore, most vortices in Table 1 obey the QG aspect ratio, so that they are quasi-spherical (see Fig. 8) when the vertical coordinate is stretched by a factor N f . This observation is consistent with numerical studies of ellipsoidal vortices, which have demonstrated that the higher stability of quasi-spherical vortices (Dritschel et al., 2005) and the tendency of aspherical vortices to relax towards sphericity in the stretched coordinate system x, y, N f z (Tsang and Dritschel, 2015).
Contrary to its aspect ratio, the magnitude of the vortex perturbation is case-dependent and does not necessarily fit in the classical QG scaling, thereby contrasting with the synoptic scale circulation. Indeed, typical vortex-averaged Rossby numbers range from 0.06 up to 0.35, a range similar to observed mesoscale and submesoscale oceanic eddies (e.g. Le Vu et al., 2017). For the intense 2020 Koobor vortex, the maximum vorticity in the vortex core is on the verge of inertial instability or even slightly beyond its threshold for linear flows, as can be seen from the positive Π values in Fig. 8 (see also the small negative values of ζ a f in Table 1). This property was maintained over one month and half from the birth of the Koobor vortex to its first breaking (see Fig. S4 of the supplement). Although it is not clear how realistic the ECMWF vortices are, their amplitude is likely underestimated since they are forced by data assimilation and the standalone model does not simulate them. We can therefore not generally conclude on the amplitude of the perturbations, in particular the smaller NH vortices. Nevertheless, the 2020 cases demonstrate that the vorticity anomaly may reach the threshold for inertial instability at which it likely saturates.
Related to their different magnitudes, the vortices also have distinct impacts on their immediate surroundings, as can be identified in Fig. 8. While the 2020 Koobor vortex was strong enough to generate a significant cyclonic PV anomaly eastward of the vortex centroid which rolls up around it on its northern edge, no such signature can be distinguished around vortex A. In a zonal plane, the 2020 negative PV patch has diameter comparable to the main vortex but reduced magnitude. The existence of this negative anomaly can be attributed to the equatorward advection of PV on the eastern flank of the main vortex and is likely responsible for the equatorward beta drift undergone by the main vortex during its ascent.

Diagnosing diabatic tendencies
By comparing the main 2020 smoke-charged vortex in ECMWF analyses and model forecasts,  have shown two systematic shortcomings of the free-running model with respect to the analysis, namely: a failure to reproduce the observed ascent of the structure at about 6 K/day in potential temperature: on the contrary, the modeled anticyclones tend to remain at constant θ; a decay of the vorticity anomaly within one week, while the observed vortex survives for more than 3 months.
This suggests that data assimilation of observed temperature and ozone profiles is a necessary ingredient to both ascent and the maintenance of the vortex in the model, while some physical or spurious processes act to dissipate the structure. In the atmosphere, it is assumed that the forcing is exerted by radiative heating through solar absorption by black carbon aerosols within the smoke (Yu et al., 2019). Increments are thus expected to play substitute in the analysis system for the diabatic tendency maintaining the vortices, as the wildfires smoke is missing in the model. Figure 9 (a, d) depicts the average composite of the ERA 5 total heating due to physics calculated over the forecast. This field is dominated by the long-wave radiative heating component since the short-wave absorption by the smoke is missing from the IFS. The dominant feature shown by Fig. 9 is a damping of the temperature anomaly T with respect to radiative equilibirum temperature which can be cast into a Newtonian relaxation as:

Temperature and vorticity tendencies
where τ rad is the radiative damping rate. Figure 9 suggests τ rad 6 − 7 days. This is consistent with the life time of the structure in the ECMWF forecast which is about one week , suggesting that radiative dissipation plays a major role in the decay of the vortex in the model but also in the real atmosphere. However, since the model itself cannot sustain the vortices, it comes to data assimilation to maintain the structure. It was demonstrated by  that the IFS extracts information from the thermal signature and the ozone anomaly detected by satellite instruments. The temperature increments are shown in Fig. 9 (b, e). In contrast with the heating rates, which are instantaneous forcing of the temperature by physical processes, the temperature increments result from the balanced dynamical response of the temperature field to the forcing induced by the discrepancy of measured radiances with respect to their estimated values in the assimilation procedure. Hence, the temperature pattern does not match the observed aerosol bubble colocated with the vortex which would be expected for the SW aerosol absorption. It rather exhibits a multipolar structure, essentially oriented vertically. In the 2017 case with limited vertical ascent rate, the T increment vertical dipole mainly cancels the radiative damping. In 2020, it is rather a tripole which enables the maintenance but also the ascent of the dipolar temperature anomaly structure.
Due to this dynamical adjustment, temperature is not the only field incremented via data assimilation. Figure 9 (c,g) shows the increments of relative vorticity ζ. They contribute again to the maintenance (2017 case) and ascent and maintenance (2020 case) of the vortices. Moreover, similarly to the composite vortex, the composite vorticity and temperature increment tendencies exhibit the structure of a balanced vortex: the 2017 anticyclonic ζ tendency monopole is sandwiched between the two extrema of temperature increment dipole whereas in 2020 the tripole of T increment alternates with the two layers of vorticity anomalies. Overall, the T and ζ patterns of the 2020 assimilation increment are qualitatively consistent with the expected effect of a localized heating in a rotating atmosphere (as sketched in fig.10 of Hoskins et al., 2003).
Together with the apparently balanced vortex structure of the increments, the similarity with the theoretical response in Hoskins et al. (2003) suggests that potential vorticity (or Lait PV) is the appropriate field to consider for a more straightforward interpretation of the increments in terms of missing diabatic processes. This is the focus of the next subsection.

Lait PV and ozone tendencies due to assimilation increments
Composite time average increments of Lait PV Π are presented in Fig. 10. As stressed above, Π increments result from a combination of both temperature and vorticity and appear due to the implicit incorporation of temperature observations into an approximately dynamically-balanced response by the 4D-Var scheme used in the ERA5. Contrary to the mean vortex structure in Fig. 8, we notice a stark contrast between the NH and SH vortices. In the 2017 case, the Π increment shows a dominant monopole structure whose extremum is slightly shifted upward (by about 500 m) with respect to the Π extremum (Fig. 8). The effect of the increment is thus mainly to reinforce the existing Π anomaly and, secondarily, to drive its ascent. In contrast, the 2020 case shows essentially a vertical dipole structure which forces an ascent of the Π center. The dipole is furthermore tilted along a north-west to south-east axis. A qualitative analysis (Appendix B) shows that the observed westward tilt of the increment dipole emerges as a consequence of the negative zonal wind shear encountered by the main 2020 Koobor vortex along its ascent. More precisely, it is the combination of this shear and the distribution of the increment tendency over the 12 hours analysis time window which is responsible for the tilt.
The different patterns of Lait PV Π (Fig. 8) and its increment (Fig. 10) around the vortex are mimicked in the ozone  Figure 9. Time average composite sections of ERA 5 total heating rate (a, c), increment-induced temperature (b, e) and vorticity tendency (c, f) following two selected smoke charged vortices, the major vortex from the 2020 Australian wildfires described in  (a, b, c) and vortex A (d, e, f). The green lines are contours of anticyclonic vertical vorticity (corresponding to 3 · 10 −5 and 5 · 10 −5 s −1 for the top panels, and −1.5 · 10 −5 and −2.5 · 10 −5 s −1 for the bottom panels). The black contours are temperature anomaly with respect to the zonal mean. Note that the displayed horizontal range is reduced by a factor of 2 for vortex A. field, as shown in Fig. 11. As described above, the ascending anticyclonic vorticity anomalies are accompanied by negative ozone anomalies (Fig. 11 a and c). Compared with Figs. 8 and 10, Fig. 11 shows that the patterns of the ozone anomaly bubble and its increment are very similar to those of Π. Overall, model increments tend to (i) counter the dissipation and (ii) drive a cross-isentropic vertical motion of PV and ozone anomalies.
At this point, it is tempting to interpret the low ozone, low absolute PV anticyclones as resulting from the vertical advection of smoke-charged tropospheric air bubbles conserving both their low absolute PV and their tropospheric tracer content during their ascent in the stratosphere. In the case of the main 2020 plume, this perception is broadly consistent with the behavior of PV and PV anomaly at the center of the vortex, which respectively remain constant or exhibit a steady increase related to the large vertical gradient of PV within the stratosphere (Fig. S4 of the supplement). It is also consistent with our analysis of the developing stage of the mother vortex O in Sec. 3.2. We note, however, that the analogy between inert tracer and PV transport arises from their Lagrangian conservation in the adiabatic inviscid case and breaks in the presence of diabatic processes (Haynes and McIntyre, 1987). The maintenance of this relationship hence calls for dedicated theoretical and numerical investigations.

Conclusions
The generation of smoke-charged vortices rising in the southern stratosphere was discovered following the Australian wilfires at the turn of 2019 . We find here that similar events took place in 2017 in the northern stratosphere after the British Columbia wilfires that culminated in early August 2017 when a large plume of smoke and of low potential vorticity air has been injected in the lower stratosphere on 12 August 2017. We show that soon after, a vortical structure developed in the plume and started to rise in the stratosphere. During the following days this smoke vortex became an isolated bubble and was transported eastward across the Atlantic while getting elongated by the wind shear  Figure 11. Same as Figure 9, but for ozone anomalies from the zonal mean ((a), (c)) and ozone increments ((b), (d)). and reaching 19 km in altitude at its top. Subsequently, that structure dissociated over western Europe into three separate vortices that could be followed until the middle of October 2017. Two of them kept moving eastward in the middle latitudes and performed a round the globe travel, rising up to a potential temperature of 530 K (21 km). The third one has transited to the subtropics across the Asian monsoon anticyclone and moved eastward until Asia. It also rose higher, to 570 K (23 km).
Our study demonstrates again the ability of the advanced assimilation system used at the ECMWF which is the basis of the ERA5 to exploit the signature left by smoke vortices in the temperature and ozone fields and reconstruct the balanced vortical structure. The balance is produced during the iterations of the assimilation and therefore the assimilation increment reflects this balance by providing an increment on the wind as well as on the temperature. The increment is consistent with the expected effect of a localized heating due to the radiative absorption of the smoke. It fights the longwave radiative dissipation and moves the vortices upward. It is quite likely, however, that the detection of vortices by the assimilation system is limited by the sensitivity of the satellite instruments to the disturbances in the temperature and ozone fields. There are also detection limits in CALIOP and limited chances to overpass small-scale structures if they are sparse enough. It is therefore possible that a number of small-scale vortices escape the direct methods used in this study. It is also possible that dynamical constrains, like the ambient shear, limit the existence of such vortices. It is thus difficult to conclude from our study what is the generality of such structures and their global impact. It is however quite clear that they provide an effective way for smoke plumes to stay compact and concentrated inducing a rise of the order of 10 km or more in the stratosphere that, in turns, increases the life time and the radiative impact mentioned in previous studies (Ditas et al., 2018) from months to years knowing the character of the Brewer Dobson circulation (e.g. Butchart, 2014).
Long lived anticyclones dubbed as "frozen-inanticyclones" (FrIAC) by Manney et al. (2006) have already been reported in the Arctic stratosphere. They share a long life span with the smoke anticyclones but they differ in many other respects. The FrIAC are deep barotropic structures extending from 550 K to 1300 K and are observed in the polar region after the breaking of the winter polar vortex. They are generated by the isentropic intrusion of low latitude air with low PV (Allen et al., 2011;Thiéblemont et al., 2011). Instead, the smoke vortices exhibit a strong baroclinic structure with a vertical temperature dipole and are propelled through the layers of the stratosphere by their internal heating bringing tropospheric air to the middle stratosphere at latitudes where the Brewer circulation does the opposite (Butchart, 2014).
As volcanic plumes in the stratospheric are usually made of secondary sulfuric acid aerosols that are considerably less absorbing than the black carbon, raising compact plumes are not expected to be seen in the stratosphere after extratropical volcanic eruptions. Such a case, however, has been reported after the 2019 eruption of the Raikoke (Chouza et al., 2020;Muser et al., 2020). This event which displays a slow rise and is not associated with a vortex in the ERA5 might fall below the detection limit but opens the question of the possible role of heating by volcanic aerosols.
The conditions to maintain the stability of the smoke vortices and those leading to their final dilution are also not yet understood and this work is only opening a path to be explored. Appendix A: Tracking of the smoke bubbles Figure A1 shows two screenshots as examples of the interactive tracking of the smoke bubbles in CALIOP sections of the scattering ratio. The first case (Fig. A1a) is taken from bubble O on 19 August when it was elongated within the Atlantic thalweg. The second case (Fig. A1b) is taken from bubble A on 24 September while it was moving westward over the Atlantic. The interactive method used to perform the tracking is displayed in a an animation in Sec. S1 of the Supplement.
Appendix B: Increment structure around the ascending vortex in a sheared background wind The average composite of Lait PV increment for the 2020 anticyclone (Fig. 10) exhibits a tilted dipole structure, whose inclination by far exceeds the slight dip observed for the vortex axis and figured by the green vorticity contours in Fig. 10 b. This tilt causes the dipole to appear not only in vertical panes but also in the horizontal one (Fig. 10 a). From a purely kinematic point of view, the dipolar structure of the increment shows that the forecast model systematically underestimates both upward and westward motions. While it is unlikely that the ECMWF forecast model underestimates zonal wind magnitude over such a long period, this missing westward displacement filled by the increments may be attributed to missing vertical ascent in an easterly sheared background wind. This can be grasped by comparing trajectories advected by a sheared flow with or without an ascent. The two situations are contrasted for a simplified case in Tab. B1. Over the finite analysis time window ∆t, the vertical velocity causes the parcel (here, the vortex) to be advected in a region of different mean wind and drift away from its non-ascending cousin. Hence, the additional motion brought by the ascent has both vertical and horizontal components, which can be expressed as: The tendency of Π associated with a translation at speed (Ẋ,Ż) incr is: If we consider the isolated bubble of maximum potential vorticity anomaly P max and length scales L x and L z = αL x , the magnitude of the tendency required to translate it along x and z are then: and their ratio γ is: In the case of the 2020 vortex, α 7.8 10 −3 , Λ 3. 10 −3 s −1 and the assimilation window is ∆t = 0.5 day so that γ 0.5, which is quantitatively consistent with Fig. 10, panel b).
When the shear is properly taken into account, the westward drift induced by the increments can thus be reconciled with the idea that the ascent is the main process lacking from the model.
Author contributions. The four authors performed together the analysis. HL, BL and AP wrote the manuscript. All authors contributed to the final version.
Competing interests. We declare that no competing interest are present. Table B1. Analytical trajectory of a bubble in an idealized sheared wind profile: U (z) = U 0 +Λz with (forecast) or without (analysis) vertical ascent, and resulting average speed.
forecast analysis (atmosphere) vertical velocityŻ = 0Ż = W horizontal velocityẊ = U 0Ẋ = U 0 + ΛZ vertical trajectory Z(t) = Z 0 Z(t) = Z 0 + W t horizontal trajectory X(t) = X 0 + U 0 t X(t) = X 0 + U 0 t + ΛW 2 t 2 averaged vertical velocity over ∆t Warning: The animations play interactively on most operating systems and pdf viewers by single or double click on the active links. In case it does not work (such as in the PDF file produced by arXiv), here are the links to be copied in a navigator -

S1 Tracking of the smoke bubble
The tracking screenshot shows how the tracking of the smoke bubble is interactively performed on a scattering ratio section.

S2 Time evolution of vortex O
The vortex O animation shows the time evolution of vortex O in six panels: a) PV at the isentropic level of the vortex according to PV tracking; b) ozone mixing ratio at the isentropic level of the vortex according to ozone tracking; c) PV altitude-latitude section at the longitude of the vortex according to PV tracking; d) ozone mixing ratio altitude-latitude section at the longitude of the vortex according to ozone tracking; e) PV altitude-latitude section at the longitude of the vortex according to PV tracking; f) ozone mixing ratio altitude-latitude section at the longitude of the vortex according to ozone tracking. In all panels the white cross indicates the location of the centroid vortex according to the relevant tracking. In panels c-f, the contours show the potential temperature.
Using a standard reader, the animation can be played in a continuous way or slide by slide to follow the evolution.

S3 Time evolution of vortices A, B1 and B2
The vortex A animation, vortex B1 animation and vortex B2 animation show, respectively, the time evolution of vortices A, B1 and B2 in two superimposed panels : a) PV at the isentropic level of the vortex according to PV tracking; b) ozone mixing ratio at the isentropic level of the vortex according to ozone tracking. The black cross indicates the location of the centroid vortex according to the relevant tracking S4 Trajectories of vortices O, B1 and B2 in 2017 S5 Potential vorticity at the core of the vortex Figure S4 show the temporal series of the analysed potential vorticity (PV) and the PV anomaly at the core of the main 2020 vortex. The PV anomaly is defined at each latitude and altitude as the deviation with respect to the zonal average at the same latitude and altitude, the altitude being the hybrid vertical coordinate. From the birth of the vortex until the beginning of March, the PV is close to zero and even slightly positive, meaning that the vortex is on the verge of inertial instability. The PV anomaly grows as the low PV bubble rises within an ambient air of increasing PV with height. This sequence ends with the breaking of the vortex by the end of February  and is followed by another sequence where the value of analysed PV suggests an im-