Estimating Precipitating Energy Flux, Average Energy, and Hall Auroral Conductance From THEMIS All-Sky-Imagers With Focus on Mesoscales

Recent attention has been given to mesoscale phenomena across geospace (∼10 s km to 500 km in the ionosphere or ∼0.5 RE to several RE in the magnetosphere), as their contributions to the system global response are important yet remain uncharacterized mostly due to limitations in data resolution and coverage as well as in computational power. As data and models improve, it becomes increasingly valuable to advance understanding of the role of mesoscale phenomena contributions—specifically, in magnetosphere-ionosphere coupling. This paper describes a new method that utilizes the 2D array of Time History of Events and Macroscale Interactions during Substorms (THEMIS) white-light all-sky-imagers (ASI), in conjunction with meridian scanning photometers, to estimate the auroral scale sizes of intense precipitating energy fluxes and the associated Hall conductances. As an example of the technique, we investigated the role of precipitated energy flux and average energy on mesoscales as contrasted to large-scales for two back-to-back substorms, finding that mesoscale aurora contributes up to ∼80% (∼60%) of the total energy flux immediately after onset during the early expansion phase of the first (second) substorm, and continues to contribute ∼30–55% throughout the remainder of the substorm. The average energy estimated from the ASI mosaic field of view also peaked during the initial expansion phase. Using the measured energy flux and tables produced from the Boltzmann Three Constituent (B3C) auroral transport code (Strickland et al., 1976; 1993), we also estimated the 2D Hall conductance and compared it to Poker Flat Incoherent Scatter Radar conductance values, finding good agreement for both discrete and diffuse aurora.


INTRODUCTION
A fundamental question at the heart of the dynamics and coupling of Earth's magnetosphere-ionosphere system is, "How is energy deposited into the ionosphere from the magnetosphere?" That itself has a complex answer comprised of multiple parts such as electromagnetic energy, particle precipitation, and waves. Great strides have been made in characterizing the global energy input due to precipitation; for example, Newell et al. [1] built an empirical precipitation model using Defense Meteorological Satellite Program (DMSP) statistics to determine the auroral power from diffuse, monoenergetic discrete, and broadband discrete aurora as a function of solar wind coupling. Results were binned by magnetic local time (MLT) and latitude. Their work has been useful for describing auroral power at large spatial scales, but does not explain contribution of specific scale sizes or provide event-based numbers that could be used in global climate models to investigate the ionosphere-thermosphere response to an intense, mesoscale (∼10 s-500 km) auroral form. DMSP statistics-based modeling has been improved over the years, leading to the construction of the Oval Variation, Assessment, Tracking, Intensity, and Online Nowcasting (OVATION) models [2][3][4] to the current state-of-the art version, OVATION Prime, all of which still excel at the more global precipitation response without capturing the spatially and temporally dynamic aurora that occur on mesoscales.
In a newly developed empirical model, Auroral Spectrum and High-Latitude Electric field and variability (ASHLEY) [5], there is no assumed spectra distribution for the precipitating particles, and the DMSP observed energetic particles in 19-energy channels have been directly binned according to the geomagnetic local time, latitude, and solar wind conditions. ASHLEY has shown some advantage to specify the soft electron precipitation, but is still limited to the large-scale phenomena. Another very recent DMSP statistics-based modeling work employed machine learning techniques to build a new total electron energy flux particle precipitation nowcast model called PrecipNet [6]. PrecipNet better captures the dynamic changes of the auroral flux than OVATION Prime, improving capabilities to reconstruct mesoscale phenomena and making great progress towards nowcasting and forecasting auroral response to solar wind driving.
Though these and other works make great strides towards resolving the role of mesoscale auroral precipitation, the field still lacks a rigorous understanding and characterization of the contribution of the dynamic and energetic mesoscale auroral forms that are embedded within the large-scale system. These are important to resolve, as studies such as Deng et al. [7] and Sheng et al. [8] have illustrated with the Global Ionosphere-Thermosphere model (GITM) how mesoscale phenomena in the ionosphere can vastly alter the ionospheric and mesospheric response.
A major reason for the continued deficiency in mesoscale characterization has been inadequate data sources to fully resolve spatial and temporal ambiguities. The DMSP examples above are perhaps one of the best-suited datasets because they provide in situ measurements of the particle flux and energy, as well as 2D views of the aurora from above as the satellite flies overhead. However, statistics must be employed, since orbiting spacecraft measurements lack continuous coverage over a given 2D area. It is not possible to investigate auroral evolution over a short time in the region of interest. DMSP orbits are also not often optimized for covering the nightside auroral oval, as most of their measurements are made near the dawn-dusk meridian.
In contrast, ground-based observations consistently observe the same geographic location over time, allowing for the temporal changes to be recorded. However, most datasets do not have sufficiently high resolution for observing the aurora's spatial evolution with a wide spatial coverage in the nightside auroral oval. For example, meridian scanning photometers (MSPs), which observe multiple wavelengths (such as blue, green, and redline aurora), are limited to a line of constant longitude. Incoherent Scatter Radar (ISR) can measure electron densities to determine the precipitated energy, energy flux, and resulting changes in conductance in 1D or 2D (e.g., [9,10]), but the measurement area is limited to ∼100 s km in horizontal spatial extent at 100 km altitude.
All-sky-imagers (ASIs) provide auroral intensity evolution in 2D space as well as over time. ASIs have been used to resolve auroral arcs and to obtain local 2D precipitation information such as energy flux (e.g., [11][12][13][14][15][16]), Pedersen conductance (e.g., [17]), and arc occurrence and distribution information (e.g., [18]). For example, Lam et al. [17] used the THEMIS ASI at Fort Yukon to demonstrate a relationship between the Pedersen conductance (as determined from Poker Flat Incoherent Scatter Radar (PFISR) data) to white light intensity within 3 mho or 40%. They assumed that the Pedersen conductance is proportional to the square root of the white light intensity, based on Kosch et al. (1998) who argued for this relationship after making the approximations that auroral optical intensity is proportional to ion production rate and that electron production is approximately equal to the electron loss.
Since 2006, NASA's THEMIS mission has provided an array of white-light ASIs across Canada and Alaska that monitor the majority of the nightside auroral oval at a 3 s cadence and up to 1 km resolution, providing a global view at temporal and spatial resolutions required to study the aurora on mesoscales [19]. This allows us to expand upon previous demonstrations of ASI capabilities to estimate precipitated energy flux and average energy in 2D, but on a continental scale. Our recent work has shown that the THEMIS ASI array can observe evolution of precipitating energy flux over the whole size of a substorm and that it is possible to estimate contributions of mesoscale precipitation over large-scale precipitation [20].
In this paper, we describe the technique that utilizes the whitelight ASI array that was not detailed in Nishimura et al. [20], including important upgrades made through the present work, supplemented by MSPs and auroral transport code calculations, to estimate precipitating energy flux, average energy, and Hall conductance, as well as what percent of the energy flux is contributed by mesoscale auroral forms during two back-toback substorms.

Background and Overview
One method to estimate the precipitating energy flux and average energy utilizes the ratio of auroral intensities at different wavelengths (colors). When electrons precipitate into the upper atmosphere, the more energetic ones reach lower altitudes before encountering high enough atmospheric densities to thermalize, whereas the less energetic electrons are stopped and thermalize at higher altitudes. Therefore, the red aurora occurring at higher altitudes are produced by less energetic populations colliding with atomic oxygen, whereas the green and blue aurora occurring at lower altitudes are caused by more energetic populations colliding with oxygen and molecular nitrogen, respectively. Methodologies that utilize the intensity ratios between the different auroral color wavelengths have been developed to determine the integrated energy flux and the average energy of the precipitating particle distribution (e.g., [21][22][23][24][25]); for example, more red than blue light indicates a less energetic population, whereas more blue than red light indicates the opposite.
Janhunen [13] used the figures that plotted intensity ratios from Rees and Luckey [21] and Strickland et al. [22] to develop an inversion method that determines precipitating energy flux from multiwavelength ASIs. They found that the energy flux reconstructed from only the green images was very good, almost indistinguishable from the reconstructions that included the red images. When comparing the reconstructed arc energy flux with energy flux measured by the FAST satellite, they found a good agreement (12 mW m −2 vs. 10 mW m −2 , respectively). Kauristie et al. [14] used the methods laid out by Janhunen [13] to estimate the energy flux measured by the All-Sky Camera in Kilpisjärvi and compared those values against those measured by the incoherent scatter radar EISCAT, finding a correlation of r 0.72 between datasets. Partamies et al. [15] further tested the technique by comparing 100 All-Sky Camera images at Longyearbyen with EISCAT energy fluxes, finding that in 35% of all of the blue and green image inversions the relative errors were less than 50%, and in 90% of the blue and green image inversions relative errors were less than 100%. Dahlgren et al. [16] also developed an inversion technique to determine both energy flux and average energy from ASIs combined with EISCAT data. Instead of relying on Rees and Luckey [21], which assumes all secondary, etc., electrons dissipate their energy locally, they determined the average energies of the precipitating electrons from Auroral Structure and Kinetics narrow angle camera emission ratios modelled with the Southampton ion chemistry and electron transport model (described in [11,26]), using atmospheric parameters for the conditions during their case study.
Hecht et al. [24,25] also presented relationships between different auroral wavelength ratios and the average energy of the precipitation using the B3C auroral electron transport code [27,28]. Strickland et al. [22,27,28] made a detailed attempt to investigate dependencies of red line emission on key parameters, unlike Rees and Luckey [21] who did not incorporate knowledge of O abundance relative to N2, which is necessary when using the red line emission. We rely upon and slightly update the methodologies developed by Strickland et al. [22] and Hecht et al. [23][24][25], which are described in Determining Energy Flux and Average Energy From Aurora Color Ratios. These methodologies, which utilize a model atmosphere to incorporate the appropriate physics and chemistry, improve calculations of the emissions and conductance as functions of the incident energy flux and average energy as compared to Rees and Luckey [21], which Nishimura et al. [20] used. We note that our techniques additionally include more physics and chemistry than the techniques employed by Kosch et al. (1998) and Lam et al. [17], who did not calculate emissions, energy flux, or average energy, but deduced the Pedersen conductance assuming that it is proportional to the square root of the white light intensity.
In order to employ the network of THEMIS ASIs, we developed log-log linear analytical functions to convert from white light counts to red (630.0 nm), green (557.7 nm), and blue (427.8 nm) intensities using overlapping meridian scanning photometer (MSP) data to calibrate. The steps we took are detailed in Converting White Light Intensities to Color Intensities. After obtaining the ratios, we estimated the precipitating electron energy flux and average energy. Using the 2D capabilities of the ASI mosaic, we determined what percentage of the energy flux is contributed by various mesoscales (Results).
From the energy flux and average energy, we estimated the Hall conductance in 2D over time (Determining the Hall Conductance) and compared to Poker Flat Incoherent Scatter Radar (PFISR) data (Results). Because the white light can be related to the OI (557.7 nm) greenline and (427.8 nm) blueline (shown in Converting White Light Intensities to Color Intensities), it does well estimating energy flux (Q). At higher average electron energies that are associated with intense aurora (>5 keV), the Hall conductivity scales very well with Q with little dependence on average electron energy (see Figure 16.5fig165 of [29] and Supplementary Figure S1). At lower energies, there is an energy dependence; but within a factor of two, the major dependent variable is still Q. Since this paper is focused on mesoscales associated with the more intense precipitation events, the Hall conductivities become a byproduct of a measurement of Q.

THEMIS All-Sky-Imagers
The THEMIS all-sky-imager (ASI) array consists of about 20 unfiltered (white light) fish-eye cameras at any given time spanning Canada and Alaska, covering a large section of the auroral oval with up to 1 km resolution near zenith, 3 km at a 40°e levation angle, 10 km at a 20°elevation angle, and 25 km at a 10°e levation angle (see Figure 6 from [19] [19]. For this study, because ASIs do not differentiate between particle species, we assume the visible light is due to precipitating electrons. Hecht et al. [30] discuss that the electron flux Q can be derived from the blue intensity, and that the proton flux is a small fraction of the electron flux. As we will show, the white light and blueline are well-correlated, so this assumption carries through. Before we continue, we highlight the caution one needs to take with ASI data at low elevation angles. For studies that require low error in energy flux, energy, or conductance magnitudes, elevation angles above ∼70°will be more trustworthy as they look approximately up the magnetic field line local to the station. The farther from zenith in the ASI field of view (low elevation angles), the more flux tubes the line-ofsight observation cuts through. Because the image is 2D, aurora at different heights may not be observed integrally, but are spread out over several latitudes (and/or longitudes) in the image. This is because all white light is mapped to 110 km altitude and then projected onto the 2D plane. To address the question of mesoscale contribution to auroral precipitation, however, we require a larger field of view; otherwise, the ASIs do not contribute additional perspectives to what we can surmise from, say, ISRs. The uncertainty in mapping can be mitigated in cases where multiple ASI field of views (FOVs) overlap by validating and/or adjusting the mapping altitude to match aurora locations viewed from different imagers. For this study, as will be shown in Results, the arcs were well-aligned between different imagers. Lastly, a future work could be to try further calibrating ASI data to DMSP data, which does not have the low elevation angle error.

NORSTAR Meridian Scanning Photometer
The NORSTAR Meridian Scanning Photometer (MSP) array consists of four ground-based photometers designed to simultaneously measure luminosity in four key auroral bands (470.9, 486.1, 557.7, and 630.0 nm). To use the color ratio relationships discussed in Determining Energy Flux and Average Energy From Aurora Color Ratios, we converted the 470.9 nm intensity to 427.8 nm by multiplying the former by 4.5 [31]. As Omholt [32] showed, the ratio is constant between these auroral emissions. Each MSP utilizes a high speed, 8-channel filter wheel (to observe auroral plus background channels), coupled to a highly sensitive photo-multiplier tube (PMT) which mechanically scans from near horizon-to-horizon in the North-South direction. Each scan consists of approximately 544 elevation steps, taking 30 s to complete a full horizon to horizon scan in high resolution mode. PMT counts from these scans are then background subtracted and converted into units of Rayleighs (R). The green 557.7 nm, blue 470.9 nm, and 486.1 nm emissions are projected to an altitude of 110 km and the red 630.0 nm emission is typically projected to 230 km. We plotted the Fort Smith (FSMI) MSP line of sight (LOS) on top of the FSMI THEMIS ASI field of view (FOV) in Figure 1A, mapped to 110 km.

Poker Flat Incoherent Scatter Radar
In order to test how well our optical methods are able to determine the Hall conductance, we compared to the conductance determined from the Poker Flat Incoherent Scatter Radar (PFISR) altitude resolved electron density data. From this perspective the PFISR observations are considered a "truth" dataset, at least for the calculation of the Hall conductance. We plotted the PFISR location on top of the Fort Yukon (FYKN) THEMIS ASI FOV in Figure 1B. The PFISR electron density observations used in this study were estimated from received power measured using the Barker Coded pulses in the 121beam "Grid11x11_11" mode on 16 February 2010. In this mode, the radar transmitted long-pulses and Barker codes, which are optimize for F region and lower E region observations, respectively, but it did not include alternating code pulses that are optimized for the whole E-region. This 121-beam mode uses a 130 µs, 13 baud Barker Code sampled at 5 µs, which provides electron density observations every 0.75 km with a range resolution of 1.5 km and achieves 48 samples per beam in 14.157 s. The electron density data was post-integrated to achieve ∼30 s integration time. This mode does not provide a vertical beam, so the highest elevation angle beam, 89°, was used to compute conductance.
The Barker code is a pulse-compression code that only decodes properly when the correlation time of the target is longer than the time length of the pulse (e.g., [33]). ISR theory (e.g., [34]) predicts that the correlation time of ionospheric plasma depends strongly on the ion temperature resulting in a decrease of the correlation time as the ion temperature increases. For a 130 µs pulse and given the ion temperature observed by PFISR during the event under study, we expect the performance of the Barker code to significantly decrease for the top of the E region, e.g., above ∼120 km. Above this altitude, the Barker code electron density observations will underestimate the electron density, but the Barker code observations should well-cover the electron densities contributing to the Hall conductance. However, it can be shown that the Pedersen conductivity has a peak at the location where the ratio of the ion cyclotron frequency to ion neutral collision frequency equals one; this altitude typically occurs at ∼120 km (e.g., Sangalli et al., 2009; Burchill et al., 2011; Oyama et al., 2012). Therefore, for this event in particular, a significant portion of the Pedersen conductivity is not well-resolved with the ISR observations, and the Barker code observations from PFISR would underestimate the Pedersen conductivities. Preliminary analysis (not shown) of the PFISR long pulse data shows that there is in fact greater density at altitudes above ∼120 km that contribute to the Pedersen conductance. Unfortunately, the long pulse data has low range resolution (∼72 km) and would require deconvolution (e.g., [35]) to be more useful for conductance calculations. We therefore focus only on the Hall conductance. This works to the strength of the ASI method, since the Hall conductance is essentially a function of Q when average energies exceed 5 keV [29]. Q is essentially a function of the blueline emission [22], which we show later is well-resolved with white light data. We calculated the Hall conductivity using Eq. 10 from [36], where Ω i corresponds to the ion cyclotron frequency, B corresponds to the magnetic field, e corresponds to the elementary charge, and z represents the altitude. The altitude resolved electron density, n e (z), is obtained using PFISR. The ion neutral collision frequency, ] in , is calculated using equation 4.88 from Schunk and Nagy [37] for the non-resonant ion neutral collision frequencies and expressions in Table 4.5tbl45 of Schunk and Nagy [37] for the resonant ion-neutral collision frequencies.
This formula has been implemented in the "flipchem" Python package. Eq. 8 is approximately valid above 80 km altitude [36].
To obtain conductance, we perform an altitude integration of conductivity using a Simpson method integration routine. These formulas have been validated with previous results presented in Kaeppler et al. [9]. Results are discussed in Results.

Determining Energy Flux and Average Energy From Aurora Color Ratios
Following Hecht et al. [24,25], we used the B3C auroral electron transport code [27,28] to calculate lookup tables that allow the integrated energy flux (Q; erg cm −2 s −1 ), the average energy (E; keV), and the oxygen scaling factor (fO) to be derived from measurements of the green, red, and blue emissions (see Figure 2). Since the publication of Hecht et al. [24], the transport code has been updated especially with respect to the green line emission. Hence, the table reproduced below in  Figure 8 in the earlier publication. The B3C code solves the coupled set of linear Boltzmann equations for electrons, H + , and H atom fluxes with full collisional processes [38] for a given set of cross sections and neutral densities for altitudes from 500 km down to 90 km. The B3C calculations are made with atmospheric neutral densities provided by the empirical MSIS90 model [39] that is parameterized by the geomagnetic Ap and solar flux F10.7 indices. (MSIS90 describes the neutral temperature and densities in Earth's atmosphere from ground to thermospheric heights.) The code accepts a choice of either a Maxwellian or a Gaussian incident precipitating electron flux spectrum and the associated Q and E for each geographic latitude and longitude of interest. We ran B3C for Gaussian and Maxwellian distributions to provide values for both discrete (Gaussian) and diffuse (Maxwellian) aurora. Figure 2A depicts the Gaussian distribution's average energies as a function of green to blue and red to blue ratios; Figure 2B depicts the Maxwellian version. A more complete discussion of the MSIS model is given in Background and Overview and Determining Energy Flux and Average Energy From Aurora Color Ratios of Hecht et al. [25].
The specific attributes of the model in this work are as follows. MSIS90 was run for 16 February 2010, the date of our event study. The location was chosen as Fort Smith, Canada at 60.0°N, 248.2°E, which is where our primary instrumentation were located. The F10.7 was 85.5, F10.7 average was 84.7, and the daily Ap was 22 (we also incorporated the other six Ap values for 16 February 2010). The reference O/N2 (an fO of 1) for a column density of 10 17 (cm −2 ) of N2, following Strickland et al. [40], is 0.53. This reference value allows a comparison of fO values for different model runs. This is a disturbed atmosphere for higher latitudes as would be suggested by comparing the reference O/N2 (0.53) with the quiet (0.91) and disturbed (0.45) values used in Hecht et al. [30]. Figure 1 of Hecht et al. [30] suggests that uncertainties in the chosen model atmosphere will affect the derived E values. However, for the E values of interest here we would expect the derived values of E might be lower by 10-20% if the actual atmosphere were somewhat more disturbed than the chosen reference model. The exact bias on the derived E and fO values depends upon the area on the plot in Figure 1 being considered.
To determine the average energies more precisely from the non-uniformly gridded look-up table, we first interpolated the 6 × 30 [fO × E] look-up table values across log-log space to a 50 × 50 array, allowing for 50 red/blue values and 50 green/blue values, each with a corresponding E and fO value. Q was determined by dividing the 427.8 nm intensity (in Rayleighs) by the yield (Rayleigh/erg cm −2 s −1 ), which is a function of E and fO [22] and is another output from the B3C run.

Converting White Light Intensities to Color Intensities
In order to use the methods described in Determining Energy Flux and Average Energy From Aurora Color Ratios, we had to first convert the ASI white light counts to red (630.0 nm), green (557.7 nm), and blue (427.8 nm) intensities. We accomplished this by fitting log-log linear functions via ordinary least-squares linear regression [41] to scatter plots of white light counts vs. green light intensity at the longitude line where the ASI data overlaps the NORSTAR MSP data at Fort Smith, NWT. We then fit scatter plots of the MSP blue light vs. green light, and the MSP red light vs. green light. We fit the data points from 07:00-10:00 UT on 16 February 2010.
Before developing the fits, we processed the ASI data by removing background light by subtracting a quiet-time snapshot from each station. The snapshot was chosen from a dark period immediately preceding any aurora observed in the ASI field of view, prior to the substorm. Starlight was removed by applying a median filter, since the starlight moves over time from the original dark sky snapshot. We ignored MSP data with <50 Rayleighs to avoid color ratio spikes that can occur when the denominator is near zero.
We performed the fitting in log-log space. In fitting the green to the white, we performed an ordinary least squares fit bisector (minimizing the sum of the squares of both white and green residuals [41]) to the median line over the white light range of 600-20,000 counts, to put the weight in the brighter aurora and because count levels below 600 show increased scatter and signalto-noise ratio (see Figure 3A). (Above 20,000 counts had too few data points per bin.) Figure 3A plots the median and quartiles in yellow (using a white count bin size of 50 counts) and the fitted line in green. The quartiles (dashed yellow lines) are included to demonstrate the spread (uncertainty) in the data. We slightly altered the fit line below 600 counts so it would lay at or above the lower quartile, which still weights the higher counts but doesn't ignore the shift in median slope below 600 counts completely. When we tried shifting the fit below 600 counts so that it matched the median line, the result was that we frequently over-estimated the green light at low emissions, making it impossible to model low energies (since the green light median does not approach zero when the white light approaches zero).
Instead of fitting the white light counts to blue and red, we fit the blue and red intensities to the green intensity to calibrate colors using the same instrument (the MSP) at the exact same location (removing cross-instrument related errors). Figure 3B plots the blue vs. green fit in blue, and the median and quartiles in yellow (using a bin size of 50 R to calculate the medians and quartiles). We applied the three main methods for OLS described in Isobe et al. [41] that minimizes the sum of the squares of the green residuals, the blue residuals, and both green and blue, finding that minimizing the green residuals resulted in the smallest variance in both the slope and the y-intercept. There is no blue light below 225 R in the fitting because, as mentioned above, we did not include MSP data under 50 R and we converted the 470.9 nm intensity to 427.8 nm by multiplying the former by 4.5.
We next show the red vs. green fitting ( Figure 3C). We split the fitting into three functions, described in Eq. 5-7. The largest deviation in the fit from the median is between the green line intensities 10 3 -10 4 R, where the redline scatter is greatest. We initially tried fitting to the median line within 10 3 -10 4 R by splitting the fitting even further, but discovered that because the slope is approximately zero in this range, the 10 ( 0.98+0.62log10 ( green intensity )) The variances in the slope and the y-intercept [41], respectively, in the log-log space are as follows: Eq. 2: 0.000336 and 0.00589. Eq. 3: 0.018 and 0.071 (when compared to the median line below 400 cts).
Eq. 6: 9.40e-5 and 0.0008. Eq. 7: 0.0040 and 0.0185. We note that in early iterations of our fitting process, we discovered that the red and green aurorae-especially at the time of the bright auroral arc seen in Figure 4 from ∼07:18-07:25 UT-were morphologically offset from one another, with the red mapping to lower latitudes than the green (see Supplementary Figure S1). The result was an unphysical situation in which the green overtook the red at a higher latitude and the red overtook the green at a lower latitude, and the artificial result was much higher average energies where the green arc mapped and much lower average energies where the red arc mapped. To align the 630.0 and 557.7 nm emissions, the effective altitude was found to be 160 km. This is below the expected emission altitude close to 200 km (e.g., Solomon, 1988; [42]). This effective altitude may appear lower than expected due to horizontal winds which are known to blow the redline emission (e.g., [43]). Ultimately, because the point of comparing red to green was to obtain an analytical relationship between intensities, we mapped the red to 160 km to match auroral morphologies and to build the scatterplot in Figure 3C.
The scatterplot in Figure 3C shows that the red and the green intensities are fairly uncorrelated, which may be due to several factors. First, the red to green ratio depends on both the average energy of the precipitating electron distribution and the energy flux. For Gaussian electron distributions that follow the Knight relation [44,45], Q and E should be broadly correlated; but for Maxwellians they are not [30,46,47]. Furthermore, the redline and greenline chemistries have different dependences on atomic oxygen density, which also decreases the correlation. This uncertainty in the average energy Frontiers in Physics | www.frontiersin.org October 2021 | Volume 9 | Article 744298 estimation could be explored further in future studies by defining when and where discrete aurora occur in order to build a functional relationship between green and red only during discrete auroral times. Figure 4, discussed below, as well as Supplementary Figures S3, S4 do however show the energy comparison between MSP and ASI datasets to be relatively consistent. See also Table 1 and Supplementary Tables S1, S2 which quantify the differences between ASI and MSP derived average energies, along with other parameters. Figure 4 plots the MSP data ( Figures 4A,C,E) next to the THEMIS ASI-derived auroral colors ( Figures 4B,D,F) for comparison and validation. Additionally, we determined the energy flux and average energy using the methods described in Background and Overview with the MSP color ratios as input ( Figures 4G,I) and the ASI color ratios as input (Figures 3H,J). Visually the comparison is good, but to quantify the differences we interpolated the THEMIS-derived products to match the MSP location and cadence and calculated the ratio and the difference. The mean, standard deviation, median, and upper and lower quartiles of the ratios and the differences are compiled in Table 1 for both the Figure 4 time window (07:10-08:00 UT) and the full time window used for the calibration (07:00-10:00 UT). We present the ratios, the differences, and a similar table that compares the ASI and MSP color intensities in the supplemental materials (Supplementary Figure S3 and Supplementary Table S1). Finally, we also present the percent differences between the ASI and MSP measured and derived quantities in Supplementary Figure S4 and Supplementary Figure S2. Table 1 Table 1 also tells us that the average of the ratio of average energies is 1.0, with a median of 0.9 and upper and lower quartiles of 1.2 and 0.6, respectively, for the Figure 4 timeframe. These and the rest of Table 1 indicate the level of uncertainty assuming the MSP data as "truth". The very slight underestimation in the differences may be due to the fact we did not remove the <600 R greenline data points from the comparison which were not heavily weighted in Figure 3 for the calibrations. As is shown in Supplementary  Figure S2, for the full calibration time (07:00-10:00 UT), the percent differences are: for Q, 47.5 ±36.1% with median 39.6% and quartiles 18.4 and 68.6%; for E, 46.0 ± 38.2% with median 35.3% and quartiles 15.8 and 67.6%. We also present the green, blue, and red percent differences in the Supplementary Material.
With these relationships to estimate the color ratios, we determined the energy flux and average energy for all pixels that comprise the ASI mosaic. The energy flux and average energy were then interpolated to a 400 × 400 uniform geographic grid with 0.1°latitude and 0.4°longitude resolution (see Figure 5 for an example). Any points at elevation angles lower than 10°were dropped We first applied 2D median filtering with a window size of 3 pixels in both directions. Where more than one ASI fields of view overlap, the pixel with the brightest value usually comes from the imager with better viewing conditions (e.g., higher elevation angle and clearer sky conditions), which are more reliable. We therefore assigned each of the 400 × 400 grid points the maximum energy flux (or maximum energy) value found in any pixel that fell within ±0.05°latitude and ±0.2°l ongitude of the grid point. The uniform grid system is much simpler than the imager-dependent coordinates that overlap with surrounding imagers.
The data in the uniform grid system is used to calculate total energy flux and average energy across the imager coverage (see Results). Conductance maps are also obtained. Data in this coordinate system are also used to provide data to users. (Note that this grid could be altered to a different desired resolution, e.g., 500 × 500 or 300 × 300.)

Determining the Hall Conductance
We calculated the height-integrated Hall conductance from the incident precipitating electron flux energy spectrum using the B3C auroral transport code for the same F10.7 and Ap values used in Determining Energy Flux and Average Energy From Aurora Color Ratios. The code was run for either a Maxwellian or a Gaussian incident precipitating electron flux spectrum with an associated integrated energy flux and an average energy for each geographic latitude and longitude of interest. A Maxwellian distribution is representative of diffuse auroral electron precipitation whereas a Gaussian distribution is descriptive of an accelerated distribution associated with discrete auroral precipitation [48].
Similar to our methods to obtain energy flux and energy in Determining Energy Flux and Average Energy From Aurora Color Ratios, we used B3C to create a look-up  1 | The mean, standard deviation, median, and upper and lower quartiles of the ratios and the differences between ASI-and MSP-determined energy flux (Q) and average energy (E). The first four lines are for the Figure 4 time window (07:10-08:00 UT). The bottom four lines are for the full time window used for the calibration (07:00-10:00 UT).

Mean
Std   The Hall conductance was then determined across the 400 × 400 ASI grid by interpolating the look-up table to the given parameters. The energy flux and energies in the B3C look-up table that correspond to a conductance value were treated as a uniform 2D grid that could be interpolated to each observed energy flux and energy stored in the ASI grid. Next, the geographic location in the look-up table closest to the given ASI point (found using the haversine between two geographic points) was used.
Although conductance is computed rigorously with the B3C auroral transport code, for contextual purposes we also calculated the Hall (Σ H conductance with the widely-used Robinson et al. (1987) formulae, where Σ P is the Pedersen conductance: 0.85 (9) In this case, we simply used the energy flux and average energy value at each of the 400 × 400 grid points as inputs to the formulae. We did this for both the Gaussian and Maxwellian versions, but recognize that the Robinson formulae were intended to provide simple approximations of the conductance. Figure 6 displays snapshots of the energy flux and average energies determined by the THEMIS ASIs throughout the first substorm lasting from ∼07:18-09:30 UT. The values are from the Gaussian distribution; we include the Maxwellian version in the supplemental material as Supplementary Figure S5. Movies of the substorm are also included in the supplemental material (Supplementary Movies S1-S4). The ASIs capture the substorm onset, expansion, and westward-traveling surge as the vortical auroral form that has elevated energy flux and propagated westward. Note the relatively good agreement in arc location, energy flux, and energy across different ASI FOVs, suggesting overall the mapping is good. We also note that the strongest energy fluxes and average energies come in mesoscales during the expansion phase (e.g., 07:28:30 UT). We next interpolated the data onto the 400 × 400 pixel grid to determine the total energy flux deposited in the ASI FOV (above 10°elevation) as well as the average of the average energy. The total energy flux at any given point in time was found by integrating the grid points over the 2D area:

RESULTS
We ran the ASI grid through a 2D median filter to determine what amount of the total energy flux was deposited as mesoscale auroral forms. For example, to determine the amount of energy flux contributed by scale sizes 500 km and less, for each pixel we found the median across 500 km in 2D (both latitude and longitude directions). We then found the total for that scale size by applying Eq. 10 to the median filtered data, and subtracted that total from the total energy flux: We weighted the average energies by the energy flux to perform a weighted average of the average energies: This puts emphasis on the energies that comprise the aurora. Figure 7 displays these results alongside the magnetograms (Figures 7B,C) and keograms (Figures 7D,E) from the two selected stations that saw the substorm onsets: FSMI and FYKN. Figure 7A displays the SML index [49], which is a measure of the global substorm evolution and strength. The first substorm starts at ∼07:18 UT as a bright, narrow arc that propagates poleward across FSMI's field of view until the westward traveling surge moves over the keogram longitude (at ∼7:30 UT), showing auroral expansion both poleward and equatorward. The tertiary expansion is seen in Figure 7E as the westward traveling surge brightens and travels over FYKN after 08:00 UT. A bright, latitudinally narrow arc forms at ∼66°m agnetic latitude (mlat) during the recovery phase and moves equatorward (seen in Figure 7D from ∼08:45-09:33 UT) until the second substorm onset at ∼09:33 UT, seen in Figure 7E. Figure 7F displays the total energy flux deposited within all available THEMIS ASIs in black, the total energy flux from scale sizes under 500 km, 300 km, and 100 km in red, green, and blue, respectively. Figure 7G displays the percentage of the total energy flux that is contributed by those three scale sizes. Figure 7H displays the average of the average energy.
The main take-aways from Figure 7 are that the mesoscales (<500 km scale sizes) are important in the deposition of energy flux to the ionosphere, particularly so immediately after substorm onset during the rapid expansion phase where they contribute up to ∼80 and ∼60% of the total energy flux in substorms 1 and 2, respectively ( Figure 7G, at 07:20 UT and 09: 36-09:38). Even the <100 km scale sizes, which only contribute ∼10% of the energy flux throughout the substorm on average, contribute >30% of the energy flux during the rapid expansion phase. Because the difference between the <500 km and <300 km total energy flux is roughly 10%, whereas the difference between the <100 km and <300 km total energy flux is roughly 20%, we can also surmise that more of the energy flux is deposited between 100-300 km scale sizes than 300-500 km scale sizes. The rapid expansion phase is also the most energetic, as the average of the average energy spikes then for both substorms at ∼6.75 and ∼5 keV, respectively ( Figure 7H). The maximum average energy is of course, higher. The mesoscale contribution begins to decrease towards the end of the expansion phase, dropping significantly by the beginning of the recovery phase. From Figure 7 and Movies 1-2mmc1-2, we see this is a combination of the mesoscale auroral forms evolving into larger-scale auroras as well as becoming less intense.
These results are related to those presented by Partamies et al. [18], who studied the substorm evolution of auroral arc structures. They found that auroral arcs occur in all conditions as the main element of the aurora, leading to the expectation that mesoscale auroral forms should play an important role in the total precipitating energy flux. They also found that auroral arcs occupy 8.7% of the total growth phase duration, 13.2% of the total expansion phase duration, and 9.2% of the total recovery phase duration. These percentages indicate the prevalence of arcs during the expansion phase-which is when we found the most contribution from mesoscales to the total energy flux. Partamies et al. [50] showed that the average auroral emission intensity is lower during the growth phases of substorms than the bright arcs during the expansion phases, which also aligns with our results showing the total energy flux is lower during the growth phase than the expansion phase.
Auroral arcs in the form of longitudinally localized poleward boundary intensifications and streamers have been correlated to mesoscale plasma sheet fast flows (on the scale of a few R E ) (e.g., [51]; [52][53]; [54][55][56]). These mesoscale fast flows were shown to play a large role in magnetotail magnetic flux transport [57,58], and now look to also play an important role in providing the precipitated energy flux during substorms. Further analysis to determine the nature of the mesoscale auroral forms-whether or not they are streamers-could be done to investigate how much of the mesoscale contribution to precipitated energy flux is due to these mesoscale plasma sheet fast flows.
We also see that the peak total energy flux of each substorm occurs in the middle of the rapid expansion phase, immediately after the peak in mesoscale contribution ( Figure 7F). The vertical dashed lines mark the start of the sharp increase in total energy flux and the end of the sharp decrease in total energy flux that occurs during the substorm expansion phase for each substorm (07:18 and 07:33:20 for the first, 09:34:52 and 09:45 UT for the second). Note that these are coincident with the sharp increase and decrease of the mesoscale contribution. They are also approximately coincident with the substorm onset and the end of the substorm expansion phase, defined by the SML index. The SML index also shows that the first substorm had two more expansions before the final recovery, during each of which there was an increase in total energy flux as well as an increase in mesoscale contribution.
Analyzing the Supplementary Movie S1, the increased mesoscale contribution at the end of the first substorm (∼09: 18-09:34:42 UT) appears to be due to the latitudinally narrow, pre-onset arc that formed and moved equatorward prior to the second substorm. The increased mesoscale contribution at the end of the second substorm (∼10:51 UT onward) appears to be from the aurora breaking up and the formation of pulsating aurora (e.g., [59]), a common occurrence after substorms due to whistler mode chorus waves in the magnetosphere. Figure 8 plots the Hall conductance as determined by PFISR, B3C with ASI inputs, and Robinson formula with ASI inputs (methods described in Poker Flat Incoherent Scatter Radarand Determining the Hall Conductance). We used the FYKN ASI to determine when a discrete auroral arc was over PFISR, during which times we used the Gaussian distribution to determine average energy and energy flux. For other times, we used the Maxwellian distribution. Figure 8A plots the FYKN ASI keogram at the PFISR longitude. The horizontal line represents the PFISR latitude. Mapping the aurora to 110 km, the PFISR measurements lie somewhere between a 31-35°elevation angle in the FYKN field of view. Note that there were gaps in the PFISR data, which is when the black line in Figures 8B,C disappears. Figures 8D-G illustrates the Hall conductance from B3C plotted in the 400 × 400 grid with the PFISR location marked. Figure 8D (09:36:30 UT) is a 2D snapshot of the Hall conductance assuming a Maxwellian distribution, during the time prior to the large increase in conductance due to the bright discrete auroral arc seen in Figure 8F. Figures 8F,G (09:40 and 09:41:30 UT) plot the Hall conductance assuming a Gaussian distribution while the discrete auroral arc is overhead. We see a decrease in conductance when the arc moves slightly away from PFISR in Figure 8G. Figure 8H (10:30:30 UT) again plots the Hall conductance assuming a Maxwellian distribution after the diffuse aurora has returned. Because PFISR measurements are at a lower elevation angle within the FYKN ASI FOV, there is a greater likelihood of introducing error. First, it is possible that an arc maps to a higher or lower altitude than 110 km, in which case a narrow arc above PFISR could be incorrectly mapped in latitude. However, comparing the ASI-informed conductance with the PFISRdetermined conductance, we see that the ASIs reproduce relative variations in conductance that is derived from PFISR very well. The fact that the conductances vary at the same time for both the ASI-and PFISR-informed values suggests the mapping is not a primary issue in this case.
Additionally, the ASI is not looking up the magnetic field line at that latitude, whereas PFISR and B3C are integrating along the field line. The result is that not all of the auroral intensity that is actually above PFISR will be mapped to the PFISR latitude and some intensity that is not above PFISR could be mapped to the PFISR latitude, which can result in an under-or over-estimation of the average energy and/or energy flux by the ASI. However, we see that the B3C model informed by the ASIs does exceptionally well at supplying the Hall conductance magnitudes, as it closely follows the PFISR Hall conductance during both discrete and diffuse aurora. Therefore, although the above cautions remain, in this case the methods do well. In contrast, the ASI-informed Robinson formula overestimated the Hall conductance during discrete aurora and underestimated during diffuse aurora. This highlights the improvements physics-based atmospheric modeling provides.
We list the mean and standard deviations, median and quartiles of the differences, ratios, and percent differences between the various B3C informed by ASI (ASI B3C ) and PFISR results in Table 2. Plots of the percent differences over time are included in Supplementary Figure S6. The diffuse aurora (Maxwellian) has a better match than the discrete, which is expected because discrete aurora is more variable and comes on smaller scale sizes, which are more likely to be mapped less precisely for the purposes of exact comparison between datasets. For the discrete aurora, the median percent difference is 24.2% whereas for the diffuse aurora the median percent difference is 18.1%.

SUMMARY
This paper presented a new methodology to utilize the mosaic of THEMIS all-sky-imagers to study auroral precipitation in 2D over time, estimating parameters such as energy flux, average energy, and Hall conductance. Using two consecutive substorms for a case study, results demonstrated the importance of mesoscale auroral forms (10 s km to 500 km scale sizes) to the overall energy flux deposition to the ionosphere, particularly during the initial expansion phase immediately after onset. Roughly 60-80% of the total energy flux immediately after onset, and ∼30-55% thereafter, is from mesoscale auroral forms. These results confirm what we may have hypothesized from the Partamies et al. [18,50] observations, since they had found that auroral arcs occupy the greatest percentage of the substorm expansion phase duration (contrasted with the growth and recovery phases) and are brighter than during the growth phase.
Because prior work highly relied upon statistics or single ASIs, the mesoscale contribution has not been so well-resolved before over a continental scale, especially for specific case studies. The results indicate the importance of including mesoscale auroral forms in global models to appropriately capture the related physics. Additionally, the paper demonstrated the 2D time series data product that can be provided to modelers when ASI data are available, allowing for case-specific data-informed modeling.
The methods showed a good comparison between the calibrated ASI data and the MSP data, as were quantified in Table 1 and Supplementary Tables S1, S2. As expected, the greatest scatter is in 630.0 nm wavelength. As a future work, the redline ASI cameras that are currently coming online as part of the TREx array of ASIs could be used to specify the 630.0 nm wavelength more accurately (https://www.ucalgary.ca/aurora/ projects/trex). Comparing and/or calibrating to DMSP energy fluxes and average energies would also provide additional information on the uncertainties and/or improve our estimates. We do point out as mentioned above that the inversion method in Janhunen [13] found little difference in the energy flux results when the red aurora was included, since the energy flux is mostly related to the blue intensity [22,30], which is a close function of the green intensity ( Figure 3).
The Hall conductance was well-captured by the B3C auroral transport code, informed by ASI energy fluxes and average energies, when compared to PFISR-determined Hall conductance.
The ASI-informed Robinson formula overestimated the Hall conductance during discrete aurora and underestimated during diffuse aurora. As expected, the conductance increased during the intense discrete aurora. These values are not well-constrained in current models that rely on global, diffuse aurora as precipitation input only. These data could therefore be useful to modelers who want to include the larger conductance values that result from mesoscale, discrete aurora. As a future work, a date should be selected when PFISR can provide alternating code observations of the whole E region for a better comparison, and also to include the Pedersen conductance.

AUTHOR CONTRIBUTIONS
CG led the work, updated techniques developed by TN, performed the analyses, and wrote the paper. TN devised the original ASI-MSP techniques and advised CG on the details. MC ran the B3C code to provide the conductance look-up table for CG to use with ASI inputs, and described the methods in the text. JHH ran the B3C code to provide Figure 1 and the look-up table to determine energy fluxes and average energies from color ratios, as well as the associated text in the manuscript. SK calculated the conductance from the PFISR data and contributed to writing the PFISR section. DMG provided processed MSP data and text describing the instruments. AR provided processed PFISR data and text describing the instrument and dataset. LL contributed discussions about mesoscale auroral phenomena. YD contributed discussions about mesoscale phenomena in the ionosphere and insights on what parameters are important to specify. ED runs THEMIS-ASI in Canada, and contributed discussions on our methodology and use of the MSP data. JE manages the B3C code and gives permission and guidance for its use. All co-authors contributed comments, edits, and text to the manuscript.

FUNDING
The work in this paper was gratefully funded by NASA grants 80NSSC20K0725, 80NSSC18K0657, 80NSSC20K0604, 80NSSC20K0195, and 80NSSC20K1786, AFOSR Grant FA9559-16-1-0364, and NSF grants AGS 1602862 and AGS-1907698. The THEMIS mission is supported by NASA contract NAS5-02099, NSF grant AGS-1004736, and CSA contract 9F007-046101. This material is based upon work supported by the Poker Flat Incoherent Scatter Radar which is a major facility funded by the National Science Foundation through cooperative agreement AGS-1840962 to SRI International.