INTERSTELLAR NEUTRAL HELIUM IN THE HELIOSPHERE FROM IBEX OBSERVATIONS. III. MACH NUMBER OF THE FLOW, VELOCITY VECTOR, AND TEMPERATURE FROM THE FIRST SIX YEARS OF MEASUREMENTS

Analyses of interstellar neutral helium (ISN He) measurements from the Interstellar Boundary Explorer (IBEX; McComas et al. 2009b; Möbius et al. 2009a, 2009b) by Bzowski et al. (2012) and Möbius et al. (2012) showed a very strong coupling of the He inflow velocity vector (speed and direction) and the upstream temperature in a four-dimensional (4D) "tube" of tightly coupled parameters (McComas et al. 2012). These studies also produced the surprising conclusion that the statistically most likely solution for the flow of ISN He observed during the 2009 and 2010 ISN observation seasons seemed to differ by ∼4° in the inflow direction and ∼3.5 km s⁻¹ in speed from the solution obtained from a coordinated analysis of previous measurements (Möbius et al. 2004), mostly based on direct sampling of ISN He by the GAS instrument on board Ulysses (Witte 2004; Witte et al. 2004). For this different inflow vector, the consistent temperature along the IBEX tube was essentially identical to the GAS value of 6300 K. The same analysis by Bzowski et al. (2012) also showed that the uncertainty range for the obtained solution was relatively wide along the parameter correlation tube in the 4D parameter space, and that the direction and inflow velocity determined by Witte et al. (2004) were not ruled out by the IBEX-Lo analysis, but were statistically less likely. For the Ulysses flow vector, however, the temperature would have to be higher than the temperature found by Witte (2004) by ∼2500 K.

These findings by Bzowski et al. (2012) and Möbius et al. (2012) were of great impact and produced some controversy. If the inflow speed was really 3–4 km s⁻¹ slower than previously thought, then that would require a change in the configuration of the boundary region of the heliosphere compared with earlier views: the heliosphere would most likely feature no fast-mode bow shock (McComas et al. 2012; Zank et al. 2013), and so the disturbance of the local ISN gas by the heliosphere would spatially reach much farther upstream. In addition, Frisch et al. (2013) considered the simplest explanation of the difference between the results of the measurements taken at different epochs, i.e., that the direction of the flow was changing with time. Their analysis of the published measurements and error bars carried out over the past 40 years suggested that the presence of a temporal change was statistically more likely than the lack of a change. In any case, Bzowski et al. (2012), Möbius et al. (2012), and McComas et al. (2012) made it clear that the combined flow vector and temperature published for the Ulysses observations were inconsistent with the 4D tube of parameters required by the IBEX observations. This 4D tube remains intact today, as recently shown by McComas et al. (2015a) and Leonard et al. (2015), and as we show in this paper and other papers in this special issue (Möbius et al. 2015b; Schwadron et al. 2015).

Still, some researchers have voiced concerns about the early IBEX flow vector. Lallement & Bertaux (2014) suggested that the IBEX results were wrong because of errors in the data treatment by the IBEX team (specifically, that the effects of an instrument dead time were ignored) and possible inaccuracies in the modeling. Thus, they thought that the original ISN He flow parameters and temperature from Ulysses were the correct ones and, consequently, considering any change in the inflow parameters with time was not needed. The misleading suggestions of a flawed data treatment by the IBEX team were rebutted by Möbius et al. (2015a, 2015b) and Frisch et al. (2015). The mechanism of reduction in data throughput on board IBEX and the methods used to account for it in the data analysis are explained in detail by Swaczyna et al. (2015).

Since the publication of the original analysis there has been an important development in our understanding of the IBEX-Lo measurements of ISN He: Kubiak et al. (2014) discovered the Warm Breeze—a new population of neutral helium in the heliosphere. The Warm Breeze was shown to flow from ecliptic longitude between 230° and 250° and latitude between 8° and 18° at a speed between 7 and 15 km s⁻¹. Its temperature is between 7000 and 21,000 K, with a most likely value of 15,000 K, and the abundance relative to the primary ISN He ranges from 4% to 10%. The parent population assumed in the analysis was given by the homogeneous Maxwell–Boltzmann distribution function at 150 AU from the Sun. The two hypotheses for the origin of the Warm Breeze which have been identified as the most likely by the discoverers are as follow: (1) it is the secondary population of ISN He arising in the outer heliosheath due to the neutralization of interstellar He⁺ ions via charge exchange with the pristine, undisturbed ISN He population flowing through the outer heliosheath, or (2) it is a non-thermalized flow of neutral He in the local interstellar matter that may originate at the hypothetic nearby boundary of the interstellar cloud which the Sun is embedded in and a medium with a much higher temperature, e.g., as in the hypothesis by Grzedzielski et al. (2010) for the origin of the IBEX Ribbon. The discovery of the Warm Breeze is important in the context of finding the parameters of the primary ISN He because, as we show in the paper, it modifies the ISN He flux distribution observed by IBEX in all orbits.

Another important development was the reanalysis of ISN He observations by the GAS instrument on Ulysses. Bzowski et al. (2014) carried out their analysis using the Warsaw Test Particle Model (WTPM), which also was used by Bzowski et al. (2012) for their analysis of IBEX observations from 2009 and 2010. The analysis method was analogous to the method used in the IBEX analysis and the research team included the personnel who carried out the IBEX analysis and researchers who previously analyzed the Ulysses observations (Banaszkiewicz et al. 1996; Witte et al. 1996, 2004; Witte 2004). The data analysis included the entire set of data fom the GAS experiment, including the subset of data taken during the last orbit of Ulysses that ended in 2007, which had not previously been analyzed. The results indicated (to within errors) that the velocity vector of ISN He from before does not feature any statistically significant changes with time, but that the temperature is much higher than had been obtained from the earlier analysis (Witte 2004). These findings were also confirmed by Wood et al. (2015), who used an analysis method for the Ulysses data analogous to the method used by Lee et al. (2012) and Möbius et al. (2012) to analyze the IBEX measurements.

Finally, Leonard et al. (2015) analyzed measurements of ISN He from 2012 to 2014, restricted to times when the spin axis tilt was less than 02 out of the ecliptic plane. These authors found that using data from times with large (07–49) spin axis pointing out of the ecliptic gave quite different results, meaning that the analytic approximations used were not adequate for such pointing intervals. At the same time, McComas et al. (2015a) combined the Leonard et al. (2015) analysis with Warsaw model calculations using data from 2013 and 2014. These authors found that the very tight coupling between the interstellar He inflow vector (longitude, latitude, speed) and the upstream temperature continued to define the same 4D tube in the parameter space as found earlier by McComas et al. (2012). They further indicated that the resolution of the apparent "Ulysses-IBEX enigma" was simply another location along the tube with a flow vector closer to the Ulysses value but a much higher temperature (by at least 1000 K). They suggested using combined IBEX/Ulysses values of v_ISN = ∼26 km s⁻¹, λ_ISN = ∼255°, β_ISN = ∼5°, and T_ISN = ∼7000–9500 K and discussed the implications of such a substantially warmer region of the interstellar medium than previously thought.

This is one of three papers presenting the latest elements of the analysis of ISN He observed by IBEX, carried out using the analysis concept originally devised by Bzowski et al. (2012). The other two papers include a presentation of the WTPM model used to interpret the observations (Sokół et al. 2015b) and of the uncertainty system for the data and parameter fitting method (Swaczyna et al. 2015). The latter paper also presents corrections for the interface throughput reduction, described by Möbius et al. (2015a), and a very accurately determined IBEX spin axis pointing for all analyzed orbits. These three papers are part of a larger set of papers on interstellar neutrals as measured by IBEX, outlined by McComas et al. (2015b) in this Astrophysical Journal Supplement Series Special Issue.

In the present paper, we analyzed data from all of the seasons of ISN He observations from 2009 through 2014. For the first time, we included the influence of the Warm Breeze on the ISN He flux in the analysis. Compared with the original analysis from 2012, we have improved and optimized WTPM, which allows us to analyze a data set more than four times larger than that in Bzowski et al. (2012). We have updated and refined the previously used model of ionization losses of ISN He by Bzowski et al. (2013) needed to correctly account for the losses of ISN He inside the heliosphere during the observations covering half of the solar cycle length. This update was presented by Bochsler et al. (2014) and Sokół & Bzowski (2014). We also reanalyzed the observational aspects of the measurements. In particular, we developed an analytical model of data throughput throttling effects due to ambient electrons and refined the IBEX spin axis pointing determination, previously analyzed by Hłond et al. (2012). In addition, we collected all known uncertainties in the data and the correlations between them, and we construct a homogeneous, state-of-the-art system of accounting for the data uncertainties in the process of fitting the parameters of ISN He gas. These aspects of the analysis are presented by Swaczyna et al. (2015). Since Kubiak et al. (2014) pointed out the potential sensitivity of the simulated ISN He signal to an energy sensitivity threshold, we also carried out an analysis of the energy sensitivity threshold of the IBEX-Lo instrument. We found that it is more than 20 eV and, consequently, does not affect the analysis of the measurements of ISN He during the spring season of observations of ISN He by IBEX (Galli et al. 2015; Sokół et al. 2015a).

Our paper is complementary to data analyses carried out using alternative methods, presented by Leonard et al. (2015), Möbius et al. (2015a), and Schwadron et al. (2015). The analysis and results that we show here are an extension of the analysis presented by McComas et al. (2015a) based on the data collected in 2013 and 2014.

In this study, we analyze data from IBEX energy step 2, histogram (HB) mode, taken during the ISN He observation seasons 2009 through 2014. An individual data point is the count rate from an individual 6° bin in spin-angle, averaged over the "good time" intervals during an orbit, adopted from the ISN List data product. Details of the observation process relevant for our analysis are presented by Swaczyna et al. (2015). The initial data selection and cleaning, including the identification of the good time intervals in the data, are presented by Möbius et al. (2012) and Leonard et al. (2015). The ISN Good Times List is a data product maintained by the IBEX Science Team. It includes the time intervals when there were no issues with the synchronization of the spin pulse with the IBEX spin period, when Earth and the Moon were far from the field of view, when the count rate from magnetospheric sources outside the ISN signal angle range was low (assessing its magnitude inside the ISN range is challenging because the magnetospheric H is not distinguishable from the ISN atoms), and when the electron count rates, which are responsible for the interface throttling effect discussed by Möbius et al. (2015a) and Swaczyna et al. (2015), are not excessive. We also require reliable spin axis pointing information from the IBEX navigation system. A graphical representation of the ISN good time interval distribution during the orbits used in the analysis is presented in Figure 1.

Zoom In Zoom Out Reset image size

Each of the first six seasons of IBEX-Lo observations of ISN He had some special features. 2009 was the first season and began directly after the commissioning of the spacecraft. Therefore, observations of the Warm Breeze, which typically start in the middle of November of the preceding year, could not be taken. In addition, the level of background due to electrons was elevated for the first few orbits (up to orbit 12), which was eliminated by a change of the on board processing program. 2010 was the first full observation season, the season from which Kubiak et al. (2014) took data to analyze the Warm Breeze. Just as in 2009, the observations were carried out using the regular mode of IBEX-Lo, with all 8 IBEX-Lo energy steps being stepped through in sequence. Data from one orbit, namely, pre-peak orbit 62, were lost due to an on board computer reset. Because IBEX was optimized for its primary measurements of Energetic Neutral Atoms for the outer heliosphere (McComas et al. 2009a, 2009b), and not for ISN viewing, the ISN List coverage of the orbits with the highest count-rate region in Earth's orbit was relatively low. This potentially results in a relatively high statistical weight for the orbits with a lower intensity of the ISN He signal, where the time coverage was better, but which likely include large contributions from the Warm Breeze or ISN H.

In 2011 and 2012, the instrument was operated in special modes during orbits 110 through 114, and 150a through 156a, respectively, and sensor-accumulated rates that are used in the throughput correction algorithm are not available at the required cadence. As a result, the ISN He count rates from the ISN seasons 2011 and 2012 cannot be precisely corrected for the interface throughput reduction. However, as shown by Swaczyna et al. (2015), the ISN He parameter determination is not seriously affected by this. Data from orbits 111 and 155 were lost because of on-orbit spacecraft issues. Additionally, the IBEX orbit was changed in the middle of 2011 to a more stable one against perturbations from the Moon (McComas et al. 2011), and with a longer, ∼9-day period. Therefore, the IBEX spin axis repointing scheme also had to be modified, and after the orbit change the spin axis has been adjusted twice per orbit, close to the perigee and close to the apogee, and consequently the data from individual orbits are split into orbital arcs a and b. In addition, the target for the spin axis latitude was also modified. During seasons 2009 through 2011, the spin axis was maintained ∼07 above the ecliptic, and in 2012 it was alternated between ∼ −01 and ∼07 (Hłond et al. 2012; Swaczyna et al. 2015).

In seasons 2013 and 2014, no data were lost and thus the coverage is approximately symmetric in ecliptic longitude around the maximum of the ISN He flux. An important change, however, was a reduction in the post-acceleration voltage of the instrument because of instrumental reasons, which resulted in a reduction of the counting efficiency, and thus a reduction in the total ISN He count rate. In 2013, the spin axis latitude was maintained between ∼ −02 and ∼04, and in 2014 it was alternated between ∼0° and −49 in an effort to further strengthen the ISN observations. The good time coverage during the 2013 season is the best obtained so far.

Hence, measurements from different seasons are not precisely comparable to each other. Each season has its own specifics, which must be appropriately taken into account in the analysis. In particular, the equivalence of orbits between the seasons is limited because of the different distribution of good time intervals during the orbits and the non-identical distribution of the High Altitude Science Operations (HASO) intervals from season to season. Also, the altitude at which HASO starts and ends, as well as the maximum distance of the spacecraft from the Earth differs, especially for the orbits from 2009 to 2011 compared to those from 2012 to 2014.

An important aspect of data selection is maintaining the homogeneity of the sampled population. Möbius et al. (2012) and Bzowski et al. (2012) noticed the presence of the Warm Breeze during the early orbits of the season and realized that the ISN H signal is expected at later times during each season (end of March/beginning of April). As a consequence, they decided to focus their analysis of ISN He on those orbits for which the contributions from these populations are expected to be minimal. After the analysis of the Warm Breeze by Kubiak et al. (2014) and the ISN H by Saul et al. (2012, 2013), and in particular by Schwadron et al. (2013), it became clear that observations taken from ecliptic longitudes larger than ∼160° include a component of ISN H (see also, e.g., Figure 3 in Kubiak et al. 2013). Leonard et al. (2015) analyzed measurements of ISN He from 2012 to 2014 restricted to orbits when the spin axis tilt out of the ecliptic was less than 02 and suggested that an important Warm Breeze contribution was present in the data from all orbits with ecliptic longitude less than ∼115° (see also Möbius et al. 2015b). They recommended analyzing the ISN He signal based on a subset of orbits restricted to the interval of ecliptic longitudes of spin axis 115°–160°, i.e., approximately from the second week of February until mid-March. We follow this recommendation, which is supported by our own analysis which follows.

IBEX spins around its own axis, which is maintained close to the Sun and near the ecliptic plane. The IBEX-Lo boresight scans a great circle of 360° in the sky, and spin-angles from 180° to 360° point to the ram hemisphere, i.e., to the hemisphere that contains the direction of the Earth motion around the Sun. The peak of the primary ISN He component enters the instrument around spin-angle 264°. The signal due to a single Maxwell–Boltzmann population of ISN He is expected to closely follow a Gaussian shape as a function of the spin-angle. Hence, it should be symmetric about the spin-angle corresponding to the peak flux. Any additional population should manifest itself by breaking this symmetry, unless it happened to be symmetric about identical spin-angles for all of the orbits in question, as in the case of the ISN population, which seems unlikely. Therefore, a good indicator for the presence of an additional population in the signal is a disturbed symmetry of the observed signal.

To test the symmetry of the signal, we analyze its central moments. A moment ${m}^{(n)}$ of the nth order (i.e., the nth moment) is defined as

$$\begin{eqnarray}&&{m}^{(n)}=\displaystyle \sum _{i}{\psi }_{i}^{n}{F}_{i}/\displaystyle \sum _{i}{F}_{i},\end{eqnarray} \tag{ 1 }$$

where ψ_i is the ith spin-angle during a given orbit and F_i is the measured count rate for this bin.

The first moment ${m}^{(1)}$ corresponds to the "center of mass" of the signal and for a symmetric distribution in the spin-angle in which it coincides with the peak value. The central moments ${M}^{(n)}$ of the distribution are defined as

$$\begin{eqnarray}&&{M}^{(n)}=\displaystyle \sum _{i}{F}_{i}{\left({\psi }_{i}-{m}^{(1)}\right)}^{n}/\displaystyle \sum _{i}{F}_{i}.\end{eqnarray} \tag{ 2 }$$

Using this definition, the first central moment is always 0. For a symmetric distribution, the square root of the second central moment is proportional to the width of the distribution, and the third central moment corresponds to the skewness of the distribution and is expected to be zero due to the symmetry.

We found that a small positive skewness occurs in the simulations carried out under the assumption of a single homogeneous Maxwell–Boltzmann population in the gas ahead of the heliosphere. Generally, however, for Earth longitudes larger than ∼115°, the skewness should disappear if only the primary He component is contributing to the observations. We took this as a guideline for the data selection: we calculated the central moments and their uncertainties from the data using the full uncertainty system developed by Swaczyna et al. (2015). We compared the skewness of the observed ISN He count rates with the Warm Breeze subtracted and without it subtracted. The expectation was that the data with the Warm Breeze contribution not subtracted could show some skewness, but subtraction of the Warm Breeze should take it out.

In reality, it turned out that some orbits at the beginning and the end of the ISN He observation seasons show a residual skewness even after subtraction of the Warm Breeze. Aside from some exceptions, subtracting the Warm Breeze, calculated precisely for the observation conditions for a given season with the parameters found by Kubiak et al. (2014), gives smaller skewness but does not eliminate it altogether within the error bars at the early and late orbits of the season, i.e., at the orbits for which the contribution from the Warm Breeze is expected to be the largest. This is illustrated in Figure 2, which shows the skewness of the data for all seasons with uncertainties calculated using the full uncertainty system from Swaczyna et al. (2015). Evidently, subtracting the Warm Breeze eliminates the skewness very well for season 2010: only orbit 67 departs from 0 within the uncertainty range after correction, whereas all orbits except 65 show a statistically significant skewness before correction. This is understandable, since Kubiak et al. (2014) determined the Warm Breeze parameters based solely on data from the 2010 season. Subtracting the Warm Breeze also reduces the skewness satisfactorily for seasons 2009, 2011, and 2014. For the 2012 season, the remnant skewness is largest in the early and late orbits during the season. An exception is season 2013, for which subtracting the Warm Breeze alleviates the skewness only a little in the early orbits and in fact increases it in the late orbits.

Zoom In Zoom Out Reset image size

This analysis suggests that the choice of the range of spin axis ecliptic longitudes from 115° to 160° made by Leonard et al. (2015) and Möbius et al. (2015b) is good for studies of pristine ISN He. It spans 45°, and to approximately match this range in the other dimension, i.e., in ecliptic latitude, we use spin-angle bins from 252° to 282°. In this way, we symmetrically straddle the peak of the count rate, which occurs between the bins corresponding to spin-angles 264° and 270°, and we cut out the far-tail data points from all of the orbits, which have the largest remnant contribution from the Warm Breeze (see Figure 9) and ISN H. We applied these criteria consistently to all six observation seasons. We verified that the data selected in all of the orbits adopted for the analysis fit a Gaussian shape (i.e., departures do not exceed 2σ).

The final data selection is presented in Figure 3, which corresponds to Figure 1 in Sokół et al. (2015a). This figure illustrates the contrast of the data points adopted for ISN He analysis relative to the background. The count rates are shown season by season, sorted in decreasing order, and scaled to the maximum value registered during the given season (so that the season maximum corresponds to 1). Data from individual seasons are presented as vertically spaced a decade apart, to facilitate viewing. The data range corresponding to the ISN signal is clearly visible as the sharply decreasing sequence of points on the left-hand side of the figure. The background corresponds to the data portion starting from approximately 0.25 on the horizontal scale. The region between ∼0.15 and ∼0.25 on the horizontal scale corresponds to a mixture of ISN He and the Warm Breeze, in varying proportions. Note that as a result of the adopted data selection, we leave out the data points with count rates less than approximately 2% of the seasonal maximum, and that the percentage of points used in the analysis is approximately equal from season to season.

Zoom In Zoom Out Reset image size

The data used in this study, along with their uncertainties and the ancillary information needed to reproduce the results, are part of the IBEX Data Release 9. Details are provided by Schwadron et al. (2015), Swaczyna et al. (2015), and Sokół et al. (2015b).

As the physical model, we adopted the distribution function of ISN He ahead of the heliosphere in the form of a Maxwell–Boltzmann function in the reference frame comoving with the unperturbed flow of interstellar gas, with the velocity vector and temperature homogeneous in space. We assume that these conditions apply at a finite distance from the Sun, which we adopted at 150 AU. The rationale for this is explained in Sokół et al. (2015a, see also further discussion in McComas et al. 2015b). This assumption is identical to that adopted by Bzowski et al. (2012). The atoms from the population enter the heliosphere and follow the Keplerian hyperbolic trajectories due to gravitational attraction from the Sun and suffer ionization losses due to their interaction with the solar EUV radiation and solar wind. Details of the calculation of the simulated signal are presented by Sokół et al. (2015b). In order to include the Warm Breeze in the model, we adopt the model presented by Kubiak et al. (2014) and calculate the count rates expected for a given data set. This simulated signal is then subtracted from the data. The Warm Breeze is taken into account throughout the entire analysis presented in this paper. The consequence of neglecting it in the analysis are discussed by Swaczyna et al. (2015), who identified it as the main contributor to the very high values of reduced χ² obtained by Bzowski et al. (2012) in their original analysis of the first two ISN He seasons.

The baseline model fitting method is presented by Swaczyna et al. (2015) in Section 4. In brief, it is a multi-tier minimization of the difference between the data and the model, scaled by the magnitude of the data uncertainty. This difference is denoted as χ² and defined in Equation (16) in their paper. The optimum parameter set corresponds to the global minimum value of χ². To find this optimum parameter set, i.e., the ecliptic longitude λ and latitude β of the inflow direction, the inflow speed v, and the ISN He temperature T ahead of the heliosphere, the model is first calculated on a carefully selected grid of 4297 points in the parameter space. The grid is deployed around a set of equally spaced ecliptic longitude values. For all of the ecliptic longitudes sampled, the grid points with the lowest χ² values are identified and the respective minimum χ² values are found for these longitudes by the minimization of a 3D paraboloid in parameter space, fitted to a limited subset of points around the minima for each λ value. In the last step, the global optimum parameter set is found from a χ² minimization for the model searching for the minimum in χ² as a function of λ. With this minimum identified, the final search for all 4 parameters is performed.

We also used a supplementary method, which approximately corresponds to the method originally developed by Lee et al. (2012) and applied by Möbius et al. (2012, 2015b) and Leonard et al. (2015). In this method, the measured signal is assumed to conform to a Gaussian function of the spin-angle, and the search of the ISN He parameters is based on the analysis of the parameters of this Gaussian function, fitted to a time series of the measurements in each orbit. Thus, the symmetry of the signal around the peak is implicitly assumed and the analysis seems to be less prone to a bias from asymmetries due to the hypothetical presence of an unaccounted population in the signal.

We have adopted a similar approach as an alternative to our direct-fitting method. For the selected orbits and selected spin-angle range 240°–294°, we take the signal summed over all of the good time intervals for a given orbit and fit the Gaussian function of the spin-angle at each orbit in the following form:

$$\begin{eqnarray}&&F\left(\psi \right)=b+f\mathrm{exp}\left[-{\left(\psi -{\psi }_{0}\right)}^{2}/{\sigma }^{2}\right],\end{eqnarray} \tag{ 3 }$$

where b corresponds to the background, f to the peak height of the signal, ψ₀ to the peak position in spin-angle, and σ to the peak width. The Gaussian function parameter fitting is done by χ² minimization using the uncertainty system described by Swaczyna et al. (2015). For each orbit, the Gaussian fit is described by four parameters, and thus the model has a number of parameters equal to four times the number of orbits. Thus, for each orbit j, we obtain a set of local ISN He beam parameters {b_j, f_j, ψ_0j, σ_j} in the IBEX reference frame. The covariance matrix form was obtained as described by Swaczyna et al. (2015), thus the correlations between bins and orbits are also included in the Gaussian representation of the signal. Subsequently, we perform global ISN He parameter fitting by χ² minimization, adopting as the data this series of Gaussian beam parameters, and as the model the series of the Gaussian parameters {f_j, ψ_0j, σ_j} with their covariance matrix, obtained from fitting Equation (3) to the simulated signal for a given ISN He parameter set. The uncertainties are obtained in a fashion similar to  the direct fit method.

We fit the ISN He parameters to the data from all six observation seasons over the spin-angle range from 252° to 288°, for all orbits in the Earth ecliptic longitude range 115° to 160°, fitting the model directly to the observed count rate (with the contribution from the Warm Breeze subtracted). We found that the inflow velocity vector of ISN He obtained from the fit is λ_ISNHe = 25575 ± 033, β_ISNHe = 516 ± 007 in the J2000 ecliptic coordinates, v_ISNHe = 25.76 ± 0.27 km s⁻¹, and temperature T_ISNHe = 7440 ± 190 K. The inflow direction in galactic coordinates is l_ISNHe = 377 ± 016, b_ISNHe = 1494 ± 030, with the correlation coefficient between the uncertainties in the two coordinates equal to −0.86. A plot of χ²(λ) for this calculation is presented in Figure 4. The correlation and covariance matrices for the best-fit parameters are given in Equations (4) and (5):

$$\begin{eqnarray}\mathrm{Cor}=\left(\begin{array}{lcccc} & \lambda & \beta & T & v\\ \lambda & 1.000 & -0.326 & -0.913 & -0.892\\ \beta & -0.326 & 1.000 & 0.285 & 0.294\\ T & -0.913 & 0.285 & 1.000 & 0.950\\ v & -0.892 & 0.294 & 0.950 & 1.000\end{array}\right)\end{eqnarray} \tag{ 4 }$$

$$\begin{eqnarray}&&\mathrm{Cov}={{\boldsymbol{\Sigma }}}_{\mathrm{ISN}\mathrm{He}}\;=\\ \quad \left(\begin{array}{lcc} & \lambda (^\circ ) & \beta (^\circ )\\ \lambda & 0.10734 & -0.00787\\ \beta & -0.00787 & 0.00544\\ T & -0.05641 & 0.00396\\ v & -0.07965 & 0.00591\end{array}\right. & \left.\begin{array}{cc}T({10}^{3}\;{\rm{K}}) & v(\mathrm{km}\;{{\rm{s}}}^{-1})\\ -0.05641 & -0.07965\\ 0.00396 & 0.00591\\ 0.03553 & 0.04879\\ 0.04879 & 0.07421\end{array}\right).\end{eqnarray} \tag{ 5 }$$

Zoom In Zoom Out Reset image size

In addition to the ISN He parameters, we calculated the scaling factors to rescale the simulated flux expressed in the physical units to count rates. These parameters, calculated separately for each observation season, were searched simultaneously with the ISN parameter fitting using the least-squares χ² fitting procedure presented by Sokół et al. (2015b) in their Section 2.7. For the seasons 2009 through 2014, we obtained scaling factors equal to 1.676 · 10⁻⁵, 1.666 · 10⁻⁵, 1.511 · 10⁻⁵, 1.450 · 10⁻⁵, 7.709 · 10⁻⁶, and 7.977 · 10⁻⁶, respectively. These quantities, albeit technical in nature, provide an important insight into the stability of the instrument sensitivity and the fidelity of the ionization loss model, as we detail in the discussion section.

The ISN parameter uncertainties reported in the first paragraph of this section are taken directly from the covariance matrix. Based on the covariance matrix, one could also split the uncertainties into uncertainty parts along and across the correlation tube, as described by Swaczyna et al. (2015, in particular, see their Equations (22) and (23)). Finally, the correlation line, superimposed on the cuts through the calculation grid in the parameter space, is shown in Figure 5. Note that the temperatures in the figures and in the matrices from Equations (4) and (5) are given in thousands of Kelvins. The comparison of the data and simulations, and the residuals in two different views, are shown in Figure 6.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

As pointed out by Swaczyna et al. (2015), there is a very high correlation between temperature and speed: the correlation coefficient is equal to ∼0.95. Also, a very strong correlation exists between the temperature and longitude, and speed and longitude. This manifests itself in the existence of a "correlation tube" in the 4D parameter space, as described above and shown in three cuts in Figure 5. The minimum value of χ² = ∼467 is very high compared with the expected value, which is equal to the number of degrees of freedom N_dof = 254. It is larger by more than 9 standard deviations (equal to $\sqrt{2{N}_{\mathrm{dof}}}$) than the expected value. This suggests that a component in the signal still remains unaccounted for in the analysis, or that the physical model adopted here is not fully adequate, or that the uncertainty system is not complete. This latter possibility makes us scale up the uncertainties obtained formally from the covariance matrix.

A formally unacceptably high value of the minimum χ² obtained from the fitting procedure is not an unusual situation in physics and astronomy. One of the reasons for the high value of χ² could be that despite all of our efforts, the uncertainties are still underestimated. In that case, the simplest procedure to bring the χ² value into an acceptable range is to scale the measurement uncertainties by a common factor equal to the square root of the obtained minimum χ² value divided by the number of degrees of freedom (Olive et al. 2014):

$$\begin{eqnarray}&&S={\left({\chi }^{2}/{N}_{\mathrm{dof}}\right)}^{1/2}.\end{eqnarray} \tag{ 6 }$$

This is equivalent to multiplying the data covariance matrix by S². This procedure does not change the position of the χ² minimum (i.e., the best-fit parameter set), but now the minimum χ² value is exactly equal to the number of degrees of freedom, i.e., to the expected value. This leads to the magnification of the obtained covariance matrix of the fitted parameters by the same factor S², and thus the final uncertainties are also increased by the factor S. All of the following uncertainties, including the uncertainties in Tables 1 and 2, are obtained following this procedure. A discussion of this error scaling is provided in the discussion section.

The very high correlation between speed and temperature is obtained directly from the fitting process without any additional assumptions. Thus, the existence of the "tube" in the 4D parameter space, shown in the original ISN He analysis by Lee et al. (2012), Möbius et al. (2012), and mostly by McComas et al. (2012) based on a simplified analytical model, is experimentally confirmed using a full model of the ISN He flow. Its presence suggests that, effectively, the most robust parameter that can be derived from this analysis is the Mach number of the ISN He flow in front of the heliosphere:

$$\begin{eqnarray}&&{M}_{\mathrm{ISN}\mathrm{He}}=\displaystyle \frac{{v}_{{\rm{B}}}}{\sqrt{\frac{5}{3}\frac{{k}_{{\rm{B}}}T}{{m}_{\mathrm{He}}}}},\end{eqnarray} \tag{ 7 }$$

where v_B is the bulk speed of ISN He in front of the heliosphere, m_He is the mass of He atoms, T the ISN He temperature, and k_B is the Boltzmann constant. The ISN He Mach number obtained from the fit to the entire data set is equal to 5.079 ± 0.028, where the uncertainty was calculated using the full covariance matrix obtained in the fitting, listed in Equation (5) and scaled by the factor S obtained from Equation (6).

To further verify the robustness of the result, we carried out a complementary analysis using the Gaussian parameter fitting method described above, which is closely related to the method used by Möbius et al. (2012, 2015b) and Leonard et al. (2015). The results obtained from this method are in close agreement with the global fit discussed above (λ = 2561 ± 0.6, β = 513 ± 0.12, T = 7330 ± 350 K, v = 25.6 ± 0.5 km s⁻¹, and M_ISNHe = 5.074 ± 0.034), and the conclusions about the presence and strength of the correlation between the parameter pairs are supported. Again, the minimum value of χ² is well above the expected value, even more so than for the direct-fitting method. The results from the direct-fitting and Gaussian parameter fitting methods are within their respective 1σ uncertainties. They are collected in Table 1. In addition to the results just discussed, this table includes results from fits to a subset of data that will be explored in Section 5.2.

In addition to fitting the data from all of the observation seasons together, we performed fits to data from the individual seasons. The results are collected in Table 2. The global fit and the fits from individual seasons yield very similar, but not identical, correlation lines, as shown in Figure 7, which illustrates cuts through the χ²(λ) lines in the 4D parameter space. The correlation lines from individual seasons are very close to each other and to the correlation lines obtained for the entire data set, as well as to the correlation line found by McComas et al. (2012). Additionally, Figure 7 presents the scaled 2σ covariance ellipsoid for the inflow longitude versus the Mach number of the ISN He flow, calculated from the inflow speed and temperature obtained from the global fit and from the fits to individual seasons. Clearly, the Mach number obtained from the fits and the fit longitude of the inflow are correlated, as suggested by the correlation matrix and as it emerges from the simplified analytical theory by Lee et al. (2012, see also McComas et al. 2012).

Zoom In Zoom Out Reset image size

The inflow longitudes from individual seasons are in agreement with each other at the 2σ level of their individual uncertainties, listed in Table 2, but outside the 1σ uncertainties of the result obtained for the combined data set. Regardless of the correlation with the fit inflow longitude, the Mach numbers obtained for individual seasons and from the global analysis are in agreement within ∼2.5%, which is much better than the spread obtained in the temperature and the components of the inflow speed vector. Thus, it can be regarded as the most tightly constrained parameter of the ISN He inflow into the heliosphere observed by IBEX.

In the discussion, we consider the following related issues: (1) does the result depend on the data selection, and (2) how do the results compare for analyses of the same data using different methods? We also discuss the question of what the high minimum χ² value suggests and point out some implications of the obtained results.

5.1. Do the Results from Individual Seasons Statistically Agree with the Result from the Global Fit?

The inflow longitudes obtained from the fits to all of the individual seasons do not form any statistical trend with time. To check whether results from individual seasons are statistically consistent with the result from the global fit, we check if differences between the parameters from a single season ${{\boldsymbol{\pi }}}_{k}=\left({\lambda }_{k},{\beta }_{k},{T}_{k},{v}_{k}\right)$ and from the global fit ${{\boldsymbol{\pi }}}_{\mathrm{ISN}\mathrm{He}}=\left({\lambda }_{\mathrm{ISN}\mathrm{He}},{\beta }_{\mathrm{ISN}\mathrm{He}},{T}_{\mathrm{ISN}\mathrm{He}},{v}_{\mathrm{ISN}\mathrm{He}}\right)$ are statistically consistent with zero. The fit results are vectors in the parameter space and we need to define the distance between two points in this space, taking into account the uncertainties of the vector positions and correlations between the vector elements. For this distance, we adopt the following definition, which takes into account the uncertainties and correlations between the parameters obtained from the fits:

$$\begin{eqnarray}{\chi }_{k}^{2} & = & {\left({{\boldsymbol{\pi }}}_{k}-{{\boldsymbol{\pi }}}_{\mathrm{ISN}\mathrm{He}}\right)}^{{\rm{T}}}\cdot {\left({{\boldsymbol{\Sigma }}}_{k}{S}_{k}^{2}+{{\boldsymbol{\Sigma }}}_{\mathrm{ISN}\mathrm{He}}{S}_{\mathrm{ISN}\mathrm{He}}^{2}\right)}^{-1}\\ & & \cdot \left({{\boldsymbol{\pi }}}_{k}-{{\boldsymbol{\pi }}}_{\mathrm{ISN}\mathrm{He}}\right),\end{eqnarray} \tag{ 8 }$$

where k enumerates the individual seasons, ${{\boldsymbol{\Sigma }}}_{k}$ is the covariance matrix obtained from the fit to the kth season, ${{\boldsymbol{\Sigma }}}_{\mathrm{ISN}\mathrm{He}}$ is the covariance matrix for the global fit (Equation (4)), and S_k, S_ISNHe are the scaling factors calculated from Equation (6) for the kth season and the global fit, respectively. The covariance matrix for the difference between the parameters obtained for an individual season and the global fits is calculated by the uncertainty propagation given by the sum of the respective covariance matrices. Note that we scale our obtained covariance matrix by the respective factors S_k, S_ISNHe, to acknowledge the increased uncertainty due to the high value of the minimum χ² for the fits, and that we assume that the covariance matrix for an individual season ${{\boldsymbol{\Sigma }}}_{k}$ is independent of the covariance matrix for the global fit ${{\boldsymbol{\Sigma }}}_{\mathrm{ISN}\mathrm{He}}.$ Since the matrix ${{\boldsymbol{\Sigma }}}_{k}{S}_{\mathrm{ISN}}^{2}$ in all entries is much smaller than ${{\boldsymbol{\Sigma }}}_{k}{S}_{k}^{2},$ this assumption should not significantly affect the result. The matrix ${{\boldsymbol{\Sigma }}}_{\mathrm{ISN}\mathrm{He}}$ is given by Equation (5), whereas the matrices ${{\boldsymbol{\Sigma }}}_{k}$ are part of the IBEX data release (Schwadron et al. 2015; Sokół et al. 2015b; Swaczyna et al. 2015).

The quantity returned by Equation (8), i.e., the "distance" between the global solution ${{\boldsymbol{\pi }}}_{\mathrm{ISN}\mathrm{He}}$ and the solution ${{\boldsymbol{\pi }}}_{k}$ from an individual season k should be distributed as the χ² probability distribution with 4 degrees of freedom because we have 4 parameters and the statistical model tested (the consistency with zero) has no free parameters. The values obtained for seasons 2009–2014 are presented in Figure 8, where we overlay them on a plot of the probability distribution function (PDF) and the cumulative distribution function (CDF) of the χ² probability distribution for 4 degrees of freedom. Note that the χ² distribution for the low number of degrees of freedom does not resemble the Gaussian (normal) distribution: it is asymmetric relative to the peak (most probable) value, and the expected value (mean), equal to the number of degrees of freedom, is not equal to the most probable value. The largest value is obtained for season 2011, for which we obtain a "distance" equal to 8.22. The smallest "distance" from the global fit was obtained for 2013, i.e., for the season with the best temporal coverage. This is because 2013 is the season that has the largest influence on the final result among the individual seasons, even though the minimum χ² obtained from fitting this season deviates by 3.8σ from the expected value (see Table 2).

Zoom In Zoom Out Reset image size

We treat the "distances" between the results for individual seasons and the global fit, obtained from Equation (8), as random draws from a population that is expected to follow the chi-squared distribution with 4 degrees of freedom. We have a series of six such draws (and thus six values of "distances" ${\chi }_{k}^{2},$ k = 2009, ..., 2014, corresponding to the χ² axis in Figure 8). We ask the following question. If ${\chi }_{k}^{2}$ is randomly drawn from a population following the χ² distribution with 4 degrees of freedom, then what is the probability α_k that such a random draw returns a value not larger than this given ${\chi }_{k}^{2}$? We follow the common scheme for finding the significance level of agreement of an experimental result with expectations. The probability α_k is given by the CDF of the χ² distribution with 4 degrees of freedom. For α = 0.68 (i.e., "1σ"), ${\chi }_{k}^{2}$ would need to be ≤4.75. For α = 0.955 ("2σ"), ${\chi }_{k}^{2}$ would need to be ≤9.7. Thus, our largest sampled ${\chi }_{2011}^{2}=8.22$ agrees with the global result to better than 2σ, but worse than 1σ.

Additionally, we calculated the sum of the ${\chi }_{k}^{2}$ values from individual seasons, which is equal to 29.03, and checked the likelihood of the hypothesis that it is randomly drawn from a χ² probability distribution with 24 degrees of freedom (6 seasons with 4 parameters each). In that case, the CDF value for an argument of 29.03 is equal to 0.78. Thus, this result also confirms the consistency of the results from individual seasons. Another simple test involves calculating the arithmetic mean value of the ${\chi }_{k}^{2}$ values from the six individual seasons and comparing it with the expected value of 4. This arithmetic mean is equal to 4.8 ± 1.2, which is in very good agreement with expectations. Thus, we conclude that the overall scatter of the ISN He parameters obtained from fits to the data from individual seasons around the parameter set obtained from the fit to the full data set is in good agreement with the statistical expectations.

5.2. Does the Global Fit Result Depend on Data Selection?

5.3. Possible Reasons and Implications of the Overly High Minimum χ² Value from the Global Fit

The minimum χ² value obtained from the global fit significantly exceeds the expected value. This is likely due to the fact that the adopted physical model is missing some important aspect and/or element. One of the highly uncertain elements in the model is the Warm Breeze. The contribution from the Warm Breeze to the signal is not negligible and its influence on the fitting result is apparently substantial. The ratio of the Warm Breeze flux to the ISN He flux, calculated for the observation conditions for all six ISN He seasons and an ISN He parameter set close to the best-fit solution, is presented in Figure 9. Since the temperature of the Warm Breeze found by Kubiak et al. (2014) is at least two times greater than the temperature of ISN He, it is not surprising that its contribution at the early and late orbits during the season is largest, and for a given orbit, the spin-angles farthest from the ISN He peak are the most heavily affected. It is interesting, however, that even though the inflow longitude of the Warm Breeze is lower than the inflow longitude of ISN He, the early orbits in each season do not appear to be more affected than the late ones.

Zoom In Zoom Out Reset image size

The relative contributions from the Warm Breeze to the simulated signal are the lowest and the least variable at the peak orbits in each season, and are the largest in the data points farthest from the orbital maximum. Beyond that, however, they do not show any systematic trend: there is no strong dependence on the ionization level, which was increasing from 2009 to 2014 (Sokół & Bzowski 2014), there is no individual season when the Breeze contribution is consistently the largest, etc. The fact that the relative contributions vary from orbit to orbit and from one spin-angle bin to another is likely due to details such as the spin axis pointing, IBEX orbital motion, and length and placement of ISN good time intervals during the orbits.

Following Swaczyna et al. (2015), we verified that not subtracting the Warm Breeze from the data significantly worsens the fit quality obtained from the direct-fitting method and slightly shifts the correlation tube in parameter space, but the fitted ISN He parameters change relatively little (the inflow longitude is still less than 256°). However, using the Gaussian parameter fitting method with the Warm Breeze not subtracted from the data results in an inflow longitude of almost 257° (with the other parameters changed accordingly), and in a χ² value of 609, while the expected value is 118.

Given the strength of the Warm Breeze signal, it is clear that the accuracy of its model in ISN He parameter fitting is critical. However, Kubiak et al. (2014) were able to provide only a wide uncertainty range for its parameters: the intensity was uncertain by ∼50%, temperature and speed by ∼40%, and the inflow direction by 12°. These uncertainties were taken into account in our uncertainty system, but if the true values depart from the nominal values suggested by Kubiak et al. (2014), even within the uncertainty range, then we expect a residual contribution left in the data with the Warm Breeze subtracted, which will affect the present analysis. The Warm Breeze had been fit only to the data from the 2010 ISN He season, with a simplified uncertainty system and without the insights now available (e.g., the exact value of the IBEX-Lo energy cutoff for neutral He, obtained by Galli et al. 2015; Sokół et al. 2015a). Also, two important assumptions for that analysis were that (1) the primary ISN He parameters were exactly those found by Bzowski et al. (2012) to be the most likely parameters, and (2) that the parent population for the Warm Breeze is a homogeneous Maxwell–Boltzmann distribution. In fact, neither of these may actually be the case. Given the insight illustrated in Figure 9, the influence of this hypothetical residual on the fitting results may vary from one observation season to another. In addition, insight provided by Fuselier et al. (2014) and Galli et al. (2014) suggests that there may be a small but non-negligible contribution to the Warm Breeze signal from a local magnetospheric foreground signal, which potentially may be eliminated by comparing the measurements from the three lowest IBEX-Lo energy channels. Therefore, we are not surprised to see evidence in the present ISN He analysis that the Warm Breeze deserves to be revisited, with all seasons analyzed simultaneously. This analysis is ongoing, but given the unknown nature of the Warm Breeze on one hand and the complex observational aspects on the other, it is even more challenging than the analysis of the primary ISN He population.

Figure 6 provides insights on the contributions of the Warm Breeze. Most of the residuals are positive, even though in a perfect case the mean residual should be 0. Such behavior is allowed in the uncertainty system from Swaczyna et al. (2015) because of the correlations between the data points. A strong correlating factor in this context is the model of the Warm Breeze that is subtracted from the data. It is striking that the largest excess in the residuals tends to occur during the early orbits in all seasons, where the contribution from the Warm Breeze is largest and, simultaneously, the signal from the Warm Breeze strongly varies from season to season, unlike in the late orbits where the Warm Breeze contribution is also strong, but the season to season variation is much weaker (see Figure 9). Perhaps the overall density of the Warm Breeze relative to the density of ISN He found by Kubiak et al. (2014) was a little too high and the fitting system "compensated" for it by choosing a solution with a slightly reduced density in the region of the data where the contribution from the Warm Breeze is the largest during the first three seasons (orbits 14 and 15 and equivalent from 2010 to 2011) and relatively less variable. During the 2013 and 2014 seasons, the Warm Breeze contribution in these orbits is lower, which apparently tipped the scales in the global fit toward lower longitudes. Even though the results from three individual seasons suggest a longitude in the range of 258°–261° and the other three suggest a longitude of ∼255° (see Table 2), due to the high number of data points in these latter orbits (with two arcs per orbit, and thus twice as many data points), the global fit returned a longitude similar to that from the three seasons of high statistical weight.

Another potential cause for the high minimum χ² value may be a departure of the parent ISN He distribution function toward a kappa function. A kappa distribution with a small κ value (i.e., a highly unequilibrated distribution) seems to be ruled out by the analysis performed by Sokół et al. (2015a). However, without a detailed analysis, it cannot be excluded that the parent ISN He distribution could be better described by a kappa function with a relatively large κ value. Sokół et al. (2015a) identify the region in the sky observed by IBEX where one should look for evidence of such a departure. We plan to verify this hypothesis once the Warm Breeze contribution to the signal is better understood. However, even if the parent distribution is indeed somewhat kappa-like, i.e., it features elevated tails, it is unlikely that the inferred bulk velocity vector would be biased because the symmetry is not changed. It is likely, however, that such a hypothetical departure could increase the minimum χ² value if a fit to the Maxwell–Boltzman function is attempted. Similarly, additional filtration in the outer heliosheath due to an increased charge exchange rate would modify the signal observed by IBEX. If relatively strong and not axially symmetrical, then it could bias the obtained flow vector and increase the minimum value of χ². However, then the bias in the flow parameters should be present in all direct sampling and pickup ion experiments as well as in remote-sensing studies, for example, observations of the heliospheric helium glow in 58.4 nm (e.g., Vallerga et al. 2004).

Yet another potential cause of the high χ² value in the fits may be data points which could potentially be affected by intermittent systematic effects, as suggested by Schwadron et al. (2015). These authors fitted the positions of the observed peak of the ISN He signal using a model similar to ours. An important difference between the data treatment by Schwadron et al. (2015) and this study is that while we fit the data averaged over the good time intervals for entire orbits, Schwadron et al. (2015 and, in fact, Leonard et al. 2015 and Möbius et al. 2015b) retain their time resolution and fit the time evolution of the latitude of the ISN He beam. Based on this, they are able to identify the intervals of time when the peak positions differ from their fitted model during the portions of individual orbits. Based on this comparison, Schwadron et al. (2015) point out a few short intervals when the observed peak positions systematically differ from their best-fit model. This may be indicative of an additional intermittent signal that is strong enough to shift the fitted peak position of the signal during a time interval of the order of hours. In the absence of this systematic effect, our minimum χ² would be very unlikely to exceed the expected value due to purely statistical fluctuations because the Poisson scatter of the ISN He signal is taken into account in our uncertainty system. Detection of the hypothetical intermittent signal is very challenging in the data product we use since we have a signal integrated over entire good time intervals. The data were selected based on the well-justified criteria described earlier in the paper. We are reluctant to reject data based on their poor fit to the model. Instead, we rely on objective data selection criteria and acknowledge that the filtering thus performed may not be perfect or that the adopted model may be too simple.

The uncertainty scaling due to an excessively high minimum χ² value should be considered in a situation such as ours according to the literature (Rosenfeld 1975; Press et al. 2007; Olive et al. 2014). The procedure assumes that the errors and uncertainties responsible for the excessively high minimum χ² value are independent and Gaussian. This need not be the case, but, on the other hand, scaling the uncertainty to larger values using a well-defined scheme seems to be a better solution than neglecting the extra unknown uncertainty altogether. In our case, we have identified the main suspected cause of the overly high minimum χ² value, i.e., the non-perfect model of the Warm Breeze, and we believe that before a more extensive analysis of this component is available, tentative adoption of the uncertainty assessment we propose in this article is the most reasonable solution.

5.4. Comparison with Other Estimates

5.5. Thermal Mach Number of the ISN He Flow as the Most Robust Parameter Obtained from IBEX Observations

We have analyzed measurements from the first six ISN He observation seasons of IBEX using a sophisticated data uncertainty and parameter fitting system developed by Swaczyna et al. (2015) and the newest version of the ISN He simulation WTPM by Sokół et al. (2015b). We identified a surprisingly strong influence of the Warm Breeze, recently discovered by Kubiak et al. (2014), and accordingly performed careful data selection to obtain the cleanest sample of ISN He observations. We used two alternative parameter fitting methods: (1) the direct flux fitting method, originally used by Bzowski et al. (2012) to analyze IBEX measurements and by Bzowski et al. (2014) to analyze Ulysses measurements of ISN He, and currently strongly refined and optimized by Swaczyna et al. (2015); and (2) the Gaussian beam fitting method, which is an adaptation of the approach originally devised by Lee et al. (2012) and Möbius et al. (2012). In addition to fitting the baseline data set, we also analyzed the data subset used by Leonard et al. (2015) and Möbius et al. (2015b) and subsets from the individual observation seasons.

We found that the ISN He inflow velocity vector seems to be consistent with the velocity vector originally obtained by Witte (2004) from analysis of Ulysses measurements, but only for a temperature that is higher by at least 1100 K. We also found that the Mach number of the neutral He flow ahead of the heliosphere is equal to 5.08. These findings agree with the results of the Ulysses analysis by Bzowski et al. (2014) and Wood et al. (2015) with a very high accuracy of 2.5%, and with the recent study of IBEX results by McComas et al. (2015a). The Mach number obtained now is in full agreement with the results obtained by Bzowski et al. (2012), Möbius et al. (2012), and McComas et al. (2012).

The ISN He speed and temperature are very highly correlated, with the correlation coefficient equal to 0.95. Also, a very high correlation exists between the temperature and inflow longitude. This is because the parameters obtained from observations carried out in a plane close to the plane containing the inflow direction of ISN He are located within a narrow tube in the 4D parameter space. Consequently, the uncertainties of the ISN He parameters are much wider along the parameter correlation than across this line, as shown explicitly by McComas et al. (2012) based on a simplified analytical theory by Lee et al. (2012).

The ISN He inflow parameters were obtained assuming that ISN He adheres to a spatially homogeneous Maxwell–Boltzmann distribution at a distance of 150 AU from the Sun. We refer to this parameter as the tracking distance and physically it can be regarded as the distance at which the test atom was placed into its final solar orbit due to a collision with the ambient member particles from the interstellar gas. This concept is introduced to the analysis to acknowledge the fact that the ISN He gas is collisionless on spatial scales limited to ∼150 AU. The rationale for the adoption of a finite tracking distance and the influence of the tracking distance magnitude on the simulated ISN He signal is presented by Sokół et al. (2015b). The parameter search result is mildly sensitive to the adopted tracking distance, e.g., the speed may change by ∼0.3 km s⁻¹, as shown by McComas et al. (2012, 2015b) and Möbius et al. (2015b).

We also found evidence that the data very likely contain an unaccounted for component that manifests itself systematically in the residuals and the high value of minimum χ² found in the fitting. We suspect that at least part of this component may be due to less than optimum knowledge and removal of the Warm Breeze. However, we do not rule out that the adopted physical model of the parent population for the observed atoms being a superposition of the homogeneous ISN He and Warm Breeze Maxwell–Boltzmann functions may also be too simple.

We obtained a relatively large spread in the ISN He parameters derived from analysis of data from individual observation seasons, but statistical analysis showed that such a spread can be expected given the parameter uncertainties and correlations between them. Therefore, we believe that currently the best estimate of the ISN He inflow velocity vector, temperature, and Mach number are those we found from our analysis of all of the observation seasons together, and that a reliable uncertainty estimate is obtained from the covariance matrix of the fit, adjusted for the excessively high value of minimum χ² found from the fitting. Thus, remembering the possibile existence of an effect unaccounted for in the model that could bias the result, the best estimate for the ISN He temperature and inflow velocity vector from this work is λ_ISNHe = 2558 ± 05, β_ISNHe = 516 ± 010, T_ISNHe = 7440 ± 260 K, v_ISNHe = 25.8 ± 0.4 km s⁻¹, and M_ISNHe = 5.08 ± 0.03, with a very strong correlation between the parameters that follow a 4D tube in the parameter space. The estimates for the ISN He Mach number show a very low relative uncertainty, and thus we consider them to be the most robust result of our analysis.

The portion of the research done in SRC PAS was supported by Polish National Science Center grant 2012-06-M-ST9-00455. A.G. and P.W. thank the Swiss National Science foundation for financial support. Work by the US authors was supported by the IBEX mission, which is part of NASA's Explorer Program.