Interstellar flow longitude from pickup ion cut-off observations at 1 AU with STEREO and ACE

A precision determination of the interstellar neutral (ISN) flow direction is important in several ways. As a cardinal axis of the heliosphere it has strong leverage on the plane subtended by the ISN velocity and interstellar magnetic field vector, which controls the shape of the heliosphere and its interaction with the interstellar medium. Continuing observations of the ISN flow through the heliosphere for several decades allow the search for potential temporal variations of the ISN flow and comparison with astronomical observations. Recent efforts to obtain a consistent ISN vector and temperature with Ulysses and IBEX neutral gas observations point to remaining uncertainties and potential systematic effects. In particular, IBEX measurements provide a very precise relation between ISN flow longitude and speed via the hyperbolic trajectory equation, but they contain larger uncertainties along the parameter tube defined by this relation. The pickup ion (PUI) cut-off variation with ecliptic longitude at 1 AU can provide a complementary determination of the ISN flow longitude with high precision. We compare STEREO PLASTIC and ACE SWICS observations with a simple analytical model of the cut-off. We perform a Pearson correlation analysis of the cut-off as a function of ecliptic longitude with its mirrored function and obtain the symmetry axis with a statistical uncertainty <0.1°. Here we test variations of this value due to Poisson fluctuations in the original data with simulations and due to systematic effects with multi year and location ACE SWICS and STEREO PLASTIC samples.


Introduction
The interstellar neutral (ISN) velocity vector V ISN∞ relative to the Sun is one of the key quantities that control the global interaction between the interstellar medium and the heliosphere. Together with the interstellar magnetic field B IS it defines the B IS -V ISN∞ plane, which determines the shape of, and the ISN flow deflection in, the outer heliosheath [1,2]. B IS is deduced from the IBEX ribbon, consistent with the orientation based on the heliospheric asymmetry [3,4,5] and TeV cosmic ray anisotropy [6]. After a long history of studying the local interstellar parameters with UV backscatter observations [7,8], pickup ions (PUIs) [9,10], and finally direct observations of ISN He with Ulysses GAS [11], the Interstellar Boundary Explorer (IBEX) provides a more complete view of the H, He, O, and Ne ISN flow [12]. A precise ISN flow vector taken over an extended time of more than one decade will also enable detection of small temporal variations in the ISN flow or place tight upper limits on them, and thus resolve a current debate [13,14,15]. In addition, such local measurements over time, in comparison with astronomical absorption line observations of nearby interstellar clouds [16], will facilitate the identification of local structures, such as turbulent eddies, in the local interstellar cloud.
IBEX returns a very precise relation between ISN flow longitude λ ISN∞ and V ISN∞ via the hyperbolic trajectory equation, but leaves a much larger uncertainty separately for λΙ SN∞ and V ISN∞ . These results define a 4-dimensional parameter tube [17] connecting flow longitude, latitude, speed, and temperature, which have been consolidated recently [18,19,20,21]. Related values and uncertainties for λ ISN∞ are 75.9±0.5 o [18] and 75.5±1.4 o [21], respectively, in agreement with most recent analyses of Ulysses GAS data that have led to λΙ SN∞ = 75.3-1.1+1.2 o [22] and 75.54±0.19 o [23]. However, these results lead to a much higher ISN temperature than deduced previously in [24]. While all results now agree within their mutual uncertainties, the best-fit values still vary noticeably. An independent determination of λ ISN∞ will substantially tighten the flow vector in combination with further narrowing the IBEX parameter tube through a growing database. The analysis of the variation of the pickup ion (PUI) energy cut-off as a function of ecliptic longitude, which reflects the ISN flow pattern at 1 AU, has recently been suggested as the appropriate precision method [25].
Reporting on this new PUI cut-off method is a very appropriate contribution for the celebration of Ed Stone's 80 th birthday. This topic touches upon two of his most remarkable achievements, the Voyager interstellar mission [26], which continues to pose challenges to our understanding of the connection between the interstellar flow and magnetic field at the heliospheric boundary, and the ACE mission [27], which has provided the legacy PUI data set with the SWICS sensor [28] that is being used together with continuing STEREO PLASTIC [29] observations to test for systematic effects and to evaluate long term variations in the ISN flow direction.
In this paper, we will take a first look at the effects that statistical and systematic uncertainties have on the ISN flow longitude as deduced from the PUI cut-off. We will start by briefly laying out the method in Section 2. Will evaluate in Section 3 how Poissonian noise added to an artificial STEREO PLASTIC-like dataset affects the statistical uncertainty of the result. In Section 4 we will analyze several individual annual data sets of ACE SWICS and STEREO A PLASTIC, to determine average values and standard errors, which we then compare with the results of the entire multi-year data set. Section 5 will contain a summary, a discussion of the implications, and suggested future work.

Determination of the ISN flow longitude from the PUI cut-off
The ISN flow longitude can be derived independently from characteristic spatial structures in the inner heliosphere, the so-called focusing cone [30,31] and crescent, both identifiable with PUI observations at 1 AU [32]. Although a careful statistical study largely eliminated the influence of solar wind structures and potential sensor efficiency variations [32], PUI flux observations likely retain considerable variations relative to the parent neutral gas structures due to solar wind structures [33], temporal and spatial variations in the ionization rates, and PUI transport effects [34,35]. A more robust method uses the fact that the PUI cut-off speed (the high-energy edge of the relatively flat PUI velocity distribution [9]) is a function of the radial ISN flow component V r and the solar wind speed V SW and thus varies with the observer longitude λ Obs [36]. In its dimensionless form v r = V r /V E , normalized to the Earth's speed at 1 AU, the radial ISN flow varies with λ as [25] where λ = λ Obsλ ISN∞ -180 o and λ ISN∞ is the ISN flow direction outside the heliosphere. v r reaches its maximum exactly upwind, i.e. for λ Obs = λ ISN∞ + 180 o where λ = 0. v r decreases symmetrically with increasing angle |λ| and v r = 0 holds at ±λ 0 , where the ISN bulk flow reaches its perihelion at 1 AU [37,38] with cos λ 0 = -(1 + v ISN∞ 2 ) −1 . As discussed in [25] the PUI cut-off speed w' CutOff in the solar wind frame (primed variable), normalized to the solar wind speed V SW , is to a good approximation given by Note in equation (2) and in Figure 1 V r is taken positive for inward ISN flow directions. Because PUI distributions are gyrotropic and mostly scatter in pitch angle it is advantageous to accumulate their spectra on phase space surfaces constant in w'. With the energy and angular resolution of STEREO PLASTIC the He + observations taken in the spacecraft frame have been transformed into the solar wind frame to facilitate accumulation along constant w'. ACE SWICS integrates He + over its entire field-of-view in constant energy per charge in the spacecraft frame, or in w = w' + 1. As we will see, this makes the ACE SWICS observations more sensitive to variations in the interplanetary magnetic field (IMF) direction and in the solar wind speed. For consistency in the presentation, we show all PUI cut-off results in terms of w' CutOff and transform the ACE SWICS observations according to w' CutOff = w CutOff -1. To obtain the most reliable cut-off determination, the initial torus of the PUI velocity distribution, which has been clearly identified with STEREO PLASTIC observations [39], needs to be in the sensor field-of-view. To achieve this condition, STEREO PLASTIC observations have been constrained so the angle between B IMF and the solar wind lies between 70 o and 110 o . Because ACE SWICS observations have to be used in the spacecraft frame, where changes in the IMF orientation lead to variations in the cut-off, these data have been constrained to angles in the range 75 o -105 o .
To obtain a proxy value for the cut-off the He + PUI counts accumulated in a 1 o cadence in ecliptic longitude by STEREO PLASTIC are used on a fixed grid in w', normalized to the maximum count value for each 1 o bin to stay as close to the observable as possible. For ACE SWICS the normalized phase space density is used to correct for differences in duty cycle and efficiencies as a function of particle energy. A hyperbolic tangent is fitted to these data sets, as described in [25], according to (3) where f n represents the normalized counts or the normalized phase space density, respectively, and A is the steepness of the cut-off. How the proxy for the cut-off determined in this way is shifted by a constant relative to the actual cut-off of the PUI distribution for the two data sets has been discussed in [40]. They also showed that this shift does not influence the determination of the ecliptic longitude of the ISN flow because of the independence of the shift in value and the longitude dependence and because of the symmetry of the longitude variation.
The cut-off w' CutOff (λ Obs ) then undergoes a correlation analysis that invokes the mirror symmetry of the longitude variation according to equation (1). The longitude with the maximum of the Pearson correlation coefficient is identified as the ISN flow direction λ ISN , with a very small statistical fit uncertainty for the combined STEREO A PLASTIC 7-year data set that includes March 2007 through March 2014 [25]. It is evident that the small uncertainty reported in this first analysis may underestimate effects due to stochastic variations in the observations and that most systematic effects have not been considered yet. In the following, we will first obtain a more realistic estimate of the statistical uncertainties involved and then look into systematic effects with individual year ACE and STEREO data sets. Most of the important systematic effects that influence PUI observations, such as solar wind speed and flux variations, changes in IMF strength and orientation, as well as in wave power, are distributed stochastically over ecliptic longitude and thus vary in the same fashion from year to year. Therefore, multiple year compilations will see substantially reduced systematic effects.

Estimation of the statistical uncertainties
To estimate the statistical uncertainties, artificial data sets were generated and then varied in a Monte Carlo scheme with Poisson distributed counts. To construct the artificial data we followed the dependence of the observed PUI distribution as a function of w' and λ Obs . For the variation in w' the decrease near the cut-off and the decrease toward lower energy were represented by different hyperbolic tangent functions. Only the hyperbolic tangent at the high-energy cut-off matters for this analysis. For the longitude variation, a cosine dependence was chosen, with the maximum in w' CutOff at the presumed upwind direction of the ISN flow. The parent function of the artificial PUI distribution data sets is shown in a color-coded spectral representation in the left panel of figure 1. Instead of w' the radial component of the PUI injection velocity V r is shown. Poissonian noise was added to each pixel in w' and λ Obs with a random generator, according to the mean total number of counts expected over the course of one year. In this way 200 trials were prepared, one of which is shown in the right panel of figure 1.
To find the cut-off w' CutOff two different methods were used. In method I, the fit of a hyperbolic tangent to the decrease toward the high-energy end was applied as described above. In method II, a typical decrease function was determined from the data set averaged over 240 o ≤ λ Obs ≤ 270 o , which was then used in a correlation analysis for each 1 o bin, with w' varying to obtain the maximum correlation coefficient.
The results from this Monte Carlo simulation with artificial PUI distribution data, assuming λ ISN = 75 o , analyzed with both methods are compiled in Table 1. The presumed ISN flow direction is found very close to the input value. The statistical uncertainties are ±0.62 o for method I and 0.40 o for method II, both noticeably larger than the mere fit uncertainty of ±0.04 o reported in [25]. As pointed out above, in this simulation the counting statistics for 1 year was used. Considering that [25] used a 7-year data set, the uncertainties will likely be reduced by √7, which would leave the uncertainty obtained with method II at 0.15 o about a factor of 3.5 higher. Barring any systematic effects, a statistical uncertainty ≤0.15 o would then be achievable. Method II appears to be superior to the fit with a hyperbolic tangent, likely because this correlation method is less sensitive to any variations in the slope or shape in the decreasing part of the PUI distribution. None of these uncertainties contain any systematic effects, the magnitude of which we tackle next.

Systematic effects on the ISN flow longitude in multiple year observations
To obtain a first estimate of systematic effects on the ISN flow longitude derived with the mirror correlation method described in [25], we apply the same method to several individual year data sets from STEREO A PLASTIC and ACE SWICS.    It is obvious that the single year data sets show much larger variations. While this is already expected due to much reduced counting statistics, as is evident from the larger statistical error bars in the plots, there also appear to be visible systematic deviations from the rather smooth variation of the cutoff with λ Obs in the STEREO A multi-year data set. Near λ Obs = 100 o as well as around 240 o and 270 o there appear to be visible systematic deviations to smaller w' CutOff values, surrounded by somewhat larger values. These deviations are visible in both data sets, but much more pronounced in the ACE SWICS observations. We will discuss potential reasons for these variations at the end of this section. Figure 3 shows the Pearson correlation coefficient obtained by mirroring the results about a chosen ecliptic longitude λ M , as described in Section 2 and in [25] for the two 2007 data sets in figure 2. Both mirror correlations vary generally according to a cosine function with λ M and show a pronounced maximum near the ISN upwind direction. However, the maximum correlation coefficient is substantially reduced compared to the result in [25], approximately 0.5 for STEREO PLASTIC and 0.3 for ACE SWICS, and the longitudes for the maximum deviate noticeably from the value found for the combine STEREO data set. Because the original data set shown in figure 2 contains systematic deviations, such a behavior is expected. Table 2 shows the ISN flow direction λ ISN along with the fit error and the maximum correlation coefficient for each annual data set used from ACE SWICS and STEREO A PLASTIC. The third line from the bottom shows the results obtained with the combined ACE SWICS and STEREO A PLASTIC data set, for which figure 4 shows the longitude variation of the mirror correlation coefficient for the combined ACE SWICS (left panel) and combined STEREO A PLASTIC (right panel) observations. The bottom line contains the weighted mean values for both data sets, with the maximum correlation coefficient used as the weight, along with the standard deviation of the individual results and the standard error for the multi-year data sets. The year-to-year results of STEREO PLASTIC show a 2.2 o standard deviation, but a weighted mean that falls close to the value obtained with the combined data set and a standard error of 0.9 o that signals that the two different results are in agreement within the stated error bars. The annual ACE SWICS results show a much larger standard deviation of 9.

Discussion of the annual and combined results
Clearly, individual yearly data appear to be affected substantially by systematic effects, which pull the ISN flow direction obtained with the PUI cut-off method. These systematic effects are much stronger for the ACE SWICS data as can be expected from the fact that so far no provisions have been taken to adjust for the difference that ACE SWICS data are used in the spacecraft frame, while the STEREO PLASTIC data are transformed into the solar wind frame. Also the statistical fit errors and their standard counterparts for the weighted mean of the yearly results substantially understate the complete set of uncertainties involved. Systematic effects are visible also in the STEREO PLASTIC results and are definitely larger than the statistical uncertainties, solely from the fits or based on the Monte Carlo simulations in Section 3, but they are much smaller than for ACE SWICS.
Although the STEREO PLASTIC observations have been constrained to solar wind speeds <450 km/s to exclude CMEs and high speed solar wind streams for now, all PUI observations have been averaged over 1 o in ecliptic longitude, or approximately one day, while the solar wind varies on time scales faster than that. Because the observations are binned in w' the observed cut-off is modulated over ecliptic longitude in V r /V sw, thus retaining the memory of V sw variations. In addition, fast solar wind compressions and rarefactions can lead to variations of the cut-off due to adiabatic changes of the PUI distribution [41]. Obviously, the ACE SWICS observations are affected more strongly. First, the SWICS data have not been constrained in solar wind speed to maintain high counting statistics and to see effects from solar changes more clearly. Second, ACE SWICS data are not transformed into the solar wind frame. CMEs and fast solar wind streams coincide with large-scale variations in the IMF orientation. Although the data have been constrained to IMF angles relative to the solar wind direction >75 o the observed cut-off is an average over the SWICS field-of-view, with smaller IMF angles reducing the observed value for otherwise identical conditions. Therefore, we will identify the parameter dependencies of these variations in the future and then build them into the model of the PUI cut-off to compensate for these systematic effects. However, these steps are beyond the scope of the current study.
As can be seen from the right panel of figure 4, the resulting best-fit ISN flow longitude based on the 7-year STEREO A PLASTIC PUI cut-off observations from DOY 100 2007 through DOY 100 year data set through DOY 100 2013 is due to an important systematic effect on the PUI cut-off analysis that is revealed by this difference. The average of the differences between the PUI cut-off values taken for the two different data sets shows a small, but noticeable, asymmetry about the ISN flow direction, which explains the direction difference by 0.5 o and goes along with a variation in the average solar wind speed because in the current analysis the relevant radial ISN flow velocity component appears normalized to V SW as w r = V r /V SW . In other words, a slow change in V SW over the year(s) used for the cut-off analysis may lead to a systematic effect on the resulting ISN flow direction. We find that the average solar wind speed mostly decreases over most of the 7 years of the STEREO PLASTIC data set. For decreasing solar wind speed, the resulting flow direction is shifted toward larger longitude and into the opposite direction for increasing average solar wind speed.
Note that the STEREO data set, with decreasing average solar wind speeds, is taken during the rise to the present solar maximum. Conversely, the ACE SWICS data set mostly covers the previous solar maximum and the declining phase, likely with the opposite behavior. Interestingly, the longitude obtained from the combined ACE SWICS data set appears noticeably smaller than the combined STEREO PLASTIC result. These findings reemphasize that the values shown in Table 2 contain significant systematic uncertainties and thus are not to be interpreted yet as new values for the ISN flow longitude. We have merely started to investigate how large the statistical and systematic uncertainties inherent in the data sets and this new analysis method are, and which effects are behind them. The influence of variations in V SW on the PUI cut-off analysis can be mitigated by a modification in the normalization of the PUI distributions with emphasis on the radial component V r of the ISN flow velocity at 1 AU rather than using w as in previous PUI investigations. In addition, the variation of the cut-off speed with V SW is known and can be fully included in the model that is used for the analysis. This is beyond the scope of this paper and will be the subject of a follow-on analysis.

Conclusions and outlook
In this paper we have built on the finding that the cut-off speed of the interstellar PUI distribution has a maximum exactly in the upwind direction, with a symmetric decrease with longitude away from the upwind direction [36]. This behavior has led to a method, which uses a Pearson correlation analysis between the longitude variation of the cut-off and its mirror function to determine the ISN flow direction with a very small fit uncertainty [25].
We have investigated the influence of Poissonian fluctuations in the observed PUI count rates on the correlation result and found: -A correlation method, in which a generic cut-off function is compared with the observed PUI spectra, produces results with substantially smaller uncertainties. -A realistic purely statistical uncertainty for a comparable data set as used in [25] is about twice the reported fit uncertainty.
-A single-year statistical uncertainty with STEREO PLASTIC is likely ≈±0.2 o . Furthermore, we have investigated the overall influence of any systematic effects on the resulting ISN longitude by applying the mirror correlation method to ACE SWICS and STEREO PLASTIC data for individual years, obtaining weighted mean values of these individual results, and comparing them with results from the multi-year data sets. We have found: -Individual year data sets include substantial systematic effects, likely due to the distribution of largescale solar wind structures over the year, as one can see by inspection from the longitude variation of the PUI cut-off values.
-ACE observations are affected much more severely than STEREO observations, with standard deviations around the multi-year weighted mean results of 9. -The multi-year results of the two data sets and the weighted means come close to each other, but still fall outside the 1-σ uncertainties.
-Continuous slow variation of the average solar wind speed, as is expected over the course of the solar cycle, also leads to systematic variation of the resulting ISN flow longitude, which can be addressed within the model for the cut-off and a modification of the normalization of the PUI distributions.
These results re-emphasize the promise that, if the systematic effects have been identified and largely corrected for, the ISN flow longitude can be obtained with very small uncertainties of ≈±0.4 o from STEREO for individual years and ≤0.15 o with the entire growing data set. Systematic effects are much more prominent for the ACE SWICS data. Likely candidates for systematic effects are solar wind speed variations for both data sets, and, in addition, IMF direction changes for ACE SWICS. Although, the PUI distributions have been obtained for the highest time resolution as a function of v/V SW , the implicit variation of V r /V SW is only taken into account for the entire data set thus far, and no adjustment for potential adiabatic changes of the cut-off due to rapid solar wind speed changes [41] have been made. The STEREO observations are largely immune to IMF directional changes as long as the initial PUI torus remains in the sensor field-of-view [39], which is one of the selection criteria for the data because there is no effect on the cut-off in the solar wind frame. The cut-off obtained with ACE SWICS data, on the other hand, is affected by such changes, because the initial torus distribution moves noticeably in v/V SW in the spacecraft frame. This makes the ACE observations much more vulnerable to solar wind structures. However, with the variation of the initial torus distribution with IMF direction known from STEREO PLASTIC observations, this effect can be built into a forward model of the ACE SWICS PUI distributions and thus compensated for. For 2007 even a direction validation is possible because co-located observations with ACE and STEREO exist for this year.
Currently, the PUI cut-off analysis is still subject to systematic effects that apparently increase the uncertainties substantially, for the STEREO multi-year result likely of the order of the standard error of the weighted mean over the yearly results of 0.9 o . It is evident from the improvements made with analyzing both multi-year data sets over individual years that most of the systematic effects are distributed over longitude in a stochastic manner from year to year, thus allowing a compensation. By taking into account solar wind speed variations on the highest time resolution and including IMF variations in the modeling of the SWICS PUI distributions the aforementioned effects will be largely compensated in the future. In addition, a full interplanetary parameter analysis, including solar wind, IMF, and wave power variations, will reveal further systematic effects, and a method used in [32], where they used dual-year data sets with all possible permutations in the combination of years, will allow a further reduction of the systematic effects.
Multi-year observations with ACE and STEREO provide for a precision determination of the ISN flow longitude, which is complementary to the precise determination of a 4-dimensional ISN parameter tube with IBEX [18,19,20,21]. In the future, the combination of these observations will allow a direct comparison with astronomical observations of the physical state of the local interstellar cloud  [16] for any potential small-scale spatial variations and/or any temporal variations of the ISN parameters over the past decades. This analysis presents another stepping stone on the journey to explore the outer boundary of the heliosphere and the surrounding interstellar medium, for which Dr. Ed Stone has paved the way and continues to do so with his effective and tireless lead of key heliophysics missions, such as the Advanced Composition Explorer and the Voyager Interstellar Mission.