© E. Robert et al., Published by EDP Sciences 2023

1 Introduction

The auroras are characterized by their visible emissions, which are generated by the physical interactions between the Sun, Earth’s magnetic field, and the upper atmosphere. By analyzing particle precipitation, we can gain a deeper understanding of the Sun-Earth system, specifically in regard to the ionosphere, thermosphere, and magnetosphere. Particle precipitation interacts with the molecules and atoms of the thermosphere and ionosphere. This interaction results in the phenomenon of fluorescence. Resulting visible emissions are a good indicator of the environment composition and dynamics. Auroras occur at altitudes between 90 km and 300 km, which is too low for direct measurement by in-situ satellites and too high for high-altitude balloons (Barthelemy et al., 2018). Although rocket-based measurements are feasible, they are not capable of providing continuous observations. Ground-based observations with radar and/or optical instruments provide a complementary way to monitor auroral emissions and retrieve precipitating electron fluxes.

Several papers have already proposed methods for electron precipitation reconstruction (Lanchester et al., 1994, 1997; Dahlgren et al., 2011; Simon Wedlund et al., 2013; Adachi et al., 2017). In the work from Lanchester et al. (1994), the European Incoherent SCATter Association (EISCAT) radar data, photometer measurements, and a forward model of the ionosphere were used to retrieve the energies of precipitating electrons. Dahlgren et al. (2011) use Auroral Structure and Kinetics (ASK) measurements, the Southampton ion chemistry, and the electron transport model as described by Lanchester et al. (1994) to compute the intensity ratio of first negative band at 562 nm, and O emission at 777.4 nm in order to retrieve the mean energy of precipitating electrons in a small field of view of 3.1° × 3.1°. They also use EISCAT data and Auroral Large Imaging System (ALIS) measurements ( 427.8 nm line) with 4 stations to provide an independent method to obtain the mean energy of the precipitating electrons. In Simon Wedlund et al. (2013) study, they implement a two-step inversion method to retrieve a two-dimensional (2D) map of precipitating electron fluxes using optical observations from 4 stations of the ALIS network. The first inversion is a tomographic-like inversion providing the 3D Volume Emission Rate (hereafter VER) of the 1NG band at 427.8 nm. A second inversion using an ionospheric model is then used to retrieve the electron fluxes. EISCAT measurements are used to confirm the validity of the reconstruction for one point along the radar beam directed along the magnetic field.

The previous studies illustrate the usefulness of forward models and inversion methods in determining the energy characteristics of precipitated electrons when combining ionospheric models with optical or radar measurements. In our study, we explore a new approach using an optimization method to compare simulations obtained with the Transsolo ionospheric transport code and VER resulting from ALIS optical observations enabling us to reconstruct electron precipitation characteristics over a large geographical region, approximately 150 km × 150 km in size. This approach complements the previous work of Simon Wedlund et al. (2013) in which a direct inversion of the VER is performed using a forward model of the ionosphere based on a 1D kinetic Monte Carlo model (Sergienko & Ivanov, 1993). Here we use a different transport code (Transsolo) and an optimization procedure to fit the outputs of the code to the observed VER. To compare the validity and accuracy of our method, we use the same data set as in Simon Wedlund et al. (2013), obtained on March 5th, 2008 at 18:41:30 UT and 18:42:40 UT. The goal is first, to test the capabilities of the Transsolo code to reconstruct the precipitation without any other simulations and thus give more consistency to the calculation by using only one kind of approach, and second, to test the differences between the Monte Carlo forward model and the kinetic simulations.

In the upcoming sections, we will introduce the physical model and the optical data, explain the optimization method, and then present and interpret the results of our study with an analysis of the reliability of our approach.

2 Physical model and optical data

In this section, we aim to provide a comprehensive overview of both the kinetic model Transsolo and the optical data utilized in the optimization process. First, we explain the Transsolo kinetic model, discussing its inputs, operation, and outputs. Then, we introduce the ALIS data, which have been used by Simon Wedlund et al. (2013) to calculate the VER using the tomographic-like method outlined in Gustavsson (2000).

In this section, we present the adjustments made in the context of our study in order to provide a clearer picture of the decisions that were taken to ensure an efficient optimization as shown in the next section.

2.1 Transsolo code

The description of Transsolo, as outlined in Figure 1, is as follows: Transsolo is a 1D multi-beam kinetic code that solves the non-conservative Boltzmann equation for suprathermal electrons (Lummerzheim & Lilensten, 1994; Lilensten & Blelly, 2002). It calculates the transport of these particles along magnetic field lines using the DISORT routine, which stands for Discrete Ordinates Radiative Transfer Program for a Multi-Layered Plane-Parallel Medium (Stamnes et al., 1988). The code calculates the interactions of these precipitated particles with atmospheric components and simulates the resultant VER, as described by Vialatte (2017).

Figure 1 Schematic view of the Transsolo code operation. Inputs are represented in blue box and outputs in red box. The yellow boxes are the kinetic transport and photoionisation calculation modules. From Barthelemy et al. (2018).

Figure 2 Geomagnetic activity on March 5th, 2008 recorded by the Tromsø Geophysical Observatory’s magnetometers at Bear Island (BJN), Nordkapp (NOR), Sørøya (SOR) and Tromsø (TRO). Time is along the horizontal axis, while the measured variation in the horizontal component, as well as geomagnetic latitude, is along the vertical axis, a scale indicating 200 nT is in the bottom left corner. The vertical blue dashed line represents the two datasets used for the study at 18:41:30 UT and 18:42:40 UT.

Figure 3 Visualization of VER on March 5th, 2008, at 18:41:30 UT for the 1NG line at 427.8 nm. (a) 2D projection at 100 km altitude centred on Skibotn corresponding to X = 0 km and Y = 0 km location (see the red point), (b) 3D mapping from 100 km to 240 km altitude, (c) VER field-aligned profile at different altitudes for a given location in the map (here X = −24.1 km and Y = 6.6 km at 100 km altitude).

Figure 4 Visualization of VER on March 5th, 2008, at 18:42:40 UT (70 seconds later than dataset 1) for the 1NG line at 427.8 nm. (a) 2D projection at 95 km altitude centred on Skibotn corresponding to X = 0 km and Y = 0 km location (see the red point), (b) 3D mapping from 95 km to 240 km altitude, (c) VER field-aligned profile at different altitudes for a given location in the map (here X = −54.6 km and Y = 3.8 km at 95 km altitude).

Inputs from the user in Transsolo are the following:

• Datetime.

• Geomagnetic index Ap.

• Solar index F _10.7.

• Geographic coordinates (latitude, longitude).

• Energetic characteristics of the electron precipitation: energy flux (E _tot) in erg·cm⁻²·s⁻¹, mean energy (E _m) in eV and type of distribution (among Maxwellian, Gaussian or quasi-monoenergetic¹). To calculate the energy distribution, we consider an exponential grid of 100 energy values (E) ranging from 0.1 eV to 55 keV.

• The angular distribution of the particle precipitation population at the top of the atmosphere is considered as Gaussian. This means that we applied Gaussian weights (Davis et al., 2014) to the flux for the different angular sectors. We consider 32 angles sectors: 16 downward and 16 upward.

To compute the thermospheric components, Transsolo uses the neutral atmospheric NRLMSISE 2000 model (Picone et al., 2002) while the ionospheric components are computed from the IRI model (Bilitza & Reinisch, 2008). These models take as inputs the parameters mentioned above, except for the energy characteristics and angular distribution of the electron precipitation which are used for the calculation of the kinetic transport. NRLMSISE00 considers the following species (and IRI the same corresponding singly charged ions:

• Hydrogen (H),

• Oxygen (O),

• Dioxygen (O₂),

• Nitrogen (N),

• Dinitrogen (N₂).

The Transsolo code produces auroral volume emission rate profiles for these species present in the upper atmosphere, covering an altitude range from 80 km to 520 km, with 100 layers distributed in thickness to ensure equal optical depth per layer (see “TRANSSOLO simulations” in Fig. 5). The VER unit is photons·cm⁻³·s⁻¹.

Figure 5 Simplified diagram of the optimization method. INPUTS – (top) Example of few Transsolo runs for mean energies of 9.5 eV, 102 eV, 497 eV, 3.3 keV, 10 keV, 22 keV and 30 keV, and energy flux set at 10.8 erg·cm−2·s−1, with a Maxwellian distribution; (bottom) Example of an VER field-aligned profile for a given location in the map of the dataset 1 (X = −46.6 km and Y = −13.4 km at 100 km altitude). MINIMIZATION – Cost function from the least-square method to minimize the difference between the VER measurements and the VER simulations. OUTPUTS – Result of the minimization. For this example, only one energy distribution is needed to find the best fit: a mean energy of 1.9 keV and an energy flux of 8.6 erg·cm−2·s−1.

Figure 6 Comparison of the optimization results with (left) a Maxwellian (top)/Quasi-monoenergetic (bottom) energy distribution at location 1 (X; Y) = (−72.9; −126.1), and (right) a Maxwellian (top)/Quasi-monoenergetic (bottom) energy distribution at location 2 (X; Y) = (−46.6; −13.4), for the dataset 1. The blue curve corresponds to the VER field-aligned profile, and the red curve to the volume emission rate of the blue line 1NG 427.8 nm simulated by the Transsolo code.

2.2 The ALIS network

The Auroral Large Imaging System (ALIS) is a network of six stations located in Northern Sweden and Norway at the following sites: Kiruna, Skibotn, Abisko, Silkkimuotka, Merasjärvi, and Tjautjas. A map of the ALIS optical network can be found in Figure 1 from Simon Wedlund et al. (2013). Stations are equipped with high-resolution CCD cameras of 1024 × 1024 pixels with an average field of view of 70°. Cameras are equipped with narrowband interference filters of 4 nm to measure the intensity of the following auroral lines²:

•  1NG 427.8 nm corresponding to the blue emission,

• O I(¹ S – ¹ D) 557.7 nm corresponding to the green line,

• O I(¹ D – ¹ P) 630.0 nm corresponding to the red line,

• O I(3p ³ P – 3s ³ S) 844.6 nm corresponding to the IR line.

The ALIS network datasets consist of pictures taken every 5 s at each station with a spatial resolution of 100 m. More details about the ALIS network can be found in the PhD thesis of Brändström (2003).

By pointing all cameras at the same region of the sky, a tomographic-like reconstruction of the VER in photons·cm⁻³·s⁻¹ can be made from ALIS images. ALIS VER ranges typically from 80 km to 260 km altitude with a spatial resolution of 2.5 km. This tomographic-like reconstruction method was developed by Gustavsson (2000). As already mentioned, the tomographic maps studied here come from the work of Simon Wedlund et al.(2013) courtesy of the Royal Belgian Institute for Space Aeronomy (BIRA-IASB). The tomographic reconstruction maps generated are oriented along the local magnetic field line rather than stacked vertically in altitude. They are displayed as a function of the East-to-West distance in km from the Skibotn station (X), the South-to-North distance in km from Skibotn (Y) and the altitude (Z), i.e. an overall grid size of 100 × 100 × 74 pixels.

The visualization of dataset 1 and dataset 2 is respectively given in Figure 3 and 4. The datasets are further described in the next section.

2.3 Our case study

On March 5th, 2008 a magnetic substorm was recorded by the Norwegian magnetometer stations (see Fig. 2). The expansion phase of the substorm, which is characterized by the negative values in the horizontal component of the magnetic field (H), took place at 18:55 UT. We noticed that the southern and westernmost magnetometer stations (TRO and SOR) were located within the eastward electrojet prior to the onset of the substorm. The strongest disturbance caused by the substorm was observed in the two northernmost magnetometers (NOR and BJN), suggesting that the study area is located slightly to the south of the main activity. The magnetic disturbances recorded were moderate, with amplitudes of around 200 nT, during and after the time frame of our interest.

Our study is based on VER maps located at Skibotn 69.35°N, 20.36°E on March 5th, 2008 at 18:41:30 UT (referred to as Dataset 1) and 18:42:40 UT (Dataset 2), and reconstructed using only Skibotn, Abisko, Silkkimuotka, and Tjaujas stations.The datasets were obtained during the growth phase and prior to the onset of the substorm (see dashed vertical lines in Fig. 2).

2.3.1 Transsolo inputs

The overall objective is to retrieve the energy characteristics of the electron precipitation population at the time of this substorm. To do this, we perform the optimization described in Section 3 by executing several runs of the Transsolo code. Here, we describe the different inputs used.

At the relevant date and time of our data, the values for Ap and F _10.7 were recorded as 7 and 68.2 solar flux units, respectively. However, it should be noted that the Ap index is a global and daily index and may not accurately reflect local geomagnetic activity at the time of data. To account for this, we arbitrarily increased the Ap value to 40, which corresponds to a Kp value of between 3 and 4, and is considered reasonable for this type of substorm (as shown in Fig. 2). This adjustment allows us to take into account the potentially higher intensity of local geomagnetic activity during the event and its effects on the upper atmosphere.

Regarding the energy inputs, Transsolo needs mean energy and energy flux values and a distribution type (see Sect. 1). We prepared different mean energies of the electron distribution, ranging from 3 eV to 35 keV, using a logarithmic grid with 60 values and different flux values from 0.1 erg·cm⁻²·s⁻¹ to 1000 erg·cm⁻²·s⁻¹, also based on a logarithmic grid with 60 values. The choice of the size of the vectors is a compromise between a reasonable computation time (in our case 10 h) and a significant sampling. Two hypotheses were made regarding the shape of the distribution: Maxwellian or quasi-monoenergetic.

As explained in Section 1, the Transsolo code produces auroral volume emission rate profiles for each species and excited energy levels. In our case, the output consists of two sets of the VER intensity cubes, showing the intensity of auroral emission as a function of mean energy (eV), energy flux erg·cm⁻²·s⁻¹, and altitude (km). Each cube has dimensions of 60 (eV) × 60 (erg·cm⁻²·s⁻¹) × 100 (km). However, in this study, we focus on the blue line emission of 1NG 427.8 nm, as it provides information about the energy deposited in the upper atmosphere, according to the research of Simon Wedlund et al. (2013). The blue line emission is exclusively produced by electron impacts on nitrogen molecules, while other complex processes such as photodissociation and recombination also contribute to the emission of the red and green lines (Witasse, 2000).

2.3.2 ALIS datasets

The datasets used for this case study are shown in Figures 3 and 4 with a) a VER map at an altitude of 100 km for dataset 1 and at an altitude of 95 km for the dataset 2 (the X-axis/Y-axis corresponds to the West-East/South-North distance from Skibotn (in km)); b) a 3D representation of VER as a function of altitude (Z) oriented along the local magnetic field line; and c) an example of a VER field-aligned profile as a function of altitude at a given pixel in the map. We chose to limit the size of the VER map data to 81 × 81 × 74 pixel to avoid uncertainties related to the pixels on the edge of the map.

In addition, we chose altitudes between 100 km and 237.5 km for dataset 1, and between 95 km and 237.5 km for dataset 2, because the accuracy of the tomographic data is less reliable beyond 240 km. Indeed this problem lies in the fact that at too high altitude, the signal strength is weaker and therefore the signal-to-noise is weaker. The lower limit is chosen because of the number of points needed for optimization which must take into account a sufficient part of the emission peak and avoid the part with a low signal-to-noise ratio.

For the sake of clarity, not all maps are shown here. Only maps of the VER measurements with maximum emission are shown in Appendix A. For dataset 1, the VER is greatest between 102.5 km and 117.5 km altitude with a maximum emission in the southwest. An auroral arc of moderate intensity is visible to the north of the map. For dataset 2, the maximum emission is located between 100 km and 117.5 km altitude with auroral arc extending from West to East in the north of the map. This arc is assumed to be the same arc as that of dataset 1, but it has intensified.

We can observe some localized artefacts on both VER maps, corresponding to pixels of constant intensity at all altitudes as we can see on the VER measurements of dataset 1 (Fig. 3), where a vertical line is observed at X = 0 on all maps. Other such pixels are visible in dataset 1 and also in dataset 2 (Fig. 4). These artefacts are due to imperfections in the tomographic inversion, caused in particular by the projection of one of the fields of view into the common volume. Nevertheless, they are not a problem for this case study, but they must be taken into account when interpreting the results.

In the next section, we describe the optimization method implemented between the dataset and the Transsolo simulations.

3 Optimization method

3.1 Method

This study aims to determine the energy properties of the particle precipitation population using VER measurements and the Transsolo code. To accomplish this, the study employs a weighted cost function that computes the difference between the experimental measurements (ALIS) and simulated outputs from Transsolo. By minimizing this function, the study aims to identify the Transsolo inputs that provide outputs that are as close as possible to the ALIS measurements, thereby obtaining the energy characteristics of the precipitated electrons. It is noteworthy that, according to the results of Simon Wedlund et al. (2013), multiple energy distributions are needed to properly reconstruct the blue VER field-aligned profiles. Indeed, the VER observed in the upper atmosphere may be attributed to the interaction between different electron populations of varying energy distributions. Since the equation is linear, it would be possible to sum the contribution of each distribution and to fit the output to ALIS measurements (hence having multiple Transsolo run to recover ALIS measures). However, in this study, we used only the lowest peak in the profile which corresponds to the distribution with the highest mean energy. An upper boundary needs to be set to avoid reaching the contributions of other peaks and to optimize the results effectively. However, setting these boundaries is not always straightforward, as it may correspond to a local minimum in some cases but not in others. To prevent any mistakes, we set the boundaries manually. The lower boundary is 100 km for dataset 1 and 95 km for dataset 2. For dataset 1, the upper boundary matrix then ranges from 105 km³ to 160 km across the entire map, while for dataset 2, it ranges from 105 km to 110 km. It is worth noting that the boundaries for the two datasets have been defined differently to balance accuracy and processing time. The lower and upper boundaries for dataset 1 were carefully determined on a pixel-by-pixel basis, while for dataset 2, a first upper boundary value was applied, and then it was adjusted only at the locations where the optimization resulted in extreme values.

The simplified algorithmic process is depicted in Figure 5 and can be summarized as follows:

• We utilize a rectangular parallelepiped of VER measurements with dimensions X (West-East distance from Skibotn) × Y (South-North distance from Skibotn) × Z (altitude). Additionally, a cube of the VER simulations is generated using the Transsolo code with dimensions mean energy × energy flux × altitude.

• The cost function between the VER field-aligned profile and the VER simulated by the Transsolo code is minimized for each geographical position (i.e., VER field-aligned profile).

• We determine the optimal pair of mean energy and energy flux of the VER profile by selecting the best fit of the VER field-aligned profile.

• This procedure is performed for the two distinct distributions separately, which are subsequently compared.

The cost function to minimize is calculated using the least squares method:(1)with z corresponding to the altitude ranging from 95 km/100 km to 240 km, and the weight W of each measurement being:(2)

The weighting W is performed on the relative emission rate of each point of the considered profile and allows us to give more importance to the emissions close to the emission peak.

To keep the calculation time reasonable we set the (E _m (eV), E _tot erg·cm⁻²·s⁻¹) grid to 60 × 60 values. In order to improve the precision of the results we perform an interpolation into the grid. A 2D spline method is used for interpolation⁴.

This optimization method is performed on all VER field-aligned profiles of each dataset. The results are presented as a map shown in Section 4 (Figs. 7–10).

Figure 7 (Left) Mean energy map (in eV) and (right) energy flux map (in erg·cm−2·s−1), resulting from the optimization method applied to the VER measurements for dataset 1, with a Maxwellian energy distribution. The white pixels correspond to the areas where the optimization is not feasible because the energy flux is less than 0.1 erg·cm−2·s−1. The geographic coordinates here correspond, arbitrarily, to the ones at a 100 km altitude.

Figure 8 (Left) Mean energy map (in eV) and (right) energy flux map (in erg·cm−2·s−1), resulting from the optimization method applied to the VER measurements for dataset 2, with a Maxwellian energy distribution. The geographic coordinates here correspond, arbitrarily, to the ones at a 100 km altitude.

Figure 9 (Left) Mean energy map (in eV) and (right) energy flux map (in erg·cm−2·s−1), resulting from the optimization method applied to the VER measurements for dataset 1, with a quasi-monoenergetic energy distribution. The white pixels correspond to the areas where the optimization is not feasible because the energy flux is less than 0.1 erg·cm−2·s−1. The geographic coordinates here correspond, arbitrarily, to the ones at a 100 km altitude.

Figure 10 (Left) Mean energy map (in eV) and (right) energy flux map (in erg·cm−2·s−1), resulting from the optimization method applied to the VER measurements for dataset 2, with a quasi-monoenergetic energy distribution. The geographic coordinates here correspond, arbitrarily, to the ones at a 100 km altitude.

Figure 11 Relative error on measured and simulated altitude of the emission peak (in %) with respect to the VER field-aligned profiles for dataset 1, with an optimization calculation based on (left) a Maxwellian energy distribution, and (right) a quasi-monoenergetic energy distribution. The white pixels correspond to the areas where the relative error is not calculable.

Figure 12 Relative error on measured and simulated altitude of the emission peak (in %) with respect to the VER field-aligned profiles for dataset 2, with an optimization calculation based on (left) a Maxwellian energy distribution, and (right) a quasi-monoenergetic energy distribution.

Figure 13 Relative error on measured and simulated integrated VER (in %) with respect to the VER field-aligned profiles for dataset 1, with an optimization calculation based on (left) a Maxwellian energy distribution, and (right) a quasi-monoenergetic energy distribution. The white pixels correspond to the areas where the relative error is not calculable.

Figure 14 Relative error on measured and simulated integrated VER (in %) with respect to the VER field-aligned profiles for dataset 2, with an optimization calculation based on (left) a Maxwellian energy distribution, and (right) a quasi-monoenergetic energy distribution.

Figure 15 Result of the optimization for dataset 2 at the EISCAT location. We obtain a mean energy of 2.7 keV and an energy flux of 17 erg·cm−2·s−1 with Maxwellian distribution.

3.2 Shape of VER field-aligned altitude profiles

In Figure 6, we show an example of the optimization method performed with both Maxwellian (left) and quasi-monoenergetic (right) energy distributions for dataset 1.

The VER altitude profiles in Figure 6 (blue curves on the left and right) exhibit both narrow and wide shapes depending on their location. As a result, two types of energy distribution are tested in the optimization process to identify the most appropriate energy distribution for the VER field-aligned profile. It is expected that narrow VER profiles are associated with quasi-monoenergetic distributions, while wider VER profiles are associated with wider distributions corresponding to Maxwellian distributions.

4 Result and discussion

The optimization method is applied to the datasets described in Section 2.3. First, we discuss the results of the optimization method, based on Maxwellian and quasi-monoenergetic energy distributions. Next, we analyze the correlation between the VER field-aligned profile’s emission altitude peak and the mean energy. Finally, we present the relative error in mean energy and energy flux associated with the optimization calculation and compare our results obtained with a Maxwellian distribution to the work of Simon Wedlund et al. (2013).

4.1 Maxwellian energy distribution

Figure 7 represents the mean energy and energy flux results for dataset 1 when considering a Maxwellian energy distribution. These maps present the final E _m and E _tot for each pixel and are arbitrarily plotted on the geographical reference frame of 100 km.

The left image in Figure 7 shows the mean energy map with a range of 1.8–5.2 keV, while the right image represents the energy flux map with a range of 0.1–44.3 erg·cm⁻²·s⁻¹. The colorbars are logarithmically scaled, and the maps can be compared to those in dataset 2. The white pixels on the map correspond to areas where the energy flux is too low for the optimization to be effective (less than 0.1 erg·cm⁻²·s⁻¹).

Regarding dataset 2, results are shown in Figure 8: the mean energy ranges from 2.2 to 9.5 keV, while the energy flux ranges from 2.1 to 136.7 erg·cm⁻²·s⁻¹.

By comparing the VER maps of the datasets (see Appendix A) with the reconstructed maps such as the mean energy and energy flux maps (Figs. 7 and 8), we observe that:

• There is no clear link between the column-integrated intensity maps and the mean energy. An arc in the northern part of the maps shows lower mean energies than in more diffuse regions. This observation seems to be counterintuitive. As we can see in Figure A1 of dataset 1, an emission intensity maximum is observable in the southwest region of the VER map, while we have low values of the order of 2–3 keV in this same region on the mean energy map. The same observation can be made on dataset 2 where we have an intensity maximum located in the North and from West to East (see Figs. A3 and A4), while the mean energy in this area is still on the order of 2–3 keV again (see Fig. 8).

• On the contrary, the column integrated intensity is clearly linked to the energy flux. It is also interesting to note that the arc in the north of the map intensifies between the two datasets, indicating higher energy fluxes. The areas of higher fluxes correspond well to the areas of higher intensity emission, and this observation is visible for both datasets 1 and 2.

• The optimization does not work well for narrow emission altitude profiles when using a Maxwellian distribution.

4.2 Quasi-monoenergetic energy distribution

Concerning the mean energy and energy flux maps from the simulations with a quasi-monoenergetic energy distribution, results are shown in Figure 9 for dataset 1 and in Figure 10 for dataset 2. The first observation is that the average energies are globally higher than those of a Maxwellian distribution. The mean energy ranges from 4.2 to 11.8 keV for dataset 1, and from 8 to 17.1 keV for dataset 2. Thus, the assumption of a quasi-monoenergetic energy distribution gives higher average energies. It is also clear that the optimization performed with a quasi-monoenergetic distribution works much better on narrow VER field-aligned altitude profiles (see Fig. 6).

Despite the differences in the datasets and the use of a different energy distribution than the previous case, the conclusions remain consistent:

• The column integrated intensity is clearly not linked to the mean energy of the distribution. It can be observed that these mean energies appear smaller in the arc than in diffuse regions.

• The column integrated intensity is however strongly linked with the energy flux. The values of the energy flux range from 0.1 to 45 erg·cm⁻²·s⁻¹ for dataset 1, and from 2.2 to 110.1 erg·cm⁻²·s⁻¹ for dataset 2. These values are of the same order of magnitude as those obtained in the case of a Maxwellian distribution suggesting that the integral of the VER profiles changes little with the choice of a Maxwellian or quasi-monoenergetic energy distribution, in contrast to the mean energy which is more sensitive to this input assumption. Moreover, we still find a global good fit between the compiled energy flux values and the experimental VER measurements.

4.3 Peak emission altitude of VER field-aligned profile: a proxy for mean energy ?

Kosch et al. (2001), after many authors, suggest that there exists a link between the monoenergetic electron precipitations and energy deposition altitudes. Based on the formula of Del Pozo et al. (1997), the study mentions that “the effective deposition altitude for precipitating electrons can be related to their mean energy E by:(3)where H is altitude in kilometers”.

This function is consistent with the mean energy calculated by the optimization with a quasi-monoenergetic distribution, but, as expected, it does not hold for the Maxwellian case. However, if less evident in the Maxwellian case, we are still able to find a link between the altitudes of the peak and the mean energy.

4.4 Comparison of the optimization results obtained with a Maxwellian energy distribution and a quasi-monoenergetic energy distribution

In this section, we compare the optimization results obtained with a Maxwellian energy distribution on the one hand and a quasi-monoenergetic energy distribution on the other. As mentioned in Section 3.2, the “altitude thickness” of the emission profile can vary and needs different types of energy distribution functions to allow accurate fits. This clearly appears in Figure 6.

The two VER field-aligned profiles in Figure 6 correspond to two different geographic locations taken arbitrarily, around the (0, 0) point corresponding to Skibotn station. Specifically, the left figures correspond to a location at (−72.9 km, −126.1 km), and the right figure corresponds to a location at (−46.6 km, −13.4 km). The top two figures correspond to the optimization performed with a Maxwellian distribution, while the bottom two figures correspond to a quasi-monoenergetic distribution.

On the left figures, it is evident that a narrow VER field-aligned profile is better fitted using a quasi-monoenergetic energy distribution, whereas the fit obtained using a Maxwellian distribution is unsatisfactory. On the right, a Maxwellian distribution is more effective in fitting a wide VER field-aligned profile, while the quasi-monoenergetic distribution is not suitable. These findings are consistent across other geographic coordinates in the map.

To assess the accuracy of the fit obtained through the optimization method using the two energy distributions, we calculate the relative percentage error of the measured and simulated altitude of the emission peak, and the relative percentage error of the measured and simulated integrated VER in order to obtain an indication of the error on the mean energy and on the energy flux. We then compare the integrals of the ALIS VER field-aligned profiles (Int_ALIS) and the corresponding values calculated using Transsolo simulations (Int_Transsolo) at each geographical location to get an estimation of the relative error on the energy flux. Similarly, we compare the peak emission altitude of the ALIS VER field-aligned profiles (Alt_ALIS) and the value obtained through Transsolo simulations (Alt_Transsolo) for each location. The consequences of the estimation of the mean energy are discussed below.

The relative error in the altitude of the emission peak RE_Altitude and the relative error RE_Integral in the integrated VER is calculated by the following expressions:(4) (5)

Maps, Figures 11 and 12, show the relative error in altitude of the emission peak, using Maxwellian or quasi-monoenergetic distribution (respectively left and right), for the two data sets. Identically, Figures 13 and 14 show the relative error in the integral value. Tables 1 and 2 quantify the trends, i.e. whether the result is considered correct or whether it is overestimated or underestimated. In the particular case of the integrated intensity, we consider the optimization result to be satisfactory to within ±30% of its experimental value. This assumption is based on the fact that the cross-sections are defined with an uncertainty that can easily reach 25%. Concerning the altitude of the emission peak, we consider that the simulated result is correct when the value of the emission peak altitude is equal to the value given by the ALIS measurements with a precision lower than the VER altitude resolution (2.5 km).

Dataset 1 – Maxwellian distribution (Fig. 11):

• (Table 1) Overall, the altitude of the emission peak is correctly estimated for 79.2% of the map, and its value is overestimated/underestimated by ±6% at most. Note that the 4% of the maps that are underestimated are mainly located in the low mean energy area.

• (Table 2) Similarly, the total intensity is correctly estimated for 86.4% of the map i.e. within ±30% of the experimental intensity, 1.4% of the map is underestimated, and finally, 12.2% of the map is overestimated mainly in the central-southeast area of the map. 3.4% of the map corresponds to an integral relative error of +70 to +250% mainly in the same region. It appears clearly that the distribution shape does not fit for these cases.

Dataset 1 – Quasi-monoenergetic distribution (Fig. 13):

• (Table 1) The altitude of the emission peak is correctly estimated for 70.7% of the map. The altitude is underestimated up to −6% of the measured value for 1% of the map in the low mean energy area. In the other part of the map, the altitude is overestimated of up to +10%.

• (Table 2) For the column integrated intensity, the estimate is assumed to be correct over the most part of the map (83.4%), with an underestimation in the northern arc (from West to East). Very few cases result in an overestimation of the intensity (0.1%). It appears clearly that in the arc area, the Maxwellian distribution is more appropriate to reproduce the measurements. However, the use of the quasi-monoenergetic distribution limits the number of aberrant intensity values on the map.

For this dataset, the optimization works quite well for both the column integrated intensity and the altitude of the emission peak, whether we use a Maxwellian or a quasi-monoenergetic distribution. However, we notice that some areas are better described by a specific type of energy distribution. This is the case for the arc zone which is better reconstructed by a Maxwellian distribution.

Dataset 2 – Maxwellian distribution (Fig. 12):

• (Table 1) For the most part of the map (84.9% of the geographic locations), the altitude of the emission peak is overestimated (mostly from +2 to +6%). The area of the arc where the mean energy is low is correctly estimated or underestimated by up to −6%.

• (Table 2) For the integrated intensity, we consider that its estimation is correct for about half of the map (49.8%), without correspondence with a precise area of the aurora. The rest is generally overestimated, with a small proportion of extreme values (4.2% beyond +70% relative error.)

Dataset 2 Quasi-monoenergetic distribution (Fig. 14):

• (Table 1) The altitude is mostly well estimated (for 67.2% of the map), while the altitude in the arc region is rather underestimated (up to −10%) and overestimated in the southern region up to +6%, especially where the average energies are rather high.

• (Table 2) The column intensity integral is also mostly well estimated (72.6%). Only the area near the arc is underestimated up to −70%. We do not obtain extreme values as in the previous case with the Maxwellian distribution.

If the quasi-monoenergetic distribution gives better results for this dataset, it should be noted that, again, the arc is slightly better described with the optimization performed with a Maxwellian distribution.

In conclusion, we can observe that whatever the dataset and the energy distribution used in the optimization calculation, the relative error in altitude is always less than ±10% and most often less than ±6%.

Concerning the intensity integral, the error distribution is much more spread out, between −70% and +250%. However, we consider that an integral estimate of ±30% is still quite acceptable. Furthermore, we can note that the error distribution with a quasi-monoenergetic energy distribution is less spread out (−60 to +40% globally) than in the case of the Maxwellian distribution, which gives a minority proportion (up to 3.4% of the areas covered by the map) of extreme relative error values of +70 to +250%.

It is important to notice that the link between the altitude of the peak emission and the mean energy is not linear meaning that an error of 6% on the altitude can produce a much larger error on the mean energy. This is especially critical at low altitudes where the variation of the altitude of the emission peak becomes small. At 100 km, an error of 2.5 km (the resolution in altitude of ALIS tomographic reconstructions) produces a shift of 28% on the mean energy. The most interesting observation is that some structures of the aurora are best described by a specific type of energy distribution, as we have seen with the Maxwellian distribution for the arc. In a more general way, we observe in Figure 6 that narrow ALIS emission maxima profiles were best approximated by a quasi-monoenergetic distribution assumption and that wider emission profiles were best described by a Maxwellian energy distribution. This implies that the choice of energy distribution should ideally be adapted locally, depending on the structures observed in the aurora.

In addition, it should be noted that the altitude thresholds were defined differently for dataset 1 and dataset 2. This may introduce a bias and partly explain the poorer results obtained in the Maxwellian case of the last dataset (dataset 2). However, since there is a real performance gap with the quasi-monoenergetic case for the same data set, it seems that the Maxwellian description of the energy deposition is less suitable for the emissions measured in dataset 2 meaning that the energy distributions are more monoenergetic in this case.

Finally, it is worth emphasizing that each parameter of the energy distribution (E _m or E _tot) seems to be linked to a unique observable quantity (column intensity integral and altitude of the emission peak). This allows us to decouple these two quantities. However, the width and shape of the distribution appear to be an important parameter since some regions do not fit with one of the two tested distribution shapes. A solution to this problem can be to test a Gaussian distribution with variable width. It adds a third input and thus necessitates extending the dimension of the input cubes which complicates the optimization problem. However, the three input parameters will continue to be decoupled. In other words, the altitude of the emission peak stays a proxy of the mean energy, the area defined by the emission distribution a proxy of the energy flux, and the altitude thickness of the emission peak a proxy for the width of the distribution.

4.5 Comparison with results from Simon Wedlund et al. (2013)

The datasets were initially chosen to allow comparison of the results of our optimization method with those of Simon Wedlund et al. (2013). Figure 14 of Simon Wedlund et al. (2013) plots the electron particle fluxes obtained by authors with the ALIS MART algorithm at three different date and hours:

• 05 March 2008 18:41:30 UT (dataset 1),

• 05 March 2008 18:42:40 UT (dataset 2),

• 05 March 2008 18:43:30 UT.

The reconstructions are done on a single location corresponding to the location of the EISCAT UHF beam at 110 km meaning (X; Y) = (−46; 1.9). Dataset 1 corresponds to the dash line and dataset 2 to the continuous line (Fig. 14 of Simon Wedlund et al., 2013). For this comparison, we consider the dataset 2.

Figure 15 represents the result of the optimization at the EISCAT location for the dataset 2. The energy flux obtained by Simon Wedlund et al. (2013) is of the order of 18 erg·cm⁻²·s⁻¹ for dataset 2 in agreement with 17 erg·cm⁻²·s⁻¹ obtained in this work. Regarding the mean energy of the particles, we obtain for dataset 2, E _m = 2.7 keV. This is slightly higher than the energy of the second peak reconstructed by Simon Wedlund et al. (2013) equal to 2 keV. Our results are consistent with those of Simon Wedlund et al. (2013) regarding the mean energy and the energy flux. However, it is important to keep in mind that the extracted distribution will be different since we imposed here a given shape which is not the case in the work of Simon Wedlund et al. (2013). The hypothesis is however mandatory in this optimization to allow to reduce the number of parameters. At the EISCAT location, the fit using a Maxwellian distribution seems much better for this VER field-aligned profile than the quasi-monoenergetic distribution suggesting a wide energy distribution. To conclude this section, this comparison allows us to verify the accuracy of our results and validate the use of this new optimization method, which is complementary to the one used in Simon Wedlund et al. (2013).

5 Conclusion and perspectives

An optimization method has been implemented between the VER measurements of the 1NG 427.8 nm and the Transsolo code simulations. The cost function computed from the least-square method is used to approximate the experimental measurements as accurately as possible. We thus obtained the distribution of the mean energy and energy flux of the precipitated electrons which are responsible for the volume emission rate measured by the ALIS instruments during the night on March 5th, 2008 at 18:41:30 UT (dataset 1) and at 18:42:40 UT (dataset 2).

Optimization calculations performed with a Maxwellian energy distribution yielded a mean energy ranging from 1.8 to 5.2 keV with energy flux from 0.1 to 44.3 erg·cm⁻²·s⁻¹ for dataset 1, and a mean energy from 2.2 to 9.5 keV with energy flux from 2.1 to 136.7 erg·cm⁻²·s⁻¹ for dataset 2. The optimization result at the EISCAT location is consistent with the work conducted by Simon Wedlund et al. (2013) in terms of mean energy and energy flux for the main emission peak. Assuming a quasi-monoenergetic energy distribution over the entire map, we find a mean energy of 4.2 to 11.8 keV with energy flux ranging from 0.1 to 45 erg·cm⁻²·s⁻¹ for dataset 1, and 8 to 17.1 keV with energy flux ranging from 2.2 to 110.1 erg·cm⁻²·s⁻¹ for dataset 2.

Detailing the maps, we observe that the column integral intensity is clearly linked to the energy flux and the altitude of the peak with the mean energy. Furthermore, it is important to point out that depending on the location on the map, the fit is better for monoenergetic distributions when its altitude thickness is thin or with Maxwellian distributions when it is large. We thus are able to differentiate between regions with mono-energetic acceleration processes such as potential drops from regions with other acceleration processes. A further update of this work could be to try to fit the emission profiles using Gaussian distributions and add the width of the Gaussian as a free parameter. This will change the optimization architecture because it will need to fit on a third parameter but will still be possible with those datasets.

Furthermore, our results show that there is linearity between the energy flux and the column-integrated VER and that the coupling between the mean energy and the energy flux is very weak. Therefore, it could be considered, in order to improve the optimization process, to decouple the mean energy and the energy flux and initially perform optimization only on the first free parameter (i.e. the mean energy). Once this parameter is fixed, the second free parameter (i.e. the energy flux) can be used to recover the characteristics of the energy distribution.

The uncertainties on the reconstructions are still important in the order of 30%. This is mainly due to the uncertainties of the atmospheric model and the cross-sections in the Transsolo simulations. Another way to improve these simulations can be to perform such optimization on other emission lines linked to the electron precipitations. In this perspective, the O I 844.6 nm line could be a good candidate since it is produced by two processes directly linked with the precipitation flux (O₂ dissociation and direct excitation of O). We will then be able to compare the results and then better calibrate the precipitation reconstruction.

Acknowledgments

This research has been funded by SpaceAble, the CNES SHM group and the PNST (CNRS-INSU). This work is realized within the framework of a CIFRE thesis between SpaceAble, the Centre Spatial Universitaire de Grenoble (CSUG) and the Institut de Planétologie et Astrophysique de Grenoble (IPAG). We acknowledge the use of the Transsolo code developed by Lummerzheim & Lilensten (1994) and thank Jean Lilensten for providing us code. We thank the ALIS team (IRF, Sweden) for providing the ALIS data and Simon Wedlund et al. (2013) for providing us with the tomographic maps courtesy of the Royal Belgian Institute for Space Aeronomy (BIRA-IASB). We thank Tromsø Geophysical Observatory, UiT – the Arctic University of Norway, for providing the magnetometer data. We would also like to thank SpaceAble for their support and expertise. The editor thanks two anonymous reviewers for their assistance in evaluating this paper.

References

Appendix

VER measurements – Datasets 1 and 2

Figure A1 VER maps for dataset 1 between 102.5 km and 110 km. We can observe an auroral arc of moderate intensity extending from West to East in the north of the map and an emission intensity maximum in the southwest region. The maximum arc intensity is between 107.5 km and 110 km of altitude.

Figure A2 VER maps for dataset 1 between 112.5 km and 117.5 km. The intensity of emissions from the auroral arc and the southwest region decreases with altitude.

Figure A3 VER maps of dataset 2 between 100 km and 107.5 km. We expected to observe the same auroral arc, extending from West to East in the north of the map, that of dataset 1. The emission intensity of the auroral arc increases with altitude. The south-western region with maximum emission intensity presented in dataset 1 seems to have dissipated in this dataset 2.

Figure A4 VER maps of dataset 2 between 110 km and 117.5 km. We can observe the maximum intensity of the auroral arc at 110 km of altitude then a decrease up to 117.5 km.

All Tables

All Figures

	Figure 1 Schematic view of the Transsolo code operation. Inputs are represented in blue box and outputs in red box. The yellow boxes are the kinetic transport and photoionisation calculation modules. From Barthelemy et al. (2018).
In the text

Figure 2 Geomagnetic activity on March 5th, 2008 recorded by the Tromsø Geophysical Observatory’s magnetometers at Bear Island (BJN), Nordkapp (NOR), Sørøya (SOR) and Tromsø (TRO). Time is along the horizontal axis, while the measured variation in the horizontal component, as well as geomagnetic latitude, is along the vertical axis, a scale indicating 200 nT is in the bottom left corner. The vertical blue dashed line represents the two datasets used for the study at 18:41:30 UT and 18:42:40 UT.

In the text

Figure 3 Visualization of VER on March 5th, 2008, at 18:41:30 UT for the 1NG line at 427.8 nm. (a) 2D projection at 100 km altitude centred on Skibotn corresponding to X = 0 km and Y = 0 km location (see the red point), (b) 3D mapping from 100 km to 240 km altitude, (c) VER field-aligned profile at different altitudes for a given location in the map (here X = −24.1 km and Y = 6.6 km at 100 km altitude).

In the text

Figure 4 Visualization of VER on March 5th, 2008, at 18:42:40 UT (70 seconds later than dataset 1) for the 1NG line at 427.8 nm. (a) 2D projection at 95 km altitude centred on Skibotn corresponding to X = 0 km and Y = 0 km location (see the red point), (b) 3D mapping from 95 km to 240 km altitude, (c) VER field-aligned profile at different altitudes for a given location in the map (here X = −54.6 km and Y = 3.8 km at 95 km altitude).

In the text

Figure 5 Simplified diagram of the optimization method. INPUTS – (top) Example of few Transsolo runs for mean energies of 9.5 eV, 102 eV, 497 eV, 3.3 keV, 10 keV, 22 keV and 30 keV, and energy flux set at 10.8 erg·cm−2·s−1, with a Maxwellian distribution; (bottom) Example of an VER field-aligned profile for a given location in the map of the dataset 1 (X = −46.6 km and Y = −13.4 km at 100 km altitude). MINIMIZATION – Cost function from the least-square method to minimize the difference between the VER measurements and the VER simulations. OUTPUTS – Result of the minimization. For this example, only one energy distribution is needed to find the best fit: a mean energy of 1.9 keV and an energy flux of 8.6 erg·cm−2·s−1.

In the text

Figure 6 Comparison of the optimization results with (left) a Maxwellian (top)/Quasi-monoenergetic (bottom) energy distribution at location 1 (X; Y) = (−72.9; −126.1), and (right) a Maxwellian (top)/Quasi-monoenergetic (bottom) energy distribution at location 2 (X; Y) = (−46.6; −13.4), for the dataset 1. The blue curve corresponds to the VER field-aligned profile, and the red curve to the volume emission rate of the blue line 1NG 427.8 nm simulated by the Transsolo code.

In the text

Figure 7 (Left) Mean energy map (in eV) and (right) energy flux map (in erg·cm−2·s−1), resulting from the optimization method applied to the VER measurements for dataset 1, with a Maxwellian energy distribution. The white pixels correspond to the areas where the optimization is not feasible because the energy flux is less than 0.1 erg·cm−2·s−1. The geographic coordinates here correspond, arbitrarily, to the ones at a 100 km altitude.

In the text

	Figure 8 (Left) Mean energy map (in eV) and (right) energy flux map (in erg·cm−2·s−1), resulting from the optimization method applied to the VER measurements for dataset 2, with a Maxwellian energy distribution. The geographic coordinates here correspond, arbitrarily, to the ones at a 100 km altitude.
In the text

Figure 9 (Left) Mean energy map (in eV) and (right) energy flux map (in erg·cm−2·s−1), resulting from the optimization method applied to the VER measurements for dataset 1, with a quasi-monoenergetic energy distribution. The white pixels correspond to the areas where the optimization is not feasible because the energy flux is less than 0.1 erg·cm−2·s−1. The geographic coordinates here correspond, arbitrarily, to the ones at a 100 km altitude.

In the text

	Figure 10 (Left) Mean energy map (in eV) and (right) energy flux map (in erg·cm−2·s−1), resulting from the optimization method applied to the VER measurements for dataset 2, with a quasi-monoenergetic energy distribution. The geographic coordinates here correspond, arbitrarily, to the ones at a 100 km altitude.
In the text

	Figure 11 Relative error on measured and simulated altitude of the emission peak (in %) with respect to the VER field-aligned profiles for dataset 1, with an optimization calculation based on (left) a Maxwellian energy distribution, and (right) a quasi-monoenergetic energy distribution. The white pixels correspond to the areas where the relative error is not calculable.
In the text

	Figure 12 Relative error on measured and simulated altitude of the emission peak (in %) with respect to the VER field-aligned profiles for dataset 2, with an optimization calculation based on (left) a Maxwellian energy distribution, and (right) a quasi-monoenergetic energy distribution.
In the text

	Figure 13 Relative error on measured and simulated integrated VER (in %) with respect to the VER field-aligned profiles for dataset 1, with an optimization calculation based on (left) a Maxwellian energy distribution, and (right) a quasi-monoenergetic energy distribution. The white pixels correspond to the areas where the relative error is not calculable.
In the text

	Figure 14 Relative error on measured and simulated integrated VER (in %) with respect to the VER field-aligned profiles for dataset 2, with an optimization calculation based on (left) a Maxwellian energy distribution, and (right) a quasi-monoenergetic energy distribution.
In the text

	Figure 15 Result of the optimization for dataset 2 at the EISCAT location. We obtain a mean energy of 2.7 keV and an energy flux of 17 erg·cm−2·s−1 with Maxwellian distribution.
In the text

	Figure A1 VER maps for dataset 1 between 102.5 km and 110 km. We can observe an auroral arc of moderate intensity extending from West to East in the north of the map and an emission intensity maximum in the southwest region. The maximum arc intensity is between 107.5 km and 110 km of altitude.
In the text

	Figure A2 VER maps for dataset 1 between 112.5 km and 117.5 km. The intensity of emissions from the auroral arc and the southwest region decreases with altitude.
In the text

	Figure A3 VER maps of dataset 2 between 100 km and 107.5 km. We expected to observe the same auroral arc, extending from West to East in the north of the map, that of dataset 1. The emission intensity of the auroral arc increases with altitude. The south-western region with maximum emission intensity presented in dataset 1 seems to have dissipated in this dataset 2.
In the text

	Figure A4 VER maps of dataset 2 between 110 km and 117.5 km. We can observe the maximum intensity of the auroral arc at 110 km of altitude then a decrease up to 117.5 km.
In the text