Determination of errors in derived magnetic field directions in geosynchronous orbit : results from a statistical approach

This study aims to statistically estimate the errors in local magnetic field directions that are derived from electron directional distributions measured by Los Alamos National Laboratory geosynchronous (LANL GEO) satellites. First, by comparing derived and measured magnetic field directions along the GEO orbit to those calculated from three selected empirical global magnetic field models (including a static Olson and Pfitzer 1977 quiet magnetic field model, a simple dynamic Tsyganenko 1989 model, and a sophisticated dynamic Tsyganenko 2001 storm model), it is shown that the errors in both derived and modeled directions are at least comparable. Second, using a newly developed proxy method as well as comparing results from empirical models, we are able to provide for the first time circumstantial evidence showing that derived magnetic field directions should statistically match the real magnetic directions better, with averaged errors <∼ 2, than those from the three empirical models with averaged errors >∼ 5. In addition, our results suggest that the errors in derived magnetic field directions do not depend much on magnetospheric activity, in contrast to the empirical field models. Finally, as applications of the above conclusions, we show examples of electron pitch angle distributions observed by LANL GEO and also take the derived magnetic field directions as the real ones so as to test the performance of empirical field models along the GEO orbits, with results suggesting dependence on solar cycles as well as satellite locations. This study demonstrates the validity and value of the method that infers local magnetic field directions from particle spin-resolved distributions.


Introduction
It is well-known that energetic electrons in the Earth's outer radiation belt -ranging from ∼ 3 to 8 Earth radii (R E ) -are highly dynamic and present storm-specific behaviors (e.g.Reeves et al., 2003;Chen et al., 2007b;Tu et al., 2014).Thus, monitoring, understanding, and forecasting the variations of outer-belt electrons are central topics for the space weather community.To address these topics, one basic imperative is to have long-term continuous observations with high quality and good coverage over key areas, particularly regions close to the low-altitude boundary (i.e., the lower thermosphere and mesosphere where originally trapped electrons precipitate), the internal plasma boundary (i.e., the plasmapause where wave-electron resonance prevails), as well as the high-altitude boundary (i.e., the magnetopause separating the enclosed drift shells from open ones).Among those regions, satellites in the geosynchronous orbit (GEO, a geo-equatorial circular orbit with geocentric distance of ∼ 6.6 R E ) play a unique role by monitoring the corridor through which substorm particles are injected into the inner magnetosphere, while radiation belt electrons can also be diffused outward towards the magnetopause.
Los Alamos National Laboratory has a long history of flying particle instruments aboard its geosynchronous satellites (LANL GEO hereinafter) to monitor the space environment since 1976.These instruments sample energetic electrons and protons from near the magnetic equator (Fig. 1a) over a wide range of energies, and the electron data used in this work are from the Synchronous Orbit Particle Analysis (SOPA) (Belian et al., 1992) as well as the Energy Spectrometer for Particles (ESP) (Meier et al., 1996) instruments.By themselves or in combination with others, LANL GEO particle data sets have been widely used in numerous studies leading to many significant discoveries, including identify-Y.Chen et al.: Determination of errors in derived magnetic field directions in geosynchronous orbit ing relativistic electrons as the cause of satellite deep dielectric charging (Baker et al., 1987), revealing the modulation of outer-belt electrons by solar cycle (Belian et al., 1996) and solar wind conditions (Li et al., 2005), and demonstrating the dominance of wave-particle resonance in accelerating outer-belt electrons (Chen et al., 2007a), among others.Nowadays, LANL GEO satellites provide critical complementary observations to the Van Allen Probes mission that operates inside of GEO; and in the foreseeable future, LANL GEO data sets will continue to play an irreplaceable role in scientific research as well as operational applications -such as the Dynamic Radiation Environment Assimilation Model (DREAM) (Reeves et al., 2012) -due to their long-term continuity, reliability, and high quality.
Besides resolving energy, SOPA and ESP instruments also measure particle directional distributions (Fig. 1b).SOPA's three telescopes are mounted to have different angles with respect to a satellite's spin axis (always pointing toward the Earth's center).This configuration allows each telescope to sweep out a band of the surrounding space within each spin period (∼ 10 s), and different pointing directions make each telescope sample different pitch angle ranges.Since the average magnetic field direction is more or less perpendicular to the spin axis, telescopes T 1 and T 3 will usually not be able to measure electrons near the loss cone (aligned with the magnetic field direction) as the example distributions in Fig. 1c show.Thus, measurements from all telescopes form a spinresolved distribution for each energy channel.For higher energies, ESP has a single telescope that points perpendicular to the spin axis and provides additional directional measurements.However, without a magnetometer on board, extra measures are needed to convert the directional distribution from SOPA and ESP into a more useful pitch angle distribution that is often used to characterize radiation belt dynamics (e.g., see the introduction and references in Chen et al., 2014).
Besides turning to empirical magnetic field models, one may also derive the local magnetic field direction using a physics-based technique that is first proposed by Thomsen et al. (1996) and applied to Magnetospheric Plasma Analyzer (MPA) data.This technique takes advantage of the fact that trapped-particle directional distributions should be gyrotropic, i.e., rotationally symmetric around the magnetic field line, as well as symmetric about the 90 • pitch angle.Thus, applying a principal-axis analysis to the MPA plasma directional distributions, one may generate three eigenvalues and eigenvectors and choose the most unique eigenvector as the one that is parallel or antiparallel to the local magnetic field direction (see Thomsen et al., 1996, or Chen et al., 2005, for detailed descriptions of the algorithm).Considering the fact that the low-energy plasma is often nearly isotropic on the nightside during substorm injections (e.g., Meredith et al., 1999), a related technique is developed to apply to the spin-resolved energetic electron distributions measured by SOPA and ESP.In the same vein, this technique searches for the symmetric direction in particle distributions, and details can be found in the Appendix.One example magnetic field direction derived using this technique is marked in Fig. 1c.
This study focuses on the latest method of using SOPA and ESP measurements, and testing the MPA method is left to the future (further discussions on this can be found in the Appendix).
Although the theoretical basis is solid for the above derivation technique, determining the errors associated with this technique is still a critical issue.This work aims to address this issue through estimating the errors in a statistical manner.For the first time, we provide answers to the following questions: -Does this technique outperform empirical magnetic field models?
-How large can its errors be?
-And do the errors depend on geomagnetic activity?
The cartoons in Fig. 1 illustrate the difficulty and our solution for this study.Ideally, for any given instant in time, if we were able to have all three magnetic field directions available, including the B d derived from particle distribution, the B m calculated from an empirical model, and the "real" magnetic direction B r from an in situ measurement, they would usually point in different directions (panel b).If plotting those directions inside a unit sphere as in panel d, the three points form a polar triangle with each of the side lengths proportional to the angles between each pair of unit vectors.This way, we may simply compare the length a of the side B d B r to the length b of the side B m B r to draw a conclusion.Unfortunately, in our case, the main barrier is the unknown position of B r due to the lack of in situ magnetic field measurements, and thus both values of a and b in panel d are undetermined.To overcome the barrier, we replace individual directions with statistical averages, assuming similar statistical distributions and average values for neighboring satellites, and use a triangulation method to determine the location of B r .That is, starting from two points with positions known, we first calculate their distances to B r using statistical averages from other resources; then, we draw a circle around each of the two points with a radius of the calculated average, and the intersection of circles will reveal the position of B r .In addition, since the angles between magnetic directions are mostly smaller than 10 • , we use planar triangle DRM to approximate the spherical triangle B d B r B m (panel e), which brings an ignorable error <∼ 0.5 %.Essentially, our primary goal in this work is to determine the position of R and then the length of DR.More details will be discussed in Sect.3. Hereinafter, B d , B r , and B m always refer to the statistically averaged directions of derived, real, and modeled magnetic field (i.e., unit vectors), respectively, unless being specified otherwise, and they are often shortened to D, R, and M in triangulation plots.
Instrument descriptions, data, and magnetic field models are presented in Sect. 2. Section 3 explains the statistical approaches to estimate errors in derived magnetic directions.
Section 4 discusses how to understand the results within context and their applications, and this report is concluded by a summary in Sect. 5.
2 Resources: instruments, data, and empirical global magnetic field models As mentioned in Sect. 1, local magnetic field directions are derived every 4 min from spin-resolved electron measurements from each LANL GEO satellite using the technique described in the Appendix.To get B d in this work, long-term LANL GEO data sets are used, ranging over 1996-2004 from seven satellites (1989-046, 1990-095, 1991-080, 1994-084, LANL-97A, LANL-01A, and LANL-02A) distributed globally with different geographic longitudes.
The only real magnetic field directions used in this work are from in situ measurements by several NOAA Geostationary Operational Environmental Satellites (GOES).The three-axis fluxgate magnetometers, located on a boom 3 m away from the main body of each GOES satellite, provide the magnitude and direction of the local magnetic field with a 0.512 s time resolution (Singer et al., 1996).To get B r in this work, GOES data are downloaded from the Coordinated Data Analysis Web (CDAWeb), including from GOES-08, 09 (in 1995 and1997) and GOES-10, -12 (2004).After removing the offsets in GOES data (Tsyganenko et al., 2003;Chen et al., 2005), the downloaded 1 min resolved GOES data are rebinned to 4 min to match LANL GEO data.Generally there are two GOES satellites in operation simultaneously: one at ∼ 285 • and the other at ∼ 225 • longitude.Occasionally data are available with longitudinal separations smaller than ∼ 60 • when a third GOES satellite is being activated or changing station.
For comparisons, we calculate local magnetic field directions from empirical models.We always use the International Geomagnetic Reference Field (IGRF) as the internal model, and for the external field we use three empirical models: a static model -the quiet Olson and Pfitzer magnetic field model (OP77) (Olson and Pfitzer, 1977); a simple dynamic and Kp-driven model -the Tsyganenko 1989 model (T89) (Tsyganenko, 1989); and a much more sophisticated dynamic model driven by the Disturbance Storm-Time Index (Dst) and solar wind parameters (including the pressure, interplanetary magnetic field y and z components, and interplanetary indices G 2 and G 3 ) -the Tsyganenko 2001 storm model (T01s, also called TSK03) (Tsyganenko et al., 2003).Our selection of these models is based upon previous studies (e.g., Chen et al., 2007b;Huang et al., 2008;McCollough et al., 2008), with the expectation of differing performance and the best performance of T01s from the model list in Chen et al. (2007b).It should be mentioned that we do recognize the existence of other magnetic field models (e.g., the more recent TS05 by Tsyganenko and Sitnov, 2005), and we will have more discussions on this in Sect. 4. 3 Error estimation in derived magnetic field directions using statistical approaches In this section, we focus on data in 2004 considering the simultaneous data coverage from a LANL GEO satellite 1991-080 and a NOAA satellite GOES-10.During this year, 1991-080 is ∼ 30 • west of GOES-10.Here we show both individual data examples and their statistical distributions.
Figure 2 presents an 8-day period with one major storm (minimum Dst ∼ −150 nT as in the last panel) for a glimpse of how the data, derivation, and model results compare.Panel a plots the time series of polar angles in the geocentric solar magnetospheric (GSM) system for magnetic field directions derived from 1991-080 particle distributions in comparison to polar angles of T01s model outputs.In the same format, panel b plots polar angles for measured magnetic field directions by GOES-10 in comparison to those from T01s model.Panel c depicts the angles between derived and model direc-tions for 1991-080, while panel d presents angles between real and model directions for GOES-10.Comparing panels a to b and c to d, one can see the similarities between LANL GEO and GOES data sets, such as the diurnal variations and large deviations in storm main phase.Clearly, angles in panel c and d are smaller in dayside than nightside for each satellite (a spatial feature), while angles increase significantly and simultaneously for both satellites during active times (a temporal feature).
Figure 3 presents statistical distributions of angles between magnetic field directions.Panels in the top row present deviation angles between derived and modeled field directions.As in panel a1, the mean deviation angle for T01s model has a value of 4.88 • that is the line segment length between D and M in Fig. 1e.Besides the mean values, distributions show that more than 90 % of the angles are below 10 • for dynamic magnetic field models (panel a1), while a small portion has large angle values as the long tail of the distribution in panel b1.In general, the mean angle values get smaller (with a minimum ∼ 2 • ) in dayside and larger in nightside (with a maximum ∼ 7 • for T01s), and the sizes of error bars determined from root mean squares are comparable to the mean values (panel c1).From panels in the middle row comparing measured and modeled field directions, we see similar distributions, while the mean deviation angles have slightly smaller values (panel a2) and higher percentages for low angle values (panel b2).Here the mean deviation angle for T01s has a value of 3.81 • that is the line segment length between R and M in Fig. 1e.We should note that a larger value of DM than RM does not necessarily indicate a large value of DR.
When further binned to magnetic indices, deviation angle values increase with increasing magnetic activity level, as shown by panels in the bottom row.It is interesting to see that the DM (black) and RM (red) curves trace each other very closely, and their separations are almost independent of the activity index, except for the highly active categories for which data sample numbers are too small (< 100) to make statistically significant.Results from all three magnetic field models show a similar closeness between DM and RM (not shown here), leading us to the hypothesis that the dependence of deviation angles on magnetic activities is merely caused by the degrading performance of each empirical field model, and the barely changing separations between DM and RM suggests small values for DR all the time.This hypothesis will be addressed next.

Determining the range of DR
First, before applying the triangulation method, we prove that relative positions between two magnetic field vectors have a weak azimuthal preference.As in Fig. 4a, all GOES-10 and GOES-12 data (we include two satellites for better statistics) in 2004 are plotted against model directions from T01s in a coordinate system, in which the z axis (pointing out of pa-  per) is always the local field direction calculated from the model (M), the x axis is in the z-x GSM plane and points to the Sun, and the y axis completes the right-handed orthogonal set.Thus, the position of each data bin is determined by its distance to the origin M, i.e., the deviation angle between real field direction (R) and modeled field direction (M), as well as the azimuthal angle of R with respect to the x axis.
The color in each bin indicates the count of data points (distributions with deviation angles > 20 • are not plotted here), the red circle plots the mean of all deviation angles, and the white curve shows the directional mean of deviation angles in each radial direction.Although data samples are highly unevenly distributed azimuthally, the directional mean values are still very close to the mean of all with an average absolute fluctuation level of ∼ 11 %.Therefore, we conclude that, given a statistically averaged distance of RM, we may draw a circle around the point M for all possible positions for the point R, whose exact location is, however, undetermined unless additional information is provided.Similarly, the distribution comparing the M and D from two LANL satellites (1991-080 and LANL-02a) in Fig. 4b also shows no significant azimuthal preference with an average absolute fluctuation level of 6 %.Thus, we assume that there is weak azimuthal preference for any pair of two directions in this study.
Then we apply the data analysis method aforementioned in Sect. 1 to both LANL GEO and NOAA GOES data in 2004 to estimate the range of deviation angles between D and R. First, based on the comparison of 1991-090 data in 2004 to T01s model, we draw a line connecting D and M in Fig. 4c with the segment length of 4.88 in between, which is the average value for deviation angles between the two directions as discussed in the beginning of this section.(Hereinafter all length values between two points have the unit of degree.)Then, based on analysis of GOES-10 data in 2004, a half circle is drawn around M (the lower half can be omitted due to symmetry) with a radius of 3.81 -the average value between M and R calculated from above.Point R can be anywhere on this circle, from which we estimate the median (minimum, maximum) angle between the derived direction D, and the real direction R is 6.19 (DR min = 1.07,DR max = 8.69).That is, the averaged deviation angles between D and R range within [0.28, 2.28] times of DM with a median value of 1.62×DM.Therefore, at least we can first conclude that the errors between derived and real magnetic directions are comparable to that between model and real directions.However, to further locate the exact position of R, an extra point (e.g., the imaginary point P in Fig. 4c) as well as its distance to R is needed for triangulation.

Locating point R using proxy magnetic field
To add an extra point to the construction diagram as in Fig. 4c, we developed a proxy method which approximates the real magnetic field direction for a satellite using measurements from a neighboring satellite.The proxy is derived using the equation B pt − B mt = B rs − B ms , where B pt is proxy magnetic vector for the target satellite, B rs is the real magnetic field from a neighboring (source) satellite, and B mt (B ms ) is the magnetic vector calculated from an empirical magnetic field model (T01s is used here) for the target (source) satellite, and all vectors vary with time.Since devi- ations in the modeled magnetic field are from both temporal and spatial features, the above equation assumes that the deviations in two neighboring satellites are homogeneous due to their proximity.Obviously, the validity of this assumption degrades with increasing longitude separation between two satellites.
Here, we validate this proxy method using a pair of GOES satellites when they are close enough and in situ magnetic

Derived T01s
Angle in b/t (deg)

Date (MMDD in 2004)
1991-080: derived magnetic directions, model directions vs. proxy directions Figure 6.Time series of deviation angles between derived and proxy magnetic field directions (red) and deviation angles between model and proxy directions (blue).This example covers the same 8-day period in 2004 as in Fig. 2, which includes an intense storm with the minimum Dst ∼ −150 nT on 22 January.The proxy magnetic field for LANL-GEO 1991-080 is derived from in situ measurements of NOAA GOES-10 with a ∼ 30 • longitude separation.
field data are available for both.As mentioned, GOES satellites generally have a large longitude separation of ∼ 60 • , but this separation can be smaller when a GOES satellite is relocated, although observation data during those periods are rarely available.We were fortunate enough to identify a short period with available data in 1995 when GOES-09 was moved from longitude 270 to 244 • .This movement makes the longitude separation between GOES-09 and GOES-08 increase from initially ∼ 15 to ∼ 40 • .Therefore, after applying the above equation to approximating GOES-09 magnetic field using GOES-08 measurements, proxy magnetic field directions are validated by GOES-09 measurements, as the green curve in Fig. 5a.For comparison, deviation angles between GOES-09 measurements and T01s model are also plotted.It is clear that the proxy outperforms the T01s model significantly when the longitude separation between satellites is <∼ 30 • , and both perform similarly even when the separation goes beyond ∼ 40 • by the end of the period.Therefore, since the longitude separation between GOES-10 and 1991-080 is ∼ 30 • in 2004, this proxy method can add the point P to the plot by using GOES-10 to derive proxy for 1991-080.First, locating point P requires knowing the lengths of DP and MP .Therefore, the proxy field directions are compared to derived and model field directions for 1991-080, and Fig. 6 presents a short interval as an example.A statistical study gives out an averaged DP value of 5.34 and an MP value of 4.11.Thus, we are able to plot the point P for proxy directions in Fig. 7b.
Then we need to derive the value of P R for the circle radius.Besides the above validation using the pair of GOES-08, -09 in December 1995, we also use data from two other periods: GOES-08, -09 between 20 November and 1 December 1997, and GOES-10, -12 between 1 March and 1 April 2004, both with a longitude separation of 60 • .Deviation angles between derived and proxy field directions for all three periods are plotted against longitude separation in Fig. 7a, overplotted by averaged angle values ranging from ∼ 2 to 5.Here we use the average value of 3.50 for the seg-ment length P R. As in Fig. 7b, the circle P with a P R radius of 3.50 intersects with the circle M T01s (with a radius of 3.81): the intersection point R 1 has a distance of 1.83 from D, and the intersection point R 2 has a distance of 7.71.So now the question is which point -R 1 or R 2 or both -is real.
To answer the question, we replace the T01s model with the OP77 model and repeat all of the steps above.As in Fig. 7c, we have different values of DM, MP , and MR due to the different model but the same values of DP and P R, and again there are two intersection R points.However, DR 1 in both panels b and c has the same values but the DR 2 values are different, which serves as the first piece of evidence that R 1 should be the real R point since we do expect the DR values to be independent of empirical magnetic field models.

Locating point R from grouping of points
Inspired by the proxy point P added in Fig. 7, we speculate that an alternative way of using two empirical models should also be able to add an extra point.As in Fig. 8a, after using the T01s model to place the baseline DM T01s and drawing the circle M T01s , an extra point M op77 from the model OP77 can be located from the segment lengths of DM op77 and M op77 M T01s using 1991-080 data.Then the second circle M OP77 is drawn with the radius of M OP77 R determined from GOES-10 data.Again the two circles have two intersection points: the R 1 point with a distance of 1.10 from point D, and R 2 with a distance of 8.73.To differentiate R 1 and R 2 , we introduced the second extra point M T89 using the model T89, whose position is located from the segment lengths of DM T89 and M T89 M T01s using 1991-080 data (panel b).And the third circle M T89 , with a radius of 3.81 from GOES-10 data, intersects with the circles M T01s ( M OP77 ) at points R 1a and R 2a (R 1b and R 2b ).In the ideal case, R 1 , R 1a , and R 1b should overlap (the same is true for R 2 , R 2a , and R 2b ), though it is natural to see that they do not do so exactly since statistically averaged values are used here.However, points R 1 , R 1a , and R 1b in panel b are indeed tightly clustered but not points R 2 , R 2a , and R 2b , which serves as the second piece of evidence that R 1 points should be very close the real position of the R point, instead of the widely spreading points R 2 , R 2a , and R 2b .
Based on this analysis, we conclude that in an average sense the derived magnetic field directions are closer to the real magnetic field than simulations from the three selected empirical field models used in this work.Although the DR 1 values are not the same in Fig. 7 (1.83) and Fig. 8 (1.10), this can be explained by the uncertainties in the numbers used here.For example, the 3.50 for P R used in Fig. 7 may have larger errors than other numbers due to the limited available data.However, from both DR 1 values, it is reasonable to state that the average deviation angle between derived and real magnetic field directions is smaller than the value between model and real directions by a factor of >∼ 2.  (from the T89 model) and its circle generate R 1a and R 1b points very close to the R 1 , as well as R 2a and R 2b points spreading away from R 2 .(c) R 1 and R 2 points are determined for different magnetic activity categories.Since the DR m1 and DR m2 basically stay constant with magnetic activity, the grouping of R 1 points should be much tighter than that of R 2 points.
How the DR value varies with magnetic activities can be learned in a similar way, by taking advantage of the fact that length difference between the DM and RM length stays almost unchanged with magnetic activity levels as discussed in the beginning of this section.A qualitative instead of quantitative method is employed in this step, which should guarantee that our conclusion is reliable.As in Fig. 8c, we draw a diagram using two magnetic field models (1 and 2) for two different activity categories (a and b).As just mentioned, since the distances from D to the circles around M 1a and M 1b along the DM 1a line have the constant value of DR m1 , and the circles around M 2a and M 2b (not drawn here) will both be at a tangent with the small circle D with a radius of DR m2 , we can see that the intersection points R 1a (between M 1a and M 2a ) and R 1b (between M 1b and M 2b ) stay very close to each other while R 2a and R 2b are well separated.Therefore, because we already know that the R 1 group is close to the real R point, we conclude that DR values are not sensitive to the magnetic activity levels.This supports our hypothesis that the observed increasing DM values with elevated activity levels in Figure 3 should be mainly due to the degrading performance of empirical models, as discussed in the beginning of this section.

Discussion and applications
One possible major error for this study comes from the statistical approach itself, that is, how representative the average points are in the construction plots, such as Figs. 4, 7,  and 8.For an individual case study, each point in those figures is definite and thus the triangulation method is valid.However, for two given distributions, the representativeness of the calculated mean deviation points may be questionable.Indeed, considering the variations in each distribution, the above method is only valid when the two distributions are relatively homogeneous, which again cannot be directly tested due to the lack of simultaneous derived and measured magnetic field data.Nevertheless, one indirect test can give us some indications and thus confidence for the representativeness of averages: in Fig. 8b, the distance between M T89 and M OP77 can be measured from the plot to be 5.49.Compared to the calculated value of 4.66, this indicates a ∼ 18 % error that should be acceptable.
To understand the averaged deviation angle of <∼ 2 • in derived magnetic field directions, we need to discuss what are the possible error sources for this method.The most likely error source is the large size of angular bins used to sort measured particle counts.In our case, the largest angular bin size can be ∼ 11 • so that the assigned pitch angles can have errors as large as ∼ 5.5 • .The second error source may be at times when particle distributions are close to isotropic.This can be significant for low-energy plasma particularly during substorm injections but should be alleviated for energetic radiation belt electrons (typically with several hundred keV to > 1 MeV energies like the SOPA E5 and ESP E1 channels selected for this work).For example, according to a recent pitch angle distribution (PAD) statistical study, PADs for ∼ 150 keV electrons at L ∼ 6 are statistically very close to isotropic during substorms as shown in Fig. S2b, panels A2 and B2, in the Supplement of Chen et al. (2014), while PADs for ∼ 1.5 MeV electrons at the same L are statistically highly anisotropic as in Fig. S2b, panels A2 and B2.Another possible error source is the intrinsic asymmetry in PAD due to either the statistical fluctuations in counts registered by instruments or some process that breaks down the particles' bounce movement.The former occurs when MeV electron fluxes drop significantly during storm main phases; the latter may also be possible for electrons close to the loss cone but can be ignored for stably trapped populations that make up the LANL GEO observations.All these could contribute to the small but existing errors we found here.
A direct application of the derived magnetic field direction is to sort LANL GEO particle directional measurements into PADs, as one such example shown in Fig. 9.During this double-dip storm period, substorm electron PADs in panel a vary differently from those of energetic electrons in panel b.For instance, substorm electron PADs are mainly pancakeshaped or close to isotropic during injections (e.g. at ∼ 40 and 125 h), while MeV electrons show intriguing sustained butterfly PADs in the early phase of radiation belt enhancements (e.g., throughout the day 19 March).This difference suggests that the two populations should have experienced different physical processes.Therefore, as discussed in the "Introduction" section, LANL GEO measurements have high energy and pitch angle resolutions and are distributed over multiple longitudes at GEO; thus, they are highly valuable for studying radiation belt dynamics, particularly together with simultaneous observations from Van Allen Probes inside GEO.
Additionally, since the deviation of B d is small, we may use the derived directions as real ones to test the performance of empirical models over the long term (1997)(1998)(1999)(2000)(2001)(2002)(2003)(2004). Figure 10 presents the distributions in the same format as in Fig. 3. Percentage distributions in panels a and b are similar to those in Fig. 3 except getting slightly flatter, which is consistent with the slightly increased mean values in the magnetic local time (MLT) distributions in panel c.The small spikes at noon are mainly from data before 2000, and how realistic they are will be left to future investigation by examining individual events.This larger data set allows better coverage with statistical significance extending to higher magnetic activity categories in panels d, e, and f.From lowto moderate-activity categories, dynamic models persistently perform better than the static model; however, an interesting reverse can be seen in distributions for which T01s model has the largest deviation for the very high activity range.
We further inspect the dependence of deviation angles on the solar cycle and satellite positions.As in Fig. 11a, the deviation has a general growing tendency in the rising phase of the solar cycle until reaching the maximum in ∼ 2002 and then declines afterwards.Also, with different geographic longitudes, LANL GEO satellites are located at different magnetic latitudes or equivalently at different L shells (Chen et al., 2005).By plotting the mean deviation values vs. the Lm  (McIlwain L shell) calculated from T01s model, we do see a general trend of increasing deviation values with increasing Lm by linearly fitting those data points.Indeed, the calculated Pearson's correlation coefficient has a nontrivial value of 0.41.All these suggest that the model T01s performance degrades with increasing L shells (or, latitudes), which is consistent with our general impression of empirical models.
Finally, as mentioned in Sect.2, we only chose three representative empirical magnetic field models without including the more recent sophisticated TS05 model.Although previous studies have demonstrated that T01s performs better than many other models (Chen et al., 2005;McCollough et al., 2008), no comprehensive study has been conducted to compare between T01s and TS05.Therefore, we cannot simply extend our conclusion to the TS05 model, although there are some clues suggesting comparable performances of T01s and TS05 at GEO: when statistically comparing to observations dominated by GEO data, TS05 has correlation coefficients of (0.92, 0.83, and 0.92) for magnetic field (x, y, z) components, while T01s has values of (0.91, 0.82, and 0.90) (Tsyganenko and Sitnov, 2005).We decide to leave the inclusion of the TS05 model to the future.

Summary
This work statistically estimates the errors in the local magnetic field directions derived from electrons' directional distributions measured by LANL GEO satellites.First, by comparing derived and measured magnetic field directions in GEO to outputs from empirical global magnetic field models (including a static Olson and Pfitzer quiet magnetic field model, a simple dynamic Tsyganenko 1989 model, and a sophisticated dynamic Tsyganenko 2001 storm model), we show that the errors in both derived and modeled directions are at least comparable.Second, using a newly developed proxy method as well as comparing results from multiple empirical models, we provide for the first time evidence showing that derived magnetic field directions should statistically outperform -with a ratio factor of >∼ 2 between magnetic field deviation angles -the three selected empirical models (including T01s) in matching the real magnetic directions.Additionally, our results suggest that errors in derived magnetic directions are not so much dependent on magnetospheric activities as the empirical field models.At last, after showing electron PADs observed by LANL GEO satellites, we further use the derived magnetic field directions for testing the performance of empirical field models, with results showing dependence on solar cycles as well as GEO satellite positions.This study for the first time demonstrates the validity and the value of using the symmetric nature of particle spin-resolved distributions for deriving local magnetic field directions.
6 Data availability LANL GEO data used in this study are available upon request by contacting the corresponding author Y. Chen (cheny@lanl.gov).Appendix A: Inferring magnetic field directions from LANL GEO SOPA and ESP measurements The algorithm applied to the SOPA and ESP data first bins each of the three SOPA telescopes and the lone ESP telescope into spin phase using accumulations over a 4 min window to flesh out the distribution as a function of spin phase.The count from each accumulation bin is either placed into one of 32 spin-phase bins for SOPA data or into one of 180 spin-phase bins for ESP data.Next, the spin-phase angle, φ, is found, about which the particle distribution measured by the ESP E1 (0.7-1.8 MeV) channel is most symmetric.This angle points parallel or antiparallel to the projection of the background magnetic field into the plane perpendicular to the spin axis, or, for certain particle distributions, points 90 • perpendicular to the magnetic field.These ambiguities are cleared up in the second stage of the analysis, wherein every angle, θ , measured from the spin axis, is tested in 2 • increments as a potential field line direction when combined with φ.The pair (φ, θ ) specifies a tested magnetic field direction, and the spin-resolved SOPA E5 (225-315 keV) elec-tron channel counts are binned into pitch angles under the assumption that this pair is the correct one.A smooth polynomial function is fitted to the pitch angle binned counts, and the root mean squared error (RMSE) of the fit is calculated.The pair (φ, θ ) that produces the lowest RMSE is chosen as the field direction.Because the three telescopes for SOPA may not be perfectly calibrated to one another, multiplicative constants for T 1 and T 3 are found that map the pitch angle binned counts for T 1 and T 3 so that they best match those from T 2 .This "calibration" of T 1 and T 3 is done separately for each 4 min time bin, each energy channel, and each hypothesized magnetic field direction (φ, θ ).SOPA Channel E5 was chosen to estimate the magnetic field direction because it had the best combination of anisotropy and count rate over the broadest range of conditions, but a better algorithm could be devised that analyzes all energy channels simultaneously, as in Thomsen et al. (1996), or selects the best energy at any given time.
A systematic comparison of the two methods using MPA and SOPA with ESP is outside the scope of this work; however, it would be informative to get a glimpse of how magnetic field directions from the two methods compare.Figure A1 presents one such example which compares the derived (φ, θ ) values from two methods along the orbit of LANL-01A during a 9-day period.It can be seen that in panels a and b, although directions from both methods agree well with each other mostly during quiet times, values from MPA experience large fluctuations when the satellite travels through the midnight sector (at ∼ 23:00 UT each day) with substorm injections indicated by high Auroral Electrojet Index (AE) values (e.g., on DOY 147-150), which is consistent with the discussion in Sect. 4. In comparison, directions from SOPA and ESP do not have as many large fluctuations (e.g., during the first small dip of Dst on DOY 147-148).Similar results have been seen for other LANL GEO satellites for different periods (not shown here).Indeed, if necessary, directions derived from MPA measurements can also be statistically studied using the same approach presented in the current work.

Figure 1 .
Figure 1.LANL GEO satellites measure electron directional distributions.(a) Side view of the GEO orbit.A LANL GEO satellite is usually close to but not exactly in the geomagnetic equator due to the tilted geomagnetic dipole field.(b) Rotation of the satellite platform allows the three SOPA telescopes (T 1 , 30 • to the spin axis z; T 2 , 90 • , and T 3 , 120 • ) to sample directional distributions of electrons; meanwhile, the unit local magnetic field from the empirical model (B m in blue), the one derived from electron distribution (B d in red), and the real direction (B r in black, if measured) can be different.The goal of this work is to determine the angle between B d and B r (indicated by the question mark).(c) One example electron distribution measured by SOPA.Count rates are sorted by the roll angle (defined as the azimuthal angle in the satellite spin plane: 0 • along the x (due east) direction and 90 • along y (due south)), and a derived magnetic field direction from symmetry of the distribution is marked by the white cross.The very low counts for T 2 are measured close to the loss cone.(d) In a unit sphere, the three magnetic vectors form a polar triangle B d B r B m , whose side lengths (a, b, and c) are proportional to the angles between each pair of unit vectors.(e) Polar triangle B d B r B m can be approximated by the planar triangle DRM in this study.

Figure 2 .
Figure 2. Sample magnetic field directions during an 8-day period in 2004.(a) Polar angles of derived magnetic field directions (red) from 1991-080 particle data are compared to those calculated from T01s model (blue), both plotted as a function of time.Polar angle is defined as the angle between a magnetic field direction and the z axis of GSM coordinate system.(b) Polar angles of observed magnetic field directions (black) by GOES-10 compared to those from T01s model (blue).(c) Angles between (b/t) derived and model magnetic vectors for 1991-080.Gray (black) symbols are for data in dayside (nightside).(d) Angles between measured and model vectors for GOES-10.(e) The Dst (black) and Kp (gray) indices.A major storm occurs on 22 January during the period.

Figure 3 .
Figure 3. Statistical studies comparing derived, real (measured), and model magnetic field directions in 2004.Panels in the top row are for 1991-080.(a1) Accumulative percentage vs.deviation angles between derived and modeled directions for the three empirical models T01s (black), T89 (green), and OP77 (gray).Mean angle values as well as satellite coordinates are also presented.(b1) Normalized percentages vs. deviation angles for T01s.(c1) Deviation angles are binned to MLT for the three models, and the vertical gray bars are the errors for T01s model.Panels in middle row are for GOES-10 in the same format except for comparing real and modeled directions.(a3) Deviation angles are binned to Kp for 1991-080 (black) and GOES-10 (red) using T01s model.Again the vertical bars are errors for each.The gray dotted line plots data sample number in each bin (read by the vertical axis on right).(b3) Deviation angles are binned to Dst. (c3) Deviation angles are binned to the Auroral Electrojet Index (AE).

Figure 4 .
Figure 4. Deviation distributions and estimating the deviation angle range between derived and real magnetic field directions.(a) Distributions of real directions (R) relative to model directions (marked by the white "M" in the origin pointing out of paper).The radial distance from any point to M is the deviation angle between a pair of model and real directions, and the azimuthal angle is determined in a modified local B-GSM coordinate system (and thus is not local time).Color in each bin indicates the count of data points.The overplotted white curve indicates the directional mean of deviation angles in each radial direction, compared to the red circle showing the mean of all deviation angles.(b) Distributions of derived directions (D) relative to model directions (M), directional means, and the mean circle in the same format.(c) Given averaged deviation angle values for DM (4.88 • ) and RM (3.81 • ) in Fig.3, we may estimate that the range of DR is between [1.07 • , 8.69 • ], that is [0.28, 2.28]×DM.The imaginary point P and circle in gray, if available, will help pinpoint the position of R.

Figure 5 .
Figure 5. Validating the proxy method of using measurements from a neighboring satellite.(a) In this 1-month period, deviation angles between the proxy magnetic field direction and in situ measurements (green) along the GOES-09 orbit are plotted as a function of time, compared to angles between T01s model and measurements (blue).(b) During the period, the relocation of GOES-09 makes its longitude separation with GOES-8 varying from ∼ 15 • to up to 40 • .(c) Dst (black) and Kp (gray) indices.Minor and moderate magnetic activity is observed during the period.

YFigure 7 .
Figure 7. Determining the position of R point(s) using proxy magnetic field.(a) Deviation angles between proxy and measured field directions, in three selected periods, are plotted against the longitude separations between each pair of GOES satellites.Overplotted data symbols are averaged angles for binned longitude separations.(b) The introduction of the point P and the P R circle generate R 1 and R 2 intersection points when T01s model is used.(c) The introduction of the point P and its circle generate another pair of R 1 and R 2 intersection points when OP77 model is used.

Figure 8 .
Figure 8. Determining the position of R point(s) using multiple empirical models.(a) A circle is drawn around each of the points M T01s and M OP77 from the two models, and the intersections give two candidate points R 1 and R 2 .(b) The introduction of another point M T89(from the T89 model) and its circle generate R 1a and R 1b points very close to the R 1 , as well as R 2a and R 2b points spreading away from R 2 .(c) R 1 and R 2 points are determined for different magnetic activity categories.Since the DR m1 and DR m2 basically stay constant with magnetic activity, the grouping of R 1 points should be much tighter than that of R 2 points.

Figure 9 .
Figure 9. Electron PADs -based upon the derived magnetic field directions -observed by LANL-01A SOPA during a geomagnetic storm period (7 days).(a) Pitch-angle-resolved fluxes for low-energy (131 keV) electrons evolve with time.(b) Pitch-angle-resolved fluxes for high-energy (1.2 MeV) electrons evolve with time.(c) Dst (black) and Auroral Electrojet (AE) (gray) indices for the period.The time bin size for each PAD is 4 min.LANL-01A reaches the noon local time position at ∼ 23:00 UT each day during this period.

Figure 10 .Figure 11 .
Figure 10.Results from comparing derived and model magnetic field directions for all available LANL GEO data within 1997-2004.Panels have the same format as in Fig. 3. (a) Accumulative percentage vs.deviation angles between derived and simulated directions for three empirical models T01s (black), T89 (green), and OP77 (gray).Mean angle values are also presented.(b) Normalized percentage vs.deviation angles for T01s model.(c) Mean deviation angles are binned to MLT for the three models, and the vertical gray bars are the error ranges for T01s model.(d) Deviation angles are binned to Kp ranges using different colors for three models.Again, the vertical bars are error bars for T01s model.The gray dotted line plots data sample numbers (read by the vertical axis in right).(e) Deviation angles are binned to Dst.(f) Deviation angles are binned to the Auroral Electrojet Index (AE).

YFigure A1 .
Figure A1.Magnetic field directions derived from MPA aboard LANL-01A compared to those from SOPA and ESP during a 9-day period.(a) The time series of magnetic field directions' polar angle theta (θ ) derived from MPA (red) are compared to those from SOPA and ESP (blue).LANL-01A crosses midnight at ∼ 23:00 UT each day in 2005.(b) Time series of field directions' azimuthal angle Phi (φ).(c) Dst (black) and Auroral Electrojet (AE) (gray) indices for the period.A double-dip storm occurred during the period with the minimum Dst ∼ −135 nT reached on 30 May (DOY 149).