Removing Orbital Variations From Low Altitude Particle Data: Method and Application

Particle flux measurements from polar orbiting low altitude satellites provide a view of the near Earth radiation environment that is extremely valuable for science as well as space weather monitoring. Unlike, geosynchronous satellites that sample only a limited region of space (L = ∼6.6), these low altitude satellites sample the extended radiation environment (L = 1 to >10) at a relatively high time cadence (tens of minutes) that captures its global evolution. While these data are clearly useful, it is also challenging to work with because the particle flux measurements have large orbital variations related to the changing geographic location of the satellites. These orbital flux variations can sometimes obscure the time variations of interest for scientific study or space weather hazard awareness. Here, we describe and evaluate a method for removing these variations that is based on Statistical Asynchronous Regression. We demonstrate the utility and accuracy of the method by applying it to electron flux measurements from the NOAA POES and EUMetSat MetOp low altitude satellites.

variations that may obscure the signature of those processes. For example, the SAR method was developed and applied to electron flux measurements made by the NOAA GOES satellites at GEO orbit. At this orbit, particle fluxes on the dayside may be 1-2 orders of magnitude higher than those on the nightside resulting in an extreme diurnal variation as the GEO satellites orbit Earth. In this case, the particle fluxes are enhanced on the dayside where the magnetic field is compressed and reduced on the nightside where the field is stretched outward (see Green et al., [2013], Onsager et al., [2002] for a more detailed explanation of this effect). The SAR method was successfully used to map measured GOES fluxes to one consistent location, noon local time, facilitating studies at sub-orbital time resolution. The mapped data were used to clearly demonstrate differences in driving between those geomagnetic storms that resulted in >2 MeV electron flux enhancements and those that did not ). Those differences occurred within the first 24 h after the start of the recovery phase of storms and would not have been evident without first removing the diurnal electron flux variations. The technique has also proved valuable for space weather applications. A simplified version of the SAR technique is now currently used by the NOAA Satellite Environmental Anomalies Expert System Real Time (SEAESRT) which indicates the hazard of an anomaly occurring on a satellite due to internal build up of high energy electrons (O'Brien, 2009). The SEAESRT system uses the simplified SAR technique to provide this hazard information at all GEO longitudes based on measurements from just one GOES satellite.
In a similar manner, our goal here is to map the low altitude polar satellite measurements to a fixed longitude and hemisphere to remove the large orbital variations and highlight changes that occur on timescales shorter than a day. Commonly, polar satellite data, when used in tools to enable satellite anomaly monitoring, is averaged to a daily cadence to reduce such variations. For example, the NOAA National Centers for Environmental Information (NCEI) provide daily belt indices from the POES/MetOp particle flux measurements that can be used for situational awareness to understand how the current radiation belt intensities compare to statistical averages (https://satdat.ngdc.noaa.gov/sem/poes/data/belt_indices/). The belt indices have been included as indicators used in decision trees to help satellite operators determine whether a satellite anomaly is likely caused by the radiation environment (O'Brien et al., 2012). Another example, and our specific interest here, is the neural network model known as Specifying High-altitude Electrons using Low-altitude LEO systems or SHELLS (Claudepierre & O'Brien, 2020). The SHELLS model was developed to predict global electron flux based on low altitude data. SHELLS uses a neural network to define the relationship between low altitude POES/MetOp electron flux measurements to those near the equator observed by the Van Allen Probes. Once that relationship is defined, the neural network can specify (nowcast) global electron flux using only the low altitude data. The model currently uses daily averages of POES/MetOp electron flux as input. Tailoring and validating the SAR technique for the POES/MetOp particle flux measurements would enable the SHELLS model and other tools based on LEO data to be output at a much shorter time cadence of tens of minutes rather than days.

Data
The primary data set used here is the electron fluxes from the Medium Energy Proton and Electron Detector (MEPED) carried by the POES satellites operated by NOAA and the MetOp satellites operated by EUMetSat. This constellation consists of satellites in highly inclined polar orbits at ∼850 km. Since July 1979, 11 POES satellites and three MetOp satellites have launched that carried a version of MEPED, making the data set an excellent source for understanding long term and near real time variations in the space environment. At present, we focus our initial analysis on data from five satellites (NOAA 15,18, and MetOp-02 [or MetOp-A]) only because our eventual goal is to compare the data with that from the Van Allen Probes (Mauk et al., 2012). These NOAA satellites were all operational during the Van Allen Probes mission from August 30, 2012 to February 12, 2019. The actual comparison of the POES and Van Allen Probes data is not presented here but will be described in greater detail in future work. Here we focus solely on the details and implementation of the SAR method that can be more generally applied for many scientific purposes.
The MEPED instruments measure high energy electron and proton flux with both telescope and omnidirectional dome detectors (Evans & Greer, 2000;Green, 2013). Our initial analysis considers only the electron telescope measurements but the method could also be applied to those from the proton telescopes.
GREEN ET AL.

10.1029/2020SW002638
Additionally, the telescope detectors measure the electron flux in two directions: roughly radially outward from the Earth (0° detector) and perpendicular to the radial direction (90° detector). We demonstrate the technique here using data from the 90° detector. (Data from the 0° detectors are more severely impacted by contaminating proton flux that should be evaluated and removed prior to applying the SAR technique [Lam et al., 2010].) The MEPED data are currently made available online by the NOAA NCEI at https://www. ngdc.noaa.gov/stp/satellite/poes/dataaccess.html. Several versions of data are available: raw, processed corrected, and processed uncorrected data. The raw data give the original instrument counts not translated into meaningful particle flux units. The processed corrected data (Peck et al., 2015) give the data translated from counts into flux with an attempt to correct or remove the known proton contamination in the electron detectors. Unfortunately, these data are only available up to 2014 and does not include the full time period of the Van Allen Probes era that we are most interested in. These data are also not suitable for real time products since it does not update as new data are collected. The processed uncorrected data used here provide the data in flux units (#/cm 2 -s-sr) at four integral energy channels: >40 keV (e1), >130 keV (e2), >287 keV (e3), and >612 keV (e4). This format of the data includes data beginning in 2012 and is also updated in real time as data are received by the ground stations.

Method
As mentioned in the introduction, the SAR method used here works by establishing the relationship between the statistical distributions of particle fluxes sampled at different locations. The underlying assumption here is that while the flux levels at two locations may be very different, their variations will be highly correlated. More specifically, the main premise of the method is that the percentile flux levels at each particle energy will be consistent at different locations. For example, if a satellite at GEO measures 2 MeV fluxes at the 90% level at noon magnetic local time, the method assumes that the 2 MeV fluxes will be at the 90% level simultaneously at all other local times. Mapping measurements from one location to the next requires knowledge of the cumulative distribution of the fluxes at each location or local time. This method is expected to hold as long as the flux changes are highly correlated, the relationship between different locations is monotonic, and the data are stationary (i.e., the cumulative distributions are not changing with time). At GEO, this method was applied by finding the cumulative distributions of the log 10 (>2 MeV electron flux) in 1 h local time regions over 8 years . Relating the percentile flux levels from one local time to another is done by plotting the two cumulative distributions against each other (also known as a "QQ" plot). The resultant curve is then fit to an assumed function or captured in a simple lookup table that relates the values.
At LEO, applying the SAR method follows the same logic but is slightly more involved because the satellites sample a larger range of L shells. In essence, the method is the same as used at GEO (L = ∼6) but it is done separately for each L shell and instead of comparing data in different local time bins, we compare data in different longitude bins. In addition, the data must be separated by hemisphere (northern and southern) and satellite direction (northbound and southbound). Separating by hemisphere is necessary because the magnetic field strength and thus the equatorial pitch angle of the particles being measured (which is in part responsible for the geographic variations) is different in the two hemispheres. Typically, the flux increases with increasing equatorial pitch angles. The dependence on the magnetic field strength is shown in the equation for the equatorial pitch angle of the particles, where  eq is the equatorial pitch angle, B equator is the magnetic field strength at the equator, B local is the magnetic field at the location of the satellite, and  localis the pitch angle measured at the satellite. Separating by northbound and southbound satellite direction is necessary because the 90° telescope is oriented at a slightly different angle relative to the magnetic field depending on the satellite direction. Even though the magnetic field strength is the same in the two directions, the particles being measured will have different local pitch angles (  local in Equation 1) and thus different equatorial pitch angles. Therefore, the northbound and southbound portions of the orbit will have different flux levels and distributions. Figure 1 demonstrates how the equatorial pitch angle varies as a function of geographic location and satellite direction at the altitude of the POES satellites. The figure shows the average expected equatorial pitch angle of the particles being measured by the 90° detector and how that varies as a function of L, longitude, hemisphere (North and South), and satellite direction (Northbound and Southbound). The values are averaged in 0.25 L and 10° longitude bins. Here the equatorial pitch angle is approximated using the IGRF magnetic field as B local and a simple dipole field for B equator . Clearly, the actual magnetic field is more complicated than a dipole and B equator will have some magnetic local times (MLT) dependence. In this case, Figure 1 is simply meant to qualitatively demonstrate the major differences in the sampled equatorial pitch angles particularly between the Northern and Southern hemispheres. In the northern hemisphere, there is very little difference in the sampled equatorial pitch angles as a function of longitude. Whereas, in the southern hemisphere the equatorial pitch angle may vary by 20° across longitudes at a fixed L shell. (Some areas in the northern hemisphere are left empty because there are negative values for L IGRF given in the POES/MetOp datafiles. These values are negative because particles here will encounter the atmosphere in the southern hemisphere.) Additionally, the northbound and southbound directions must be separated because they occur at different MLT for each satellite which can also affect the flux levels and overall distributions (see Figure 2). Lastly, we separate the data by geomagnetic activity because several studies indicate improvement in mapping after inclusion of such indicators (Burin des Roziers et al., 2006;Su et al., 2014). We expect the distributions to change with geomagnetic activity due to changes in the topology of the global magnetic and electric fields.
Here we use a single indicator of geomagnetic activity, the Kp index. In order to, maintain meaningful statistics, data with Kp≥4 are combined into a single bin.
To apply the method, we create cumulative distributions of the log 10 (electron flux) for each of the four integral electron energy channels from each satellite in 0.25 L and 10° East longitude bins. We do so for northern and southern hemispheres and the northbound and southbound satellite directions over the 5 year time period from 2014 through the end of 2018. We chose a 5 year time span in order to, include enough data to create a meaningful cumulative distribution but limit changes that result from the satellite orbits drifting in MLT. For example, Figure 2 shows that the MetOp satellites stay fixed in MLT but the NOAA POES satellites drift in MLT and at different rates. We expect this MLT drift to create a non-stationary data GREEN ET AL. set and cause changes in the cumulative distributions over time as shown by Asikainen and Ruopsa (2019).
In their analysis, fluxes closer to the dayside were up to two orders of magnitude higher than those on the nightside. In the 5-year period from 2014 to 2018, the N18 satellite drifts the most, from ∼15 to ∼20 MLT on its northbound section of the orbit. In order to account for the MLT drift when applying the method to a very long-term data set, cumulative distributions should be created for each year from data within a surrounding 5 year time window. This will effectively create cumulative distributions that vary with MLT as the satellites drift in time.
An example of the flux distributions and how they vary with longitude and L is shown in Figure 3 for the lowest electron energy channel on the MetOp-2 satellite and Kp = 0. The figure shows the median flux in each bin for each hemisphere and satellite direction. The expected two-belt structure is evident with an inner peak near L∼3 and an outer peak near L∼5. Superimposed on the belt structure is a longitude variation that is most evident in the southern hemisphere. Comparing Figure 3 with the equatorial pitch angles in Figure 1 shows some correlation. At any given L shell the fluxes are consistently higher at longitudes where larger equatorial pitch angles are sampled. For example, at L = 5 in the southern hemisphere (Figures 1c  and 1d) the largest equatorial pitch angles and highest fluxes (Figures 3c and 3d) are observed near 0° longitude +/−30. Also, comparing the median fluxes for northbound and southbound directions shows that they both follow the equatorial pitch angles but with small differences. One obvious discrepancy is clear in the northern hemisphere fluxes (Figures 3a and 3b). In the northern hemisphere, there is very little variation in the equatorial pitch angles as a function of longitude but there is a very large flux decrease near 0° longitude +/−30. The low fluxes at these longitudes are due to the fact that particles here will mirror at such low altitudes when they bounce down to the southern hemisphere that they are absorbed by the atmosphere and lost from the magnetosphere. Figure 4 shows how the median fluxes vary with energy. The flux levels decrease as energy increases and the peak fluxes move inward to lower L shells. Figure 5 demonstrates how the median fluxes change with Kp. Here the fluxes decrease at large L shells and intensify at lower L shells as Kp increases likely due to changes in both the global magnetic topology and convection electric field.
Once the cumulative distributions for each satellite are created, they are used to map the fluxes as follows. drifting orbit. For our reference longitude and hemisphere, we use 20-30° in the southern hemisphere for two reasons. One is that we are interested in comparing the POES/MetOp electron fluxes to the near equatorial fluxes measured by the Van Allen probes. At this longitude in the southern hemisphere, the magnetic field strength is low so the satellite is more likely to be sampling trapped particles with larger equatorial pitch angles. The other reason is simply that mapping to this longitude gives good results when analyzing GREEN ET AL.
10.1029/2020SW002638 6 of 15  how well the diurnal variations are removed as discussed further in the next section. Using either satellite direction gives similar results so we arbitrarily chose the southbound satellite direction that corresponds to ∼9 MLT. Next, to map the fluxes, we find the percentile flux level at each time step based on the compiled distributions for that satellite, energy channel, hemisphere, direction, longitude, and the measured Kp.
The mapped value will be the flux at that same percentile level for the MetOp-2 satellite and our chosen reference hemisphere (southern), direction (southbound), longitude (20-30°) and the measured Kp. Rather than creating fits between cumulative distributions for each satellite and each bin using QQ plots, we do the mapping using tabular representations of the cumulative distributions. The main reason for this choice is that the QQ plots (plots showing the quantiles of the first data set against the quantiles of the second data set) do not always follow a simple functional form such as a power law. As an example, Figure 6 shows QQ plots comparing the flux quantiles measured by the MetOp2 satellite at the reference longitude to the flux quantiles measured by the NOAA-15 satellite at six longitudes and 3 L values.
One important note here is that our application of the SAR technique not only removes the orbital variations but also effectively accounts for any differences in intercalibrations (Asikainen, 2019)    between fluxes measured at different energies as long as the fluxes are still highly correlated and the relationship between the distributions is monotonic (see Appendix A for more details on these requirements). In our application of SAR, we effectively map the fluxes from the precise energy range measured at the observing satellites to that measured by the reference satellite, MetOp-2. Of course, to ensure a high degree of correlation, we map between similar energy channels on each satellite.
The code developed to apply and analyze the SAR method and demonstrate the results is written in python. It relies heavily on the scipy , numpy (Harris et al., 2020), matplotlib (Hunter, 2007), and xarray libraries (Hoyer & Hamman, 2017). The work likely would not have been feasible to produce in such a relatively short amount of time without these fundamental tools. Figure 7 demonstrates the results of the method by showing 1 year (2014) of electron flux from the lowest energy channel (e1, >40 keV) on the NOAA 15 satellite with and without the mapping as described in the previous section. One notable difference between the mapped and original data here is that the mapped fluxes are significantly higher. The increase is primarily due to the fact that the fluxes are mapped to the southern hemisphere at longitudes where the overall flux levels are high because of the low magnetic field strength. Some orbital variations are clear in the unmapped data (especially at low L < 4) with noticeable improvement after mapping. However, over a year time period differences are difficult to discern at all L shells. In order to highlight these variations and their reduction after mapping, Figure 8 shows only the first 2 weeks of 2014 and Figure 9 shows the same type of plot but for the higher energy e3 channel (e3 GREEN ET AL.   >287 keV). At higher energies the orbital variations and change after SAR mapping is even more striking and highlights the value of applying the technique. Looking at panel c, one could potentially wrongly interpret some of the daily order of magnitude increases in the flux near L∼4 as being caused by a physical mechanism such as interaction with electromagnetic waves. In fact, the mapped data in panel d reveal a more gradual and steady flux increase.

Results: SAR Mapping and Its Accuracy
GREEN ET AL.

10.1029/2020SW002638
10 of 15  To further illustrate the orbital variation in the data prior to SAR mapping and highlight how well it is reduced after SAR mapping we use Fourier analysis. Prior to the mapping we expect to see daily peaks (and harmonics) in the Fourier spectrum of the data, since the satellites revisit the same longitude and L shell once per day (Northbound and Southbound) as Earth rotates beneath. If the SAR mapping is effective, it should completely eliminate or reduce those peaks in the spectrum. Figure 10 shows the amplitudes from Fourier analysis before and after mapping at different L shells for the e1 channel and Figure 11 shows the same for the e3 channel. All data sets are detrended first to remove longer variations using a Savitzky-Golay filter (Here we use the python scipy savgol_filter function with a window length of 25 points and a polynomial order of (1). As expected, obvious peaks in the Fourier amplitudes are apparent in the unmapped data (blue trace) at 24, 12, and 8 h. At all L shells these orbital beating harmonic amplitudes are significantly reduced or even completely eliminated by SAR mapping demonstrating that the technique works as expected.
Lastly, we use one final comparison to demonstrate the overall accuracy of the SAR mapping. Our goal here is to directly compare mapped data with observed data and quantify the difference between the two.
To do so, we compare NOAA-15 data with MetOp-2 data when the satellites are at the same L shell and the passes are nearly simultaneous (within 10 min of each other). More explicitly, we first find the percentile of each NOAA-15 flux measurement for the satellite's location (L, longitude, hemisphere, and direction). Next, applying the SAR technique, we predict the flux expected to be measured by the MetOp-2 satellite when it is at the same L shell but at its measured longitude, hemisphere and direction. We make GREEN ET AL. the prediction based on the percentiles obtained from NOAA-15 and the cumulative distributions of Me-tOp-2. If the SAR mapping is perfect, the two values will be identical. The difference between the mapped NOAA 15 data and measured MetOp-2 data is shown in Figure 12. In this plot, all the data would lie along a straight line if the mapping had no uncertainty. We use the median symmetric accuracy (Morley et al., 2018; shown in the legend as "acc") to quantify the differences. Here, the median symmetric accuracy ranges from a low of 12% at L = 1 to a high of 25% at L = 8. In this case, lower accuracy values indicate a smaller difference between the measured and mapped data. We emphasize that, across all L values shown, the percent errors are low and there is a strong 1-to-1 correspondence between the fluxes measured at Me-tOp-2 and the values mapped based on NOAA 15.
To give these plots and accuracies context, we consider these same differences between the observations from the two satellites but now with no SAR mapping (Figure 13). At all L shells, the median symmetric accuracy is improved by using the SAR mapping. The most significant improvement is at low (L = 2) and high (L = 8) L shells likely due to the fact that the longitudinal and hemispheric flux differences are most significant at these L shells.

Summary
The POES/MetOp satellite constellation and its particle detectors have been in operation for over 40 years and provide an invaluable long-term record of the space environment. The data have been and are still widely used by the science community. They are less commonly used to support near real time space weather GREEN ET AL.
10.1029/2020SW002638 12 of 15 monitoring and forecasting, in part, because the data are difficult to interpret on short suborbital (<∼90 min) timescales that are important for space weather. The large orbital flux variations that occur as the satellites move through different geographic locations may obscure physical variations of interest for space weather applications. The work presented here demonstrates and validates the use of the SAR technique for removing these orbital variations allowing the data to be used at a higher time cadence. Visual inspection as well as more detailed analysis of the resultant data set after applying the SAR mapping technique shows that the method effectively removes orbital flux variations as expected. Fourier analysis of the data before and after the technique is applied shows that the amplitude of the daily oscillations in the flux induced by the orbital variations are completely eliminated or significantly reduced. Direct comparison of SAR mapped data with observed data at similar locations and times shows that the method is highly accurate. In this work, we focused on the high energy electron flux measured by the 90° MEPED detector to demonstrate the value of the technique. The data from this detector are only a limited portion of the available measurements. The technique could also be applied to the high energy proton fluxes as well as those measured with the 0° detector. Incorporating this new form of the data into space weather tools such as the SHELLS model will be the focus of future work.

Appendix A
One requirement for the SAR method to be valid is that the relationship between the fluxes at different mapping locations is monotonic. In order to confirm that this requirement is satisfied we analyzed scatter plots of the fluxes measured by the NOAA 15 and MetOp-2 satellites when they are at different locations but similar times (within 1 day). Here we show a sample of those plots to demonstrate the monotonic relationship. Figure A1 shows a series of scatter plots comparing the >40 keV fluxes measured by the NOAA 15 and MetOp-2 satellites at the same L shell (L = 4) and Kp level (0-2) but different longitudes. Here the MetOp-2 data are measured at the reference longitude (20-30) and the NOAA-15 data are measured at eight different longitudes. One challenge to making such scatter plots is that the satellites pass through these longitudes only once per day. Thus, we compare data from the satellites measured on the same day from 2014-2016. Some scatter may be introduced because the measurements are not simultaneous.
GREEN ET AL.

10.1029/2020SW002638
13 of 15 Figure A1. Comparison between the >40 keV electron flux measured by the NOAA-15 and MetOp-2 satellites at different locations. The MetOp-2 data are measured at the reference longitude (20-30), hemisphere (Southern) and direction (Southbound). The NOAA-15 data are measured at eight longitudes in the Southern hemisphere and the Southbound direction as given in the title above each plot. Figure A2. Comparison between the >40 keV electron flux measured by the NOAA-15 and MetOp-2 satellites at the same L shells (L = 2, 4, and 6) but different longitudes. The MetOp-2 data are measured at the reference longitude (20-30), hemisphere (Southern) and direction (Southbound). The NOAA-15 data are shown at three longitudes (0,120, and 140) in the Southern hemisphere and the Southbound direction. The longitudes and L shells are given in the title of each plot. Figure A3. Comparison between the electron flux measured by the NOAA-15 and MetOp-2 satellites at the same L shell (L = 4) for three energies (e1=>40 keV, e2=>130, and e3=>287) and different longitudes. The MetOp-2 data are measured at the reference longitude (20-30), hemisphere (Southern) and direction (Southbound). The NOAA-15 data are shown at three longitudes (0, 120, and 140) in the Southern hemisphere and the Southbound direction. The longitudes and energy channel are given in the title of each plot. Figure A2 shows a similar comparison but now at different L shells. Here the scatter increases at larger L shells likely due to time changes that occur on shorter timescales than 1 day.
Lastly, Figure A3 shows how the monotonic relationship holds for different energies.

Data Availability Statements
The original POES/MetOp particle data are available at https://www.ngdc.noaa.gov/stp/satellite/poes/dataaccess.html and is described by Green (2013). The python code used to process that data and its documentation is available at https://github.com/poes-metop-sem/poes_metop_sem.