  Interaction of the solar wind with interstellar matter involves, among other,
charge exchange between interstellar neutral atoms and plasma, which results in
the creation of secondary population of interstellar neutral (ISN) atoms. The
secondary population of interstellar He was detected by Interstellar Boundary
Explorer (IBEX), but interpretation of these measurements was mostly based on
an approximation that the primary interstellar neutral population and the
secondary population were non-interacting homogeneous Maxwell-Boltzmann
functions in the outer heliosheath. We simulate the distribution function in
the outer heliosheath and inside the heliopause using method of characteristics
with statistical weights obtained from solutions of the production and loss
equations for the secondary atoms due to charge-exchange collisions in the
outer heliosheath. We show that the two-Maxwellian approximation for the
distribution function of neutral He is not a good approximation within the
outer heliosheath but a reasonable one inside the termination shock. This is
due to a strong selection effect: the He atoms able to penetrate inside the
termination shock are a small, peculiar subset of the entire secondary He
population. Nevertheless, the two-Maxwellian approximation reproduces the
density distribution of ISN He inside the termination shock well and enables a
realistic reproduction of the orientation of the plane defined by the Sun's
velocity vector through the local interstellar matter and the vector of
unperturbed interstellar magnetic field.
