<figure><img src="image_1.png" /></figure>

Cosmic-ray diffusion in two local filamentary clouds

F.R Kamal Youssef, Isabelle A Grenier

► To cite this version:

F.R Kamal Youssef, Isabelle A Grenier. Cosmic-ray diffusion in two local filamentary clouds. Astronomy and Astrophysics - A&A, 2024, 685, pp.A102. 10.1051/0004-6361/202348299. hal-04510892

HAL Id: hal-04510892

https://hal.science/hal-04510892

Submitted on 17 May 2024

**HAL** is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_1.png" /></figure>

# Cosmic-ray diffusion in two local filamentary clouds

F.R Kamal Youssef, Isabelle A Grenier

► To cite this version:

F.R Kamal Youssef, Isabelle A Grenier. Cosmic-ray diffusion in two local filamentary clouds. Astronomy & Astrophysics - A&A, 2024, 685, pp.A102. 10.1051/0004-6361/202348299. hal-04510892

HAL Id: hal-04510892

https://hal.science/hal-04510892v1

Submitted on 17 May 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Cosmic-ray diffusion in two local filamentary clouds

F. R. Kamal Youssef and I. A. Grenier

Université Paris Cité and Université Paris-Saclay, CEA, CNRS, AIM, 91190 Gif-sur-Yvette, France
e-mail: isabelle.grenier@cea.fr

Received 17 October 2023 / Accepted 14 February 2024

## ABSTRACT

**Context.** Hadronic interactions between cosmic rays (CRs) and interstellar gas have been probed in γ rays across the Galaxy. A fairly uniform CR distribution is observed up to a few hundred parsecs from the Sun, except in the Eridu cloud, which shows an unexplained 30–50% deficit in GeV to TeV CR flux.

**Aims.** To explore the origin of this deficit, we studied the Reticulum cloud, which shares notable traits with Eridu: a comparable distance in the low-density region of the Local Valley and a filamentary structure of atomic hydrogen extending along a bundle of ordered magnetic-field lines that are steeply inclined to the Galactic plane.

**Methods.** We measured the γ-ray emissivity per gas nucleon in the Reticulum cloud in the 0.16–63 GeV energy band using 14 years of Fermi-LAT data. We also derived interstellar properties that are important for CR propagation in both the Eridu and Reticulum clouds, at the same parsec scale.

**Results.** The γ-ray emissivity in the Reticulum cloud is fully consistent with the average spectrum measured in the solar neighbourhood, but this emissivity, and therefore the CR flux, is 1.57±0.09 times larger than in Eridu across the whole energy band. The difference cannot be attributed to uncertainties in gas mass. Nevertheless, we find that the two clouds are similar in many respects: both have magnetic-field strengths of a few micro-Gauss in the plane of the sky; both are in approximate equilibrium between magnetic and thermal pressures; they have similar turbulent velocities and sonic Mach numbers; and both show magnetic-field regularity with a dispersion in orientation lower than 10⁻⁷–15⁻⁷ over large zones. The gas in Reticulum is colder and denser than in Eridu, but we find similar parallel diffusion coefficients around a few times 10²⁸ cm² s⁻¹ in both clouds if CRs above 1 GV in rigidity diffuse on resonant, self-excited Alfvén waves that are damped by ion-neutral interactions.

**Conclusions.** The loss of CRs in Eridu remains unexplained, but these two clouds provide important test cases to further study how magnetic turbulence, line tangling, and ion-neutral damping regulate CR diffusion in the dominant gas phase of the interstellar medium.

**Key words.** ISM: clouds – cosmic rays – ISM: magnetic fields – gamma rays: ISM

## 1. Introduction

The Sun is located in a valley of low gas density surrounded by two walls of dense molecular clouds (see Fig. 2 of Leike et al. 2020). Part of the valley, out to~150 pc or ~250 pc depending on the direction, corresponds to the Local Bubble (see Fig. 1 of Pelgrims et al. 2020), but the elongated valley extends further away in opposite directions roughly parallel to the local spiral arm. Interactions of cosmic rays (CRs) with the interstellar gas and radiation fields present inside and around the valley produce γ rays that have been recorded at energies above 30 MeV by the Large Area Telescope (LAT) on board the Fermi Gamma-Ray Space Telescope. As the hadronic contribution to this γ-ray intensity integrates the product of the CR and gas volume densities along the lines of sight, the γ-ray map of nearby clouds can be used to probe the CR flux in different locations and in different types of environment (Lebrun et al. 1982). This method requires an accurate measurement of the gas mass in a cloud.

In the warm atomic phase of clouds, often referred to as the warm neutral medium (WNM), one can directly infer the gas mass from observations of the Hi 21 cm lines (e.g. Dickey et al. 2000; Heiles & Troland 2003). In the denser, cold atomic phase of the clouds, where the Hi brightness temperature probes the excitation (spin) temperature of the gas rather than the column density, the lines must be corrected for partial self-absorption in

order to estimate the gas mass (Murray et al. 2015, and references therein). This cold atomic phase, referred to as the cold neutral medium (CNM), generally contributes a fraction between 5% and 28% of the total atomic mass (Murray et al. 2018, 2020). H₂ molecules radiate inefficiently in the cold conditions of molecular clouds, so we must rely on observations of surrogate lines, such as the rotational lines of ¹²CO and its isotopologues, in order to quantify the H₂ mass. This is only possible where the abundance and excitation temperature of CO molecules are sufficient. Environmental variations of the XCO conversion factor between the ¹²CO(1–0) line intensity and the H₂ column density have been studied theoretically and observationally (Bolatto et al. 2013; Gong et al. 2018; Remy et al. 2017, and references therein), but the variations are too large and too uncertain to reliably use the CO data to measure the CR flux in the dense molecular phase of clouds.

Complementary information can be obtained from the thermal emission of large dust grains mixed with the gas. The emission is observed at submillimetre and infrared wavelengths (Planck Collaboration XI 2014; Planck Collaboration Int. XVII 2014). Dust grains and CRs are both tracers of the total gas in all chemical and thermodynamical forms (Grenier et al. 2005; Planck Collaboration Int. XXVIII 2015). The γ rays reveal the gas mass convolved by the ambient CR flux, whereas the dust optical depth shows the gas mass convolved by the ambient
# Cosmic-ray diffusion in two local filamentary clouds

F. R. Kamal Youssef and I. A. Grenier

Université Paris Cité and Université Paris-Saclay, CEA, CNRS, AIM, 91190 Gif-sur-Yvette, France
e-mail: isabelle.grenier@cea.fr

Received 17 October 2023 / Accepted 14 February 2024

## ABSTRACT

**Context.** Hadronic interactions between cosmic rays (CRs) and interstellar gas have been probed in γ rays across the Galaxy. A fairly uniform CR distribution is observed up to a few hundred parsecs from the Sun, except in the Eridu cloud, which shows an unexplained 30–50% deficit in GeV to TeV CR flux.

**Aims.** To explore the origin of this deficit, we studied the Reticulum cloud, which shares notable traits with Eridu: a comparable distance in the low-density region of the Local Valley and a filamentary structure of atomic hydrogen extending along a bundle of ordered magnetic-field lines that are steeply inclined to the Galactic plane.

**Methods.** We measured the γ-ray emissivity per gas nucleon in the Reticulum cloud in the 0.16–63 GeV energy band using 14 years of Fermi-LAT data. We also derived interstellar properties that are important for CR propagation in both the Eridu and Reticulum clouds, at the same parsec scale.

**Results.** The γ-ray emissivity in the Reticulum cloud is fully consistent with the average spectrum measured in the solar neighbourhood, but this emissivity, and therefore the CR flux, is 1.57±0.09 times larger than in Eridu across the whole energy band. The difference cannot be attributed to uncertainties in gas mass. Nevertheless, we find that the two clouds are similar in many respects: both have magnetic-field strengths of a few micro-Gauss in the plane of the sky; both are in approximate equilibrium between magnetic and thermal pressures; they have similar turbulent velocities and sonic Mach numbers; and both show magnetic-field regularity with a dispersion in orientation lower than 10⁻⁷–15⁻⁷ over large zones. The gas in Reticulum is colder and denser than in Eridu, but we find similar parallel diffusion coefficients around a few times 10²⁸ cm² s⁻¹ in both clouds if CRs above 1 GV in rigidity diffuse on resonant, self-excited Alfvén waves that are damped by ion-neutral interactions.

**Conclusions.** The loss of CRs in Eridu remains unexplained, but these two clouds provide important test cases to further study how magnetic turbulence, line tangling, and ion-neutral damping regulate CR diffusion in the dominant gas phase of the interstellar medium.

**Key words.** ISM: clouds – cosmic rays – ISM: magnetic fields – gamma rays: ISM

## 1. Introduction

The Sun is located in a valley of low gas density surrounded by two walls of dense molecular clouds (see Fig. 2 of Leike et al. 2020). Part of the valley, out to~150 pc or ~250 pc depending on the direction, corresponds to the Local Bubble (see Fig. 1 of Pelgrims et al. 2020), but the elongated valley extends further away in opposite directions roughly parallel to the local spiral arm. Interactions of cosmic rays (CRs) with the interstellar gas and radiation fields present inside and around the valley produce γ rays that have been recorded at energies above 30 MeV by the Large Area Telescope (LAT) on board the Fermi Gamma-Ray Space Telescope. As the hadronic contribution to this γ-ray intensity integrates the product of the CR and gas volume densities along the lines of sight, the γ-ray map of nearby clouds can be used to probe the CR flux in different locations and in different types of environment (Lebrun et al. 1982). This method requires an accurate measurement of the gas mass in a cloud.

In the warm atomic phase of clouds, often referred to as the warm neutral medium (WNM), one can directly infer the gas mass from observations of the Hi 21 cm lines (e.g. Dickey et al. 2000; Heiles & Troland 2003). In the denser, cold atomic phase of the clouds, where the Hi brightness temperature probes the excitation (spin) temperature of the gas rather than the column density, the lines must be corrected for partial self-absorption in

order to estimate the gas mass (Murray et al. 2015, and references therein). This cold atomic phase, referred to as the cold neutral medium (CNM), generally contributes a fraction between 5% and 28% of the total atomic mass (Murray et al. 2018, 2020). H₂ molecules radiate inefficiently in the cold conditions of molecular clouds, so we must rely on observations of surrogate lines, such as the rotational lines of ¹²CO and its isotopologues, in order to quantify the H₂ mass. This is only possible where the abundance and excitation temperature of CO molecules are sufficient. Environmental variations of the XCO conversion factor between the ¹²CO(1–0) line intensity and the H₂ column density have been studied theoretically and observationally (Bolatto et al. 2013; Gong et al. 2018; Remy et al. 2017, and references therein), but the variations are too large and too uncertain to reliably use the CO data to measure the CR flux in the dense molecular phase of clouds.

Complementary information can be obtained from the thermal emission of large dust grains mixed with the gas. The emission is observed at submillimetre and infrared wavelengths (Planck Collaboration XI 2014; Planck Collaboration Int. XVII 2014). Dust grains and CRs are both tracers of the total gas in all chemical and thermodynamical forms (Grenier et al. 2005; Planck Collaboration Int. XXVIII 2015). The γ rays reveal the gas mass convolved by the ambient CR flux, whereas the dust optical depth shows the gas mass convolved by the ambient
dust-to-gas mass ratio and by the ambient grain emissivity. Both the CR and dust properties vary from place to place, but the spatial correlation between the two tracers in a given cloud can be used to trace the dark gas that escapes HI and CO line observations (Grenier et al. 2005). This dark neutral medium (DINM) accumulates at the $HI-HI$ transition in the clouds (Liszt et al. 2019). It is not seen in HI or CO emission because of self-absorption in the dense HI clumps and because CO molecules are efficiently photo-dissociated or too weakly excited in the diffuse H₂ gas.

Given the uncertainty in gas mass estimates in the dark and molecular phases of clouds, one can only infer the CR spectrum with adequate precision in the atomic phase of the interstellar medium (ISM). The derivation of the γ-ray emissivity, $q_{HI}$, per gas nucleon in the atomic gas relies on the spatial correlation between the γ-ray intensity and HI column density ($N_{HI}$), but the analysis must include maps of the other gas phases in order to avoid biasing the HI emissivity spectrum (Planck Collaboration Int. XXVIII 2015). These spectra are typically obtained at energies above a few hundred MeV in order to take advantage of the better angular resolution of the LAT in order to resolve the clouds and their different phases.

Such emissivity measurements have been performed in different clouds of the local ISM (e.g. Abdo et al. 2010; Ackermann et al. 2012a; Planck Collaboration Int. XXVIII 2015; Casandjian 2015; Tibaldo et al. 2015; Remy et al. 2017; Joubaud et al. 2019) and across the Milky Way (e.g. Ackermann et al. 2013; Acero et al. 2016). The local emissivity spectra appear to be often consistent with the local average, $q_{LIS}(E)$, which has been obtained for the whole atomic gas seen at Galactic latitudes ranging between 7° and 70° (Casandjian 2015). This consistency implies that the CR flux is rather uniform within a few hundred parsecs of the Sun (Grenier et al. 2015).

At variance with this uniformity, a significantly lower γ-ray emissivity spectrum has been found in a nearby cloud named Eridu. The spectrum is 34% lower than the $q_{LIS}$ average (Joubaud et al. 2020). The deviation is statistically significant (14 σ) and there is no confusion with other clouds along the lines of sight. The deviation cannot be attributed to uncertainties in the HI mass because the minimum amount of gas inferred from the HI lines (optically thin case) gives an upper limit to the emissivity that is still 25% below $q_{LIS}$. More CRs would outshine the observed γ rays. The γ-ray analysis has probed CRs with energies above several GeV. Their penetration in the cloud is likely unhampered given the modest column densities (1–5×10²⁰ cm⁻²) and low volume densities (around 7 cm⁻³) of the diffuse gas (Joubaud et al. 2019). The Eridu cloud therefore appears to be pervaded by a ~30% lower CR flux, but with the same particle energy distribution as in the local ISM.

The Eridu cloud is a diffuse, atomic, and filamentary cirrus. It is located outside the edge of the Orion-Eridanus superbubble, at high Galactic latitudes of around 50° and at an altitude of about 200–250 pc below the Galactic plane (Joubaud et al. 2019). Figure 1 shows that the filament is inclined with respect to the Galactic plane and the dust polarisation data at 353 GHz from Planck (Planck Collaboration Int. XIX 2015) show that the filament extends along a bundle of well-ordered magnetic-field lines, pointing within 45° of the southern Galactic pole.

Magnetic loops are pushed out of the Galactic disc by gas fountains. In order to investigate their impact on CR transport to the halo, we searched for another filament with gaseous and magnetic characteristics comparable to those of Eridu. The Reticulum cloud is also an atomic and filamentary cloud that lies at the edge of the Local Valley, at a distance we estimate to be

<figure><img src="image_3.png" /><figcaption>**Fig. 1.** Dust optical depth $\tau_{353}$ at 353 GHz in Galactic coordinates towards the Eridu (top) and Reticulum (bottom) filaments. The black lines give the orientation of the plane-of-the-sky $B_{sky}$ magnetic field inferred from the dust polarisation observations at 353 GHz and at an angular resolution of 14' (FWHM). The red contours delineate the filaments at HI gas column densities of 3 × 10²⁰ cm⁻² and 1.2×10²⁰ cm⁻², respectively.</figcaption></figure>

about 240 pc (Appendix A). Figure 1 shows that, like Eridu, it is largely inclined with respect to the Galactic plane. It points towards the southern Galactic halo and has well-ordered magnetic field lines aligned with the filament axis. We therefore aim to measure the γ-ray emissivity spectrum per gas nucleon in the Reticulum filament and to compare it with the average in the local ISM and with the emissivity spectrum in the Eridu cloud.

The paper is structured as follows: The data are presented in Sect. 2 and the models are described in Sect. 3. The results and the detection of dark gas are detailed in Sect. 4. The CR content of the clouds is discussed in Sects. 5.1–5.5.

## 2. Data

We selected an analysis region around the Reticulum filament extending from 251.5 to 302.5 in Galactic longitude and from -54.5 to -30.5 in Galactic latitude, as shown in Fig. 2. We masked regions with large column densities from the Large Magellanic Cloud (LMC) and Small Magellanic Cloud (SMC), and from the edge of the Gum Nebula.

### 2.1. Dust data

In order to trace the dust column density, we used the all-sky map of the dust optical depth, $\tau_{353}$, at 353 GHz (Planck Collaboration XI 2014, updated by Irfan et al. 2019). It was constructed at an angular resolution of 5' by modelling the intensity, $I_v$, of the dust thermal emission recorded in frequency $\nu$ by Planck and the Infrared Astronomical Satellite (IRAS). The dust optical depth relates to the modified black-body intensity $B_v(T_d)$ for a dust temperature $T_d$ and spectral index, $\beta$, according to
dust-to-gas mass ratio and by the ambient grain emissivity. Both the CR and dust properties vary from place to place, but the spatial correlation between the two tracers in a given cloud can be used to trace the dark gas that escapes HI and CO line observations (Grenier et al. 2005). This dark neutral medium (DINM) accumulates at the $HI-HI$ transition in the clouds (Liszt et al. 2019). It is not seen in HI or CO emission because of self-absorption in the dense HI clumps and because CO molecules are efficiently photo-dissociated or too weakly excited in the diffuse $H_2$ gas.

Given the uncertainty in gas mass estimates in the dark and molecular phases of clouds, one can only infer the CR spectrum with adequate precision in the atomic phase of the interstellar medium (ISM). The derivation of the γ-ray emissivity, $q_{HI}$, per gas nucleon in the atomic gas relies on the spatial correlation between the γ-ray intensity and HI column density ($N_{HI}$), but the analysis must include maps of the other gas phases in order to avoid biasing the HI emissivity spectrum (Planck Collaboration Int. XXVIII 2015). These spectra are typically obtained at energies above a few hundred MeV in order to take advantage of the better angular resolution of the LAT in order to resolve the clouds and their different phases.

Such emissivity measurements have been performed in different clouds of the local ISM (e.g. Abdo et al. 2010; Ackermann et al. 2012a; Planck Collaboration Int. XXVIII 2015; Casandjian 2015; Tibaldo et al. 2015; Remy et al. 2017; Joubaud et al. 2019) and across the Milky Way (e.g. Ackermann et al. 2013; Acero et al. 2016). The local emissivity spectra appear to be often consistent with the local average, $q_{LIS}(E)$, which has been obtained for the whole atomic gas seen at Galactic latitudes ranging between 7° and 70° (Casandjian 2015). This consistency implies that the CR flux is rather uniform within a few hundred parsecs of the Sun (Grenier et al. 2015).

At variance with this uniformity, a significantly lower γ-ray emissivity spectrum has been found in a nearby cloud named Eridu. The spectrum is 34% lower than the $q_{LIS}$ average (Joubaud et al. 2020). The deviation is statistically significant ($14 \sigma$) and there is no confusion with other clouds along the lines of sight. The deviation cannot be attributed to uncertainties in the HI mass because the minimum amount of gas inferred from the HI lines (optically thin case) gives an upper limit to the emissivity that is still 25% below $q_{LIS}$. More CRs would outshine the observed γ rays. The γ-ray analysis has probed CRs with energies above several GeV. Their penetration in the cloud is likely unhampered given the modest column densities ($1-5 \times 10^{20} \text{ cm}^{-2}$) and low volume densities (around $7 \text{ cm}^{-3}$) of the diffuse gas (Joubaud et al. 2019). The Eridu cloud therefore appears to be pervaded by a ~30% lower CR flux, but with the same particle energy distribution as in the local ISM.

The Eridu cloud is a diffuse, atomic, and filamentary cirrus. It is located outside the edge of the Orion-Eridanus superbubble, at high Galactic latitudes of around 50° and at an altitude of about 200–250 pc below the Galactic plane (Joubaud et al. 2019). Figure 1 shows that the filament is inclined with respect to the Galactic plane and the dust polarisation data at 353 GHz from Planck (Planck Collaboration Int. XIX 2015) show that the filament extends along a bundle of well-ordered magnetic-field lines, pointing within 45° of the southern Galactic pole.

Magnetic loops are pushed out of the Galactic disc by gas fountains. In order to investigate their impact on CR transport to the halo, we searched for another filament with gaseous and magnetic characteristics comparable to those of Eridu. The Reticulum cloud is also an atomic and filamentary cloud that lies at the edge of the Local Valley, at a distance we estimate to be

<figure><img src="image_3.png" /><figcaption>**Fig. 1.** Dust optical depth $\tau_{353}$ at 353 GHz in Galactic coordinates towards the Eridu (top) and Reticulum (bottom) filaments. The black lines give the orientation of the plane-of-the-sky $B_{sky}$ magnetic field inferred from the dust polarisation observations at 353 GHz and at an angular resolution of 14' (FWHM). The red contours delineate the filaments at HI gas column densities of $3 \times 10^{20} \text{ cm}^{-2}$ and $1.2 \times 10^{20} \text{ cm}^{-2}$, respectively.</figcaption></figure>

about 240 pc (Appendix A). Figure 1 shows that, like Eridu, it is largely inclined with respect to the Galactic plane. It points towards the southern Galactic halo and has well-ordered magnetic field lines aligned with the filament axis. We therefore aim to measure the γ-ray emissivity spectrum per gas nucleon in the Reticulum filament and to compare it with the average in the local ISM and with the emissivity spectrum in the Eridu cloud.

The paper is structured as follows: The data are presented in Sect. 2 and the models are described in Sect. 3. The results and the detection of dark gas are detailed in Sect. 4. The CR content of the clouds is discussed in Sects. 5.1–5.5.

## 2. Data

We selected an analysis region around the Reticulum filament extending from 251.5 to 302.5 in Galactic longitude and from -54.5 to -30.5 in Galactic latitude, as shown in Fig. 2. We masked regions with large column densities from the Large Magellanic Cloud (LMC) and Small Magellanic Cloud (SMC), and from the edge of the Gum Nebula.

### 2.1. Dust data

In order to trace the dust column density, we used the all-sky map of the dust optical depth, $\tau_{353}$, at 353 GHz (Planck Collaboration XI 2014, updated by Irfan et al. 2019). It was constructed at an angular resolution of 5' by modelling the intensity, $I_v$, of the dust thermal emission recorded in frequency $\nu$ by Planck and the Infrared Astronomical Satellite (IRAS). The dust optical depth relates to the modified black-body intensity $B_v(T_d)$ for a dust temperature $T_d$ and spectral index, $\beta$, according to
<figure><img src="image_4.png" /><figcaption>**Fig. 2.** Dust optical depth measured at 353 GHz from the *Planck* and IRAS data and displayed on a 10.'x7.5 grid at the angular resolution of *Planck* (top panel) and of *Fermi*-LAT (middle). The bottom panel shows the γ-ray photon counts produced by CR interactions with gas and recorded by *Fermi*-LAT in 0.25×0.25 bins in the 0.16–63 GeV energy band (other emissions of non-gaseous origin are subtracted). Three zones are excluded from the analysis region: the SMC on the left, LMC in the middle, and the edge of the Gum nebula in the upper right. The Reticulum cloud is the long filament crossing the region below the LMC.</figcaption></figure>

$I_v = \tau_{353} B_v(T_a)(\nu/353 \text{ GHz})^\beta$. The dust opacity at 353 GHz is defined for a hydrogen column density $N_H$ as the ratio $\tau_{353}/N_H$.

Figure 2 shows the distribution of the dust optical depth in the analysis region at the original *Planck* resolution in the upper panel and at the angular resolution of *Fermi*-LAT for the 0.16–63 GeV band in the middle one.

## 22 Gamma-ray data

We used 14 yr of *Fermi*-LAT Pass 8 P8R3 survey data for energies between 0.16 and 63 GeV (Bruel et al. 2018) and we analysed the data in eight energy bands that are bounded by $10^{2.2}$, $10^{2.4}$, $10^{2.6}$, $10^{2.8}$, $10^{3.0}$, $10^{3.2}$, $10^{3.6}$, $10^{4.0}$, and $10^{4.8}$ MeV.

The *Fermi*-LAT point spread function (PSF) strongly varies with energy. In order to retain sufficient photon statistics in this faint region and to preserve a good contrast in photon counts across small cloud structures, we combined photons of different PSF types in each energy band. Photons are labelled from the worst (PSF0) to best (PSF3) quality in angle reconstruction. We used PSF 3 in the three lowest-energy bands, PSF 2 and 3 in the next band, PSF 1 to 3 in the next two bands, and PSF 0 to 3 in the two highest energy bands. The 68% containment radius of the resulting PSF varies from 0.9 above 630 MeV to 1.4 at 250 MeV and to 2.1 at 160 MeV. The resulting photon map is shown in the lower panel of Fig. 2 where we subtracted the point sources and the very smooth isotropic and inverse-Compton (IC) components to display the γ rays produced in the clouds.

In order to reduce to less than 10% the contamination of the photon data by residual CRs and by Earth atmospheric γ rays, we applied tight selection criteria such as the SOURCE class and photon arrival directions within an energy-dependent and PSF-type-dependent angle varying from 90° to 105° from the Earth zenith, namely 90° in the first energy band, 100° in the second, and 105° at higher energies. We used the associated instrument response functions P8R3_SOURCE_V3.

We modelled the interstellar emission in a broader region than the analysis one in order to take into account the extent of the *Fermi*-LAT PSF. The 5'-wide peripheral band around the analysis perimeter corresponds to the PSF >95% containment radius for our selections. It is displayed in Fig. 3.

The positions and flux spectra of the γ-ray sources in the field were provided by the *Fermi*-LAT 4FGL-DR3 12-yr source catalogue (Abdollahi et al. 2022). The analysis region and the 5'-wide peripheral band contain 220 point sources. Their flux spectra were calculated with the spectral characteristics given in the catalogue.

The observation also includes IC γ rays produced by CR electrons up-scattering the interstellar radiation field in the Milky Way. This was modelled by GALPROP<fnref n="1" /> and we used the LRYusifovXCO5z6R30_Ts150_Rs8 model version (Ackermann et al. 2012b). The isotropic spectrum for the extragalactic and residual instrumental backgrounds is given by the *Fermi*-LAT science support centre<fnref n="2" />. We used the isotropic spectrum provided for each PSF type of the SOURCE class.

## 23. Gas data

In order to map the neutral atomic hydrogen (Hα) in the region, we used the 21-cm Galactic All Sky Survey GASS III data (Kalberla & Haud 2015). The survey has an angular resolution of 16.2, a velocity resolution of 0.82 km s⁻<fnref n="1" /> in the local standard of rest, and an rms sensitivity of 57 mK. Lines in the analysis region span velocities from -125 to +495 km s⁻<fnref n="1" />.

Following the method developed by Planck Collaboration Int. XXVIII (2015), we decomposed the Hα spectra into sets of individual lines with Pseudo-Voigt profiles. Based on the (l, b, v) distribution of the fitted lines, we identified seven nearby cloud complexes that are coherent in velocity and in position across the region. We named these clouds according to the constellation in their direction, namely Reticulum (Ret), Far Hydrus (HyiF), North Hydrus (HyiN), South Hydrus (HyiS), East Dorado (DorE), West Dorado (DorW), and Horologium (Hor). They are displayed in Fig. 3 and they span velocities from 3.3 to 20.6 km s⁻<fnref n="1" />, -15.7 to -1.6 km s⁻<fnref n="1" />, -0.8 to 13.2 km s⁻<fnref n="1" />, -14.8 to 2.5 km s⁻<fnref n="1" />, 2.5 to 4.9 km s⁻<fnref n="1" />, 5.8 to 33.8 km s⁻<fnref n="1" />, and -14.8 to 5.8 km s⁻<fnref n="1" />, respectively. These clouds are local. We gathered lines from the Magellanic clouds in a single map. We added a ninth component gathering the very weak lines found at velocities below -8 km s⁻<fnref n="1" />.

We produced $N_{HI}$ column-density maps for each complex by integrating the fitted lines that have a central velocity within the chosen velocity interval for each (l, b) direction. This method corrects for potential line spill-over from one velocity interval (i.e. cloud) to the next. The line fits produce small residuals between the modelled and observed Hα spectra. In order to preserve the total Hα intensity, we distributed these residuals among

<footnote n="1">¹ http://galprop.stanford.edu</footnote>

<footnote n="2">² http://fermi.gsfc.nasa.gov/ssc/data/access/lat/BackgroundModels.html</footnote>
<figure><img src="image_4.png" /><figcaption>**Fig. 2.** Dust optical depth measured at 353 GHz from the *Planck* and IRAS data and displayed on a 10.'x7.5 grid at the angular resolution of *Planck* (top panel) and of *Fermi*-LAT (middle). The bottom panel shows the γ-ray photon counts produced by CR interactions with gas and recorded by *Fermi*-LAT in 0.25×0.25 bins in the 0.16–63 GeV energy band (other emissions of non-gaseous origin are subtracted). Three zones are excluded from the analysis region: the SMC on the left, LMC in the middle, and the edge of the Gum nebula in the upper right. The Reticulum cloud is the long filament crossing the region below the LMC.</figcaption></figure>

$I_v = T_{353} B_v(T_a)(\nu/353 \text{ GHz})^\beta$. The dust opacity at 353 GHz is defined for a hydrogen column density $N_H$ as the ratio $T_{353}/N_H$.

Figure 2 shows the distribution of the dust optical depth in the analysis region at the original *Planck* resolution in the upper panel and at the angular resolution of *Fermi*-LAT for the 0.16–63 GeV band in the middle one.

## 22 Gamma-ray data

We used 14 yr of *Fermi*-LAT Pass 8 P8R3 survey data for energies between 0.16 and 63 GeV (Bruel et al. 2018) and we analysed the data in eight energy bands that are bounded by $10^{2.2}$, $10^{2.4}$, $10^{2.6}$, $10^{2.8}$, $10^{3.0}$, $10^{3.2}$, $10^{3.6}$, $10^{4.0}$, and $10^{4.8}$ MeV.

The *Fermi*-LAT point spread function (PSF) strongly varies with energy. In order to retain sufficient photon statistics in this faint region and to preserve a good contrast in photon counts across small cloud structures, we combined photons of different PSF types in each energy band. Photons are labelled from the worst (PSF0) to best (PSF3) quality in angle reconstruction. We used PSF 3 in the three lowest-energy bands, PSF 2 and 3 in the next band, PSF 1 to 3 in the next two bands, and PSF 0 to 3 in the two highest energy bands. The 68% containment radius of the resulting PSF varies from 0.9 above 630 MeV to 1.4 at 250 MeV and to 2.1 at 160 MeV. The resulting photon map is shown in the lower panel of Fig. 2 where we subtracted the point sources and the very smooth isotropic and inverse-Compton (IC) components to display the γ rays produced in the clouds.

In order to reduce to less than 10% the contamination of the photon data by residual CRs and by Earth atmospheric γ rays, we applied tight selection criteria such as the SOURCE class and photon arrival directions within an energy-dependent and PSF-type-dependent angle varying from 90° to 105° from the Earth zenith, namely 90° in the first energy band, 100° in the second, and 105° at higher energies. We used the associated instrument response functions P8R3_SOURCE_V3.

We modelled the interstellar emission in a broader region than the analysis one in order to take into account the extent of the *Fermi*-LAT PSF. The 5'-wide peripheral band around the analysis perimeter corresponds to the PSF >95% containment radius for our selections. It is displayed in Fig. 3.

The positions and flux spectra of the γ-ray sources in the field were provided by the *Fermi*-LAT 4FGL-DR3 12-yr source catalogue (Abdollahi et al. 2022). The analysis region and the 5'-wide peripheral band contain 220 point sources. Their flux spectra were calculated with the spectral characteristics given in the catalogue.

The observation also includes IC γ rays produced by CR electrons up-scattering the interstellar radiation field in the Milky Way. This was modelled by GALPROP<fnref n="1" /> and we used the LRYusifovXCO5z6R30_Ts150_Rs8 model version (Ackermann et al. 2012b). The isotropic spectrum for the extragalactic and residual instrumental backgrounds is given by the *Fermi*-LAT science support centre<fnref n="2" />. We used the isotropic spectrum provided for each PSF type of the SOURCE class.

## 23. Gas data

In order to map the neutral atomic hydrogen (Hα) in the region, we used the 21-cm Galactic All Sky Survey GASS III data (Kalberla & Haud 2015). The survey has an angular resolution of 16.2, a velocity resolution of 0.82 km s⁻<fnref n="1" /> in the local standard of rest, and an rms sensitivity of 57 mK. Lines in the analysis region span velocities from -125 to +495 km s⁻<fnref n="1" />.

Following the method developed by Planck Collaboration Int. XXVIII (2015), we decomposed the Hα spectra into sets of individual lines with Pseudo-Voigt profiles. Based on the (l, b, v) distribution of the fitted lines, we identified seven nearby cloud complexes that are coherent in velocity and in position across the region. We named these clouds according to the constellation in their direction, namely Reticulum (Ret), Far Hydrus (HyiF), North Hydrus (HyiN), South Hydrus (HyiS), East Dorado (DorE), West Dorado (DorW), and Horologium (Hor). They are displayed in Fig. 3 and they span velocities from 3.3 to 20.6 km s⁻<fnref n="1" />, -15.7 to -1.6 km s⁻<fnref n="1" />, -0.8 to 13.2 km s⁻<fnref n="1" />, -14.8 to 2.5 km s⁻<fnref n="1" />, 2.5 to 4.9 km s⁻<fnref n="1" />, 5.8 to 33.8 km s⁻<fnref n="1" />, and -14.8 to 5.8 km s⁻<fnref n="1" />, respectively. These clouds are local. We gathered lines from the Magellanic clouds in a single map. We added a ninth component gathering the very weak lines found at velocities below -8 km s⁻<fnref n="1" />.

We produced $N_{HI}$ column-density maps for each complex by integrating the fitted lines that have a central velocity within the chosen velocity interval for each (l, b) direction. This method corrects for potential line spill-over from one velocity interval (i.e. cloud) to the next. The line fits produce small residuals between the modelled and observed Hα spectra. In order to preserve the total Hα intensity, we distributed these residuals among

<footnote n="1">¹ http://galprop.stanford.edu</footnote>

<footnote n="2">² http://fermi.gsfc.nasa.gov/ssc/data/access/lat/BackgroundModels.html</footnote>
<figure><img src="image_5.png" /><figcaption>**Fig. 3.** $N_{HI}$ column densities in the Reticulum (Ret), Far Hydrus (HyiF), North Hydrus (HyiN), South Hydrus (HyiS), East Dorado (DorE), West Dorado (DorW), Horologium (Hor), and the Magellanic clouds for optically thin HI. The faint column densities in the negative-velocity map (NegVel) are multiplied by 4 for display. The maps include a 5°-wide peripheral band around the analysis region (black perimeter) to account for the wings of the *Fermi*-LAT PSF in the γ-ray model.</figcaption></figure>

the fitted lines according to their relative strength in each velocity channel.

In order to correct the HI intensity for potential self-absorption, we produced the $N_{HI}$ maps for a set of uniform spin temperatures, $T_s$, of 100, 125, 150, 200, 300, 400, 500, 600, 700, and 800 K and also for the optically thin case. Nguyen et al. (2018) have shown that a simple isothermal correction of the emission spectra with a uniform spin temperature can well reproduce the more precise $N_{HI}$ column densities that are inferred from the combination of emission and absorption spectra. Given the wide range of potential spin temperatures we could not test individual temperatures for each cloud, but we used the same temperature for all clouds and we let it vary when modelling the dust and γ-ray data.

The nearby clouds in the region are mainly composed of atomic gas. We looked for CO(1–0) emission at 115 GHz towards them, but did not find published observations in that area. We checked that Planck detected no CO(1–0) intensity up to 1 K km s⁻¹ in the region (Planck Collaboration XIII 2014).

## 24. Cloud distances

Leike et al. (2020) used the stellar data from the Gaia-DR2, ALLWISE, PANSTARRS, and 2MASS surveys to model the dust extinction and to reconstruct the 3D map of dust clouds within ~400 pc from the Sun with a spatial resolution of 2 pc. Distances to compact clouds and their dust densities are expected to be well constrained by the precise locations of the stars. In Figs. A.1 and A.3, we display cuts of the 3D dust distribution in 20 pc wide intervals, spanning distances from 150 pc to about 350 pc. The shapes and locations of the clouds isolated in the

HI data in Fig. 3 can be roughly recognised in the dust distribution. We also probed the dust distribution along lines of sight exhibiting large gas column densities in the different clouds (see Fig. A.2). From the dust and gas correspondence we infer average distances of 240 ± 30 pc, 220 ± 30 pc, 220 ± 30 pc, 170 ± 30 pc, 310 ± 30 pc, 340 ± 30 pc, and 280 ± 30 pc to the Reticulum, North Hydrus, South Hydrus, Far Hydrus, East Dorado, West Dorado, and Horologium clouds, respectively. The 20 pc uncertainty reflects the size of dust overdensities that can be isolated in the 3D dust cube.

We applied the same method to revisit the distance to the Eridu cloud. The rim of the Orion-Eridanus superbubble is clearly detected at a distance around 200 pc in the dust 3D data. The Eridu cloud lies beyond this edge. Figure A.3 shows a second dust front at a distance of 300 ± 30 pc with an inclination and extent in the sky that is compatible with that of the elongated Eridu cloud seen in HI. This distance places the cloud at an altitude of 260 pc below the Galactic plane, near the edge of the 220 pc thick HI layer of the local ISM (Nakanishi & Sofue 2016). A much larger distance would be disfavoured. The HI line velocities of the cloud are also typical of the local ISM and not of intermediate- or high-velocity clouds outside the disc. Extending the 3D dust analysis to larger distances than available in the Leike et al. (2020) cube would help confirm the distance to Eridu.

## 3. Models and analyses

Figure 2 shows that the dust optical depth and the γ-ray intensity produced by CRs in the different clouds are strongly correlated. The dust and γ rays both trace, to first order, the total gas column density $N_{HI}$. In order to detect neutral gas not observed in the HI
<figure><img src="image_5.png" /></figure>

**Fig. 3.** $N_{HI}$ column densities in the Reticulum (Ret), Far Hydrus (HyiF), North Hydrus (HyiN), South Hydrus (HyiS), East Dorado (DorE), West Dorado (DorW), Horologium (Hor), and the Magellanic clouds for optically thin HI. The faint column densities in the negative-velocity map (NegVel) are multiplied by 4 for display. The maps include a 5°-wide peripheral band around the analysis region (black perimeter) to account for the wings of the *Fermi*-LAT PSF in the γ-ray model.

the fitted lines according to their relative strength in each velocity channel.

In order to correct the HI intensity for potential self-absorption, we produced the $N_{HI}$ maps for a set of uniform spin temperatures, $T_s$, of 100, 125, 150, 200, 300, 400, 500, 600, 700, and 800 K and also for the optically thin case. Nguyen et al. (2018) have shown that a simple isothermal correction of the emission spectra with a uniform spin temperature can well reproduce the more precise $N_{HI}$ column densities that are inferred from the combination of emission and absorption spectra. Given the wide range of potential spin temperatures we could not test individual temperatures for each cloud, but we used the same temperature for all clouds and we let it vary when modelling the dust and γ-ray data.

The nearby clouds in the region are mainly composed of atomic gas. We looked for CO(1–0) emission at 115 GHz towards them, but did not find published observations in that area. We checked that Planck detected no CO(1–0) intensity up to 1 K km s⁻¹ in the region (Planck Collaboration XIII 2014).

## 24. Cloud distances

Leike et al. (2020) used the stellar data from the Gaia-DR2, ALLWISE, PANSTARRS, and 2MASS surveys to model the dust extinction and to reconstruct the 3D map of dust clouds within ~400 pc from the Sun with a spatial resolution of 2 pc. Distances to compact clouds and their dust densities are expected to be well constrained by the precise locations of the stars. In Figs. A.1 and A.3, we display cuts of the 3D dust distribution in 20 pc wide intervals, spanning distances from 150 pc to about 350 pc. The shapes and locations of the clouds isolated in the

HI data in Fig. 3 can be roughly recognised in the dust distribution. We also probed the dust distribution along lines of sight exhibiting large gas column densities in the different clouds (see Fig. A.2). From the dust and gas correspondence we infer average distances of 240 ± 30 pc, 220 ± 30 pc, 220 ± 30 pc, 170 ± 30 pc, 310 ± 30 pc, 340 ± 30 pc, and 280 ± 30 pc to the Reticulum, North Hydrus, South Hydrus, Far Hydrus, East Dorado, West Dorado, and Horologium clouds, respectively. The 20 pc uncertainty reflects the size of dust overdensities that can be isolated in the 3D dust cube.

We applied the same method to revisit the distance to the Eridu cloud. The rim of the Orion-Eridanus superbubble is clearly detected at a distance around 200 pc in the dust 3D data. The Eridu cloud lies beyond this edge. Figure A.3 shows a second dust front at a distance of 300 ± 30 pc with an inclination and extent in the sky that is compatible with that of the elongated Eridu cloud seen in HI. This distance places the cloud at an altitude of 260 pc below the Galactic plane, near the edge of the 220 pc thick HI layer of the local ISM (Nakanishi & Sofue 2016). A much larger distance would be disfavoured. The HI line velocities of the cloud are also typical of the local ISM and not of intermediate- or high-velocity clouds outside the disc. Extending the 3D dust analysis to larger distances than available in the Leike et al. (2020) cube would help confirm the distance to Eridu.

## 3. Models and analyses

Figure 2 shows that the dust optical depth and the γ-ray intensity produced by CRs in the different clouds are strongly correlated. The dust and γ rays both trace, to first order, the total gas column density $N_{HI}$. In order to detect neutral gas not observed in the HI
data, we used the fact that it is permeated by both CRs and dust grains. We therefore extracted γ-ray and dust residuals between the observation and the model expected from the $N_{HI}$ distribution. We then used the spatial correlation between the positive residuals found in dust and in γ rays to extract the additional gas column densities (see Sect. 3.3). This additional gas corresponds to the DNM at the HI-H$_2$ transition where the dense HI becomes optically thick to 21-cm line emission and CO emission is absent or too faint to trace the H$_2$ gas. The analysis results yield small amounts of dark gas in the Reticulum and North Hydrus clouds (see Sect. 4.2).

### 3.1. Dust model

Provided that the dust-to-gas mass ratio and the mass emission coefficient of the grains are uniform (but not necessarily the same) in each gas phase of a particular cloud, the dust optical depth $\tau_{353}(l, b)$ can be modelled as a linear combination of the different gas column densities ($N_{HI}$ and $N_{H_{DNM}}$), each with free γ normalisation to be fitted to the data, as in Planck Collaboration Int. XXVIII (2015). A free isotropic term, $\gamma_{iso}$, was added to the model to account for the residual noise and the uncertainty in the zero level of the dust maps (Planck Collaboration XI 2014; Casandjian et al. 2022).

The $\tau_{353}(l, b)$ model can be expressed as:

$$ \tau_{353}(l, b) = \sum_{i=1}^{9} \chi_{\gamma_{HI},i} N_{\gamma_{HI},i}(l, b) + \sum_{i=1}^{2} \chi_{\gamma_{DNM},i} N_{\gamma_{DNM},i}(l, b) + \gamma_{iso}, \quad (1) $$

where $N_{\gamma_{HI}}(l, b)$ denotes the nine $N_{\gamma_{HI}}$ maps depicted in Fig. 3 and $N_{\gamma_{DNM}}(l, b)$ stands for the column density in the two DNM components constructed from the γ-ray data (see Sect. 3.3). We degraded the optical-depth map from its original resolution of 5' down to 16.2' to match that of the GASS HI survey.

The $\gamma_{HI,i}$ coefficients give the average values of the dust $\tau_{353}/N_H$ opacity in the different clouds. The $\gamma_{DNM,i}$ parameters can probe opacity changes in the denser DNM phases. They were fitted to the data using a least-square $\chi^2$ minimisation. We expect the uncertainties of the model to exceed those of the dust map because of our assumption of uniform grain distributions through the clouds and because of the limited capability of the HI data to trace the total gas. As we cannot precisely estimate the model uncertainties, we set a fractional error of 11% that yields a reduced $\chi^2$ value of one in the dust fit. This fraction is slightly larger than the 1% to 9% measurement uncertainties in the $\tau_{353}$ values across this region.

### 3.2. Gamma-ray model

Earlier studies have indicated that the Galactic CRs radiating at 0.1–100 GeV have diffusion lengths much larger than typical cloud dimensions and that they propagate through the different gas phases (Grenier et al. 2015). The observed interstellar γ-ray emission can therefore be modelled, to first order, by a linear combination of the gas column densities present in the different clouds seen along the lines of sight, in both the HI and DNM gas phases. The model also includes a contribution from the Galactic IC intensity, $I_{IC}(l, b, E)$, point sources of non-interstellar origin with individual flux spectra $S_j(E)$, and an isotropic intensity, $I_{iso}(E)$, accounting for the extragalactic γ-ray background and for residual CRs misclassified as γ rays. We verified that the soft emission from the Earth limb is not detected in the present energy range for the choice of maximum zenith angles.

The γ-ray intensity, $I(l, b, E)$, in γ cm$^{-2}$ s$^{-1}$ sr$^{-1}$ MeV$^{-1}$ in each $(l, b)$ direction, is modelled in each energy band E as:

$$ \begin{aligned} I(l, b, E) = & q_{\text{LIS}}(E) \times \sum_{i=1}^{9} \chi_{\gamma_{HI},i}^{2}(E) N_{\gamma_{HI},i}(l, b) \times \\ & + q_{\text{LIS}}(E) \times \sum_{i=1}^{2} \chi_{\gamma_{DNM},i}^{2}(E) \tau_{353}^{\text{DNM},i}(l, b) \times \\ & + q_{\text{IC}}(E) I_{\text{IC}}(l, b, E) + q_{\text{iso}}(E) I_{\text{iso}}(E) \times \\ & + \sum_{j} q_{S_j}(E) S_j(E) \delta(l - l_j, b - b_j), \end{aligned} \quad (2) $$

where $N_{\gamma_{HI},i}(l, b)$ denotes the HI column density map of each cloud and the extraction of the dust DNM $\tau_{353}^{\text{DNM},i}$ maps is explained in Sect. 3.3. The $q$ coefficients of the model are determined from fits to the Fermi-LAT data.

The $q_{\gamma_{HI},i}$, $q_{\text{IC}}$, and $q_{\text{iso}}$ parameters are simple normalisation factors to account for possible deviations from the input spectra taken for the local gas γ-ray emissivity $q_{\text{LIS}}(E)$ and for the isotropic and IC intensities. The input $q_{\text{LIS}}$ spectrum was based on the correlation found between the γ radiation and the $N_{\gamma_{HI}}$ column densities derived from the Leiden-Argentine-Bonn (LAB) survey for a spin temperature of 140 K, at latitudes $10^\circ \le l \le 70^\circ$ (Casandjian 2015). The $q$ parameter in the model can therefore compensate for differences in the HI data

(calibration, angular resolution, spin temperature) and for cloud-to-cloud variations in CR flux. Such differences will affect the normalisation equally in all energy bands whereas a change in CR penetration inside a specific cloud will show as an energy-dependent correction, primarily as a loss at energies below one GeV (Schlickeiser et al. 2016).

The source flux spectra $S_j(E)$ were calculated with the spectral characteristics given in the 4FGL-DR3 catalogue. The $q_{S_j}(E)$ normalisation factors in the model allow for possible changes due to the longer exposure used here and to the use of a different interstellar background for source detection in the 4FGL-DR3 catalogue. The $q_{S_j}(E)$ normalisation factors of the sources inside the analysis region or in the peripheral band were fitted individually, except for those inside the SMC and LMC regions and for those found in the lowest dust zones ($\tau_{353} < 3 \times 10^{-6}$), which were treated collectively in three separate sets with a single normalisation factor per set.

In order to compare the model with the LAT photon data in the different energy bands, we multiplied each component map by the energy-dependent exposure and we convolved it with the energy-dependent PSF and energy resolution. The convolution with the response functions was done for each PSF event type independently in small energy bands. We also calculated the isotropic map for each PSF event type in the same small energy bands. We then summed the photon count maps in energy and for the different combinations of PSF types to get the modelled count maps in the eight energy bands of the analysis for the same photon selection as in the data.

Figure 4 gives the photon yields that were obtained in the entire energy band. It shows that the emission originating from the gas can be distinguished from the much smoother emission backgrounds and that the LAT angular resolution is sufficient to separate the various clouds, as well as the different gas phases within the clouds. The distinct spatial distributions of the clouds and of their phases allow rather independent measurements of their individual γ-ray emissivities. On the contrary, we do not expect a firm separation of the isotropic and IC spectra across
data, we used the fact that it is permeated by both CRs and dust grains. We therefore extracted γ-ray and dust residuals between the observation and the model expected from the $N_{HI}$ distribution. We then used the spatial correlation between the positive residuals found in dust and in γ rays to extract the additional gas column densities (see Sect. 3.3). This additional gas corresponds to the DNM at the HI-H$_2$ transition where the dense HI becomes optically thick to 21-cm line emission and CO emission is absent or too faint to trace the H$_2$ gas. The analysis results yield small amounts of dark gas in the Reticulum and North Hydrus clouds (see Sect. 4.2).

### 3.1. Dust model

Provided that the dust-to-gas mass ratio and the mass emission coefficient of the grains are uniform (but not necessarily the same) in each gas phase of a particular cloud, the dust optical depth $\tau_{353}(l, b)$ can be modelled as a linear combination of the different gas column densities ($N_{HI}$ and $N_{H_{DNM}}$), each with free γ normalisation to be fitted to the data, as in Planck Collaboration Int. XXVIII (2015). A free isotropic term, $\gamma_{iso}$, was added to the model to account for the residual noise and the uncertainty in the zero level of the dust maps (Planck Collaboration XI 2014; Casandjian et al. 2022).

The $\tau_{353}(l, b)$ model can be expressed as:

$$ \tau_{353}(l, b) = \sum_{i=1}^{9} \chi_{\gamma_{HI},i} N_{\gamma_{HI},i}(l, b) + \sum_{i=1}^{2} \chi_{\gamma_{DNM},i} N_{\gamma_{DNM},i}(l, b) + \gamma_{iso}, \quad (1) $$

where $N_{\gamma_{HI}}(l, b)$ denotes the nine $N_{\gamma_{HI}}$ maps depicted in Fig. 3 and $N_{\gamma_{DNM}}(l, b)$ stands for the column density in the two DNM components constructed from the γ-ray data (see Sect. 3.3). We degraded the optical-depth map from its original resolution of 5' down to 16.2' to match that of the GASS HI survey.

The $\gamma_{HI,i}$ coefficients give the average values of the dust $\tau_{353}/N_H$ opacity in the different clouds. The $\gamma_{DNM,i}$ parameters can probe opacity changes in the denser DNM phases. They were fitted to the data using a least-square $\chi^2$ minimisation. We expect the uncertainties of the model to exceed those of the dust map because of our assumption of uniform grain distributions through the clouds and because of the limited capability of the HI data to trace the total gas. As we cannot precisely estimate the model uncertainties, we set a fractional error of 11% that yields a reduced $\chi^2$ value of one in the dust fit. This fraction is slightly larger than the 1% to 9% measurement uncertainties in the $\tau_{353}$ values across this region.

### 3.2. Gamma-ray model

Earlier studies have indicated that the Galactic CRs radiating at 0.1–100 GeV have diffusion lengths much larger than typical cloud dimensions and that they propagate through the different gas phases (Grenier et al. 2015). The observed interstellar γ-ray emission can therefore be modelled, to first order, by a linear combination of the gas column densities present in the different clouds seen along the lines of sight, in both the HI and DNM gas phases. The model also includes a contribution from the Galactic IC intensity, $I_{IC}(l, b, E)$, point sources of non-interstellar origin with individual flux spectra $S_j(E)$, and an isotropic intensity, $I_{iso}(E)$, accounting for the extragalactic γ-ray background and for residual CRs misclassified as γ rays. We verified that the soft emission from the Earth limb is not detected in the present energy range for the choice of maximum zenith angles.

The γ-ray intensity, $I(l, b, E)$, in γ cm$^{-2}$ s$^{-1}$ sr$^{-1}$ MeV$^{-1}$ in each $(l, b)$ direction, is modelled in each energy band E as:

$$ \begin{aligned} I(l, b, E) = & q_{\text{LIS}}(E) \times \sum_{i=1}^{9} \chi_{\gamma_{HI},i}^{2}(E) N_{\gamma_{HI},i}(l, b) \times \\ & + q_{\text{LIS}}(E) \times \sum_{i=1}^{2} \chi_{\gamma_{DNM},i}^{2}(E) \tau_{353}^{\text{DNM},i}(l, b) \times \\ & + q_{\text{IC}}(E) I_{\text{IC}}(l, b, E) + q_{\text{iso}}(E) I_{\text{iso}}(E) \times \\ & + \sum_{j} q_{S_j}(E) S_j(E) \delta(l - l_j, b - b_j), \end{aligned} \quad (2) $$

where $N_{\gamma_{HI},i}(l, b)$ denotes the HI column density map of each cloud and the extraction of the dust DNM $\tau_{353}^{\text{DNM},i}$ maps is explained in Sect. 3.3. The $q$ coefficients of the model are determined from fits to the Fermi-LAT data.

The $q_{\gamma_{HI},i}$, $q_{\text{IC}}$, and $q_{\text{iso}}$ parameters are simple normalisation factors to account for possible deviations from the input spectra taken for the local gas γ-ray emissivity $q_{\text{LIS}}(E)$ and for the isotropic and IC intensities. The input $q_{\text{LIS}}$ spectrum was based on the correlation found between the γ radiation and the $N_{\gamma_{HI}}$ column densities derived from the Leiden-Argentine-Bonn (LAB) survey for a spin temperature of 140 K, at latitudes $10^\circ \le l \le 70^\circ$ (Casandjian 2015). The $q$ parameter in the model can therefore compensate for differences in the HI data

(calibration, angular resolution, spin temperature) and for cloud-to-cloud variations in CR flux. Such differences will affect the normalisation equally in all energy bands whereas a change in CR penetration inside a specific cloud will show as an energy-dependent correction, primarily as a loss at energies below one GeV (Schlickeiser et al. 2016).

The source flux spectra $S_j(E)$ were calculated with the spectral characteristics given in the 4FGL-DR3 catalogue. The $q_{S_j}(E)$ normalisation factors in the model allow for possible changes due to the longer exposure used here and to the use of a different interstellar background for source detection in the 4FGL-DR3 catalogue. The $q_{S_j}(E)$ normalisation factors of the sources inside the analysis region or in the peripheral band were fitted individually, except for those inside the SMC and LMC regions and for those found in the lowest dust zones ($\tau_{353} < 3 \times 10^{-6}$), which were treated collectively in three separate sets with a single normalisation factor per set.

In order to compare the model with the LAT photon data in the different energy bands, we multiplied each component map by the energy-dependent exposure and we convolved it with the energy-dependent PSF and energy resolution. The convolution with the response functions was done for each PSF event type independently in small energy bands. We also calculated the isotropic map for each PSF event type in the same small energy bands. We then summed the photon count maps in energy and for the different combinations of PSF types to get the modelled count maps in the eight energy bands of the analysis for the same photon selection as in the data.

Figure 4 gives the photon yields that were obtained in the entire energy band. It shows that the emission originating from the gas can be distinguished from the much smoother emission backgrounds and that the LAT angular resolution is sufficient to separate the various clouds, as well as the different gas phases within the clouds. The distinct spatial distributions of the clouds and of their phases allow rather independent measurements of their individual γ-ray emissivities. On the contrary, we do not expect a firm separation of the isotropic and IC spectra across
<figure><img src="image_6.png" /><figcaption>**Fig. 4.** Photon yields on a 0.25 pixel grid from the main components of the γ-ray model in the 0.16–63 GeV band. The yields come from the Hi and DNM column densities derived for the optically thin case in the different clouds, and from the isotropic background, Galactic IC emission, and γ-ray point sources.</figcaption></figure>

this region. We used a binned maximum-likelihood with Poisson statistics to fit the model coefficients $q$ to the LAT data in each of the eight energy bands and in the overall one.

### 3.3. Analysis iterations

In order to extract the DNM gas present in both the dust and γ-ray data, we iteratively coupled the γ-ray and dust models. We started by fitting the dust optical depth and the γ-ray data independently according to Eqs. (1) and (2) with only the Hi maps as gaseous components and with ancillary (other than gas) components. We then built maps of positive residuals (see how below) between the data (dust or γ rays) and the best-fit contributions of the model components. We separated these excesses into two DNM maps according to their location, one in the Reticulum cloud and the other one in the North Hydrus area.

The next loops in the iterative process provided the DNM templates estimated from the dust to the γ-ray model ($\tau_{353}^{\text{DNM}}$ in Eq. (2)) and conversely provided the DNM column density maps derived from the γ rays to the dust model ($N_{\text{HI}}^{\text{DNM}}(l, b)$ in Eq. (1)). We repeated the iteration between the γ-ray and dust models until the log-likelihood of the fit to γ-ray data saturated.

The estimates of the $q$ and γ model coefficients and of the DNM templates improve at each iteration since there is less and less need for the other components, in particular the Hi ones, to compensate for the missing gas. They can still miss gas at some level because the DNM maps provided by the γ rays or dust emission have limitations (e.g. dust emissivity variations, limited γ-ray sensitivity).

Care must be taken in the extraction of the positive residuals because of the noise around zero. A simple cut at zero would induce a positive offset bias, and so we denoised the residual maps using the multi-resolution support method implemented in the MR filter software (Starck & Pierre 1998). For the dust residuals, we used a multi-resolution thresholding filter, with six scales in the B-spline-wavelet transform (à trous algorithm), Gaussian noise and a 2-σ threshold to filter. For the γ-ray residuals, we used the photon residual map obtained in the 0.16–63 GeV band and we filtered it with six wavelet scales, Poisson noise, and a 2-σ threshold. The stability of the iterative analysis and the results of the dust and γ-ray fits are discussed in the following section.

## 4. Results

Before presenting the γ-ray emissivity spectra obtained in the Reticulum cloud in Sect. 4.4, we discuss how the fits to the γ-ray and dust data improve with optical depth corrections to the $N_{\text{HI}}$ column densities (Sect. 4.1) and with the addition of DNM gas templates (Sect. 4.2). We also assess the robustness of the best-fit models with jackknife tests in Sect. 4.3.

### 4.1. Hi optical depth correction

We do not know the exact level of Hi opacity in the different clouds. This opacity induces systematic uncertainties in the $N_{\text{HI}}$ column densities, and therefore also in the Hi contributions to the γ-ray and dust models. Nguyen et al. (2019) showed that a
<figure><img src="image_6.png" /><figcaption>**Fig. 4.** Photon yields on a 0.25 pixel grid from the main components of the γ-ray model in the 0.16–63 GeV band. The yields come from the Hi and DNM column densities derived for the optically thin case in the different clouds, and from the isotropic background, Galactic IC emission, and γ-ray point sources.</figcaption></figure>

this region. We used a binned maximum-likelihood with Poisson statistics to fit the model coefficients $q$ to the LAT data in each of the eight energy bands and in the overall one.

### 3.3. Analysis iterations

In order to extract the DNM gas present in both the dust and γ-ray data, we iteratively coupled the γ-ray and dust models. We started by fitting the dust optical depth and the γ-ray data independently according to Eqs. (1) and (2) with only the Hi maps as gaseous components and with ancillary (other than gas) components. We then built maps of positive residuals (see how below) between the data (dust or γ rays) and the best-fit contributions of the model components. We separated these excesses into two DNM maps according to their location, one in the Reticulum cloud and the other one in the North Hydrus area.

The next loops in the iterative process provided the DNM templates estimated from the dust to the γ-ray model ($\tau_{353}^{\text{DNM}}$ in Eq. (2)) and conversely provided the DNM column density maps derived from the γ rays to the dust model ($N_{\text{HI}}^{\text{DNM}}(l, b)$ in Eq. (1)). We repeated the iteration between the γ-ray and dust models until the log-likelihood of the fit to γ-ray data saturated.

The estimates of the $q$ and γ model coefficients and of the DNM templates improve at each iteration since there is less and less need for the other components, in particular the Hi ones, to compensate for the missing gas. They can still miss gas at some level because the DNM maps provided by the γ rays or dust emission have limitations (e.g. dust emissivity variations, limited γ-ray sensitivity).

Care must be taken in the extraction of the positive residuals because of the noise around zero. A simple cut at zero would induce a positive offset bias, and so we denoised the residual maps using the multi-resolution support method implemented in the MR filter software (Starck & Pierre 1998). For the dust residuals, we used a multi-resolution thresholding filter, with six scales in the B-spline-wavelet transform (à trous algorithm), Gaussian noise and a 2-σ threshold to filter. For the γ-ray residuals, we used the photon residual map obtained in the 0.16–63 GeV band and we filtered it with six wavelet scales, Poisson noise, and a 2-σ threshold. The stability of the iterative analysis and the results of the dust and γ-ray fits are discussed in the following section.

## 4. Results

Before presenting the γ-ray emissivity spectra obtained in the Reticulum cloud in Sect. 4.4, we discuss how the fits to the γ-ray and dust data improve with optical depth corrections to the $N_{\text{HI}}$ column densities (Sect. 4.1) and with the addition of DNM gas templates (Sect. 4.2). We also assess the robustness of the best-fit models with jackknife tests in Sect. 4.3.

### 4.1. Hi optical depth correction

We do not know the exact level of Hi opacity in the different clouds. This opacity induces systematic uncertainties in the $N_{\text{HI}}$ column densities, and therefore also in the Hi contributions to the γ-ray and dust models. Nguyen et al. (2019) showed that a
<figure><img src="image_7.png" /><figcaption>**Fig. 5.** Evolution of the log-likelihood ratio of the γ-ray fit in the 0.16–63 GeV band as a function of the HI spin temperature used to calculate the $N_{HI}$ column densities in the clouds.</figcaption></figure>

uniform isothermal correction of the HI emission spectra across a cloud provides a good approximation to more precise estimates based on emission and absorption HI spectra. We therefore produced the dust and γ-ray models for a wide range of optical depth corrections, from a rather optically thick case ($T_s = 100$ K) to the optically thin case that yields the minimum amount of gas in the cloud. Given the moderate $N_{HI}$ column densities of the present clouds, the difference between those extreme corrections reaches 25% on a few pixels and less than 1% on the mean column density of a cloud.

Even though the spin temperature is uniform, the HI opacity correction depends on the HI intensity recorded in each direction and velocity channel, so it varies across the cloud and slightly modifies the $N_{HI}$ spatial distribution of a cloud. Figure 5 shows that the goodness of the γ-ray fit (after the DNM iteration) significantly improves as we increase the HI spin temperature. The $N_{HI}$ maps best match the spatial structure of the γ-ray flux emerging from the atomic gas for $T_s \ge 300$ K; therefore for small HI opacities. Given the saturation of the log-likelihood ratio towards the optically thin case, we present the analysis results obtained for this thin case in the rest of the paper, unless otherwise mentioned. The residual maps between the γ-ray data and best-fit model show that using the same $T_s$ for all the clouds can reproduce the γ-ray data without any significantly deviant cloud.

As the spin temperature $T_s$ increases, the optical depth correction and $N_{HI}$ column densities decrease, and the γ-ray emissivity $q_{HI}$ and dust opacity $\gamma_{HI}$ of the HI clouds increase. They do so respectively by 4–11% and 3–7% as the spin temperature rises over the entire range. HI optical depth corrections have therefore a small effect in this region because the clouds are not massive, so they are fairly transparent to HI line radiation.

## 4.2. DNM detection

In order to check for the presence of DNM gas not traced by HI line intensities, we first performed the dust and γ-ray fits with only the $N_{HI}$ maps as gaseous components. The residuals are displayed in the upper panels of Fig. 6. Positive residuals indicate that the data exceed the model predictions. Coincident excesses are seen in dust and in γ rays even though the two fits are independent. They cross the region in diagonal, from the North Hydrus cloud (upper left corner) down along the Reticulum filament. A coincident increase in CR density and in dust opacity is

very improbable. Hadronic CR interactions with dust grains produce undetectable γ-ray emission (Grenier et al. 2005), so this correlated excess is likely due to additional gas. The distribution of positive γ-ray excesses along the Reticulum cloud is sparse, but the large extent of negative (blue) excesses across the whole analysis region indicates that the γ-ray fit has unduly pushed the HI emissivities of the clouds upward to compensate for some missing gas.

The DNM gas column density is built by correlating and filtering the dust and γ-ray excess signals in the iteration between both model fits (see Sect. 3.3). We built separate DNM maps for the North Hydrus and Reticulum clouds. Adding them improves the fits quality (see the middle maps of Fig. 6). A log-likelihood ratio increase by 161 between the γ-ray fits with and without the DNM templates implies a very significant DNM detection. The spectrum of the γ-ray emission associated with both DNM components is consistent with that found in the HI gas phase over the whole energy range (see below). This gives further support to a gaseous origin of the γ-ray emission in the DNM structures.

Positive residuals remain in the final dust fit, even after the addition of the DNM maps. Those in spatial coincidence with the DNM phase are partially due to the non-linear increase in dust opacity as the gas density rises towards the dense molecular phase (Remy et al. 2017). They also result from the limited γ-ray photon statistics that limit the accuracy of the DNM template provided to the dust fit. The negative and positive dust residuals seen at longitudes $l < 270^\circ$ reflect limitations in the $T_{353}$ data map due to the difficult separation of the dust emission from other submillimetre emission in these very faint zones (Planck Collaboration XI 2014; Irfan et al. 2019, and references therein). These regions were excluded from the DNM analysis.

## 4.3. Best fits and jackknife tests

The best dust and γ-ray fits with the models described by Eqs. (1) and (2) include the nine HI and two DNM maps for the gaseous components. Tables 1 and 2 give the best-fit coefficients obtained for the HI and DNM phases of the main Reticulum and North Hydrus clouds. The results found for the other small local clouds will be discussed in a subsequent paper (their emissivity spectra are rather close to the $q_{LIS}$ expectation). We focus here on the Reticulum filament to compare it with the Eridu one and we include the emissivity measurement in the North Hydrus cloud as it partially overlaps the northern end of Reticulum. The γ-ray emissivities are given in each of the eight energy bands and for the whole 0.16–63 GeV band.

The final iteration of the dust and γ-ray fits leads to the residual maps presented in the middle panels of Fig. 6. The γ-ray residuals are small and randomly distributed, so our linear model provides an excellent statistical fit to the γ-ray data in the 0.16–63 GeV energy band across the whole region. The same conclusion is reached in the eight separate energy bands that are not shown here.

The statistical errors on the best-fit coefficients are inferred from the information matrix of the fits (Strong 1985). They include the correlation between parameters. They are generally small (2–9% for the HI gas γ-ray emissivities, 15% and 10% for DNM emissivities for Ret and HyiN, 13% for $q_{IC}$, 9% for $q_{ISO}$, 0.4–1.3% for the HI dust opacities, 8% and 3% for $\gamma_{DNM,Ret}$ and $\gamma_{DNM,HyiN}$ and 2% for $\gamma_{ISO}$) because the clouds are spatially well resolved in the analysis and because tight correlations exist between the gas, dust and γ-ray maps (see Fig. 2).

We also checked the magnitude of potential systematic uncertainties due to the linear-combination approximation
<figure><img src="image_7.png" /><figcaption>**Fig. 5.** Evolution of the log-likelihood ratio of the γ-ray fit in the 0.16–63 GeV band as a function of the HI spin temperature used to calculate the $N_{HI}$ column densities in the clouds.</figcaption></figure>

uniform isothermal correction of the HI emission spectra across a cloud provides a good approximation to more precise estimates based on emission and absorption HI spectra. We therefore produced the dust and γ-ray models for a wide range of optical depth corrections, from a rather optically thick case ($T_s = 100$ K) to the optically thin case that yields the minimum amount of gas in the cloud. Given the moderate $N_{HI}$ column densities of the present clouds, the difference between those extreme corrections reaches 25% on a few pixels and less than 1% on the mean column density of a cloud.

Even though the spin temperature is uniform, the HI opacity correction depends on the HI intensity recorded in each direction and velocity channel, so it varies across the cloud and slightly modifies the $N_{HI}$ spatial distribution of a cloud. Figure 5 shows that the goodness of the γ-ray fit (after the DNM iteration) significantly improves as we increase the HI spin temperature. The $N_{HI}$ maps best match the spatial structure of the γ-ray flux emerging from the atomic gas for $T_s \ge 300$ K; therefore for small HI opacities. Given the saturation of the log-likelihood ratio towards the optically thin case, we present the analysis results obtained for this thin case in the rest of the paper, unless otherwise mentioned. The residual maps between the γ-ray data and best-fit model show that using the same $T_s$ for all the clouds can reproduce the γ-ray data without any significantly deviant cloud.

As the spin temperature $T_s$ increases, the optical depth correction and $N_{HI}$ column densities decrease, and the γ-ray emissivity $q_{HI}$ and dust opacity $\gamma_{HI}$ of the HI clouds increase. They do so respectively by 4–11% and 3–7% as the spin temperature rises over the entire range. HI optical depth corrections have therefore a small effect in this region because the clouds are not massive, so they are fairly transparent to HI line radiation.

## 4.2. DNM detection

In order to check for the presence of DNM gas not traced by HI line intensities, we first performed the dust and γ-ray fits with only the $N_{HI}$ maps as gaseous components. The residuals are displayed in the upper panels of Fig. 6. Positive residuals indicate that the data exceed the model predictions. Coincident excesses are seen in dust and in γ rays even though the two fits are independent. They cross the region in diagonal, from the North Hydrus cloud (upper left corner) down along the Reticulum filament. A coincident increase in CR density and in dust opacity is

very improbable. Hadronic CR interactions with dust grains produce undetectable γ-ray emission (Grenier et al. 2005), so this correlated excess is likely due to additional gas. The distribution of positive γ-ray excesses along the Reticulum cloud is sparse, but the large extent of negative (blue) excesses across the whole analysis region indicates that the γ-ray fit has unduly pushed the HI emissivities of the clouds upward to compensate for some missing gas.

The DNM gas column density is built by correlating and filtering the dust and γ-ray excess signals in the iteration between both model fits (see Sect. 3.3). We built separate DNM maps for the North Hydrus and Reticulum clouds. Adding them improves the fits quality (see the middle maps of Fig. 6). A log-likelihood ratio increase by 161 between the γ-ray fits with and without the DNM templates implies a very significant DNM detection. The spectrum of the γ-ray emission associated with both DNM components is consistent with that found in the HI gas phase over the whole energy range (see below). This gives further support to a gaseous origin of the γ-ray emission in the DNM structures.

Positive residuals remain in the final dust fit, even after the addition of the DNM maps. Those in spatial coincidence with the DNM phase are partially due to the non-linear increase in dust opacity as the gas density rises towards the dense molecular phase (Remy et al. 2017). They also result from the limited γ-ray photon statistics that limit the accuracy of the DNM template provided to the dust fit. The negative and positive dust residuals seen at longitudes $l < 270^\circ$ reflect limitations in the $T_{353}$ data map due to the difficult separation of the dust emission from other submillimetre emission in these very faint zones (Planck Collaboration XI 2014; Irfan et al. 2019, and references therein). These regions were excluded from the DNM analysis.

## 4.3. Best fits and jackknife tests

The best dust and γ-ray fits with the models described by Eqs. (1) and (2) include the nine HI and two DNM maps for the gaseous components. Tables 1 and 2 give the best-fit coefficients obtained for the HI and DNM phases of the main Reticulum and North Hydrus clouds. The results found for the other small local clouds will be discussed in a subsequent paper (their emissivity spectra are rather close to the $q_{LIS}$ expectation). We focus here on the Reticulum filament to compare it with the Eridu one and we include the emissivity measurement in the North Hydrus cloud as it partially overlaps the northern end of Reticulum. The γ-ray emissivities are given in each of the eight energy bands and for the whole 0.16–63 GeV band.

The final iteration of the dust and γ-ray fits leads to the residual maps presented in the middle panels of Fig. 6. The γ-ray residuals are small and randomly distributed, so our linear model provides an excellent statistical fit to the γ-ray data in the 0.16–63 GeV energy band across the whole region. The same conclusion is reached in the eight separate energy bands that are not shown here.

The statistical errors on the best-fit coefficients are inferred from the information matrix of the fits (Strong 1985). They include the correlation between parameters. They are generally small (2–9% for the HI gas γ-ray emissivities, 15% and 10% for DNM emissivities for Ret and HyiN, 13% for $q_{IC}$, 9% for $q_{ISO}$, 0.4–1.3% for the HI dust opacities, 8% and 3% for $\gamma_{DNM,Ret}$ and $\gamma_{DNM,HyiN}$ and 2% for $\gamma_{ISO}$) because the clouds are spatially well resolved in the analysis and because tight correlations exist between the gas, dust and γ-ray maps (see Fig. 2).

We also checked the magnitude of potential systematic uncertainties due to the linear-combination approximation
<figure><img src="image_8.png" /><figcaption>**Fig. 6.** Residual maps (data minus best-fit model, in sigma units) in γ-ray counts in the 0.16–63 GeV band (σ<sub>γ</sub>) and in dust optical depth (σ<sub>T353</sub>), obtained when including different sets of gaseous components in the models: without (top) and with (middle) the DNM components in addition to the HI templates. The N<sub>H_DNM</sub> column-density maps (bottom) extracted from the γ-ray residuals (left) and from dust residuals (right) are shown in the bottom row. The models are derived for an HI at the optically thin case.</figcaption></figure>

<table><caption>**Table 1.** Best-fit dust opacities in the HI and DNM phases of the Reticulum and North Hydrus clouds, and isotropic optical depth.</caption><thead><tr><th>γ<sub>HI,Ret</sub> (10<sup>-27</sup> cm<sup>2</sup>)</th><th>γ<sub>DNM,Ret</sub> (10<sup>-27</sup> cm<sup>2</sup>)</th><th>γ<sub>HI,HyiN</sub> (10<sup>-27</sup> cm<sup>2</sup>)</th><th>γ<sub>DNM,HyiN</sub> (10<sup>-27</sup> cm<sup>2</sup>)</th><th>γ<sub>iso</sub></th></tr></thead><tbody><tr><td>8.91 ± 0.07 ± 0.12</td><td>5.05 ± 0.40 ± 0.58</td><td>8.68 ± 0.03 ± 0.05</td><td>14.90 ± 0.42 ± 0.44</td><td>(2.79 ± 0.06 ± 0.12) 10<sup>-7</sup></td></tr></tbody></table>

**Notes.** The first uncertainties are statistical (1σ) and the second give the standard deviations obtained in the jackknife fits.

<table><caption>**Table 2.** Best-fit γ-ray emissivities spectra per gas nucleon in the HI and DNM phases of the Reticulum and North Hydrus clouds.</caption><thead><tr><th>Energy (MeV)</th><th>q<sub>HI,Ret</sub></th><th>q<sub>DNM,Ret</sub> (10<sup>26</sup> cm<sup>-2</sup>)</th><th>q<sub>HI,HyiN</sub></th><th>q<sub>DNM,HyiN</sub> (10<sup>26</sup> cm<sup>-2</sup>)</th></tr></thead><tbody><tr><td>10<sup>2.2</sup>–10<sup>2.4</sup></td><td>1.56 ± 0.23 ± 0.09</td><td>0.21 ± 1.22 ± 0.28</td><td>1.17 ± 0.15 ± 0.04</td><td>0.78 ± 0.52 ± 0.18</td></tr><tr><td>10<sup>2.4</sup>–10<sup>2.6</sup></td><td>1.27 ± 0.13 ± 0.07</td><td>0.74 ± 0.61 ± 0.26</td><td>1.10 ± 0.08 ± 0.03</td><td>0.61 ± 0.26 ± 0.09</td></tr><tr><td>10<sup>2.6</sup>–10<sup>2.8</sup></td><td>0.98 ± 0.10 ± 0.06</td><td>0.85 ± 0.41 ± 0.20</td><td>1.02 ± 0.06 ± 0.02</td><td>0.65 ± 0.17 ± 0.06</td></tr><tr><td>10<sup>2.8</sup>–10<sup>3.0</sup></td><td>1.16 ± 0.08 ± 0.04</td><td>0.85 ± 0.29 ± 0.12</td><td>1.04 ± 0.04 ± 0.02</td><td>0.77 ± 0.12 ± 0.04</td></tr><tr><td>10<sup>3.0</sup>–10<sup>3.2</sup></td><td>0.98 ± 0.07 ± 0.04</td><td>0.85 ± 0.27 ± 0.14</td><td>1.03 ± 0.04 ± 0.01</td><td>0.64 ± 0.11 ± 0.05</td></tr><tr><td>10<sup>3.2</sup>–10<sup>3.6</sup></td><td>0.98 ± 0.07 ± 0.04</td><td>1.46 ± 0.26 ± 0.11</td><td>1.00 ± 0.04 ± 0.02</td><td>0.70 ± 0.11 ± 0.04</td></tr><tr><td>10<sup>3.6</sup>–10<sup>4.0</sup></td><td>1.02 ± 0.16 ± 0.09</td><td>1.20 ± 0.55 ± 0.32</td><td>1.12 ± 0.09 ± 0.03</td><td>0.98 ± 0.23 ± 0.09</td></tr><tr><td>10<sup>4.0</sup>–10<sup>4.8</sup></td><td>0.93 ± 0.38 ± 0.18</td><td>0.05 ± 1.14 ± 0.25</td><td>1.13 ± 0.21 ± 0.06</td><td>1.18 ± 0.55 ± 0.17</td></tr><tr><td>10<sup>2.2</sup>–10<sup>4.8</sup></td><td>1.11 ± 0.04 ± 0.03</td><td>0.99 ± 0.15 ± 0.09</td><td>1.05 ± 0.02 ± 0.01</td><td>0.75 ± 0.07 ± 0.03</td></tr></tbody></table>

**Notes.** The first uncertainties are statistical (1σ) and the second give the standard deviations obtained in the jackknife fits.

underlying the models. In order to evaluate the impact of possible spatial changes in the model parameters and/or in the mean level of HI self-absorption, we repeated the last iteration of the dust and γ-ray fits a thousand times over random subsets of the analysis region, each time masking out 20% of the pixels with 1°-wide randomly selected squares. In γ rays, the jackknife tests were performed for the total 0.16–63 GeV energy band and for each of

the eight bands. Figure 7 shows the distribution obtained for the best-fit coefficients relating to the Reticulum and North Hydrus clouds and to the large-scale IC and isotropic components. The Gaussian-like distributions indicate that the emissivities and opacities presented in Tables 1 and 2 are statistically stable and that they are not driven by subset regions in each cloud complex. The standard deviations found in the jackknife distributions
<figure><img src="image_8.png" /><figcaption>**Fig. 6.** Residual maps (data minus best-fit model, in sigma units) in γ-ray counts in the 0.16–63 GeV band (σ<sub>γ</sub>) and in dust optical depth (σ<sub>T353</sub>), obtained when including different sets of gaseous components in the models: without (top) and with (middle) the DNM components in addition to the HI templates. The N<sub>H_DNM</sub> column-density maps (bottom) extracted from the γ-ray residuals (left) and from dust residuals (right) are shown in the bottom row. The models are derived for an HI at the optically thin case.</figcaption></figure>

<table><caption>**Table 1.** Best-fit dust opacities in the HI and DNM phases of the Reticulum and North Hydrus clouds, and isotropic optical depth.</caption><thead><tr><th>γ<sub>HI,Ret</sub> (10<sup>-27</sup> cm<sup>2</sup>)</th><th>γ<sub>DNM,Ret</sub> (10<sup>-27</sup> cm<sup>2</sup>)</th><th>γ<sub>HI,HyiN</sub> (10<sup>-27</sup> cm<sup>2</sup>)</th><th>γ<sub>DNM,HyiN</sub> (10<sup>-27</sup> cm<sup>2</sup>)</th><th>γ<sub>iso</sub></th></tr></thead><tbody><tr><td>8.91 ± 0.07 ± 0.12</td><td>5.05 ± 0.40 ± 0.58</td><td>8.68 ± 0.03 ± 0.05</td><td>14.90 ± 0.42 ± 0.44</td><td>(2.79 ± 0.06 ± 0.12) 10<sup>-7</sup></td></tr></tbody></table>

**Notes.** The first uncertainties are statistical (1σ) and the second give the standard deviations obtained in the jackknife fits.

<table><caption>**Table 2.** Best-fit γ-ray emissivities spectra per gas nucleon in the HI and DNM phases of the Reticulum and North Hydrus clouds.</caption><thead><tr><th>Energy (MeV)</th><th>q<sub>HI,Ret</sub></th><th>q<sub>DNM,Ret</sub> (10<sup>26</sup> cm<sup>-2</sup>)</th><th>q<sub>HI,HyiN</sub></th><th>q<sub>DNM,HyiN</sub> (10<sup>26</sup> cm<sup>-2</sup>)</th></tr></thead><tbody><tr><td>10<sup>2.2</sup>–10<sup>2.4</sup></td><td>1.56 ± 0.23 ± 0.09</td><td>0.21 ± 1.22 ± 0.28</td><td>1.17 ± 0.15 ± 0.04</td><td>0.78 ± 0.52 ± 0.18</td></tr><tr><td>10<sup>2.4</sup>–10<sup>2.6</sup></td><td>1.27 ± 0.13 ± 0.07</td><td>0.74 ± 0.61 ± 0.26</td><td>1.10 ± 0.08 ± 0.03</td><td>0.61 ± 0.26 ± 0.09</td></tr><tr><td>10<sup>2.6</sup>–10<sup>2.8</sup></td><td>0.98 ± 0.10 ± 0.06</td><td>0.85 ± 0.41 ± 0.20</td><td>1.02 ± 0.06 ± 0.02</td><td>0.65 ± 0.17 ± 0.06</td></tr><tr><td>10<sup>2.8</sup>–10<sup>3.0</sup></td><td>1.16 ± 0.08 ± 0.04</td><td>0.85 ± 0.29 ± 0.12</td><td>1.04 ± 0.04 ± 0.02</td><td>0.77 ± 0.12 ± 0.04</td></tr><tr><td>10<sup>3.0</sup>–10<sup>3.2</sup></td><td>0.98 ± 0.07 ± 0.04</td><td>0.85 ± 0.27 ± 0.14</td><td>1.03 ± 0.04 ± 0.01</td><td>0.64 ± 0.11 ± 0.05</td></tr><tr><td>10<sup>3.2</sup>–10<sup>3.6</sup></td><td>0.98 ± 0.07 ± 0.04</td><td>1.46 ± 0.26 ± 0.11</td><td>1.00 ± 0.04 ± 0.02</td><td>0.70 ± 0.11 ± 0.04</td></tr><tr><td>10<sup>3.6</sup>–10<sup>4.0</sup></td><td>1.02 ± 0.16 ± 0.09</td><td>1.20 ± 0.55 ± 0.32</td><td>1.12 ± 0.09 ± 0.03</td><td>0.98 ± 0.23 ± 0.09</td></tr><tr><td>10<sup>4.0</sup>–10<sup>4.8</sup></td><td>0.93 ± 0.38 ± 0.18</td><td>0.05 ± 1.14 ± 0.25</td><td>1.13 ± 0.21 ± 0.06</td><td>1.18 ± 0.55 ± 0.17</td></tr><tr><td>10<sup>2.2</sup>–10<sup>4.8</sup></td><td>1.11 ± 0.04 ± 0.03</td><td>0.99 ± 0.15 ± 0.09</td><td>1.05 ± 0.02 ± 0.01</td><td>0.75 ± 0.07 ± 0.03</td></tr></tbody></table>

**Notes.** The first uncertainties are statistical (1σ) and the second give the standard deviations obtained in the jackknife fits.

underlying the models. In order to evaluate the impact of possible spatial changes in the model parameters and/or in the mean level of HI self-absorption, we repeated the last iteration of the dust and γ-ray fits a thousand times over random subsets of the analysis region, each time masking out 20% of the pixels with 1°-wide randomly selected squares. In γ rays, the jackknife tests were performed for the total 0.16–63 GeV energy band and for each of

the eight bands. Figure 7 shows the distribution obtained for the best-fit coefficients relating to the Reticulum and North Hydrus clouds and to the large-scale IC and isotropic components. The Gaussian-like distributions indicate that the emissivities and opacities presented in Tables 1 and 2 are statistically stable and that they are not driven by subset regions in each cloud complex. The standard deviations found in the jackknife distributions
<figure><img src="image_9.png" /><figcaption>**Fig. 7.** Distributions of the best-fit coefficients obtained in the thousand jackknife fits to the 0.16–63 GeV γ-ray data (two left columns) and dust data (two right columns) for the optically thin HI case. The vertical lines and overlaid numbers give the mean values of the distributions. The units are $10^{26}$ cm$^{-2}$ for $q_{DNM}$, $10^{-27}$ cm$^2$ for $\gamma_{HI}$ and $\gamma_{DNM}$, and $10^{-7}$ for $\gamma_{ISO}$. $q_{HI,i}$, $q_{IC}$ and $q_{ISO}$ are dimensionless normalisation factors.</figcaption></figure>

amount to 1–6% for the HI gas γ-ray emissivities, 9% and 4% for the Ret and HyiN DNM emissivities, 8% for $q_{IC}$, 5% for $q_{ISO}$, 0.6–3% for the HI dust opacities, 12% and 3% for $\gamma_{DNM,Ret}$ and $\gamma_{DNM,HyiN}$ and 4% for $\gamma_{ISO}$. Both statistical errors ($1\sigma$) and systematic uncertainties (jackknife standard deviations) are given in Tables 1 and 2. We quadratically added them to construct the final uncertainties on the $q$ and $\gamma$ coefficients.

Figure 8 shows the stability with γ-ray energy of the best-fit coefficients relating to the Reticulum and North Hydrus clouds and to the IC and isotropic intensities. We find no significant deviation with respect to the average coefficient fitted for the entire energy band (dashed lines). This implies that the HI and DNM gas emissivity spectra in the two clouds have the same shape as the average $q_{LIS}(E)$ spectrum measured in the local ISM. We further discuss their relative normalisation with respect to $q_{LIS}$ in the next section.

Figure 8 and the value of the relative intensity ratio $q_{IC} = 1.18 \pm 0.15$ in the whole energy band further show that the IC spectrum detected in these directions of the sky is fully consistent in shape and in intensity with the GALPROP prediction based on the complex evolution of the interstellar radiation field and of the CR lepton spectra across the Milky Way. The isotropic spectrum found in this small region of the sky is also close to the larger-scale estimate obtained between 10° and 60° in latitude<fnref n="3" />. The relative intensity ratio is $q_{ISO} = 0.84 \pm 0.08$ in the 0.16–63 GeV energy band. The interstellar foregrounds are faint enough to separate the IC and isotropic components despite rather similar spatial distributions across the analysis region (see Fig. 3).

Planck Collaboration XI (2014) and Remy et al. (2017) found marked gradients in dust opacity across the sky. These variations reflect a systematic opacity increase as the gas becomes denser and the grains evolve in size and chemistry (Remy et al. 2017; Ysard et al. 2019, and references therein). Diffuse atomic cirri at high latitudes exhibit low opacities around a mean value of $(7.0 \pm 2.0) \times 10^{-27}$ cm$^2$ (Planck Collaboration XI 2014) or $(7.1 \pm 0.6) \times 10^{-27}$ cm$^2$ (Planck Collaboration Int. XVII 2014) whereas larger values ranging from $10 \times 10^{-27}$ to $16 \times 10^{-27}$ cm$^2$ are found in the denser atomic envelopes of nearby clouds (Remy et al. 2017), and values from $15 \times 10^{-27}$ to $17 \times 10^{-27}$ cm$^2$ characterise the DNM phase (Planck Collaboration Int. XXVIII 2015; Remy et al. 2017; Joubaud et al. 2020). Dust opacities further increase in the dense, CO-bright, molecular phase, reaching values as high as $60 \times 10^{-27}$ cm$^2$ (Remy et al. 2017). The dust opacities we found in the HI phase of the Reticulum and North Hydrus clouds are consistent with this trend. They are close to the low diffuse cirrus values and, in particular, to the dust opacity of $(6.4 \pm 0.1) \times 10^{-27}$ cm$^2$ measured in the Eridu filament (Joubaud et al. 2020). The dust opacity we found in the DNM phase of the North Hydrus cloud fully agrees with other estimates. We caution the reader that the unusually low $\gamma_{DNM}$ value found for the Reticulum cloud may not be reliable because of the difficulty discussed in Sect. 4.2 to fully extract the $N_{H_{DNM}}$ column density along Reticulum from the γ-ray data. The final dust residual map in Fig. 6 shows a leftover positive excess along the Reticulum filament. Taking this excess into account would increase the $\gamma_{DNM}$ value.

## 44. Cosmic rays in the Reticulum and Eridu filaments

Figure 9 compares the γ-ray emissivity spectrum we found per gas atom in the Reticulum filament to that in the similar filament of Eridu and in the more massive North Hydrus cloud. We also compare it to the average $q_{LIS}$ spectrum measured in the local atomic gas at latitudes between 10° and 70° across the sky (Casandjian 2015). We compare them for the spin temperatures that best fitted their respective γ-ray data; therefore at $T_S = 100$ K for Eridu, 140 K for $q_{LIS}$, and for optically thin HI for Reticulum and North Hydrus. The results indicate that the CR population that permeates the HI and DNM phases of the Reticulum and North Hydrus clouds has the same energy distribution and same flux as the average population in the local ISM. This population applies to the 100–150 pc long sector subtended by the Reticulum and North Hydrus clouds along the rim of the Local Valley. As expected for the low $N_{HI}$ column densities characterising these clouds, we find no significant spectral

<footnote n="3">³ http://fermi.gsfc.nasa.gov/ssc/data/access/lat/BackgroundModels.html</footnote>
<figure><img src="image_9.png" /><figcaption>**Fig. 7.** Distributions of the best-fit coefficients obtained in the thousand jackknife fits to the 0.16–63 GeV γ-ray data (two left columns) and dust data (two right columns) for the optically thin HI case. The vertical lines and overlaid numbers give the mean values of the distributions. The units are $10^{26}$ cm$^{-2}$ for $q_{DNM}$, $10^{-27}$ cm$^2$ for $\gamma_{HI}$ and $\gamma_{DNM}$, and $10^{-7}$ for $\gamma_{ISO}$. $q_{HI,i}$, $q_{IC}$ and $q_{ISO}$ are dimensionless normalisation factors.</figcaption></figure>

amount to 1–6% for the HI gas γ-ray emissivities, 9% and 4% for the Ret and HyiN DNM emissivities, 8% for $q_{IC}$, 5% for $q_{ISO}$, 0.6–3% for the HI dust opacities, 12% and 3% for $\gamma_{DNM,Ret}$ and $\gamma_{DNM,HyiN}$ and 4% for $\gamma_{ISO}$. Both statistical errors ($1\sigma$) and systematic uncertainties (jackknife standard deviations) are given in Tables 1 and 2. We quadratically added them to construct the final uncertainties on the $q$ and $\gamma$ coefficients.

Figure 8 shows the stability with γ-ray energy of the best-fit coefficients relating to the Reticulum and North Hydrus clouds and to the IC and isotropic intensities. We find no significant deviation with respect to the average coefficient fitted for the entire energy band (dashed lines). This implies that the HI and DNM gas emissivity spectra in the two clouds have the same shape as the average $q_{LIS}(E)$ spectrum measured in the local ISM. We further discuss their relative normalisation with respect to $q_{LIS}$ in the next section.

Figure 8 and the value of the relative intensity ratio $q_{IC} = 1.18 \pm 0.15$ in the whole energy band further show that the IC spectrum detected in these directions of the sky is fully consistent in shape and in intensity with the GALPROP prediction based on the complex evolution of the interstellar radiation field and of the CR lepton spectra across the Milky Way. The isotropic spectrum found in this small region of the sky is also close to the larger-scale estimate obtained between 10° and 60° in latitude<fnref n="3" />. The relative intensity ratio is $q_{ISO} = 0.84 \pm 0.08$ in the 0.16–63 GeV energy band. The interstellar foregrounds are faint enough to separate the IC and isotropic components despite rather similar spatial distributions across the analysis region (see Fig. 3).

Planck Collaboration XI (2014) and Remy et al. (2017) found marked gradients in dust opacity across the sky. These variations reflect a systematic opacity increase as the gas becomes denser and the grains evolve in size and chemistry (Remy et al. 2017; Ysard et al. 2019, and references therein). Diffuse atomic cirri at high latitudes exhibit low opacities around a mean value of $(7.0 \pm 2.0) \times 10^{-27}$ cm$^2$ (Planck Collaboration XI 2014) or $(7.1 \pm 0.6) \times 10^{-27}$ cm$^2$ (Planck Collaboration Int. XVII 2014) whereas larger values ranging from $10 \times 10^{-27}$ to $16 \times 10^{-27}$ cm$^2$ are found in the denser atomic envelopes of nearby clouds (Remy et al. 2017), and values from $15 \times 10^{-27}$ to $17 \times 10^{-27}$ cm$^2$ characterise the DNM phase (Planck Collaboration Int. XXVIII 2015; Remy et al. 2017; Joubaud et al. 2020). Dust opacities further increase in the dense, CO-bright, molecular phase, reaching values as high as $60 \times 10^{-27}$ cm$^2$ (Remy et al. 2017). The dust opacities we found in the HI phase of the Reticulum and North Hydrus clouds are consistent with this trend. They are close to the low diffuse cirrus values and, in particular, to the dust opacity of $(6.4 \pm 0.1) \times 10^{-27}$ cm$^2$ measured in the Eridu filament (Joubaud et al. 2020). The dust opacity we found in the DNM phase of the North Hydrus cloud fully agrees with other estimates. We caution the reader that the unusually low $\gamma_{DNM}$ value found for the Reticulum cloud may not be reliable because of the difficulty discussed in Sect. 4.2 to fully extract the $N_{H_{DNM}}$ column density along Reticulum from the γ-ray data. The final dust residual map in Fig. 6 shows a leftover positive excess along the Reticulum filament. Taking this excess into account would increase the $\gamma_{DNM}$ value.

## 44. Cosmic rays in the Reticulum and Eridu filaments

Figure 9 compares the γ-ray emissivity spectrum we found per gas atom in the Reticulum filament to that in the similar filament of Eridu and in the more massive North Hydrus cloud. We also compare it to the average $q_{LIS}$ spectrum measured in the local atomic gas at latitudes between 10° and 70° across the sky (Casandjian 2015). We compare them for the spin temperatures that best fitted their respective γ-ray data; therefore at $T_S = 100$ K for Eridu, 140 K for $q_{LIS}$, and for optically thin HI for Reticulum and North Hydrus. The results indicate that the CR population that permeates the HI and DNM phases of the Reticulum and North Hydrus clouds has the same energy distribution and same flux as the average population in the local ISM. This population applies to the 100–150 pc long sector subtended by the Reticulum and North Hydrus clouds along the rim of the Local Valley. As expected for the low $N_{HI}$ column densities characterising these clouds, we find no significant spectral

<footnote n="3">³ http://fermi.gsfc.nasa.gov/ssc/data/access/lat/BackgroundModels.html</footnote>
<figure><img src="image_10.png" /><figcaption>**Fig. 8.** Spectral evolution relative to the local interstellar spectrum $q_{\text{LIS}}$ of the γ-ray emissivities found for the HI and DNM phases of the Reticulum and North Hydrus clouds (top and middle panels). Spectral evolution of the isotropic and IC intensities found in the region relative to the input spectra (bottom panel). $q_{\text{DNM}}/q_{\text{LIS}}$ is given in $10^{26}$ cm$^{-2}$. The results are shown for the optically thin HI case. The error bars combine the statistical and systematic uncertainties.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>**Fig. 9.** Spectral evolution of the γ-ray emissivities per gas nucleon in the Eridu (green), Reticulum (red), and North Hydrus (cyan) clouds. The orange curve gives the average $q_{\text{LIS}}$ emissivity measured in the atomic gas at latitudes $10^{\circ} \le |b| \le 30^{\circ}$ (Casandjian 2015). The emissivities are shown for the HI optical-depth correction that best fits the γ-ray data ($T_{\text{S}} = 100$ K for Eridu, 140 K for $q_{\text{LIS}}$, and optically thin for the others). The shaded bands include the statistical and systematic uncertainties.</figcaption></figure>

deviation at low energies that would mark a loss or concentration of CRs when penetrating into the denser DNM clumps (Schlickeiser et al. 2016; Bustard & Zweibel 2021).

Figure 9 illustrates that the CR flux in the Reticulum filament is $1.57 \pm 0.09$ times that pervading the similar Eridu filament. This discrepancy decreases to $1.51 \pm 0.09$ if we compare them for the same HI spin temperature of 100 K, to $1.38 \pm 0.08$ if we assume them both optically thin. The minimum ratio is $1.32 \pm 0.08$ if Eridu is optically thin and Reticulum thicker with $T_{\text{S}} = 100$ K, a case yielding a worse γ-ray fit to the data for both clouds. The loss mechanism that affects CRs in Eridu should be energy independent to preserve the same emissivity spectral shape in the three clouds across the entire energy band.

Hα hydrogen recombination line emission has been recorded along the Reticulum cloud in the Southern Alpha Sky Survey Atlas (Gaustad et al. 2001), suggesting that the cloud be a cooling filament or cooling shell surrounded by recombining warm ionised gas. The spatial correlation between the dust and γ rays detects any gas in addition to the HI gas, regardless of its chemical or ionisation state. The inferred DNM column densities are dominated by the dense CNM and diffuse H₂, but they also include the much more dilute ionised gas. The γ-ray emissivity we obtained in the HI phase of Reticulum is therefore not strongly biased by untraced ionised gas surrounding it. As a further check, we used the Hα map reprocessed by Finkbeiner (2003) as an additional gas template into the dust and γ-ray fits, even though the Hα emission measure does not scale with the $N_{\text{H}^+}$ column density, but it integrates the product of the electron and ion densities along the line of sight, $\text{EM}_\alpha \propto \text{LOS} \, \text{nen}_{\text{H}^+} \, \text{dl}$. With this crude additional template, the Reticulum γ-ray emissivity remains consistent with the local-ISM $q_{\text{LIS}}$ average and is still $1.51 \pm 0.11$ times larger than in
<figure><img src="image_10.png" /><figcaption>**Fig. 8.** Spectral evolution relative to the local interstellar spectrum $q_{\text{LIS}}$ of the γ-ray emissivities found for the HI and DNM phases of the Reticulum and North Hydrus clouds (top and middle panels). Spectral evolution of the isotropic and IC intensities found in the region relative to the input spectra (bottom panel). $q_{\text{DNM}}/q_{\text{LIS}}$ is given in $10^{26}$ cm$^{-2}$. The results are shown for the optically thin HI case. The error bars combine the statistical and systematic uncertainties.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>**Fig. 9.** Spectral evolution of the γ-ray emissivities per gas nucleon in the Eridu (green), Reticulum (red), and North Hydrus (cyan) clouds. The orange curve gives the average $q_{\text{LIS}}$ emissivity measured in the atomic gas at latitudes $10^{\circ} \le |b| \le 30^{\circ}$ (Casandjian 2015). The emissivities are shown for the HI optical-depth correction that best fits the γ-ray data ($T_{\text{S}} = 100$ K for Eridu, 140 K for $q_{\text{LIS}}$, and optically thin for the others). The shaded bands include the statistical and systematic uncertainties.</figcaption></figure>

deviation at low energies that would mark a loss or concentration of CRs when penetrating into the denser DNM clumps (Schlickeiser et al. 2016; Bustard & Zweibel 2021).

Figure 9 illustrates that the CR flux in the Reticulum filament is $1.57 \pm 0.09$ times that pervading the similar Eridu filament. This discrepancy decreases to $1.51 \pm 0.09$ if we compare them for the same HI spin temperature of 100 K, to $1.38 \pm 0.08$ if we assume them both optically thin. The minimum ratio is $1.32 \pm 0.08$ if Eridu is optically thin and Reticulum thicker with $T_{\text{S}} = 100$ K, a case yielding a worse γ-ray fit to the data for both clouds. The loss mechanism that affects CRs in Eridu should be energy independent to preserve the same emissivity spectral shape in the three clouds across the entire energy band.

Hα hydrogen recombination line emission has been recorded along the Reticulum cloud in the Southern Alpha Sky Survey Atlas (Gaustad et al. 2001), suggesting that the cloud be a cooling filament or cooling shell surrounded by recombining warm ionised gas. The spatial correlation between the dust and γ rays detects any gas in addition to the HI gas, regardless of its chemical or ionisation state. The inferred DNM column densities are dominated by the dense CNM and diffuse H₂, but they also include the much more dilute ionised gas. The γ-ray emissivity we obtained in the HI phase of Reticulum is therefore not strongly biased by untraced ionised gas surrounding it. As a further check, we used the Hα map reprocessed by Finkbeiner (2003) as an additional gas template into the dust and γ-ray fits, even though the Hα emission measure does not scale with the $N_{\text{H}^+}$ column density, but it integrates the product of the electron and ion densities along the line of sight, $\text{EM}_\alpha \propto \text{LOS} \, \text{nen}_{\text{H}^+} \, \text{dl}$. With this crude additional template, the Reticulum γ-ray emissivity remains consistent with the local-ISM $q_{\text{LIS}}$ average and is still $1.51 \pm 0.11$ times larger than in
Eridu over the whole energy band for the choice of best fitting spin temperatures.

Assuming that the same CR flux pervades the recombining and neutral hydrogen, the γ-ray signal related to the Hα template yields $N_{HI}$ column densities around $10^{19}$ cm$^{-2}$ in Reticulum, an order of magnitude below those in the neutral phase. Column densities of the Galactic warm ionised medium, integrated to the Reticulum distance, are also an order of magnitude below the $N_{HI}$ ones and they would not spatially correlate with the HI clouds to impact the values of the $q_{HI}$ emissivities (Gaensler et al. 2008).

# 5. Discussion

The Reticulum and Eridu clouds are both light and filamentary neutral clouds, with ordered magnetic fields on a few parsec scale, and a large inclination with respect to the Galactic plane. Yet, the first has $57 \pm 9\%$ more CRs than the second. In the next sections, we try to characterise the dynamical and magnetic states of those clouds at a parsec scale in order to search for differences that can influence CR transport inside them at the much smaller scales of their gyro-motion. CR momenta parallel to the local magnetic field randomly vary by interaction with magnetic perturbations. The Alfvén and slow modes of interstellar magnetohydrodynamic (MHD) turbulence are thought to be inefficient scatterers because of their large anisotropy at small scales (Chandran 2000; Xu et al. 2016). We therefore first focus on gas properties that can impact CR diffusion in the self-confinement scenario where the particles scatter on Alfvén waves they generate by the streaming instability (Kulsrud & Pearce 1969). Both interstellar and self-excited MHD turbulence should be severely damped by ion-neutral interactions in these neutral clouds, but CR gyration along tangled magnetic-field lines as well as the perpendicular (super)diffusive wandering of the lines also induce an effective CR diffusion in space, with a mean free path close to the magnetic-field coherence length (Lazarian 2016; Hu et al. 2022; Sampson et al. 2023). We therefore compare the dispersion in gas velocity and in magnetic bending in the two clouds. All CR transport modes strongly depend of the state of the MHD turbulence, so we derive the $B_{sky}$ magnetic-field strength from dust polarisation observations and we compare the sonic Mach numbers and gas-to-magnetic pressure ratios along the cores of both filaments.

## 5.1. Gas densities in the atomic phase

Given the elongated aspect of the Reticulum and Eridu filaments, we used two different assumptions for the cloud geometry, namely a cylindrical fibre or a shell seen edge on. These geometries bracket the unknown depth along the line-of-sight (LOS) in order to estimate HI gas volume densities, $n_{HI}$. In the fibre configuration, the LOS depth is close to the observed width in the plane of the sky. In order to estimate the latter, we used the Gaussian full widths at half max (FWHM) of the $N_{HI}$ column-density profiles observed along cuts perpendicular to the filament axis. Figure 10 shows the eight and seven cuts respectively chosen along the Reticulum and Eridu filaments. The resulting depths are shown in the middle panel of Fig. 11 with an uncertainty reflecting that in cloud distance. The mean depth and rms dispersion are $5.0 \pm 0.6$ pc in Reticulum and $10.5 \pm 1.0$ pc in Eridu. In the shell shaped configuration, the LOS depth is closer to half the chord length through the thick shell. Given the mean widths measured above and a shell radius of curvature of about 100 pc,

<figure><img src="image_12.png" /><figcaption>**Fig. 10.** $N_{HI}$ maps of the Reticulum (top) and Eridu (bottom) filaments, smoothed by a Gaussian kernel of 0.3' and 0.2', respectively, to highlight the cloud edges. The white lines show the perpendicular cuts used to estimate the cloud widths. The cyan squares show 1'-wide squares adopted to estimate gas volume densities, $B_{sky}$ field strengths and parallel diffusion coefficients $\kappa_{\|}$ along the filament core.</figcaption></figure>

we estimated the LOS depths to be roughly 20 pc for Reticulum and 30 pc for Eridu. These limits are represented as triangles in the middle plot of Fig. 11. They provide reasonable upper limits to the cloud depth.

We derived the mean gas volume density in a 1'-wide square centred on the peak of each $N_{HI}$ profile cutting perpendicularly through the Reticulum and Eridu filaments. The results are shown in the bottom plot of Fig. 11. The shell geometry provides lower limits to the gas densities.

Figure 11 illustrates that the larger $N_{HI}$ column densities recorded in the Eridu filament are likely due to longer LOS depths and that the Reticulum filament has denser substructures than Eridu. In both geometries, the density estimates along Reticulum are typical of the CNM and they are consistent with the detection of comparable $N_{H_{DNM}}$ column densities in the DNM. For both geometries, the Eridu volume densities are more typical of the unstable lukewarm neutral medium (LNM) that transitions between the stable WNM and CNM phases. They indicate that CNM clumps are likely present inside both clouds, but should be more numerous and more densely packed in Reticulum than in Eridu.
Eridu over the whole energy band for the choice of best fitting spin temperatures.

Assuming that the same CR flux pervades the recombining and neutral hydrogen, the γ-ray signal related to the Hα template yields $N_{HI}$ column densities around $10^{19}$ cm$^{-2}$ in Reticulum, an order of magnitude below those in the neutral phase. Column densities of the Galactic warm ionised medium, integrated to the Reticulum distance, are also an order of magnitude below the $N_{HI}$ ones and they would not spatially correlate with the HI clouds to impact the values of the $q_{HI}$ emissivities (Gaensler et al. 2008).

# 5. Discussion

The Reticulum and Eridu clouds are both light and filamentary neutral clouds, with ordered magnetic fields on a few parsec scale, and a large inclination with respect to the Galactic plane. Yet, the first has $57 \pm 9\%$ more CRs than the second. In the next sections, we try to characterise the dynamical and magnetic states of those clouds at a parsec scale in order to search for differences that can influence CR transport inside them at the much smaller scales of their gyro-motion. CR momenta parallel to the local magnetic field randomly vary by interaction with magnetic perturbations. The Alfvén and slow modes of interstellar magnetohydrodynamic (MHD) turbulence are thought to be inefficient scatterers because of their large anisotropy at small scales (Chandran 2000; Xu et al. 2016). We therefore first focus on gas properties that can impact CR diffusion in the self-confinement scenario where the particles scatter on Alfvén waves they generate by the streaming instability (Kulsrud & Pearce 1969). Both interstellar and self-excited MHD turbulence should be severely damped by ion-neutral interactions in these neutral clouds, but CR gyration along tangled magnetic-field lines as well as the perpendicular (super)diffusive wandering of the lines also induce an effective CR diffusion in space, with a mean free path close to the magnetic-field coherence length (Lazarian 2016; Hu et al. 2022; Sampson et al. 2023). We therefore compare the dispersion in gas velocity and in magnetic bending in the two clouds. All CR transport modes strongly depend of the state of the MHD turbulence, so we derive the $B_{sky}$ magnetic-field strength from dust polarisation observations and we compare the sonic Mach numbers and gas-to-magnetic pressure ratios along the cores of both filaments.

## 5.1. Gas densities in the atomic phase

Given the elongated aspect of the Reticulum and Eridu filaments, we used two different assumptions for the cloud geometry, namely a cylindrical fibre or a shell seen edge on. These geometries bracket the unknown depth along the line-of-sight (LOS) in order to estimate HI gas volume densities, $n_{HI}$. In the fibre configuration, the LOS depth is close to the observed width in the plane of the sky. In order to estimate the latter, we used the Gaussian full widths at half max (FWHM) of the $N_{HI}$ column-density profiles observed along cuts perpendicular to the filament axis. Figure 10 shows the eight and seven cuts respectively chosen along the Reticulum and Eridu filaments. The resulting depths are shown in the middle panel of Fig. 11 with an uncertainty reflecting that in cloud distance. The mean depth and rms dispersion are $5.0 \pm 0.6$ pc in Reticulum and $10.5 \pm 1.0$ pc in Eridu. In the shell shaped configuration, the LOS depth is closer to half the chord length through the thick shell. Given the mean widths measured above and a shell radius of curvature of about 100 pc,

<figure><img src="image_12.png" /><figcaption>**Fig. 10.** $N_{HI}$ maps of the Reticulum (top) and Eridu (bottom) filaments, smoothed by a Gaussian kernel of 0.3' and 0.2', respectively, to highlight the cloud edges. The white lines show the perpendicular cuts used to estimate the cloud widths. The cyan squares show 1'-wide squares adopted to estimate gas volume densities, $B_{sky}$ field strengths and parallel diffusion coefficients $\kappa_{\|}$ along the filament core.</figcaption></figure>

we estimated the LOS depths to be roughly 20 pc for Reticulum and 30 pc for Eridu. These limits are represented as triangles in the middle plot of Fig. 11. They provide reasonable upper limits to the cloud depth.

We derived the mean gas volume density in a 1'-wide square centred on the peak of each $N_{HI}$ profile cutting perpendicularly through the Reticulum and Eridu filaments. The results are shown in the bottom plot of Fig. 11. The shell geometry provides lower limits to the gas densities.

Figure 11 illustrates that the larger $N_{HI}$ column densities recorded in the Eridu filament are likely due to longer LOS depths and that the Reticulum filament has denser substructures than Eridu. In both geometries, the density estimates along Reticulum are typical of the CNM and they are consistent with the detection of comparable $N_{H_{DNM}}$ column densities in the DNM. For both geometries, the Eridu volume densities are more typical of the unstable lukewarm neutral medium (LNM) that transitions between the stable WNM and CNM phases. They indicate that CNM clumps are likely present inside both clouds, but should be more numerous and more densely packed in Reticulum than in Eridu.
<figure><img src="image_13.png" /><figcaption>Fig. 11. Distributions of HI column densities (top), LOS depths (middle), and mean HI volume densities (bottom) in the Reticulum (magenta) and Eridu (green) clouds. The offset distance is the separation between a 1'-wide square along a filament and the first one (of largest latitude or declination). The different symbols refer to the fibre or shell cloud geometry. The latter provides upper/lower limits to the depth/density. The LOS depth shading is due to the distance uncertainty. The data points and error bars in volume density respectively give the mean value in each ~5 pc-wide square and the rms dispersion at 0.5–0.7 pc scale inside each square.</figcaption></figure>

## 52. Self-confined diffusion

CRs gyrate and flow along magnetic field lines. They scatter off magnetic fluctuations with sizes close to their gyro-radius (e.g. 2 µpc for a 10 GeV CR in a $B = 5$ µG field). For instance, in the frame moving with Alfvén waves along the field, the Lorentz force acting on a CR during a gyro-resonant interaction changes the particle’s pitch angle ($\tan \theta = \mathbf{p} \wedge \mathbf{B}/\mathbf{p} \cdot \mathbf{B}$), but not its energy, so random phase shifts between the gyro-motion and the encountered waves cause a random walk in pitch angle and reduce the CR effective speed along the field line. See Ruszkowski & Pfrommer (2023) for a nice pictural description of gyro-resonant scatterings in pitch angle. Whereas very high energy CRs are likely scattered by interstellar MHD turbulence, the lower-energy CRs that we sample here in γ rays are more likely scattered by Alfvén waves that they self-generate via the resonant streaming instability (Kulsrud & Pearce 1969). The anisotropy of interstellar MHD turbulence increases with decreasing scale in the cascade (Goldreich & Sridhar 1995). The

eddies become elongated along the field because of turbulent reconnection. As a result, during the resonant encounter of a CR with a gyro-radius-long eddy, the gyro-orbit perpendicular to the field encompasses other uncorrelated eddies that cancel out the net change in pitch angle (Chandran 2000). Unlike Alfvén and slow modes, fast compressive turbulent modes remain much more isotropic, and are therefore efficient scatterers on small scales. In the neutral atomic gas, however, they mostly affect CRs of much larger energies (>>100 TeV) than those of interest here (Xu et al. 2016).

In the self-confinement scenario, a slight anisotropy in the CR pitch-angle distribution can excite Alfvén waves with the right geometry to induce further pitch-angle scattering (Kulsrud & Pearce 1969). The required level of anisotropy to trigger the instability corresponds to a bulk drift velocity $v_D$ larger than the ion Alfvén speed $v_{ion} = B/4\pi\rho_i$ (the drift velocity is that of the frame where the CR distribution is isotropic). The wave growth rate by this streaming instability depends on the integrated number density of exciting CRs and on the mass density of interstellar gas ions, $\rho_i$, as:

$$ \Gamma_{\text{grow}}(k) = \frac{\pi^{3/2} e}{2 c} (v_D - v_A^{\text{ion}}) \frac{1}{\sqrt{\rho_i}} \times \frac{a_j - 3}{a_j - 2} Z_j n_j (R_j > \frac{B}{k}), \quad (3) $$

where $B$ is the mean field strength and we add the contributions of different CR nuclei with charge $Z_j e$, rigidity $R_j$, and volume number density spectrum $n_j(R)$ with a power-law spectral index $a_j$ in the range of interest (Kulsrud & Pearce 1969; Kulsrud & Cesarsky 1971). The CR density is integrated in rigidity above the threshold corresponding to the resonance condition $R_j = B/k$ for the wave number $k$. We considered the spectra of CRs with $Z_j \le 28$ that reach the heliopause (Boschini et al. 2020). They correspond to the $q_{\text{LIS}}$ emissivity discussed above. We integrated the spectra over 1 GV in rigidity and fitted the $a_j$ spectral indices in the 1–20 GV range that contributes most to the growth rate.

The actual level of CR anisotropy that causes the streaming instability and the growth of scattering waves is extremely difficult to estimate from observations and from theory. We therefore follow Jiang & Oh (2018) and Bustard & Zweibel (2021) in expressing the net drift velocity $v_D - v_A^{\text{ion}}$ with respect to the Alfvén wave frame as a diffusive flux:

$$ v_D - v_A^{\text{ion}} = \frac{\kappa_{\parallel}}{\ell_{\text{CR}}}, \quad (4) $$

with the characteristic length for a change in CR flux, $\ell_{\text{CR}}$, and the diffusion coefficient $\kappa_{\parallel}$ parallel to the mean B field.

Wave amplitudes and scattering rates are reduced by wave damping processes that depend on the properties of the background MHD turbulence and on ambient plasma conditions. The main damping mechanism in the neutral gas phases that we have probed in γ rays is the ion-neutral one. Following Hopkins et al. (2021) and Xu & Lazarian (2022), we checked that the linear Landau, non-linear Landau, and turbulent damping rates are two to three orders of magnitude smaller than the ion-neutral one in the typical environments of the Reticulum and Eridu clouds.

Momentum exchange between the ions and hydrogen or helium neutral atoms, all with the same equilibrium temperature $T$, leads to the damping rate given by Soler et al. (2016):

$$ \Gamma_{\text{in,S16}} = \frac{32 k_B T}{9 \pi A_i m_H} n \frac{\Box}{\Box \Box \Box \Box \Box \Box \Box \Box} \frac{\Box}{\Box} + \frac{1}{1 + A_i} + \sigma \frac{x_{\text{He}}}{i\text{He}} \frac{\Box}{x_{\text{HI}}} \frac{A_{\text{He}}}{A_i + A_{\text{He}}} \quad (5) $$
<figure><img src="image_14.png" /><figcaption>Fig. 12. Parallel diffusion coefficients obtained for the Reticulum (top) and Eridu (bottom) clouds in the fibre (left) and shell (right) geometries, assuming the damping rate of Eq. (5), a characteristic CR flux variation length $l_{CR} = 400$ pc, and a uniform LOS depth equal to the shell one or to the average of the fibre ones. For each geometry, the Eridu plots show $\kappa_{\parallel}$ values for the same CR flux exciting the Alfvén waves as in Reticulum (left) or for the reduced CR flux corresponding to the cloud lower γ-ray emissivity (right).</figcaption></figure>

where $k_B$ is the Boltzmann constant, $m_H$ is the H atom mass, $x_{HI}$ and $x_{He}$ are the neutral hydrogen and helium abundance by number, relative to the total number of hydrogen nucleons, and $A_{He}=4$ is the atomic weight of helium. The ion-neutral momentum-transfer cross sections are $\sigma_{HI} = 10^{-14}$ cm$^2$ and $\sigma_{iHe} = 3 \times 10^{-15}$ cm$^2$ from Vranjes & Krstic (2013).

The ion abundance, $x_i$, and mean ion atomic weight, $A_i$, vary with the gas state (Draine 2011). The predominant form of ions gradually varies from $H^+$ in the WNM to about half $H^+$ and half $C^+$ in the CNM, and their total abundance decreases from about 2% in the WNM down by two orders of magnitude in the CNM. We used the heating, cooling, chemistry, and ionisation model of Kim et al. (2023) to calculate the equilibrium gas temperature $T$, ion abundance $x_i$, and average ion weight $A_i$ as a function of the gas densities observed in the neutral phase ($n_{HI}$) of the two clouds, and for the interstellar conditions of the local ISM ($x_{He} = 0.1$, $x_C = 1.6 \times 10^{-4}$ in the gas phase, solar metallicity, primary CR ionisation rate of $2 \times 10^{-16}$ s$^{-1}$, and local FUV intensity). This allowed us to estimate the damping rates of Alfvén waves in different locations along both clouds.

Another formalism for ion-neutral damping is given by Kulsrud & Pearce (1969):

$$ \Gamma_{\text{in,}K69} = \frac{1}{2} \frac{n_i m_i}{m_n + m_i} \langle \sigma \nu \rangle_{\text{in},i} \quad (6) $$

where $n_n$ is the neutral number density, $m_n$ is the mean mass of neutrals, and $m_i = A_i m_H$ is the mean mass of ions. The ion-neutral momentum transfer rates are $\langle \sigma \nu \rangle_{\text{in}} = 3.25 \times 10^{-9}$ cm$^3$ s$^{-1}$ for collisions between H and $H^+$, $2.39 \times 10^{-9}$ cm$^3$ s$^{-1}$ between H and $C^+$, and $1.42 \times 10^{-9}$ cm$^3$ s$^{-1}$ between He and

$H^+$ (from Table 2.1 of Draine 2011). This formulation does not depend on the gas temperature.

We adopted the Soler et al. (2016) damping rate as the baseline for this work. Assuming a uniform LOS depth equal to the mean value of the estimates shown in Fig. 11 in each geometry, we find rather uniform damping rates close to $6 \times 10^{-9}$ s$^{-1}$ across both clouds, in both geometries. The two damping rate formulations yield comparable values across the clouds in the shell geometry, but not in the denser fiber configuration where the Kulsrud & Pearce (1969) formulation yields a few times larger rates than the Soler et al. (2016) ones (one to three times larger along Eridu and one to five times larger along Reticulum). This difference is due to the drop in temperature and increase in ion atomic weight in the denser fibre gas in the Soler et al. (2016) formulation.

We further note that these rates likely over-predict the actual damping rate in the partially ionised plasma because they do not take the charged dust grains into account even though each grain carries a charge that varies from 0.1e (small grains) to 20e or 30e (large grains) in the atomic gas phases of interest here (Draine 2004). The presence of dust charges reduces the dissipation of Alfvén waves by at least an order of magnitude (Hennebelle & Lebreuilly 2023). Charged grains also lower the drift velocity required for CR excitation of the waves. The gas-to-dust mass ratio being of order 100, the ion Alfvén speed is limited to about ten times the gas one ($v_{A}^{\text{ion}}/v_{A}^{\text{gas}} \approx \rho/\rho_{\text{dust}}$) whereas, if dust is ignored, one expects much larger $v_{A}^{\text{ion}}$ speeds in all gas densities above $1 \text{ cm}^{-3}$ where the abundance of gaseous ions drastically falls below 1% (Kim et al. 2023). Charged dust therefore allows larger wave intensities for gyro-resonant scattering in the CNM and denser gas phases (Hennebelle & Lebreuilly 2023).
<figure><img src="image_13.png" /><figcaption>Fig. 12. Parallel diffusion coefficients obtained for the Reticulum (top) and Eridu (bottom) clouds in the fibre (left) and shell (right) geometries, assuming the damping rate of Eq. (5), a characteristic CR flux variation length $l_{CR} = 400$ pc, and a uniform LOS depth equal to the shell one or to the average of the fibre ones. For each geometry, the Eridu plots show $\kappa_{\parallel}$ values for the same CR flux exciting the Alfvén waves as in Reticulum (left) or for the reduced CR flux corresponding to the cloud lower γ-ray emissivity (right).</figcaption></figure>

where $k_B$ is the Boltzmann constant, $m_H$ is the H atom mass, $x_{HI}$ and $x_{He}$ are the neutral hydrogen and helium abundance by number, relative to the total number of hydrogen nucleons, and $A_{He}=4$ is the atomic weight of helium. The ion-neutral momentum-transfer cross sections are $\sigma_{HI} = 10^{-14}$ cm$^2$ and $\sigma_{iHe} = 3 \times 10^{-15}$ cm$^2$ from Vranjes & Krstic (2013).

The ion abundance, $x_i$, and mean ion atomic weight, $A_i$, vary with the gas state (Draine 2011). The predominant form of ions gradually varies from $H^+$ in the WNM to about half $H^+$ and half $C^+$ in the CNM, and their total abundance decreases from about 2% in the WNM down by two orders of magnitude in the CNM. We used the heating, cooling, chemistry, and ionisation model of Kim et al. (2023) to calculate the equilibrium gas temperature $T$, ion abundance $x_i$, and average ion weight $A_i$ as a function of the gas densities observed in the neutral phase ($n_{HI}$) of the two clouds, and for the interstellar conditions of the local ISM ($x_{He} = 0.1$, $x_C = 1.6 \times 10^{-4}$ in the gas phase, solar metallicity, primary CR ionisation rate of $2 \times 10^{-16}$ s$^{-1}$, and local FUV intensity). This allowed us to estimate the damping rates of Alfvén waves in different locations along both clouds.

Another formalism for ion-neutral damping is given by Kulsrud & Pearce (1969):

$$ \Gamma_{\text{in,K69}} = \frac{1}{2} \frac{n_i m_i}{m_n + m_i} \langle \sigma \nu \rangle_{\text{in},i} \quad (6) $$

where $n_n$ is the neutral number density, $m_n$ is the mean mass of neutrals, and $m_i = A_i m_H$ is the mean mass of ions. The ion-neutral momentum transfer rates are $\langle \sigma \nu \rangle_{\text{in}} = 3.25 \times 10^{-9}$ cm$^3$ s$^{-1}$ for collisions between H and $H^+$, $2.39 \times 10^{-9}$ cm$^3$ s$^{-1}$ between H and $C^+$, and $1.42 \times 10^{-9}$ cm$^3$ s$^{-1}$ between He and

$H^+$ (from Table 2.1 of Draine 2011). This formulation does not depend on the gas temperature.

We adopted the Soler et al. (2016) damping rate as the baseline for this work. Assuming a uniform LOS depth equal to the mean value of the estimates shown in Fig. 11 in each geometry, we find rather uniform damping rates close to $6 \times 10^{-9}$ s$^{-1}$ across both clouds, in both geometries. The two damping rate formulations yield comparable values across the clouds in the shell geometry, but not in the denser fiber configuration where the Kulsrud & Pearce (1969) formulation yields a few times larger rates than the Soler et al. (2016) ones (one to three times larger along Eridu and one to five times larger along Reticulum). This difference is due to the drop in temperature and increase in ion atomic weight in the denser fibre gas in the Soler et al. (2016) formulation.

We further note that these rates likely over-predict the actual damping rate in the partially ionised plasma because they do not take the charged dust grains into account even though each grain carries a charge that varies from 0.1e (small grains) to 20e or 30e (large grains) in the atomic gas phases of interest here (Draine 2004). The presence of dust charges reduces the dissipation of Alfvén waves by at least an order of magnitude (Hennebelle & Lebreuilly 2023). Charged grains also lower the drift velocity required for CR excitation of the waves. The gas-to-dust mass ratio being of order 100, the ion Alfvén speed is limited to about ten times the gas one ($v_A^{\text{ion}}/v_A^{\text{gas}} \approx \rho/\rho_{\text{dust}}$) whereas, if dust is ignored, one expects much larger $v_A^{\text{ion}}$ speeds in all gas densities above $1 \text{ cm}^{-3}$ where the abundance of gaseous ions drastically falls below 1% (Kim et al. 2023). Charged dust therefore allows larger wave intensities for gyro-resonant scattering in the CNM and denser gas phases (Hennebelle & Lebreuilly 2023).
<figure><img src="image_15.png" /><figcaption>Fig. 13. Estimates of the parallel diffusion coefficient in the 1⁻-wide squares chosen along the core of the Reticulum (magenta) and Eridu (green) filaments, for the fibre (top) and shell (bottom) geometries, for the same CR flux pervading the two clouds (filled symbols) and for the reduced CR flux in Eridu (open symbols). The data points and error bars respectively give the mean value in each ~5 pc-wide square and the rms dispersion at 0.5–0.7 pc scale inside each square. Different n<sub>H</sub> and κ<sub>||</sub> scales are used for each geometry.</figcaption></figure>

A quasi-steady state is reached when the growth rate Γ<sub>grow</sub> of the gyro-resonant waves is balanced by the damping rate Γ<sub>in</sub>, giving a parallel diffusion coefficient:

$$ \kappa_{\parallel} = \frac{2 c \Gamma_{\text{in}} \ell_{\text{CR}}}{\pi^{3/2} e} \frac{\sqrt{\rho_i}}{-\frac{\alpha_{j-3}}{\alpha_{j-2}} Z_j n_j(R_j > \frac{B}{k})} \propto T^{1/2} x_i^{1/2} n_{\text{HI}}^{3/2} \ell_{\text{CR}}. \quad (7) $$

We note that, under the steady-state self-confinement transport assumption, the diffusion coefficient depends on two opposite trends that limit its potential range inside a cloud since both the gas temperature and ionisation fraction rapidly drop as the gas density increases. We also note that the $κ_{||}/ℓ_{CR}$ formulation for the net drift velocity cancels the $κ_{||}$ dependence on the ion Alfvén speed, and therefore also the dependence on the mean magnetic field strength.

The distance that separates the Reticulum and Eridu filaments ranges from 170 pc to 270 pc given their respective length and distance uncertainty to the Sun. They are also about 150 pc away from other clouds exhibiting γ-ray emissivities close to the average $q_{\text{LIS}}$. Given the 40–60% change in CR flux between Eridu and the other clouds, we have adopted $ℓ_{CR}=400$ pc as a typical CR gradient length scale for the CR energies probed with the LAT. This macroscopic length scale is motivated by differences in CR flux inside various clouds and not by the actual CR gradients that prevail in the immediate vicinity of the Eridu and Reticulum clouds, or those that are induced by magnetic bottlenecks, $v_A^{\text{ion}}$ changes, or CR losses in the cloud envelopes (Schlickeiser et al. 2016; Bustard & Zweibel 2021). Given its uncertainty, we explicitly quote the $ℓ_{CR}$ dependence in our results.

Figure 12 presents maps of the parallel diffusion coefficient obtained for both clouds when using a uniform LOS depth across the whole cloud (the mean one in the fibre configuration or the

upper limit in the shell geometry). We further relaxed the uniformity of the LOS depth by calculating the diffusion coefficient in each of the 1°×1° squares chosen along both filaments, using the fibre LOS depth estimated in each square. The results are shown in Tables B.1 and B.2 and as a function of the mean HI gas density in each square in Fig. 13. For a characteristic CR gradient scale length $ℓ_{CR}$ of 400 pc, we find $κ_{||}$ values that typically range from 6 ×10<sup>28</sup> to 14 ×10<sup>28</sup> cm<sup>2</sup> s<sup>-1</sup> across the Reticulum cloud in both geometries, and over a smaller range from 5 ×10<sup>28</sup> to 8 ×10<sup>28</sup> cm<sup>2</sup> s<sup>-1</sup> for Eridu because of its lower contrast in gas column density. The median $κ_{||}$ value in all the pixels composing the squares along the core of each filament differs by less than 20% between the two clouds in each geometry.

The maps show that the diffusion coefficient does not vary substantially as the gas density increases from the cloud envelopes to the cores because the increase in gas density is tempered by the loss in ion abundance and because the damping rise is tempered by the drop in gas temperature and increase in ion mass. The overall variation from edge to core is of the order of two to three. We stress that these $κ_{||}$ estimates are based on the cloud properties sampled at parsec scales rather than at the actual micro-parsec scales of gyro-resonance, so in-situ variations may well be much more severe.

Figure 12 illustrates that, when using the same CR flux to excite Alfvén waves in the Reticulum and Eridu clouds, both filaments exhibit very similar diffusion properties over most of their area in each of the fibre or shell configurations. The larger $κ_{||}$ values found along the spine of Reticulum in Fig. 13 correspond to the densest cores where CRs can flow more rapidly due to the enhanced damping. In the self-confinement theory, CRs create the magnetic fluctuations they interact with, so fewer CRs have less power to generate perturbations to modify their pitch angle. The CR flux reduction in Eridu with respect to Reticulum leads to an increase in diffusion length by the flux ratio of 1.57 ± 0.09. Taking this ratio into account, Fig. 12 shows that CRs should diffuse at about the same speed in the densest parts of Eridu and Reticulum, but twice as fast at the periphery of Eridu than in the envelope of Reticulum, in both geometries. However, the diffusivity difference between the two clouds is too small and the uncertainties in the growth and damping rates are too large to relate the loss of CRs in Eridu to different self-confinement efficiencies.

Cosmic rays transfer energy to the magnetic fluctuations if they have an anisotropic distribution with a degree of anisotropy larger than $v_A^{\text{ion}}/c$, therefore if their drift velocity $v_D$ along the mean magnetic field exceeds $v_A^{\text{ion}}$ (Kulsrud & Pearce 1969). If the CRs are well scattered, there is a direct relationship between their anisotropy and spatial density gradient (Zweibel 2013), so a difference in CR flux inside Eridu and Reticulum may originate from different CR gradients in their immediate environment (different $ℓ_{CR}$ values). Both clouds lie near the dense walls of the Local Valley, but Eridu stands in the vicinity of the Orion-Eridanus superbubble where steep CR gradients are possible even though we have found no evidence for fresh acceleration or reacceleration of CRs in the compressed outer rim of the superbubble (Joubaud et al. 2020). Steeper gradients, however, would excite more streaming instabilities and decrease the diffusion coefficient rather than increase it compared to Reticulum.

On the other hand, Reticulum may be a cooling fragment of an old shell (possibly from a supernova remnant) surrounded by steep CR gradients to efficiently confine the particles inside the shell, but a local boost of the wave growth rate inside this cloud is hard to reconcile with the fact that its CR flux matches that found
<figure><img src="image_14.png" /><figcaption>Fig. 13. Estimates of the parallel diffusion coefficient in the 1⁻-wide squares chosen along the core of the Reticulum (magenta) and Eridu (green) filaments, for the fibre (top) and shell (bottom) geometries, for the same CR flux pervading the two clouds (filled symbols) and for the reduced CR flux in Eridu (open symbols). The data points and error bars respectively give the mean value in each ~5 pc-wide square and the rms dispersion at 0.5–0.7 pc scale inside each square. Different n<sub>H</sub> and κ<sub>||</sub> scales are used for each geometry.</figcaption></figure>

A quasi-steady state is reached when the growth rate Γ<sub>grow</sub> of the gyro-resonant waves is balanced by the damping rate Γ<sub>in</sub>, giving a parallel diffusion coefficient:

$$ \kappa_{\parallel} = \frac{2 c \Gamma_{\text{in}} \ell_{\text{CR}}}{\pi^{3/2} e} \frac{\sqrt{\rho_i}}{-\frac{\alpha_{j-3}}{\alpha_{j-2}} Z_j n_j(R_j > \frac{B}{k})} \propto T^{1/2} x_i^{1/2} n_{\text{HI}}^{3/2} \ell_{\text{CR}}. \quad (7) $$

We note that, under the steady-state self-confinement transport assumption, the diffusion coefficient depends on two opposite trends that limit its potential range inside a cloud since both the gas temperature and ionisation fraction rapidly drop as the gas density increases. We also note that the $κ_{||}/ℓ_{CR}$ formulation for the net drift velocity cancels the $κ_{||}$ dependence on the ion Alfvén speed, and therefore also the dependence on the mean magnetic field strength.

The distance that separates the Reticulum and Eridu filaments ranges from 170 pc to 270 pc given their respective length and distance uncertainty to the Sun. They are also about 150 pc away from other clouds exhibiting γ-ray emissivities close to the average $q_{\text{LIS}}$. Given the 40–60% change in CR flux between Eridu and the other clouds, we have adopted $ℓ_{CR}=400$ pc as a typical CR gradient length scale for the CR energies probed with the LAT. This macroscopic length scale is motivated by differences in CR flux inside various clouds and not by the actual CR gradients that prevail in the immediate vicinity of the Eridu and Reticulum clouds, or those that are induced by magnetic bottlenecks, $v_A^{\text{ion}}$ changes, or CR losses in the cloud envelopes (Schlickeiser et al. 2016; Bustard & Zweibel 2021). Given its uncertainty, we explicitly quote the $ℓ_{CR}$ dependence in our results.

Figure 12 presents maps of the parallel diffusion coefficient obtained for both clouds when using a uniform LOS depth across the whole cloud (the mean one in the fibre configuration or the

upper limit in the shell geometry). We further relaxed the uniformity of the LOS depth by calculating the diffusion coefficient in each of the 1°×1° squares chosen along both filaments, using the fibre LOS depth estimated in each square. The results are shown in Tables B.1 and B.2 and as a function of the mean HI gas density in each square in Fig. 13. For a characteristic CR gradient scale length $ℓ_{CR}$ of 400 pc, we find $κ_{||}$ values that typically range from 6 ×10<sup>28</sup> to 14 ×10<sup>28</sup> cm<sup>2</sup> s<sup>-1</sup> across the Reticulum cloud in both geometries, and over a smaller range from 5 ×10<sup>28</sup> to 8 ×10<sup>28</sup> cm<sup>2</sup> s<sup>-1</sup> for Eridu because of its lower contrast in gas column density. The median $κ_{||}$ value in all the pixels composing the squares along the core of each filament differs by less than 20% between the two clouds in each geometry.

The maps show that the diffusion coefficient does not vary substantially as the gas density increases from the cloud envelopes to the cores because the increase in gas density is tempered by the loss in ion abundance and because the damping rise is tempered by the drop in gas temperature and increase in ion mass. The overall variation from edge to core is of the order of two to three. We stress that these $κ_{||}$ estimates are based on the cloud properties sampled at parsec scales rather than at the actual micro-parsec scales of gyro-resonance, so in-situ variations may well be much more severe.

Figure 12 illustrates that, when using the same CR flux to excite Alfvén waves in the Reticulum and Eridu clouds, both filaments exhibit very similar diffusion properties over most of their area in each of the fibre or shell configurations. The larger $κ_{||}$ values found along the spine of Reticulum in Fig. 13 correspond to the densest cores where CRs can flow more rapidly due to the enhanced damping. In the self-confinement theory, CRs create the magnetic fluctuations they interact with, so fewer CRs have less power to generate perturbations to modify their pitch angle. The CR flux reduction in Eridu with respect to Reticulum leads to an increase in diffusion length by the flux ratio of 1.57 ± 0.09. Taking this ratio into account, Fig. 12 shows that CRs should diffuse at about the same speed in the densest parts of Eridu and Reticulum, but twice as fast at the periphery of Eridu than in the envelope of Reticulum, in both geometries. However, the diffusivity difference between the two clouds is too small and the uncertainties in the growth and damping rates are too large to relate the loss of CRs in Eridu to different self-confinement efficiencies.

Cosmic rays transfer energy to the magnetic fluctuations if they have an anisotropic distribution with a degree of anisotropy larger than $v_A^{\text{ion}}/c$, therefore if their drift velocity $v_D$ along the mean magnetic field exceeds $v_A^{\text{ion}}$ (Kulsrud & Pearce 1969). If the CRs are well scattered, there is a direct relationship between their anisotropy and spatial density gradient (Zweibel 2013), so a difference in CR flux inside Eridu and Reticulum may originate from different CR gradients in their immediate environment (different $ℓ_{CR}$ values). Both clouds lie near the dense walls of the Local Valley, but Eridu stands in the vicinity of the Orion-Eridanus superbubble where steep CR gradients are possible even though we have found no evidence for fresh acceleration or reacceleration of CRs in the compressed outer rim of the superbubble (Joubaud et al. 2020). Steeper gradients, however, would excite more streaming instabilities and decrease the diffusion coefficient rather than increase it compared to Reticulum.

On the other hand, Reticulum may be a cooling fragment of an old shell (possibly from a supernova remnant) surrounded by steep CR gradients to efficiently confine the particles inside the shell, but a local boost of the wave growth rate inside this cloud is hard to reconcile with the fact that its CR flux matches that found
in other local clouds with different histories and environments. The detection of Hα recombination emission towards Reticulum reveals the presence of ionised gas close to the cloud. Efficient streaming and gyro-scattering in the ionised gas may have reduced the CR anisotropy before entering the cloud, thereby preventing the further excitation of scattering waves inside the cloud. But this scenario would lead to faster diffusion in Reticulum than in Eridu since we find no ionised shell around the latter.

We note that $κ_{||}$ linearly scales with $ℓ_{CR}$. This is the reason why the present estimates drastically differ from those found by Armillotta et al. (2021) for the same gas density. They find very smooth CR densities in the CNM with about two orders of magnitude larger gradient scale lengths in their simulations (Lucia Armillotta, priv. comm.).

The $κ_{||}$ values we find in these individual clouds in the self-confinement scenario compare well with the range of recent estimates obtained for the entire Milky Way for isotropic diffusion. The latter vary from $4 \times 10^{28}$ to $20 \times 10^{28}$ cm$^2$ s$^{-1}$ for 10 GeV/n CRs (Jóhannesson et al. 2016, 2019; Weinrich et al. 2020; De La Torre Luque et al. 2022). They are based on direct measurements of CR nuclear spectra measured inside and near the heliosphere (by Voyager, ACE-CRIS, HEAO-3, CREAM, AMS-02, and PAMELA) and on different propagation codes to model spallation reactions (namely GALPROP, USINE, and DRAGON-2). The different estimates vary by a factor five even though the models merely try to constrain an isotropic, homogeneous, environment-independent coefficient in the Galaxy. The convergence between the Galactic averages and the values presented in Figs. 12 and 13 is nice, but it primarily relies on the choice of 400 pc for the $ℓ_{CR}$ gradient scale and on the assumption of CR self-confinement in these clouds. As discussed above, variations in the cloud properties (temperature, ionisation fraction, and gas density) tend to compensate each other to reduce the range of $κ_{||}$ values in this transport scenario. The choice of $ℓ_{CR}$ was not tuned to match Milky Way values, but motivated by the 3D separation between clouds exhibiting different γ-ray emissivities. Nevertheless, while we can compare the diffusion coefficients presented in Figs. 12 and 13 for the two clouds because they rely on the same assumptions and on the same spatial scale and methods to measure the relevant interstellar parameters, we should consider the absolute value of these coefficients with great care because of the unknown state of the plasma and unknown level of CR anisotropy at the gyro-resonance scale.

## 5.3. Velocity dispersion and sonic Mach numbers

As indicated when opening the discussion section, other CR transport modes relate to the turbulent wandering and bending of magnetic-field lines, in particular in super-Alfvénic turbulence (Lazarian & Xu 2021; Hu et al. 2022; Xu & Lazarian 2022; Sampson et al. 2023, and references therein). We therefore explored the turbulence level in gas velocity and in $B_{sky}$ orientation in each cloud, and we estimated the sonic Mach numbers of the turbulence to gauge the plasma compressibility.

We took advantage of the pseudo-Voigt decomposition of the HI spectra we used for the cloud separation to follow changes in LOS velocity for the main HI line across each cloud. Figure 14 presents the maps of the $σ_v$ line widths of the main line across Eridu (top panel) and Reticulum (bottom panel). We also mapped the position-to-position velocity dispersion at the smallest scale resolved by the GASS using a two-point structure

<figure><img src="image_16.png" /><figcaption>Fig. 14. Maps of the velocity line width, $σ_v$, and of the velocity dispersion function, $S_v(δ_{GASS})$, of the main HI emission line across the Eridu (top) and Reticulum (bottom) clouds. We note the different colour scale for each cloud.</figcaption></figure>

function in LOS velocity:

$$S_v(x, \delta) = \frac{\sum_{i=1}^{N} \gamma_v [v_{LOS}(x) - v_{LOS}(x + \delta_i)]^2}{\gamma_v} \quad (8)$$

where the velocity dispersion compares the line central velocity in a given direction $x$ and in a displaced direction $v(x + \delta)$. The average is taken over displaced directions spanning an annulus $δ/2 \le \delta_i \le 3δ/2$ around the $x$ direction, with $δ$ equal to the effective GASS FWHM angular resolution of 16.2 (Kalberla et al. 2010). The resulting $S_v$ maps are shown for the Reticulum and Eridu clouds in the top and bottom panels of Fig. 14, respectively.

The $σ_v$ line widths are commensurate with the thermal broadening expected from the temperature map inferred in the previous section. The velocity resolution of GASS and the uncertainty in the modelled gas temperature do not allow a reliable subtraction of the thermal broadening to measure the dispersion in turbulent velocities. The $σ_v$ line widths therefore provide upper limits to the turbulent motions along the lines of sight. On the other hand, because of the average over the $δ$ ring, the $S_v$ structure function is closer to a lower limit to the turbulent motions along the line of sight. Both $S_v$ and $σ_v$ are one-dimensional measures of the turbulent velocity field.
in other local clouds with different histories and environments. The detection of Hα recombination emission towards Reticulum reveals the presence of ionised gas close to the cloud. Efficient streaming and gyro-scattering in the ionised gas may have reduced the CR anisotropy before entering the cloud, thereby preventing the further excitation of scattering waves inside the cloud. But this scenario would lead to faster diffusion in Reticulum than in Eridu since we find no ionised shell around the latter.

We note that $κ_{||}$ linearly scales with $ℓ_{CR}$. This is the reason why the present estimates drastically differ from those found by Armillotta et al. (2021) for the same gas density. They find very smooth CR densities in the CNM with about two orders of magnitude larger gradient scale lengths in their simulations (Lucia Armillotta, priv. comm.).

The $κ_{||}$ values we find in these individual clouds in the self-confinement scenario compare well with the range of recent estimates obtained for the entire Milky Way for isotropic diffusion. The latter vary from $4 \times 10^{28}$ to $20 \times 10^{28}$ cm$^2$ s$^{-1}$ for 10 GeV/n CRs (Jóhannesson et al. 2016, 2019; Weinrich et al. 2020; De La Torre Luque et al. 2022). They are based on direct measurements of CR nuclear spectra measured inside and near the heliosphere (by Voyager, ACE-CRIS, HEAO-3, CREAM, AMS-02, and PAMELA) and on different propagation codes to model spallation reactions (namely GALPROP, USINE, and DRAGON-2). The different estimates vary by a factor five even though the models merely try to constrain an isotropic, homogeneous, environment-independent coefficient in the Galaxy. The convergence between the Galactic averages and the values presented in Figs. 12 and 13 is nice, but it primarily relies on the choice of 400 pc for the $ℓ_{CR}$ gradient scale and on the assumption of CR self-confinement in these clouds. As discussed above, variations in the cloud properties (temperature, ionisation fraction, and gas density) tend to compensate each other to reduce the range of $κ_{||}$ values in this transport scenario. The choice of $ℓ_{CR}$ was not tuned to match Milky Way values, but motivated by the 3D separation between clouds exhibiting different γ-ray emissivities. Nevertheless, while we can compare the diffusion coefficients presented in Figs. 12 and 13 for the two clouds because they rely on the same assumptions and on the same spatial scale and methods to measure the relevant interstellar parameters, we should consider the absolute value of these coefficients with great care because of the unknown state of the plasma and unknown level of CR anisotropy at the gyro-resonance scale.

## 5.3. Velocity dispersion and sonic Mach numbers

As indicated when opening the discussion section, other CR transport modes relate to the turbulent wandering and bending of magnetic-field lines, in particular in super-Alfvénic turbulence (Lazarian & Xu 2021; Hu et al. 2022; Xu & Lazarian 2022; Sampson et al. 2023, and references therein). We therefore explored the turbulence level in gas velocity and in $B_{sky}$ orientation in each cloud, and we estimated the sonic Mach numbers of the turbulence to gauge the plasma compressibility.

We took advantage of the pseudo-Voigt decomposition of the Hα spectra we used for the cloud separation to follow changes in LOS velocity for the main Hα line across each cloud. Figure 14 presents the maps of the $σ_v$ line widths of the main line across Eridu (top panel) and Reticulum (bottom panel). We also mapped the position-to-position velocity dispersion at the smallest scale resolved by the GASS using a two-point structure

<figure><img src="image_15.png" /><figcaption>Fig. 14. Maps of the velocity line width, $σ_v$, and of the velocity dispersion function, $S_v(δ_{GASS})$, of the main Hα emission line across the Eridu (top) and Reticulum (bottom) clouds. We note the different colour scale for each cloud.</figcaption></figure>

function in LOS velocity:

$$S_v(x, \delta) = \frac{\sum_{i=1}^{N} \gamma_i [v_{LOS}(x) - v_{LOS}(x + \delta_i)]^2}{N} \quad (8)$$

where the velocity dispersion compares the line central velocity in a given direction $x$ and in a displaced direction $v(x + \delta)$. The average is taken over displaced directions spanning an annulus $δ/2 \le \delta_i \le 3δ/2$ around the $x$ direction, with $δ$ equal to the effective GASS FWHM angular resolution of 16.2 (Kalberla et al. 2010). The resulting $S_v$ maps are shown for the Reticulum and Eridu clouds in the top and bottom panels of Fig. 14, respectively.

The $σ_v$ line widths are commensurate with the thermal broadening expected from the temperature map inferred in the previous section. The velocity resolution of GASS and the uncertainty in the modelled gas temperature do not allow a reliable subtraction of the thermal broadening to measure the dispersion in turbulent velocities. The $σ_v$ line widths therefore provide upper limits to the turbulent motions along the lines of sight. On the other hand, because of the average over the $δ$ ring, the $S_v$ structure function is closer to a lower limit to the turbulent motions along the line of sight. Both $S_v$ and $σ_v$ are one-dimensional measures of the turbulent velocity field.
<figure><img src="image_17.png" /><figcaption>**Fig. 15.** Sonic Mach numbers of turbulent motions in the 1°-wide squares chosen along the Eridu (green) and Reticulum (magenta) filaments, in the fibre (top) and shell (bottom) geometries, assuming isotropic turbulence and using the $S_v$ dispersion function (dark filled symbols) or the $\sigma_v$ line width (light open symbols) to trace the one-dimensional turbulent velocity dispersion. The data points and error bars respectively give the mean value in each ~5 pc-wide square and the rms dispersion at 0.5–0.7 pc scale inside each square. The histogram shows the Mach number distribution in cold CNM clumps from Heiles & Troland (2003).</figcaption></figure>

Figure 14 shows that the line width and velocity dispersion values decrease from the cloud edges to the cores of both clouds, as expected from CNM regions. It also indicates that the dispersion in the velocity field of the Eridu cloud is about twice larger than in the more compact Reticulum filament (both in $S_v$ and $\sigma_v$). This is expected for a more diffuse cirrus (Audit & Hennebelle 2010).

We calculated the sonic Mach numbers of the turbulence in the 1°-wide squares sampling the core of each filament in the case of isotropic turbulence ($\delta v_{3D} = 3\delta v_{1D}$). Figure 15 shows that these estimates are consistent with the distribution of Mach numbers obtained in cold CNM clumps by Heiles & Troland (2003) when comparing the H1 line width with the thermal broadening inferred from spin temperatures. Figure 15 indicates that the turbulence is mildly supersonic along both Eridu and Reticulum in the fibre configuration and it is trans-sonic in the warmer shell geometry. The Mach numbers compare well in both clouds and the turbulence appears to be moderately compressible in both of them.

## 5.4. Dispersion in the $B_{sky}$ magnetic field orientation

In order to study the geometrical complexity of the magnetic field lines in the Reticulum and Eridu clouds we used the dust polarisation data recorded by Planck at 353 GHz. We degraded the original 4.8 resolution to an angular resolution of 14(FWHM) in order to increase the signal to noise ratio. The data come from the third public data release (PR3), which is described in the Planck 2018 Release Explanatory Supplement<fnref n="4" />. The Q and U Stokes parameters provided in the Planck database are defined with respect to the Galactic coordinates, but in a direction opposite to the IAU convention along the longitude

axis, so we define the polarisation angle:

$$ \psi_{IAU} = \frac{1}{2} \arctan(-U, Q), \qquad (9) $$

where $\psi$ is zero at the north Galactic pole and it follows the IAU convention by increasing counter-clockwise with increasing longitude. The arctan(sin, cos) is the arc-tangent function that solves the $\pi$ ambiguity taking into account the sign of the cosine.

The nearby Reticulum and Eridu clouds are the dominant dust clouds along their respective lines of sight, so the dust polarisation measurements give us access to the magnetic field in the plane of the sky, $B_{sky}$, in those clouds. In order to quantify the regularity of the $B_{sky}$ direction, we used the polarisation angle dispersion function (Planck Collaboration Int. XIX 2015):

$$ S_{\psi}(x, \delta) = \frac{1}{N} \sum_{i=1}^{N} \psi(x) - \psi(x + \delta_i), \qquad (10) $$

where $\Delta\psi_{xi} = \psi(x) - \psi(x + \delta_i)$ is the difference between the polarisation angle in a given direction $x$ and in a displaced direction $x + \delta_i$. In practice, we calculated $\Delta\psi_{xi}$ from the Stokes parameters as (Planck Collaboration Int. XIX 2015):

$$ \Delta\psi_{xi} = \frac{1}{2} \arctan(Q_i U_x - Q_x U_i, Q_i Q_x + U_i U_x), \qquad (11) $$

where the $x$ and $i$ indices stand for the central and displaced directions, respectively. The average in $S_\psi$ is taken over displaced directions spanning an annulus around the $x$ direction with an inner radius $\delta/2$ and an outer radius $3\delta/2$. We took 15° for the $\delta$ lag.

The $S_\psi$ dispersion values range from 0° for a well ordered $B_{sky}$ magnetic field to $S_\psi = 90°$ for a disordered one. When the polarisation signal is dominated by noise, $S_\psi$ converges to $\pi/12$ ($\sim 52°$). The results on the $B_{sky}$ field orientation and on its dispersion are displayed as field lines in Fig. 1, as $\psi$ angle maps in Fig. 16, and as $S_\psi$ angle dispersion maps in Fig. 17.

Atomic filamentary clouds have near-magnetic fields often aligned with the filament axis (Clark et al. 2014; Clark & Hensley 2019) whereas the field orientation turns to perpendicular or to disorder in the dense molecular parts (Soler et al. 2013; Planck Collaboration Int. XXXV 2016). The atomic Eridu and Reticulum clouds follow the expected trend. Figures 1, 16, and 17 jointly show that the $B_{sky}$ magnetic field is well ordered along the axis of both filaments, with a stable direction along nearly their whole length, and with a dispersion lower than 10° on a 1-pc scale over large parts of the filament extent. Moreover, these ordered fields are almost perpendicular to the orientation of the Galactic magnetic field inferred in the solar neighbourhood from dust polarisation in local clouds (Planck Collaboration Int. XLIV 2016) or by the field draping around the heliosphere (Zirnstein et al. 2016).

The Reticulum field appears to be slightly more ordered in its dense part than in Eridu, likely because of the compression that originally formed the compact fibre or shell. The dispersion in field orientation appears to be rather uniform across the whole of Eridu and is a little more contrasted between the filament axis and edges of Reticulum. The parsec-scale data, however, do not suggest that CRs have slowed down in Reticulum because of a much more tangled field. To the contrary, compared to Eridu, they should flow along more ordered fields over a large fraction of the filament length.

<footnote n="4">⁴ http://wiki.cosmos.esa.int/planck-legacy-archive</footnote>
<figure><img src="image_16.png" /><figcaption>Fig. 15. Sonic Mach numbers of turbulent motions in the 1'-wide squares chosen along the Eridu (green) and Reticulum (magenta) filaments, in the fibre (top) and shell (bottom) geometries, assuming isotropic turbulence and using the $S_v$ dispersion function (dark filled symbols) or the $\sigma_v$ line width (light open symbols) to trace the one-dimensional turbulent velocity dispersion. The data points and error bars respectively give the mean value in each ~5 pc-wide square and the rms dispersion at 0.5–0.7 pc scale inside each square. The histogram shows the Mach number distribution in cold CNM clumps from Heiles & Troland (2003).</figcaption></figure>

Figure 14 shows that the line width and velocity dispersion values decrease from the cloud edges to the cores of both clouds, as expected from CNM regions. It also indicates that the dispersion in the velocity field of the Eridu cloud is about twice larger than in the more compact Reticulum filament (both in $S_v$ and $\sigma_v$). This is expected for a more diffuse cirrus (Audit & Hennebelle 2010).

We calculated the sonic Mach numbers of the turbulence in the 1'-wide squares sampling the core of each filament in the case of isotropic turbulence ($\delta v_{3D} = 3\delta v_{1D}$). Figure 15 shows that these estimates are consistent with the distribution of Mach numbers obtained in cold CNM clumps by Heiles & Troland (2003) when comparing the H1 line width with the thermal broadening inferred from spin temperatures. Figure 15 indicates that the turbulence is mildly supersonic along both Eridu and Reticulum in the fibre configuration and it is trans-sonic in the warmer shell geometry. The Mach numbers compare well in both clouds and the turbulence appears to be moderately compressible in both of them.

## 5.4. Dispersion in the $B_{sky}$ magnetic field orientation

In order to study the geometrical complexity of the magnetic field lines in the Reticulum and Eridu clouds we used the dust polarisation data recorded by Planck at 353 GHz. We degraded the original 4.8 resolution to an angular resolution of 14(FWHM) in order to increase the signal to noise ratio. The data come from the third public data release (PR3), which is described in the Planck 2018 Release Explanatory Supplement<fnref n="4" />. The Q and U Stokes parameters provided in the Planck database are defined with respect to the Galactic coordinates, but in a direction opposite to the IAU convention along the longitude

axis, so we define the polarisation angle:

$$ \psi_{IAU} = \frac{1}{2} \arctan(-U, Q), \qquad (9) $$

where $\psi$ is zero at the north Galactic pole and it follows the IAU convention by increasing counter-clockwise with increasing longitude. The arctan(sin, cos) is the arc-tangent function that solves the $\pi$ ambiguity taking into account the sign of the cosine.

The nearby Reticulum and Eridu clouds are the dominant dust clouds along their respective lines of sight, so the dust polarisation measurements give us access to the magnetic field in the plane of the sky, $B_{sky}$, in those clouds. In order to quantify the regularity of the $B_{sky}$ direction, we used the polarisation angle dispersion function (Planck Collaboration Int. XIX 2015):

$$ S_{\psi}(x, \delta) = \frac{1}{N} \sum_{i=1}^{N} \psi(x) - \psi(x + \delta_i), \qquad (10) $$

where $\Delta\psi_{xi} = \psi(x) - \psi(x + \delta_i)$ is the difference between the polarisation angle in a given direction $x$ and in a displaced direction $x + \delta_i$. In practice, we calculated $\Delta\psi_{xi}$ from the Stokes parameters as (Planck Collaboration Int. XIX 2015):

$$ \Delta\psi_{xi} = \frac{1}{2} \arctan(Q_i U_x - Q_x U_i, Q_i Q_x + U_i U_x), \qquad (11) $$

where the $x$ and $i$ indices stand for the central and displaced directions, respectively. The average in $S_\psi$ is taken over displaced directions spanning an annulus around the $x$ direction with an inner radius $\delta/2$ and an outer radius $3\delta/2$. We took 15' for the $\delta$ lag.

The $S_\psi$ dispersion values range from 0' for a well ordered $B_{sky}$ magnetic field to $S_\psi = 90'$ for a disordered one. When the polarisation signal is dominated by noise, $S_\psi$ converges to $\pi/12$ ($\sim 52'$). The results on the $B_{sky}$ field orientation and on its dispersion are displayed as field lines in Fig. 1, as $\psi$ angle maps in Fig. 16, and as $\pi$ angle dispersion maps in Fig. 17.

Atomic filamentary clouds have near-magnetic fields often aligned with the filament axis (Clark et al. 2014; Clark & Hensley 2019) whereas the field orientation turns to perpendicular or to disorder in the dense molecular parts (Soler et al. 2013; Planck Collaboration Int. XXXV 2016). The atomic Eridu and Reticulum clouds follow the expected trend. Figures 1, 16, and 17 jointly show that the $B_{sky}$ magnetic field is well ordered along the axis of both filaments, with a stable direction along nearly their whole length, and with a dispersion lower than 10' on a 1-pc scale over large parts of the filament extent. Moreover, these ordered fields are almost perpendicular to the orientation of the Galactic magnetic field inferred in the solar neighbourhood from dust polarisation in local clouds (Planck Collaboration Int. XLIV 2016) or by the field draping around the heliosphere (Zirnstein et al. 2016).

The Reticulum field appears to be slightly more ordered in its dense part than in Eridu, likely because of the compression that originally formed the compact fibre or shell. The dispersion in field orientation appears to be rather uniform across the whole of Eridu and is a little more contrasted between the filament axis and edges of Reticulum. The parsec-scale data, however, do not suggest that CRs have slowed down in Reticulum because of a much more tangled field. To the contrary, compared to Eridu, they should flow along more ordered fields over a large fraction of the filament length.

<footnote n="4">⁴ http://wiki.cosmos.esa.int/planck-legacy-archive</footnote>
<figure><img src="image_18.png" /><figcaption>**Fig. 16.** Map of the polarisation angle $\psi$ at 353 GHz with 14' resolution in the Eridu (top) and Reticulum (bottom) clouds. The angle gives the orientation of the $B_{sky}$ field. The black contours delineate the filaments at $N_{HI}$ column densities of $3 \times 10^{20}$ cm$^{-2}$ and $1.2 \times 10^{20}$ cm$^{-2}$, respectively.</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>**Fig. 17.** Map in the Eridu (top) and Reticulum (bottom) clouds of the polarisation angle dispersion function $S_{\psi}$ at 353 GHz with 14' resolution and $\delta=15'$ lag. The black contours delineate the filaments at $N_{HI}$ column densities of $3 \times 10^{20}$ cm$^{-2}$ and $1.2 \times 10^{20}$ cm$^{-2}$, respectively.</figcaption></figure>

Because of the lateral super-diffusion of turbulent magnetic field lines, CRs with large enough pitch angles can bounce back on magnetic mirrors and leap onto adjacent lines (see Fig. 1 of Lazarian & Xu 2021) whereas the low pitch-angle CRs

<figure><img src="image_20.png" /><figcaption>**Fig. 18.** Plane-of-sky magnetic-field strengths in the 1'-wide squares chosen along the Reticulum (magenta) and Eridu (green) filaments, in the fibre (top) and shell (bottom) geometries, using the $S_{\psi}$ dispersion function in the ST (dark filled symbols) or DCF method (light open symbols). The data points and error bars respectively give the mean value in each ~5 pc-wide square and the rms dispersion at 0.5-0.7 pc scale inside each square. The offset distance of the DCF data points has been slightly enlarged for clarity. The grey band shows the mean CNM $B_{sky}$ value from Heiles & Troland (2005) assuming $B_{sky}^2 = \frac{2}{3}B^2$.</figcaption></figure>

(in the loss cone for mirroring) keep flowing on and undergo gyro-resonant pitch-angle scattering. Mirror diffusion efficiently slows down the effective CR speed along the mean field, but only for the small fraction of particles with large pitch angles, so the overall reduction in $K_{\parallel}$ is small. It is typically less than a factor of 2 and fades out for CR energies below 100 GeV (Lazarian & Xu 2021). The fact that we find comparable dispersion levels in $B_{sky}$ orientation in the two clouds does not suggest widely different levels of mirror diffusion, unless the field gets much more tangled in the more numerous CNM clumps of Reticulum at unresolved scales below one parsec. We are investigating this possibility with interstellar MHD simulations.

## 5.5. Magnetic-field strength

The total magnetic field is composed of a mean component, $B_0$, and a small fluctuating one, $|\delta B| \ll |B_0|$, so that the magnetic energy density is $\frac{B^2}{8\pi} = \frac{1}{8\pi} [B^2 + \delta B^2 + 2\delta B \cdot B_0]$. To infer the mean field strength from dust polarisation, the Davis-Chandrasekhar-Fermi (DCF) method assumes that isotropic turbulent motions initiate incompressible, transverse, Alfvén waves ($\delta B \cdot B_0 = 0$) in the frozen field (Davis 1951; Chandrasekhar & Fermi 1953). They also assume that the kinetic energy of the turbulence is equal to the fluctuating magnetic energy density, and that the dispersion in polarisation angle $\psi$ traces the dispersion in magnetic-field orientation. From the equipartition condition, $\frac{1}{2} \rho \dot{v}^2 = \frac{\delta B^2}{8\pi}$ and from $S \approx \delta B/B$ one obtains:

$$B_{sky}^{\text{DCF}} = f \cdot \frac{4\pi\rho}{S_{\psi}} \frac{\delta v}{S_{\psi}}, \qquad (12)$$

where the f correction factor $v_{\parallel}$ ries with authors between 0.3 and 0.7, for instance $f = 1/3$ for isotropic turbulence (Chandrasekhar & Fermi 1953). Other values are reviewed by Skalidis & Tassis (2021, Sect. 2.2.3). We adopted their $f = 0.5$ recommendation for this work.
<figure><img src="image_17.png" /><figcaption>Fig. 16. Map of the polarisation angle $\psi$ at 353 GHz with 14' resolution in the Eridu (top) and Reticulum (bottom) clouds. The angle gives the orientation of the $B_{sky}$ field. The black contours delineate the filaments at $N_{HI}$ column densities of $3 \times 10^{20}$ cm$^{-2}$ and $1.2 \times 10^{20}$ cm$^{-2}$, respectively.</figcaption></figure>

<figure><img src="image_18.png" /><figcaption>Fig. 17. Map in the Eridu (top) and Reticulum (bottom) clouds of the polarisation angle dispersion function $S_{\psi}$ at 353 GHz with 14' resolution and $\delta=15'$ lag. The black contours delineate the filaments at $N_{HI}$ column densities of $3 \times 10^{20}$ cm$^{-2}$ and $1.2 \times 10^{20}$ cm$^{-2}$, respectively.</figcaption></figure>

Because of the lateral super-diffusion of turbulent magnetic field lines, CRs with large enough pitch angles can bounce back on magnetic mirrors and leap onto adjacent lines (see Fig. 1 of Lazarian & Xu 2021) whereas the low pitch-angle CRs

<figure><img src="image_19.png" /><figcaption>Fig. 18. Plane-of-sky magnetic-field strengths in the 1'-wide squares chosen along the Reticulum (magenta) and Eridu (green) filaments, in the fibre (top) and shell (bottom) geometries, using the $S_{\psi}$ dispersion function in the ST (dark filled symbols) or DCF method (light open symbols). The data points and error bars respectively give the mean value in each ~5 pc-wide square and the rms dispersion at 0.5-0.7 pc scale inside each square. The offset distance of the DCF data points has been slightly enlarged for clarity. The grey band shows the mean CNM $B_{sky}$ value from Heiles & Troland (2005) assuming $B_{sky}^2 = \frac{2}{3}B^2$.</figcaption></figure>

(in the loss cone for mirroring) keep flowing on and undergo gyro-resonant pitch-angle scattering. Mirror diffusion efficiently slows down the effective CR speed along the mean field, but only for the small fraction of particles with large pitch angles, so the overall reduction in $K_{\parallel}$ is small. It is typically less than a factor of 2 and fades out for CR energies below 100 GeV (Lazarian & Xu 2021). The fact that we find comparable dispersion levels in $B_{sky}$ orientation in the two clouds does not suggest widely different levels of mirror diffusion, unless the field gets much more tangled in the more numerous CNM clumps of Reticulum at unresolved scales below one parsec. We are investigating this possibility with interstellar MHD simulations.

## 5.5. Magnetic-field strength

The total magnetic field is composed of a mean component, $B_0$, and a small fluctuating one, $|\delta B| \ll |B_0|$, so that the magnetic energy density is $\frac{B^2}{8\pi} = \frac{1}{8\pi} [B^2 + \delta B^2 + 2\delta B \cdot B_0]$. To infer the mean field strength from dust polarisation, the Davis-Chandrasekhar-Fermi (DCF) method assumes that isotropic turbulent motions initiate incompressible, transverse, Alfvén waves ($\delta B \cdot B_0 = 0$) in the frozen field (Davis 1951; Chandrasekhar & Fermi 1953). They also assume that the kinetic energy of the turbulence is equal to the fluctuating magnetic energy density, and that the dispersion in polarisation angle $\psi$ traces the dispersion in magnetic-field orientation. From the equipartition condition, $\frac{1}{2} \rho \dot{v}^2 = \frac{\delta B^2}{8\pi}$ and from $S \approx \delta B/B$ one obtains:

$$B_{sky}^{\text{DCF}} = f \cdot \frac{4\pi\rho}{S_{\psi}} \frac{\delta v}{S_{\psi}}, \qquad (12)$$

where the f correction factor $v_{\psi}/\sqrt{\rho}$ with authors between 0.3 and 0.7, for instance $f = 1/3$ for isotropic turbulence (Chandrasekhar & Fermi 1953). Other values are reviewed by Skalidis & Tassis (2021, Sect. 2.2.3). We adopted their $f = 0.5$ recommendation for this work.
<figure><img src="image_21.png" /><figcaption>Fig. 19. Histograms of the gas magnetisation levels in the 1°-wide squares chosen along the Eridu (green) and Reticulum (magenta) filaments, for each geometry, and using the $S_v$ velocity dispersion (filled bars) or the $\sigma_v$ line width (markers) in the ST formulation to infer $B_{sky}$ and the total magnetic-field strength.</figcaption></figure>

The DCF method tends to overestimate the actual field strength because of line-of-sight and beam averaging in the data and because the method overlooks the strong anisotropy in ISM turbulence and the role of compressible modes in the cross-term $\delta B \cdot B_0 \neq 0$ (Skalidis & Tassis 2021). To include the compressible modes, these authors propose an alternative formulation that agrees within 20–30% with MHD simulation data:

$$B_{\text{sky}}^{\text{ST}} = \frac{2\pi\rho}{S_{\psi}} \frac{\delta v}{S_{\psi}}. \qquad (13)$$

We applied both methods to calculate the plane-of-sky magnetic-field strengths in the 1°-wide squares chosen along the Reticulum and Eridu filaments, with the velocity dispersion function $S_v$ as a proxy for the turbulent velocity $\delta v$. The results are shown in Fig. 18 and in Tables B.1 and B.2. Following the discussion by Skalidis & Tassis (2021), the $B_{sky}^{ST}$ and $B_{sky}^{DCF}$ estimates likely bracket the true $B_{sky}$ value, the latter being about 40% (50%) larger than the former in the case of Reticulum (Eridu). We also calculated upper limits to $B_{sky}$ using the $\sigma_v$ line width as an upper limit to the turbulent velocity. The results are provided in Tables B.1 and B.2. They are about three times larger than the $S_v$-based estimates for each method.

Figure 18 shows that the $B_{sky}$ field strength is rather uniform along the core of each cloud and that it is about one and a half times larger along Eridu than along Reticulum: around 3.1 (1.7) µG for Eridu and 1.9 (0.9) µG for Reticulum in the fibre (shell) configuration if we use the $S_v$ velocity dispersion to trace the turbulent motions. Three times larger upper limits to the fields are found when using the $\sigma_v$ line width instead. The difference between Eridu and Reticulum does not depend much on the cloud geometry, but one should keep in mind that the magnetic-field component along the line of sight may significantly differ in both clouds. Figure 18 also illustrates that the $B_{sky}$ strengths in Eridu and Reticulum compare reasonably with the median total field $B_{tot} = 6 \pm 1.8$ µG inferred from Zeeman splitting of Hα lines and Monte Carlo simulations (Heiles & Troland 2005) if one assumes $B_{sky}^{2} = (2/3) B^{2}$ for a random field orientation with respect to the line of sight (Heiles & Troland 2005).

We cannot estimate whether the MHD turbulence is super-Alfvénic or not since both $B_{sky}$ derivations assume equipartition between the fluctuating kinetic and magnetic energy densities, whether the latter is dominated by $\delta B^2/8\pi$ in the case of incompressible turbulence (DCF) or by $\delta BB/4\pi$ when compressible modes are present (ST). In both cases, the Alfvén Mach numbers will be close to one. The large-scale order we observe in the $B_{sky}$ orientation in both clouds suggests sub- or trans-Alfvénic turbulence as often found for the atomic gas phases in interstellar MHD simulations.

We can estimate the gas magnetisation level from the squared ratio of the sound to Alfvén speed, $c_s/v^2$ (Soler et al. 2013), which is close to the gas to magnetic pressure ratio. To derive the Alfvén velocity of the total gas, we assumed $B^2 = (3/2)B_{sky}^2$ and the ST formulation for $B_{sky}$. Figure 19 shows the normalised distributions of the ratios found in the 1°-wide squares sampling the gas along the core of both filaments. If we use the $B_{sky}(S_v)$ estimate, the gas and magnetic pressures compare well across both clouds in their fibre geometry and the gas pressure slightly exceeds the magnetic one in the shell configuration. The magnetisation is skewed to values below one (i.e. to dominant magnetic pressure) when using the $B_{sky}(\sigma_v)$ upper limits to the field, but not severely, so equivalent gas and magnetic pressures may prevail over large parts of the clouds. Therefore, whether in terms of $B_{sky}$ orientation dispersion or in terms of magnetic stiffness against turbulent gas motions, we find no important difference between the two clouds at a parsec scale that would suggest more magnetic-field tangling and slower CR diffusion in one of them.

## 6. Conclusions

We compared the CR flux in the atomic phase of two nearby clouds that have many properties in common. They are both elongated filaments, which cross the same hollow sector of the Local interstellar Valley at comparable distances from the Sun. Both are threaded by ordered magnetic fields that are mostly aligned with the long axis of the filament and that are largely
<figure><img src="image_20.png" /><figcaption>Fig. 19. Histograms of the gas magnetisation levels in the 1°-wide squares chosen along the Eridu (green) and Reticulum (magenta) filaments, for each geometry, and using the $S_v$ velocity dispersion (filled bars) or the $\sigma_v$ line width (markers) in the ST formulation to infer $B_{sky}$ and the total magnetic-field strength.</figcaption></figure>

The DCF method tends to overestimate the actual field strength because of line-of-sight and beam averaging in the data and because the method overlooks the strong anisotropy in ISM turbulence and the role of compressible modes in the cross-term $\delta B \cdot B_0 \neq 0$ (Skalidis & Tassis 2021). To include the compressible modes, these authors propose an alternative formulation that agrees within 20–30% with MHD simulation data:

$$B_{\text{sky}}^{\text{ST}} = \frac{2\pi\rho}{S_{\psi}} \frac{\delta v}{S_{\psi}}. \qquad (13)$$

We applied both methods to calculate the plane-of-sky magnetic-field strengths in the 1°-wide squares chosen along the Reticulum and Eridu filaments, with the velocity dispersion function $S_v$ as a proxy for the turbulent velocity $\delta v$. The results are shown in Fig. 18 and in Tables B.1 and B.2. Following the discussion by Skalidis & Tassis (2021), the $B_{sky}^{ST}$ and $B_{sky}^{DCF}$ estimates likely bracket the true $B_{sky}$ value, the latter being about 40% (50%) larger than the former in the case of Reticulum (Eridu). We also calculated upper limits to $B_{sky}$ using the $\sigma_v$ line width as an upper limit to the turbulent velocity. The results are provided in Tables B.1 and B.2. They are about three times larger than the $S_v$-based estimates for each method.

Figure 18 shows that the $B_{sky}$ field strength is rather uniform along the core of each cloud and that it is about one and a half times larger along Eridu than along Reticulum: around 3.1 (1.7) µG for Eridu and 1.9 (0.9) µG for Reticulum in the fibre (shell) configuration if we use the $S_v$ velocity dispersion to trace the turbulent motions. Three times larger upper limits to the fields are found when using the $\sigma_v$ line width instead. The difference between Eridu and Reticulum does not depend much on the cloud geometry, but one should keep in mind that the magnetic-field component along the line of sight may significantly differ in both clouds. Figure 18 also illustrates that the $B_{sky}$ strengths in Eridu and Reticulum compare reasonably with the median total field $B_{tot} = 6 \pm 1.8$ µG inferred from Zeeman splitting of Hα lines and Monte Carlo simulations (Heiles & Troland 2005) if one assumes $B_{sky}^{2} = (2/3) B^{2}$ for a random field orientation with respect to the line of sight (Heiles & Troland 2005).

We cannot estimate whether the MHD turbulence is super-Alfvénic or not since both $B_{sky}$ derivations assume equipartition between the fluctuating kinetic and magnetic energy densities, whether the latter is dominated by $\delta B^2/8\pi$ in the case of incompressible turbulence (DCF) or by $\delta BB/4\pi$ when compressible modes are present (ST). In both cases, the Alfvén Mach numbers will be close to one. The large-scale order we observe in the $B_{sky}$ orientation in both clouds suggests sub- or trans-Alfvénic turbulence as often found for the atomic gas phases in interstellar MHD simulations.

We can estimate the gas magnetisation level from the squared ratio of the sound to Alfvén speed, $c_s/V^2$ (Soler et al. 2013), which is close to the gas to magnetic pressure ratio. To derive the Alfvén velocity of the total gas, we assumed $B^2 = (3/2)B_{sky}^2$ and the ST formulation for $B_{sky}$. Figure 19 shows the normalised distributions of the ratios found in the 1°-wide squares sampling the gas along the core of both filaments. If we use the $B_{sky}(S_v)$ estimate, the gas and magnetic pressures compare well across both clouds in their fibre geometry and the gas pressure slightly exceeds the magnetic one in the shell configuration. The magnetisation is skewed to values below one (i.e. to dominant magnetic pressure) when using the $B_{sky}(\sigma_v)$ upper limits to the field, but not severely, so equivalent gas and magnetic pressures may prevail over large parts of the clouds. Therefore, whether in terms of $B_{sky}$ orientation dispersion or in terms of magnetic stiffness against turbulent gas motions, we find no important difference between the two clouds at a parsec scale that would suggest more magnetic-field tangling and slower CR diffusion in one of them.

## 6. Conclusions

We compared the CR flux in the atomic phase of two nearby clouds that have many properties in common. They are both elongated filaments, which cross the same hollow sector of the Local interstellar Valley at comparable distances from the Sun. Both are threaded by ordered magnetic fields that are mostly aligned with the long axis of the filament and that are largely
inclined with respect to the Galactic plane and to the large-scale magnetic-field direction in the local ISM. Both clouds have magnetic pressures of the same magnitude as the thermal ones. Nevertheless, we find that the γ-ray emissivity per gas nucleon in one cloud, Reticulum, is $1.57 \pm 0.09$ times larger than in the other cloud, Eridu, over the whole energy range of analysis, from 0.16 to 63 GeV. The difference cannot be attributed to uncertainties in gas mass because the correlation between the dust and γ-ray data allows an efficient detection of any additional gas that is not traced by H I line emission (as we find in Reticulum) and the γ-ray contribution of this additional gas is included in the model. Furthermore, we cannot reduce the CR flux ratio below $1.32 \pm 0.08$ when allowing for strong differences in H I spin temperature between the two clouds (at the expense of a poorer γ-ray fit).

Beyond characterising the global magnetic orientation and average field strength of a few microGauss in the clouds, we studied turbulence characteristics at a parsec scale that might hint at different environmental conditions for CR diffusion, namely the dispersion in gas velocity and in magnetic-field orientation, differences in the sonic Mach numbers (compressibility) of the turbulence, and differences in magnetisation level (magnetic stiffness against gas turbulent motions). We find that the two clouds compare well in all these aspects. The only notable difference is that the core of the Reticulum cloud has collapsed (or has been compressed) to larger densities typical of the CNM and DNM phases, whereas a large fraction of the more diffuse gas in Eridu may be transitioning through the unstable, lukewarm phase. MHD simulations will help in studying the level of magnetic-field tangling in the different atomic and DNM gas phases, or in shell compression behind a slow shock, in order to search for changes that may impact the collective drift of the CR population, such as changes in the turbulent energy fraction of fast compressive modes or in the density of magnetic mirrors (Lazarian et al. 2023).

Alternatively, the 1–500 GeV/n CRs we probed in γ rays may primarily scatter off MHD waves they have excited through streaming instability. If so, one expects a stronger damping in the denser Reticulum cloud, but we show that the difference in expected steady-state diffusion coefficient between the two clouds can be small (less than a factor of two) because the increase in diffusion length with the ambient gas density is moderated by the corresponding drop in gas temperature and ionisation fraction. We plan to further investigate the impact of charged dust grains, which should reduce the damping rates along the core of the two clouds (Hennebelle & Lebreuilly 2023). However, we only probed the cloud properties at a parsec scale, whereas the growth and damping of the relevant MHD waves occur at wavelengths of a few micro-parsecs where the interstellar gas properties cannot be explored. Moreover, in the self-confinement mode, the wave growth rate in each cloud relates to the degree of anisotropy in CR pitch angle or, equivalently, to the steepness of the local CR spatial density gradient (Zweibel 2013), which could differ in the vicinity of the clouds or in bottlenecks along their edge (different $l_{CR}$).

The damping rates we estimate suggest that the magnetic perturbations propagating along the field lines should dissipate over typical lengths of tens to hundreds of micro-parsecs. ($l \approx 10^{10}$ in $\Gamma$-ray, $\Gamma$-ray, Spangler et al. 2011). The median of their distribution is $A$ in

twice as large in Eridu as in Reticulum. The dissipation lengths correspond to proton gyro-radii for energies in the~10 GeV to ~1 TeV band, with a median energy three times larger in Eridu than in Reticulum. The waves must be constantly replenished or we would observe a spectral change in γ-ray emissivity in

the corresponding ~1 to ~100 GeV band as the higher-energy particles would leak out more rapidly. The fact that we find a change in CR flux but not in spectrum between the two clouds is therefore an important constraint. The same overall spectrum has been observed in several clouds in the local ISM (Grenier et al. 2015), but with improved photon statistics the Fermi-LAT data can now serve to precisely test the uniformity of this spectrum in different types of atomic clouds and independently of the CR flux pervading them.

**Acknowledgements.** We thank Jeong-Gyu Kim for sharing his interstellar model for the solar neighbourhood conditions, Lucia Armillotta for sharing information from her work, the members of the Interstellar Institute (I2) CNRS International Research Network for enlightening discussions, and the referee for helping us clarify several points in the discussion. We acknowledge the financial support by the LabEx UnivEarthS (ANR-10-LABX-0023 and ANR-18-IDEX-0001) for this work. The Fermi LAT Collaboration acknowledges generous ongoing support from a number of agencies and institutes that have supported both the development and the operation of the LAT as well as scientific data analysis. These include the National Aeronautics and Space Administration and the Department of Energy in the United States, the Commissariat à l'Energie Atomique and the Centre National de la Recherche Scientifique / Institut National de Physique Nucléaire et de Physique des Particules in France, the Agenzia Spaziale Italiana and the Istituto Nazionale di Fisica Nucleare in Italy, the Ministry of Education, Culture, Sports, Science and Technology (MEXT), High Energy Accelerator Research Organization (KEK) and Japan Aerospace Exploration Agency (JAXA) in Japan, and the K.A. Wallenberg Foundation, the Swedish Research Council and the Swedish National Space Board in Sweden. Additional support for science analysis during the operations phase is gratefully acknowledged from the Istituto Nazionale di Astrofisica in Italy and the Centre National d'Etudes Spatiales in France.

## References

Abdo, A. A., Ackermann, M., Ajello, M., et al. 2010, ApJ, 710, 133

Abdollahi, S., Acero, F., Baldini, L., et al. 2022, ApJS, 260, 53

Acero, F., Ackermann, M., Ajello, M., et al. 2016, ApJS, 223, 26

Ackermann, M., Ajello, M., Allafort, A., et al. 2012a, ApJ, 756, 4

Ackermann, M., Ajello, M., Atwood, W. B., et al. 2012b, ApJ, 750, 3

Ackermann, M., Ajello, M., Baldini, L., et al. 2013, ApJ, 772, 154

Armillotta, L., Ostriker, E. C., & Jiang, Y.-F. 2021, ApJ, 922, 11

Audit, E., & Hennebelle, P. 2010, A&A, 511, A76

Bolatto, A. D., Wolfire, M., & Leroy, A. K. 2013, ARA&A, 51, 207

Boschini, M. J., Della Torre, S., Gervasi, M., et al. 2020, ApJS, 250, 27

Briel, P., Burnett, T. H., Digel, S. W., et al. 2018, arXiv e-prints, [arXiv:1810.11394]

Bustard, C., & Zweibel, E. G. 2021, ApJ, 913, 106

Casandjian, J.-M. 2015, ApJ, 806, 240

Casandjian, J.-M., Ballet, J., Grenier, I., & Remy, Q. 2022, ApJ, 940, 116

Chandran, B. D. G. 2000, Phys. Rev. Lett., 85, 4656

Chandrasekhar, S., & Fermi, E. 1953, ApJ, 118, 113

Clark, S. E., & Hensley, B. S. 2019, ApJ, 887, 136

Clark, S. E., Peek, J. E. G., & Putman, M. E. 2014, ApJ, 789, 82

Davis, L. 1951, Phys. Rev., 81, 890

De La Torre Luque, P., Mazziotta, M. N., Ferrari, A., et al. 2022, J. Cosmology & Astropart. Phys., 2022, 008

Dickey, J. M., Mebold, U., Stanimirovic, S., & Staveley-Smith, L. 2000, ApJ, 536, 756

Draine, B. T. 2004, in The Cold Universe, 32, 213

Draine, B. T. 2011, Physics of the Interstellar and Intergalactic Medium (Princeton University Press)

Finkbeiner, D. P. 2003, ApJS, 146, 407

Gaensler, B. M., Madsen, G. J., Chatterjee, S., & Mao, S. A. 2008, PASA, 25, 184

Gaustad, J. E., McCullough, P. R., Rosing, W., & Van Buren, D. 2001, PASP, 113, 1326

Goldreich, P., & Sridhar, S. 1995, ApJ, 438, 763

Gong, M., Ostriker, E. C., & Kim, C.-G. 2018, ApJ, 858, 16

Grenier, I. A., Casandjian, J.-M., & Terrier, R. 2005, Science, 307, 1292

Grenier, I. A., Black, J. H., & Strong, A. W. 2015, ARA&A, 53, 199

Heiles, C., & Troland, T. H. 2003, ApJ, 586, 1067

Heiles, C., & Troland, T. H. 2005, ApJ, 624, 773

Hennebelle, P., & Lebreuilly, U. 2023, A&A, 674, A149

Hopkins, P. F., Squire, J., Chan, T. K., et al. 2021, MNRAS, 501, 4184

Hu, Y., Lazarian, A., & Xu, S. 2022, MNRAS, 512, 2111