A semi-empirical model for mesospheric and stratospheric NO y produced by energetic particle precipitation

The MIPAS Fourier transform spectrometer on board Envisat has measured global distributions of the six principal reactive nitrogen (NOy) compounds (HNO3, NO2, NO, N2O5, ClONO2, and HNO4) during 2002–2012. These observations were used in a previous study to detect regular polar winter descent of reactive nitrogen produced by energetic electron precipitation (EPP) down to the lower stratosphere, often called EPP indirect effect. It has further been shown that the observed fraction of NOy produced by EPP (EPP-NOy) has a nearly linear relationship with the geomagnetic Ap index when taking 5 into account the time lag introduced by transport. Here we exploit these results in a semi-empirical model for computation of EPP-modulated NOy densities and wintertime downward fluxes through stratospheric and mesospheric pressure levels. Since the Ap dependence of EPP-NOy is distorted during episodes of strong descent in Artic winters associated with elevated stratopause events, a specific parameterization has been developed for these episodes. This model accurately reproduces the observations from MIPAS and is also consistent with estimates from other satellite instruments. Since stratospheric EPP-NOy 10 depositions lead to changes in stratospheric ozone with possible implications for climate, the model presented here can be utilized in climate simulations without the need to incorporate many thermospheric and upper mesospheric processes. By employing historical geomagnetic indices, the model also allows for reconstruction of the EPP indirect effect since 1850. We found secular variations of solar cycle-averaged stratospheric EPP-NOy depositions in the order of 1 Gigamole. In particular, we model a reduction of the EPP-NOy deposition rate during the last three decades, related to the coincident decline of geo15 magnetic activity, that corresponds to 1.8% of the NOy production rate by N2O oxidation. As the decline of the geomagnetic activity level is expected to continue in the coming decades, this is likely to affect the longterm NOy trend by counteracting the expected increase caused by growing N2O emissions.

Abstract.The MIPAS Fourier transform spectrometer on board Envisat has measured global distributions of the six principal reactive nitrogen (NO y ) compounds (HNO 3 , NO 2 , NO, N 2 O 5 , ClONO 2 , and HNO 4 ) during 2002-2012.These observations were used previously to detect regular polar winter descent of reactive nitrogen produced by energetic particle precipitation (EPP) down to the lower stratosphere, often called the EPP indirect effect.It has further been shown that the observed fraction of NO y produced by EPP (EPP-NO y ) has a nearly linear relationship with the geomagnetic A p index when taking into account the time lag introduced by transport.Here we exploit these results in a semiempirical model for computation of EPP-modulated NO y densities and wintertime downward fluxes through stratospheric and mesospheric pressure levels.Since the A p dependence of EPP-NO y is distorted during episodes of strong descent in Arctic winters associated with elevated stratopause events, a specific parameterization has been developed for these episodes.This model accurately reproduces the observations from MIPAS and is also consistent with estimates from other satellite instruments.Since stratospheric EPP-NO y depositions lead to changes in stratospheric ozone with possible implications for climate, the model presented here can be utilized in climate simulations without the need to incorporate many thermospheric and upper mesospheric processes.By employing historical geomagnetic indices, the model also allows for reconstruction of the EPP indirect effect since 1850.We found secular variations of solar cycleaveraged stratospheric EPP-NO y depositions on the order of 1 GM.In particular, we model a reduction of the EPP-NO y deposition rate during the last 3 decades, related to the coincident decline of geomagnetic activity that corresponds to 1.8 % of the NO y production rate by N 2 O oxidation.As the decline of the geomagnetic activity level is expected to continue in the coming decades, this is likely to affect the long-term NO y trend by counteracting the expected increase caused by growing N 2 O emissions.

Introduction
Both solar protons and energetic magnetospheric electrons affect the chemistry in the stratosphere and mesosphere.These energetic particles can alter atmospheric composition either via in situ production of reactive nitrogen and hydrogen species, or by subsidence of air rich in odd nitrogen from its source region, the upper mesosphere and lower thermosphere.In the stratosphere these reactive species gain importance by participating in the catalytic ozone destruction cycles.The in situ production is also called the "direct" effect of energetic particle precipitation (EPP).It is particularly important in the context of solar protons, because only these have been shown to penetrate deep enough into the stratosphere (Jackman et al., 2008, and references therein).
Subsidence of air rich in reactive nitrogen from above is called the EPP indirect effect (EPP IE) (Randall et al., 2007); it is particularly important in the context of auroral electrons.The indirect effect is limited to polar night regions, first, because this strong subsidence can only take place in the downwelling branch of the overturning circulation, i.e., in the po-B.Funke et al.: EPP-NO y model lar winter mesosphere, and second, because odd nitrogen can survive its transport through the mesosphere only in the absence of sunlight.
Based on measurements with the Michelson Interferometer for Passive Atmospheric Sounding (MIPAS), Funke et al. (2005a) estimated the total amount of NO x ([NO x ] = [NO] + [NO 2 ]) descending from the EPP source region into the stratosphere at 2.4 Gigamole (GM) for the southern polar winter 2003.However, since below 50 km NO x is partly converted into its reservoirs (Stiller et al., 2005), the assessment of the indirect EPP effect requires consideration of the entire NO y family ( In a recent paper, Funke et al. (2014a) provide quantitative estimates of the total amount of EPP-NO y for the years 2002-2012, also inferred from MIPAS measurements.In a subsequent paper, Funke et al. (2014b) showed that the EPP IE, i.e., the descended EPP-NO y , is highly correlated with geomagnetic activity, as indicated by the A p index, in Southern Hemisphere (SH) winters and dynamically unperturbed Northern Hemisphere (NH) winters.This suggests that the indirect effect during those winters is driven by the EPP source strength rather than by variations of subsidence.Similar tight correlations with the A p index have been found in seasonally averaged upper stratospheric polar winter NO 2 column density in both hemispheres observed by the Global Ozone Monitoring By Stars (GOMOS) instrument taken during 2002-2006(Seppälä et al., 2007) ) and in estimates of SH EPP-NO x depositions from Halogen Occultation Experiment (HALOE) observations during 1992-2005 (Randall et al., 2007).
However, in NH winters with perturbed dynamics, characterized by episodes of sudden stratospheric warmings (SSW) and associated elevated stratopause (ES) events, accelerated descent in the reformed polar vortex leads to much stronger odd nitrogen descent than in quiescent winters with a similar geomagnetic activity level.Holt et al. (2013) investigated the influence of SSW/ES events on the transport of odd nitrogen produced by EPP from the mesosphere-lower thermosphere to the stratosphere using the Whole Atmosphere Community Climate Model (WACCM).They found that the NO x amount that descends to the stratosphere is strongly affected by the timing of the event, resulting in higher amounts for mid-winter SSW/ES events compared to those occurring in late winter.This behavior was attributed by these authors to the pronounced seasonal dependence of the strength of the vertical winds following an event.
In recent years, the potential impact of particle precipitation on regional climate has been gaining the attention of the climate modeling community.Solar forcing recommendations for the recently launched Climate Model Intercomparison Project Phase 6 (CMIP6) (Eyring et al., 2016) include, for the first time, the consideration of energetic particle effects (Matthes et al., 2016).EPP is strongly linked to solar activity and hence to the solar cycle, either directly by coronal mass ejections producing solar energetic particles or indirectly by the impact of the solar wind on the Earth's magnetosphere.EPP-induced ozone changes are thought to modify the thermal structure and winds in the stratosphere that, in turn, modulate the strength of the polar vortex.The introduced signal could then propagate down to the surface, introducing significant solar-like variations of regional climate (Baumgaertner et al., 2009;Rozanov et al., 2012;Seppälä et al., 2014;Maliniemi et al., 2014).Today, there is an increasing number of chemistry climate models capable of dealing with EPP effects; however, not all of them extend up into the upper mesosphere/lower thermosphere, where a large fraction of EPP-induced odd nitrogen production occurs.Those models with their upper lid in the mesosphere, i.e., which do not represent the entire EPP source region, require an odd nitrogen upper boundary condition, accounting for EPP productions higher up, in order to allow for simulation of the introduced EPP IE in the model domain (Baumgaertner et al., 2009;Rozanov et al., 2012).
In this paper, we provide a detailed semi-empirical model for retrodiction/prediction of the indirect EPP-NO y as a function of the geomagnetic A p index that has been adjusted to the decadal MIPAS EPP-NO y record.In order to account for the pronounced EPP-NO y increases during ES events, a specific parameterization has been included for these episodes.The aim of this model is to provide the stratospheric and lower mesospheric NO y to chemistry climate models that do not explicitly model upper mesospheric and thermospheric EPP effects.A further application of this model is the reconstruction of the EPP indirect effect on secular timescales by employing historical geomagnetic indices.
This paper is organized as follows: Sect. 2 describes the MIPAS EPP-NO y data used for adjusting the semi-empirical model.Section 3 provides a detailed description of the semiempirical model for hemispheric EPP-NO y amounts and fluxes during polar winters, excluding ES episodes in the NH.The extension of the latter model with respect to the consideration of ES events is provided in Sect.4, and the detection of such events is discussed in Sect. 5.The modeled EPP indirect effect is compared to available observational estimates in Sect.6.The application of the model as an odd nitrogen upper boundary condition for chemistry climate models with their upper lid in the mesosphere is discussed in Sect.7, and Sect.8 deals with the reconstruction of the EPP indirect effect during 1850-2015.

Observations
The MIPAS instrument (Fischer et al., 2008) on the polarorbiting Envisat satellite provided global stratospheric and mesospheric measurements of temperature (von Clarmann et al., 2003), NO x (Funke et al., 2005b), NO y (Mengistu Tsidu et al., 2004) and numerous other trace species (e.g., von Clarmann et al., 2009(e.g., von Clarmann et al., ) during 2002(e.g., von Clarmann et al., -2012.From these data, the contribution of NO y produced by EPP has been distinguished from that produced by N 2 O oxidation using a tracer correlation method that is based on coincident CH 4 and CO observations (Funke et al., 2014a).The latter tracer is used to restrict the EPP-NO y detection to observations containing mesospheric air.The EPP-NO y uncertainty is dominated by the multiplicative component of the NO y systematic retrieval error that is about 10 %.The scatter in the tracer correlation results in a precision of inferred EEP-NO y of 0.5 ppbv, which can be considered to be the 1-σ detection limit, particularly at lower altitudes.Other uncertainties act systematically upon the estimated EEP-NO y and lead to a possible underestimation.This is particularly true for the end of the winter and during stratospheric warming episodes.For further details on the error analysis of this method, see Funke et al. (2014a).The altitude resolution is given by that of the MIPAS NO y data used to derive the EPP-NO y record and ranges from 4-6 km in the stratosphere to 6-9 km in the mesosphere.EPP-NO y VMR profiles have been converted to number density profiles using temperature and pressure information inferred from the same MIPAS spectra from which the NO y data were also retrieved.
Here, we use the daily zonal mean climatology of EPP-NO y densities, available on latitudinal bins of 10 • with global coverage.From this EPP-NO y record, we determine the hemispheric EPP-NO y total amounts N t (z, t) by first vertically integrating the NO y densities from z 0 = 40 hPa to pressure level z.The amount in GM within each latitude bin φ is then calculated as the product of the respective zonal mean column density and the area A(φ) covered by the bin.In a second step these individual contributions are summed up for each hemisphere, i.e., where [EPP-NO y ] is the density of the EPP-related NO y contribution.In a similar way, the differential EPP-NO y amount N d (z, t) in units of GM km −1 , i.e., the vertical differential of N t (z, t), is calculated by This quantity is proportional to the hemispherically averaged mean density of EPP-NO y .Finally, we derive the hemispherically integrated EPP-NO y flux F (z, t) through z from where L(z, t) is the hemispheric photochemical loss rate of EPP-NO y below z (in units of GM day −1 ).The latter was obtained from box model calculations that have been constrained by observed fields of temperature, O 3 , and NO x (see Funke et al., 2014b, for more details).Equation ( 3) is only valid if there are no local EPP-NO y productions below z.
For this reason, we exclude episodes of solar proton forcing from the calculated EPP-NO y flux data.In principle, precipitating electrons from the radiation belts, depositing their energy primarily in the middle and upper mesosphere, could also induce local productions in the altitude range of interest, although recent studies (e.g., Sinnhuber et al., 2014) have indicated that their contributions are on average negligibly small in the polar winter upper stratosphere and lower mesosphere when comparing to the EPP indirect effect.It should be noted, however, that for isolated cases (e.g., large geomagnetic storms) the productions from radiation belt electrons can be important even in the lower mesosphere (Andersson et al., 2014;Arsenovic et al., 2016).
3 Semi-empirical model for EPP-NO y in SH and NH winters (excluding ES episodes) In this section, we develop an empirical model of hemispheric EPP-NO y differential amounts and fluxes in SH and NH winters (excluding ES episodes) as a function of the geomagnetic A p index, altitude, and time, based on the EPP-NO y distributions inferred from MIPAS during 2002-2012.Note that the model is developed for hemispherically integrated quantities, where the total EPP-NO y is conserved.These quantities, however, can be converted into zonal mean densities and fluxes by imposing the observed latitudinal distribution of EPP-NO y (see Sect. 7).
Our model is based on the linear dependence of the observed stratospheric and mesospheric EPP-NO y on the A p index as demonstrated in Funke et al. (2014b).They performed a multi-linear regression of monthly averaged EPP-NO y amounts to the average A p indices of the current and 3 preceding months in order to empirically account for time lags introduced by transport and its dispersion.In a more theoretical approach, the EPP-NO y differential amount N d (z, t) at pressure level z and time t can be described by with the Green function G depending on the mean transport time , its dispersion , and the photochemical loss rates l experienced during the descent from the source region, the latter depending on altitude and time.We further assume that the temporal variation of the photochemical modulation during this descent is slow compared to the dispersion of transport times such that where N d (z, t) is the spatio-temporal distribution of the EPP-NO y amounts for a constant A p index of unity, and A p (z, t) It describes the propagation of the A p modulation from the source region down to the stratosphere.This Green function has the same mathematical structure as that describing the transport of a passive tracer (Andrews et al., 1999) and can be approximated by an inverse Gaussian function, i.e., with the mean transport time (z) from the source region to the pressure level z and the width of the distribution (z).
For N d (z, t) we use the following empirical function: where t is the number of days passed since 1 July for the NH and since 1 January for the SH.N m (z) is the maximum EPP-NO y differential amount encountered at the pressure level z during the course of the winter, t N m (z) the occurrence time of this maximum, and w N (z) a parameter determining the temporal width of the distribution.We have chosen this function among several candidate analytical functions because it describes the observed distribution very closely and it allows us to express the temporal evolution by a few, physically meaningful parameters.
The parameters , , N m , t N m , and w N m have been adjusted by performing a nonlinear least squares fit of Eq. (5) to the observed daily vertical distributions of EPP-NO y amounts for each pressure level z, excluding periods of SPE events (orange-shaded areas in Fig. 3a and b).Also, episodes of SSW/ES events in NH winters (grey-shaded areas in Fig. 3a) are excluded.These events require a specific parameterization that is discussed in Sect. 4.
Since, due to the complex temporal structure of the A p evolution, multiple maxima of the objective function of this optimization problem are expected for the fit parameters and , we use a quasi-global rather than local minimization strategy.This is, we scan, within reasonable bounds, thespace and adjust for each pair of and the corresponding quantities N m , t N m , and w N m .An additional constraint for has been introduced by assuming 2 ∼ , which would be the case for a linear relationship between vertical advection and diffusion.Waugh and Hall (2002) have reported a dependence of 2 = 0.7 for stratospheric transport.We obtain the smallest χ 2 values for (z) = 0.35 (z) + 4.24. (10) A smaller scaling factor of 0.35, compared to the value of 0.7 derived by Waugh and Hall (2002) for tracer transport from the tropical transition layer into the stratosphere, is reasonable since less eddy diffusion is expected for vertical transport within the polar vortex.The empirically determined "offset" of 4.24 days indicates that dispersion is more pronounced for short transport times (i.e., in the mesosphere) compared to longer transport times (i.e., in the stratosphere).
As a consequence, the maximum of the distribution function G(t , , ) is shifted to shorter transport times compared to its mean value in the mesosphere.Physically, this could be caused by the predominantly diffusive entry of EPP-NO y from the auroral source region in the upper mesosphere and lower thermosphere, where the mean circulation is upward (Smith et al., 2011).On the other hand, local odd nitrogen productions by radiation belt electron precipitation in the mesosphere would cause a similarly dispersed spectrum of transport times.Table 1 lists the derived parameters N m (z), t N m (z), and w N m (z) for pressure levels between 30 and 0.01 hPa for both hemispheres.The best fitting values for (z) are shown in Fig. 1 (diamonds).Although transport times may vary over the winter season in dependence on the strength of the vertical winds, our adjusted values of are time-independent and thus represent seasonal averages, implicitly weighted with the actual EPP-NO y amount by the fitting algorithm.Therefore, the obtained values are most representative of the period of the EPP-NO y maximum occurrence.Note that the fit of (z) becomes unstable below 0.5 hPa for NH winters due to small signal-to-noise ratios, caused by the low EPP-NO y amounts together with the large dynamical variability (not shown).In the SH, increases steadily towards lower pressure levels (higher pressures), as expected.The fitted mean transport times are in very good agreement with those derived from the SH mid-winter descent velocities estimated in Funke et al. (2014b) (indicated by dashed lines).
The increase in towards lower pressure levels (higher pressures) is even more pronounced in the NH, above 0.3 hPa, where the fitted values exceed the transport times derived from the mid-winter descent velocities.Longer mesospheric transport times at the time of the EPP-NO y maximum occurrence are expected in the NH due to the deceleration of mesospheric descent in the second half of the winter (Funke et al., 2014b).Below, the fitted NH mean transport times become shorter again (and closer to those derived from the midwinter descent velocities) since the EPP-NO y reaching those pressure levels has been transported through the mesosphere primarily during the first half of the winter.
Mean transport times from the EPP source region to the indicated pressure level in the SH (blue diamonds) and NH (red diamonds, only above 0.5 hPa) as derived from the best-fit mean value (z) of the inverse Gaussian used as an A p weighting function.Transport times derived from the mid-winter vertical velocities of Funke et al. (2014b) are shown by dashed lines.Solid lines correspond to the transport times expressed as a function of t N m (z), used in the semi-empirical model.
Mean transport times and occurrence times of the observed EPP-NO y maximum t N m (z) are closely linked in the SH; however, the latter is shorter than the former.Such a time lag is likely related to the seasonal dependence of meso-spheric downward velocities (being larger around solstice), introducing a distortion of the temporal evolution of the differential EPP-NO y amounts.(z) can be reasonably well expressed by as indicated by the solid blue line in Fig. 1.
In the NH, the corresponding parameterization in terms of t N,NH m (z) reproduces (z) above 0.3 hPa.Below, however, the resulting (z) would be significantly underestimated.This is expected because only EPP-NO y descending during the first part of the winter reaches the stratosphere due to the deceleration of mesospheric descent around mid-winter.As a consequence, the stratospheric NH EPP-NO y maximum occurs much earlier than in the SH despite the longer transport times as derived from the NH mid-winter descent velocities, the latter providing an estimate of in the vertical range where no fitted values are available.On the other hand, NH mid-winter descent velocities below 0.3 hPa can be expressed reasonably well as a function of t N,SH m (z), and we obtain   Figure 2 shows the seasonal evolution of EPP-NO y differential amounts for a constant A p index of 10 corresponding to the average A p during the Envisat mission lifetime in both hemispheres, i.e., 10 × N d (z, t).As expected, stratospheric NH differential amounts are considerably smaller than those in the SH, the latter exceeding the former by a factor of 8 around 10 hPa.Both distributions reflect the decrease in descent rates from the mesosphere to the stratosphere, leading to a change in the vertical gradient of the EPP-NO y tongue with a "knee" at around 1 hPa.The deceleration of vertical transport below the stratopause is also responsible for the increased amounts, there, due to compression of EPP-NO y .The latter occurs because the temporal increase in NO y is proportional to the vertical gradient of the descent rate, which follows from mass conservation.
In contrast to the SH, minimum differential amounts are found around 0.2 hPa in Arctic winters.The minimum is caused by the deceleration with time of the vertical velocity occurring in the mid-winter at pressure levels above 0.2 hPa and the acceleration below that pressure level, causing a local depletion of EPP-NO y .This sudden deceleration of mesospheric descent is also responsible for the "splitting" of the EPP-NO y tongue into a slowly descending mesospheric branch, reaching the 0.2 hPa level around April, and a rapidly descending stratospheric branch, comparable to the typical SH pattern.
Figure 3a and b show the observed and modeled temporal evolutions of EPP-NO y differential amounts at pressure levels of 0.03, 0.3, 2, and 10 hPa in both hemispheres.There is generally good agreement, indicating that most of the interannual variability encountered in the observed EPP-NO y amounts can be reproduced by the semi-empirical model, particularly in the SH.As expected, the agreement in NH winters is not as good due to the more pronounced dynamical modulation.This is particularly the case for winters with SSW and ES events, which are not accounted for in the model for quiescent NH winters.During and after these events, EPP-NO y amounts are underestimated by the model by up to an order of magnitude, highlighting the need for specific parameterizations as presented in Sect. 4. Also, the typically lower EPP-NO y amounts in the NH, which in some winters are close to the detection limit, are necessarily more dispersed relative to the model results.However, the ability to reproduce singular features, such as the "peaky" evolution of N d (z, t) during January 2007 in the NH related to a shortterm increase in geomagnetic activity, provides confidence in the model.
The semi-empirical model for EPP-NO y fluxes through given pressure levels has been constructed in a similar way to the model for the EPP-NO y amounts.Since the flux F (z, t) can be expressed as the product of the NO y differential amount N d (z, t) at z and the EPP-NO y descent rate (Funke et al., 2014b), the latter being independent of A p , we can assume the same A p dependence A p (z, t) as for the EPP-NO y differential amounts, i.e., F (z, t) = F (z, t)A p (z, t). (13) For F (z, t) we use the same type of function as for A d (z, t) in Eq. ( 9), i.e., This empirical function is then adjusted to the "observed" fluxes F obs (z, t) through the vertical level z, divided by the A p propagation function at z, i.e., As discussed in Sect.2, F obs (z, t) is derived from the temporal changes in the sum of the observed EPP-NO y total amounts N t (z, t) and the accumulated photochemical losses L(z, t).We assume that any reduction of the total amount is caused by vertical mixing due to vortex rupture and that under these conditions vertical velocities tend to be zero, such that we can limit the "observed" fluxes to non-negative values (see Funke et al., 2014b).Again, we exclude periods of SPE and SSW/ES events (shaded areas in Fig. 3a).Table 2 lists the derived parameters F m (z), t F m (z), and w F m (z) for pressure levels between 30 and 0.02 hPa for both hemispheres.It also provides the modeled seasonal EPP-NO y deposition, T 10 , below z corresponding to the 2002-2012 average geomagnetic forcing (A p = 10).Below 0.02 hPa, the modeled depositions are 0.48 GM in the NH and 1.21 GM in the SH, that is, nearly 3 times more in the former than in the latter.
Figure 4 shows the seasonal evolution of the modeled EPP-NO y fluxes in both hemispheres, again for A p = 10.Maximum fluxes of 0.07 GM day −1 in the NH and 0.22 GM day −1 in the SH are found at the uppermost pressure levels during the winter solstice.Towards lower altitudes, both NH and SH fluxes are decreasing.In the mesosphere, this decrease is mainly related to photochemical losses.At lower altitudes, dynamical loss due to mixing out of the polar vortex is responsible for the flux decrease.Assuming that the fraction of EPP-NO y mixed out of the vortex is not being transported further downwards, the vertical gradient of the seasonally integrated fluxes (T 10 of Table 2) hence represents the deposition profile of EPP-NO y at the end of the winter.

Parameterization for elevated stratopause events
The challenge of parameterizing EPP-NO y amounts and fluxes during ES events resides mainly in the scarcity of observational data during these events.MIPAS has recorded NO y data with sufficient temporal coverage only during two events occurring in January 2004 and February 2009.Holt et al. (2013) have shown, using WACCM simulations with constant geomagnetic forcing, that, besides the geomagnetic activity level, the event timing is a crucial driver of the strength of odd nitrogen descent because of the seasonal dependence of residual vertical wind speeds.In addition, it is also likely that the EPP-NO y amount in the source region is modulated by photochemical losses, again resulting in  smaller EPP-NO y depositions during events occurring later in the winter.The two events observed by MIPAS have rather different characteristics regarding the geomagnetic activity level and timing and hence cover a large range of the expected variability.
Our approach to provide a general parameterization of odd nitrogen descent during ES events is, first, to parameterize the EPP-NO y amounts and fluxes individually for each of the two observed events, and then to exploit dependencies of the obtained parameters on the event timing.The time evolution of the differential EPP-NO y amounts N d (z, t) during ES winters is finally calculated by adding the EPP-NO y residual amounts during the ES event to the "quiescent" differential amounts.
The modeling of the individual 2004 and 2009 ES events is performed in a very similar way as for the quiescent winters (see Eqs. 5 and 13 of Sect.3).First, we adjust the parameters of the A p propagation function and the spatio-temporal term to the EPP-NO y residual amounts, that is, the difference between the observed differential amounts and those modeled for quiescent conditions, after the onset of the event at the pressure levels reached by the descending NO y tongue.The onset time t ES 0 is defined as the time (days since 1 July) when the ES-related differential amount increases at 0.02 hPa, resulting in t 2004 0 = 196 and t 2009 0 = 221.We further define t ES m (z) as the time lag between t ES 0 and the occurrence of the EPP-NO y maximum at pressure level z, and we assume that In order to account for the fact that there are no residual amounts before the event, we apply the following correction to the spatio-temporal term, resulting in a stretching of the temporal distribution before the occurrence of the EPP-NO y maximum, while leaving it unchanged afterwards.The adjusted values of t ES m (z) do not differ significantly between the two events in most of the vertical range where the EPP-NO y maximum occurs.In the vicinity of the lowermost pressure level reached around equinox, however, t ES m (z) increases drastically due to the deceleration of descent.We parameterize therefore t ES m (z) as a function of t ES 0 and z, with t ES m (z) being a fourth-order polynomial of ln(z) (see Table 3 for coefficients).The parameter w ES (z), describing the width of the EPP-NO y peak after the ES event, is in first order height-independent, and has been adjusted to the common value of 0.15 for both events, which corresponds to a full width at half maximum of 13 days.Finally, the adjusted values of the A p -normalized maximum amounts N ES m (z) for both events are shown in Fig. 5a (blue and red squares for 2004 and 2009, respectively).The observed maximum amounts peak around the stratopause, being a factor of 3 higher in 2004 compared to 2009.This difference becomes smaller with height and decreases at the uppermost levels to a factor of 2.
Similarly, we fit the spatio-temporal term in Eq. ( 14) to the observed EPP-NO y fluxes during both events.As a good approximation, we use the same width parameter w ES (z) = 0.15 and the same lag time t ES m (z) as for the residual amounts.The adjusted values of the maximum A p -normalized fluxes F ES m (z) for both events are shown as blue and red squares in Fig. 5b.For both events, the maximum fluxes decrease monotonically towards lower pressure levels, indicating photochemical and/or dynamical losses.The ratio of the 2004 and 2009 fluxes is nearly constant, with a value of 4.2 in the mesosphere.Since the downwelling flux at 0.01 hPa is in first order the product of the EPP-NO y amount and vertical velocity during the ES event in the source region, we expect this ratio to be influenced by the seasonal variations of both quantities.Around the stratopause and below, fluxes decrease faster in 2009 because the maximum occurrence time is closer to equinox at these pressure levels.The observed A p -normalized maximum flux can be parameterized as a function of t ES 0 and z by with F ES m (z) being a fourth-order polynomial of ln(z) (see Table 3 for coefficients).takes the values of 0.00567 and 0.00125 for the 2004 and 2009 events, respectively.
We calculate descent velocities ω ES m (z), dividing the adjusted values of The obtained values are shown with blue and red squares for 2004 and 2009, respectively, in Fig. 5c.Mesospheric descent rates are a factor of 2 higher in 2009 (2 km per day) compared to 2009 (1 km per day), in qualitative agreement with the model results of Holt et al. (2013) for events with similar timing.The adjusted descent rates for the 2004 event are also in good agreement with previous estimates of the vertical component of the meridional circulations using MIPAS temperatures and diabatic heating rates (Randall et al., 2015).During both events, the descent rates decrease rapidly towards the stratopause, where they take values of around 200 m day −1 .The velocity ratio between both events is rather constant down to the stratopause.We thus use a similar parameterization as given by Eq. ( 18) for the vertical velocities, with ω ES m (z) being a sixth-order polynomial of ln(z) (see Table 3   in order to obtain an event time-dependent quantity related to the EPP-NO y amount in the source region.Both (t ES 0 ) and (t ES 0 ) are expected to maximize at solstice (day 173) and to reach the zero level close to equinox.We use therefore a similar expression as provided in Eq. ( 9) to parameterize their dependence on the event timing: Figure 6 shows the resulting values of F ES m at 0.4 hPa and ω ES at 0.08 hPa as a function of t ES 0 .Both time dependencies can be qualitatively compared to the total EPP-NO x amounts crossing the 0.41 hPa level (i.e., the integral flux through this level) and maximum descent at the 0.08 hPa level for a large number of ES events simulated by WACCM as presented by Holt et al. (2013) (their Figs.6a and 9c, respectively).In the parameterization, as well as in the WACCM simulations, descent rates decay almost linearly with t ES 0 , reaching minimum values around the end of March (note that the day of the event is defined in Holt et al. (2013) as the central date of the preceding SSW, typically 8 days before t ES 0 ).ES events starting around 1 February are characterized by half of the maximum descent rate for solstice ES events.Both parameterized and WACCM-simulated fluxes decrease with time more nonlinearly and reach the background level in mid-February (50 % of the solstice flux around mid-January).In our semiempirical model, the more pronounced nonlinearity of the flux decay compared to descent is introduced by Eqs. ( 20) and ( 21), in consonance with our hypothesis that the EPP-NO y flux depends on the temporal evolutions of both the descent rates and the EPP-NO y amounts in the source region.
As mentioned above, the time evolution of the differential EPP-NO y amounts N d (z, t) during ES winters is finally calculated by adding the EPP-NO y residual amounts during the ES event (as a function of ES (z), N ES m (z), t ES m (z), and w ES m (z)) to the "quiescent" differential amounts (Eq.5).In order to illustrate the impact of the ES event timing on the differential amounts, Fig. 7 shows the seasonal evolution of the resulting differential amounts for a constant geomagnetic level A p = 10 in NH winters with ES events occurring on 15 December, 15 January, 10 February and 5 March.A maximum differential amount of 0.14 GM km −1 is predicted for a December event, exceeding the maximum amounts in SH winters (encountered around 10 hPa) by a factor of nearly 2. On the other hand, the EPP-NO y amounts after March events are hardly distinguishable from the background.Also, the lowermost pressure level reached by the EPP-NO y tongue differs significantly between the events, being 7, 3, 1, and 0.2 hPa, respectively, for these events.
Figure 8 compares the observed and modeled temporal evolutions of EPP-NO y differential amounts at pressure levels 0.03, 0.1, 0.3, and 1 hPa during the ES winters 2004 and 2009.The generally good agreement demonstrates the capability of the extended semi-empirical model to reproduce the observed EPP-NO y also in NH winters with ES events, in contrast to the model for quiescent winters (compare to Fig. 3a).

Determination of ES onsets relevant for EPP-NO y
While the detection of ES events and the determination of their onset date, t ES 0 , is straightforward for the winters observed by MIPAS, a specific criterion is required in general in order to model the EPP-NO y distribution in longer time periods.Further, if the semi-empirical model is used in chemistry climate models to provide an odd nitrogen upper boundary condition (see Sect. 7), the detection of ES events needs to be performed online.An obvious quantity for the detection would be the stratopause height derived from zonally averaged polar temperatures as suggested by de la Torre et al. ( 2012) and Chandran et al. (2013).However, the un-equivocal detection of ES events from the polar zonal mean temperature profile suffers from short-term excursions of the stratopause height that result from transient wave forcing.The temporal smoothing, needed to reduce this effect, disables the online detection on a daily basis, as required to account for the fast increase in EEP-NO y in the mesosphere at the beginning of an event.Also, if the semi-empirical model is used for reconstruction of the EPP indirect effect over historical time periods, and reanalyzed meteorological data need to be employed for ES detection, the latter would suffer from the poor representation of mesospheric temperatures in the reanalysis data.
An alternative approach to detection of ES events and determination of t ES 0 from reanalysis or model temperature fields takes advantage of the tight anti-correlation of the upper stratospheric and mesospheric meridional temperature gradients during ES events.In particular, the difference between the zonal mean temperature averaged over 0-30 • N and that averaged over 70-90 • N at 1 hPa, in the following referred to as T 30-70 , from MIPAS observations during 2002-2012 shows pronounced increases of up to 55 K during the 2004 and 2009 events, not reached during quiescent winters.We have tested this criterion using 1 hPa temperatures from a transient model simulation with the ECHAM/MESSy Atmospheric Chemistry (EMAC) model, covering the period 1979-2014.The EMAC model is a numerical chemistry and climate simulation system that includes sub-models describing tropospheric and middle atmosphere processes and their interaction with oceans, land and human influences (Jöckel et al., 2010).It uses the second version of the Modular Earth Submodel System (MESSy2) to link multi-institutional computer codes.The core atmospheric model is the 5th gener- ation European Centre Hamburg general circulation model (Roeckner et al., 2006).For the present study we applied EMAC (ECHAM5 version 5.3.02,MESSy version 2.50) at the T42L90MA resolution, i.e., with a spherical truncation of T42 (corresponding to a quadratic Gaussian grid of approx.2.8 by 2.8 • in latitude and longitude) with 90 vertical hybrid pressure levels of up to 0.01 hPa.Vorticity, divergence, and temperature fields have been relaxed to ERA-Interim reanalysis data (Dee et al., 2011) below 1 hPa, facilitating the simulation of dynamic events that have occurred during 1979-2014.
Figure 9 shows the temporal evolution of T 30-70 over the whole simulation period.The 53 K level (indicated by the dotted red line) is exceeded during all reported ES events in the present and last decades, namely 2004, 2006, 2009, and 2013 (indicated by vertical dashed lines).Elevated stratopause events are also detected in 1985 and 1987.Table 4 lists the dates of the first day exceeding the 53 K threshold during all events.These dates coincide very precisely (within 1 day) with t ES 0 , as determined from the MIPAS data for the 2004 and 2009 events.Also, the onset dates for the 2006 and 2013 events are consistent with the analysis per- formed in previous works (Manney et al., 2008;Pérot et al., 2014).An inspection of the modeled stratopause evolution in 1985 and 1987 confirms elevated stratopauses after SSWs in these winters (not shown).Also, winters with T 30-70 below the 53 K threshold show no elevated stratopause, despite the occurrence of several SSWs.Since EPP was not considered in the EMAC simulation, we look at the simulated CO evolution in the lower mesosphere in order to prove whether enhanced descent indeed occurred during the detected ES events.CO is an adequate tracer of mesospheric air due to increasing concentrations towards the upper mesosphere/lower thermosphere and relatively long photochemical lifetimes in polar winter.The solid blue line in Fig. 9 represents the 70-90 • N CO anomalies (with respect to the climatological seasonal mean) at 1 hPa.Noticeable increases in CO are found after all ES events, although the magnitude of the increase after the 1987 event is rather small, most likely related to the late onset.Chandran et al. (2013) investigated the occurrence of ES events, as detected from polar stratopause jumps, in MERRA reanalysis data (Rienecker et al., 2011(Rienecker et al., ) covering 1979(Rienecker et al., -2011. .They identified the events detected by our T 30-70 criterion (see Table 4), but also additional events in winters 1980/81, 1983/84, 1989/90, 1994/95, and 2009/10.In these additional ES events, the maximum values of T 30-70 of the nudged EMAC simulation remained well below 53 K. Further, no significant CO increases at 0.5 hPa were simulated with the EMAC model, except for 1983/84.In this particular winter, however, the CO enhancements have occurred already before the event.The EPP-NO y evolution in the 2009/10 NH winter, which has been observed by MIPAS, does not show indications for ES-related odd nitrogen intrusions.Most of the additional ES events detected by Chandran et al. (2013) were accompanied by minor stratospheric warmings, in contrast to the events detected by the T 30-70 criterion that were preceded by major SSWs.
We thus conclude that our criterion based on T 30-70 allows us to detect the ES events with strong descent of mesospheric air and associated efficient deposition of EPP-NO y in the stratosphere.Also, we found that the first crossing time of the T 30-70 =53 K threshold provides a reasonable estimate of the onset time, t ES 0 .

EPP indirect effect during 1978-2014 and comparison with previous estimates
Figure 10a and b show the semi-empirical model estimates of the EPP-NO y depositions in the SH and NH winters during 1978-2014 together with previous estimations.First, we observe a generally good agreement between the results of the semi-empirical model and the estimates of the EPP indirect effect provided by Funke et al. (2014b).This is not surprising, since both are based on the same MIPAS observations, but this comparison gives us a good measure of the quality of the fitting of the model to the actual measurements from which it has been derived.Figure 10a also shows the estimates on the EPP indirect effect of Randall et al. (2007) for the SH winters 1992-2005 from HALOE NO x solar occultation observations in the upper stratosphere (note that NO x is nearly equivalent to NO y at these altitudes and hence comparable to our results).Like Funke et al. (2014b), they also used a tracer correlation method to extract the EPP-NO x contribution, but, in contrast to the MIPAS-derived depositions, they derived it from the accumulated NO x flux through the 45 km altitude level (∼1 hPa; see Fig. 10a).The flux was calculated from the observed NO x density at that level assuming a constant SH polar winter descent rate of 400 m day −1 .Due to the sparse sampling of HALOE they made important assumptions about the latitudinal distribution of the EPP-NO x inside the vortex www.atmos-chem-phys.net/16/8667/2016/Atmos.Chem.Phys., 16, 8667-8693, 2016 that led to uncertainties as large as 100 % in their estimates (see the grey-shaded area in Fig. 10a).Further, a rather conservative threshold was used for discriminating the EPP-NO x from the background NO x , which might have offset their resulting estimates.We are interested in evaluating the consistency of the MIPAS and HALOE estimates, particularly in terms of inter-annual variability.For that purpose, and in order to account for possible biases related to the different measurements and estimation methods, we adjusted an offset and a scale factor to the HALOE estimates.The determined scale factor of 0.7 is well within the range encompassed by the "average" and "maximum excess NO x " estimates of Randall et al. (2007).The offset of 0.5 GM is rather high but plausible, since comparable EPP-NO y depositions are expected during SH winters with similarly low geomagnetic activity levels such as 1996 (HALOE-derived estimate of 0.1) and 2007 (MIPAS-derived estimate of 0.6).Note also that a negative bias of 0.5 GM in the EPP-NO y depositions would be introduced by an underestimation of about 2 × 10 8 cm −3 in the EPP contribution to the NO x densities at 1 hPa, which is comparable to the conservative threshold for EPP-NO x discrimination from background values used by Randall et al. (2007) (see their Fig. 6).The inter-annual variations of our modeled EPP-NO y depositions below 1 hPa are highly consistent with the adjusted estimates from HALOE in the 1992-2005 period.Particularly during 1993-1998, the agreement is excellent.Holt et al. (2012) provided observational EPP-NO x deposition estimates for NH winters by employing the same method as Randall et al. (2007) but using MERRA-derived vertical velocities instead of a fixed 400 m day −1 descent rate in the flux calculation.They used NO x observations from the Limb Infrared Monitor of the Stratosphere (LIMS) for the winter 1978/79, MIPAS for the winters 2002/03 and 2003/04, as well as Atmospheric Chemistry Experiment Fourier transform spectrometer (ACE-FTS) data for the Arctic winters in 2004-2009.They reported depositions below the 2000 and 3000 K potential temperature surfaces, corresponding roughly to the 0.7 and 0.1 hPa pressure levels, respectively.Again, we adjust these estimates to those derived by Funke et al. (2014b) using a time-independent offset and scale.
For both vertical levels, we determine an offset of 0.2 GM, being considerably lower than for the SH estimates of Randall et al. (2007).This might be partly related to the latitude coverage of the employed instruments (global sampling in the case of LIMS and MIPAS, and 60-85 • N for ACE-FTS in the NH), resulting in a better polar coverage compared to HALOE (< 55 • S during May-August).These sampling differences might also explain why no scaling (derived scale factor of 0.995) needs to be applied in order to fit the 3000 K deposition estimates to those of Funke et al. (2014b) for the 0.1 hPa level.The modeled seasonally integrated fluxes through this pressure level during the NH winters 2004/05, 2005/06 and 2006/07, not available from MI-PAS data, show very good agreement with the estimates from ACE-FTS shifted by 0.2 GM.
However, the EPP-NO x depositions below 2000 K need to be scaled by a factor of 0.77 in order to achieve consistency with those of Funke et al. (2014b) for the 0.7 hPa level.A possible explanation for this mismatch is the use of MERRAderived vertical velocities in Holt et al. (2012) to convert the EPP-NO x densities in fluxes.These velocities might be overestimated by up to 40 % at this pressure level (Funke et al., 2014b).However, differences might also be introduced by comparing depositions below pressure levels and depositions below potential temperature surfaces, because they are characterized by different latitudinal and temporal variations and the EPP-NO y fluxes have strong gradients, particularly in this vertical region.After applying the adjustment to the estimates of Holt et al. (2012), the agreement of the inter-annual variations with those of the semi-empirical model is reasonably good.Bailey et al. (2014) employed the same method as Holt et al. (2012) in SOFIE NO observation.In this case, we apply an offset of 0.12 GM and a scale factor of 0.77 to compare their estimates of the seasonally integrated EPP-NO y fluxes through 0.7 hPa in the winters 2008/09 and 20012/13 to the semi-empirical model, and again find reasonable agreement.
The modeled SH seasonal EPP-NO y depositions during the 1978-2014 period, covering three solar cycles, are on average 1.26 GM below the 0.1 hPa level and 0.99 GM below 1 hPa.The large EPP indirect effect in 2003 -the strongest during the MIPAS observation period -is only exceeded by that of the Antarctic winter 1991 with about 30 % higher depositions.The average NH depositions in 1978-2014 are 0.50 (0.25) GM below the 0.1 (1) hPa level.The EPP-NO y deposition of the extraordinary ES winter 2003/04 is at 3 GM below 0.1 hPa the strongest of the whole period, followed by the 1984/85 ES winter with 1.9 GM.The average contribution of the ES events to the EPP-NO y depositions in the 1978-2014 period is only 4 % (0.02 GM at 0.1 hPa and 0.01 GM at 1 hPa).This indicates that strong descent episodes related to ES events, while being of high relevance for the EPP-NO y evolution during individual ES winters, seem to play only a minor role on longer timescales.However, the average EPP-NO y contributions due to ES events increase noticeably when considering only the last decade.The question of whether the clustering of ES events during the latter period is part of the natural variability or indicative of a tendency, however, still remains open.

EPP-NO y upper boundary conditions for atmospheric models
A major purpose of this semi-empirical model is to provide an upper boundary condition (UBC) for chemistry climate models with an upper lid in the mesosphere.These models leave a large fraction of the EPP source region (extending to the lower thermosphere) uncovered and hence do not allow for a detailed simulation of the EPP indirect effect.However, EPP can still be taken properly into account by prescribing NO y at the upper model lid.This can either be done by specifying a flux of NO x into the top of the model domain (e.g., Baumgaertner et al., 2009) or by specifying a NO x concentration at the uppermost model layer(s) (e.g., Reddmann et al., 2010).Taking into account that NO y NO x at pressure levels higher than approximately 1 hPa, our semiempirical model allows for both types of NO x UBCs, although the prescription of fluxes should be restricted to model levels at 0.02 hPa and lower altitudes in order to minimize contaminations by local productions related to radiation belt electrons (see above).
Typically, the NO y (or NO x ) flux or concentration is assumed to have a zonally homogeneous distribution.Baumgaertner et al. (2009) also assumed a homogeneously distributed flux within 55-90 • latitude, roughly coinciding with the mesospheric polar vortex.We have analyzed the latitudinal distribution of the MIPAS-derived EPP-NO y averaged over the 2002-2012 period in order to come up with a more realistic distribution.The derived meridional dependency is then used to distribute the differential amounts and fluxes in the respective hemisphere.Figure 11 shows the fraction of the hemispheric differential amount polewards of a given latitude separately for SH winters, quiescent NH winters, and ES episodes as a function of pressure.SH and quiescent NH distributions are very similar, with around 60 % of the EPP-NO y at latitudes > 65 • and 90 % at > 50 • .The distributions tend to widen below 0.1 hPa by about 5 • .No pronounced variation of the distributions along the winter season has been encountered.During ES events, the distribution is more confined over the pole, with 60 % of the EPP-NO y at latitudes > 75 • and 90 % at > 60 • .No significant variation with height is found during ES episodes.Normalized latitudinal distributions (φ, z) of EPP-NO y concentrations in the vertical range 1-0.01 hPa are provided in Tables A1, A2, and A3 for SH winters, NH winters, and ES episodes, respectively.Following Eq. ( 2), the UBC for prescribing EPP-NO y concentrations (in units of cm −3 ) is then given by (in units of cm −2 s−1) is given by Background NO y concentrations, i.e., the NO y contributions not related to the EPP indirect effect, are not negligible in the lower/middle mesosphere and need to be considered when prescribing concentrations at the model's top layer(s).When specifying fluxes, this step is not required because the NO y entering the model domain in polar winters is in good approximation, exclusively originating from the EPP source.We model the background NO y concentrations by fitting the following regression function to the seasonal composite of the background [NO with t being here the day of the year.The regression coefficients a n (φ, z) and b n (φ, z) for pressure levels within 1-0.01 hPa are listed in Tables B1-B7.Note that this parameterization of NO bg y does not provide a full description of the observations since inter-annual variations (e.g., introduced by the QBO) are not considered.For prescription of NO y in models with upper lids above 1 hPa, the consideration of merely seasonal variations of the background NO y is a good approximation.Figure 12 compares the resulting latitudetime NO y distribution of the semi-empirical model to the MIPAS observations at 0.5 and 0.02 hPa.While the background contribution at the latter pressure level is nearly 2 orders of magnitude smaller than the EPP contribution, this is not the case at 0.5 hPa.Here, the background is comparable to the EPP contribution in many NH and SH winters.Overall, the modeled NO y densities reproduce very well the observed latitude distribution and time evolution, except for episodes of large solar proton events (e.g., October/November 2003 and January/March 2012).This is expected since the semiempirical model does not account for the EPP direct effect.
The consideration of ES-related enhanced odd nitrogen descent in the UBC for chemistry climate models is relatively straightforward in nudged model simulations since the ES onset dates, needed to drive the semi-empirical model, are known beforehand.However, its consideration in freerunning simulations would require one to diagnose ES events "online" in a quasi-instantaneous manner, e.g., by analyzing the modeled temperature fields averaged over a narrow time window covering the past model day.Then, in case of the detection of an event onset, it would be required to update the UBC, accounting for the ES event as described in Sect.4, from the ES onset date to at least the end of the actual NH winter season.The main task is hence the implementation of the online ES detection scheme in the model system.We have proposed in Sect. 5 a detection criterion, T 30-70 > 53 K, based on the difference of 0-30 and 70-90 • N temperature averages at 1 hPa, which allowed for quasi-instantaneous detection of ES event onsets (associated with enhanced odd nitrogen descent) in EMAC simulations.Its application in other model systems, however, might require an adjustment of the detection threshold (which could be achieved by calculating T 30-70 from a nudged simulation covering 1980-2014 and tuning the threshold such that it is exceeded only for ES events listed in Table 4).

Historical reconstruction of the EPP indirect effect
The semi-empirical model also allows for a historical reconstruction of EPP-NO y depositions for the period covered by the A p record (i.e., since 1932).This period can by extended by use of the aa index (available since 1868) and the Helsinki A k index (available for .Both aa Atmos.Chem.Phys., 16,2016 www.atmos-chem-phys.net/16/8667/2016/and A k indices provide a similar proxy of geomagnetic activity as A p , however, based on observation from only one (or two) stations.Therefore, both data sets can be combined, although biases have to be accounted for.Such a combined, de-biased data set, expressed as a homogenized A p index, has been generated as part of the solar forcing recommendations for CMIP6, available at http://solarisheppa.geomar.de/solarisheppa/cmip6.A detailed description of the methodology for homogenization of these three indices can be found in Matthes et al. (2016).We use here the extended A p index for the EPP-NO y reconstruction in the period 1850-2014, corresponding to the historical simulation as part of the CMIP6 DECK experiments (Eyring et al., 2016).
Figure 13 shows the reconstruction of seasonal stratospheric EPP-NO y depositions below 0.5 hPa (∼ 50 km) in both hemispheres in the period 1850-2014 covering solar cycles 9-24.The temporal evolution of the EPP indirect effect follows closely that of solar variability on multi-decadal timescales, as expected due to the strong link of solar and geomagnetic activity.On shorter timescales, EPP-NO y depositions also show a solar cycle modulation, however, with maxima shifted in tendency towards the declining phase of the cycle, as expected due to the correlation of A p with the solar wind.A long-term variation of EPP-NO y depositions with highest amounts in cycles 19-22 (Modern maximum) is clearly visible.On average, depositions have increased by a factor of 3 from the Gleissberg minimum around 1900 to the recent Modern maximum.The highest EPP-NO y amounts since 1850 were deposited into the stratosphere during the 1991 SH winter, but the 2003 SH and 2004 NH winters are also among the four strongest EPP winters since 1850.On the other hand, the prolonged solar minimum around 2008 led to exceptionally small EPP-NO y depositions that are as small as during the Gleissberg minimum.In this sense, solar cycle 23 had one of the largest amplitudes of EPP variability in this 164-year period.
Looking at the variability during the last three solar cycles covered by the "satellite era", we model a reduction of the average global EPP-NO y deposition rate of 0.8 GM year −1 that corresponds to 1.7 % of the global production rate by N 2 O oxidation.This is likely to affect the long-term NO y trend by counteracting the expected increase caused by growing N 2 O emissions (about 6 % in the same period), although the impact is most probably limited to mid to high latitudes.

Conclusions
We have presented a semi-empirical model for computation of hemispheric energetic particle precipitation (EPP)-NO y amounts transported to stratospheric and mesospheric pressure levels, as well as the associated vertical fluxes, during Antarctic and Arctic winters.The model has been trained with the EPP-NO y record inferred from MIPAS observations during 2002-2012 (Funke et al., 2014a).Inter-annual variations of the EPP indirect effect at a given time of the winter are related to variations of the EPP source strength, the latter being considered to depend linearly on the geomagnetic A p index.A finite impulse response approach is employed to describe the impact of vertical transport on this modulation at given pressure levels.The seasonal dependence of the EPP-NO y vertical distribution, driven by variations of chemical losses and transport patterns, is assumed to be independent of inter-annual dynamical variability.This assumption is shown to be a reasonably good approximation for SH winter and dynamically quiescent NH winters.For episodes of accelerated descent associated with elevated stratopause (ES) events in Arctic winters, however, this assumption does not hold and a dedicated parameterization of the spatio-temporal EPP-NO y distribution needs to be employed.This parameterization takes into account the dependence of the EPP-NO y amounts and fluxes during ES-related descent episodes on the event timing in accordance with results from the model study of Holt et al. (2013).
In order to consider accelerated descent during ES events in the semi-empirical model, a criterion for ES detection is required.de la Torre et al. ( 2012) identified ES events by looking at abrupt increases in the polar cap stratopause height, the latter defined as the altitude between 20 and 100 km where temperature maximizes.We have shown, by analyzing temperature and CO from EMAC simulations nudged to ERA-Interim reanalysis data, that the upper stratospheric meridional temperature gradient, expressed as the difference of 0-30 and 70-90 • N temperature averages at 1 hPa, provides a reliable alternative criterion for detection of strong descent episodes.Further, the proposed criterion allows for a quasi-instantaneous detection of a commencing ES event since no temporal smoothing of the temperature time series is required.It is therefore well suited to drive the semi-empirical model when used to prescribe NO y concentrations or fluxes in transient simulations with models lacking a detailed representation of the EPP source region.However, due to its definition as a discrete threshold, our criterion may need to be adjusted when applied to other models.
We have quantified the EPP indirect effect in both hemispheres during 1978-2014 with the semi-empirical model, considering the ES events as detected from the EMAC simulations.The resulting wintertime EPP-NO y depositions have been compared to observational estimates from satellite instruments, including LIMS, HALOE, MIPAS, ACE-FTS, and SOFIE.In order to account for multiplicative and additive biases between the different estimates, related to instrumental uncertainties and/or differences in the employed estimation methods, we have adjusted the other instrument's estimates to the MIPAS data by applying an offset and a scale factor.The resulting homogenized time series of observational EPP-NO y deposition estimates is in very good agreement with the results of our semi-empirical model.The simulated average EPP-NO y deposition per year into the stratosphere during 1978-2014 was found to be 1.26 GM in the SH and 0.5 GM in the NH.Strong descent associated with ES events in Arctic winters, while being of high relevance during individual events, led to an increase of only 4 % in the average NH deposition during these 3 decades.
A major purpose of the semi-empirical model is to provide an odd nitrogen upper boundary condition (UBC) for chemistry climate models with their upper lid in the mesosphere and, thus, missing the EPP-NO y production occurring above.This is achieved by distributing the hemispheric EPP-NO y amounts and fluxes at given pressure levels in latitude bands, using the MIPAS average meridional distribution during 2002-2012, and expressing them either as concentrations (in units of cm −3 ) or as molecular fluxes (in units of cm −2 s−1).In order to avoid top boundary artifacts in the models when specifying NO y concentrations at latitudes not dominated by EPP, we also provide a background NO y contribution (from N 2 O oxidation in the stratosphere) obtained from a simple regression model adjusted to the MIPAS seasonal 2002-2012 composite.The resulting UBC model hence provides global zonal mean NO y concentrations and EPP-NO y molecular fluxes on an adaptable pressure level and latitude grid as a function of time for upper model lids within 1-0.01 hPa.
Odd nitrogen UBCs have previously been used in chemistry climate models not extending into the EPP source region for representation of the EPP indirect effect.In some model studies, the UBC was taken directly from NO x observations (e.g., Reddmann et al., 2012;Päivärinta et al., 2013), which, however, implies the restriction to the relatively short time period spanned by the observations.In other cases, a simple parameterization in dependence of the seasonally averaged A p index (Baumgaertner et al., 2009) was employed (e.g., Baumgaertner et al., 2011;Rozanov et al., 2012), enabling extended simulations over multi-decadal time periods.
Our UBC model is designed for the latter application and represents an improved parameterization due to its more detailed representation of geomagnetic modulations, latitudinal distribution, and seasonal evolution, as well as the ability to reproduce odd nitrogen enhancements due to ES events in Arctic winters.It has been successfully tested in simulations carried out with the EMAC model (Matthes et al., 2016).
By employing historical geomagnetic indices, as provided with the CMIP6 solar forcing, we also estimated the EPP indirect effect since 1850.We found long-term changes in solar cycle-averaged stratospheric EPP-NO y depositions on the order of 1 GM, which can be attributed to secular variations of geomagnetic and solar activity.Inter-annual variations along the solar cycle were particularly pronounced during solar cycles 16, 22, and 23, with cycle amplitudes of up to 2.5 GM.We also found a reduction in the EPP-NO y deposition rate during the last 3 decades related to a decline of geomagnetic activity that corresponds to 1.8 % of the NO y production rate by N 2 O oxidation.The negative trend in the geomagnetic activity level is closely related to the reduction of solar cycle amplitudes encountered for cycles 23 and 24.As the decline of solar activity is expected to continue in the coming decades (Steinhilber and Beer, 2013), this is likely to affect the long-term NO y trend by counteracting the expected increase caused by growing N 2 O emissions (Ravishankara et al., 2009).
A limitation of our semi-empirical model for reconstructions on multi-decadal timescales, however, is related to potential secular variations of meridional circulation patterns in the mesosphere (Baumgaertner et al., 2010).A deviation from the dynamical mean state characteristic for the 2002-2012 period could lead to modifications of the EPP indirect effect not considered in our model.However, such dynamically induced variations are expected to be small compared to the geomagnetically induced variations.In particular, it seems unlikely that mesospheric circulation changes could outweigh the simulated reduction of stratospheric EPP-NO y depositions in the last decades related to the decline of solar variability.

Data availability
The semi-empirical UBC model is available as IDL and MATLAB routines at http://solarisheppa.geomar.de/solarisheppa/cmip6 (Funke, 2016) for its use with geomagnetic proxy data provided with the CMIP6 solar forcing (available on the same webpage).

Figure 2 .
Figure 2. Temporal evolution of the vertical distribution of EPP-NO y differential amounts N d (z) during September-May in the NH (left) and March-November in the SH (right) as a result of the empirical model for a constant A p index of 10 (2002-2012 average).Note the different color scales for the NH and SH.

Figure 3 .
Figure 3. (a) Observed (red diamonds) and fitted (solid black) temporal evolution of hemispheric EPP-NO y amounts during 2002-2012 at the pressure levels 0.03, 0.3, 2, and 10 hPa (top to bottom) in the NH.A p (z, t) of the empirical model for quiescent dynamical conditions is shown with blue lines (dotted blue lines indicate A p levels with a spacing of 5 A p units).Shaded areas have been excluded from the fit due to perturbed dynamics (i.e., SSW and ES events, grey) and large SPE events (orange).Note that the y axis range does not cover the very high NO y amounts encountered in the NH 2003/04 winter.(b) As for panel (a) but for the SH.Note the variable y axis range.

Figure 4 .
Figure 4. Temporal evolution of EPP-NO y fluxes F EEP (z) through the vertical level z indicated by the y axis during September-May in the NH (left) and March-November in the SH (right) as a result of the empirical model for the 2002-2012 average A p index of 10.

Figure 6 .
Figure 6.Dependence of the maximum EPP-NO y flux F ES m through the 0.4 hPa level (black) and the descent velocity ω ES m at 0.08 hPa (red) on the ES event timing.

Figure 7 .
Figure 7. Temporal evolution of the vertical distribution of EPP-NO y differential amounts N d (z) during September-May in NH winter with ES events occurring on 15 December (upper left), 15 January (upper right), 10 February (bottom left), and 5 March (bottom right), resulting from the empirical model for a constant A p index of 10 (2002-2012 average).

Figure 8 .
Figure 8. Observed (red diamonds) and modeled (solid black) temporal evolution of NH EPP-NO y amounts during ES winters 2004 (left) and 2009 (right) at pressure levels 0.03, 0.1, 0.3, and 1 hPa (top to bottom).The orange-shaded area indicates the period affected by the large SPE event of October/November 2003.

Figure 9 .
Figure 9. Temporal evolution of T 30-70 (in K) = T (0-30 • N) −T (70-90 • N) at 1 hPa (solid red) and 70-90 • N zonal mean CO anomalies (with respect to the climatological seasonal mean, in ppmv) at 0.5 hPa (solid blue).The threshold of T 30-70 = 53 K for ES detection is indicated by the dotted red line.Detected event onsets are marked by vertical dashed lines.

Figure 10 .
Figure 10.(a) Inter-annual variation of seasonal EPP-NO y depositions during 1978-2014 calculated with the semi-empirical model in the SH below the pressure levels of 1 hPa (top) and 0.1 hPa(bottom).EPP-NO y deposition estimates from MIPAS observations(Funke et al., 2014b) are indicated by filled red diamonds.HALOE-derived estimates ofRandall et al. (2007) are shown by the grey-shaded area (limited by their "average" and "maximum excess NO x " estimates) and are shifted by 0.5 GM (dashed blue line) in order to facilitate comparisons to the MIPAS estimates.Open blue symbols represent the adjusted "maximum excess" estimates (scaled by a factor of 0.7 to fit the MIPAS estimates).(b) As for panel (a) but for the NH and for the 0.7 hPa (top) and 0.1 hPa (bottom) levels.Modeled EPP-NO y depositions without consideration of ES events are also shown by dashed lines.Observational deposition estimates from satellite data are also shown: MIPAS(Funke et al., 2014b) (filled red diamonds);MIPAS (2003MIPAS ( -2004)),ACE- FTS (2005-2009), and LIMS (1979)(Holt et al., 2012) (filled blue diamonds); SOFIE (2009-2013)(Bailey et al., 2014) (filled green diamonds).Open blue and green symbols represent the adjusted estimates ofHolt et al. (2012) andBailey et al. (2014), respectively, after applying a constant offset (colored dotted lines) and an scale factor in order to facilitate comparisons with the MIPAS-derived depositions ofFunke et al. (2014b) (see text for more details).Note that years indicated on the x axis correspond to the second year of the season; e.g.,"2003" means "winter 2002/2003".
Figure 11.Latitudinal distribution of EPP-NO y from MIPAS averaged over 2002-2012, shown as the fraction of the hemispheric differential amount polewards of the indicated latitude, for SH winters (top), NH winters (middle), and ES winters (bottom).Averages for individual months are shown by colored lines.The amountweighted seasonal average is indicated by the thick black line.

Figure 12 .
Figure 12.Latitude-time sections of NO y densities observed by MIPAS (left) and from the UBC model (right) at 0.02 hPa (top) and 0.5 hPa (bottom).

Figure 13 .
Figure 13.Reconstructed stratospheric EPP-NO y deposition below 0.5 hPa (∼ 50 km) in the SH (red) and NH (without ES events, dark blue) during 1850-2015.NH depositions with consideration of ES events (only 1979-2014) are shown with the light blue line.Solar cycle average depositions are indicated with the dashed lines.

Table 1 .
Parameters N m , t N m , and w N of the empirical model for the vertical and temporal distributions of EPP-NO y differential amounts (see Eq. 5).

Table 2 .
Parameters F m , t F m , and w F of the empirical model for the vertical and temporal distributions of EPP-NO y fluxes (see Eq. 13) through a given pressure level in both hemispheres.The seasonally accumulated EPP-NO y amounts T 10 for a constant A p index of 10 (2002-2012 average) are also listed.

Table 3 .
Coefficients of the polynomial expansion n i=0 a i ln(p) i used for parameters of the EPP-NO y model for ES events.Values in parentheses should be read as powers of 10.

Table 4 .
Start and end dates of ES events detected by the T 30-70 > 53 K criterion.