Free Access
Issue
A&A
Volume 601, May 2017
Article Number A78
Number of page(s) 27
Section Interstellar and circumstellar matter
DOI https://doi.org/10.1051/0004-6361/201629632
Published online 08 May 2017

© ESO, 2017

1. Introduction

The structural, dynamical and thermodynamical evolution of interstellar clouds is often probed with observations of the H i 21-cm line for the atomic gas (Kalberla et al. 2005, 2010; Peek et al. 2011) and of CO rotational emission lines for the molecular gas (Heyer & Dame 2015). However, measuring gas column densities in the dense media, both atomic and molecular, is hampered by the opacity of the interstellar medium (ISM) to these radiations and we have no efficient means to correct for line opacities in the absence of absorption data against background point sources. Complementary information can be gained from the sub-millimeter and infrared thermal emission of large dust grains mixed with the gas (Planck Collaboration XI 2014), and from γ rays produced in interactions of high-energy cosmic rays (CRs) with interstellar gas nucleons (Lebrun et al. 1982). Dust grains and CRs trace all chemical and thermodynamical forms of the gas to large depths into the clouds (Planck Collaboration Int. XXVIII 2015), but they trace only integrals of the gas column densities, NH, along sight lines, and they bear no kinematic information. To use them, we rely on critical assumptions, namely a uniform dust-to-gas mass ratio and uniform grain emissivity across a cloud and, likewise, a uniform CR flux through a cloud. Those assumptions still need to be tested in a variety of clouds and across complex phase changes, especially on consideration that dust grains and their radiation properties vary from the diffuse ISM to molecular clouds (Stepnik et al. 2003; Martin et al. 2012; Roy et al. 2013a; Ysard et al. 2013; Flagey et al. 2009; Planck Collaboration XI 2014; Planck Collaboration Int. XVII 2014; Planck Collaboration Int. XXVIII 2015). One can also test the uniformity of the γ-ray emissivity spectrum of the gas in well-sampled nearby clouds to verify the smooth penetration of CRs with energies above ~ 1 GeV (Skilling & Strong 1976; Cesarsky & Volk 1978; Morlino & Gabici 2015).

The recent progress in H i, CO, dust, and γ-ray observations provides excellent opportunities to confront these observations at a resolution of a few parsecs inside nearby clouds to probe biases in the different tracers and to search for missing gas in the census. This comparison can shed light on the highly opaque atomic gas forming the cold neutral medium (CNM), where the H ibrightness temperature probes the excitation (spin) temperature of the gas rather than the column density (Murray et al. 2015, and references therein). One can also probe how the 12CO, J = 1 → 0, emission rapidly saturates at high column density and how it fades out in the diffuse molecular gas where CO lines are too weakly excited and/or CO molecules are too efficiently photodissociated to map the corresponding molecular hydrogen (van Dishoeck & Black 1988; Wolfire et al. 2010). Large H iopacities and the lack of CO emission both conspire to hide large gas masses at the H–H2 interface. This dark neutral medium (DNM) can be jointly revealed by cosmic rays and dust (Grenier et al. 2005).

A more precise census of the gas column densities is key to establishing and quantifying variations in dust opacity (optical depth per gas nucleon) across cloud structures. Mapping those variations as a function of grain temperature, T, and opacity spectral index, β gives important diagnostics for grain evolution through mantle accretion, coagulation, and ice coating (Köhler et al. 2015).

With an approach similar to the multi-tracer analysis of the Chamaeleon complex (Planck Collaboration Int. XXVIII 2015), we present an analysis of the anticentre region, which exhibits more massive clouds, with higher levels of star-formation activity, so we can test the properties of the gas tracers to larger column densities and molecular fractions. The analysed region spans a few hundred parsecs in distance in the local ISM (see below). It encompasses the well-known Taurus, Auriga, and Perseus molecular clouds (Ungerechts & Thaddeus 1987; Lombardi et al. 2010; Schnee et al. 2008), the massive complex associated with the California nebula (Lada et al. 2009), and a nearby, diffuse, hardly molecular cloud, which we refer to as Cetus.

Our goal is to test the dust and γ-ray emission as tracers of the total gas against linear combinations of the emission line intensities recorded in the H iand CO observations, a new free-free emission map of the ionised gas, plus the DNM gas which escapes the radio and millimeter spectroscopic surveys at the atomic-molecular interface. We have separated the H iand CO clouds in space (i.e., direction, velocity, and distance) to study their individual contributions to the total dust and γ-ray signals. We have used the model results to derive scaling factors such as the CO-to-H2 conversion factor, XCO, in each cloud, and the average dust opacity in the different gas phases for each cloud. We have then used the multi-wavelength data to follow dust evolution with increasing NH in several molecular clouds with masses spanning nearly two orders of magnitude. A companion paper will present the results on the transitions between gas phases, and on the mass fractions and contribution of CO-dark H2 to the molecular phase in individual regions.

The paper is structured as follows: the data are presented in Sect. 2 and the models are described in Sect. 3. The model results and the detection of gas in the DNM, CO-saturated, and ionised phases are detailed in Sect. 4. The CR content of the clouds, the XCO ratios, and the dust opacities are discussed in Sects. 5, 6, and 7, respectively. Additionally, we present the H iand CO component separation into nearby cloud complexes in Appendix A, the derivation of the free-free intensities in Appendix B, the dependence of the goodness of fit on the H ispin temperature in Appendix C, the best-fit coefficients of the dust and γ-ray models in Appendix D, and a table with historical γ-ray measurements of XCO factors in E.1.

thumbnail Fig. 1

Maps of the Cetus (Cet), Main Taurus (TauM), South Taurus (TauS), North Taurus (TauN), Perseus (Per), and California (Cal) clouds, and of the Galactic disc background (Gal) in H icolumn density, NHI, for a spin temperature of 400 K, and in CO intensity, WCO. The NHI and WCO maps are respectively labelled HI and CO. The map labelled ff shows the intensity Iff of free-free emission toward this region. The last subplot shows the coverage of the GALFA (dark grey) and EBHIS (light blue) H i data.

Open with DEXTER

2. Data

The analysis region extends from 139° to 191° in Galactic longitude and from − 56° to − 3° in Galactic latitude. We have masked regions with large column densities from the background Galactic disc, at − 9°<b< − 3° for l> 183° and l< 152° and regions with small column densities at declinations below − 5°where the H imaps have lower angular resolution. All maps in Fig. 1 have been projected onto the same 0.̊125-spaced Cartesian grid as that of the CfA CO survey (Dame et al. 2001). We have checked that the use of the Cartesian projection, instead of equal-area bins, does not significantly bias the results. The largest differences are less than a fifth of the uncertainties found in the jackknife tests discussed in Sect. 4.2.

2.1. HI and CO data

For this analysis we have used the 4-resolution GALFA-HI survey where available (Peek et al. 2011) and we have complemented it with the 10.8-resolution EBHIS survey elsewhere (Winkel et al. 2016). As shown in Fig. 1, the GALFA survey covers 74% of the region. GALFA data-cubes1 were re-sampled into the 0.̊125-spaced Cartesian grid. We have used the narrow-band cubes with their original velocity resolution of 0.18 km s-1 in the local standard of rest (LSR). The EBHIS frequency sampling is coarser, with a velocity resolution of 1.44 km s-1. In the analysed region covered by both survey we have found a tight correlation in column densities between the two: (1)In order to trace the molecular gas we have used 12CO J = 1 → 0 observations at 115 GHz from the moment-masked CfA CO survey of the Galactic plane (Dame et al. 2001; Dame & Thaddeus 2004). Most of the survey is based on a 0.̊125-spaced Cartesian grid except for the high-latitude clouds, at b ≲ − 50°, which have been interpolated from 0.̊25 to 0.̊125.

2.2. HI and CO velocity separation

We have decomposed the H iand CO spectra into individual lines and we have used this information to identify six nearby entities that are coherent in velocity, position (in Galactic coordinates), and distance, in addition to the emission from the background Galactic disc. The H ivelocity resolution in the region covered by the EBHIS survey is coarser, but still adequate to separate the main structures because of the limited confusion along sight lines at medium-to-high latitudes. The separation process is three-fold and is described by Planck Collaboration Int. XXVIII (2015) and in Appendix A. The choice of six main entities has been guided by the structure in (l, b, v) density and by distance information obtained from stellar reddening (Schlafly et al. 2014). Details on the (l, b, v) cuts are given in Table A.1. The main structures found in the local ISM, namely Cetus, Main Taurus, South Taurus, North Taurus, Perseus, and California, together with the background Galactic structures are depicted in NHI column densities and in WCO intensities in Fig. 1. The set of H imaps includes all of the H iemission observed in this region except for the high velocity clouds (Wakker et al. 2008). In order to investigate the effect of the unknown H ioptical depth, we have derived all the NHI maps for a sample of uniform spin temperatures (150, 200, 300, 400, 500, 600, 700, and 800 K) and for the optically thin case.

Distance estimates to the local clouds have been compiled from the photometric measurements of PanSTARRS-1 and the detection of reddening fronts toward stellar groups (Schlafly et al. 2014). We have selected the stars along lines of sight intersecting each (and only one) cloud within its WCO contour at 0.5 K km s-1 and we have calculated the average distance in each sample. For the reddening fronts to the stars toward Cetus, Main Taurus, South Taurus, North Taurus, Perseus, and California we find average distances of 190 ± 30 pc, 140 ± 30 pc, 160 ± 10 pc, 190 ± 50 pc, 270 ± 20 pc, and 410 ± 20 pc, respectively. They compare reasonably with the geometry inferred from the modelling of the far-UV continuum emission (Lim et al. 2013). The North Taurus component is dominated by the Auriga cloud near the northern edge of our analysis region; it also includes MBM 12 molecular cloud from Magnani et al. (1985) catalogue. The South Taurus component includes MBM 6, MBM 12 and MBM 18.

thumbnail Fig. 2

Maps of dust optical depth at 353 GHz (top left), radiance (top right), spectral β index (bottom left), and colour temperature (lower right), from Planck Collaboration XI (2014).

Open with DEXTER

2.3. Dust data

In order to trace the dust column density, we have used the optical depth, τ353, obtained at 353 GHz by Planck Collaboration XI (2014). They have modelled the spectral energy distributions of the thermal emission of the large grains with modified black-body curves in specific intensity, Iν = τν0Bν(Td)(ν/ν0)β. The optical depth τν0 at frequency ν0, the dust colour temperature Td, and the spectral index β were fitted to the data from Planck and the Infrared Astronomical Satellite (IRAS) across the sky at 30 resolution. The fits were then repeated at 5 resolution by fixing β as obtained in the first step. This procedure limited the noise impact on the Tdβ degeneracy in faint regions. The modelled intensities were checked to be consistent with the data at all frequencies (see Fig. 11 in Planck Collaboration XI 2014). The contamination from CO line emission in the 353 GHz filter band amounts to a few per cent of the total signal. It was not removed so as to avoid increasing the noise in directions away from CO clouds.

The optical depth, τν, and opacity, σν, of the dust at frequency ν follow the relations: (2)We have also calculated the total power radiated by the dust, also referred to as radiance: (3)Figure 2 shows the dust properties thus measured across the analysis region.

The corresponding analysis of the Chamaeleon region (Planck Collaboration Int. XXVIII 2015) compared different dust tracers, namely the optical depth and an estimate of the dust extinction, AVQ, empirically corrected for starlight intensity. The latter distribution best followed the gas column densities in the Chamaeleon region, but we have checked that this is not true in the anticentre region in particular toward regions of dust temperatures exceeding 23 K, so we present the results of our analyses using only the τ353 optical depth. Understanding why the empirical correction used for AVQ is less efficient in the anticentre region will be investigated separately.

2.4. Ionised gas

Ionised gas is visible in Hα emission (Finkbeiner 2003) in the California nebula around l = 160.̊1 and b = − 12.̊4, in the G159.6-18.5 H iiregion, and along the Eridanus loop. The Califonia Nebula, alias NGC 1499 or S 220, is ionised by ξ Persei, an O7.5III star located at a distance of 380 ± 70 pc (van Leeuwen 2007) that has run away from the Perseus OB 2 association (Shinohara & Ishida 1976). The 1.̊2-diameter region G159.6-18.5 has been seen by IRAS and the Wide-Field Infrared Survey Explorer (WISE) (Andersson et al. 2000; Anderson et al. 2014). It is excited by another run-away O9.5V star, HD 278942, at a distance of 190 ± 40 pc (van Leeuwen 2007), which would place it in front rather than behind the Perseus cloud (Bally et al. 2008).

The composite Hα map varies in angular resolution from 6 to 1° across the region of interest and Hα emission can be efficiently absorbed behind dense clouds (see for instance the dusty absorption lane across G159.6-18.5 in Fig. 8 of Bally et al. 2008). We have therefore preferred measures of the free-free intensity at mm wavelengths to trace the ionised gas.

The Planck LFI map at 70 GHz (release 2.01 data, Planck Collaboration I 2016) shows free-free emission in spatial coincidence with the Hα emission. The angular resolution of 14 of LFI better suits our analysis than the 1°-resolution of the free-free map inferred at 22 GHz from the 9-year WMAP data (Bennett et al. 2013). To separate the free-free intensity in the 70 GHz data, we have successively removed the contributions from the Cosmological Microwave Background (CMB), from dust emission, and from point sources. The method we have used to extract the free-free emission is described in Appendix B.

The resulting free-free dominated map is displayed in Fig. 1. It closely resembles the Hα intensity map of the region (Finkbeiner 2003), but for a uniform resolution of 14 and a reduced contrast, by a factor 2 to 8, between the intensities of G159.6-18.5 and NGC 1499. because the Perseus clouds are more transparent to mm waves.

The integral of the free-free intensity along a sight line through ionised gas with temperature Te and electron and ion number densities, ne and ni, scales at 70 GHz as: (4)where the emission measure is expressed in cm-6 pc (Karzas & Latter 1961; Carson 1988). If the observed intensity results from discrete and homogeneous nebulae with comparable uniform electron densities, , Iff roughly scales with the ion column density, (5)so we have used the Iff map as a template for the NHII distribution in the gas model and we discuss the scaling ratio obtained in the fit between NHII and Iff in Sect. 4.4.

2.5. Gamma-ray data

We have used six years of the Pass 8 photon data from Fermi-LAT, the associated instrument response functions (IRFs, P8R2_CLEAN_V6), and the corresponding isotropic spectrum for the extragalactic and residual instrumental backgrounds (Atwood et al. 2013; Acero et al. 2015). We have applied tight rejection criteria (CLEAN class selection, photon arrival directions within <100° of the Earth zenith and in time intervals when the LAT rocking angle was inferior to 52°) in order to reduce the contamination by residual CRs and Earth atmospheric γ rays in the photon data (see Nolan et al. 2012, for details). We have sampled the IRFs, the exposure map, the γ-ray emissivity spectrum, qLIS, of the local interstellar gas (Casandjian 2015), and the spectrum of the isotropic background in 14 energy bins, 0.2 dex in width and centred from 102.3 to 104.9 MeV.

To ensure photon statistics robust enough to follow details in the spatial distributions of the different interstellar components, we have analysed the data in broad and independent energy bands, bounded by 102.6, 102.8, 103.2, 103.6, and 105 MeV. We have also analysed the data in the integrated 102.6–105 MeV band. The low-energy threshold was set to take advantage of the best LAT point-spread function (PSF), which strongly degrades at low energy. The analysis of the data in two lower-energy bands (102.2–102.4 and 102.4–102.6 MeV) proved to be less reliable to separate the gas phases in compact clouds. To optimise the angular resolution, we have selected all detected photons above 1.6 GeV, but only those at lower energy that converted to pairs in the front section of the tracker (Atwood et al. 2009; Ackermann et al. 2012a). Given the broad ranges of the energy bands, we have not corrected the fluxes for the energy resolution, which is better than 14% at these energies (Atwood et al. 2013). To account for the spill-over of emission produced outside the analysis region, but reconstructed inside it, we have modelled point sources and interstellar contributions in a region 4° wider than the analysis region. This choice corresponds to the 99.5% containment radius of the PSF in the lowest energy band. For the H iemission outside the analysis region not covered by EBHIS and GALFA surveys we have used the Leiden/Argentine/Bonn (LAB) survey (Kalberla et al. 2005).

The positions and the flux spectra of the γ-ray sources in the field are provided by the Fermi-LAT Third Source Catalog (Acero et al. 2015). The observed γ-ray emission also includes a contribution from the large-scale Galactic inverse Compton (IC) emission emanating from the interactions of CR electrons with the ISRF. The GALPROP2 parameter file 54-LRYusifovXCO4z6R30-Ts150-mag2 has been used to generate an energy-dependent template of the Galactic IC emission across the analysis region (Ackermann et al. 2012e).

thumbnail Fig. 3

Top: γ-ray counts of gaseous origin recorded in the 0.4–100 GeV energy band in a 0.̊125 pixel grid. γ-ray emissions other than due to cosmic-ray interactions in the gas have been subtracted. The map has been smoothed with a Gaussian kernel of 0.̊14 dispersion for display. Bottom: dust optical depth measured at 353 GHz and displayed at the Fermi-LAT angular resolution for comparison.

Open with DEXTER

3. Models and analyses

To first order, the dust and γ-ray emissions should both trace the total gas column density, NH. In order to compare them in Fig. 3, we have convolved the dust optical depth with the LAT PSF on the one hand, and we have subtracted non-gaseous γ-ray emissions obtained by the fitting (Sect. 4) from the γ-ray data on the other hand. Figure 3 shows strong similarities in the spatial distributions of both tracers, but it also reveals differences in their dynamical range in several places.

To detect neutral gas unaccounted for in NHI and WCO, we have used the fact that it is permeated by both cosmic rays and dust grains. We have therefore extracted the significant γ-ray and dust residuals over the NHI, WCO, and NHII expectations and we have used the spatial correlation between those residuals to infer the additional gas column densities (see Sect. 3.3). We have separated the residuals for two types of environments: in regions of weak or no CO intensity, below 7 K km s-1, which correspond to the DNM at the H-H2 transitions; and in regions toward large CO intensities, above 7 K km s-1, to capture additional H2 gas where 12CO emission saturates and rarer isotopologues such as 13CO or C18O should be used. We refer to this saturated-CO molecular component as “COsat”.

3.1. Dust model

In the case of a uniform dust-to-gas mass ratio and uniform mass emission coefficient of the grains, the dust optical depth linearly scales with the total NH. We have therefore modelled τ353(l,b) in each direction as a linear combination of the gaseous contributions from the different phases (H ii , H i, DNM, CO-bright, COsat), with free normalisations to be fitted to the data, as in Planck Collaboration Int. XXVIII (2015). We have added a free isotropic term, yiso, to account for the residual noise and the uncertainty in the zero level of the dust data (Planck Collaboration XI 2014). The τ353(l,b) model can be expressed as: (6)where NHI,i(l,b), WCO,i(l,b), and Iff(l,b) respectively denote the NHI, WCO, and free-free maps of the clouds depicted in Fig. 1. and stand for the column densities in the DNM and COsat components deduced from the coupled analyses of the γ-ray and dust data (see Sect. 3.3).

The y model parameters have been estimated using a χ2 minimization. We expect the model uncertainties to exceed the measurement errors in τ353(l,b) because of potential variations in grain properties through the clouds and because of the limitations of the gas tracers (survey sensitivities, emission saturation, self-absorption, etc.). As we cannot precisely determine the model uncertainties, we have set them to a fractional value of the data and we have determined this fraction to be 19% by reaching a reduced χ2 of unity. This fraction is larger than the 3% to 9% error in the measurement of τ353 across this region (Planck Collaboration XI 2014).

3.2. Gamma-ray model

Earlier studies have indicated that the bulk of the Galactic CRs radiating at 0.4–100 GeV have diffusion lengths far exceeding typical cloud dimensions and that they permeate all the H i-bright, DNM, and CO-bright gas phases. The observed γ-ray emission can therefore be modelled, to first order, by a linear combination of the same gaseous components as in the dust model. We have assumed that the emissivity spectrum of the gas follows the average one obtained in the local ISM (qLIS(E), Casandjian 2015), but we have left a free normalisation in each energy band to account for possible deviations in CR density and spectrum. The model includes other radiation components such as the Galactic IC radiation, IIC(l,b,E), the isotropic intensity mentioned above, Iiso(E), and point sources with individual flux spectra Sj(E). We have verified that the soft emission from the Earth limb is not detected in the present energy range for the choice of maximum zenith angle. The soft and transient emission from Sun and Moon is not expected to be detected as the number of γ-ray photons they emit over 6 years is negligible compared to those of the ISM components in the energy range studied. The γ-ray intensity I(l,b,E), expressed in cm-2 s-1 sr-1 MeV-1, can thus be modelled as: (7)with the and maps extracted from the coupled dust and γ-ray analyses (see Sect. 3.3).

thumbnail Fig. 4

Photon yields arising in the γ-ray model in the 0.4–100 GeV band from the H i(HI label) and CO-bright (CO label) phases of the Cetus, Main Taurus, South Taurus, North Taurus, Perseus, California, and Galactic disc clouds, from the ionised gas (ff label), from the DNM and COsat gas column densities, from the Galactic IC emission, and from the isotropic background (iso label) and γ-ray point sources.

Open with DEXTER

The input qLIS spectrum was based on four years of LAT data and on the correlation between the γ radiation and the NHI column densities derived from the LAB survey, for a spin temperature of 140 K, at latitudes between 7° and 70° (Casandjian 2015). The qHI,i scale factors in the model can therefore compensate for differences in the H idata (calibration, angular resolution, spin temperature) and potentially for cloud-to-cloud variations in CR flux. Such differences will affect the normalizations equally in all energy bands whereas a change in CR penetration in a specific cloud will show as an energy-dependent correction. For each cloud, the average γ-ray emissivity spectrum per H atom in the atomic phase is estimated from the product of the qLIS spectrum and the best-fit qHI,i normalization. This emissivity can be used to estimate the gas mass present in the other DNM, CO, and COsat parts of the cloud if one assumes a uniform CR flux across the whole structure.

The model includes 126 points sources from the 3FGL catalogue (Acero et al. 2015) inside the analysis region. Their flux spectra, Sj(E), have been computed with the spectral characteristics provided in the catalogue. Their qSj(E) normalization in the model allows for possible changes due to the longer exposure and different photon reconstruction dataset used in the present analysis (six years of Pass 8 data instead of four years of Pass 7 reprocessed data for 3FGL) and due to the use of a different interstellar background for source detection in 3FGL. The sources have been fitted simultaneously. The sources present in the 4°-wide peripheral band around the analysis region have been merged into a single map, Sext(l,b,E), and its normalization, qSext, has been left free in each energy band. The IC emission has been derived from a GALPROP model (Ackermann et al. 2012e) and isotropic emission has been determined over the whole sky with a different interstellar background model (Casandjian 2015), so we have left their scaling free in each energy band.

In order to compare with the LAT photon data in the different energy bands, we have multiplied by the exposure and convolved by the LAT PSF each model component. Figure 4 gives the photon yields in the total energy band obtained by the fitting (Sect. 4) from those components. It shows that the emission originating from the gas dominates over other types of emission and that the LAT angular resolution allows the spatial separation of the various clouds, and of the different gas phases within the clouds. The q model parameters have been fitted to the LAT photon data in each energy band using a binned maximum-likelihood with Poisson statistics.

3.3. Analysis iterations

An important aspect of our analysis is the iterative coupling of the dust and γ-ray models in order to extract the common DNM and COsat gas components which are present in both datasets, but for which we have no a priori templates. Both components show up as positive residuals over the expectations from the H i-bright, CO-bright, and free-free-bright gas components. To extract them, we have built maps of the positive residuals between the data (dust or γ rays) and the best-fit contributions from the NHI, WCO, free-free, and ancillary (other than gas) components. We have kept only the positive residuals above the noise (see below). We have separated the residuals into DNM and COsat maps according to the WCO intensity in each direction (below and above 7 K km s-1 for the DNM and COsat, respectively). The DNM and COsat templates estimated from the dust are provided to the γ-ray model ( and in Eq. (7)); conversely, the DNM and COsat column densities derived from the γ rays are provided to the dust model ( and in Eq. (6)). We have started the iterative process by fitting the dust optical depth with the H i, CO, free-free, and isotropic components to build the first DNM and COsat maps for the γ-ray model. We have then iterated between the γ-ray and dust models until a saturation in the log-likelihood of the fit to γ-ray data is reached (from the third to the fourth iteration).

The estimates of the q and y model coefficients and the DNM and COsat templates improve at each iteration since there is less and less need for the other components, in particular the H iand CO ones, to compensate for the missing gas. They still do at some level because the DNM and COsat templates 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 Gaussian noise around zero. A simple cut at zero would induce an offset bias, so we have de-noised the residual maps using the multi-resolution support method implemented in the MR filter software (Starck & Pierre 1998). The stability of the iterative analysis and the results of the dust and γ-ray fits are discussed in the following sections.

thumbnail Fig. 5

Number distribution of the dust model coefficients over 1000 jackknife fits for an H ispin temperature of 400 K. yHI,i, yCOsat, and yDNM are in units of 10-26 cm2, yCO,i in 10-6 K-1 km-1 s, yff in 3.8 10-11 Jy-1 sr and yiso in 10-6.

Open with DEXTER

thumbnail Fig. 6

Number distribution of the γ-ray model coefficients obtained for 1000 jackknife fits for an H ispin temperature of 400 K. qCO,i are in units of 1020 cm-2 K-1 km-1 s, qff in 3.8 × 1015 cm-2 Jy-1 sr, qDNM and qCOsat in 1025 cm-2. qHI,i, qIC and qiso are normalisation factors.

Open with DEXTER

thumbnail Fig. 7

Residual maps (data minus best-fit model, in sigma units) in dust optical depth (στ353) and in γ-ray counts in the 0.4–100 GeV band (σγ) obtained when including different sets of gaseous components in the models: without (left) and with (middle) the DNM and COsat components in addition to the H i, CO, and free-free templates; a close-up view (right) of the γ-ray residuals without (top) and with (bottom) the free-free map. The grey circles in the left-hand column highlight a DNM shell. The thin black contours outline the shape of the CO clouds at the 7 K km s-1 level chosen to separate DNM and COsat components. The thick black contours in the right-hand column outline the shape of the NGC 1499 and G159.6-18.5 H iiregions at the 2 × 104 Jy sr-1 level in free-free emission at 70 GHz. The models are derived for an H ispin temperature of 400 K.

Open with DEXTER

4. Results

4.1. H I optical depth correction

The unknown level of H iopacity in the different maps induces systematic uncertainties on the NHI column densities, therefore on the H icontributions to the models. The γ rays can help constrain the average level of H ioptical-depth correction applicable to the whole region by comparing the TS-dependent contrast of the NHI maps with the structure of the γ-ray flux emerging from the H igas. We have not tested different spin temperatures for each cloud complex. The maximum log-likelihood value of the γ-ray fit peaks for an H ispin temperature near 400 K (see Fig. C.1). In the following, we present the results obtained for this temperature, unless otherwise mentioned.

Spin temperatures are generally measured from pairs of H iemission and absorption spectra against bright radio sources (Heiles & Troland 2003; Mohan et al. 2004; Kanekar et al. 2011; Roy et al. 2013b; Murray et al. 2015). Line-of-sight harmonic means assume a single spin temperature along a sight line over the whole H iline or per velocity channel, as in the present work (monophasic TS). The average value near 400 K that we find is consistent with the sparse distribution of previous monophasic spin measurements in this region (Mohan et al. 2004; Kanekar et al. 2011) and with the temperature span of 200–600 K found at high Galactic latitudes for the range of NHI column densities dominating the present maps. But the monophasic assumption is known to bias the spin temperatures to higher values and to characterise the mixture of the cold and warm neutral mediums (CNM and WNM) along the line of sight, rather than the physical temperatures of individual structures (Murray et al. 2015). Monophasic temperatures tend to decrease with increasing NHI column densities because of the increasing proportion of CNM along the line of sight (Kanekar et al. 2011).

4.2. Best fits and jackknife tests

The best dust and γ-ray fits that could be achieved with the models described by Eqs. (6)and (7)include all the free-free, H i, DNM, CO, and COsat templates for the gaseous components. Table D.1 gives the corresponding best-fit coefficients and Sect. 4.5 discuss the goodness of fit with these models. Sections 4.3 and 4.4 further discuss how the models have improved by adding gas templates other than H iand CO. We focus here first on the determination of the uncertainties associated with each parameter in the models.

The large spatial extents of the maps and the existence of tight correlations between the gas, dust, and γ-ray distributions yield small statistical errors on the best-fit coefficients. They have been inferred from the information matrix of the fit (Strong 1985) and they include the effect of correlations between parameters. The gas parameters of the local clouds have precisions of 1–5% and 1–4% in dust and γ rays, respectively. Only the small contributions of the H iiregions and of the Cetus CO clouds to the γ-ray data have respective uncertainties of 10% and 14%.

We have checked the magnitude of systematic uncertainties in our linear modelling approximations, hence of spatial changes in the model and/or in the mean level of H iand CO self-absorption. We did so by repeating the last dust and γ-ray fits a thousand times over subsets of the analysis region, namely after masking out 20% of the pixels with a sum of 2.625°-wide, randomly selected squares. In γ rays, the jackknife tests have been performed only for the total 0.4–100 GeV energy band.

Figures 5 and 6 show the distributions thus obtained for the best-fit coefficients that are significantly detected. Most of them show Gaussian-like distributions, thereby indicating that the results presented in Table D.1 are statistically stable and that the average coefficients that describe our models are not driven by subset regions in each cloud complex. Several distributions exhibit long, non-Gaussian tails when the corresponding clouds subtend small solid angles (e.g., the H iiregions or the CO clouds from the Galactic disc or Perseus). The long tails reflect the indeterminacy of the parameter when a large fraction of their maps are masked.

The standard deviations found in the jackknife distributions amount to 1–3% and 1–4% for the extended H iand DNM components in dust and γ rays, respectively. The deviations are slightly larger, respectively 1–5% and 2–10%, for the more compact CO clouds in the local ISM.

We have quadratically added the 1σ fitting error and the standard deviation of the jackknife distributions to give the statistical errors listed in Table C.1. The jackknife tests have been carried out for the total energy band. We have assumed the same relative deviations, σq/q, for the individual energy bands as for the total one since the jackknife samples test potential non-uniformities in the models at larger scales than the angular resolution of the data.

The results indicate that the small contribution of the faint Cetus CO clouds to the total γ-ray photon counts is only detected below 4 GeV. Similarly, the faint CO clouds from the edge of the Galactic disc are not detected against the much larger γ-ray contribution from the Galactic atomic gas (see Fig. 4). In the dust fit, all components are significantly detected.

We have checked the convergence of the iterative scheme described in Sect. 3.3 by monitoring the increase in the maximum log-likelihood value of the γ-ray fit as we progressed into the iterations until saturation. The log-likelihood ratio 2Δln(L) of 1180 between the last and first iterations indicates a large improvement in the quality of the fits because of a better separation of the different gas structures. The initial fits yielded too-large H iand CO parameters to compensate for the missing gas. As the construction of the DNM and COsat maps improved through the iterations, the bias on the H iand CO parameters relaxed. Between the first and last iterations, the H iand CO γ-ray emissivities dropped by 5–21% and 4–22%, respectively, depending on the cloud complex and its spatial overlap with DNM and COsat structures. Similarly, the dust opacities of the H iand CO clouds decreased by 10–42% and 36–87%, respectively.

thumbnail Fig. 8

Hydrogen column density maps in the DNM and COsat components, derived from the dust (top) and γ-ray (bottom) data. The demarcation between the two components is outlined by the white contour at 7 K km s-1 in CO line intensity.

Open with DEXTER

4.3. Detection of DNM and COsat gas

Prior to building the complete models, we have checked for the presence of substantial amounts of gas that are not traced by the H iand CO line intensities or free-free emission, but that are perceptible in parallel in the independent γ-ray and dust data sets. To do so, we have performed the dust and γ-ray fits with only the H i, CO, and free-free templates for the gaseous components. We find that the data exceed the model predictions in several extended regions (see the positive residuals in the left column of Fig. 7) and that the excesses correlate in space in the independent dust and γ-ray data. Adding the DNM and COsat templates to the models significantly improves the quality of the fits (see the middle column of Fig. 7). The log-likelihood ratio of 6427 between the γ-ray fits with and without the DNM and COsat components implies that those components are detected at a formal confidence level greater than 80σ (Neyman & Pearson 1933). Hadronic interactions between cosmic rays and dust grains, or IC up-scattering of the dust thermal radiation would yield too few γ rays to explain such a correlation (Grenier et al. 2005). A coincident variation both in CR density and in dust-to-gas ratio is very improbable, so the CRs and dust jointly reveal gas in excess of the H i, CO, and ionised expectations in those directions. The spectrum of the γ-ray emission associated with both components is consistent with that found in the other gas phases over the 0.4–100 GeV energy range (see Sect. 5). This gives further support to a gaseous origin of the γ-ray emission in the DNM and COsat structures and to a smooth penetration of the CRs from the atomic envelopes to the dense parts of the clouds with saturated CO lines.

A log-likelihood ratio of 1404 between the models with different or with equal emissivities for the DNM and COsat components strongly supports their separation into distinct components. The DNM map traces column densities of dense H iand/or diffuse H2 at the atomic-molecular transition of the clouds whereas the COsat map reveals dense H2 in excess of that linearly traced by WCO because of the saturation of the CO lines (WCO remains constant as NH2 keeps increasing). Without kinematic information, we cannot separate the DNM and COsat contributions from the different cloud complexes when they overlap in direction. To extract the column densities displayed in Fig. 8, we have assumed that the DNM and COsat gas is pervaded by the average CR density measured in the atomic phase of the different clouds in the region (Grenier et al. 2005; Planck Collaboration Int. XXVIII 2015). We show in Sect. 5 that the small cloud-to-cloud dispersion in γ-ray emissivity per gas nucleon justifies the use of the average. We defer the discussion of the DNM and COsat gas masses and their relation to the H i-bright, 12CO-bright, and 13CO-bright phases to a companion paper.

4.4. Detection of ionised gas

We have also checked the significance of the addition of the free-free template to the final γ-ray and dust models, and we have tested the gaseous origin of the corresponding γ-ray signal. We develop both points in this section.

The log-likelihood ratio of 116 between the models with and without the free-free template supports a γ-ray detection of the H iiregions at a formal confidence level greater than 10σ. The close-up view in Fig. 7 shows the clear detection of a γ-ray excess toward the bright NGC 1499 region and the more marginal gain in adding the small column densities of ionised hydrogen from the fainter G159.6-18.5 region. By replacing the free-free map by a broad Gaussian source located toward the centre of NGC 1499, we have checked that the fit significantly requires an extended excess (7σ), with a FWHM comprised between 1.4° and 3.6° at the 95% confidence level. The Gaussian source, however, yields a poorer fit than the more elongated free-free emission.

The two H iiregions, NGC 1499 and G159.6-18.5, are respectively excited by a giant and a main-sequence O star. Part of their intense stellar light is upscattered to γ rays by the local CR electrons (Orlando & Strong 2007). Given the distance, luminosity, and effective temperature of each star, we have calculated the IC γ-ray flux produced in its vicinity and we have integrated the flux over the angular sizes of the H iiregion, in the 0.4–100 GeV energy band. We find IC fluxes of ~ 2.9 × 10-11 cm-2 s-1 and ~ 5.9 × 10-12 cm-2 s-1, respectively, that cannot account for the γ-ray fluxes of (2.5 ± 0.4) × 10-8 cm-2 s-1 and (2.3 ± 0.4) × 10-9 cm-2 s-1 seen in correlation with the free-free emission. Hence, the γ-ray excesses seen toward NGC 1499 and partially toward G159.6-18.5 are more likely due to hadronic interactions of the local CRs in the ionised gas.

The weak correlation between the dust optical depth τ353 and the free-free template highlighted by the jackknife test probably reflects the difficulty of estimating τ353 toward regions of strong spectral gradients (see Fig. 2). In particular because the β index and the grain temperature have been derived at different angular resolutions (Planck Collaboration XI 2014).

Assuming a uniform CR flux in the H iand H iigas phases, we have used the average γ-ray emissivity of the atomic gas in the region (see Sect. 5) and the best-fit qff parameter obtained for the 0.4–100 GeV band to convert the free-free intensities to hydrogen column densities and compare them with the values expected from Eq. (5). The electron temperature of H iiregions is known to vary with Galactocentric radius because of the metallicity gradient (Alves et al. 2012). Adopting a mean value in the solar neighbourhood between 7000 K and 8000 K, we find an average electron density of (4.3 ± 0.6) cm-3 in the H iiregions sampled here.

At the distance d ≈ 380 pc to ξ Persei inside NGC 1499 (van Leeuwen 2007), the 80  radius of the ionised region is L = 9 pc. The fractional thickness of the shell, relative to L, is estimated to be l = 0.38 (Shinohara & Ishida 1976) and the nebula subtends a solid angle 4πσ with at the exciting star. The Hα flux is corrected for the interstellar absorption as Fα0 = 100.292AVFα, with AV the visual extinction up to 500 pc derived from the 3D extinction map of Green et al. (2015). According to Shinohara & Ishida (1976) for equal electron and ion volume densities, the mean electron density of a shell-like nebula such as NGC 1499 can be expressed as with the effective volume (8)and the excitation parameter (9)The electron temperature Te of the nebula can be estimated by comparing the free-free and Hα emissions. The Hα intensity in Rayleigh units is expressed as a function of the emission measure as: (10)Combining equations (4)and (10), we obtain (11)Given this temperature, we can estimate the mean electron density.

thumbnail Fig. 9

Spectral evolution, relative to the local interstellar spectrum qLIS, of the γ-ray emissivities of the local gas components. qCO is given in units of 1020 cm-2 (K km s-1)-1 and qDNM in 1025 cm-2. The qHI normalizations are given for a spin temperature of 400 K. The dashed lines compare the results obtained for the whole 0.4–100 GeV band to the results of the four independent bands. All error bars are symmetrical; horizontal error bars at the boundaries are truncated.

Open with DEXTER

If we apply those calculations to NGC 1499 within a low contour at Iff = 2 × 104 Jy sr-1, we find a mean temperature of (6100 ± 2700) K and a mean electron density cm-3. If we adopt a more restrictive contour where the γ-ray residuals exceed 0.8σ in the best-fit model without the free-free template, we find K and a mean density cm-3 in good agreement with the gas density estimate inferred from the correlation between the γ rays and the nebular free-free emission. This agreement gives further weight to a hadronic origin of the γ rays associated with NGC 1499.

4.5. Final model residuals

The residuals obtained between the dust and γ-ray data and the best fits that include all gaseous components are presented in the middle column of Fig. 7. They show that the linear model provides an excellent fit to the γ-ray data in the overall energy band, as well as in the separate energy bands which are not shown. The residuals are consistent with noise at all angular scales except, marginally, toward the brightest CO peaks where the model tends to over-predict the data. Significant positive residuals remain at small angular scales in the dust fit. They closely follow the spatial distribution of the DNM and they are likely due to the limitations in angular resolution and in sensitivity of the γ-ray DNM template compared to its dust homologue. Small-scale clumps in the residual structure can also reflect localised variations in dust properties per gas nucleon that are not accounted for in the linear models, in particular toward the brightest CO regions. These effects are discussed in Sect. 7.

A wide, 17°-diameter shell, centred on l = 158° and b = − 24° is apparent in the initial dust residuals, and also marginally in the γ-ray residuals (see the left column of Fig. 7). The shell is rich in DNM gas (see Fig. 8). We have found no coincident structure at other wavelengths that could explain this shell in terms of a nearby H iiregion or old supernova remnant.

5. γ-ray emissivity of the gas

The γ-ray emissivity spectra of the local gaseous components are presented in Fig. 9. We note that all the qHI normalizations exceed one because the present fits preferred a larger H ispin temperature than the value of 140 K adopted for the derivation of qLIS across the whole sky off the Galactic plane (Casandjian 2015). We find no significant spectral variations except in Perseus for which we observe opposite trends for qHI decreasing and qCO increasing with energy. This behaviour can be attributed to the cross talk between the very compact set of H iand CO clouds in Perseus as the LAT PSF degrades at low energies.

The results indicate that the CR population permeating the various phases of the different clouds has the same energy distribution as the average in the local ISM (Casandjian 2015). We find no spectral signature that would betray exclusion or concentration processes in a particular cloud or with increasing gas density, from the diffuse atomic gas in Cetus up to the densities of 103 − 4 cm-3 sampled by 2.6-mm CO line emission. Together with similar findings in the Chamaeleon clouds (Planck Collaboration Int. XXVIII 2015), these results support the theoretical predictions that CRs with energies between a few hundreds of MeV and a few tens of GeV smoothly diffuse through the H i-bright and CO-bright parts of the ISM (Skilling & Strong 1976; Cesarsky & Volk 1978; Padovani & Galli 2011), up to column densities of 1021 − 22 cm-2 (Morlino & Gabici 2015). The volume densities of the gas covered by H iand CO observations are modest, but they hold the bulk of the interstellar gas mass, γ rays should reliably trace most of the interstellar gas.

thumbnail Fig. 10

Distribution with Galactocentric radius (top) and altitude above the Galactic disk (bottom) of the 0.4–10 GeV emissivities measured in the atomic gas of nearby clouds for H ispin temperatures between 125 and 150 K. The solid and grey band respectively give the average emissivity and ±1 rms dispersion in the sample. The dashed line marks the average emissivity measured across the sky at | b | ≥ 7° (Casandjian 2015).

Open with DEXTER

We can compare the γ-ray emissivities measured in the anticentre clouds for a spin temperature of 150 K with the average value found in the local ISM for Tspin = 140 K (Casandjian 2015), as shown in Fig. 10. The 35 ± 1% lower emissivity that we find for the atomic gas in the Galactic-disc background is not useful for comparison as it integrates all distances beyond the local ISM, through regions of decreasing CR density along lines of sight outward through the Galaxy (Abdo et al. 2010; Ackermann et al. 2011) and to high altitudes above the Galactic disc (Tibaldo et al. 2015). Within a few hundred parsecs, we find small cloud-to-cloud differences in γ-ray emissivity. The Cetus, South Taurus, and North Taurus clouds compare well with the average in the local ISM whereas the California, Main Taurus and Perseus clouds appear to be 15 ± 3%, 1 ± 2%, and 36 ± 6% more emissive than the average, respectively. Except for Perseus, these variations are commensurate with previous measurements. Figure 10 shows a cloud-to-cloud dispersion of ± 9% that does not relate to the cloud altitude above the Galactic disc, nor to the cloud location with respect to the local spiral arm (Galactocentric distance). This dispersion largely stems from uncertainties in the derivation of the NHI column densities. We note indeed that the qHI emissivity is 10% to 20% larger in the Cepheus-Polaris, Main Taurus, California, RCrA and Chamaeleon clouds where NHI column densities in the 80th percentile exceed 1021 cm-2, a level at which H ioptical depth corrections become significant.

The emissivity in Perseus is at variance with the other nearby clouds, but less so at high energy where the LAT PSF allows a firmer separation between the compact atomic and molecular phases of the cloud. The origin of the high emissivity will be investigated using photon selections with improved angular resolution as soon as the photon statistics permit the analysis. For the anticentre region under study, we derive very consistent averages regardless of whether we include the Perseus data to the other five measurements or not. We have used this average emission rate in the 0.4–100 GeV band, s-1 sr-1 H-1, in the atomic gas to infer column densities in the other gas phases.

6. XCO factors

The XCO factor is often applied to convert WCO intensities to hydrogen-equivalent NH2 column densities in the molecular phase. In γ rays, under the assumption of a uniform CR flux in the H i-bright and CO-bright phases, it is derived as XCOγ = qCO/ (2qHI), independently in each energy band. This ratio also assumes that helium is uniformly mixed by number with hydrogen in the ISM. Likewise, assuming a uniform dust-to-gas mass ratio and a uniform emission coefficient κ353 of the grains at 353 GHz, the XCO factor is derived in the dust fit as XCOτ = yCO/ (2yHI). The results of the two methods for the six CO clouds selected in our analysis are given in Table 1. The γ-ray results in the individual energy bands are displayed in Fig. 11.

Table 1

XCO factors in 1020 cm-2 K-1 km-1 s for the dust and γ-ray fits.

thumbnail Fig. 11

Evolution of the XCO factor, as measured in γ rays for various linear resolutions in the different clouds. The open circle marks the value obtained in the overall 0.4–100 GeV band, in close agreement with the weighted average of the four independent energy bands (black line) and its 1σ errors (dashed lines).

Open with DEXTER

6.1. XCO measurements in γ rays

We find XCO values for the present clouds that compare well with previous γ-ray measurements in other nearby clouds. For XCO factors obtained with Fermi-LAT at comparable angular resolution, we can cite the values of in the Cepheus-Polaris complex (Abdo et al. 2010), of in the RCrA cloud (Ackermann et al. 2012c), of 1.21 ± 0.02 in the Orion clouds (Ackermann et al. 2012b), and of in the Chameleon clouds (Planck Collaboration Int. XXVIII 2015), in units of 1020 cm-2 K-1 km-1 s. We note that the XCO factors in those well-resolved clouds are all close to or lower than 1020 cm-2 K-1 km-1 s. Table E.1 lists historical γ-ray estimates of XCO. It shows that the estimates in nearby clouds have remained stable or only slightly decreased as the angular resolutions of the H i, CO, and γ-ray surveys have improved. The lower level of cross-correlation between the H iand CO structures on the one hand, a better separation of the γ radiation produced in the atomic and molecular phases on the other hand, and a lower contamination from unresolved γ-ray point sources has primarily improved the precision of the measurements. The latest analyses with <10.8  FWHM resolution in H i, 8.5  resolution in CO, and Fermi-LAT data consistently yield low XCO values, ranging between 0.5 ×1020 cm-2 K-1 km-1 s and 1.2 ×1020 cm-2 K-1 km-1 s. Systematic uncertainties on the NHI column densities due to H ioptical depth corrections impact the estimation of the γ-ray emissivities per gas nucleon. Analyses for different H ispin temperatures show that the resulting systematic uncertainty on the XCO factors amount to (0.1–0.2) ×1020 cm-2 K-1 km-1 s.

The values listed in Table E.1 highlight a recurrent discrepancy between the XCO measurements at small scale in nearby clouds and the two-to-three times larger values found at kiloparsec scales when averaging over a population of molecular complexes in spiral arms or in Galactocentric rings. We return to this point at the end of this section.

In Fig. 11 we take advantage of the energy-dependent PSF of the LAT to probe various linear scales within the present nearby clouds. The scale is derived using the half-width at the half-maximum of the PSF integrated over the energy band for the qLIS spectral shape. We find no evidence of XCO changes at parsec scales except in Perseus.

thumbnail Fig. 12

Evolution of the XCO factor measured in γ rays as a function of the average WCO intensity, , the surface fraction of dense gas, SFdense, the mean visual extinction in the CO-bright phase, , the average CO line width, , the H2 mass in the CO-bright phase, MH2, and the velocity dispersion of the CO line. Dashed lines give the best (χ2) linear regressions.

Open with DEXTER

Because the qHI and qCO variations in Perseus are inversely coupled in energy, the resulting change in XCO is likely due to an increased level of cross-correlation between the compact H iand CO phases as the LAT PSF degrades (Figs. 4 and 9). If we force the γ-ray emissivity of the H igas to be the same as the average found among the other clouds, we obtain a larger XCO value of (0.68 ± 0.04)×1020 cm-2 K-1 km-1 s, close to that measured above 4 GeV with the best LAT angular resolution. However, the use of this value implies a significant (4.6σ) degradation of the fit quality at lower energies. The likelihood analysis significantly supports a lower XCO factor in Perseus. A low value is also indicated by the average of 0.3 ×1020 cm-2 K-1 km-1 s found at 0.4-pc resolution in the dust, H i, and CO study of Lee et al. (2014).

Several observations and numerical models have suggested that the XCO factor changes with the density, turbulent velocity, and metallicity of the gas, and with the penetration of the CO-dissociating radiation inside a cloud (Polk et al. 1988; Bell et al. 2006; Sheffer et al. 2008; Liszt et al. 2010; Glover & Mac Low 2011; Shetty et al. 2011a; Liszt & Pety 2012; Lee et al. 2014; Bertram et al. 2016). One therefore expects intrinsic XCO gradients in clouds. Calculations by Bell et al. (2006) show that XCO should dramatically drop inward from the outskirts of a molecular cloud as the photo-dissociating radiation gets absorbed. The XCO factor drops down to a minimum value at AV around 2 mag. The XCO ratio then rises again in the denser cores because of the saturation of the optically thick CO lines. One goal of our analysis is to explore whether these intrinsic variations can change the mean XCO value averaged over a cloud, depending on how diffuse it is. We cannot yet measure the steepest XCO gradients at the periphery of clouds (AV≪ 1 mag) near the sensitivity threshold of CO surveys because of the complex 3D interface with the abundant atomic phase. We cannot map XCO gradients pixel by pixel in γ rays because of the photon statistics, but our sample includes more or less translucent clouds of the local ISM and we can explore changes in the mean value of XCO for clouds in different states.

Figure 12 compares γ-ray measurements of XCO in nearby clouds from this work and Planck Collaboration Int. XXVIII (2015). We restrict the sample to estimates obtained in the Chamaeleon and anti-centre clouds because they are based on the same analysis method, the same γ-ray energy bands to ensure the same angular resolution of the LAT, and the same sampling resolution in the H iand CO data. They both include the DNM gas in the models. These clouds subtend large solid angles in the sky and their column density structures are well resolved.

The distributions in Fig. 12 show that the average XCO factor per cloud does not depend on the H2 mass mapped in CO, nor on the overall dynamics of the cloud characterised by the velocity dispersion of the CO lines. But XCO appears to depend on the cloud structure and in particular on its diffuseness. To reflect changes in the latter, we have explored several characteristics:

  • the mean WCO intensity, , in a cloud, taken above1 K km s-1 to avoid noise fluctuations;

  • the surface fraction of dense regions with large WCO intensity within a cloud, SFdense = SWCO> 7 K km s-1/SWCO> 1 K km s-1 (with S a solid angle), which gauges the relative weight of diffuse and dense molecular regions in the determination of XCO;

  • the mean visual extinction toward a cloud, , taken from the NICER M2a 12-resolution AJ map (Juvela & Montillaud 2016), translated into AV using a colour ratio of 3.55 according to the extinction law of Cardelli et al. (1989). Only values in directions with WCO> 1 K km s-1 and AV> 0.8 mag. have been retained in the average to avoid noise contamination;

  • the average CO line width, , in a cloud since more CO line photons can escape the cloud and ultimately be detected as the line width grows. The resulting increase in WCO intensity implies a decrease in XCO with (Shetty et al. 2011b).

Figure 12 shows that, for all these diagnostics, XCO tends to decrease from diffuse to more compact clouds. The change in XCO is approximately a factor of two. The decrease remains significant if we use the Perseus XCO value found when forcing the H iγ-ray emissivity to the local average, or if we exclude the Perseus data point. A uniform XCO factor is strongly rejected (> 30σ) for all the diffuse diagnostics, whether we modify or exclude the Perseus data point or not.

The trend is consistent with the expectation that XCO decreases from the diffuse envelopes, exposed to the photodissociating radiation field, to the shielded cores of the CO clouds (Bell et al. 2006; Glover & Mac Low 2011; Shetty et al. 2011a; Bertram et al. 2016). Hence, the average XCO factor per cloud in our sample is seen to vary with the relative solid angle subtended by the diffuse and dense parts of the clouds, reflected in SFdense and . The turbulence level and the degree of virialization of the cloud (e.g., the ratio of kinetic to gravitational energy) also influences the average XCO factor by allowing the ISRF to penetrate more or less deeply into the turbulent layers of the cloud (Bell et al. 2006; Shetty et al. 2011b; Bertram et al. 2016). The mean XCO factor per cloud in Fig. 12 does decrease with the average line width of the CO lines. The present measurements also give weight to the idea that the mass or size of a cloud does not play an important role in determining the XCO factor (Shetty et al. 2011b).

According to XCO calculations as a function of visual extinction AV (Bell et al. 2006), XCO should drop by more than two orders of magnitude with increasing AV, down to a minimum at 1 ≲ AV ≲ 3 mag. For long-lived clouds, exposed to the local cosmic-ray ionization rate of 1.4 × 10-17 s-1 (Grenier et al. 2015) and to the local ISRF, the minimum XCO value reaches below 1020 cm-2 K-1 km-1 s for hydrogen densities larger than 103 cm-3 and/or for turbulent velocities larger than 1 km s-1. Yet, the measured CO lines have widths below 0.7 km s-1 and the clouds in our sample are not massive enough to sustain such large average densities through most of their volume. Moreover, models of dense clouds show a steady increase in XCO at all NH2 column densities (Shetty et al. 2011a), at variance with the declining trends found in Fig. 12.

We have also compared the mean visual extinction, , mean CO intensity, , and XCO factor in each cloud with the simulation results of Glover & Mac Low (2011). The trend seen in Fig. 12 roughly compares in slope with the simulation expectations, but not in absolute values, even though the efficiency of the radiation screening compares well in the observations and simulations. We show in Sect. 7.1 that the AV/NH ratio measured in the clouds compares well with the Bohlin et al. (1978) value adopted in the simulations. For a given in a cloud, the simulations under-predict the mean intensity by a factor 20, and overpredict XCO by a factor 40. The comparison should be taken with care because both the observation and simulation samples are sparse, because the observations span a smaller dynamical range in than the simulations, and because the and averages depend on the threshold applied to define CO cloud edges, which differs in both datasets. Nevertheless, the comparison highlights that, for initial volume densities ≤ 300 cm-3 applicable to the set of observed clouds, the simulations produce underluminous CO clouds at low densities compared to the observations. This was noted by Levrier et al. (2012). The CO deficit in their simulation reached a factor of ten in column density in the peripheral regions exposed to the UV radiation, where NH2≲ 2 × 1020 cm-2. At larger column densities, the simulations and data compare much better. So, there seems to be a general problem with our understanding of chemistry at low gas density in the photo-dominated regions (PDR), either in the treatment of the UV attenuation or in the absence of warm chemistry driven by intermittent energy deposition in regions where the interstellar turbulence is dissipated (Godard et al. 2014).

The list of historical XCO measurements given in Table E.1 and the values we are adding with the present sample of anti-centre clouds confirm a recurrent discrepancy that exists between the mean XCO factors measured with parsec resolution in nearby clouds and the XCO averages obtained at the scale of spiral arms (Local and Perseus arms) or of the Galactic disk. The former are two to three times lower than the large-scale values, which are close to 2 ×1020 cm-2 K-1 km-1 s (see Table E.1). Simulations of molecular clouds have warned that XCO factors can vary by up to a factor of five in the Milky Way because of the complexity of the overlap and blending in position and velocity of the CO clouds emitting along the lines of sight (Shetty et al. 2011b). Smith et al. (2014) have estimated an average value of 2.2 ×1020 cm-2 K-1 km-1 s, close to the large-scale γ-ray measurements, when including all regions where CO intensities are potentially observable (i.e. wherever >0.1 K km s-1). But the spatial distribution of H2 column densities and CO intensities in their Galactic simulation implies significant deviations from this average, depending for instance on the location of cloud complexes with respect to the spiral arms.

Because of intrinsic XCO gradients across clouds, XCO averages in large-scale observations respond to the mix of clouds in different states present in the region under study. At least three potential reasons can explain why the large-scale XCO factor can exceed the local values. A first possibility is that a large fraction of the molecular mass lies in diffuse, gravitationally unbound clouds where XCO factors can reach very large values of 1021−22 cm-2 K-1 km-1 s (Bell et al. 2006; Glover & Mac Low 2011). A second possibility is that the average is driven by a large number of giant molecular clouds which are expected to exhibit XCO values of a few 1020 cm-2 K-1 km-1 s because the gas is predominantly in the optically-thick CO regime. A third influence is the level of shearing in the clouds. For instance, clouds sheared out by Galactic differential rotation are more exposed to photo-dissociation and should exhibit larger XCO values (Smith et al. 2014).

In parallel, the change in XCO estimates with scale may be due to determination biases induced by the sampling resolution in the gas maps, such as the increased difficulty at large distance to separate the clumpy CNM clouds and the DNM envelopes from their CO-bright phase. Misattributing a small part of the dense atomic gas or diffuse H2 would significantly bias the XCO factor upward. For instance, confusing 10–20% of the CNM with the CO-related signal could explain the observed change in XCO by a factor of two between the local ISM and spiral arms (Planck Collaboration Int. XXVIII 2015). In nearby galaxies, Sandstrom et al. (2013) have indeed measured a systematically larger XCO factor in highly inclined galaxies than in face-on ones where the pile-up along sight lines is reduced. More tests are needed to disentangle the origin of the XCO changes with scale. With larger photon statistics at high γ-ray energies, we should be able to investigate how the angular resolution and confusion between gas phases affect the calibration of XCO beyond the solar neighbourhood.

thumbnail Fig. 13

Ratio of XCO measurements with dust and γ rays in clouds of different average dust optical depth .

Open with DEXTER

Table 2

Average dust opacities, τ353/NH, and NH/E(BV) ratios in the gas phases of the different clouds.

6.2. XCO measurements with dust

The results of the dust model provide an independent measure of XCO. Figure 13 shows that the dust-derived XCO factors are systematically 30% to 130% larger than the γ-ray values. Systematic discrepancies had already been noted in previous studies (Planck Collaboration XIX 2011; Bolatto et al. 2013), but on the basis of different H iand CO surveys, different correlation methods, and different angular resolutions. In the present study, as in Planck Collaboration Int. XXVIII (2015), the dust and γ-ray models are fitted to the same NHI and WCO maps, at the same resolution. Bolatto et al. (2013) proposed that the inclusion of the DNM gas in γ-ray analyses and its absence in dust analyses could explain the lower XCO values obtained in γ rays, but, by construction, the inclusion of the DNM map cannot bias the XCO estimate (Grenier et al. 2005; Planck Collaboration Int. XXVIII 2015). Furthermore, the recent analyses include the DNM in both the dust and γ-ray fits and we have verified in this work that the XCO values change by less than 20% when the DNM phase is added or not in the fits.

Alternatively, the standard derivation of the XCO factor as XCOτ = yCO/ (2yHI) is based on a uniform dust opacity, τ353/NH, across the gas phases. A change in dust emission properties, in particular an increase in τ353/NH at large gas column densities, can significantly bias the XCO factor upward (Planck Collaboration Int. XXVIII 2015). We show in Sect. 7.2 evidence for such a departure from linearity in dense CO regions where τ353/NH can increase by a factor of up to six due to grain evolution compared to their properties in the atomic gas. Outside the densest molecular parts (brightest in CO), the molecular hydrogen is more diffuse and the increase in τ353/NH is more modest, by a factor ranging between 1.2 and about 2, so the resulting bias on the derivation of the average value of XCO across a cloud is less than a factor of two.

Dust reddening and extinction maps have also been used to evaluate XCO in nearby clouds. We find values in close agreement with the extinction study of Chen et al. (2015) in the main Taurus and Perseus clouds (their Taurus and Perseus E1 values). The present Taurus South complex overlaps their Taurus E1, E2, and E3 structures with respective XCO factors of (0.84 ± 0.01), (1.41 ± 0.02), and (1.69 ± 0.04)×1020 cm-2 K-1 km-1 s bracketing our result in Taurus South.

Using the 2MASS extinction map of Dobashi et al. (2013), Paradis et al. (2012) found XCO= (2.27 ± 0.9)×1020 cm-2 K-1 km-1 s in the Taurus-Perseus region, without decomposing the different clouds. Their AV/NH estimate of 5.33 ×10-22 mag cm2 in this region is 20% lower than the all-sky | b | > 10° average, but this difference is not sufficient to reconcile their measurement based on dust reddening with our values based on dust emission. The determination of AV/NH in the H iphase directly impacts the derivation of XCO in the embedded molecular phase. The different methods used to infer AV from stellar data yield systematic differences as large as 1 mag at low AV, toward the atomic and diffuse ISM (Dobashi et al. 2013; Green et al. 2015; Juvela & Montillaud 2016), so further analyses are necessary to evaluate their impact on XCO measurements.

Pineda et al. (2010) have used 2MASS data and the NICER method to map the dust reddening and to derive XCO near 2.1 ×1020 cm-2 K-1 km-1 s in the Taurus Main complex. In this case, the difference with our result can be explained because they did not separate the extinction related to H ialong the CO sight lines, but they attributed all the dust to the molecular phase. In our analyses, we separate the dust or γ rays associated with each phase along all sight lines.

7. Dust evolution

7.1. Average dust opacities per gas phase

Our analyses yield average dust properties per H atom in each gas phase. The best-fit yHI coefficients directly give the average dust opacities, , in the H istructures. We have used the γ-ray estimates of the XCO factors to derive the average dust opacities in the CO-bright phase of each complex, = yCO/ (2XCOγ). The results are listed in Table 2 together with similar measurements in the Chamaeleon clouds (Planck Collaboration Int. XXVIII 2015). To ease the comparison with other dust reddening measurements, we have converted the dust opacities to NH/E(BV) ratios using the E(BV) /τ353 slope of (1.49 ± 0.03) × 104 found in a large correlation study between τ353 and reddening measurements toward 53 399 quasars (Planck Collaboration XI 2014). We did not attempt to estimate dust opacities in the ionised gas because of the large temperature gradients which impact the τ353 derivation toward the H iiregions, with uncertainties in τ353 doubling near the hotspots.

The correlation existing between the γ-ray intensity and dust optical depth in the DNM yields an important measure of the average dust opacity in this phase. The value = qHI/qDNM assumes a uniform CR flux across the H iand DNM phases, an assumption that is corroborated by the fact that both phases exhibit comparable γ-ray emissivity spectra and moderate gas volume densities. According to Fig. 14, we find no opacity changes at parsec scales in the extended DNM structures. The average, = (16.5 ± 0.3) × 10-27 cm2, compares very well with the value of (17.2 ± 0.5) × 10-27 cm2 obtained in the DNM surrounding the Chamaeleon clouds (Planck Collaboration Int. XXVIII 2015).

Using optically thin H i, Planck Collaboration XI (2014) found marked gradients in dust opacity across the sky, with values ranging from 6.6 to about 11 ×10-27 cm2in the atomic gas, and a high-latitude value of (7.0 ± 2.0)×10-27 cm2in good agreement with the opacity of (7.1 ± 0.6)×10-27 cm2derived in high-latitude cirruses with a different analysis method (Planck Collaboration Int. XVII 2014). The opacities continued to increase in the CO-bright phase, up to 18 ×10-27 cm2for an XCO factor of 1020 cm-2 K-1 km-1 s. Our results confirm that dust opacities vary from cloud to cloud and they can reach values well in excess of the high-latitude value, both in the atomic and molecular phase.

thumbnail Fig. 14

Average dust opacities in the DNM measured in γ rays for different linear resolutions in the clouds. The open circle marks the measurement in the 0.4–100 GeV energy band, in close agreement with the weighted average of the four independent energy bands (black line) and its ± 1σ errors (dashed lines).

Open with DEXTER

We find 50% to 135% larger opacities in the anticentre H iclouds than in the more diffuse H icirruses, even though we have corrected the NHI column densities for a spin temperature of 400 K. The systematic uncertainties on τ353/NH due to the H ioptical depth correction are listed in Table 2 as the second set of errors. Given the likelihood variations shown in Fig. C.1, we have taken the standard deviation in the set of measurements as a function of spin temperature for each parameter to represent such an uncertainty. Considerations of those uncertainties cannot reconcile the individual cloud measurements and the cirrus value, even for the rather tenuous H iclouds in Cetus and Taurus North. The cloud-to-cloud dispersion indicates environmental effects beyond the latitude dependence noted by Liszt (2014) in the atomic gas. The opacities listed in Table 2 do not relate to the H imass of the corresponding clouds. The 15% dispersion in this sample rather relates to the mean visual extinction of the H icloud and to the surface fraction of large and small NHI column densities across the cloud. We present in the next section evidence for gradual opacity changes as the gas density increases that explain why the average opacity in an H icloud depends on its structure. The observation of opacity variations within the atomic phase at column densities below 1021 cm-2 challenges the current dust evolution models discussed in the next section.

We find systematically larger mean opacities in the CO-bright phase than in the H iphase, thereby confirming the continued increase in τ353/NH at large H2 column densities. The latter has been evaluated independently in each cloud (independent XCO factors), so the cloud-to-cloud changes in that Planck Collaboration XI (2014) have noted across the sky are not due to their use of a uniform XCO factor. In Table 2, the mean opacity increases by 30% to 100% between the H iand CO-bright phases, independent of the H2 mass locked in the cloud. The dispersion in rather relates to the surface fraction of dense regions with large WCO intensities, SFdense, in the averaging. To study environmental changes in τ353/NH within molecular clouds, we must take into account the added complication of spatial gradients in XCO.

7.2. Environmental changes in dust opacity

In order to follow changes in dust opacity across the clouds, we have built two maps of the total NH column density in the region. The first one, NH γ, takes advantage of the CR interactions with gas in all chemical forms and thermodynamical states. The γ-ray intensity from the gas has been obtained from the LAT data in the overall energy band after subtraction of the γ-ray counts unrelated to gas in the best-fit model. In order to reduce the Poisson noise and get more robust photon statistics, we have degraded the resolution to 0.̊375. We have converted this γ-ray intensity into NH using the average emission rate, , obtained in the atomic phase in the region. The second map, NH mλ, uses the higher-resolution information from the H i, CO, and dust data, and the mass scaling provided by the best-fit γ-ray model (uniform H ispin temperature of 400 K, XCOγ factors in the clouds, γ-ray emissivity in the DNM, ionised, and COsat gas as in the H i). The opacity distribution at the smallest angular scales should be considered with care because the scaling factors used to construct NH mλ have been derived at the scale of a whole cloud.

Figure 15 shows the 2D histogram of the correlation between the sightline integrals in dust optical depth and in NH γ column density. The distribution deviates from linearity above about 4 × 1021 cm-2. The opacity maps, shown in Fig. 16, indicate that τ353/NH evolves progressively within the clouds, especially when the medium becomes molecular, in the diffuse DNM and even more so in the CO-bright parts. The increase in opacity does not relate to an increase in the specific power, 4πR/NH, radiated by the grains, so the rise in opacity cannot be attributed to a larger heating rate (if we consider that the large grains are in thermal equilibrium to equate the radiated and absorbed powers). The largest opacities rather correspond to a ~ 30% lower heating rate in the well-shielded CO regions.

Thanks to the total gas tracing capability of the CRs and to the complementarity of the NH γ and NH mλ maps, it is very unlikely that opacity variations as large as those seen in Fig. 16 are caused by large deficits in the gas column densities. Our results give strong support to earlier indications that the big grains in the CO filaments of the Taurus cloud have larger sub-mm emissivities than in the more diffuse media. The previous studies noted opacity changes by a factor ranging from ~ 2 (Flagey et al. 2009; Planck Collaboration XXV 2011; Ysard et al. 2013) to (Stepnik et al. 2003) above the diffuse-ISM value. In the present analysis, we find opacity enhancements exceeding a factor of three and reaching a factor of six toward the CO clouds, even though we provide a measure of the additional H2 gas that is not linearly traced by WCO in the COsat directions. The largest opacities do not rise as high in the Taurus clouds as they do in the Chamaeleon ones (Planck Collaboration Int. XXVIII 2015), but a quantitative comparison would require the separation of the DNM and COsat components in a new analysis of the Chamaeleon region. We also find that the opacity starts to increase in the DNM, therefore at lower gas densities than the few thousand per cm3 sampled in CO. The grain opacity at the H i–H2 transition compares very well with that in the DNM surrounding the Chamaeleon clouds (Planck Collaboration Int. XXVIII 2015).

thumbnail Fig. 15

Upper: 2D histogram of the correlation between the total gas column density, NH γ, measured with the 0.4–100 GeV interstellar γ rays, and the dust optical depth at 353 GHz, convolved with the LAT response for an interstellar spectrum. The maps were sampled at a 0.̊375 resolution. Lower: evolution of the dust opacity in NH bins. The error bars give the standard errors of the means and the grey band gives the standard deviation of the opacities in each bin.

Open with DEXTER

thumbnail Fig. 16

Spatial variations of the dust opacities (left) and specific power (right) with the total gas measured by NH γ at 0.̊375 resolution (top) and by NH mλ at 0.̊125 resolution (bottom). The tilded quantities are convolved with the LAT response for an interstellar spectrum. The black contours outline the shape of the CO clouds at the 7 K km s-1 level chosen to separate DNM and COsat components.

Open with DEXTER

thumbnail Fig. 17

Evolution of the dust opacities as a function of the molecular fraction in the gas column (left), the total gas column density NH (middle) and the dust color temperature (right). To estimate the molecular fraction, the DNM is assumed to be 50% molecular. The error bars give the standard errors of the means and the grey band gives the standard deviation of the opacities in each bin.

Open with DEXTER

In order to determine how the opacity gradient depends steeply on the ambient ISM, Fig. 17 shows how the opacity varies with the molecular fraction in the gas column, fH2, with the total NH mλ column density as traced by the multi-wavelength data, and with the dust colour temperature. We observe a gradual increase as the gas becomes molecular and a steeper change at NH above 5 × 1021 cm-2 due to two factors: change in molecular fraction and in grain properties in the molecular gas. The amplitude of the gradient is best traced by the drop in grain temperature. Figure 16 further shows that the ambient molecular gas density is not the only factor of evolution because we detect notable cloud-to-cloud differences in the opacity gradients. Compare for instance the large opacities that often exceed 25 ×10-27 cm2 along the Main Taurus and Perseus CO clouds with the 50% lower values obtained toward comparable WCO intensities near both ends (in longitude) of the California cloud.

Several studies have noted that the opacity enhancement occurred together with a reduction of the IR emission from polycyclic aromatic hydrocarbons (PAHs) and from very small grains stochastically heated by the ambient ISRF (e.g., Stepnik et al. 2003; Flagey et al. 2009; Ysard et al. 2013). Grain evolution in dense media was invoked to explain the radiation changes (see the previous references and Martin et al. 2012; Köhler et al. 2015). It involves accretion of aliphatic-rich carbonaceous mantles from the gas phase, grain coagulation into fluffy aggregates, and the formation of an ice mantle on the surface of the aggregates in the coldest regions. Recent theoretical studies such as Köhler et al. (2015) suggest that the spectral changes in β relate to gas density-dependent processes and that the variations in opacity, spectral index β, and colour temperature are closely connected. Figure 18 shows the quantitative relation between those three parameters in the observations. The data corroborate the prediction that the opacity should increase for β above 1.5 and temperature below 17 K. According to this model, the directions with temperatures near 16–18 K, β ≥ 1.5, and τ353/NHmλ values about twice as large as the diffuse-ISM value of ~ 7×10-27 cm2would be rich in grains with carbonaceous mantles. The mantle accretion would therefore occur primarily in the DNM. Aggregates would dominate the directions below 16 K, at 1.6 ≲ β ≲ 1.8, and where τ353/NH is at least four times above the diffuse-ISM value. Such conditions are found inside the sampled CO clouds. Ice coating of the aggregates would explain the largest opacity values in Fig. 18, at the lowest temperatures and largest spectral indices, β.

8. Conclusions

We have analysed the gas, dust and CR content of several nearby anti-centre clouds including the Cetus, Taurus, Auriga, Perseus, and California clouds. We have performed an iterative fit of the total gas column density as traced by γ rays and the dust optical depth at 353 GHz. We have modelled both tracers as linear combinations of the gas column densities traced by H i, 12CO, and free-free emissions, and we have extracted the morphology and column density of the DNM from the joint γ-ray and dust data. We have verified the robustness of our set of parameters through the use of jackknife tests; these parameters have enabled us to constrain the properties of the ISM in these clouds.

thumbnail Fig. 18

Evolution of the dust opacities with the total gas measured by NH γ at 0.̊375 resolution (top) and by NH mλ at 0.̊125 resolution (bottom) as a function of the dust color temperature (x-axis) and spectral β index (y-axis)

Open with DEXTER

The main results are summarised as follows.

  • On the cosmic rays. The γ-ray emissivity of the gas in the analysedclouds has the same energy spectrum as in other clouds of the localISM and we find no dependence of their average γ-ray emissivitywith Galactocentric radius (i.e., distance from the local spiralarm), nor with height above the Galactic plane. In the0.4–100 GeV energy band and at the precisionlevel of the current LAT observations, we find no evidence of CRexclusion or CR concentration in the clouds, up to the 12CO-bright molecular regions.

  • On the dark gas. The γ-ray and dust data jointly reveal significant amounts of gas in addition to that seen in H i, free-free, and 12CO emissions. The diffuse large-scale structures are associated with the DNM at the transition between the atomic and molecular phases. They gather dense atomic hydrogen and diffuse H2 in unknown proportions, but with column densities equivalent to those found in the H iand CO emitting parts. In the molecular phase, the γ rays and dust reveal filaments of dense gas in addition to that proportionality traced by WCO where the 12CO line intensities saturate.

  • On the H iiregions. The H iiregions, NGC 1499 and G156.6-18.5, are jointly detected in dust and γ-ray emission. The corresponding γ-ray flux cannot be attributed to the IC up-scattering of the bright stellar radiation by the local CR electrons. It is more likely caused by the hadronic interactions of the local CR nuclei in the ionised gas. Hence for NGC 1499, we have used the average γ-ray emissivity of the atomic gas in the region to infer a mean electron density of 4.3 ± 0.6 cm-3 for an electron temperature of 8000 K, in good agreement with the modelling of the free-free and Hα emissions of the H iiregion.

  • On the XCO factors. We provide independent measurements of the XCO factor from the dust and γ-ray analyses in six different clouds, against the same H iand CO data. As in the Chamaeleon complex, we find that the dust-derived values are systematically larger, by 30% to 130%, than the γ-ray estimates. The difference is likely due to chemical and structural evolution of the dust grains with increasing NH. The XCO factors measured in γ rays range from (1.04 ± 0.05)×1020 cm-2 K-1 km-1 s in Taurus South to (0.44 ± 0.04)×1020 cm-2 K-1 km-1 s in Perseus. Together with the estimate found in the Chamaeleon clouds (Planck Collaboration Int. XXVIII 2015), these measurements indicate that XCO values tend to decrease with the average WCO intensity of a cloud, with the average CO line width, and with the surface fraction subtended by the brightest CO clumps. The more diffuse CO clouds therefore tend to have larger average XCO factors. Models of the formation and photodissociation of H2 and CO molecules predict a marked decline in XCO from the diffuse envelopes of molecular clouds to their dense cores. Hence the XCO average should qualitatively vary with the surface fraction of dense regions in the cloud structure and with the CO line width, as in our data. We find a modest change by typically a factor of two in XCO with cloud state, but our sample lacks a giant molecular cloud to extend the state range. The amplitude of these variations already limit the precision of molecular mass estimates based on the use of CO intensities and of a mean XCO conversion factor. We further note that the low XCO values measured here are at variance with the theoretical predictions from PDR chemistry for the moderate gas densities filling most of the volume of the observed clouds. The present measurements confirm a recurrent discrepancy by a factor of two to three between the XCO averages obtained in nearby clouds and at large scale in the Galaxy. Simulations predict XCO variations of this amplitude in the Galactic disc owing to the variety of CO-line blending situations occurring along the lines of sight (Shetty et al. 2011b; Smith et al. 2014). A dominance of diffuse, unbound or sheared clouds, or of giant molecular clouds, can drive the XCO factor of the ensemble to large values. Additionally, the increased level of cross-correlation at large distance between the CNM, DNM, and CO-bright phases of the ISM could bias the large-scale XCO determination upward.

  • On the dust opacities. The dust opacity at 353 GHz, τ353/NH, appears to rise by a factor of three from low column densities in the atomic gas (about 10-26 cm2) to cold grains in molecular gas at NH5 × 1021 cm-2. The rise can reach a factor of six at very low dust temperatures below 14 K. The amplitude of the rise is comparable to the variations observed in the Chamaeleon clouds (Planck Collaboration Int. XXVIII 2015). As the observed specific power radiated by the grains decreases in the cold molecular regions the changes cannot be attributed to a larger heating rate, but they are more likely caused by a chemical or structural change in the grains. The magnitude of the rise severely limits the use of the thermal emission of the large dust grains to trace the total gas. The linear regime is limited to NH<3 × 1021 cm-2 in these clouds. The confirmation of large opacity variations across clouds directly impacts the gas mass estimates inferred from dust emission at sub-mm and mm wavelengths to derive star-forming efficiencies in the Galaxy and in external galaxies. We quantify the coupled changes in grain temperature, opacity, and β index of the thermal radiation to help model the grain evolution.

In a second paper we will study how the DNM and COsat gas relate to the H i-bright, 12CO-bright, and 13CO-bright phases. We will derive the masses of the DNM and of the COsat phase for a sample of substructures in the Cetus, Taurus, California and Chameleon clouds, while taking care to avoid regions of cloud confusion along sightlines. We will follow the evolution of the CO-dark H2 fraction both within a cloud and from cloud to cloud and we will place our results in the context of other observational and theoretical studies.

In a subsequent study, we will confront the gas, γ-ray, and dust distributions in the nearby clouds using dust column densities inferred from stellar reddening. The joint comparison of dust emission and extinction data to other gas tracers will provide additional information on the relation between the total gas and dust, including the abundant dark neutral medium, and will help to better constrain the evolution of dust properties, in particular with regard to potential variations of the extinction curve (Schlafly et al. 2016).

Adding more data points to the trends shown in Fig. 12 is essential to confirm that the average XCO factor in a cloud tends to decrease for increasingly dense CO molecular clouds, in response to different XCO gradients spanning each cloud (Bell et al. 2006; Glover & Mac Low 2011). Particular attention should be paid to the use of the same method to derive a consistent set of XCO values, and to the addition of clouds with little confusion between the different gas phases to ensure a precise separation of the total gas tracers in each phase. The current γ-ray analyses cannot yield XCO maps across individual clouds. The variations in dust opacity with increasing gas column density prevent the use of dust emission maps to probe XCO gradients in nearby clouds. In this context, it is essential to increase the number of clouds with measured average XCO values and to confront the observed XCO trends with PDR models in order to interpret the slope and dispersion of the trends presented here. Confirming that XCO factors vary from cloud to cloud within the local ISM because of their different dynamical and chemical structures stresses that large-scale XCO averages can vary because of a different mix of cloud states in different Galactic regions, in addition to the Galaxy-wide XCO gradients expected from metallicity and ISRF gradients. For instance, it opens the way for observable differences in XCO between ensembles of clouds inside and outside of the spiral arms since the inter-arm clouds that have been sheared after their passage through an arm are more susceptible to photo-dissociation, thus more prone to having large XCO averages (Smith et al. 2014).


Acknowledgments

We thank the referee for an interesting discussion. 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’Énergie 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’Études Spatiales in France. The authors acknowledge the support of the Agence Nationale de la Recherche (ANR) under award number STILISM ANR-12-BS05-0016.

References

Appendix A: Component separation

Table A.1

Limits of the H iand CO clouds in longitude, latitude, and velocity.

We have decomposed the H iand CO spectra into sets of individual lines as in Planck Collaboration Int. XXVIII (2015). The method starts with the detection and localization in velocity of isolated lines or significant peaks in the spectra. Each HI and CO spectrum is then fitted by a sum of pseudo-Voigt line profiles centred on the centroids of the detected lines or peaks, within ± 2.7 to 3 km s-1. The fits are generally good, with less than 15% difference between the observed and fitted integrals over all velocities. In order to preserve the observed photometry, we calculate the residuals between the observed and fitted spectra in each channel and we redistribute them among the fitted lines, proportionally to their intensity in this channel. We finally construct NHI column-density and WCO intensity maps of a specific cloud by integrating the fitted lines that have their velocity centroid within a 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. To determine boundaries between clouds in the position-velocity space (l, b, v), we have used the velocity centroids of all the fitted lines and their density in the (l, b, v) space to search for line clusters (i.e., clouds). To do so, we have used the Density-Based Spatial Clustering of Applications with Noise (DBSCAN) algorithm (Ester et al. 1996) to isolate arbitrarily-shaped clusters and to draw contours of minimum density around them. The six H iand CO regions used in this analysis have been delimited by following those curves of minimal number density in (l,b,v) space. The cuts follow piecewise continuous lines in (b,v) over specific longitude intervals. The latter and the points used to define the broken lines are given in Table A.1.

Appendix B: Free-free emission in the 70 GHz data

In order to separate the free-free intensity in the Planck LFI map at 70 GHz from non-gaseous emission components, we have proceeded as follows. We have first subtracted the CMB intensity as modelled at 70 GHz by Bobin et al. (2016) from the joint analysis of WMAP and Planck observations. We have then taken advantage of the sparse spatial correlation between the thermal dust intensity I70 at 70 GHz and the dust optical depth τ353 at 353 GHz in the wavelet domain to filter out the dust component. To do so, we have transformed both maps with eight scales of the B-spline wavelet with the mr-transform tool of Starck & Pierre (1998). For each scale, we have obtained the best (χ2) linear regression between the wavelet distributions of I70 and τ353 in the regions of both significant dust optical depth (τ353 > 5 × 10-6) and low free-free intensity ( < 0.3 mKCMB at 22 GHz). We have used this regression to remove the part linearly correlated with τ353 from the 70 GHz transform at each scale. After this correction, we have applied the inverse transform to reconstruct the filtered map. We have then used the mr-detect tool with the same eight wavelet scales to detect and remove point sources. Only the first scale was used for source detection. The final map has been translated into units of specific intensity (Jy/sr) using the conversion coefficients3 for the 70 GHz band. We have applied the colour correction interpolated for the − 0.14 spectral index expected for free-free emission at 70 GHz.

Appendix C: H I spin temperature

thumbnail Fig. C.1

Evolution of the log-likelihood ratio of the γ-ray fit in the 0.4–100 GeV band as a function of the H ispin temperature used to calculate the NHI column densities in the clouds.

Open with DEXTER

The evolution of the maximum log-likelihood of the γ-ray fit in the 102.6−5 MeV band can be used to constrain the choice of average H ispin temperature that best matches the NHI column density maps of the different clouds to the structure of the γ rays produced by CR interactions in the atomic gas. We find that the fit quality peaks near 400 K (Fig. C.1). We have not attempted

to use different spin temperatures for the different clouds present in the analysis region.

Appendix D: Best fit coefficients of the γ-ray and dust models

The best fit corresponds to an HI spin temperature of 400 K. For each parameter we first give the statistical uncertainty from the fit and jackknife tests, then the standard deviation obtained by varying the H ispin temperature.

Table D.1

Best-fit coefficients of the γ-ray model in each energy band.

Appendix E: Previous γ-ray measurements of the XCO conversion factor

Table E.1

Estimates of the average XCO factor in units of 1020 cm-2 (K km s-1)-1.

All Tables

Table 1

XCO factors in 1020 cm-2 K-1 km-1 s for the dust and γ-ray fits.

Table 2

Average dust opacities, τ353/NH, and NH/E(BV) ratios in the gas phases of the different clouds.

Table A.1

Limits of the H iand CO clouds in longitude, latitude, and velocity.

Table D.1

Best-fit coefficients of the γ-ray model in each energy band.

Table E.1

Estimates of the average XCO factor in units of 1020 cm-2 (K km s-1)-1.

All Figures

thumbnail Fig. 1

Maps of the Cetus (Cet), Main Taurus (TauM), South Taurus (TauS), North Taurus (TauN), Perseus (Per), and California (Cal) clouds, and of the Galactic disc background (Gal) in H icolumn density, NHI, for a spin temperature of 400 K, and in CO intensity, WCO. The NHI and WCO maps are respectively labelled HI and CO. The map labelled ff shows the intensity Iff of free-free emission toward this region. The last subplot shows the coverage of the GALFA (dark grey) and EBHIS (light blue) H i data.

Open with DEXTER
In the text
thumbnail Fig. 2

Maps of dust optical depth at 353 GHz (top left), radiance (top right), spectral β index (bottom left), and colour temperature (lower right), from Planck Collaboration XI (2014).

Open with DEXTER
In the text
thumbnail Fig. 3

Top: γ-ray counts of gaseous origin recorded in the 0.4–100 GeV energy band in a 0.̊125 pixel grid. γ-ray emissions other than due to cosmic-ray interactions in the gas have been subtracted. The map has been smoothed with a Gaussian kernel of 0.̊14 dispersion for display. Bottom: dust optical depth measured at 353 GHz and displayed at the Fermi-LAT angular resolution for comparison.

Open with DEXTER
In the text
thumbnail Fig. 4

Photon yields arising in the γ-ray model in the 0.4–100 GeV band from the H i(HI label) and CO-bright (CO label) phases of the Cetus, Main Taurus, South Taurus, North Taurus, Perseus, California, and Galactic disc clouds, from the ionised gas (ff label), from the DNM and COsat gas column densities, from the Galactic IC emission, and from the isotropic background (iso label) and γ-ray point sources.

Open with DEXTER
In the text
thumbnail Fig. 5

Number distribution of the dust model coefficients over 1000 jackknife fits for an H ispin temperature of 400 K. yHI,i, yCOsat, and yDNM are in units of 10-26 cm2, yCO,i in 10-6 K-1 km-1 s, yff in 3.8 10-11 Jy-1 sr and yiso in 10-6.

Open with DEXTER
In the text
thumbnail Fig. 6

Number distribution of the γ-ray model coefficients obtained for 1000 jackknife fits for an H ispin temperature of 400 K. qCO,i are in units of 1020 cm-2 K-1 km-1 s, qff in 3.8 × 1015 cm-2 Jy-1 sr, qDNM and qCOsat in 1025 cm-2. qHI,i, qIC and qiso are normalisation factors.

Open with DEXTER
In the text
thumbnail Fig. 7

Residual maps (data minus best-fit model, in sigma units) in dust optical depth (στ353) and in γ-ray counts in the 0.4–100 GeV band (σγ) obtained when including different sets of gaseous components in the models: without (left) and with (middle) the DNM and COsat components in addition to the H i, CO, and free-free templates; a close-up view (right) of the γ-ray residuals without (top) and with (bottom) the free-free map. The grey circles in the left-hand column highlight a DNM shell. The thin black contours outline the shape of the CO clouds at the 7 K km s-1 level chosen to separate DNM and COsat components. The thick black contours in the right-hand column outline the shape of the NGC 1499 and G159.6-18.5 H iiregions at the 2 × 104 Jy sr-1 level in free-free emission at 70 GHz. The models are derived for an H ispin temperature of 400 K.

Open with DEXTER
In the text
thumbnail Fig. 8

Hydrogen column density maps in the DNM and COsat components, derived from the dust (top) and γ-ray (bottom) data. The demarcation between the two components is outlined by the white contour at 7 K km s-1 in CO line intensity.

Open with DEXTER
In the text
thumbnail Fig. 9

Spectral evolution, relative to the local interstellar spectrum qLIS, of the γ-ray emissivities of the local gas components. qCO is given in units of 1020 cm-2 (K km s-1)-1 and qDNM in 1025 cm-2. The qHI normalizations are given for a spin temperature of 400 K. The dashed lines compare the results obtained for the whole 0.4–100 GeV band to the results of the four independent bands. All error bars are symmetrical; horizontal error bars at the boundaries are truncated.

Open with DEXTER
In the text
thumbnail Fig. 10

Distribution with Galactocentric radius (top) and altitude above the Galactic disk (bottom) of the 0.4–10 GeV emissivities measured in the atomic gas of nearby clouds for H ispin temperatures between 125 and 150 K. The solid and grey band respectively give the average emissivity and ±1 rms dispersion in the sample. The dashed line marks the average emissivity measured across the sky at | b | ≥ 7° (Casandjian 2015).

Open with DEXTER
In the text
thumbnail Fig. 11

Evolution of the XCO factor, as measured in γ rays for various linear resolutions in the different clouds. The open circle marks the value obtained in the overall 0.4–100 GeV band, in close agreement with the weighted average of the four independent energy bands (black line) and its 1σ errors (dashed lines).

Open with DEXTER
In the text
thumbnail Fig. 12

Evolution of the XCO factor measured in γ rays as a function of the average WCO intensity, , the surface fraction of dense gas, SFdense, the mean visual extinction in the CO-bright phase, , the average CO line width, , the H2 mass in the CO-bright phase, MH2, and the velocity dispersion of the CO line. Dashed lines give the best (χ2) linear regressions.

Open with DEXTER
In the text
thumbnail Fig. 13

Ratio of XCO measurements with dust and γ rays in clouds of different average dust optical depth .

Open with DEXTER
In the text
thumbnail Fig. 14

Average dust opacities in the DNM measured in γ rays for different linear resolutions in the clouds. The open circle marks the measurement in the 0.4–100 GeV energy band, in close agreement with the weighted average of the four independent energy bands (black line) and its ± 1σ errors (dashed lines).

Open with DEXTER
In the text
thumbnail Fig. 15

Upper: 2D histogram of the correlation between the total gas column density, NH γ, measured with the 0.4–100 GeV interstellar γ rays, and the dust optical depth at 353 GHz, convolved with the LAT response for an interstellar spectrum. The maps were sampled at a 0.̊375 resolution. Lower: evolution of the dust opacity in NH bins. The error bars give the standard errors of the means and the grey band gives the standard deviation of the opacities in each bin.

Open with DEXTER
In the text
thumbnail Fig. 16

Spatial variations of the dust opacities (left) and specific power (right) with the total gas measured by NH γ at 0.̊375 resolution (top) and by NH mλ at 0.̊125 resolution (bottom). The tilded quantities are convolved with the LAT response for an interstellar spectrum. The black contours outline the shape of the CO clouds at the 7 K km s-1 level chosen to separate DNM and COsat components.

Open with DEXTER
In the text
thumbnail Fig. 17

Evolution of the dust opacities as a function of the molecular fraction in the gas column (left), the total gas column density NH (middle) and the dust color temperature (right). To estimate the molecular fraction, the DNM is assumed to be 50% molecular. The error bars give the standard errors of the means and the grey band gives the standard deviation of the opacities in each bin.

Open with DEXTER
In the text
thumbnail Fig. 18

Evolution of the dust opacities with the total gas measured by NH γ at 0.̊375 resolution (top) and by NH mλ at 0.̊125 resolution (bottom) as a function of the dust color temperature (x-axis) and spectral β index (y-axis)

Open with DEXTER
In the text
thumbnail Fig. C.1

Evolution of the log-likelihood ratio of the γ-ray fit in the 0.4–100 GeV band as a function of the H ispin temperature used to calculate the NHI column densities in the clouds.

Open with DEXTER
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.