# Cometary ion drift energy and temperature at comet 67P-Churyumov/Gerasimeko

Hayley N. Williamson,1 Annie Johansson,2 Romain Canu-Blot,1,3 Gabriella
Stenberg Wieser,1,3 Hans Nilsson,1,3 Fredrik L. Johansson,4 and Anja
Moeslinger 1,3
1Swedish Institute of Space Physics, Kiruna 981 92, Sweden
2Luleå University of Technology
3Umeå University
4ESA/ESTEC, Noordwijk, Netherlands E-mail<EMAIL_ADDRESS>
(Accepted XXX. Received YYY; in original form ZZZ)

###### Abstract

The Ion Composition Analyzer (ICA) on the Rosetta spacecraft observed both the
solar wind and the cometary ionosphere around comet 67P/Churyumov-Gerasimenko
for nearly two years. However, observations of low energy cometary ions were
affected by a highly negative spacecraft potential, and the ICA ion density
estimates were often much lower than plasma densities found by other
instruments. Since the low energy cometary ions are often the highest density
population in the plasma environment, it is nonetheless desirable to
understand their properties. To do so, we select ICA data with densities
comparable to those of Rosetta’s Langmuir Probe (LAP)/Mutual Impedance Probe
throughout the mission. We then correct the cometary ion energy distribution
of each energy-angle scan for spacecraft potential and fit a drifting Maxwell-
Boltzmann distribution, which gives an estimate of the drift energy and
temperature for 3521 scans. The resulting drift energy is generally between
11–18 eV and the temperature between 0.5–1 eV. The drift energy shows good
agreement with published ion flow speeds from LAP during the same time period
and is much higher than the cometary neutral speed. We see additional higher
energy cometary ions in the spectra closest to perihelion, which can either be
a second Maxwellian or a kappa distribution. The energy and temperature are
negatively correlated with heliocentric distance, but the slope of the change
is small. It cannot be quantitatively determined whether this trend is
primarily due to heliocentric distance or spacecraft distance to the comet,
which increased with decreasing heliocentric distance.

###### keywords:

methods: data analysis – comets: individual: 67P – plasmas – space vehicles:
instruments

††pubyear: 2024††pagerange: Cometary ion drift energy and temperature at comet
67P-Churyumov/Gerasimeko–B

##  1 Introduction

The interaction of the solar wind with a comet has been observed at multiple
comets, primarily flybys of comets 21P/Giacobini-Zimmer (Brandt et al., 1985),
1P/Halley (Gringauz et al., 1986), and 26P/Grigg-Skjellerup (Grensemann &
Schwehm, 1993). However, a better understanding of comet magnetospheres
required long-term observations. Subsequently, the Rosetta mission, designed
to rendezvous with and orbit a comet, was launched in 2004. It arrived at
comet 67P/Churyumov-Gerasimenko on August 8, 2014 and orbited until September
30, 2016, when it was deorbited into the comet nucleus. While comet 67P was
less active than the comets of previous flyby missions, Rosetta’s escort of
the comet from 3.6 AU to perihelion at 1.2 AU and back out to 3.8 AU offered
the first opportunity to observe a changing comet. In particular, the plasma
environment of comet 67P changed drastically with heliocentric distance,
beginning with a solar wind-dominated environment that became increasingly
cometary ion-dominated as higher insolation led to higher outgassing rates
(Nilsson et al., 2015, 2017). Neutral particles outgassed from the nucleus are
subsequently ionized by photoionization and electron impact ionization. These
newborn cometary ions radially expand from the nucleus and are accelerated by
electric fields (Nilsson et al., 2018; Moeslinger et al., 2023a), to
eventually be picked up by and mass load the solar wind (e.g., Coates, 2004;
Szegö et al., 2000).

The plasma environment around the comet is largely dominated by three
different types of electric fields, all of which vary with both distance from
the comet and comet distance from the sun. When ions are born from the neutral
cometary atmosphere, they originally expand radially away from the nucleus.
The electrons have a much higher thermal speed than the ions, creating a
pressure gradient that leads to a charge separation and hence an ambipolar
electric field. This ambipolar electric field accelerates the ions such that
their flow velocity is significantly greater than that of the neutrals (Vigren
& Eriksson, 2017) and is strongest close to the comet nucleus (Deca et al.,
2019). The ambipolar field magnitude is proportional to the electron number
density gradient, $\sim 1/r$ (Vigren & Eriksson, 2017). Thus, ions are
continuously accelerated as they travel radially outwards from the comet.
Higher energy ions are therefore assumed to have been born further away from
the spacecraft than lower energy ions.

However, the ambipolar field is not the only electric field the cometary ions
experience. Because of the small size of the comet ionosphere and weak
interplanetary magnetic field strength, the ion gyroradius (particularly for
the cometary ions) is much larger than the ionosphere size, even close to
perihelion, making them essentially unmagnetized (Williamson et al., 2022).
However, electrons have a much smaller gyroradius and as such experience
$\vec{E}\times\vec{B}$ drift. The different motions of the ions and electrons
thus create a polarisation electric field in an anti-sunward direction. The
polarisation electric field is strongest at intermediate heliocentric
distances: low column densities far from perihelion will cause a small
electric field magnitude, while at high column densities near perihelion, the
cometary ionosphere modifies the solar wind electric field significantly
enough that the ionosphere does not experience the tailward electric field
(Nilsson et al., 2018). Eventually, at the furthest distances from the
nucleus, the ions are picked up by the solar wind convective electric field
and begin to gyrate in the interplanetary magnetic field. Both low energy
newborn ions, in the process of being accelerated by the ambipolar electric
field, and higher energy pickup ions were observed by Rosetta at comet 67P
(Berčič et al., 2018; Masunaga et al., 2019; Nicolaou et al., 2017; Nilsson et
al., 2017; Stenberg Wieser et al., 2017).

While Rosetta orbited relatively close to the comet nucleus, the expansion of
the comet magnetosphere with decreasing heliocentric distance allowed for
measurements of various solar wind-cometary ion interaction regions
(Williamson et al., 2020). Near perihelion, the cometary magnetosphere was
extensive enough to exclude the solar wind entirely, called the solar wind ion
cavity, and at times even exclude the interplanetary magnetic field, forming a
diamagnetic cavity (Goetz et al., 2016a, b). The time period where the
cometary ionosphere reaches a comparable pressure to the solar wind and as
such is able to significantly deflect the solar wind is of particular
interest, as it indicates a fully-formed induced magnetosphere. Observations
just outside of the solar wind ion cavity, where cometary ion momentum flux is
greater than the solar wind momentum flux, are analogous to being just inside
a magnetopause. Solar wind ion distributions during this time are similar to
those found in a magnetosheath, suggesting the presence of a bow shock
(Williamson et al., 2022). This highly dynamic region is likely to have formed
various types of plasma instabilities (e.g. Kelvin-Helmholtz instabilities
seen at Venus and Mars, Poh et al. (2021); Ruhunusiri et al. (2016); Wolff et
al. (1980)); however, thorough knowledge of the plasma pressure tensor is
useful for analysis of such instabilities. Thus far, such analysis has been
complicated by incomplete knowledge of ion temperatures at comet 67P due to
instrumental effects, as described below.

Rosetta routinely experienced highly negative spacecraft potential
($V_{s/c}$), sometimes exceeding -15 V (Odelstad et al., 2015; Odelstad et
al., 2017), a result of exposed positively charged connectors on the large
solar panels (Johansson et al., 2020). The negative spacecraft charging has
the effect of repelling electrons and attracting lower energy positive ions.
This affects plasma measurements, such as those made by the Ion Composition
Analyzer (ICA) (Nilsson et al., 2007), part of the Rosetta Plasma Consortium
(Carr et al., 2007). ICA was designed to study the interaction between the
solar wind and cometary ions by measuring positive ions through a wide range
of energies and masses. As such, it is capable of measuring the low energy
ions that contribute the bulk of the pressure close to perihelion. However,
the spacecraft potential both accelerates and distorts the trajectories of
these ions, changing their energy distribution. This has the largest effect on
ions with energies $\lessapprox|2V_{s/c}|$ (Bergman et al., 2020a; Bergman et
al., 2021a; Bergman et al., 2021b). In addition, ion densities measured by ICA
were often much lower than those measured by the RPC Langmuir Probe (LAP)
(Eriksson et al., 2007) and Mutual Impedance Probe (MIP) (Trotignon et al.,
2007). It is unclear why this is the case, as the ICA densities can be up to
two orders of magnitude lower than those found by LAP/MIP. One likely
possibility is that the distortion of low energy ion trajectories by the
spacecraft potential, combined with the limited ICA field of view, prevents
ICA from measuring the full cometary ion distribution. If too few low energy
ions are measured, a distribution cannot be inferred, and so no temperature
can be calculated. Thus, finding the temperature necessary for the thermal
pressure requires data with a reasonably complete distribution.

Here, we select ICA data that contains a low energy cometary ion distribution
representative of the surrounding plasma by selecting time periods with
densities comparable to that of LAP/MIP and correct the data for the
spacecraft potential. The resulting distribution is then fit to find both the
drift energy and temperature of the cometary ions across a wide span of the
comet’s orbit, allowing us to draw conclusions on how both parameters change
with heliocentric distance. This can then be used in future work to help
calculate the ion thermal pressure, while also giving insight into the
changing distributions of the cometary ions.

##  2 Instrument description

In this work, we use data from three instruments in the Rosetta Plasma
Consortium (RPC): the Ion Composition Analyzer (ICA), the Langmuir Probe
(LAP), and the Mutual Impedance Probe (MIP). The bulk of the data used is from
ICA, with cross-calibrated LAP-MIP data used for ICA data selection and
spacecraft potential.

###  2.1 Ion Composition Analyzer

The Ion Composition Analyzer (ICA) (Nilsson et al., 2007) was an ion mass
spectrometer in the Rosetta Plasma Consortion (RPC) package (Carr et al.,
2007). The instrument design consisted of an electrostatic entrance deflection
system, used to determine ion arrival angle, a hemispherical electrostatic
analyzer for energy analysis, a magnetic mass analyzer, and a microchannel
plate detector. ICA measured ions in the energy range of a few eV/q to 40
keV/q in 96 energy steps, with $\Delta E/E=0.07$ for $E>30\mathrm{eV}$ and
0.30 for $E<30\mathrm{eV}$, and had a mass resolution of $M/\Delta M\sim 2$.
This mass resolution can distinguish between masses 1, 2, 4, 8, 16, and 32
amu/q. In this paper, we treat all ions with masses > 16 amu/q as water group
ions originating from the comet. Each full energy-angle scan took 192 s.

ICA had a nominal angular field of view of $2\pi$ sr, with
$90^{\circ}\,\times\,360^{\circ}$ coverage with a resolution of
$5.625^{\circ}\,\times\,22.5^{\circ}$ in elevation and azimuth respectively.
However, this field of view is not applicable to the lowest energy ions, such
as those discussed here. Elevation bins are measured using electrostatic
deflection. Thus, the elevation range is dependent on the energy bin being
measured, and so for low energies ($\lessapprox$ 100 eV), only the central
elevation bins are scanned over due to difficulties in setting deflector
voltages accurately. Additionally, the negative spacecraft potential highly
distorts ion trajectories at these energies (Bergman et al., 2020a). As such,
all data used in this study was summed over all angular directions to produce
a 1D distribution in energy. For the spectra we use the heavy ion data from
the level 4 mass-separated Phys-Mass data (Nilsson, 2021).

####  2.1.1 Energy table correction

The energy table required updating upon arrival at comet 67P, when it was
discovered that the ICA high voltage system had an offset compared to on-
ground calibration results (Nilsson et al., 2017). After update, the
uncertainty in ion energies is roughly 3 eV. However, there is an additional
high voltage drift related to the instrument temperature, which occurs when
ICA cooled to below $13.5^{\circ}$C. This can be corrected with the formula:

$E_{c}=\begin{cases}E_{0}+(13.5-T_{s})\times 0.7&\text{if
}T_{s}<13.5^{\circ}\text{C}\\\ E_{0}&\text{otherwise}\end{cases}$ (1)

where $E_{c}$ is the corrected energy, $T_{s}$ is the sensor temperature, and
$E_{0}$ is the original energy. This formula has been applied to all data used
in this study.

###  2.2 Langmuir Probe and Mutual Impedance Probe

The Rosetta Plasma Consortium also contained a Langmuir probe instrument, LAP
(Eriksson et al., 2007), and a Mutual Impedance Probe, MIP (Trotignon et al.,
2007). Data from these two instruments provided spacecraft potential and
electron density for this study, which was compared to the ICA ion density.
LAP consisted of two spherical Langmuir probes, called LAP1 and LAP2, located
on two booms mounted to the spacecraft. The boom lengths were 2.2 m and 1.6 m
for LAP1 and LAP2 respectively. The booms were positioned such that at least
one of the probes was always sunlit and capable of measuring a radial outflow
from the comet (Edberg et al., 2015). LAP primarily measured the current to
the probes as a function of bias voltage, which can then be used to calculate
plasma parameters such as spacecraft potential, electron density, and ion flow
speed.

The Mutual Impedance Probe had two pairs of receiving and transmitting
electrodes mounted to the LAP1 boom 1 m apart. MIP then could retrieve the
plasma density by finding the plasma resonance frequency. If densities were
low, MIP could also use the LAP2 probe as a transmitting electrode, thus
increasing the baseline to 4 m. MIP density retrievals, unlike LAP, are
independent of electron temperature (as long as the plasma frequency is within
its measurement range) and not sensitive to the spacecraft potential. Thus,
because the two measurement techniques are so different, they have been cross-
calibrated to ensure high accuracy (Johansson et al., 2021). This cross-
calibrated electron density dataset (Eriksson et al., 2020) has been used in
this study and is shown in Fig. 1. We additionally use the derived LAP ion
speeds for comparison with our fit parameters. The effective ion speed (i.e.,
a combination of the thermal and drift speeds) can be calculated from the ion
current slope of the LAP bias voltage sweep using the equation:

$\frac{dI}{dV}=\frac{q^{2}n_{i}A_{LAP}}{2m_{i}u_{i}}$ (2)

where $dI/dV$ is the slope of the ion current region, $q$ is the ion charge
(assumed to be singly charged), $n_{i}$ the plasma density from MIP, $A_{LAP}$
is the surface area of the Langmuir probe, $m_{i}$ the ion mass (assumed to be
19 in the LAP dataset), and $u_{i}$ is the effective ion speed (Vigren et al.,
2017).

##  3 Methods

###  3.1 Data selection

We use a few criteria to select ICA data that is both inside the cometopause
and include most of the low energy ion distribution. First, we average the
solar wind, cometary ion, and electron densities and pressures for every 12
hours, approximately the rotation period of the comet, to smooth out
variations from asymmetric outgassing. To ensure the data is inside the
cometopause, we select only data where the cometary ion momentum flux is
larger than the solar wind momentum flux (calculated as in Williamson et al.
(2020)). ICA ion densities are routinely lower than LAP/MIP electron
densities, assumed to be due to unobserved low energy ions. Selecting ICA data
where the cometary ions are a significant percentage of the electron densities
should help to mitigate the issue of the missing ions. We found that setting a
threshold of $n_{c}\geq 0.15n_{e}$ where $n_{c}$ is the cometary ion density
and $n_{e}$ is the LAP/MIP electron density allowed for both a significant
number of time periods that meet the condition as well as ensure ICA was
observing the majority of the local plasma distribution.

Applying both of these conditions to the ICA data results in 259 12-hour time
periods suitable for analysis, 148 of which are in the solar wind ion cavity.
We then take a random sample of 32 of these time periods for use in the
fitting procedure, giving 4676 suitable ICA scans to fit as ICA was not on for
the entirety of every 12 hour period. Figure 1 shows the 12-hour averaged
electron, cometary ion, and solar wind densities, as well as heliocentric
distance, for the entirety of the main Rosetta mission, from August 2014 to
September 2016. The solar wind ion cavity is clearly visible as a gap around
perihelion where the solar wind densities are non-existent, and times when
Rosetta was inside the cometopause are roughly a few months on either side of
the solar wind ion cavity. We show the 12-hour time periods selected for
fitting with orange X markers.

Figure 1: Electron (green), cometary ion (blue), and solar wind (purple)
densities, averaged over 12 hours, for the entire Rosetta mission. Electron
densities are from the LAP/MIP cross-calibrated data and cometary and solar
wind ion densities are from ICA. The black dashed line shows heliocentric
distance. The orange crosses mark the 32 randomly chosen time periods for
fitting.

The densities of the 32 time periods are also shown in Figure 2. We note that,
even though our conditions ensure they are all inside the cometopause, the
time periods still cover a wide range of conditions, with the local solar wind
density varying by several orders of magnitude. Eleven of the time periods are
fully inside the solar wind ion cavity, with solar wind densities equal to
zero, while the other 21 have both cometary and solar wind ions present. The
comparison between the cometary ion densities and electron densities, showing
that they are of a similar order of magnitude, can be seen in this figure.

Figure 2: Electron (green), cometary ion (blue), and solar wind (purple)
densities, averaged over 12 hours, for the 32 chosen time periods.

###  3.2 Fitting method

To determine the bulk speed (drift energy) and 1D temperature of the low
energy ions, we fit a drifting Maxwell-Boltzmann distribution to the data. The
same model was used previously in Bergman et al. (2021a) and Moeslinger et al.
(2023b). The fit is done in energy space and the fit function is given by

$f_{E}(E)=n\,\sqrt{\frac{1}{\pi
kT_{i}E_{i}}}\,\exp\left(-\frac{E+E_{i}}{kT_{i}}\right)\,\sinh\left(\frac{2\sqrt{EE_{i}}}{kT_{i}}\right)$
(3)

where $n$ is the ion density, $k$ the Boltzmann constant and hence $kT_{i}$ is
the ion temperature in eV, and $E_{i}=m_{i}v_{i}^{2}/2$ is the ion drift
energy, equivalent to the bulk speed in energy space. The full derivation of
this equation can be found in Appendix B of Bergman et al. (2021a). For some
scans, the energy distribution has two peaks. For these cases both peaks are
fit with a two-peak version of Equation 3 where the total fit is a sum of the
fits for each peak. As will be shown in Section 4, some scan distributions
have an additional higher energy tail, corresponding to pickup ions. These
tails would be best fit with a kappa distribution; however, here we choose to
focus only on the low energy Maxwellian peaks, as these ions comprise the bulk
of the plasma. Indeed, the value of the distribution function at low energies
is several orders of magnitude higher than that of higher energies, meaning
the kappa distribution tail contributes little to the bulk plasma properties.

The Maxwellian model fit to the data from each scan is inspected manually and
discarded if the fit quality is poor. We additionally calculate a modified
index of agreement (Legates & McCabe Jr., 1999) for each fit to confirm the
fit quality. Slightly less than 25% of the 4676 scans could not be fit
successfully, leaving 3521 scans from the 32 time periods shown in Fig. 2 with
good quality fits.

###  3.3 Spacecraft potential correction

As described in section 1, spacecraft potential affects ions with energies
roughly $<|2V_{s/c}|$. As we are interested in the low energy newborn cometary
ions for this work, we must correct for the acceleration due to spacecraft
potential. While Bergman et al. (2020a) and Bergman et al. (2020b) show that
the spacecraft potential distorts the incoming direction of the low energy
ions, only the one dimensional energy distribution is needed to obtain bulk
speed and temperature. Thus, we only correct for the energy added to the
positive ions. This is done simply by adding the average (negative) spacecraft
potential to the corresponding ICA ion energies. This method was confirmed
with modeling in Bergman et al. (2021a). The new ion energies per charge will
then be:

$\frac{E_{new}}{q}=-V_{s/c}+\frac{E_{ICA}}{q}$ (4)

where $q$ is the ion charge (assumed to be one for the cometary ions). An
example where a drifting Maxwell-Boltzmann distribution has been fit to the
data as described in Section 3.2 is shown in Figure 3. For this sample the
mean spacecraft potential was -13.7 V. The original data is shown in purple,
with stars marking the measurements, while the data corrected for spacecraft
potential is shown in blue, with X’s marking the measurements. In the figure,
it is evident that accounting for the spacecraft potential not only changes
the peak energy as expected, but also changes the fitted temperature. The peak
energy is always decreased and the temperature increased after the spacecraft
potential correction, since the negative spacecraft potential both accelerates
and adiabatically cools the ions.

Figure 3: A single ICA scan from February 14, 2016 fit with a drifting
Maxwell-Boltzmann distribution. Energy per charge in eV is shown on the x axis
and distribution function ($\mathrm{s^{3}m^{-6}}$) on the y axis. The
original, uncorrected ICA data and corresponding fit are shown in purple,
while the data corrected for a spacecraft potential of -13.7 V is shown in
blue. Resulting fit parameters are shown in the legend.

##  4 Results

###  4.1 Fit examples

An example of a simple one peak fit is shown in Fig. 3 and panel (c) of Fig.
4. A higher energy tail not fit by the drifting Maxwell-Boltzmann function can
be seen in Fig. 3; however, the distribution function is typically several
orders of magnitude lower than that of the data fit by the function. Many
scans, particularly closer to perihelion, contain more prominent high energy
tails that appear to fit an exponential or kappa distribution such as that in
panel (a) of Fig. 4; we leave the study of these distributions to a future
work. Other scans include a second clearly Maxwellian peak of smaller
magnitude at higher energies. An example of this is shown in panel (b) Fig. 4.
As these form a minority of the scans, we only include the lower energy peak
in our analysis, as the second peak is always lower amplitude than the primary
low energy peak, but fit the second peak when it appears to ensure the correct
fit for the lower energy peak.

As mentioned previously, approximately 25% of the scans were unable to be
successfully fit. This was primarily due to two reasons: either the
distribution lacked a low energy peak, or there was not enough data for
fitting. Particularly close to perihelion, scans could have a broad, highly
asymmetric distribution that could not be fit with the drifting Maxwell-
Boltzmann, as they did not have a distinct peak. In other samples, scans
occasionally lack a low energy peak entirely, indicating ICA likely did not
observe the full cometary ion distribution in that scan. We also set the
requirement that a scan must have at least 5 data points in order to be fit,
and that the resulting parameters must have $E_{i}>5$ eV and $kT_{i}>0.4$ eV
as arbitrary lower limits to approximate the lowest values resolvable by the
instrument based on examination of the data. The temperature lower limit was
the result of the instrument energy resolution and the energy lower limit
ensured there were measurements below the peak for a better fit. Scans
excluded for this reason tended to be the samples farthest from perihelion,
later in the mission. Thus, the ratio of skipped scans to total scans in a
12-hour sample as a function of heliocentric distance has two peaks, with one
peak close to perihelion from the unsuitable distributions and the second far
from perihelion due to either an unresolvable distribution or too few data
points.

Figure 4: Examples of the three most typical energy distributions in the
samples. Energy per charge in eV is shown on the x axis and distribution
function ($\mathrm{s^{3}m^{-6}}$) on the y axis. Data is given in black and
the fits in purple. Panel (a) shows a fit with one low energy peak followed by
a high energy tail, panel (b) shows a distribution with two Maxwell-Boltzmann
peaks, with the additional second peak fit in green and the total fit in blue,
and panel (c) shows a single low energy peak. The panels are ordered by
increasing heliocentric distance.

###  4.2 Parameters from all fits

Figure 5: The peak drift energy (x axis) and temperature (y axis) for all 3521
scans with a successful fit. If a scan distribution included more than one
peak, we have taken the energy and temperature of the lowest energy peak.
Color indicates heliocentric distance, as shown by the colorbar on the right.

In Fig. 5 we show the resulting fit parameters, drift energy (bulk velocity)
and temperature, for all successful fits, with color indicating heliocentric
distance in AU. Broadly, while there is a significant amount of scatter, there
is also a relatively small range of values that contain most of the data,
visible in the histogram in Fig. 6. 50% of the data has drift energies between
11 and 18 eV (approximately 11–14 km/s flow speed assuming
$\mathrm{H_{2}O^{+}}$) and temperatures between 0.5 and 1 eV. Both
$\mathrm{E_{i}}$ and $\mathrm{kT_{i}}$ are log-normally distributed across the
whole dataset (not shown), with their individual distributions shown in the
top and right panels of Fig. 6. The individual histograms for $E_{i}$ and
$kT_{i}$ indicate that the lower limit cutoffs of 5 eV and 0.5 eV respectively
do not drastically bias the data, as both show peaks well above the cutoff
values.

Figure 6: A histogram of all the fitted peak drift energy and temperature
shown in Figure 5. The top panel shows a histogram of the energy in purple,
the middle shows a 2D histogram of both energy and temperature, and the right
panel shows temperature in green.

###  4.3 Variations with heliocentric distance

One area of interest is the variation of the drift energy and temperature with
heliocentric distance, as this controls the size and complexity of the comet
environment. Additionally, finding such a trend would allow for generalization
across the broader ICA dataset. However, determining the effect of
heliocentric distance is complicated by the varying spacecraft distance to the
comet, which exponentially increased with decreasing heliocentric distance.
The cometocentric distance in our dataset thus ranged from a maximum of 397 km
at a heliocentric distance of 1.25 AU to a minimum of 9.5 km at 3 AU. The
increase on spacecraft cometocentric distance near perihelion was an effect of
the increasing comet activity; the spacecraft endeavored to remain in roughly
the same density corridor. In theory, the changing spacecraft distance should
also affect the ion distributions observed by Rosetta, and so it would be
desirable to disentangle the two effects.

Figure 7: The distance from the comet to the Sun in AU vs the spacecraft
distance to the comet in km for all the data used in this study.

Fig. 7 shows the strong correlation between heliocentric distance,
$\mathrm{R_{sun}}$, and cometocentric distance, $\mathrm{r_{com}}$. The data
examined for this study, while covering a wide portion of the Rosetta mission,
simply does not cover enough parameter space where one distance varies much
more than the other in order to fully disentangle the two effects. Therefore
here we only show trends with respect to heliocentric distance, with
discussion on the cometocentric dependence.

To ascertain trends in the data, we analyze the datasets using the Mann-
Kendall test (Mann, 1945; Kendall, 1975; Gilbert, 1987) and Theil-Sen
regression (Theil, 1992; Sen, 1968). The Mann-Kendall test is a non-parametric
test for monotonic trends. To assess the significance of the trend in the ICA
data, we employed bootstrapping to determine the statistical distribution of
the Mann-Kendall test under the null hypothesis of no trend. This involved
repeated random sampling with replacement from the original dataset. The
significance of the trend is indicated by the observed positive or negative
deviation of the Mann-Kendall test result (of the observed data) from the
null-hypothesis distribution. The test is two-tailed, meaning it is capable of
determining the presence of both positive and negative trends. We evaluate the
test at a 95% confidence level. We use a Theil-Sen regression to estimate a
slope for the parameter trends. Theil-Sen regression was chosen because it is
also a non-parametric test, and therefore no assumptions about the
distribution of the data are required. Additionally, it is less sensitive to
outliers than a traditional least-squares regression, which is useful for our
dataset.

From Fig. 5, it appears that there is more scatter, and hence a larger scan-
to-scan variation, in drift energy and temperature for lower heliocentric
distances (darker colors). To test this assumption, we calculate the
interquartile range of the variation for each sample period and validate the
trend using the Mann-Kendall statistical test. The interquartile range shows a
significant decreasing trend with heliocentric distance for both
$\mathrm{E_{i}}$ and $\mathrm{kT_{i}}$. Thus, the plasma does indeed vary more
from scan to scan at lower heliocentric distances, indicating a more dynamic
plasma near perihelion.

We additionally use the Mann-Kendall test to investigate if the drift energy
and temperature decrease with increasing heliocentric distance, followed by
Theil-Sen regression to estimate a slope. While Fig. 8 shows the mean
parameters per 12-hour sample to minimize clutter, both the Mann-Kendall
significance test and Theil-Sen estimate were done using the full dataset.
Both drift energy and temperature have a highly significant decreasing trend
with increasing heliocentric distance. However, the decrease is small, with
the Theil-Sen regression giving trends of $-2.05R_{sun}+19.7$ for drift energy
and $-0.63R_{sun}+2.5$ for temperature with $R_{sun}$ in AU. Thus, while
distance from the Sun does have a statistically significant effect on the low
energy ion distributions, it is not a strong one.

Figure 8: Variations of the fitted parameters with heliocentric distance in
AU. Panel (a) shows the mean drift energy per sample (eV) in dark purple, the
interquartile range (25% - 75%) for a sample in surrounding light blue, and
the Theil-Sen slope fit as a black dashed line. Panel (b) shows the same for
temperature (eV) in green. The Theil-Sen slope estimates are given in the
panel legends with the 95% confidence intervals. Panel (c) shows the Pearson
correlation coefficient between E and T for each sample, with a dotted line
indicating 0. Gray markers indicate a correlation coefficient with a p-value
above 0.05 and are not statistically significant.

There is no overall strong correlation between $\mathrm{E_{i}}$ and
$\mathrm{kT_{i}}$ in the full dataset, although this is not always the case
for individual 12 hr samples. As seen in panel (c) of Fig. 8, many of the
samples have a negative correlation, i.e. $\mathrm{kT_{i}}$ decreases as
$\mathrm{E_{i}}$ increases. However, there are four samples with a
statistically significant positive correlation, and other samples with no
significant correlation at all. It is unclear why both positive and negative
correlations appear, and the type of correlation does not appear to relate to
heliocentric distance. Possible causes of these varying correlations will be
explored in Section 5.

Unsurprisingly, given the negative correlation between $\mathrm{R_{sun}}$ and
$\mathrm{r_{com}}$, the drift energy and temperature do show a positive
correlation with cometocentric distance (Fig. 9 in the appendix shows the
cometocentric equivalent to Fig. 8). While the significance of this
correlation is similar to those for heliocentric distance, it is highly
dependent on the subset of data chosen. For example, the temperature no longer
shows any trend with cometocentric distance if only data with
$\mathrm{R_{sun}}<1.6$ AU (which equates to $\mathrm{r_{com}}>170$ km) are
selected for the Theil-Sen regression. Even so, attempting to separate the
effects of heliocentric vs cometocentric distance by choosing portions of data
does not remove the inherent bias in the dataset, as the two distances are
much more strongly correlated with each other than with either fit parameter.
With the data used in this study it is therefore not possible to determine
which distance has a stronger effect on the energy and temperature, although
the correlation with heliocentric distance remains even if a subset of the
data is chosen. The relative importance of heliocentric vs cometocentric
distance in the trends in energy and temperature are discussed below.

##  5 Discussion

As stated previously, all data used in this study was selected based on the
average cometary ion density being at least 15% that of the LAP electron
density. However, while this is true for the 12 hr averages, it does vary on a
scan-to-scan basis and helps to illuminate one possible reason ICA densities
are generally lower than LAP densities. In some ICA scans, the low energy
Maxwellian disappears and the distribution function overall has much lower
amplitude. As the density is the integral of the distribution function, this
subsequently drops as well. Thus, during times of comparatively low densities,
it would seem ICA is not detecting the low energy cometary ions. Likely this
is at least partially due to a field-of-view effect, as the spacecraft
potential severely distorts the trajectories of the low energy ions. Near
perihelion counterstreaming cold cometary ions were observed in the region of
the diamagnetic cavity (Bergman et al., 2021b), indicating it is indeed
possible for multiple cold ion beams to exist, which may not all be in the ICA
field-of-view. Because of the trajectory distortion from the spacecraft
potential, the field-of-view is increased, although directional information is
no longer accurate. It is therefore possible that scans with a relatively
complete energy distribution include ions that would normally not be detected
by the instrument, allowing us to make a better assessment of the plasma
properties. Comparison of the fitted drift energies with the LAP ion speeds
shows generally good agreement to within 1–2 sigma of each other, indicating
ICA is likely measuring the majority of the plasma. Some examples of the LAP
ion speeds compared to those found here can be seen in Appendix B. However, it
should be noted that despite the consistency of the ICA and LAP-MIP data,
drift energies in this study are generally higher than previously reported LAP
values of 2–8 km/s, approximately 0.4–6 eV for $\mathrm{H_{2}O^{+}}$ (e.g.,
Vigren et al., 2017). It is possible that ICA can only detect the full ion
distribution when flow speeds are relatively high, which could explain the
difference between the bulk speeds found here and in previous studies. Despite
this caveat, the LAP-MIP comparison indicates that the ion distribution data
used in this work are reasonably complete and useful for scientific
interpretation.

As shown in section 4.3, while we show trends with heliocentric distance, it
is likely that cometocentric distance is also a contributor to the changes in
drift energy and temperature. Certainly, a larger cometocentric distance can
contribute to the increase in E and T close to perihelion for several reasons.
For example, at farther cometocentric distances, the same instrument field of
view covers a wider spatial area with significant ion production. Thus, ICA
will observe ions produced at different distances that have experienced
differing amounts of acceleration. Such observations of ions originating from
different distances can increase the measured temperature. In general, the
longer ion travel time from the ionization point to the spacecraft allows for
more processes to affect the plasma, many of which will increase either the
drift energy or temperature. However, some of these processes are also
correlated with cometary activity, as the comet ionosphere is larger and more
fully developed close to perihelion. Here, we focus on electric fields and
wave activity as drivers of the increase in drift energy and temperature.

As explained in the introduction, the comet environment is largely governed by
three different electric fields, in order of increasing distance from the
comet: the ambipolar field, the polarisation electric field, and finally the
solar wind convective electric field. The ambipolar field, with its magnitude
of $\sim 1/r_{com}$, may be one reason the mean drift energy is higher closer
to perihelion when both the comet atmosphere size and spacecraft distance from
the nucleus increase. However, effects from all three of these fields can be
observed in the ion energy spectra. Although flow direction is not shown here,
the additional higher energy ion populations seen in panels (a) and (b) of
Fig. 4 have a strong antisunward velocity component, indicating acceleration
by the polarisation or solar wind electric fields. The second Maxwellian peak
is likely ions born farther from the spacecraft that have been accelerated by
the polarisation electric field or begun the pickup process. Likewise, the
high energy ions have been seen as a signature of pickup ion detection (Behar
et al., 2016; Berčič et al., 2018; Broiles et al., 2015; Coates et al., 2015;
Moeslinger et al., 2023a; Nicolaou et al., 2017; Williamson et al., 2022),
with the highest energy ions likely coming from far upstream of the
spacecraft. This mixed population is thus no longer in thermal equilibrium and
so the distribution evolves from the drifting Maxwell-Boltzmann to a more
kappa-like distribution. Thus, from the bottom to the top of Fig. 4 (i.e.,
decreasing heliocentric distance), we see the evolution of the ion spectra
away from equilibrium, although the cold Maxwellian population remains.

Waves also cannot be ignored as a cause for the inverse correlation with
heliocentric distance, particularly for temperature. Indeed, the presence of a
Maxwell-Boltzmann distribution itself is indication of wave activity, as the
ions are born with the temperature of the neutral gas, with a distribution
function that then evolves towards a Maxwellian as it expands (see e.g. Fig. 5
in Vigren & Eriksson (2017)). Many types of plasma waves were observed at
comet 67P by Rosetta, including lower hybrid waves at low cometocentric
distances (André et al., 2017; Karlsson et al., 2017; Stenberg Wieser et al.,
2017), singing comet waves generated by pickup ions (Goetz et al., 2020), ion
acoustic waves (Gunell et al., 2017a, b), possible ion Bernstein waves
(Odelstad et al., 2020), and mirror mode waves (Tello Fallau et al., 2023;
Volwerk et al., 2016). With a larger comet ionosphere and longer travel
distance to the spacecraft closer to perihelion, there is a higher likelihood
that the ions measured by ICA have experienced wave heating. Thus, wave
activity is a likely explanation for higher temperatures at lower heliocentric
distances. Because wave activity is highly transient, it is also a logical
reason for the large scan-to-scan variability in both energy and temperature.
This is the case throughout all the data, as seen in Fig. 8, but is more
prominent closer to perihelion. In general, the comet ionosphere has the
highest energy input close to perihelion, and subsequently the largest and
most complex structure. This then is likely the reason for the larger
interquartile range of both E and T close to perihelion.

It can be seen in panel (c) of Fig. 8 that, while the majority of the samples
have a negative correlation between E and T, this is not consistent, and shows
no trend with heliocentric distance. At times, indicated by the green dots in
the figure, there is no statistically significant correlation at all. As with
the high variability, this is likely due to the transient and dynamic nature
of the comet ionosphere. When examining E vs T for individual samples, some
cases show a clear shift in time, moving from e.g. a cluster of scans with
relatively low E and high T to later scans with higher E and lower T. Such a
change is a clear indication of an event occurring, although whether this is
due to a change in the upstream solar wind which then affects the ionosphere
or due to a more localized event is not possible to determine from the data.
One such shift in time occurs for the most positively correlated sample seen
in Fig. 8(c), with E and T both being originally low, then suddenly
increasing, where they remain high for the rest of the time period. The
cluster of data points from the beginning of the sample is only weakly
positively correlated, but both a stronger correlation in the second cluster
of data points and the shift to higher E and T create a strong positive
correlation for the whole sample. Thus, this shifting in time can affect the
correlation for a single sample, meaning that the changing plasma conditions
also affect the correlation between E and T.

##  6 Conclusions

We fit more than 3500 energy scans from the ICA instrument with a drifting
Maxwell-Boltzmann function. By correcting for spacecraft potential and only
choosing time periods with densities comparable to other instruments, we are
able to focus on the low energy cometary ions that would otherwise be
difficult to analyze. By doing so, we see that the drift energy and
temperature are generally within a relatively narrow range of 11–18 eV and
0.5–1 eV, respectively, for much of the Rosetta mission. This is despite large
changes in the comet ionosphere during the time periods analyzed, such as the
formation of the solar wind ion cavity. The drift energy (equivalent to flow
speed) is higher than the mission-wide medians reported from the Rosetta LAP
data (Johansson et al., 2021); this may be due to selection bias, as the flow
may need to be higher energy in order to be observed in full by ICA. However,
for the scans analyzed, the fitted drift energy is consistent with the
coinciding LAP data.

While the drift energy and temperature do not change drastically throughout
our samples, we do see two statistically significant changes: both energy and
temperature decrease with increasing heliocentric distance, as does the
interquartile range or scan-to-scan variability. This can largely be
attributed to the comet ionosphere being both larger in size and more complex
at the lowest heliocentric distance, although whether the correlations are
strictly due to heliocentric distance or spacecraft distance from the comet,
which increased towards perihelion, cannot be determined. We attribute the
higher energy and temperature to both the accelerating electric fields and
wave-particle interactions, which have more opportunities to accelerate or
heat the ions at low heliocentric/large spacecraft distances. As the comet
ionosphere increases in complexity, we can observe an evolution of the
cometary ion distribution from a single cold Maxwellian population to include
higher energy accelerated pickup ions.

## Acknowledgements

The work of H. Williamson and A. Johansson was supported by Swedish National
Space Agency grant 2021-000105. A. Moeslinger was supported by SNSA grant
132/19. The authors thank the members of the Rosetta Plasma Consortium team
for their helpful discussions.

## Data Availability

The Rosetta RPC-ICA, RPC-LAP, and RPC-MIP data are publicly available through
the ESA Planetary Science Archives at https://psa.esa.int/. Fit parameters for
all samples used in the paper are available at
https://data.irf.se/data/williamson2024mnras/. Mann-Kendall tests were
performed with the pyMannKendall software package (Hussain & Mahmud, 2019).
Scipy v1.12.0 was used for Theil-Sen regression (Virtanen et al., 2020).
Colormaps were used courtesy of CMasher (van der Velden, 2020).

## References

  * André et al. (2017) André M., et al., 2017, Monthly Notices of the Royal Astronomical Society, 469, S29
  * Behar et al. (2016) Behar E., Nilsson H., Wieser G. S., Nemeth Z., Broiles T. W., Richter I., 2016, Geophysical Research Letters, 43, 1411
  * Bergman et al. (2020a) Bergman S., Stenberg Wieser G., Wieser M., Johansson F. L., Eriksson A., 2020a, Journal of Geophysical Research: Space Physics
  * Bergman et al. (2020b) Bergman S., Stenberg Wieser G., Wieser M., Johansson F. L., Eriksson A., 2020b, Journal of Geophysical Research: Space Physics, 125
  * Bergman et al. (2021a) Bergman S., et al., 2021a, Monthly Notices of the Royal Astronomical Society
  * Bergman et al. (2021b) Bergman S., Wieser G. S., Wieser M., Nilsson H., Vigren E., Beth A., Masunaga K., Eriksson A., 2021b, Monthly Notices of the Royal Astronomical Society, 14, 1
  * Berčič et al. (2018) Berčič L., Behar E., Nilsson H., Nicolaou G., Stenberg Wieser G., Wieser M., Goetz C., 2018, Astronomy and Astrophysics, 613, 1
  * Brandt et al. (1985) Brandt J., Farquhar R., Maran S., Niedner M., von Rosenvinge T., 1985, in Carusi A., Valsecchi G., eds, , Vol. 115, IAU Colloq. 83: Dynamics of Comets: Their Origin and Evolution. Springer, Dordrecht, p. 405, doi:10.1007/978-94-009-5400-7_35
  * Broiles et al. (2015) Broiles T. W., et al., 2015, Astronomy and Astrophysics, 583
  * Carr et al. (2007) Carr C., et al., 2007, Space Science Reviews, 128, 629
  * Coates (2004) Coates A. J., 2004, Advances in Space Research, 33, 1977
  * Coates et al. (2015) Coates A. J., Burch J. L., Goldstein R., Nilsson H., Wieser G. S., Behar E., 2015, Journal of Physics: Conference Series, 642
  * Deca et al. (2019) Deca J., Henri P., Divin A., Eriksson A., Galand M., Beth A., Ostaszewski K., Horányi M., 2019, Physical Review Letters, 123, 55101
  * Edberg et al. (2015) Edberg N. J., et al., 2015, Geophysical Research Letters, 42, 4263
  * Eriksson et al. (2007) Eriksson A. I., et al., 2007, Space Science Reviews, 128, 729
  * Eriksson et al. (2020) Eriksson A., Gill R., Johansson E. P., Johansson F. L., 2020, Rosetta RPC-LAP archive of derived plasma parameters from the COMET ESCORT 4 mission phase, RO-C-RPCLAP-5-ESC4-DERIV2-V1.0, doi:https://doi.org/10.57780/esa-vcc50n4, https://archives.esac.esa.int/psa/ftp//INTERNATIONAL-ROSETTA-MISSION/RPCLAP/RO-C-RPCLAP-5-ESC4-DERIV2-V1.0
  * Gilbert (1987) Gilbert R. O., 1987, Technical Report PNL-4660, Statistical Methods for Environmental Pollution Monitoring, https://www.osti.gov/biblio/7037501. Pacific Northwest National Lab. (PNNL), Richland, WA (United States), https://www.osti.gov/biblio/7037501
  * Goetz et al. (2016a) Goetz C., et al., 2016a, Monthly Notices of the Royal Astronomical Society, 462, S459
  * Goetz et al. (2016b) Goetz C., et al., 2016b, Astronomy and Astrophysics, 588, 1
  * Goetz et al. (2020) Goetz C., Plaschke F., Taylor M. G. G. T., 2020, Geophysical Research Letters, 47, e2020GL087418
  * Grensemann & Schwehm (1993) Grensemann M. G., Schwehm G., 1993, Journal of Geophysical Research, 98, 20907
  * Gringauz et al. (1986) Gringauz K. I., et al., 1986, Nature, 321
  * Gunell et al. (2017a) Gunell H., et al., 2017a, Monthly Notices of the Royal Astronomical Society, 469, S84
  * Gunell et al. (2017b) Gunell H., et al., 2017b, Astronomy and Astrophysics, 600, 1
  * Hussain & Mahmud (2019) Hussain M. M., Mahmud I., 2019, Journal of Open Source Software, 4, 1556
  * Johansson et al. (2020) Johansson F. L., Eriksson A. I., Gilet N., Henri P., Wattieaux G., Taylor M. G. G. T., Imhof C., Cipriani F., 2020, Astronomy & Astrophysics, 642, A43
  * Johansson et al. (2021) Johansson F. L., et al., 2021, Astronomy & Astrophysics, 653, A128
  * Karlsson et al. (2017) Karlsson T., et al., 2017, Geophysical Research Letters, 44, 1641
  * Kendall (1975) Kendall M. G., 1975, Rank correlation methods, 4th ed., 2d impression edn. Griffin London, London
  * Legates & McCabe Jr. (1999) Legates D. R., McCabe Jr. G. J., 1999, Water Resources Research, 35, 233
  * Mann (1945) Mann H. B., 1945, Econometrica, 13, 245
  * Masunaga et al. (2019) Masunaga K., Nilsson H., Behar E., Stenberg Wieser G., Wieser M., Goetz C., 2019, Astronomy & Astrophysics, 630, A43
  * Moeslinger et al. (2023a) Moeslinger A., Nilsson H., Stenberg Wieser G., Gunell H., Goetz C., 2023a, Journal of Geophysical Research: Space Physics, 128
  * Moeslinger et al. (2023b) Moeslinger A., Wieser G. S., Nilsson H., Gunell H., Williamson H. N., LLera K., Odelstad E., Richter I., 2023b, Journal of Geophysical Research: Space Physics, 128, e2022JA031082
  * Nicolaou et al. (2017) Nicolaou G., Behar E., Nilsson H., Wieser M., Yamauchi M., Berčič L., Wieser G. S., 2017, Monthly Notices of the Royal Astronomical Society, 469, S339
  * Nilsson (2021) Nilsson H., 2021, ROSETTA-ORBITER 67P RPCICA 4 EXT3 RESAMPLED AND CALIBRATED V1.0, RO-C-RPCICA-4-EXT3-PHYS-MASS-V1.0, https://pdssbn.astro.umd.edu/holdings/ro-c-rpcica-4-ext3-phys-mass-v1.0/dataset.shtml
  * Nilsson et al. (2007) Nilsson H., et al., 2007, Space Science Reviews, 128, 671
  * Nilsson et al. (2015) Nilsson H., et al., 2015, Astronomy & Astrophysics, 583, 1
  * Nilsson et al. (2017) Nilsson H., et al., 2017, Monthly Notices of the Royal Astronomical Society, 469, S252
  * Nilsson et al. (2018) Nilsson H., et al., 2018, Astronomy and Astrophysics, 616, 1
  * Odelstad et al. (2015) Odelstad E., et al., 2015, Geophysical Research Letters, 42, 10126
  * Odelstad et al. (2017) Odelstad E., Stenberg-Wieser G., Wieser M., Eriksson A. I., Nilsson H., Johansson F. L., 2017, Monthly Notices of the Royal Astronomical Society, 469, S568
  * Odelstad et al. (2020) Odelstad E., et al., 2020, Journal of Geophysical Research: Space Physics, 125, 1
  * Poh et al. (2021) Poh G., et al., 2021, Journal of Geophysical Research: Space Physics, 126
  * Ruhunusiri et al. (2016) Ruhunusiri S., et al., 2016, Geophysical Research Letters, 43, 4763
  * Sen (1968) Sen P. K., 1968, Journal of the American Statistical Association, 63, 1379
  * Stenberg Wieser et al. (2017) Stenberg Wieser G., et al., 2017, Monthly Notices of the Royal Astronomical Society, 469, S522
  * Szegö et al. (2000) Szegö K., et al., 2000, Space Science Reviews, 94, 429
  * Tello Fallau et al. (2023) Tello Fallau A., Goetz C., Simon Wedlund C., Volwerk M., Moeslinger A., 2023, Annales Geophysicae, 41, 569
  * Theil (1992) Theil H., 1992, in Raj B., Koerts J., eds, , Henri Theil’s Contributions to Economics and Econometrics: Econometric Theory and Methodology. Springer Netherlands, Dordrecht, pp 345–381, doi:10.1007/978-94-011-2546-8_20, https://doi.org/10.1007/978-94-011-2546-8_20
  * Trotignon et al. (2007) Trotignon J. G., et al., 2007, Space Science Reviews, 128, 713
  * Vigren & Eriksson (2017) Vigren E., Eriksson A. I., 2017, The Astronomical Journal, 153, 150
  * Vigren et al. (2017) Vigren E., et al., 2017, Monthly Notices of the Royal Astronomical Society, 469, S142
  * Virtanen et al. (2020) Virtanen P., et al., 2020, Nature Methods, 17, 261
  * Volwerk et al. (2016) Volwerk M., et al., 2016, Annales Geophysicae, 34, 1
  * Williamson et al. (2020) Williamson H. N., Nilsson H., Stenberg Wieser G., Eriksson A. I., Richter I., Goetz C., 2020, Geophysical Research Letters, 47
  * Williamson et al. (2022) Williamson H. N., Nilsson H., Stenberg Wieser G., Moeslinger A., Goetz C., 2022, Astronomy & Astrophysics, 660, A103
  * Wolff et al. (1980) Wolff R. S., Goldstein B. E., Yeates C. M., 1980, Journal of Geophysical Research: Space Physics, 85, 7697
  * van der Velden (2020) van der Velden E., 2020, The Journal of Open Source Software, 5, 2004

##  Appendix A Variations with cometocentric distance

Figure 9 shows the fit parameters as a function of cometocentric distance in
the same format as Fig. 8. As with the previous figure, only the mean values
for each sample are plotted, but the Theil-Sen slope estimate and the Mann-
Kendall trend and significance were calculated using the full dataset. Both
energy and temperature have a statistically significant increasing trend with
cometocentric distance.

Figure 9: Variations of the fitted parameters with cometocentric distance in
km. Panel (a) shows the mean drift energy per sample (eV) in dark purple, the
interquartile range (25% - 75%) for a sample in surrounding light blue, and
the Theil-Sen slope fit as a black dashed line. Panel (b) shows the same for
temperature (eV) in green. The Theil-Sen slope estimates are given in the
panel legends with the 95% confidence intervals. Panel (c) shows the Pearson
correlation coefficient between E and T for each sample, with a dotted line
indicating 0. Gray markers indicate $p>0.05$ and are not statistically
significant.

##  Appendix B LAP speed vs ICA fitted speed

We show some examples of the LAP effective ion flow speeds, derived from
$dI/dV$ for each bias voltage sweep, compared to the $E_{i}$ from the fitted
ICA data. The ICA fitted drift energy has been converted to bulk speed
assuming a particle mass of 19 amu, consistent with the LAP dataset. The LAP
effective ion flow speed data was taken from the analyzed sweep (ASW) data in
Eriksson et al. (2020), filtered for quality values (a measure that depends on
the goodness of the fit on the ion current slope and the MIP spectrum quality)
above 0.7. The figure shows that the LAP and ICA derived speeds generally
overlap, although they are not typically the exact same. In many individual
samples, such as that in Fig. 10(b), they also show similar trends in time.

Figure 10: The LAP dI/dV effective ion flow speed in purple compared to the
results from this study in green. Panel (a) shows data for two months that
contain multiple fitted sample periods, while panel (b) shows a zoom in on one
of the samples, November 24, 2015. The ICA fitted drift energy has been
converted to bulk speed in km/s.


# On-Site Production of Quasi-Continuous Ultra-High Vacuum Pipes

Matthias Angerhausen  Guido Buchholz  Jef Hoste  Marion Purrio  Achim Stahl
Lars Stein  Patrick Toussaint

###### Abstract

We present a design study for a new production technology for ultra-high
vacuum pipes. The pipes are produced in a fully automatised process in
sections of hundreds of meters directly in the later location of usage. We
estimate the effort for such a production and show that it might be
substantially lower than the effort for an off-site production of
transportable sections.

††journal: Vacuum

[FEF] organization=Forschungs- und Entwicklungsgesellschaft Fügetechnik GmbH,
addressline=
Driescher Gässchen 5, city=Aachen, postcode=52062, country=Germany

[ISF] organization=Welding and Joining Institute, RWTH Aachen University,
addressline=
Pontstrasse 49, city=Aachen, postcode=52062, country=Germany

[Werkhuizen Hengelhoef] organization=Werkhuizen Hengelhoef,
addressline=Hengelhoefstraat 162, city=Genk, postcode=3600, country=Belgium

[aperam] organization=Aperam Stainless Europe, addressline=Swinnenwijerweg 5,
city=Genk, postcode=3600, country=Belgium

[RWTH] organization=RWTH Aachen University, addressline=Templergraben 55,
city=Aachen, postcode=52062, country=Germany

##  1 Motivation

The Einstein Telescope is a gravitational wave detector of the third
generation in preparation in Europe [1]-[4]. It is based on six Michelson
interferometers enhanced by Fabry-Perot arm cavities, with an arm length of
$10~{}\mathrm{km}$ each. The interferometers will be located in underground
tunnels in the shape of a $10~{}\mathrm{km}$ equilateral triangle. The laser
beams will travel in ultra-high-vacuum (UHV) pipes. A diameter of roughly
$1~{}\mathrm{m}$ is required to enclose the beams. Plain pipe sections of
$500~{}\mathrm{m}$ length are envisioned, connected to each other by pumping
stations. In total $120~{}\mathrm{km}$ of these UHV pipes will be required
plus an additional $10$ to $20~{}\mathrm{km}$ of smaller pipes. It will be the
largest UHV system ever built.

The Einstein Telescope will measure gravitational waves through the changes in
the length of the arms caused by the passing waves. The telescope will be able
to detect changes down to $10^{-20}~{}\mathrm{m}$. Fluctuations in the
pressure inside the UHV pipes along the path of the lasers induce fluctuations
in the index of refraction. This might lead to variations in the phase of the
laser waves, mimicking a passing wave. A vacuum level of
$10^{-10}~{}\mathrm{mbar}$ at room temperature is required to limit the impact
of this effect. Furthermore, a vacuum at this level reduces the probability of
laser photons scattered from the mirror surfaces to be rescattered from the
rest gas back into the laser beam. These photons would be added incoherently
to the beam and would create noise, too.

The cost of the Einstein Telescope is driven by the cost of the
infrastructure, mainly the tunnels and the vacuum system. The observatory is
seen as an infrastructure hosting the Einstein Telescope as its initial
instrument. It might be replaced after a decade of operation by a new
instrument with even higher sensitivity. This instrument will use the same
tunnels and vacuum system. Its specifications on the vacuum are still unknown.
To serve the next generation, the vacuum system should be able to reach
$10^{-11}~{}\mathrm{mbar}$, if this is achievable with a reasonable effort.

Conventionally, beam pipes for gravitational wave detectors and similarly for
particle accelerators are built from pipe sections prefabricated off-site and
welded together in the tunnels. The length of the prefabricated sections are
limited to about $20~{}\mathrm{m}$ by transportation issues. We investigated
the possibility to transport coils of sheet metal into the tunnels and to
fabricate $500~{}\mathrm{m}$-long pipe sections in a continuous process
directly in the tunnels. In this paper we summarise our findings. We present a
highly reliable, fully automated production procedure. With a single robot the
full length of pipes of $120~{}\mathrm{km}$ can be produced in about two years
with potentially large savings in cost.

##  2 Pipe Production Technology

###  2.1 Requirements

The exact diameter of the pipes is not known, yet, but it will be
approximately $1~{}\mathrm{m}$. It is determined by the diameter of the laser
beams. We plan for $1~{}\mathrm{m}$ pipes. The interferometers will include
arm-cavities of $10~{}\mathrm{km}$ length plus around $100~{}\mathrm{m}$ of
distance between the beam splitter and the cavities. The mirrors of the three
low-frequency interferometers will be cooled to cryogenic temperature
($\approx 10~{}\mathrm{K}$), while the arms will be kept at room-temperature.
The sections around the cold mirror will need special attention. Here we
consider only the room-temperature pipes. Approximately $120~{}\mathrm{km}$
are required. In addition tens of kilometres of UHV pipes of smaller diameter
will be needed for auxiliary optics, which we ignore here.

To reach a vacuum of $10^{-10}~{}\mathrm{mbar}$ or better, hydrocarbons have
to be removed from the inner surfaces by thorough cleaning. The pipes have to
be baked under vacuum at a temperature of at least $120^{\circ}~{}\mathrm{C}$
to detach water from the walls. For the quality of the surface, the pipes are
built from a stainless steel. For example 304L would be suitable. Some
investigations about new materials are ongoing, but the baseline is 304L. We
stick to the baseline.

Vacuum firing [5]-[8] might be necessary to remove hydrogen from the bulk of
the pipe walls. Vacuum firing would have to be done prior to the production of
the pipes. It has no direct impact on the production technology.

The pipes have to withstand the environmental pressure. Including a reasonable
safety margin, a wall thickness of $7~{}\mathrm{mm}$ to $8~{}\mathrm{mm}$ is
necessary. We follow an approach, which is used by the current gravitational
wave detectors LIGO and VIRGO [9, 10]. We adopt walls of $3~{}\mathrm{mm}$ to
$4~{}\mathrm{mm}$ thickness reinforced by stiffener rings. This approach eases
the weight and therefore simplifies the handling, reduces the cost of
material, and facilitates the bending of the initial material into pipes.
Those stiffener rings are not included in our production scheme. They will
have to be added later.

Reliability is probably the most critical requirement for the pipe production.
The occurence of a leak will require manual intervention for repair, with a
correspondingly large effort for the localisation of the leak and its repair.
We try to limit these interventions to a few for the whole project.

There will be four pipes running in parallel in each tunnel. Fig. 1 shows a
typical tunnel cross section.

Figure 1: A typical layout of the vacuum pipes in a tunnel with
$6.5~{}\mathrm{m}$ diameter. The yellow area around the pipes indicates free
space for access to the pipes and for insulation.

The UHV-pipes will be built in sections of $500~{}\mathrm{m}$ length. They
will be connected to each other through pumping stations. The pumping stations
include the pumps, but also instrumentation for diagnosis and valves.

###  2.2 Production Concept

####  2.2.1 State of the art in tube and pipeline construction

Conventionally the vacuum pipes would be built from prefabricated sections of
$15~{}\mathrm{m}$ to $20~{}\mathrm{m}$ of length. The sections would be bent
from sheet metal, welded, leak checked, and cleaned in a factory, then
transported on site, lowered through one of the access shafts into the
tunnels, carefully aligned to each other, clamped and butt welded into the
final pipe. The vacuum system of the current generation of gravitational wave
detectors was built in this way. The length of the sections is limited by
transportation, if they are produced in a remote factory, or by the space
restrictions in lowering the pipes into the tunnels for an on-site production
facility. For the Einstein-Telescope we would need between $6\,000$ and
$8\,000$ sections. The cost of the vacuum system will be dominated by the
personnel required for the handling, especially for the welding in the
tunnels. We wanted to avoid the thousands of welding joints between pipe
sections and improve the reliability of the welding of the pipes by a single
seam weld along the length of the pipes. This was the starting point of our
project.

####  2.2.2 An alternative approach

We propose to transport rolls of sheet metal into the tunnels and to develop a
fully automatised robot, that produces pipe sections of up to
$500~{}\mathrm{m}$ length in the tunnels in the position were the pipes will
be installed eventually.

##### Raw material

For UHV applications the standard material is stainless steel of type 304L or
316L. We base our concept on this material. The wall thickness will be
$3~{}\mathrm{mm}$ to $4~{}\mathrm{mm}$, just sufficient to withstand the
atmospheric pressure. Stiffener rings around the outer circumference of the
pipes will be added later to stabilise the pipes. The sheet metal will be
produced by a steel company in coils (see fig. 2). Typical dimension of the
coils are an inner diameter of $610~{}\mathrm{mm}$ and an outer diameter of
$1\,715~{}\mathrm{mm}$. At a thickness of $4~{}\mathrm{mm}$ a length of the
sheets of $510~{}\mathrm{m}$ is possible. In standard production processes
these coils are produced in widths of up to $2~{}\mathrm{m}$. For a
$1~{}\mathrm{m}$ pipe we will need $3.14~{}\mathrm{m}$ of circumference. We
assume to purchase coils of half of this width and weld them together to the
full width prior to the underground production. The weight of a
$1.57~{}\mathrm{m}$ wide coil is approximately $25~{}\mathrm{t}$. About $500$
of these coils are needed.

Figure 2: Roll of sheet metal as produced by a steel company.

##### Preparation on Surface

Upon arrival on site, a few steps of preparation need to be done. These
require an un-coiling and re-coiling of the material. The edges of the coil
need to be prepared for the welding seam. We might want to vacuum-fire the
material to reduce the outgasing of hydrogen. Then, two coils will be welded
together to create a coil with the full width of $3.14~{}\mathrm{m}$. Finally
the material is precleaned, bagged, and transported underground with dedicated
transport machinery. The elevators are easily capable of taking a
$50~{}\mathrm{t}$ coil.

##### Pipe Fabrication Underground

In the tunnel material from a full coil is fed into a fully automatised robot.
There are two basic options for welding pipes: a spiral weld or a longitudinal
seam. We decided for a longitudinal seam. It simplifies the welding and avoids
having to rotate the bulky pipe sections around their axis. Fig. 3 shows a
sketch of the bending scheme. We will use dry bending to avoid contamination
of the pipe sections with hydrocarbons. The cross section is determined by the
size of the coils. The machinery for bending and welding will fit behind this
cross section. We estimated a length of the robot of about $30~{}\mathrm{m}$.

The pipe is welded from inside out. This has two advantages: The cleaner
surface will face the vacuum while the root of the weld is on the outside,
where defects, if they do occur, mostly will be located. The seam is located
at the bottom of the pipe. It is welded in gravity position, which is the most
stable position for welding. To arrive at the position of the weld, the
welding head is attached to the end of a cantilever, which reaches into the
forming pipe from the direction of the coil. The seam between the original
half-coils will be located on the top of the pipe. It will have its root on
the outside, too.

Figure 3: Bending of the metal sheet into a pipe.

We investigated a number of different welding technologies. We think gas metal
arc welding is a viable option, but we decided for laser welding under vacuum.
There are several arguments in favour of welding in a vacuum: A very stable
weld pool can be achieved in a vacuum, leading to a significantly improved
seam quality. Lowering the atmospheric pressure during the weld introduces
less hydrogen into the weld and avoids oxidation of the weld area. With an
effective energy input and lower vacuum requirements than, for example,
electron beam welding, a very stable process can be achieved.

For the application in our concept a mobile vacuum is required. Fig. 4 shows
the welding head of such a system. Strong pumps reduce the pressure to a few
$\mathrm{mbar}$ underneath the metal cup surrounding the laser head. In our
case the rim of the cup will have to be adapted to the shape of the pipe and
the cup will have to slide along the seam as the welding progresses. We
believe that a pressure of 30 $\mathrm{mbar}$ can be reached in this
configuration. But, neither the non-flat geometry nor the sliding of the cup
on a pipe has been tested, yet. The production speed will be limited by the
welding. A speed of $1~{}\mathrm{m}/\mathrm{min}$ should be possible, It
allows to manufacture a $500~{}\mathrm{m}$ section in a single day.

Figure 4: Example of a head for laser welding in vacuum (source: U. Reisgen et
al., ”Laser beam in mobile vacuum”. Proceedings: Lasers in Manufacturing - LIM
2017, Munich, Germany).

There are three potential options for the position of the robot.

  * 1.

The robot can be mounted stationary in a cavern at the end of the tunnel. The
pipe section will be pushed from the machine into the tunnel as its length
increases. A support of the pipe is required that allows the pipe to slide
along the tunnel with sufficiently low friction. Once a section is completed,
it will be pushed further along the tunnel into its final position.

  * 2.

The robot is located on rails in the tunnel. It produces the pipe sections
directly in their envisaged position along the tunnel. It starts at the end of
the section and slowly moves away from this position as the pipe section grows
in length. The pipe section is put on a hydraulic stand, which allows to lift
it into its final position in the tunnel.

  * 3.

A combination of the first two options: The robot is movable on rails between
productions, but stationary during the production of a single section. It is
moved to the end of the section and pushes the section on sliding fixtures
into its final position as it is produced. This is the most complicated
option, but it avoids sliding pipe sections over kilometres and has less
constraints on space. It produces all four pipe sections from a position in
the empty part of the tunnel while it fills the tunnel behind with four pipes.

The second option is our preferred solution, but the other two options are
viable, too. After the individual sections are completed they have to be
connected into the final pipe and cleaned from inside with a detergent-
spraying robot.

A $500~{}\mathrm{m}$ pipe section will extent during bake-out by roughly
$1.5~{}\mathrm{m}$. We assume that each section will be fixed in the tunnel at
its centre and extent towards the ends. Every position except for the centre
must be able to slide in its fixtures. In between sections a bellow must be
installed, capable of absorbing $1.5~{}\mathrm{m}$. To install pumps, gauges,
etc. wholes will be extruded from the pipe walls in certain positions and
flanges welded to it. We expect, pumping will be necessary every
$500~{}\mathrm{m}$.

This is a preliminary concept. A few aspects are not taken into account, yet.
These are:

  * 1.

The addition of stiffener rings on the outer circumference of the pipes.

  * 2.

Several methods are available to control the quality of the seam weld, for
example ultrasonic or X-ray imaging to show inner defects, seam tracking to
detect geometrical intolerances or eddy current measurements. We did not
implement any in our concept, yet. Eventually a leak-test of the completed
pipe sections is required.

  * 3.

We did not work out the details of the cleaning robot, yet.

  * 4.

Baffles need to be installed along the pipes to absorb scattered light and to
prevent its back-scattering into the laser beam. In the order of $100$ to
$200$ of these baffles will be required in each $10~{}\mathrm{km}$ arm. The
installation procedure is still unclear. We envisage a manual installation
after cleaning of the completed pipes.

  * 5.

The pipes need to be baked at temperatures between $80$ and
$120^{\circ}~{}\mathrm{C}$ to remove water molecules from the inner surfaces.
We assume that the same procedures can be applied as for conventionally
produced pipes, but no details are worked out, yet.

  * 6.

Despite the high resistance of stainless steel against corrosion, some
corrosion protection of the outside of the pipes will be needed. They have to
survive for fifty years in the humid atmosphere of the tunnels. The pipes will
have to be painted or coated. The procedure still needs to be worked out.

####  2.2.3 Evaluation

##### Risk Assessment

Our concept is based on a new production technology. Extensive prototyping
will be necessary to verify the validity of the new concept. Especially as it
includes a welding technology which is state-of-the-art, but has not been used
for pipe welding before. Should it turn out, that laser welding under vacuum
is not feasible, GMA-Hybrid welding or gas metal arc welding would be
alternatives. The partners of the project are preparing for the prototyping.

The largest risk in production are vacuum leaks. Let us first remember that we
are minimising the amount of welding with our concept. The highest risk in a
conventional system are the connection welds, where pipes, that never are
perfectly round and seldom have a constant diameter have to be aligned
carefully for each connection. Specialised clamping devices are used to force
the pipe ends into the same shape and to a precise distance between the weld
faces. Deviations from the ideal weld preparation (which are not unlikely due
to the complicated alignment process) can lead to weld defects, that must be
detected and repaired. With our concept we reduce the number of these welds
from around $6\,000$ to 8 000 to only $240$. The longitudinal weld along the
seam of the pipe sections reduces the length of the welding further, compared
to the commonly used spiral welding.

But the largest risk still are welding defects. We will use a fully
automatised continuous welding process. Most defects appear at the beginning
of the welding process when parameters are not stable yet. The long welds
reduce the number of these ramp-up phases and therefor also the risk of
defects. State-of-the-art monitoring of the quality will be applied to
identify potential welding problems. Smaller defects can be repaired in situ
with the usual techniques. Unfortunately, if defective sections must be
replaced, this procedure is difficult and needs large effort. Leak detection
can only be done in the tunnel, after a section is completed.

##### Production Time

The speed of production will be limited by the continuous welding of the
longitudinal seam. At least $0.5~{}\mathrm{m}/\mathrm{min}$ can be achieved,
but also $1~{}\mathrm{m}/\mathrm{min}$ is plausible, so that a
$500~{}\mathrm{m}$-section can be produced in one day within two shifts. Prior
to the production, the machine has to be prepared. For example, it has to be
loaded with all the materials. After the production of a section, some service
and maintenance will be required on the machine. All together, we assume that
it will take three days for each section. We estimate that a team of five
people is needed for the operations.

In total, 240 sections of $500~{}\mathrm{m}$ are needed. They can be produced
in two to three years with a single machine. Parallel production in different
tunnels is possible, if more than one machine and the corresponding personnel
is made available.

##### Interference with other Installations

During the production of pipe sections in a tunnel, this tunnel will be mostly
blocked for other installations. Installation of the completed sections in
their final position and joining of the sections should be possible though.
Some space must be reserved in an adjacent cavern for logistics. The tunnel
must be kept free for transportation of the materials.

Above ground a space will be needed for the reception of the coils and other
materials, their uncoiling, pre-welding, cleaning, etc. Space will also be
needed for the storage of materials.

##### Cost

The density of the material is $7.93~{}\mathrm{g}/\mathrm{cm}^{3}$ (304L or
316L). With a wall thickness of $4~{}\mathrm{mm}$ and a diameter of
$1~{}\mathrm{m}$, the weight of the pipes is $99.6~{}\mathrm{kg}/\mathrm{m}$
or about $50~{}\mathrm{t}$ per section. The cost of the raw material in
October 2022 was about $3.5$ €$/\mathrm{kg}$ for 304L and $5$ €$/\mathrm{kg}$
for 316L. This gives a total cost of the raw material between $42$ mio. € and
$61$ mio. €. It should be noted, that prices are subject to fluctuations and
that the estimates are only up-to-date on a daily basis.

For the whole system $480$ coils of $25~{}\mathrm{t}$ each need to be
transported to the site. The cost of transportation depends on the distance
from the steel manufacturer to the site. The coils could be transported by
train or by trucks. For any realistic distance the cost of transportation will
be small compared to the cost of the material. We will ignore it in the
following.

We estimate that the cost of the main machine (coil forming and welding) is
about $5$ mio. €. The cost of the machine for the pre-fabrication above ground
is estimated at another $2.5$ mio. €. Special machines for the transportation
of the coils will be needed. We estimate $1.25$ mio. € for five of those
machines. The elevator will have to be modified for the transportation of the
coils and the machinery. This might cost another $1$ mio. €. Overall the basic
investment will be around $10$ mio. €.

Teams of five people will operate the production. For a continuous production
three of these teams will be needed at a cost of about $1.1$ mio. € anually. A
back office at a cost of $0.75$ mio. € per year has to be added. This sums up
to around $5$ mio. € for personnel for a production time of two to three
years.

Please keep in mind, that this is a preliminary concept. We only laid out the
most important production steps. The cost of these are dominated by the cost
of the material. We believe that the missing cost for stiffeners, cleaning,
etc. will not change this statement significantly.

##  3 Conclusions and Next Steps

We developed a new production technology for long vacuum pipes. It is based on
a fully automatised machine moving through the tunnels, producing continuous
vacuum pipes from a coil of raw material at a speed of about $1$ meter per
minute. We propose to use laser welding under vacuum to close the longitudinal
seam of the pipe sections.

The new production concept looks promising. Effort and cost for transportation
and production are potentially lower compared to an off-site production of
transportable sections. The fully automatised long welding seams promise an
improved reliability of the production.

Here we described the core of the new production technology. Many details,
such as the mounting of the stiffeners, the cleaning of the pipes, or the
installation of baffles are not worked out, yet. The partners of this project
are now preparing to work out all the details and to produce a first
prototype.

## Acknowledgement

This project was executed as a cross-border cooperation of SMEs within the
interreg project ET2SMEs (https://et2smes.eu/). The ET2SMEs project is carried
out under the Interreg V-A Euregio Meuse-Rhine Programme, with financial
support from the European Regional Development Fund (ERDF).

## References

  * [1] M. Punturo et al., The Einstein Telescope: A third-generation gravitational wave observatory, Class. Quant. Grav. 27 (2010) 194002.
  * [2] ET steering committee, Einstein Telescope: Science Case, Design Study and Feasibility Report, available from https://apps.et-gw.eu/tds/?content=3&r=17196.
  * [3] ET steering committee, Design Report Update 2020 for the Einstein Telescope, available from https://apps.et-gw.eu/tds/?content=3&r=17245.
  * [4] M. Maggiore et al., Science Case for the Einstein Telescope, arXiv:1912.02622 [astro-ph.CO]; https://arxiv.org/abs/1912.02622.
  * [5] R. Calder and G. Lewin, Reduction of stainless-steel outgassing in ultra-high vacuum, Br. J. Appl. Phys. 18 (1967) 1459, https://iopscience.iop.org/article/10.1088/0508-3443/18/10/313/.
  * [6] G. Grosse and G. Messer, Extreme High Vacuum, Proc. of the 7th Int. Vac. Cong., Vienna, 1977, p. 223 https://cds.cern.ch/record/455565/files/p213.pdf.
  * [7] G. Grosse and G. Messer, Proc. of the 8th Int. Vac. Cong., Cannes, 1980, p. 399.
  * [8] J. Kamiya et al., Evaluation of Vacuum Firing Effect on Stainless Steel from Vacuum and Surface Point of View, e-Journal of Surface Science and Nanotechnology 21 (2023), https://doi.org/10.1380/ejssnt.2023-015.
  * [9] M. Zucker and F. Dylla, LIGO Vacuum System Design, Construction, Features and Faults, LIGO public note, 2019, https://dcc.ligo.org/LIGO-G1900137/public.
  * [10] V. Brisson et al., Ultra-high vacuum qualification of the prototype module for the 2×3 km arms of the VIRGO interferometer, Vacuum 60 (2001) 9, https://www.sciencedirect.com/science/article/pii/S0042207X00002384.