© The Authors 2024

1 Introduction

Observations of Jupiter in the UV range have revealed the existence of extremely bright polar auroral phenomena. The first observational evidence of these phenomena was provided by observations of UV emissions from atomic hydrogen (H Lyman-α emission) and molecular hydrogen (Lyman band emissions: ${\text{B}}^1{\text{Σ}}_{\text{u}}^{+}\to {\text{X}}^1{\text{Σ}}_{\text{g}}^{+}$ and Werner bands: ${\text{C}}^1{\text{Π}}_{\text{u}}^{}\to {\text{X}}^1{\text{Σ}}_{\text{g}}^{+}$) during the flyby of the Voyager 1 spacecraft (Broadfoot et al. 1979). In addition, Jovian aurorae have been extensively studied by the International Ultraviolet Explorer (IUE) for approximately two decades (Clarke et al. 1980; Livengood et al. 1992; Gladstone & Skinner 1989; Harris et al. 1996). These observations, in the mid-UV range (between 120 nm and 170 nm), have allowed the characterization of the power of the aurorae, and also the study of their structure, variability, and intensity. Thanks to these spectral measurements, the first models of UV auroral emissions were developed (Yung et al. 1982; Gladstone & Skinner 1989). This demonstrated that the direct excitation of molecular hydrogen by electrons and absorption of CH₄ below 140 nm could, overall, reproduce UV auroral spectra well in the range of 120–170 nm.

The Hubble Space Telescope (HST) has dramatically contributed to the study of Jovian aurorae, thanks to observations of auroral structures by the Faint Object Camera (FOC; Dols et al. 1992; Gérard et al. 1993, 1994; Prangé et al. 1998) and Wide Field Planetary Camera 2 (WFPC2; Clarke et al. 1996, 1998; Grodent et al. 1997), as well as UV spectral observations by the Space Telescope Imaging Spectrograph (STIS) on board HST (HST/STIS; Gustin et al. 2002). Despite these numerous studies, the morphology of the Jovian aurora remains very difficult to describe in an exhaustive way as their structure is complex and includes spatially and temporally variable substructures (see Grodent 2015). However, these various observations have made it possible to characterize, in a simple way, the morphology of Jupiter’s aurorae: the main auroral emissions form a partially closed oval, with highly variable structures in the polar region inside the main emissions and equatorward emission structures outside, including the footprints of Io, Europa, and Ganymede.

The auroral emission on Jupiter occurs due to the interaction between its magnetosphere and atmosphere. This interaction leads to the precipitation of energetically charged particles along magnetic field lines. In these regions, electrons are the primary species that precipitate, as stated by Rego et al. (2001). In the UV domain [80 nm, 180 nm], Jupiter’s auroral spectral emission is dominated by H Lyman-α emission and the de-excitation of H₂ molecules by electronic transitions: ${\text{B}}^1{\text{Σ}}_{\text{u}}^{+}\to {\text{X}}^1{\text{Σ}}_{\text{g}}^{+}$ (Lyman bands) and the R¹, P, and Q branches of the ${\text{C}}^1{\text{Π}}_{\text{u}}^{}\to {\text{X}}^1{\text{Σ}}_{\text{g}}^{+}$ transition (Werner bands). Moreover, this spectral emission is strongly influenced by the presence of hydrocarbons in Jupiter’s atmosphere, including methane, which mainly absorbs UV photons below 140 nm.

The main studies resulting from the various spectral observations of Jupiter’s auroral regions are those concerning the energy characterization of the electrons precipitating in these regions. Although numerous probes have studied Jupiter’s magnetosphere at the spacecraft altitude, such as Voyager (Russell 1993), Ulysses (Zarka 1998), Galileo, Cassini-Huygens (Hansen et al. 2004), and the New Horizon flyby (Krupp 2007), measurements and observations of magnetospheric plasma have never characterized the energy distributions of the electrons that precipitate immediately above the polar aurora’s altitudes. There are complex processes, such as acceleration by inertial Alfvén waves (Hess et al. 2010, 2013; Saur et al. 2018) and by whistler waves (Elliott et al. 2018), and ion and electron inverted-V structures (Mauk et al. 2017, 2018; Clark et al. 2017), that alter the energy flux distribution of electrons between the altitude of measurement and the high-latitude ionospheric regions. Thus, to measure the shape of the energy flux distributions of electrons precipitating in auroral regions, only low altitude measurements can be effective.

For the present study, we used a combination of spectral observations and modeling to study these characteristics. The method we adopted is phenomenological; it was first proposed by Yung et al. (1982) and is based on the far-UV (FUV) color ratio (CR) of the auroral H₂ emission spectrum. Our aim is to take advantage of the wavelength-dependent absorption of auroral emission by hydrocarbons such as CH₄, which is the third most abundant molecule in Jupiter’s stratosphere (e.g., Moses et al. 2005; Hue et al. 2018). Thus, by considering the unab-sorbed part of the H₂ emission spectrum, we can derive a CR that allows us to characterize the energy distributions of the electron fluxes that precipitate in these regions, as well as their characteristic energies. In the case of absorption by methane, this ratio is defined by $\text{CR}=\frac{I(155\text{\hspace{0.17em}nm}-162\text{\hspace{0.17em}nm})}{I(123\text{\hspace{0.17em}nm}-130\text{\hspace{0.17em}nm})}$ (Gustin et al. 2013), where $I\left({\lambda }_{\min }-{\lambda }_{\max }\right)={\displaystyle {\int }_{{\lambda }_{\min }}^{{\lambda }_{\max }}I} \text{d}\lambda$ and I is the spectrum flux intensity and where the ranges [123–130 nm] and [155–162 nm] represent, respectively, a range absorbed by CH₄ and the unabsorbed part of the spectrum. This allows us to infer, for a fixed emission angle, that over the range [123–130 nm] increasing spectral absorption means that electrons penetrate deeper into the atmosphere before being thermalized, and are therefore more energetic. This method is very advantageous since it does not require the use of absolute spectra of auroral emission as CR can be measured using only arbitrary units within the same spectra.

The relationship between the CR and the characteristic electron energy E₀ (CR(E₀)) is monotonic, and is modeled using a combination of an electronic transport model and a H₂ UV emission model in the auroral regions. This supposes modeling the excitation of H₂ molecules by electron collisions, before calculating their de-excitation from rovibrational levels producing UV emissions. Each modeled spectrum allows the CR to be linked to the characteristic energy of the electrons precipitated in the atmospheric model. Finally, by varying the characteristic energy of the electrons in the transport model, we build the relationship CR(E₀), which is compared with the observed CR to estimate E₀.

Several studies (e.g., Trafton et al. 1994, 1998; Grodent et al. 2001; Gustin et al. 2002, 2016; Ajello et al. 2005; Gérard et al. 2014) have used this technique to characterize the electron energy in auroral regions using HST observations. In most of these studies, the auroral electron transport was modeled by initial mono-energetic, Maxwellian, or kappa phenomenologi-cal flux distributions (see the example by Gustin et al. 2016). However, despite the high spectral resolution of the HST/STIS observations exploited in previous studies, the signal-to-noise ratio (S/N) was limited and the spatial coverage of the aurora was partial and highly dependent on the planet’s tilt axis. Thus, it is only since Juno’s arrival (Bolton et al. 2017) and the Ultraviolet Spectrograph (UVS) observations (Gladstone et al. 2017; Bonfond et al. 2017) that we have had full access to Jupiter’s northern and southern local time. In addition, UVS observations are highly spatially resolved near Jupiter’s closest approach, with a spectral resolution of around 0.6–2.4 nm. This allows correct spectral sampling with a better S/N than the HST/STIS observations.

In this study, we map for the first time the characteristic energy of electrons precipitating in Jupiter’s auroral regions using Juno/UVS observations. We have developed a new UV emission model of H₂, inspired by the models of Dols et al. (2000), Gustin et al. (2002) and Menager (2011), in a more optimized version that takes into account nine H₂ electronic states including cascade excitation and auto-absorption in the Lyman and Werner bands. This model is now available for the community, and can be used in every electron transport model. The excited states of H₂ are calculated through the outputs of our TransPlanet electronic transport model (Stamnes & Rees 1983b; Gronoff 2009; Menager 2011; Benmahi 2022). Additionally, for the CR modeling, we modeled the electronic transport using mono-energetic initial electron flux distributions and a kappa distribution (Coumans et al. 2002; Scherer et al. 2018) derived from Jupiter Energetic Particle Detector Instrument (JEDI) measurements obtained during the first 20 perijoves (PJs; Salveter et al. 2022) of the Juno mission.

The outline of our study is as follows. We first describe the electronic transport model and the UV emission model. In the second step, we describe the Juno/UVS observations and explain the mapping of the characteristic energy method. Finally, we present our results and discussions before concluding.

2 Models

2.1 Electron transport model: TransPlanet

To simulate electron precipitation in Jupiter’s atmosphere, we used the TransPlanet model developed in collaboration with the Institut de Planétologie et d'Astrophysique de Grenoble (IPAG). This transport code was first created by Lilensten et al. (1989) and was modified and improved by Blelly et al. (1996) for application to terrestrial cases (Simon et al. 2007). The model was diversified and adapted to several planets over the years. The core of this algorithm was used in Trans-Mars (Witasse et al. 2002, 2003; Simon et al. 2009; Nicholson et al. 2009), Trans-Venus (Gronoff et al. 2007, 2008), Trans-Titan (Lilensten et al. 2005a,b; Gronoff et al. 2009a,b), Trans-Uranus, Trans-Jupiter (Menager et al. 2010), Aeroplanets (Gronoff et al. 2012a,b, 2014), and recently Trans-Planet (Benmahi 2022). The Trans* code core is divided into two parts: a kinetic part that calculates the interaction of precipitating electrons with atmospheric particles and a fluid part that is implemented only in the TRANSCAR and TRANS4 versions (Lathuillere et al. 1997; Simon et al. 2007), which uses a 13-moment fluid closure description, calculating, among other parameters, the number density, velocity, heat flux, and plasma temperature (of electrons and ions). Thus, in the Trans* version we used for this study, there is no fluid part. Compared with the various existing electron transport codes, Trans* allows multi-stream modeling of electron transport with electron scattering over a wide range of electron energies and pitch angles. The electronic transport model we use in this study is detailed in Appendix A. The uncertainties in that class of models has been studied in Gronoff et al. (2012a,b) and highlighted that the cross sections are one of the major sources of uncertainties. Efforts have been made to improve cross section datasets through the Atomic and Molecular Cross section for Ionization and Aurora Database (ATMOCIAD; Gronoff et al. 2021) and also through comparison with the experimental data (Wedlund et al. 2011).

For this study, we modeled electron transport taking into account only magnetospheric electron precipitation. Secondary electrons resulting from ionization by solar UV radiation are neglected as their penetration capacity in Jupiter’s atmosphere is low and considerably above the homopause of the hydrocarbons considered, such as CH₄, C₂H₂, or C₂H₆.

The atmospheric model of the auroral region we used to model electronic transport is described in Grodent et al. (2001) and presented in Fig. 1. This model is 1D, and takes into account the majority of neutral species (H, H₂, He, and CH₄) that predominate in Jupiter’s atmosphere. It extends from the tropopause at a pressure of ~100 mbar (altitude ~100 km above the cloud level) to the upper thermosphere (altitude ~2300 km above the cloud level), corresponding to a pressure of ~10⁻⁹ mbar. The initial density of thermalized electrons considered in the model is that obtained by Hinson et al. (1998) from radio occultations during the Voyager 2 flyby. Because of the limited data available, the initial electron temperature is thought to be similar to the temperature of the neutral atmosphere. In addition, as the atmospheric model used is 1D, we do not take into account any spatial or temporal variability in the abundance of neutral species in auroral regions, particularly the variability of CH₄. Hence, since methane is the only tracer used in this study to model CR, any variability in its abundance can influence the CR. Thus, in this study we consider a homogeneous and steady chemical composition in Jupiter’s whole aurora, which probably represents a significant approximation.

The electron-matter interactions considered in the physics of electron transport are the elastic interactions in which total kinetic energy is conserved, and the inelastic interactions of electrons with atmospheric particles illustrated in Table 1. The elastic cross sections e + H (Kingston & Walters 1980), e + He (Porter et al. 1987), e+H₂ (Muse et al. 2008), and e +CH₄ (Davies et al. 1989) considered are measured in different energy ranges, which may not entirely cover the energy grid² ranges needed for the electron transport modeling. For cross sections that do not cover the entire energy grid considered, we used power-law extrapolations to fill the gap. Above 400 eV the cross section σ(E) is considered to be proportional to E^−0.65 (Wedde & Strand 1974), and above 2.2 keV the cross section is proportional to E⁻¹ (Rees 1989). For inelastic cross sections we used the same approach as for elastic collisions. The power-law evolution of these cross sections at high energies makes it possible to use a decreasing logarithmic extrapolation to cover the entire range of the modeling energy grid.

As described in Appendix A, by solving the Boltzmann equation we calculate the electron flux F(z, E) as a function of altitude and electron energy. This flux results from the interaction between magnetospheric electrons precipitating in auroral regions and the neutral atmospheric particles considered in our model. To model auroral emission by H₂, we used F(z, E) as an initial condition in our UV emission model.

Fig. 1 Atmospheric model described by Grodent et al. (2001), which considers only the neutral compounds (H, H2, He, CH4, and C2H2) that predominate in Jupiter’s atmosphere. For electronic transport modeling, only H, H2, He, and CH4 compounds were considered. For the UV emission model, spectral absorption by CH4 and C2H2 was taken into account.

2.2 H₂ UV emission model

The R and P branches of the Lyman band of ${\text{H}}_2({\text{B}}^1{\text{Σ}}_{\text{u}}^{+}\to {\text{X}}^1{\text{Σ}}_{\text{g}}^{+})$ correspond to the group of rovibrational transitions that produce spectral lines with wavelengths in the range [80– 190 nm]. For the Werner band, in addition to the R and P branches, there is a third branch, the Q branch, corresponding to the rovibrational transitions that produce spectral lines in the spectral range [80–160 nm]. There are also other transitions in the UV spectrum of H₂, at shorter wavelengths, arising from the excited levels ${{\text{B}}^{\prime }}^1{\text{Σ}}_{\text{u}}^{+}$, ${B^{″}}^1{\text{Σ}}_u^{+}$, ${\text{D}}^1{\text{Π}}_{\text{u}}^{-}$, ${\text{D}}^1{\text{Π}}_{\text{u}}^{+}$, ${D^{\prime }}^1{\text{Π}}_u^{-}$ and ${{\text{D}}^{\prime }}^1{\text{Π}}_{\text{u}}^{+}$, and whose spectral emissions are less intense compared to Lyman and Werner band emissions and lie respectively in the wavelength ranges [85–125 nm], [79–110 nm], [75–110 nm], [75–110 nm], [78–107 nm], and [78–107 nm]. For this study, the H₂ UV emission model in auroral regions we developed takes into account the excited states B, C, B′, B″, D, and D′, as illustrated here.

According to the atmospheric model used (molecular abundances and thermal profile), we begin by calculating the number density n(z, X, v_i, J_i) at altitude z of the H₂ ground state levels. Thus, assuming that neutral species are thermalized in the atmosphere, the population of the ground state of H₂ follows the Boltzmann distribution given by $$n\left(z,X,v_i,J_i\right)=n_{{\text{H}}_2}(z)\frac{g_I(i)\left(2J_i+1\right)e^{-\frac{E_i}{k_{B}T(z)}}}{{\displaystyle {\sum }_k{g}_I(k)\left(2J_k+1\right)e^{-\frac{E_k}{k_{B}T(z)}}}},$$(1)

where $n_{{\text{H}}_2}(z)$ (cm⁻³) is the density of H₂ at altitude z; X, v_i, and J_i are respectively the ground electronic level n_i, and the vibrational and rotational quantum numbers; ɡ_I(i) and E_i are respectively the degree of degeneracy of the i state and its energy; k_B is the Boltzmann constant; and T(z) the temperature at altitude z. In the denominator, the sum is made over all the rovibrational ground state levels of H₂.

An H₂ molecule can be excited into a n_j, v_j, and J_j state by various processes. It can be excited directly by absorbing a photon or by collision with an electron or other atmospheric particles. It can also be excited in this state by cascade de-excitation from higher states. Unlike models of H₂ UV auroral emission that use the Born³ approximation to calculate the excitation rates of the different excited states of H₂ (e.g., Waite et al. 1983), in our model we calculate the excitation rates of the considered electronic levels through electronic transport by modeling ${\text{e}}^{-}+{\text{H}}_2\to {\text{H}}_2^{\ast }+{\text{e}}^{-}$ collisions. Transitions from the EF, GK, and $\text{H}\overline{\text{H}}$ states to the ${\text{X}}^1{\text{Σ}}_{\text{g}}^{+}$ ground state are forbidden due to the selection rule on g → g transitions, and thus a non-negligible part of the B and C states are populated by these transitions (Liu et al. 2002). In this model we also take into account the excitation of H₂ to the EF, GK, and $\text{H}\overline{\text{H}}$ states, as well as the cascade populating of the B and C states. Excitation by other collisional processes with neutral particles is neglected because the atmospheric temperature is not high enough to produce UV emission from collisions of H₂ molecules with neutral particles (e.g., ${\text{H}}_2+{\text{H}}_2={\text{H}}_2^{\ast }+{\text{H}}_2$).

Thus, the volume excitation rate (cm⁻³ s⁻¹) of a rovibrational state j is a linear combination of the direct excitation rate g_direct and the cascade excitation rate g_cascad and is given by $$g\left(z,n_j,v_j,J_j\right)=g_{\text{direct }}\left(z,n_j,v_j,J_j\right)+g_{\text{cascad }}\left(z,n_j,v_j,J_j\right).$$(2)

2.2.1 Direct excitation rate

The direct excitation rate of H₂ by electron collisions is described by the following formula: $$g_{\text{direct }}\left(z,n_j,v_j,J_j\right)={\displaystyle \sum _{i}n}\left(z,X,v_i,J_i\right){\displaystyle \int {\sigma }_{ij}}(E)F(z,E)\text{d}E.$$(3)

Here the indices j and i are used to identify the upper quantum state and the ground state, respectively; F(z, E) (cm⁻² s⁻¹ eV⁻¹) is the electron flux at altitude z with energy between E and E + dE and is modeled by electron transport; σ_ij(E) (cm²) is the excitation cross section of the jlevel from the i level by collision with an electron of energy E (see Liu et al. 1998), and is given by $$\begin{array}{l}{\sigma }_{ij}(E)=4\pi a_0^2{f}_{ij}\frac{R_y^2}{E{E}_{ij}}[{\text{C}}_0\left(\frac{1}{x^2}-\frac{1}{x^3}\right)\hfill \\ +{\displaystyle \sum _{k=1}^4{\text{C}}_k(x-1){\text{e}}^{-k\alpha x}+{\text{C}}_5+\frac{{\text{C}}_6}{x}+\ln (x)],}\hfill \end{array}$$(4)

where a₀ is the Bohr radius, $f_{ij}={1.4992}^{-16}{A}_{ji}{E}_{ij}^2\frac{2J_i+1}{2J_j+1}$ (dimen-sionless) is the oscillator strength of the transition between the i and j levels, and where A_ij (s⁻¹) is the Einstein factor of the j → i transition; $R_y=\frac{m_e{\text{e}}^2}{8h^2{ϵ}_0^2}$ is the Rydberg constant; E_ij is the energy of the transition from the i state to the j state; E is the energy of the incident electron; and $x=\frac{E}{E_{ij}}$. The coefficients C_k and α were obtained experimentally by Liu et al. (1998, 2003) by fitting the excitation functions of the transitions ${\text{X}}^1{\text{Σ}}_{\text{g}}^{+}\to {\text{B}}^1{\text{Σ}}_{\text{u}}^{+}$ and $X^1{\text{Σ}}_{ɡ}^{+}\to C^1{\text{Π}}_u$ (see Table 2). This parameter was measured for the ungerade (odd) and gerade (even) levels of the excited states B and C (Liu et al. 1998, 2003). However, these factors were not measured for the B′, B″, D, D′ excited levels, and in this study, following Menager (2011), we consider that these coefficients are also valid for all ungerade states.

2.2.2 Cascade excitation rate

Electron collisions populate the g states including the EF, GK, and $\text{H}\overline{\text{H}}$ states from the ${\text{X}}^1{\text{Σ}}_{\text{g}}^{+}$ ground state. However, since g → g are forbidden dipolar transitions, the process of populating the g levels differs from that of the u levels. This populating is described by the same relationship (see formula (3)) as that for the u states, but with a different cross section given by $${\sigma }_{ik}(E)=F(x)FC\left(v_i,J_i,n_k,v_k,J_k\right)S_r\left(J_i,J_k\right),$$(5)

where i always refers to the ground state and k to the upper g state (EF, GK or $\text{H}\overline{\text{H}}$), $x=\frac{E}{E_{ik}}$, E is the energy of the incident electron, and E_ik is the excitation threshold of the i → k transition.

The function F(x) describes the excitation of the EF state from the ground state. This function is given by the formula $$\begin{array}{l}F(x)=\pi a_0^2\frac{R_y}{E}{\text{C}}_5[\frac{{\text{C}}_0}{{\text{C}}_5}\left(\frac{1}{x^2}-\frac{1}{x^3}\right)+{\displaystyle \sum _{m=1}^4\frac{{\text{C}}_m}{{\text{C}}_5}}(x-1){\text{e}}^{-m\alpha x}\\ +\left(x-\frac{1}{x}\right)],\end{array}$$(6)

where C_m and a are described in Table 2 for transitions to g states.

F(x) was measured only for the excited EF state. For the GK and $\text{H}\overline{\text{H}}$ states, Liu et al. (2002) suggested using the same excitation function, but multiplying it by a scaling factor to take into account the excitation efficiency of the different g states. This factor is 0.8 for GK and 0.35 for $\text{H}\overline{\text{H}}$.

FC(v_i, J_i, n_k, v_k, J_k) are the Franck-Condon factors that describe the overlap of the wave functions of the i and k states and depend on the quantum numbers v_i and J_i for the ground state and on n_k, v_k, and J_k for the upper g state. The Franck-Condon factors follow the selection rules ΔJ = 0, ±2 and were calculated by Hervé Abgrall and Evelyne Roueff for Liu et al. (2002).

Finally, the function S_r(J_i, J_k) represents the rotational terms that were calculated by Abgrall et al. (1999) and adopted by Liu et al. (2003), and are given by the following formula: $$\begin{array}{l}{S}_r\left(J_i,J_k\right)=\beta {\delta }_{J_i,J_k}+(1-\beta )[\frac{3\left(J_k+1\right)\left(J_k+2\right)}{2\left(2J_k+3\right)\left(2J_k+5\right)}\delta J_{i_i,J_{k+2}}\\ +\frac{J_k\left(J_k+1\right)}{\left(2J_k-1\right)\left(2J_k+3\right)}\delta J_i,J_k+\frac{3J_k\left(J_k-1\right)}{2\left(2J_k-1\right)\left(2J_k-3\right)}\delta J_i,J_{k-2}].\end{array}$$(7)

Here δ is the Kronecker parameter and β is an anisotropy parameter for which Liu et al. (2003) recommended a value of 0.6.

Cascade excitation of the B and C states mainly increases emission from the low vibrational levels of B and accentuates emission from the C state by a smaller proportion. The cascade excitations of the B′, B″, D, and D′ states are not taken into account due to the lack of appropriate data.

2.2.3 Volume emission rate

The discrete volume emission rate η (cm⁻³ s⁻¹) of a transition from a state j to a state i at altitude z is given by $$\eta \left(z,n_j,v_j,J_j\to X,v_i,J_i\right)=g\left(z,n_j,v_j,J_j\right)\frac{A_{j\to i}}{A_j^{\text{tot }}},$$(8)

where $A_j^{\text{tot}}\left({\text{s}}^{-1}\right)$ is the total Einstein factor of the upper level (j) and is given by $$A_j^{\text{tot}}=A_j^{\text{cont}}+{\displaystyle \sum _j{A}_{j\to i}},$$(9)

with $A_j^{\text{cont }}$ the probability that level j transmits into the continuum. This results in the dissociation of H₂ into two fragments with kinetic energy E_c, whose expression is given by $$A_j^{\text{cont }}={\displaystyle \sum _{J_i=J_j-1}^{J_j+1}{\displaystyle {\int }_0^{\infty }{A}_{j\to X,E_c,J_i}^{E_c}}}\left(E_c\right)\text{d}{E}_c,$$(10)

where $A_{j\to X,E_c,J_i}^{E_c}\left(E_c\right)\left({\text{s}}^{-1}{\text{eV}}^{-1}\right)$ is the differential probability of dissociation of the j state into two fragments of kinetic energy E_c and quantum number J_i.

The Einstein factors for the B → X, C → X, B′ → X, and D → X transitions were calculated by Abgrall et al. (1994) and are available in the MOLAT⁴ database. Those for transitions B″ → X and D′ → X were obtained by personal communication from Abgrall and Roueff (Menager et al. 2010). The differential probability of dissociation $A_{j\to X,E_c,J_i}^{E_c}\left(E_c\right)$ are also available in the MOLAT database and were calculated by Abgrall et al. (1997) only for transitions from the B, C, B′, and D states.

The differential volume emission rate in the continuum ${\eta }_{\lambda }^{\text{cont }}$ (cm⁻³ s⁻¹ nm⁻¹) comes from excited states above the dissociation threshold of H₂. Its intensity at wavelength λ and altitude z after the dissociation of molecules from excited state j to dissociated state i is given by $${\eta }_{\lambda }^{\text{cont }}\left(z,n_j,v_j,J_j\to v_i,J_i\right)\text{d}\lambda =g\left(z,n_j,v_j,J_j\right)\frac{A_{\lambda }\left(n_j,v_j,J_j\to v_i,J_i\right)}{A_j^{\text{tot }}}\text{d}\lambda ,$$(11)

where A_λ(n_j,v_j,J_j v_i,J_i) is the differential probability ofdis-sociation of state j, as a function of wavelength λ, and is obtained directly from $A_{j\to X,E_c,J_i}^{E_c}\left(E_c\right)$.

Quantum transitions n_j,v_j,J_j → n_i, v_i, J_i can be represented by the wavelength λ which corresponds to the energy of each transition. Thus, by substitution, for discrete transitions we can write η(z, n_j,v_j,J_j → X, v_i, J_i) = η^discr(z, λ), and for continuum transitions we can write ${\eta }_{\lambda }^{\text{cont }}\left(z,n_j,v_j,J_j\to v_i,J_i\right)={\eta }_{\lambda }^{\text{cont }}(z,\lambda )$.

2.2.4 Auto-absorption

When an emitted photon from an excited state of H₂ is energetic enough to excite an H₂ molecule initially in the ground state, it can be reabsorbed by H₂ to emit another photon of lower energy. This self-absorption therefore tends to attenuate the UV emission spectrum toward short wavelengths (below 120 nm) and amplify it toward low-energy wavelengths. To take account of this phenomenon in the model we describe here, we used the results of Jonin et al. (2000), who experimentally studied the UV spectrum of H₂ in the wavelength range [90 nm; 120 nm]. Thus, considering the volume emission rate η_j→i(z₀) at altitude z₀ by the transition between quantum levels i and j, the volume emission rate reaching altitude z is given by $${\eta }_{j\to i}^{\text{transmitted }}(z)={\eta }_{j\to i}\left(z_0\right){\kappa }_{j\Rightarrow i}(z),$$(12)

where κ_j→i(Z) is a dimensionless attenuation factor depending on the extinction coefficient ϵ_ji, the column density ζ_i(z, z₀) of the quantum state i, and is given by $${\kappa }_{j\to i}(z)=1-\frac{1}{1+0.9948{\left({ϵ}_{ji}{\zeta }_i\left(z,z_0\right)\right)}^{1.44}}$$(13)

with $${ϵ}_{ji}=A_{ij}\frac{2J_i+1}{2J_j+1}\frac{2472\times {10}^{-6}{\lambda }_{ij}^3}{T{(z)}^{0.5}},$$(14)

where λ_ij is the wavelength of the j → i transition given in centimeters.

Thanks to this approach, self-absorbed photons are redistributed to lower-energy transitions using branching ratios calculated with the appropriate Einstein factors. In addition, as most of the continuum emission of the excited states of H₂ that we consider in this study occurs above 120 nm, continuum self-absorption is not taken into account in this model.

2.2.5 Synthetic spectrum

The flux I(λ) (cm⁻² s⁻¹ nm⁻¹ sr⁻¹) of UV emission from the atmosphere, in the θ direction, at infinite spectral resolution and without taking into account absorption by hydrocarbons (e.g.; CH₄, C₂H₂, C₂H₆, …) is given by $$I(\lambda )=\frac{1}{4\pi \cos (\theta )}\underset{z_0}{\overset{\infty }{{\displaystyle \int }}}{\eta }_{\lambda }^{\text{tot}}(z,\lambda )\text{d}z,$$(15)

where ${\eta }_{\lambda }^{\text{tot}}(z,\lambda )$ (cm⁻³ s⁻¹ nm⁻¹) is the total differential volume emission rate at altitude z, which is the linear combination of the differential volume emission rate in the continuum $({\eta }_{\lambda }^{\text{cont }}(z,\lambda )$ and the differential volume emission rate of discrete transitions $\left({\eta }_{\lambda }^{\text{discr}}(z,\lambda )\right)$ in (cm⁻³ s⁻¹ nm⁻¹)).

The thermal agitation of H₂ molecules having a mass $m_{{\text{H}}_2}$ leads to the discrete emission lines with a spectral broadening defined by $$\frac{\text{Δ}\lambda }{\lambda }=\sqrt{\frac{2k_{\text{B}}T}{m_{{\text{H}}_2}{c}^2}}.$$(16)

In our atmospheric model, at altitudes of around 400 km above the cloud level (average altitude of the auroral emission peak according to Bonfond et al. 2015) where the average temperature is around 600 K (Grodent et al. 2001), the mean spectral broadening of H₂ emission at 140 nm is around Δλ ~ 0.002 nm. This Doppler broadening is well below the broadening of the instrumental spectral resolution. In addition, since optical depth is very low above 400 km altitude, all spectral lines are optically thin. Thus, to calculate the differential volume emission rate of discrete transitions ${\eta }_{\lambda }^{\text{discr }}(z,\lambda )$, we used this mean Doppler broadening for all the spectral lines considered in our model at all altitudes. Thus, considering the Doppler profile f_Δλ(λ) (nm⁻¹) with a full width at half maximum Δλ representing the mean Doppler broadening in the atmospheric model used, and using the volume emission rate of discrete transitions η^discr(z, λ), we have ${\eta }_{\lambda }^{\text{discr}}(z,\lambda )={\displaystyle {\sum }_{{\lambda }^{\prime }}{\eta }^{\text{discr}}\left(z,{\lambda }^{\prime }\right)f_{\text{Δ}\lambda }\left(\lambda -{\lambda }^{\prime }\right)}$. Finally, the synthetic spectrum I_synth(λ) is calculated, taking into account the instrumental resolution Δλ′ with $$I_{\text{synth}}(\lambda )=I\left({\lambda }^{\prime }\right)\ast g_{\text{Δ}{\lambda }^{\prime }}\left(\lambda -{\lambda }^{\prime }\right),$$(17)

where ɡ_Δλ′(λ) is a Gaussian function with full width at half maximum Δλ′ and λ′ is a dummy variable.

2.2.6 Hydrocarbon absorption

When magnetospheric electrons penetrate deeply enough below the homopause, the emission produced by the de-excitation of H2 is attenuated at certain wavelengths by absorption from hydrocarbons. The attenuated spectral emission may be given by $$I(\lambda )=\frac{1}{4\pi \cos (\theta )}\underset{z_0}{\overset{\infty }{{\displaystyle \int }}}{\eta }_{\lambda }^{\text{tot}}(z,\lambda ){\text{e}}^{-{\tau }^{\text{tot}}(\lambda ,z)}\text{d}z,$$(18)

where ${\tau }^{\text{tot }}={\tau }_{{\text{CH}}_4}+{\tau }_{{\text{C}}_2{\text{H}}_2}+{\tau }_{{\text{C}}_2{\text{H}}_6}+\dots$ is the total optical depth of hydrocarbons considered in the atmospheric model used; z₀ is the minimum altitude; and ${\tau }_m(z,\lambda )={{{\displaystyle \int }}^{\text{}}}_z^{\infty }{n}_m\left(z^{\prime }\right){\sigma }_m(\lambda )\text{d}{z}^{\prime }$ is the optical depth of species m, whose number density is given by n_m(z), whose absorption cross section is given by σ_m(λ), and where z′ is a dummy variable. In this study we used the cross sections of CH₄ (Au et al. 1993; Kameta et al. 2002; Lee et al. 2001) and C₂H₂ (Cooper et al. 1995; Nakayama & Watanabe 2004; Wu et al. 2001) measured experimentally in the UV range.

2.2.7 Comparison with laboratory spectra and validation

To validate the H₂ UV emission model, we compared our simulations with experimental results from the study of Liu et al. (1995). In this study, the authors used a 100 eV mono-energetic electron beam for bombarding H₂ molecules and measured the UV emission. The spectral resolution was of about 0.0125 nm in the spectral range (114–170 nm). In addition, H₂ molecules with a density of $n_{{\text{H}}_2}=4.55\times {10}^{12}{\text{cm}}^{-3}$ were placed in a cell at a pressure of around 0.4 µbar and a temperature of 300 K.

To validate our calculations, we compared relative spectra and neglected absolute intensities between the modeled and observed spectra in order to avoid the quantity of H₂ used in the Liu et al. (1995) experiment. To do this, we considered a thin atmosphere with a temperature of 300 K, in which we precipitated 100 eV mono-energetic electrons (best fit A) in order to match the experimental conditions of Liu et al. (1995) as faithfully as possible. We then followed the same approach by precipitating a Maxwellian electron (best fit B) beam with an average energy of 100 eV. In Fig. 2, we represent this comparison graphically in four spectral ranges (panels i to iv), each spanning 1 nm in order to evaluate the differences at very high spectral resolution.

For the rest of the spectrum, our model is in very good agreement with the experimental results. We note that the amplitudes of some spectral lines are not perfectly reproduced by the model. These differences are minimal and are generally caused by the inhomogeneity of the electron energy spectrum exciting H₂ molecules. In reality, the incident electron beam is not perfectly mono-energetic, and not perfectly Maxwellian either. As shown in the study by Dols et al. (2000), the variation in the energy of the electrons exciting the H₂ molecules may have a strong effect on some particular spectral lines. This influences their widths and amplitudes.

We also note a spectral shift of around 0.005 nm at some wavelengths between model and observation (e.g., between 123.8 and 124 nm in panel i, between 121.7 and 122 nm in panel ii, and between 120.4 and 121 nm in panel iii). According to Dols et al. (2000), this shift is caused by thermal expansion of the structure of the spectrometer’s sensor during the measurements by Liu et al. (1995).

Fig. 2 Examples of comparisons of modeled synthetic spectra with experimental spectra (red line) obtained by Liu et al. (1995) with a spectral resolution of 0.0125 nm. The blue line and the dashed green line represent the best fits obtained by precipitating a 100 eV mono-energetic (best fit A) and Maxwellian (best fit B) electron flux distributions. The four spectral windows shown have a width of 1 nm. In the ii window around 121.567 nm, the Lyman-α line was filtered out.

3 Juno/UVS observations

The aim of the Juno mission, launched in August 2011, is to study the planet Jupiter and its environment (Bolton et al. 2017). Its insertion into a highly elliptical polar orbit was achieved on July 5, 2016, and its first PJ was carried out on August 27, 2016. Since then the spacecraft has made several dozen PJs, with a periodicity of around 53.5 days during the nominal mission (until PJ37), leading to close flybys of the polar regions allowing us to study the interaction of the magnetosphere with the Jovian atmosphere. Juno hosts several scientific instruments, including the UltraViolet Spectrograph (UVS; Gladstone et al. 2017). The UVS is specifically designed to study Jupiter’s atmosphere and auroral emissions in the extreme UV (EUV) and FUV domains. The wavelengths ranging from 68 nm to 210 nm are dispersed over a 256 spatial channel × 2048 spectral channel sensor (Davis et al. 2011; Greathouse et al. 2013; Gladstone et al. 2017). The spectrometer’s slit has a dog-bone shape, and is oriented parallel to the axis of the spacecraft’s rotation. This slit has a field of view at the edges of 2.55° × 0.2° and a spectral resolution of around 1.9–3.0 nm, and a field of view at the center of 2° × 0.025° with a spectral resolution of ~1.3 nm (Greathouse et al. 2013). Juno is a spin-stabilized probe with a period of about 30 s. As a result, the UVS slit is scanned across the sky to measure the UV emission spectrum in its field of view including the emission spectrum from Jupiter’s poles (Bonfond et al. 2017). Each detected photon is associated with ancillary information including (latitude, longitude), x and y coordinates on the UVS detector, wavelength, emission angle from the planet. Counts recorded by UVS are converted into physical flux units using the instrument effective area derived from thousands of stellar observations during regular calibration phases (Hue et al. 2019, 2021). This photon list is rearranged in latitude–longitude-wavelength data cubes for each hemisphere. Latitude and longitude is sampled every 1° and we used a 0.1 nm spectral sampling that fulfills the Nyquist criterion based on UVS spectral PSF. In addition, to increase the S/N of the UV emission spectra, we only map photons measured by the two UVS wide slits, for which the spectral resolution is around 2.1 nm (Greathouse et al. 2013), and discarded the photons coming from the narrow slit.

For the present study, we used spectral data obtained during PJ32, from 2021-Feb.-21 16:23:09 UTC to 2021-Feb.-21 22:38:45 UTC, for a total acquisition time of 6.25 h. This dataset includes 1.17 h acquisition time for the northern polar region and 4.04 h for the southern polar region, with about 1 h between acquisitions at the two regions (Juno’s passage over the equator).

The acquisition time at the south pole is naturally longer, due to the inclination of the semimajor axis of Juno’s orbit to Jupiter’s equatorial plane. During each 30 s spin, the UVS field of view intercepts Jupiter. Thus, each point on Jupiter or on the sky is looked at with an exposure of ~18 ms during one spacecraft spin.

Figure 3 shows the integrated non-absorbed (N.A.) UV emission over the northern and southern polar regions. In order to isolate the auroral photons from the solar emission backscattered by the Jovian atmosphere, we established a selection criterion for pixels within the aurora (see Fig. C.1). We only selected pixels corresponding to UV emission spectra with a signal-to-noise ratio (S/N) ≥ 3. To evaluate the S/N, we consider the average of the unabsorbed part of the UV emission spectrum of H₂ in the spectral range [155 nm; 162 nm] to define the signal. For the noise, we estimated it within the same spectral range by subtracting the average signal value and calculating the standard deviation.

The UVS observations are co-added over the acquisition time into a large spectral datacube (latitude vs. longitude vs. wavelength), from which the $\text{CR\hspace{0.17em}CR}=\frac{I(155 \text{nm}-162 \text{nm})}{I(125 \text{nm}-130 \text{nm})}$ is then calculated, characterizing the absorption of the UV emission spectrum by CH₄. It is important to note that we cannot use the initially defined wavelength range at the denominator (i.e., 123130 nm) because of the uncertain calibration due to the detector degradation, which is due to gain sag on and around Lyman-α in the wide slit region of the detector. Hence, we start at 125 nm instead. Polar maps of the CR in the northern and southern hemispheres are shown in Fig. 4. As a consequence of this different wavelength range, the minimum CR in these maps is about 1.8, which is higher than the minimum CR ~ 1.1 observed by Gustin et al. (2013) and corresponds to an unabsorbed UV emission spectrum. Regarding the maximum CR value, we have an overall CR_max ≈ 30 for both poles.

Fig. 3 Integrated non-absorbed (N.A.) UV emission from Jupiter’s auroral regions observed by Juno/UVS during PJ32 in the SIII jovicentric reference frame. The acquisition time over the northern polar region is 4282 s. At the southern polar region the acquisition time is 14 561 s. The plus sign (+) in red represents the average solar longitude during the selected acquisition times for the northern and southern polar regions individually. Thus, the Sun’s longitudinal path between its mean position during acquisition at the northern hemisphere and the southern hemisphere is approximately 130° westward.

Fig. 4 Color ratio of Jupiter’s auroral regions (left: northern hemisphere; right: southern hemisphere) observed during PJ32. The panels show the CR calculated for each pixel of the UV emission map as defined by Gustin et al. (2016) with $\text{CR}=\frac{I(155 \text{nm}-162 \text{nm})}{I(125 \text{nm}-130 \text{nm})}$.

Fig. 5 Examples of kappa distributions with different characteristic energies. The three examples use κ = 2.5 and Q0 = 10 erg cm–2 s–1, and three different characteristic energies: E0 = 1 keV, E0 = 10 keV, and E0 = 65 keV.

4 Method

To map the characteristic energy of primary electrons precipitating in auroral regions, we modeled the CR using the TransPlanet electronic transport model combined with the H2 UV emission model. For these simulations, we considered the Grodent et al. (2001) atmospheric model (an atmosphere of H, H2, He, and CH₄) with an altitude range from 100 km (~1 mbar) to 2300 km (~5.3×10^–12 bar) above the cloud level (see Fig. 1). To simulate the electron transport, we used two types of initial electron flux distribution. First, we used a mono-energetic distribution Φ(E, z_max) characterized only by a characteristic energy E₀ and given by Φ(E, z_max) ~ δ(E – E₀), where z_max represents the altitude at which the initial electron flux is injected into the atmosphere. In a second step, we used a kappa-type distribution Φ(E, z_max) ~ ƒ_κ(E, 〈E〉) (Coumans et al. 2002) characterized by an average energy 〈E〉 and a κ parameter governing the logarithmic gradient of the distribution toward high energies. The kappa distribution used in this study is given by $$f_{\kappa }(E,\langle E\rangle )=Q_0\frac{4}{\pi }\frac{\kappa (\kappa -1)}{{(\kappa -2)}^2}\frac{E}{\langle E\rangle }\frac{{\langle E\rangle }^{\kappa -1}}{{\left(\frac{2E}{\kappa -2}+\langle E\rangle \right)}^{\kappa +1}},$$(19)

where Q₀ is the total energy flux and 〈E〉 is given as a function of the characteristic energy E₀ by the expression $\langle E\rangle =2E_0\frac{\kappa }{\kappa -2}$, with E0 representing the energy of the maximum amplitude of the distribution.

In Fig. 5, we give some examples of kappa electron flux distributions. Unlike the Maxwellian distributions used in previous electron transport models (e.g., Gustin et al. 2016), this distribution extends to higher energies. The value of kappa (κ = 2.5) that controls the amplitude of the distribution toward high energies and that we use in this study was derived from observed electron fluxes by Juno/JEDI during the first 20 PJs (Salveter et al. 2022). Hence, this kappa value gives rise to a realistic distribution of electron energy flux precipitating in Jupiter’s auroral regions.

After modeling the electron transport, the resulting Φ(E, z) electron flux is used to calculate the excitation rates of H₂ before simulating the atmospheric spectral emission in a given θ direction. Depending on the characteristic energy of the initial precipitated electron flux distributions, the electrons penetrate to varying depths into the atmosphere. As a result, the UV spectral emission of H₂ is absorbed to varying degrees by hydrocarbons below the homopause.

In the UV emission model, we considered absorption by CH₄, which absorbs mainly below 140 nm (see Fig. 6 for the optical depth of CH₄ and C₂H₂). We also included absorption by C₂H₂, whose distribution profile is shown in Fig. 1. This second hydrocarbon absorbs mainly in the spectral range [150–153 nm] and in other small ranges below 140 nm (see Fig. 6). This means that longward of 145 nm, only C2H2 can attenuate the UV emission spectrum. However, short of this wavelength, CH₄ is the major absorber and C₂H₂ has a weak influence on the amplitude of the UV spectrum.

During PJ32, the close flyby of Jupiter’s poles is such that each point in both polar zones (north or south) is observed with a different viewing angle. Accordingly, each spectrum is measured with a different median θ emission angle (see Fig. C.2).

In the next step, for each characteristic energy, we modeled the spectral emission with a spectral resolution of 2.1 nm (comparable to UVS spectra), varying the emission angle between 0° and 80°. This allowed us to establish the relationship CR(E₀, θ), which links the CR, the emission angle, and the characteristic energy of the initial electron flux distribution injected at the top of the modeled atmosphere. With the assumption that the CR(E₀, θ) function is monotonic, we could then invert this relationship to calculate E₀ using the CR observed by UVS and the emission angle maps.

In this study we distinguished the CR(E₀, θ) relationship in the northern and southern hemispheres according to the magnetic dip angle in these auroral regions (see Appendix D.1). Thus, for each type of electron flux distribution and for each polar region, the relationship CR(E₀, θ) is modeled and fitted⁵ using the formula D.2 (see Appendix D.2). For the case of mono-energetic initial electron flux distribution, the CR(E₀, θ) relationship we modeled is shown in gold in the left panel of Fig. D.2 for the north pole, and separately in the left panel of Fig. D.3 for the south pole. In these same figures, we plot the fit of the relationship CR(E₀, θ) in green grid lines⁶, and we evaluated the uncertainty of this fit by plotting the absolute differences between the fit and the model for the north and south poles. The mean absolute difference is around 0.6, which is comparable to the uncertainty of the CR(E₀, θ) relationship that we modeled for the case of the mono-energetic distribution (see Table D.1).

For the case of the kappa distribution, the results are shown in the right panel of Fig. D.2 for the north pole and in the right panel of Fig. D.3 for the south pole, using the same conventions. The average absolute difference between fits and models is around 0.5, which is comparable to the uncertainty of the modeled CR(E₀, θ) relationship in this case. Whatever the initial electron flux distribution, our results shows that neglecting the emission angle leads to an underestimation of the characteristic energy.

Fig. 6 Optical depth calculated over the atmospheric column for CH4 and C2H2. The transparent green and cyan bands represent the absorption spectral ranges used for the CR calculations. For CH4 the absorption spectral range is considered between 125 nm and 130 nm and for C2H2 between 150 and 153 nm. The transparent magenta band is the non-absorbed (N.A.) spectral range over which hydrocarbon absorption was assumed to be negligible.

Fig. 7 Characteristic energy maps E0 obtained from the CR(E0, θ) relationships (see Figs. D.2 and D.3) modeled for the case of an initial mono-energetic electron flux distribution, and from the CR observed during PJ32, at the north (left panel) and south (right panel) poles. Iso-energy lines are defined for values of 1, 5, 10 keV, in steps of 20 keV between 10 and 300 keV, and then in steps of 100 keV between 300 and 900 keV.

5 Results and discussion

In comparison with the CR(E₀) relationship modeled by Gustin et al. (2016) by precipitating a Maxwellian electron flux, our results show that the CR increases about two times faster in the case of a mono-energetic distribution. In their model, Gérard et al. (2014) also modeled the CR(E₀) relationship when precipitating mono-energetic electron fluxes. For a fixed θ emission angle, our results show that the CR increases about three times faster than that of Gérard et al. (2014). This can be explained by the fact that the modeling of auroral spectral emission in these previous studies is different from our case. In our UV emission model we calculate the volume emission rate (VER) by considering all the excitation rates of all the rovibrational levels (v, J) of the H₂ electronic states that we take into account. However, in the studies by Gérard et al. (2014) and Gustin et al. (2016), the UV emission of H₂ molecules is modeled using the volume emission rate (VER) obtained directly from the cross sections of the interactions ${\text{e}}^{-}+{\text{H}}_2\to {\text{e}}^{-}+{\text{H}}_2^{*}\left({\text{B}}^1{\text{Σ}}_{\text{u}}^{+}\right)$ and ${\text{e}}^{-}+{\text{H}}_2\to {\text{e}}^{-}+{\text{H}}_2^{*}\left({\text{C}}^1{\text{Π}}_{\text{u}}\right)$ measured by Dalgarno et al. (1999). The unabsorbed spectrum is obtained by multiplying the VER by a synthetic spectrum of H₂ (see Eq. (6) in Gustin et al. 2016). The absorbed spectrum is calculated taking into account the abundance of hydrocarbons in the atmosphere in the same way as explained previously. Thus, the differences in calculated CRs are mainly due to very small differences in the VER peak altitude obtained from these different studies. In addition, the electronic transport model used in these previous studies is based on Monte Carlo simulations. Compared with our transport model, which uses a radiative transfer solver, this could introduce additional discrepancies in the results. In comparison with the studies that modeled the CR(E₀) relationship, we still have a monotonically increasing CR as a function of characteristic energy.

5.1 Mapping the characteristic energy using the mono-energetic initial distribution of the electron flux

Using the modeled CR(E₀, θ) relationship (shown in Figs. D.2 and D.3), we inverted Eq. (D.2) to determine E₀(CR, θ). We obtained characteristic energy maps E₀ from the observed PJ32 CR and emission angle maps shown in Figs. 4 and C.2, respectively.

Figure 7 shows the resulting characteristic energy maps for the mono-energetic case in the northern (left panel of Fig. 7) and southern (right panel of Fig. 7) auroral regions. The average estimated uncertainty on the characteristic energy is derived from the average uncertainty on the CR(E₀, θ) modeling presented in Table D.1. For the mono-energetic distribution, it is around 6 keV for both auroral regions. At the north pole the estimated maximum characteristic energy is around (150 ± 6) keV in the polar emission region. There are also other peaks with characteristic energies of up to (120 ± 6) keV in the injection zone of the main emission region near the pole. For the south pole the estimated maximum characteristic energy is around (170 ± 6) keV in the main emission region, with secondary peaks ranging from (90 ± 6) keV to (150 ± 6) keV and spread throughout the auroral arc.

In the polar emission zone the characteristic energy is low, with a peak of around (60 ± 6) keV. During PJ32 the region of polar emission in the south was not very bright. This may be due to exceptionally low-energy electron precipitation or it could be due to exceptionally low electron flux in general, which does not produce strong UV emission. For other PJs (e.g., PJ4, PJ5, PJ8, PJ13, PJ14), the polar emission region in the south pole was relatively bright, and its spectral emission was largely detectable by UVS (Greathouse et al. 2021).

We also observed that there is no similarity between the characteristic energies of the auroral regions connected by the magnetic field lines between the north and south poles. This seems, at first order, to indicate a different electronic precipitation between the two auroral regions. However, the overall range of energies along the auroral ovals appears to be similar.

In our modeling of the CR(E₀, θ) relationship we did not take into account the width of the horizontal extension of the auroral emission resulting from main emissions and Io’s magnetic footprint (Bonfond et al. 2009). Thus, the energies inferred in the very narrow regions bordering the auroral ovals are not usable. In consequence, the results presented here are not valid for the Io footprint region. In addition, we did not take into account the SIII longitudinal motion of Io’s footprint during the acquisition times in the northern and southern auroral regions. This could significantly decrease the emission and the CR of Io’s footprint. Furthermore, according to Bonfond (2010), the average energies of electrons precipitating in this region are around 1–2 keV, which, according to our model, cannot be detected. For electrons with energies below 5 keV, the CR ratio produced is minimal because they cannot penetrate deeply enough below the CH₄ homopause.

Fig. 8 Characteristic energy maps E0 obtained from the CR(E0, θ) relations (see Figs. D.2 and D.3) modeled for the case of an initial kappa electron flux distribution, and from the CR observed during PJ32, at the north (left panel) and south (right panel) poles. Iso-energy lines are defined in the same way as in Fig. 7.

5.2 Mapping the characteristic energy using kappa initial distribution of the electron flux

In the same way as for the case of the mono-energetic distribution, we also derived characteristic energy maps for kappa distributions of electrons precipitating in the auroral regions. Figure 8 displays the characteristic energy maps for the northern and southern polar regions (left and right panels of Fig. 8, respectively). The mean energy uncertainty of these maps is estimated around 4 keV for both auroral regions. At the north pole the maximum characteristic energy is around (40 ± 4) keV in the polar emission region. In the main emission regions we obtained a maximum characteristic energy of around (10 ± 4) keV. At the south pole we obtained peaks between 30 and 50 keV in the main emission oval.

The characteristic energy maps in Figs. 7 and 8 cannot be directly compared. In the case of a mono-energetic distribution, the characteristic energy is identical to the average energy of the distribution, whereas for the kappa distribution E₀ represents the energy of the distribution peak. In the energy maps shown in Fig. 8, each pixel represents a kappa distribution described by a characteristic energy E₀ and a parameter κ = 2.5. Therefore, only the average energy of the kappa distribution can be used for comparison with a mono-energetic distribution. For a Kappa distribution, the conversion from characteristic energy to mean energy can be obtained from the relation $\langle E\rangle =2E_0\frac{\kappa }{\kappa -2}$.

Figure 9 represents the corresponding mean energy maps of electrons precipitating in the auroral regions for the case of a kappa distribution in the north (left panel of Fig. 8) and in the south (right panel of Fig. 8). The mean energy 〈E〉 being proportional to $2\frac{\kappa }{\kappa -2}$, for a parameter κ = 2.5 its value is therefore ten times larger than the characteristic energy. The uncertainty on the energies is also ten times larger. Comparison with the maps in Fig. 7 shows that the spatially averaged mean energy of electrons precipitating in auroral regions is significantly underestimated when the CR(E₀, θ) relationship is modeled by mono-energetic distributions. However, for a fixed emission angle (θ = 0°) we compared the CR(〈E〉) relationship, at low energies, between the two cases of electron flux distribution precipitating in auroral regions (see panel c in Fig. D.4). We found that the CR(〈E〉) relationship obtained for the case of a mono-energetic distribution overestimates the mean electron energy below 〈E〉 = 90 keV. Above 90 keV, the energy of precipitating electrons is underestimated by the mono-energetic distribution assumption. This value of 〈E〉 = 90 keV, which represents the intersection between the two CR(〈E〉) relationships (obtained for the case of a mono-energetic electron flux distribution and for the case of a kappa distribution), seems to be linked to the atmospheric model used and particularly to the CH₄ homopause considered (see panel c in Fig. D.4).

Previous studies that modeled the relationship between the mean energy of precipitating electrons and CR in Jupiter’s auroral regions obtained mean energy maps that differ from our results. In the case of Gérard et al. (2014), the mean-energy maps were based on HST observations from January 2014 and considered mono-energetic distributions of electrons precipitating in the atmospheric model described by Grodent et al. (2001). Gérard et al. (2014) obtained mean energy peaks of up to 500 keV in the northern polar emission region, which is comparable to the energy peaks we obtained in the case of a kappa distribution (see Fig. 9). However, this agreement should be qualified by the fact that these observations were obtained almost 10 yr apart and may be very different. Moreover, auroral fluctuations in UV-auroral emission brightness can be significant even over relatively short periods of time. Similarly, the Gustin et al. (2016) energy measurements cannot be directly compared with our results since they used observations obtained several years after PJ32.

Furthermore, we compared the average energies obtained in our study with the in situ measurements conducted by Mauk et al. (2020) using the JEDI instrument during PJs 4, 6, 7, and 10 of the Juno mission. However, the JEDI measurements by Mauk et al. (2020) for each PJ cover only a minimal fraction of the polar regions associated with the intersection of Juno’s footprint and the auroral emission oval in the northern and southern auroral regions. We evaluated that the average energies of electrons precipitating in these auroral regions, measured by JEDI at Juno’s altitudes, range between 150 and 300 keV. In comparison with the results of our study, where we consider kappa distributions for precipitating electrons, we observe energies distributed between 100 and 200 keV in approximately 70% of the auroral region (either north or south), with peaks reaching up to 600 keV in the remainder of these zones. Thus, the mean energies derived at auroral altitudes are generally of the same order as the results obtained by Mauk et al. (2020). However, at this stage, it is premature to draw a definitive conclusion regarding the comparison between our results and the measurements of electron energy distributions conducted by Mauk et al. (2020). The auroral emission maps used to derive the mean energies of precipitating electrons are integrated over several hours, whereas JEDI measurements are almost instantaneous. Moreover, given the rapid dynamics of polar regions, including phenomena like short-lived bright flares and local temporal variations, it is challenging to go beyond a comparison of orders of magnitude. This comparison thus requires a more in-depth analysis, which will be the subject of further study.

In the present study, and in previous similar studies (e.g., Trafton et al. 1994, 1998; Gustin et al. 2002, 2016; Ajello et al. 2005; Gérard et al. 2014), the CR(〈E〉) relationship is modeled by considering a 1D atmosphere model assuming a constant homopause throughout the auroral regions, which is a very approximate hypothesis. Hydrocarbon abundances in auroral regions are expected to be influenced by the precipitation of magnetospheric charged particles. Recent observations of these abundances (e.g., Sinclair et al. 2018) demonstrate that the spatial distribution of the main hydrocarbons in the auroral region is inhomogeneous. The maximum molar fraction variability factor in the auroral region is around 1.2 for C₂H₂, 1.1 for C₂H₆, and up to 1.3 for C₂H₄ at 0.01 mbar pressure only in the north pole (see Sinclair et al. 2018). For the southern auroral region the variability is also significant, and differs from that of the north. As CH₄ photolysis is one of the sources of production of these different hydrocarbons, this suggests that the altitude of its homopause is also variable in the auroral regions. This may impact the estimates of the average energy of electrons precipitating in these regions. To evaluate this impact, we modeled the CR(〈E〉) relationship using two different CH₄ abundance profiles based on the A and C eddy diffusion models of Moses et al. (2005) and Hue et al. (2018). These profiles are not representative of CH₄ auroral abundance, but have a higher homopause compared with the atmospheric model of Grodent et al. (2001, see panel a in Fig. D.4). Using this approach, the results obtained (Fig. D.4 shown in Appendix D.3) suggest that, depending on the CH₄ abundance profile used in our atmospheric model, the CR(〈E〉) relationship increases more rapidly or less rapidly as a function of 〈E〉, as shown in Fig. D.4. The consequences of this variability for the average energy map determination may also be significant. However, this issue is beyond the scope of this article, and is left to a future investigation.

Fig. 9 Maps of the average energy 〈E〉 defined in Eq. (19) and given by the formula $\langle E\rangle =2E_0\frac{\kappa }{\kappa -2}$. These maps are calculated for the north (left panel) and south (right panel) poles directly from the maps in Fig. 8. The iso-energy lines are defined as for Figs. 7 and 8 and with the same color bar normalization.

6 Conclusions

The present study is a further step forward in the investigation of electron energies precipitating in Jupiter’s auroral regions. Inspired by Dols et al. (2000), Gustin et al. (2002) and Menager (2011), we developed a new model of UV emission from H₂ in Jupiter’s auroral regions and adapted it to Juno/UVS observations. This model is more complete than previous studies. It takes into account nine electronic states of the H₂ molecule, and we combined it with the TransPlanet electronic transport model (Stamnes & Rees 1983b; Simon et al. 2007; Gronoff 2009; Menager 2011; Benmahi 2022). The H₂ auroral emission model was validated by simulating UV emission spectra at very high spectral resolution and comparing it with the results of Liu et al. (1995). The modeled UV spectra between 125 and 170 nm are compared with some of the emission spectra observed by UVS. The fit results are in good agreement with the observations (see Figs. B.1 and B.2) except in the wavelength range between 140 nm and 150 nm where the modeled spectra are more intense because we only take into account absorption by CH₄ and C₂H₂. Between 150 and 153 nm, absorption by C₂H₂ is not strong enough for a good fit of the UVS spectra shown in Figs. B.1 and B.2. This suggests that the C₂H₂ abundance profile we used in this study is underestimated in Jupiter’s auroral regions. We did not include absorption by C₂H₆ because we do not have an auroral abundance profile for this chemical species. However, as demonstrated in Gustin et al. (2016), C₂H₆ absorbs in the interval [140 nm–150 nm] and influences the amplitude of the spectrum when electrons reach the homopause, typically for average energies above 20 keV.

Thanks to Juno/UVS observations during PJ32, we mapped the CR at the north and south polar regions and used them to map the energy of precipitating primary electrons. The relationship CR(E₀,θ) was modeled taking into account the emission angle at each observed point in both auroral regions. We used the atmosphere model of Grodent et al. (2001) and considered a constant homopause throughout the auroral regions.

In the northern and southern auroral regions, using the JRM33 magnetic field model (Connerney et al. 2022), we found that the magnetic dip angle varies between ψ ~ 60° and ψ ~ 80°. As the penetration depth of electrons precipitating in these regions is influenced by the magnetic dip angle, this has an impact on the CR. We modeled the CR(E₀) relationship for small variations in ψ and found a small influence. This allowed us to consider only the median value of ψ in each auroral region for our modeling, in order to distinguish between the modeling of the CR(E₀) relationship in the north and south.

Modeling of the CR(E₀, θ) relationship was carried out by considering mono-energetic initial electron flux distribution and kappa distribution (κ = 2.5). This allowed us to compare results from previous works to those obtained with a more realistic broadband population. We found that when considering a mono-energetic distribution, the average energy of electrons precipitating in auroral regions is globally underestimated by a factor of 3–4. We also found that at low energies (below ~ 100 keV), the CR(〈E〉) relationships intersect (see panel c in Fig. D.4). Thus, below 〈E〉 = 90 keV, the average electron energy is overestimated if the mono-energetic hypothesis is used to infer the mean energy of a broadband population. Above 90 keV, the energy of electrons precipitating into auroral regions is underestimated.

Our results clearly demonstrate the importance of considering broadband distributions representative of the actual particle observations (e.g., Salveter et al. 2022) and modeled here as a kappa distribution when modeling the CR relationship as a function of the average energy of electrons precipitating in auroral regions. The average energies inferred by this method under the above-mentioned hypotheses for the atmospheric composition profile lie in the 300–500 keV range in the polar emission region of the north. In the main emission zones we found average energies up to 550 keV, with peaks along the auroral oval. In the outer emission regions the average energies lie between 5 and 50 keV. In the south the polar emissions are much fainter (see also Greathouse et al. 2021) with a mean energy peak of about 100 keV. In the main emission regions we found several average energy peaks from 150 keV to 600 keV.

Finally, thanks to this study, our work can be readily applied to mapping the average energy of auroral electrons for all Juno mission PJs. This will allow us to establish a temporal map of the electron energies precipitating in these regions. We also aim to improve this study by taking into account the meridional and latitudinal variabilities of the CH₄ homopause.

Acknowledgements

This work was supported by the University of Liege under Special Funds for Research, IPD-STEMA Programme. C.S.W. thanks the Austrian Science Fund (FWF) and project P35954-N. V.H. acknowledges support from the French government under the France 2030 investment plan, as part of the Initiative d’Excellence d’Aix-Marseille Université – A*MIDEX AMX-22-CPJ-04.

Appendix A Description of the electron transport model (TransPlanet)

In the context of the modeling work we performed for this study, we focused on the kinetic part of the code. This part calculates the ionization rates, atmospheric particle excited states, and emission rates caused by solar ultraviolet flux and the precipitation of magnetospheric electrons. The energy deposited by electrons when they interact with atmospheric particles is often modeled using two different approaches: the continuous loss approximation (Yung et al. 1982; Gérard & Singh 1982; Waite et al. 1983; Singhal et al. 1992) and the discrete loss approximation (Kim et al. 1992; Perry et al. 1999). Contrary to other models that use these types of approximation, the Trans* code describes the interaction between suprathermal electrons and neutral atmospheric particles by self-consistently solving the dissipative Boltzmann equation. Likewise, the originality of the Trans-* solvers is that they are based on a radiative transfer solver, called DISORT (Stamnes et al. 1988). This is in sharp contrast with Monte Carlo transport models, and results in an accrued computation speed.

Boltzmann’s equation describes the interactions between particles in a gas. This equation models binary collisions between solid spheres interacting at short distance. It also takes into account the discrete aspect of the energy loss that occurs with each collision. The secondary electrons produced during ionization by collisions between precipitating primary electrons and atmospheric particles are also taken into account using this equation, and are included in the suprathermal flux. All these interactions are therefore governed by the elastic and inelastic scattering cross sections of the electrons, thanks to the Boltzmann equation, extensively described and detailed by Stamnes & Rees (1983b), Gronoff (2009), Menager (2011), and Benmahi (2022).

Electrons precipitating into the atmosphere are represented by the distribution function f(r, v, t) in phase space, where r, v, and t represent spatial position, velocity, and time, respectively, and f is given in cm⁻⁶s². Thus, the evolution of the distribution f is given by the nonconservative Boltzmann equation $$\frac{\partial f}{\partial t}+v\cdot \frac{\partial f}{\partial r}+\frac{\partial }{\partial v}\left(\frac{fX}{m_e}\right)=Q,$$(A.1)

where Q represents a source function describing the electrons produced at position r, velocity v, and time t; m_e is the electron mass; and X is an external force applied to electrons in state (r, v). In this electron transport model, the X function is described by $$X=-n_{e}L(E)\frac{v}{v},$$(A.2)

where the function L(E) describes the interaction by friction of thermalized electrons with suprathermal electrons, n_e is the number density of thermalized electrons, E is the energy of a suprathermal electron, and v = ||v|| is its velocity. Thus, the nonconservative nature of the Boltzmann equation arises from the nonconservative force X. Trans* codes use the continuous friction function L(E) proposed by Swartz et al. (1971) and established by Schunk et al. (1971) and Schunk & Hays (1971) to describe Coulomb interactions and Cerenkov emission. This function was therefore recommended by Stamnes & Rees (1983a) for the Trans* models.

In order to solve this equation, it is useful to reduce it to an equation relating to the flux I (which is given in cm⁻²s⁻¹ eV⁻¹sr⁻¹) by replacing the suprathermal electron distribution function f by the variable change $I(r,E,u,t)=\frac{v^2}{m_e}f(r,v,t)$, where $E=\frac{1}{2}m{v}^2$ mv² is the kinetic energy of the electrons and $u=\frac{v}{v}$ is their direction. Thus, the equation becomes $$\frac{1}{v}\frac{\partial I}{\partial t}+\frac{v}{v}\cdot \frac{\partial I}{\partial r}-n_e\frac{\partial }{\partial E}(L(E)I)=\frac{v^2}{m_e}Q.$$(A.3)

In the case of a plane-parallel geometry and assuming a stationary state, equation A.3, which represents the flux I along a magnetic field line B, becomes $$\mu \frac{\partial I(\tau ,\mu ,E)}{\partial \tau (z,E)}=-I(\tau ,\mu ,E)+\frac{n_e(z)}{{\displaystyle {\sum }_k{n}_k}(z){\sigma }_k^{tot}(E)}\frac{\partial }{\partial E}(L(E)I(\tau ,\mu ,E))+D(z,\mu ,E)+P(z,\mu ,E),$$(A.4)

where µ is the cosine of the angle between the magnetic field line and the direction of electron propagation, n_k(z) [cm⁻³] is the concentration of the atmospheric species k at altitude z, ${\sigma }_k^{tot}(E)$ [cm²] is the total collision cross section between an electron and the species k at energy E, ${\sigma }_k^{tot}(E)$ is also the sum of the elastic and inelastic collision cross sections, and τ is a dimensionless quantity representing the electron scattering depth (similar to the optical depth in radiative transfer) defined by $\tau (z)={\displaystyle {\int }_z^{z_{\max }}{\displaystyle {\sum }_k{n}_k}}(z){\sigma }_k^{tot}(E)\frac{dz}{\mu }$.P (z, µ,E) [cm⁻²s⁻¹eV⁻¹sr⁻¹] is a source term for the primary electron flux introduced into the atmosphere in the (µ, E) state at altitude z; this term includes incident magnetospheric electrons as well as photoelectrons produced by ionization caused by solar UV using a Beer-Lambert law of radiation absorption in the atmosphere. D(z,µ,E) [cm⁻²s⁻¹eV⁻¹ sr⁻¹ ] is a scattering term representing secondary electrons resulting from inelastic electron collisions between a primary electron and an atmospheric particle, as well as primary electrons whose energy has been dissipated by their interactions with atmospheric particles.

Appendix B Examples of UVS spectra fitting

From the PJ32 spectral cube, we selected two small auroral zones in the main emission region between 50°N and 60°N and in the outer emission region at around 70°N and 125°W. In these two regions, the median emission angles are respectively about 20^° and 45^°. In Figures B.1 and B.2, we plot in red the UV emission spectra averaged over each of the selected regions. The measured CR ratio is around 5.5 in the main emission region and 3.7 in the outer emission region.

Using an initial mono-energetic electron flux distribution, we modeled the emission spectrum to obtain the same CR as the observed spectrum for each selected region. For the main emission region, we used a mono-energetic distribution with an average energy of 90 keV to reproduce the emission spectrum in blue (Fig. B.1) with a CR of around 5.51, which represents the best fit. For the outer emission region, we obtained the best fit by modeling the UV emission spectrum using a mono-energetic distribution with an average energy of 75 keV and a CR of around 3.52.

However, depending on the spectral resolution of an observation, the emission spectrum can be fitted using any initial flux distribution in the electron transport model. On the other hand, at very high spectral resolution (typically ~0.0125 nm), the energy signature of electrons on spectral lines below 110 nm is visible, and the choice of initial electron flux distribution will result in drastically different emission spectra. This implies that, in this case, there are extra free parameters that we are unable to control in order to constrain the type of initial electron flux distribution in auroral zones. Thus, only in situ measurements of electron energy spectra can provide an answer to this problem.

Fig. B.1 Examples of spectral fit (in blue) of a spectrum observed by Juno/UVS (in red) in a small region of the main emission in the southern arc of the northern auroral oval between 60°N and 50°N latitude. The median emission angle is around 20°, the observed CR is 5.50 and the spectral resolution is around 2.3 nm.

Fig. B.2 Examples of spectral fit of a spectrum observed by Juno/UVS in a small region of the outer emission in the northern auroral oval around position 70°N and 125°W SIII. The median emission angle is around 45°, the observed CR is 3.73 and the spectral resolution is estimated at around 2.3 nm.

Appendix C Regions of auroral emission and viewing angle maps

In Fig. C.1 the magenta dots represent the pixels selected by our S/N criterion, and thus represent UV emission from the auroral region at the north pole. In the same way as above, we also used this selection criterion to isolate the auroral emission from the southern region.

Figure C.2 displays maps of the emission angles observed by UVS in the north and south polar regions. These emission angles range from 0° to 50° in the south pole, and up to 80° in the north pole.

Fig. C.1 Isolated UV emission from Jupiter’s auroral region, observed during PJ32 at the north pole by Juno/UVS. The magenta dots represent UV emission spectra where the unabsorbed part of the spectrum has a S/N ≥ 3.

Fig. C.2 Maps of emission angles in the polar regions (left: northern hemisphere; right: southern hemisphere) observed during PJ32 when Juno was flying over the polar regions. Each point was observed with different emission angles along the probe’s trajectory. Thus, in these maps, each pixel represents the median value of the emission angles of the observed point.

Appendix D Modeling of CR relationship

D.1 The dip magnetic angle and its impact on CR modeling

The modeled CR as a function of the initial energy distribution of the precipitating electrons is shown in Fig. D.1 for a fixed emisson angle θ = 0°. These results assume a mono-energetic initial electron flux distribution (Fig. D.1 left panel) and a kappa distribution (Fig. D.1 right panel). In the case of a mono-energetic distribution, the CR is modeled using 14 characteristic energy points ranging from 1 keV to 220 keV. In the case of a kappa distribution, we modeled the CR with only ten characteristic energy points ranging from 1 keV to 85 keV. We did this because the CR increases more quickly than in the mono-energetic case due to the broadening of the kappa distribution at high energies. The total precipitation flux for both distributions was set to Q₀ = 1 erg cm⁻²s⁻¹ for all electron transport simulations.

Electrons precipitating into the atmosphere are guided by magnetic field lines. Depending on the magnetic dip angle ψ of a given field line, these electrons will penetrate more deeply or less deeply in the atmosphere. By varying ψ for a given characteristic energy, the CR also varies. These variations in CR as a function of angle ψ are less than 0.5 for ψ e [60°, 75°].

In the present study we modeled CR(E₀ ,θ) as a function of the median magnetic dip angle ψ in the auroral region for each hemisphere. To this end, we used the JRM33 magnetic field model of Jupiter (Connerney et al. 2022) to calculate the magnetic dip angle at each point of the selected northern and southern auroral regions (see the example in Fig. C.1). For each electron energy distribution, the CR is modeled at the north and south poles separately, and is represented by red and blue dots, respectively (see Fig. D.1). For the north pole ψ = 65.7°, and for the south pole ψ = 74.4°.

Fig. D.1 Example of modeled CR(E0) relationship. In the left panel the modeled relationship CR(E0) corresponds to the case of a mono-energetic initial electron flux distribution. Similarly, in the right panel the relationship CR(E0) corresponds to the case of a kappa initial distribution as a function of caracteristic energy E0. The red and blue dots represent the modeled CR(E0) relationship, respectively, for the north pole with ψ = 65.7° and for the south pole with ψ = 74.4°. The dotted red and blue lines represent the best fit for the north pole and south pole, respectively, using equation D.1.

Fig. D.2 3D representation of the CR(E0, θ) modeled relationship for the case of a mono-energetic distribution (left panel) and for the case of a kappa distribution (right panel). Both relationships are calculated for the north pole (i.e., for ψ = 65.7°). The green grid lines represent the fit of the modeled CR(E0, θ) relationship. The colored surface (color bar from blue to red) represents the absolute difference between the modeled surface and the fit surface, with an average CR uncertainty of around 0.6 for the case of the mono-energetic distribution and an uncertainty of around 0.5 for the case of the kappa distribution.

For the case of mono-energetic electron flux distribution, these results demonstrate the importance of taking into account the geometry of the magnetic field lines at high energies because, typically above 120 keV (see Fig. D.1), the relationship CR(E₀) becomes dependent on the magnetic dip angle. For the kappa distribution case, the evolution of the CR(E₀) relationship is less influenced by the magnetic dip angle.

D.2 Analytic formula and fit of the CR(E₀, θ) relationship

Following the study of Gustin et al. (2016), we considered that CR(E0) follows a hyperbolic law at low energy and increases as a logarithmic law at high energy. Using our result from Fig. D.1, we derive a phenomenological relation CR(E₀) in order to fit the CR modeling as a function of characteristic energy. This step allows us to obtain an analytical form of CR as a function of energy. In addition, thanks to a MCMC fit, it makes it possible to estimate the uncertainty on the modeled CR. Our analytical formulation of CR(E₀) for a fixed θ emission angle is given by $$\begin{array}{ll}\hfill & \text{CR}\left(E_0\right)=\hfill \\ \hfill & A\cdot C\cdot \left(\tanh \left(\frac{E_0-E_c}{B}+1\right)\right)\cdot \ln {\left({\left(\frac{E_0}{D}\right)}^{\alpha }+e\right)}^{\beta },\hfill \end{array}$$(D.1)

where A is the minimum amplitude of the modeled CR; E_c is a threshold energy; and B, C, D, α, and β are fit parameters that constrain the shape of the curve throughout the energy range.

To adjust these fitting parameters, we used the Python emcee package developed by Foreman-Mackey et al. (2013), which implements the MCMC method using the Metropolis-Hastings algorithm. The fitting configuration is characterized by 250 Markov chains and 2500 iterations. These two parameters were determined after several runs of the burn-in size determination. We found that the Markov chains converge after 500–1000 iterations on average for all fit parameters for the mono-energetic and kappa distributions. The choice of 2500 iterations ensures convergence of the Markov chains in all cases.

In Fig. D.1, each CR-energy relationship that we modeled (red dots for the north and blue dots for the south) was fitted using the formula D.1. The magenta and cyan colored envelopes represent the 1σ confidence band for each fit. This allowed us to estimate the uncertainty of the CR modeling as a function of energy. We found a mean CR uncertainty of around 0.5 for the case of a mono-energetic distribution and 0.2 for the case of a kappa distribution.

For a variable emission angle θ and a fixed characteristic energy, we found that the CR(θ) relationship follows a sinusoidal law. By taking into account the variability of the emission angle and the characteristic energy simultaneously, the CR(E₀, θ) relationship is two-dimensional and is given by $$\begin{array}{ll}\hfill & \text{CR}\left(E_0,\theta \right)=\hfill \\ \hfill & A\cdot C\cdot \left(\tanh \left(\frac{E_0-E_c}{B}+1\right)\right)\cdot \ln {\left({\left(\frac{E_0}{D}\right)}^{\alpha }+e\right)}^{\beta }\hfill \\ \hfill & \left(1+\delta \cdot \sin {(\theta )}^{\gamma }\right)\hfill \end{array}$$(D.2)

where δ and γ are additional fit parameters.

As in the one-dimensional case (i.e., CR(E₀)), we modeled the CR(E₀, θ) relationship in a 2D (E₀, θ) grid. This grid is defined by 14 characteristic energy points ranging from 1 keV to 220 keV for the case of a mono-energetic distribution, and 10 characteristic energy points ranging from 1 keV to 85 keV for the case of a kappa distribution. For each characteristic energy E₀, we also considered an emission angle grid θ of ten points ranging from 0° to 80°. This corresponds to modeling 140 UV emission spectra, for each auroral region, to obtain a map of the CR(E₀, θ) relationship for the case of mono-energetic distribution. For the case of kappa distribution, this corresponds to modeling 100 UV emission spectra, for each auroral region, to obtain the maps of the CR(E₀, θ) relationship.

The electron transport modeling accuracy at each of the (E₀, θ) grid points (illustrated above) is evaluated by a total energy conservation rate of the electrons precipitated in the model. Thus, for each electron transport simulation, the average rate of conservation of the total energy precipitated is about 99% for the case of a mono-energetic distribution and 99.5% for the case of a kappa distribution, which represents an energy loss, respectively, of 1% and 0.5%. All parameters were fitted in the same way as for the 1D case, using the same burn-in size for the MCMC method.

Fig. D.3 3D representation of the CR(E0, θ) modeled relationship for the case of mono-energetic distribution (left panel) and kappa distribution, south pole (i.e., for ψ = 65.7°) and with the same conventions as Fig. D.2. The modeled CR uncertainty is estimated respectively at around 0.6 and 0.5 for both cases.

In Figures D.2 and D.3, we represent in yellow the results of the CR(E₀, θ) modeled relationship. To invert this phenomenological relationship, we used the formula D.2 and fit its parameters using the MCMC method explained in the text. In Table D.1, we presented the fit values of the parameters of the formula D.2 for each of the cases considered in our study. These parameters can be used directly in Eq. D.2 to map the characteristic energy of electrons precipitating in Jupiter’s auroral regions for any PJ. Δ_CR is the fitted noise of the modeled CR and corresponds to the uncertainty of CR(E₀, θ). The average absolute difference between the fit and the model result is comparable to the uncertainty on the CR obtained by the MCMC method.

D.3 Evolution of the CR(E₀,θ) relationship using different CH₄ abundance profiles

To evaluate the impact of methane distribution on the CR, we modeled the CR(〈E〉) relationship for an emission angle θ = 0° and using two different CH₄ abundance profiles from the A and C eddy diffusion models of Moses et al. (2005) and Hue et al. (2018). These profiles are not representative of CH₄ auroral abundance, but have different homopauses compared to that of the Grodent et al. (2001) atmospheric model (see panel a in Fig. D.4). The CR(〈E〉) relationship was modeled for the cases of a kappa distribution and a mono-energetic distribution for electrons precipitating in auroral regions.

The results obtained are shown in panels b and c in Fig. D.4. This shows us that the CR depends on the CH₄ abundance profile used in our atmospheric model. Thus, depending on the CH₄ homopause, the CR(〈E〉) relationship increases as a function of 〈E〉 more rapidly or less rapidly, as shown in Fig. D.4. The impact of this variability on the determination of the average energy map is also significant.

The results shown in panel c in Fig. D.4 show us the differences in the CR(〈E〉) relationships obtained for the case of a kappa distribution and the case of a mono-energetic distribution of precipitating electrons. Thus, using the Grodent et al. (2001) atmospheric model, we found that below 90 keV the CR(〈E〉) relationship obtained for the case of a mono-energetic distribution overestimates the average energy of precipitating electrons. Above 90 keV, the same relationship underestimates the average energy of electrons precipitating in auroral regions.

By analyzing the ratio of the maps of the mean energies derived from two different electron distributions, a kappa distribution and a mono-energetic one, we obtain the results presented in Fig. D.5. This result clarifies Fig. D.4 (panel c). It is notable that, in the main emission and polar emission zones, the mean energies derived from the kappa distribution are on average 3 to 5 times higher than those derived from the mono-energetic distribution for precipitating electrons in the northern and southern auroral regions. In Fig. D.5, in the auroral regions where the ratio of mean energies approaches ~1, the mean energies derived from the kappa distribution are around ~90 keV. In areas where the ratio is less than 1, the mean energies associated with the kappa distribution are strictly less than 90 keV.

Fig. D.4 Model of the CR(〈E〉), θ) relationship for θ = 0° using different CH4 abundance profiles. In panel a, the red, blue, and green solid lines represent the CH4 abundance profile from, respectively, the A and C eddy diffusion models of Moses et al. (2005) and Hue et al. (2018), and the Grodent et al. (2001) atmospheric model. In panel b is represented the modeled CR(〈E〉) relationship. The solid lines represent the CR modeled using kappa electron flux distribution and the dotted lines by using a mono-energetic flux distribution. The color conventions are the same as described for panel a. In panel c a zoomed-in version of panel b is shown to distinguish between the various CR(〈E〉) results at low energy.

Fig. D.5 Ratio of average energy map for the kappa-distribution case (Fig. 9) to that of the monoenergetic case (Fig. 7). Left: Northern auroral region. Right: Southern auroral region.

References

All Tables

All Figures

	Fig. 1 Atmospheric model described by Grodent et al. (2001), which considers only the neutral compounds (H, H2, He, CH4, and C2H2) that predominate in Jupiter’s atmosphere. For electronic transport modeling, only H, H2, He, and CH4 compounds were considered. For the UV emission model, spectral absorption by CH4 and C2H2 was taken into account.
In the text

Fig. 2 Examples of comparisons of modeled synthetic spectra with experimental spectra (red line) obtained by Liu et al. (1995) with a spectral resolution of 0.0125 nm. The blue line and the dashed green line represent the best fits obtained by precipitating a 100 eV mono-energetic (best fit A) and Maxwellian (best fit B) electron flux distributions. The four spectral windows shown have a width of 1 nm. In the ii window around 121.567 nm, the Lyman-α line was filtered out.

In the text

Fig. 3 Integrated non-absorbed (N.A.) UV emission from Jupiter’s auroral regions observed by Juno/UVS during PJ32 in the SIII jovicentric reference frame. The acquisition time over the northern polar region is 4282 s. At the southern polar region the acquisition time is 14 561 s. The plus sign (+) in red represents the average solar longitude during the selected acquisition times for the northern and southern polar regions individually. Thus, the Sun’s longitudinal path between its mean position during acquisition at the northern hemisphere and the southern hemisphere is approximately 130° westward.

In the text

	Fig. 4 Color ratio of Jupiter’s auroral regions (left: northern hemisphere; right: southern hemisphere) observed during PJ32. The panels show the CR calculated for each pixel of the UV emission map as defined by Gustin et al. (2016) with $\text{CR}=\frac{I(155 \text{nm}-162 \text{nm})}{I(125 \text{nm}-130 \text{nm})}$.
In the text

	Fig. 5 Examples of kappa distributions with different characteristic energies. The three examples use κ = 2.5 and Q0 = 10 erg cm–2 s–1, and three different characteristic energies: E0 = 1 keV, E0 = 10 keV, and E0 = 65 keV.
In the text

Fig. 6 Optical depth calculated over the atmospheric column for CH4 and C2H2. The transparent green and cyan bands represent the absorption spectral ranges used for the CR calculations. For CH4 the absorption spectral range is considered between 125 nm and 130 nm and for C2H2 between 150 and 153 nm. The transparent magenta band is the non-absorbed (N.A.) spectral range over which hydrocarbon absorption was assumed to be negligible.

In the text

Fig. 7 Characteristic energy maps E0 obtained from the CR(E0, θ) relationships (see Figs. D.2 and D.3) modeled for the case of an initial mono-energetic electron flux distribution, and from the CR observed during PJ32, at the north (left panel) and south (right panel) poles. Iso-energy lines are defined for values of 1, 5, 10 keV, in steps of 20 keV between 10 and 300 keV, and then in steps of 100 keV between 300 and 900 keV.

In the text

	Fig. 8 Characteristic energy maps E0 obtained from the CR(E0, θ) relations (see Figs. D.2 and D.3) modeled for the case of an initial kappa electron flux distribution, and from the CR observed during PJ32, at the north (left panel) and south (right panel) poles. Iso-energy lines are defined in the same way as in Fig. 7.
In the text

	Fig. 9 Maps of the average energy 〈E〉 defined in Eq. (19) and given by the formula $\langle E\rangle =2E_0\frac{\kappa }{\kappa -2}$. These maps are calculated for the north (left panel) and south (right panel) poles directly from the maps in Fig. 8. The iso-energy lines are defined as for Figs. 7 and 8 and with the same color bar normalization.
In the text

	Fig. B.1 Examples of spectral fit (in blue) of a spectrum observed by Juno/UVS (in red) in a small region of the main emission in the southern arc of the northern auroral oval between 60°N and 50°N latitude. The median emission angle is around 20°, the observed CR is 5.50 and the spectral resolution is around 2.3 nm.
In the text

	Fig. B.2 Examples of spectral fit of a spectrum observed by Juno/UVS in a small region of the outer emission in the northern auroral oval around position 70°N and 125°W SIII. The median emission angle is around 45°, the observed CR is 3.73 and the spectral resolution is estimated at around 2.3 nm.
In the text

	Fig. C.1 Isolated UV emission from Jupiter’s auroral region, observed during PJ32 at the north pole by Juno/UVS. The magenta dots represent UV emission spectra where the unabsorbed part of the spectrum has a S/N ≥ 3.
In the text

	Fig. C.2 Maps of emission angles in the polar regions (left: northern hemisphere; right: southern hemisphere) observed during PJ32 when Juno was flying over the polar regions. Each point was observed with different emission angles along the probe’s trajectory. Thus, in these maps, each pixel represents the median value of the emission angles of the observed point.
In the text

Fig. D.1 Example of modeled CR(E0) relationship. In the left panel the modeled relationship CR(E0) corresponds to the case of a mono-energetic initial electron flux distribution. Similarly, in the right panel the relationship CR(E0) corresponds to the case of a kappa initial distribution as a function of caracteristic energy E0. The red and blue dots represent the modeled CR(E0) relationship, respectively, for the north pole with ψ = 65.7° and for the south pole with ψ = 74.4°. The dotted red and blue lines represent the best fit for the north pole and south pole, respectively, using equation D.1.

In the text

Fig. D.2 3D representation of the CR(E0, θ) modeled relationship for the case of a mono-energetic distribution (left panel) and for the case of a kappa distribution (right panel). Both relationships are calculated for the north pole (i.e., for ψ = 65.7°). The green grid lines represent the fit of the modeled CR(E0, θ) relationship. The colored surface (color bar from blue to red) represents the absolute difference between the modeled surface and the fit surface, with an average CR uncertainty of around 0.6 for the case of the mono-energetic distribution and an uncertainty of around 0.5 for the case of the kappa distribution.

In the text

	Fig. D.3 3D representation of the CR(E0, θ) modeled relationship for the case of mono-energetic distribution (left panel) and kappa distribution, south pole (i.e., for ψ = 65.7°) and with the same conventions as Fig. D.2. The modeled CR uncertainty is estimated respectively at around 0.6 and 0.5 for both cases.
In the text

Fig. D.4 Model of the CR(〈E〉), θ) relationship for θ = 0° using different CH4 abundance profiles. In panel a, the red, blue, and green solid lines represent the CH4 abundance profile from, respectively, the A and C eddy diffusion models of Moses et al. (2005) and Hue et al. (2018), and the Grodent et al. (2001) atmospheric model. In panel b is represented the modeled CR(〈E〉) relationship. The solid lines represent the CR modeled using kappa electron flux distribution and the dotted lines by using a mono-energetic flux distribution. The color conventions are the same as described for panel a. In panel c a zoomed-in version of panel b is shown to distinguish between the various CR(〈E〉) results at low energy.

In the text

	Fig. D.5 Ratio of average energy map for the kappa-distribution case (Fig. 9) to that of the monoenergetic case (Fig. 7). Left: Northern auroral region. Right: Southern auroral region.
In the text