Free Access
Issue
A&A
Volume 587, March 2016
Article Number A73
Number of page(s) 26
Section Extragalactic astronomy
DOI https://doi.org/10.1051/0004-6361/201527746
Published online 19 February 2016

© ESO, 2016

1. Introduction

Tightly connected, dust and gas are the key ingredients of star formation, which is governed by their complex, mutual interplay as in a cyclic dance. Stars form in cold, dense molecular clouds and dust works as catalyst in transforming atomic hydrogen into molecular hydrogen (e.g., Wolfire et al. 1995). Gas is again expelled from stars during their lifetime in the form of winds and finally at the supernova (SN) stage. Dust is believed to be mainly produced in the envelopes of asymptotic giant branch (AGB) stars and at the end of the life of massive stars during the explosive SN phase (see Dunne et al. 2011, for a short summary). Supernovae shocks, on the other hand, destroy dust grains, which can form again in the interstellar medium (ISM) by an accretion process. Dust absorbs the ultraviolet (UV) emission of young stars, allowing gas to cool and condense to form new stars.

The dust and gas content of galaxies are linked to each other through metallicity, as shown in the local Universe (Leroy et al. 2011; Draine et al. 2007): the gas-to-dust mass ratio in the ISM increases as a function of metallicity. Dust has therefore been often adopted as a proxy to gas in absence of time-expensive, sub-millimeter spectroscopic observations (e.g., Magdis et al. 2012; Santini et al. 2014; Scoville et al. 2014; Genzel et al. 2015) at high redshift also, assuming that local relations hold. The energy absorbed by dust at short wavelengths is re-emitted in the infrared (IR) and sub-millimeter (sub-mm) regimes, where the thermal emission of grains dominates the spectral energy distribution (SED) of galaxies (~81000 μm). Dust emission is frequently used to trace the ongoing rate of star formation (SFR; e.g., Nordon et al. 2010; 2012; Elbaz et al. 2011; Kennicutt 1998a).

It is thus no wonder that the interest in the study of galaxy dust and gas properties has grown in the past two decades and that an always increasing number of studies based on Herschel data are dedicated to the dust emission of galaxies at all redshifts. The Herschel satellite (Pilbratt et al. 2010) provides far-infrared (FIR) photometric observations of local and distant galaxies thanks to its two onboard cameras, the Photodetector Array Camera and Spectrometer (PACS, 70160 μm, Poglitsch et al. 2010) and the Spectral and Photometric Imaging REceiver (SPIRE, 250500 μm, Griffin et al. 2010). The availability of a good SED coverage from 24 to 500 μm, obtained by combining Herschel and Spitzer observations, allows us to apply different types of galaxy emission models with the aim of deriving the main properties of dust, such as its mass, temperature, and, to some degree, composition. Different families of models, with different levels of complexity, can be broadly recognized: a) full radiation transfer treatments (e.g., Silva et al. 1998; Efstathiou et al. 2000; Piovan et al. 2006; Siebenmorgen & Krügel 2007) requiring geometric assumptions, very detailed observed SEDs, and large amounts of computational power; b) evolutionary and mixed stellar population synthesis, including FIR emission (e.g., da Cunha et al. 2008); c) physically motivated FIR dust emission models (e.g., Draine & Li 2007; Galliano et al. 2011); d) simple or multiple modified blackbody (MBB) SED fitting; e) template families following locally calibrated relations (e.g., the LT relation, Chary & Elbaz 2001; Dale & Helou 2002); and f) semiempirical template fitting (e.g., Polletta et al. 2007; Wuyts et al. 2011a; Berta et al. 2013b).

Focusing on dust emission, of particular relevance and very popular are types c) and d), i.e., SED fitting with physically motivated models, comprising several dust components or MBB fitting under the assumption that dust emissivity can be described as a simple frequency power law and dust is in thermal equilibrium.

Dunne et al. (2011) applied MBB SED fitting to z< 0.5 H-ATLAS (Eales et al. 2010) star-forming galaxies selected at 250 μm and studied the evolution of the dust mass function over the last 5 billion years, showing that the dust/stellar mass ratio was three to four times larger at z = 0.4−0.5 than today. At these redshifts, PACS and SPIRE data cover the FIR peak of the SED and part of its Rayleigh-Jeans tail (RJ), thus leading to reliable Mdust estimates.

Cortese et al. (2012) combine Herschel photometry and multiwavelength data to build detailed SEDs of nearby galaxies from the Herschel Reference Survey (HRS; Boselli et al. 2010) and fit them with MBB and Draine & Li (2007, DL07) models. In combination with radio data, these authors study the dust-to-Hi and dust-to-stellar mass ratios of their targets, finding that the former increases as a function of stellar mass, while the latter tends to decrease as a function of stellar mass. Ciesla et al. (2014) apply DL07 modeling to HRS SEDs and have built a new set of reference, local SED templates.

Dale et al. (2012) and Bianchi (2013) apply DL07 and MBB fitting to the FIR-submm SEDs of KINGFISH (Kennicutt et al. 2011) local galaxies. The former authors show evidence of sub-mm (500 μm) flux excess in dwarf galaxies with respect to the expectation based on the extrapolation from Spitzer data. They also find a factor ~2 difference between DL07- and MBB-based dust masses, while Bianchi (2013) shows that this difference is avoided if using consistent assumptions for the two approaches.

Rémy-Ruyer et al. (2013) combine KINGFISH and other samples of nearby galaxies observed by Herschel, including direct measurements of molecular gas mass, to study the gas/dust vs. metallicity relation over a 2 dex metallicity range. Dust masses are in this case estimated with Galliano et al. (2011) models. Hunt et al. (2015) come back to the KINGFISH sample and present the surface brightness profiles of these nearby galaxies. These authors again study the spatially resolved SEDs and resolved properties of dust with DL07 and MBB modeling.

Eales et al. (2012) foresee the possibility of applying a single conversion factor to derive gas masses from sub-mm fluxes, based on a Milky Way (MW) observation and recalibrated on HRS nearby galaxies with CO and HI observations. They point out that this method suffers from the limitation that their sample is not necessarily representative of all galaxy populations and metallicity dependencies play an important role both for dust- and CO-based masses.

Successfully applied to local galaxies, the SED fitting techniques and scaling relations mentioned above have also been adopted by several authors to study the IR emission of distant galaxies. Broadly speaking, two main approaches can be identified: using a single sub-mm continuum observation to scale it to dust or gas mass using known relations (e.g., Eales et al. 2012; Scoville et al. 2014); or characterizing the galaxy SED including the FIR peak (e.g. Dale et al. 2012; Galametz et al. 2013). Nevertheless, when dealing with high-z galaxies, the use of these techniques is not as straightforward as in the nearby Universe. It is necessary not only to keep in mind that the properties of dust and gas might be evolving as a function of time, but also to face the limitations of the available data.

In the simplest case of all, Eales et al. (2012) discuss the limitations of SED coverage at high redshift, and show that small dust temperature gradients (see also Hunt et al. 2015) within a galaxy can lead to incorrect measurements when using a MBB fit to the global SED. Galametz et al. (2014) point out the effects of limited angular resolution on the derivation of model parameters.

Scoville et al. (2014) apply a similar concept to the RJ side of the SED of local and high-z galaxies. They suggest that gas masses can be measured using single-band sub-mm observations, adopting dust emissivities calibrated on Planck observations of the MW. They apply this method to early ALMA observations of four small samples of galaxies at four different cosmic epochs. Using simple MBB simulations, Genzel et al. (2015) discuss the systematic effects induced by this method and show that an approach involving two sub-mm bands and a proper MBB fit would be preferable.

Magdis et al. (2011) and Magnelli et al. (2012a) study the SEDs of star-forming and sub-mm galaxies that benefit from Spitzer, Herschel, and sub-mm photometry. Modeling of MBB and DL07 is applied to their SEDs with the aim of studying the dust properties of these objects and comparing these to CO observations. In this case, Magdis et al. (2011) combine dust and CO observations to obtain a measurement of the αCO CO-to-gas conversion factor. Their method is then further refined by Magdis et al. (2012) on a larger sample of sources via stacking on FIR and sub-mm maps.

Stacking is also employed by Santini et al. (2014), who bin 24 μm-selected galaxies in a grid of z-M-SFR and study the evolution of dust/stellar mass ratio and star formation efficiency (SFE) as a function of redshift (see also Genzel et al. 2015). These authors find no evidence that the Mdust evolves with redshift at a given M or SFR. They find that the SFE at z ~ 2 was a factor ~5 higher than at z ~ 0. Béthermin et al. (2015) also include sub-mm data in their stacking analysis and thus improve on the application of DL07 models to the average SEDs of distant galaxies.

Despite the wide use of SED fitting to derive dust and gas properties of distant galaxies, the underlying assumptions and implications of instrumental limitations are rarely highlighted nor discussed. In this work, we aim to show how real-world limitations affect dust and gas mass determinations based on DL07 and MBB SED fitting. To achieve our goal, we use the deepest PACS data available to date, namely the GOODS fields observations belonging to the PACS Evolutionary Probe (PEP, Lutz et al. 2011) and the GOODS-Herschel (Elbaz et al. 2011) surveys combined with SPIRE photometry from the Herschel Multi-tiered Extragalactic Survey (HerMES; Oliver et al. 2012). The analysis consists of measurements carried out on real sources (individual galaxies and stacked photometry) and Monte Carlo simulations.

The paper is organized as follows. Sections 2 and 3 present the methods adopted to derive gas and dust masses, including the details of DL07 and MBB SED fitting. In Sect. 4 the different collections of data used in the analysis are described. Monte Carlo simulations and their results are presented in Sects. 5 and 6 in two flavors: analyzing Mdust uncertainties and systematics in the z-M-SFR space and as a function of rest-frame spectral coverage. Real-world limitations and prospects for ALMA observations are discussed in Sect. 7. Finally, Sect. 8 summarizes our findings.

Throughout this manuscript, we adopt a Λ-CDM cosmology with (h0, Ωm, ΩΛ) = (0.71, 0.27, 0.73); a Chabrier (2003) initial mass function (IMF); and the Pettini & Pagel (2004) metallicity scale. Moreover, as described in the next sections, we make use of the parameterization of the M-Z relation by Genzel et al. (2015) and of the definition of reference (main) sequence of star formation (MS) by Whitaker et al. (2014). By the term gas mass, we mean the molecular gas component of galaxies unless otherwise specified.

2. Deriving gas mass

Molecular gas is generally quantified using measurements of the intensity of rotational transitions of the CO molecule and adopting a CO-to-H2 conversion factor (see Bolatto et al. 2013, for a review). Local samples include several hundred targets (e.g., Saintonge et al. 2011), but CO detections of star-forming galaxies at intermediate- or high-redshift are still limited to modest statistics and mostly to luminous galaxies (Carilli & Walter 2013, and references therein), although they increasingly reach so-called main-sequence objects (Tacconi et al. 2013). Genzel et al. (2015) collected CO line emission measurements for 484 star-forming galaxies at z = 0−3, including more than 200 z ≤ 0.05 galaxies from the COLDGASS (Saintonge et al. 2011) survey and more than 80 at z = 0.5−1.5, z = 2.0−2.5 from the PHIBSS and PHIBSS2 (Tacconi et al. 2013, and in prep.) surveys.

Alternative, faster ways to gain gas masses for distant galaxies involve the use of locally calibrated relations linking gas mass to dust mass through gas/dust ratio scaling as a function of metallicity or to SFR through the scaling of depletion times. In this work, we extensively employ the former approach, based on the gas/dust ratio, and compare the results to CO- and to tdepl-based results.

2.1. Scaling of depletion times

Based on current CO observations of local and z ~ 1 objects (e.g., Saintonge et al. 2011, 2012; Tacconi et al. 2013), and integrating the relation found by Schmidt (1959) and Kennicutt (1998b) between volume (surface) density of star formation and of molecular gas mass, the relation between the total Mgas and SFR of MS galaxies can be described as a simple scaling with a depletion timescale mildly dependent on redshift (1)referred to as the integrated Kennicutt-Schmidt relation (KS) in what follows. Based on CO observations of z ~ 1.2 and 2.2 galaxies, compared to local galaxies, Tacconi et al. (2013) derived a dependence of τdep on redshift of (2)Based on their collection of ~500 CO measurements and on dust mass derivations, Genzel et al. (2015) derived a three-parameter expression of τdep, depending on stellar mass, M, SFR, and redshift. We use the latest update of this relation (R. Genzel, priv. comm.): (3)adopting the Whitaker et al. (2014) definition of MS and the Milky Way value of the conversion factor between CO luminosity and molecular gas mass.

Here, star formation rates include both IR and UV contributions. The SFRIR,UV are computed with the Kennicutt (1998a) calibrations, scaled to the Chabrier (2003) IMF. Infrared luminosities between 8 and 1000 μm, L(IR), are derived by integrating the results of FIR SED fitting (see Sect. 3.3). An alternative L(IR) estimate was obtained by fitting SEDs with the Berta et al. (2013b) templates library, leading to equivalent results.

2.2. The dust method

As described by Magdis et al. (2011, 2012) and Eales et al. (2012), and as widely found in the literature, it is possible to derive gas masses of (distant) galaxies by relying on measurements of dust mass, Mdust, to be transformed into Mgas, once the gas-to-dust ratio (GDR or δGDR) of the galaxy is known, i.e., (4)A number of assumptions can be adopted to estimate the GDR, from simply adopting the MW value (which turns out to work well for local Virgo cluster galaxies; Eales et al. 2012), to adopting observed relations, such as the Leroy et al. (2011) local dependence of GDR on metallicity. Here we use a revisited form of the Leroy et al. (2011) relation, consistently recalibrated to the Pettini & Pagel (2004, PP04) metallicity scale by Magdis et al. (2012) as follows: (5)with a scatter of 0.15 dex. In absence of direct spectroscopic measurements, galaxy metallicity can be derived from the known stellar mass-metallicity relation and its redshift dependence. Recently, Genzel et al. (2015) combined the M-Z parameterizations of Erb et al. (2006), Maiolino et al. (2008), Zahid et al. (2014), and Wuyts et al. (2014) into a single relation (6)where a = 8.74 and b = 10.4 + 4.46 × log (1 + z)−1.78 × [log (1 + z)]2. In what follows, we adopt this relation, which provides metallicities in the PP04 scale. Alternative parameterizations or the so-called M-Z-SFR fundamental relation (Mannucci et al. 2010, 2011) give similar results, for z ≤ 2 massive star-forming galaxies with Z ~ Z, which dominate the samples discussed below, modulo the adopted metallicity scale.

We stress that it is very important to use a self consistent set of relations, i.e., the δGDR vs. Z scaling and the M vs. Z relation must be calibrated to the same metallicity scale to produce meaningful results.

3. Deriving dust mass through FIR SED fitting

Dust mass can be estimated by fitting the FIR SED of galaxies in several ways, e.g., by means of a MBB function or with dust models such as those prepared by Draine & Li (2007, DL07). In this section we provide more details about these two approaches.

3.1. Modified blackbody

The FIR SED of galaxies can be simplistically reproduced with a single-temperature modified blackbody, assuming that emission comes from a single-temperature δ distribution of dust. In this case, the rest-frame SED is modeled as (e.g., Blain et al. 2002) (7)where Ω is the solid angle of emission, ϵν is the emissivity coefficient, and Bν(Tdust) is the Planck function (8)in units of [erg s-1 Hz-1 m-2 sr-1].

For a uniform medium (e.g., a dust cloud in our case) of optical depth, τν, radiative transfer implies that radiation is reduced by a factor exp(τν) (e.g., Rybicki & Lightman 1979). Therefore, we can also write ϵν = 1−exp(τν) (e.g., Benford et al. 1999; Omont et al. 2001). We allow ν0 to be the frequency where τ = τ0 = 1. The optical depth at the generic frequency ν can be approximated as a power law (9)Therefore, (10)On the other hand, the frequency-dependent optical depth can be written as (e.g., Beelen et al. 2006) (11)where κν is the mass absorption coefficient of dust at rest frequency ν and ρ is the total mass density. If the emission is optically thin (i.e., τ ≪ 1), we can expand ϵν in its Taylor series expression1 and thus obtain (12)Therefore, the emergent luminosity from a given dust mass is given by (13)where we have assumed that emission is isotropic over a spherical surface and κν = κ0(ν/ν0)β. For a galaxy FIR dust emission, β is in the range 1.5−2.0 (e.g., Magnelli et al. 2012a, and references therein). Here Mdust is in units of [kg], κν is given in [m2 kg-1], and therefore Lν turns out to have units of [erg s-1 Hz-1].

Making the dependence of frequency explicit, Eq. (13) can be rewritten as (14)which is the usual simplified form of the modified blackbody function.

Inverting Eq. (13), it is possible to obtain an estimate of dust mass, once the rest-frame luminosity of a galaxy at a given FIR wavelength is known, a MBB fit to the data is possible, and a value of κν (or κ0 and ν0) for dust is adopted, as (15)This holds when working at rest-frame frequencies and needs to be transformed into an observed frame in the usual way. Finally, combined with the knowledge of the gas-to-dust ratio, δGDR, Eq. (15) leads to an estimate of gas mass.

3.2. The use of κν

To compute dust masses through a modified blackbody fit, one needs to assume a value of κν at a given wavelength, possibly within the rest-frame wavelength range covered by actual photometric data, or assume a value of κ0 at the reference frequency ν0, where τν = τ0 = 1.

One viable approach is to adopt the set of κν computed by Li & Draine (2001, see their Table 6), and either use the value appropriate for the given rest-frame wavelength from their Table (or interpolated) or apply Eq. (9).

It is common (e.g., Magdis et al. 2012; Magnelli et al. 2012b) to follow the second approach, using the values of κν tabulated by Li & Draine (2001) and applying κν = κ0(ν/ν0)β. In this procedure, it is important to properly apply k correction to all terms in Eq. (15) and use the correct values of κ0, ν0, and β. In fact, Eq. (9) assumes that the slope β used while fitting the observed SEDs (see Eq. (14)) is the same β adopted in the computation of κν.

Often, β is let free to vary while fitting, if enough datapoints are available, or if there are not enough available, a fixed value (e.g., β = 1.5) is adopted. On the other hand, the value of β adopted by Draine & Li (2007) is ~2 (see Cortese et al. 2012) and the values of κν in Table 6 of Li & Draine (2001; in its electronic version with the old 2001 computation and for the MW dust mixture with qPAH = 0.47) turn out to show βLD01,a = 2.07 between λ = 100 and 600 μm, and βLD01,b = 1.68 at λ> 600μm.

The direct consequence of the mismatch between the adopted value of β and the β ~ 2 describing the actual Li & Draine (2001)κν is that one should apply an additional correction factor when applying Eq. (9) (i.e., Eq. (6) in Magdis’ paper and Eq. (4) in Magnelli’s paper) or there will be a dependence of Mdust on the value of ν0 used. This dependence on frequency has a power of roughly −0.5 when adopting β = 1.5.

Bianchi (2013) points out this effect and also shows that using the κν values tabulated by Draine (2003), i.e., those actually used in the Draine & Li (2007) models, is more appropriate and leads to consistent results when comparing MBB-based and DL07-based dust masses. Also, Bianchi points out that adopting a different value of β ≠ 2.08 leads to a dependence of model normalizations and, hence, of Mdust, on the reference wavelength adopted when computing κν starting from κ0. Finally, Bianchi (2013) shows that after accounting for the needed corrections on κν, a mild dependence on blackbody shapes is still left because that the best-fit dust temperature also depends on β (see also Bianchi et al. 1999).

In practice, these subtle differences in the treatment of κν and β have led to MBB-based dust masses changing by a factor 35 for the same SED and the same basic reference dust opacity (see also Sect. 7.2).

3.3. Draine & Li (2007) models

A more sophisticated and physically motivated approach to FIR SED fitting and derivation of dust masses is to adopt the Draine & Li (2007, DL07) dust models, which are an upgrade of those originally developed by Draine & Li (2001), Li & Draine (2001), and Weingartner & Draine (2001).

In short, these models describe interstellar dust as a mixture of carbonaceous and amorphous silicate grains, whose size distributions are chosen to mimic the observed extinction law in the MW, Large Magellanic Cloud (LMC), or Small Magellanic Cloud (SMC) bar region. As described by Draine & Li (2007), carbonaceous grains have the properties of polycyclic aromatic hydrocarbons (PAH) molecules and ions when the effective radius is a< 5.0 [nm], the properties of graphite spheres when a 10 [nm], and optical properties intermediate between those of PAH particles and graphite particles for 5 <a< 10[nm]. The ionization fraction xion(a) of the PAH particles is assumed to be the average for the diffuse ISM. The properties of grains are parameterized by the PAH index qPAH, which is defined as the fraction of dust mass in the form of PAH molecules.

The dust distribution is divided in two components: the diffuse interstellar medium (ISM), usually constituting the majority of the dust budget and dust enclosed in photodissociation regions (PDRs). The former is heated by a radiation field of constant intensity Umin. The latter, representing a fraction γ of the total amount of dust, is exposed to starlight with intensities ranging from Umin to Umax. Although PDRs usually provide a small fraction of the total dust mass, they can contribute to the majority of dust emission in mid-IR SEDs. The mass dM of dust exposed to starlight intensities in the range [U,U + dU] is given by the following power-law distribution: (16)The diffuse ISM is modeled by setting Umax = Umin.

Summing the diffuse and PDR components, we obtain (17)for α ≠ 1; here δ is the delta function representing the diffuse interstellar radiation field of intensity U = Umin = Umax.

Using the numerical methods described in Draine & Li (2001), Li & Draine (2001), and Draine & Li (2007), temperature distribution functions are computed for all particles small enough for quantized heating to be important; large grains are treated as having steady-state temperatures determined by starlight heating and radiative cooling equilibrium. With temperature distributions and dust absorption cross sections, the time-averaged IR emission for a given grain size and type is computed, and finally the power radiated per unit frequency and unit mass for the given dust mix exposed to the radiation field U, is computed summing over all grain types and sizes. Integrating between Umin and Umax, this gives the specific power per unit mass pν(qPAH,Umin,Umax) (see Eq. (9) in Draine et al. 2007).

Draine & Li (2007) thus propose to fit the emission of galaxies through a linear combination of emission from diffuse ISM dust and PDRs emission. The emission spectrum, expressed as emissivity per hydrogen nucleon, , can then be simply written as (18)in units of [erg s-1 Hz-1 sr-1 H-1]; jν is the quantity contained in the model files available online2.

As shown by Draine et al. (2007), the total dust luminosity is given by (19)where the dust-weighted mean starlight intensity, or mean radiation field, scale factor is (20)if α = 2 as in our assumptions and P0 is the power absorbed per unit dust mass in a radiation field U = 1.

In principle, the model includes six free parameters: qPAH; Umin; Umax; α; γ; Mdust. Dust mass is basically the model’s normalization. In practice, studying local galaxies in the Spitzer Nearby Galaxy Survey (SINGS) and fitting their FIR SEDs, Draine et al. (2007) demonstrate that some of the parameters can be limited to a restricted range of values or even fixed to a single value. The overall fit is not very sensitive to the dust model adopted (MW, LMC, SMC), and can be limited to MW models alone, implying only seven values of qPAH. Also, the actual values of α and Umax do not significantly influence the performance of models, and fixing them to α = 2 and Umax = 106 is a reasonable choice to successfully reproduce SINGS galaxies. Finally, since small values of Umin correspond to dust temperatures below 15 K that cannot be constrained by FIR photometry alone, in the absence of rest-frame sub-mm data, they suggest limiting Umin to the range 0.7 ≤ Umin ≤ 25. A direct consequence is a possible underestimate of dust mass, if a significant amount is stored in a very cold component, but Draine et al. (2007) conclude that omitting sub-mm data increases the scatter on Mdust by 50%, but does not induce any significant systematic offset. Most observational papers in the literature (e.g., Magdis et al. 2012; Magnelli et al. 2012a; Santini et al. 2014) adopt this optimized choice of parameters, relying on a setup that has only been thoroughly verified for local galaxies near solar metallicity (Draine et al. 2007).

thumbnail Fig. 1

Effect of varying qpah (left), Umin (center), γ (right), in Draine & Li (2007) models, independently. On top of each panel, the values of the frozen parameters are quoted. The orange squares at the bottom indicate the position of the chosen rest-frame bands used to generate artificial catalogs.

Open with DEXTER

Figure 1 shows how varying each parameter (qpah, Umin, γ) independently, and fixing the others, induces modifications on the shape of the modeled SED. We also indicate also a set of IR rest-frame bands, which are available in surveys of nearby or distant galaxies.

The choice of the actual PAH abundance, parameterized by qPAH and limited to MW models, mostly influences the short-wavelength SED (120 μm), and is not making a difference at λ> 20μm. Therefore FIR surveys of galaxies are not very sensitive to this parameter and cannot put tight constraints on its value.

The parameter Umin regulates the range of starlight intensities that are heating dust. Since the amount of dust exposed to starlight intensities in between U and U + dU is modeled as a power law Uα, the smaller Umin, the larger is the amount of dust subject to low-energy radiation. As a consequence, the dust component dominating the SED is colder for smaller values of Umin and the FIR peak shifts to longer wavelengths; at the same time, at a given wavelength and for a given total dust mass, the model intensity decreases. Thus, we expect degeneracies and possibly some amount of systematics between Mdust and Umin.

The fraction γ of dust locked in PDRs influences the repartition of the emitted energy between the Umax = Umin component (ISM) and the rest. When fitting the SEDs of local SINGS galaxies, Draine et al. (2007, see their Table 4) found that the SEDs of vast majority of galaxies were successfully reproduced with γ ≤ 0.15 with only a couple of cases extending to γ ~ 0.3−0.4.

Once the FIR SED is fitted with these models, its dust mass can be computed as (21)where the right-hand term is evaluated at a given rest-frame frequency ν of choice, which is at best covered by observations; Lν is the rest-frame luminosity at frequency ν; mH is the mass of the hydrogen nucleon; and Mdust/MH is the dust-hydrogen ratio characteristic of the adopted dust model (and tabulated together with qPAH). Making sure to use the right units for all terms, we use mH = 1.67262178 × 10-27[kg], Lν in units of [erg s-1 Hz-1], and therefore dust mass turns out to have units of [kg].

4. Available data

The PEP survey (Lutz et al. 2011) has covered the most popular blank fields with the PACS (Poglitsch et al. 2010) instrument aboard Herschel (Pilbratt et al. 2010) at 100 and 160 μm. In addition, the same fields have been observed by the HerMES survey (Oliver et al. 2012) with SPIRE (Griffin et al. 2010) at 250, 350, 500 μm. Furthermore, the GOODS-Herschel (Elbaz et al. 2011) survey has provided deeper PACS coverage of the GOODS fields, reaching the confusion limit at 100 and 160 μm. Finally, the PEP team has combined all the available PACS data in the GOODS-N and GOODS-S fields, including science verification (SV) and ECDFS observations, producing the deepest FIR maps obtained so far (Magnelli et al. 2013).

Here we make use of PEP first data release (DR13, Lutz et al. 2011), the HerMES DR2 and DR34 (Oliver et al. 2012; Roseboom et al. 2010, 2012), and combined PEP + GOODS-H data (included in PEP DR1, Magnelli et al. 2013).

4.1. Individual sources

The deepest PACS fields, GOODS-N and GOODS-S benefit from an extensive coverage at all wavelengths. We use the PACS 70, 100, 160 μm data presented by Magnelli et al. (2013); the SPIRE HerMES data (Oliver et al. 2012; Roseboom et al. 2010); the multiwavelength catalogs built by Berta et al. (2011) and by the MUSIC team (Grazian et al. 2006); the collection of spectroscopic redshift by Barger et al. (2008), Balestra et al. (2010), and Berta et al. (2011); and the photometric redshifts by Berta et al. (2011) and Wuyts et al. (2011a). We refer to the dedicated publications for more details on each dataset and to Berta et al. (2011, 2013b) for an overview.

We select objects with at least a 3σ detection in the PACS 160 μm band, which for the redshift range of interest turns out to be the best single-band proxy of star formation rate (Elbaz et al. 2011; Nordon et al. 2013).

To this generic sample of ~1400, 160 μm sources, we add the list of 61 sub-mm galaxies (SMGs) compiled by Magnelli et al. (2012a) and distributed in the main PEP fields. For these sources, a rich multiwavelength dataset is available, ranging from the optical to the sub-mm as compiled by Magnelli et al. (2012a).

4.2. Stacked photometry

Magnelli et al. (2014) studied the dust temperature, Tdust, in star-forming galaxies as a function of their position in the z-M-SFR space.

These authors used the M, SFR, z estimates by Wuyts et al. (2011a,b) in the GOODS-N, GOODS-S, and COSMOS fields. Stellar masses are based on optical-NIR SED fitting carried out adopting Bruzual & Charlot (2003) templates. Star formation rates are based on a ladder of indicators (SED fitting, mid-IR photometry, FIR photometry), calibrated on Herschel observations of ~7000 galaxies in PEP fields. Redshifts are the combination of a collection of spectroscopic measurements and photometric estimates (see above).

Magnelli et al. (2014) produced a grid in the M-SFR-z space, binned a Ks-selected sample accordingly, and finally produced an average SED (70, 100, 160, 250, 350, 500 μm) for each bin by means of stacking on Herschel maps.

We make use of the stacked photometry by Magnelli et al. (2014), plus similar Spitzer/MIPS 24 μm data, taking care to limit our analysis to those regions of the M-SFR-z space, where good stacked SEDs are available (see Figs. 4 and 5 in Magnelli et al. 2014). Table 1 summarizes the z, M, sSFR ranges of interest.

4.3. CO-detected sources

Table 1

Regions of the M-SFR-z space where good stacked SEDs are available (see Figs. 4 and 5 in Magnelli et al. 2014).

A direct measurement of Mgas comes from the intensity of spectral lines of molecular gas tracers, such as CO, modulo the conversion factor from the given molecular species to total gas. Here we collect the FIR and sub-mm photometry of galaxies that were observed in CO-line spectroscopy to compare CO-based and dust-based Mgas estimates directly.

The following samples of CO-detected objects are taken into account:

Appendix A describes the publicly available data for each case. In summary, the number of CO-detected objects that have enough FIR photometry to allow a DL07 SED fitting and that have a M estimate is: 23 sources in EGS from the PHIBSS catalog, 5 BzK galaxies by Daddi et al. (2010), six sources by Magnelli et al. (2013), ten lensed galaxies by Saintonge et al. (2013), and seven SMGs by Bothwell et al. (2013); this is a total of 51 objects.

Stellar masses of these CO-detected sources come from the respective publications dealing with each sample. They were all computed with Chabrier (2003) IMF and BC03 models or are consistent with these assumptions. For SMGs, the adopted stellar masses are similar to the scale of Hainline et al. (2011). The adopted cosmological parameters are all in line with those assumed in this work.

5. Accuracy of Mdust and Tdust, as allowed by Herschel surveys

thumbnail Fig. 2

Relative uncertainty of dust masses (defined as σ(Mdust) /Mdust), based on DL07 SED modeling, as a function of position in the z-M-SFR space. The two diagrams belong to two independent simulations, obtained noise levels of individual detections in GOODS-S and COSMOS, applied to artificial photometry based on the SEDs by Magnelli et al. (2014). For reference, black lines denote the position of the main sequence (MS, solid line) of star formation (Elbaz et al. 2011) and ± 4, ± 10 MS levels (dashed).

Open with DEXTER

Table 2

Herschel PACS and SPIRE 3σ depths in the GOODS-S and COSMOS fields, adopted in Monte Carlo simulations.

Magnelli et al. (2014) provide a tool to derive Tdust of a galaxy once its redshifts, M and SFR are known, as well as to obtain its expected SED based on the Dale & Helou (2002) template library. Genzel et al. (2015) re-analyze the data by Magnelli et al. (2014) and produce new scalings linking Tdust and Mgas to M, sSFR, and z.

In this section, we would like to study how the observational limitations of Herschel photometry affect the derivation of Mgas and Tdust. To this aim, we adopt the signal-to-noise ratios (S/Ns) of two case studies: the deepest field of the PEP survey, i.e., GOODS-S; and the widest, but ~8 times shallower, field of PEP, i.e., COSMOS. In what follows, we produce a list of artificial sources, characterized by Tdust, Mgas, and SEDs given by the relations by Magnelli et al. (2014) and Genzel et al. (2015), and fit them with MBB and DL07 models. We limit the analysis to the parameter ranges listed in Table 1. Our simulation is structured as follows:

  • 1)

    We adopt the z-M-sSFR grid by Magnelli et al. (2014) and we limit the analysis to the range of parameters over which the relation by Magnelli et al. (2014) holds (see Sect. 4.2).

  • 2)

    The recipe by Magnelli et al. (2014) produces a value of Tdust given the position in the z-M-sSFR grid along with the typical FIR SEDs of a galaxy in that position. These SEDs are based on (and limited to) the Dale & Helou (2002) templates library.

  • 3)

    These SEDs are convolved with photometric filters: MIPS 24 μm; PACS 70, 100, and 160 μm; and SPIRE 250, 350, and 500 μm.

  • 4a)

    In the case of the DL07 simulation, the convolved SEDs are fit with the DL07 models, adopting a 1% photometric uncertainty in all bands. In this way, a so-called input catalog is produced, defining the DL07 parameters and a reference Mdust to be associated with each z-M-SFR bin.

  • 4b)

    In the case of MBB simulation, Genzel et al. (2015) provide scaling relations to define the input value of Tdust and Mgas (to be transformed into Mdust adopting a gas/dust mass ratio) at any position in the M-sSFR-z space. These relations are calibrated on Magnelli et al. (2014) data, which hold for β = 1.5. Thus, using the relation by Magnelli et al. (2014) leads to similar results. In case a value β = 2.0 was adopted, temperatures need to be increased by 4 K (Magnelli et al. 2014) and dust masses corrected as discussed in Sect. 3.2.

  • 5)

    Real noise levels are associated with each band (see Magnelli et al. 2013; Berta et al. 2013b). We use the PEP/GOODS-H/HerMES noise values for individual detections of two different fields: GOODS-S (PEP plus GOODS-H and HerMES depths; see Magnelli et al. 2013; Oliver et al. 2012) and COSMOS (PEP and HerMES depths; see Lutz et al. 2011; Oliver et al. 2012). Table 2 lists the adopted depths. Two independent simulations are run: one for each set of depths. Only bands with S/N ≥ 3 are taken into account.

  • 6)

    We fit the noisified catalog with DL07 or MBB models. In the DL07 case, photometric points at λ ≥ 8μm (rest frame) are considered; in the MBB case, only bands at λ> 50μm (rest frame) are used.

  • 7)

    Solutions are found both through χ2 minimization and through 1000 Monte Carlo (MC) realizations for each entry in the synthetic catalog and are obtained by letting the photometry vary within the “observed” noise. In the two cases, uncertainties are computed based on Δχ2 or as the dispersion of all MC realizations, respectively. The two approaches lead to comparable results. The analysis and figures presented here are based on the MC approach.

5.1. Results of DL07: Relative uncertainties on Mdust

We first verify down to what precision our procedure is able to determine dust masses as a function of position in the redshift, M, SFR space, and at the depths reached by Herschel for individual detections in the GOODS-S and COSMOS fields.

Figure 2 shows the M-SFR plane in different redshift slices, color coding each bin on the basis of its average σ(Mdust) /Mdust value. Bins indicated in red have an average Mdust/σ(Mdust) smaller than 3. It was not possible to run SED fits for bins missing at the low SFR side because too few photometric bands are available. Figure 3 collapses this diagram along the M axis for the GOODS-S case, thus allowing for a more straightforward view as a function of redshift and SFR.

As expected, the relative uncertainty on Mdust becomes worse, as redshift increases (due to Malmquist bias) and SFR decreases. Redshift plays a double role by dimming fluxes and thus raising the luminosity threshold; and by pushing the rest-frame SED coverage to shorter wavelengths (see also Sect. 6). At the sensitivity of PEP and HerMES, dust masses are retrieved with a S/N ≥ 3 for galaxies on the MS of star formation down to M ~ 1010[M] up to z ~ 1. It is possible to obtain an estimate of Mdust at comparable stellar masses up to z ~ 2 for objects lying increasingly above the MS.

5.2. Results of DL07: Systematics on Mdust

It is now possible to study how well Mdust is recovered by comparing input and output dust masses. Figure 4 presents the mere comparison of input and output Mdust for all bins in the considered z-M-SFR space. Scatter increases and the incidence of outliers becomes larger as the depth of the available bands become less balanced, i.e., in the GOODS-S simulation, where MIPS and PACS are much deeper than SPIRE.

thumbnail Fig. 3

Relative uncertainty of dust masses (defined as σ(Mdust) /Mdust) at the GOODS-S depth, based on DL07 SED modeling, as a function of redshift and SFR. Color coding is as in Fig. 2. The dotted, solid, and dashed lines indicate the position of the MS of star formation at different values of stellar mass (Elbaz et al. 2011). These values were obtained from the dependence of the sSFRMS on redshift by adopting log (M/M) = 10.0, 10.5, and 11.0, respectively.

Open with DEXTER

thumbnail Fig. 4

Comparison of input and Monte Carlo (MC) output Mdust in DL07 simulations. Color coding is based on redshift, stellar mass, specific SFR, and Mdust relative uncertainty. As each dot belongs to an individual object, there is a general overlap of colors (green possibly hiding pink). The two diagrams belong to two independent simulations obtained with GOODS-S and COSMOS depths. Data points indicated in red have Mdust/σ(Mdust) < 3.

Open with DEXTER

The color coding in Fig. 4 is based on redshift, stellar mass, specific SFR, and on Mdust relative uncertainty. Critical outliers lie at intermediate to high redshift and are characterized by low M. They turn out to be sources with few photometric points available and thus with limited wavelength coverage. At the adopted depths, it is not possible to derive Mdust with an accuracy better than 3σ for them, therefore, they are indicated in red in Fig. 4.

Ignoring these critical cases, milder systematic offsets seem still to occur for those bins with larger relative uncertainties on Mdust, which tend to lie preferentially below the 1:1 line. However, for these sources with Mdust measured at >3σ, such systematic offsets are well below a factor 2.

Figure 5 puts systematic offsets in the context of the z-M-SFR space in the case of GOODS-S. With PACS and SPIRE depths strongly unbalanced, there is a trend to underestimate Mdust at intermediate sSFR. This effect is more critical at z ≃ 0.5−1.5, where SPIRE bands are most important to constrain the peak of the FIR SED. At lower redshift, deep PACS data provide a good constraint already, and at higher redshift the peak and Mdust are likely equally poorly constrained both in the “IN” and “OUT” cases.

Working on a small sample of high-redshift galaxies, Magdis et al. (2012) showed that the presence of photometric data at rest-frame wavelengths larger that 100160 μm should avoid large systematics in the estimate of dust masses. Similar results were obtained locally by Ciesla et al. (2014). Imposing the condition to have at least one 3σ detection at λrest ≥ 160μm slightly mitigates these systematics detected in our simulations, but does not solve the problem.

thumbnail Fig. 5

Systematics on dust masses (defined as (OUT-IN)/OUT) as a function of position in the z-M-SFR space, in the simulation with GOODS-S depth. Black diagonal lines indicate the position of the MS (Elbaz et al. 2011) and ± 4, ± 10 MS levels (dashed).

Open with DEXTER

5.2.1. The choice of templates and the role of β

The results presented above in Sects. 5.1 and 5.2 were obtained by fitting, with DL07 templates, the synthetic photometry computed by convolving the Magnelli et al. (2014) SEDs with photometric filters. It is worth recalling that Magnelli et al. (2014) assigned to each Dale & Helou (2002, DH02) template a value of dust temperature, Tdust, based on an MBB fit to the template itself. This Tdust was obtained fixing β = 1.5. Then they used DH02 templates to fit the Herschel stacked photometry of each z-M-SFR bin, thus linking each position in this space to Tdust. In short, we are fitting the DH02-based synthetic photometry of each z-M-SFR bin with DL07 templates.

The DL07 models were computed by adopting a fixed slope of the dust IR emissivity, β = 2.08. On the other hand, in DH02 models, β varies as a function of the intensity of the radiation field, U, following the relation β = 2.5−0.4log U at λ> 100 μm. This intensity varies in the range 0.3 <U< 1.e5, thus one obtains β = 2.5, 2.1, 1.3 for U = 1, 100, 1000.

We now test whether the systematic trends seen in Fig. 5 and described in Sect. 5.2 could be, at least partially, driven by the different adopted models in input and output. To this aim, the DL07 simulation (Sect. 5) is now repeated with a new photometry based on DL07 models themselves. This is carried out by convolving DL07 models with photometric filters at step number 4a) in the previous scheme, and adopting this new photometry during the SED fitting process instead of the photometry computed at step number 3).

Figure 6 briefly shows the results for GOODS-S (our worst case, see Sect. 5.2), i.e., reproduces Fig. 4 with the new setup. As expected, the majority of systematics are no longer found with the exception of the low S/N region of the parameters space, where dust masses tend to be overestimated. Their relative uncertainty on Mdust exceeds 33% (red points). Similar results are obtained at the COSMOS depth.

The lesson to be learned is that the value of these kind of simulations and their ability to describe systematics is limited to the ideal case that the SED shape adopted in the simulation, and to fit the real sources is the same as that of real-life galaxies. In our case, this is exemplified by the value of β. Small differences can lead to misleading results or, more worrying, to physical overinterpretations of numerical effects. Vice versa, the absence of systematics in these kind of simulations might not always be indicative of absence of systematics on real sources if the real SEDs differ from those adopted to model them.

thumbnail Fig. 6

Same as Fig. 4, for GOODS-S only and for MC results obtained using DL07-based photometry instead of DH02-based photometry.

Open with DEXTER

5.3. Results of MBB: Constraints on Tdust

Following the scheme introduced for the DL07 simulations, we now analyze the results of MBB-based SED fitting of synthetic catalogs. The procedure followed is similar to the DL07 case and is described in Sect. 5. We fix the value of β to ease the comparison of our simulations to recent literature works and later to the analysis of individual real sources for which few photometric points are available.

thumbnail Fig. 7

Comparison of input and Monte Carlo (MC) output Tdust in MBB simulations at the COSMOS depth, with fixed β. Color coding is based on redshift, stellar mass, specific SFR, and Mdust relative uncertainty.

Open with DEXTER

thumbnail Fig. 8

Study of possible biases on Tdust in a MBB fit and dependence of DH02 templates (e.g., Magnelli et al. 2014). DH02 templates are characterized by their α parameter. Black filled squares and open circles represent the values of Tdust associated with each DH02 template by fitting it with a MBB model with β = 1.5 or 2.0, respectively. The light blue grid maps the variation of Tdust as a function of redshift for each given template. The inset includes a zoom on a single template, aimed at showing the details of the dependence of the Tdust derivation on redshift. Each redshift point was artificially shifted by a small amount along the x-axis to avoid overlapping. Black crosses refer to the result obtained by applying Wien’s displacements law to DH02 templates.

Open with DEXTER

We now deal with the case at COSMOS depth, which represents the shallowest and worst case scenario in this analysis. Figure 7 compares input and output dust temperatures. Data points in the four panels are color coded as a function of redshift, M, sSFR, and Tdust relative error. Dust temperature is retrieved within few percents in most cases. Only few catastrophic failures are recorded: they are limited to bins with very poor S/N and, thus, with large uncertainties on Tdust. Only ~4% of the cases turn out to have Tdust overestimated by more than 10%. These statistics become even better at the GOODS-S depth.

When Tdust is not correct, Mdust can also be overestimated, up to a factor of ~5. Nevertheless, in these cases, the photometric S/N is not good enough to guarantee a relative error on Mdust smaller than 3 σ.

5.3.1. Biases in MBB Tdust determination

The SEDs of real galaxies are indeed not single-temperature modified blackbodies: as seen before, the FIR emission is produced by a mixture of dust grains, of different sizes and shapes, which turn into a nontrivial distribution of temperatures. Although convenient, the MBB approximation is an over-simplification in most of the cases. More complex descriptions include using multitemperature blackbodies, DL07 modeling (see previous sections), or mixed approaches as that by Magnelli et al. (2014, described in Sect. 5.2.1. In this section we study possible biases in MBB-based Tdust due to changes in β and to redshift.

As in Magnelli et al. (2014), each template in the DH02 library is convolved with PACS/SPIRE passbands (70 to 500 μm) and the broadband photometry thus obtained is fitted with MBB models with varying dust temperature, Tdust. The dust emissivity index, β, is fixed to 1.5 and 2.0 to test different scenarios. In this way, a value of Tdust is associated with each DH02 template. This is the reference dust temperature.

Each template is then shifted to increasing redshifts, in the range z = 0−3, and convolved again with filters. The z> 0 synthetic photometry is then refit with a MBB, taking care to only use those bands with a rest-frame wavelength λrest> 50μm.

For each template, characterized by the value of the parameters α (Dale & Helou 2002), Figure 8 plots the values of Tdust obtained with β = 1.5 or 2.0 at z = 0 (black squares and open circles) and those derived at z> 0 with β = 1.5 (light blue grid). The inset includes a zoom on a single template, which is aimed at showing the details of the derivation of Tdust on redshift. Each redshift point was artificially shifted by a small amount along the x-axis to avoid overlapping.

For reference, the results of applying Wien’s displacements law to DH02 templates, thus deriving dust temperatures from the wavelength of the SED’s peak (Casey et al. 2012), are shown (black crosses).

Using different values of β produces a systematic change in the inferred value of Tdust, a smaller β implying a larger Tdust. At higher redshift, Tdust further increases, but the change in temperature does not follow a linear trend, showing a saw-toothed pattern instead. This reflects the fact that as redshift increases, the rest-frame wavelength sampled by each photometric band decreases, and the fitting procedure tends to be biased toward higher dust temperatures because the SED is constrained at shorter wavelengths. When one band falls shorter than 50 μm rest frame, it is discarded and Tdust shows a jump toward colder values because suddenly a short-λ band is missing. This effect happens periodically, every time one band is discarded as a result of the λrest> 50μm requirement.

In other words, for sources lying at different redshift, the discreteness of the sampled wavelengths and their shift to the rest frame cause a few K systematic difference in the determination of Tdust.

6. Rest-frame simulations

thumbnail Fig. 9

Relative uncertainty on DL07 free parameters, as obtained with the rest-frame simulation. For simplicity, we show only the case with α = 2, Umax = 106, Umin = 0.7, qPAH = 0.47, and γ = 0.01−0.11. The left diagram is obtained using ten photometric bands from 8 to 850 μm; the right diagrams uses only six bands between 8 and 250 μm. Left/right columns belong to the case obtained by removing long- and short-wavelength bands. Color coding is based on the average photometric uncertainty as computed over all available bands. See Sect. 6.1 for more details.

Open with DEXTER

To understand the results presented in Sect. 5, we ran an additional simulation, which is limited to rest-frame SEDs, aimed at studying the effects of losing one or more photometric bands on the short- or long-wavelength side of the SED. At the same time, we also study the effects of having different amounts of noise in the data by randomly assigning S/N values to each band, independently. This means that, for each entry in the synthetic catalog, the distribution of relative photometric uncertainties among bands is neither flat nor based on the ratio of noise levels in selected Herschel surveys, but is a random combination.

The new simulation is structured as follows:

  1. First of all, a library of models (DL07 or MBB) is generated based on a grid of input parameters (qPAH, Umin, and γ in one case, and Tdust in the other). We limit the analysis for DL07 models to the setup suggested by Draine et al. (2007). In the MBB case, the parameter β is fixed to a value of 1.5 and Tdust spans the range 1050 K. All models are renormalized to a fixed total dust mass of 108[M] and then convolved with mid- and FIR photometric bands. We use ten bands for the DL07 models: Spitzer IRAC 8 μm, IRS 16 μm, MIPS 24 μm, Herschel PACS 70, 100, 160 μm, SPIRE 250, 350, 500 μm, and SCUBA 850 μm. We also performed an additional run, limited to six bands between 8 and 250 μm. Only eight bands are used for MBB models: a square box filter centered at 40 μm and with a half width of 10 μm, the 70, 100, 160, 250, 350, 500 μm PACS bands, and the SCUBA 850 μm passband. All entries in the catalog are at redshift z = 0. This is the so-called input catalog.

  2. The input catalog is degraded in two ways:

    • An increasing number of bands is removed, progressively one by one, until only two are left. This procedure is repeated twice: first photometric bands are removed from the long-wavelength side of SEDs5; then the number of bands is reset to eight (i.e., the full SED is re-installed) and bands are removed one by one from the short-wavelength side of SEDs6.

    • The uncertainty on the photometry is increased. Relative uncertainties are increased randomly and independently for each band, spanning the range between 5% and 50% in relative errors.

In practice, the degraded catalog contains all entries of the input catalog, each one modified eight times getting rid of bands on the short-λ side and eight times getting rid of bands on the long-λ side. Each of the N × 16 entries is modified ten times with random S/N values; each time, a random level of noise is assigned to each band, independently from the other bands.

By removing bands from the long-wavelength side of the SED, we are simulating the case of deep surveys, such as the deepest fields in PEP/HerMES (Lutz et al. 2011; Oliver et al. 2012), benefiting from deep PACS photometry for the majority of sources, and progressively missing SPIRE detections because of depth, confusion, and blending effects. Vice versa, when removing bands from the short-wavelength side, we are simulating the case of shallow surveys, such as H-ATLAS (Eales et al. 2010): SPIRE quickly reaches the relatively high confusion limit at these wavelengths, but a fast PACS observation cannot fully take advantage of the smaller beam and deeper confusion limit, and remains limited to bright and mostly lower redshift sources.

6.1. Reading key

The results of the rest-frame simulations are presented with two main flavors of diagrams.

First the relative uncertainty on parameters, as obtained with Monte Carlo runs, is shown. In these diagrams, each displayed dot belongs to a full MC run per one artificial object, i.e., the value and error associated with each dot are computed as the average and standard deviation of 1000 evaluations of the same entry of the synthetic catalog. Color coding is based on the average relative photometric uncertainty of all available bands, if not otherwise specified. Darker colors refer to smaller average photometric uncertainties. The x-axis shows the maximum or minimum wavelength covered by the data, λmax (or λmin), depending on whether photometric bands were removed at the long- or short-wavelength side. In correspondence of the wavelength of each photometric band, a column of dots is plotted. These dots represent all those entries in the artificial catalog with that value λmax (or λmin), therefore, the given object in the given column benefits from available photometry in all bands shortward (longward) of the λmax (λmin) it belongs to. This is exemplified by the left arrows (right arrows) at the bottom of each diagram.

The second type of diagram (see also Appendix B.1) probes for systematic effects. This diagram presents the comparison of output derived quantities to their known input values as a function of λmax,min. Each dot again belongs to a given object in the synthetic catalogs, i.e., the dot represents the average of 1000 MC evaluations.

6.2. Results of DL07: Relative uncertainties

thumbnail Fig. 10

Relative error on dust temperature (top) and dust mass (bottom) obtained with MBB simulations, when erasing long-wavelength bands (left) and short-wavelength bands (right). Color coding is based on the average relative photometric uncertainty computed on the available bands, and is explained in the legend. See text for more details on how to read these diagrams.

Open with DEXTER

Figure 9 presents the relative uncertainty on the free DL07 parameters and Mdust as a function of λmax,min for cases obtained removing long-wavelength and short-wavelength bands. We limit the diagrams describing relative uncertainties to specific values of the input parameters to avoid the piling up of too many cases.

Trivially, objects with larger photometric uncertainties have larger relative uncertainties on the derived quantities. However, since the photometric S/N of each band is independent from the others and the color coding based on an average photometric relative uncertainty is computed over all bands, in Fig. 9 there is not a smooth gradient of colors, but a mixing of cases is present.

Unless the rest-frame mid-infrared SED is sampled by the available photometry (824 μm in the specific case), and with good S/N, constraints on the value qPAH and γ are weak. As long as the peak of the FIR SED is sampled, and if the average S/N is good (darker dots and lower envelope of the relative error distribution in Fig. 9), it is possible to have an estimate of Umin within a ~30% uncertainty. A broader peak (induced by larger input values of Umin) is more difficult to constrain and the consequent uncertainty on Umin is larger.

The large uncertainties on the γ and Umin parameters have strong consequences on the determination of U. Unless the SED is robustly constrained (see below), the propagation of γ and Umin errors onto Eq. (20) significantly hinders the estimate of U. On the other hand, U can also be derived from Ldust and Mdust (see Eq. (19)) under the assumption that the value of P0 is known (e.g., Magdis et al. 2012) and that the two quantities can be measured independently enough.

The best constrained and most stable parameter turns out to be dust mass, Mdust. When several bands are available and the λ> 250μm SED is sampled, dust mass is constrained within a 30% uncertainty even with poor photometry, but as soon as sub-mm data are missing, the maximal uncertainty (poor S/N, light dots) on Mdust immediately reaches 7080%, even with the full 8−250μm wavelength range covered. When removing short-wavelength bands, the uncertainty on Mdust remains more stable than in the MBB case (see Sect. 6.4) on both the high and low S/N sides. This might be because of the choices made on the basis of SINGS results (Draine et al. 2007), effectively limiting the freedom of choosing extreme models that would cause drastic changes of normalization. In the best case of high quality photometry, the uncertainty on Mdust can be as low as 1020% even with only three FIR bands (rest-frame 100, 160, 250 μm) available.

Upper limits are simulated with relative photometric uncertainties larger than 0.33 (see Fig. 9). If we compare the lightest blue dots in a given column to those in the subsequent column that has the given band fully removed, we find that the use of upper limits can produce an improvement in the determination of Mdust, but this improvement is nevertheless modest.

Expectations for a typical PEP galaxy at z = 1.5−2.0, and observed photometry covering the rest-frame 8160 μm (i.e., λobs = 24−350 μm) can be read from the results of simulations carried out with six photometric bands and removing long-wavelength data (Fig. 9, right). The column of dots at λ = 160μm represents results for sources with rest-frame 8160 μm photometry available. In such cases, the relative uncertainty on Mdust can be as low as ~1020% for good S/N data (darker dots and lower envelope of the dots distribution), the exact value depending on the value of input parameters. Nevertheless, it easily grows beyond 50% for 3σ-only detections. Results only get worse by a factor ~2 for λrest up to 100 μm.

6.3. Results of DL07: Systematics

Possible systematic offsets are studied by comparing output and input values. We focus on dust mass estimates as the other parameters are much more poorly constrained (see Sect. 6.2). Moreover, only a short compendium of results is presented here, while a thorough description of details is given in Appendix B.1.

thumbnail Fig. 11

Relative uncertainty on Mdust for MBB (left) and DL07 (right) fits to the FIR SEDs of GOODS-N (blue) and GOODS-S (red) sources. Color coding is based on the maximum S/N of the available photometry (at λobs ≥ 100μm) for each object, ranging from a value of 3.0 to 1000.

Open with DEXTER

When removing short wavelength bands, Mdust is usually easily retrieved and there are no very significant trends related to γin in over- or underestimating Mdust. This also holds in varying qPAH,in and Umin,in. When dealing with a limited number of bands, i.e., with MC runs using only six filters in total, systematics as function of Umin,in can be triggered. This is because not only bands at the short wavelength side are progressively missing, but also the long wavelength SED is sampled only up to 250 μm rest frame.

When removing long-wavelength photometric bands, systematics show trends as a function of the position in the qPAH,in, Umin,in, and γin parameter space. At low Umin,in there is a tendency to underestimate Mdust if the band coverage is poor, while at high Umin,in there is a tendency to overestimate it. A larger value of γin produces an increased chance to overestimate Mdust. Similarly, the larger qPAH,in, the more Mdust can be overestimated.

The three effects can add up or compensate each other. If there is a general tendency to under-estimate Mdust (e.g., because Umin,in is small), then the larger values of γin or qPAH,in mitigate it. In contrast, large values of Umin,in combined to larger γin or qPAH,in can induce a tendency to significantly overestimate Mdust when the SED coverage is poor at long wavelengths.

Generally speaking, if the maximum covered rest-frame wavelength is λmax< 160−200μm systematics on Mdust can become significant, especially when the S/N of the available FIR photometry is poor.

6.4. Results of MBB: Relative uncertainties

We now analyze the case of MBB models with fixed β and how relative errors on Tdust and Mdust depend on the available photometry and on its S/N. Figure 10 shows the behavior of relative uncertainties on Tdust and Mdust, for all entries in the artificial catalog built on MBB models without distinction of input dust temperature.

The left panel shows that, as the long-wavelength part of the SED is progressively less sampled, the uncertainty on Tdust and Mdust increases, as expected. The parameter Tdust is always constrained relatively well within a 25% uncertainty. As the average photometric uncertainty increases (brighter color tones), naturally the uncertainty on the derived quantities increases. This is particularly true for Mdust. When rest-frame sub-mm data are available, with the exception of a few outliers dust masses can generally be retrieved with relative uncertainties, within 2030% in case of poor photometry. As the maximum covered wavelength decreases, the uncertainty on Mdust explodes when the available bands do not sample the SED at λ ≥ 350 μm (restframe) anymore and the photometry is poor. This turns out to be slightly worse than DL07 modeling, which produces more stable Mdust uncertainties even with SEDs only limited to ~160 μm rest frame (see Sect. 6.2).

The right-hand panel shows instead what happens when short-wavelength bands are progressively removed. The uncertainty on Tdust is now larger than in the previous case, reflecting the fact that while the long-wavelength side of the MBB has a weak dependence on Tdust, there is a stronger dependence on the short-wavelength side. Uncertainties on Mdust are not significantly enhanced with respect to the previous case.

As the temperature increases the minimum uncertainty on Tdust tends to become larger, reflecting the fact that the peak of the MBB moves to shorter wavelengths. A similar trend is also detected for Mdust, although with smaller amplitude.

The case of a PEP z = 1.5−2.0 galaxy can be studied in Fig. 10, focusing on the left panel, at λ = 40−160μm. The observed photometry (e.g., at 100, 160, 250, 350 μm) samples the rest-frame 40200 μm wavelength range. The simulation shows that it is possible to constrain Tdust within a <20% accuracy, but the uncertainty on Mdust could easily grow up to 80% with poor S/N. In the best case, the smallest uncertainty possible (~15%) is reached if photometric accuracies of ~15−20% are available in almost all bands (350 μm included).

6.5. Results of MBB: Systematics

Here the possible systematics on dust mass are tested for the MBB modeling along with Tdust reliability. As in the DL07 case, we only summarize the main results, while we defer to Appendix B.2 for fine details.

When the photometric accuracy is poor, Mdust can be systematically overestimated if bands are missing on the long-wavelength side of the SED; the problem is relatively milder when removing bands on the short-wavelength side.

The estimate of Tdust is not very sensitive to the lack of long wavelength datapoints. On the contrary, it can be systematically underestimated or overestimated when removing data on the blue side of the SED. This happens because the peak of the MBB emission shifts to longer wavelengths as dust temperature decreases, and the effect of poor SED sampling (due to the removed datapoints) is thus amplified. Overall, only ~10% of cases turn out to have Mdust,out> 2 × Mdust,in.

thumbnail Fig. 12

Comparison of MBB and DL07 dust estimates. The MBB fitting adopts β = 2.08 and a minimum rest-frame wavelength of 50 μm in the fit. The DL07 estimates are obtained with the whole 8500 μm photometric set for a maximum of nine bands, and DL07 parameters are limited as prescribed by Draine et al. (2007). Left: individual GOODS-N/S sources detected by Herschel. Right: stacked photometry by Magnelli et al. (2014) color coded on the basis of the maximum available rest-frame wavelength (see also inset).

Open with DEXTER

7. Discussion

In the previous sections, we built expectations on DL07 and MBB SED fitting applied to Herschel-detected sources. It is now time to derive dust and gas masses for real sources and stacked data (see Sect. 4) and to study their actual uncertainties and the properties of Herschel galaxies in this context.

7.1. The real world and its limitations

We apply the methods described in Sect. 3 to the GOODS-N and GOODS-S data, taking care to limit the MBB fit to λrest ≥ 50μm to avoid contamination from warm dust. We adopted a value of β = 2.08 to simplify the comparison of MBB results to DL07 dust masses, and we used the Draine (2003) revisitation of κν (see also Li & Draine 2001).

The errors analysis is only focused on Mdust. In fact, simulations have already shown that via MBB fitting, Tdust can be determined within a 30% uncertainty as long as the short-wavelength side of the SED is constrained; the values of the γ, qPAH, Umin, and DL07 parameters can hardly be constrained by the available photometry. These results are also shared in the SED fitting of real SEDs discussed here and therefore are not covered further in this section.

The left panel of Fig. 11 shows the relative uncertainty on Mdust for a MBB fit to all 160 μm-selected GOODS-N/S sources as a function of the maximum available rest-frame wavelength, λmax, of each object. The right panel of Fig. 11 deals with DL07 models. In both cases, color coding is based on the maximum S/N of the available photometry for each object, with darker symbols indicating a higher S/N.

At larger rest-frame λmax, the number of available bands is also larger, mostly because of redshift effects. In the case of the MBB fit, we emphasize that sources with smaller tend to benefit from a higher quality photometry, i.e., a higher (maximum) photometric S/N at λobs ≥ 100μm. On the other hand, no such trend is seen in the DL07 case.

The relative uncertainty on Mdust increases as the maximum covered rest-frame wavelength decreases, thus confirming the findings of Monte Carlo simulations (see Sect. 6). In the case of GOODS fields, the use of photometric upper limits on the long-wavelength side does not provide a significant advantage in SED fitting because SPIRE data are much shallower than PACS data (up to a factor of ~10; see Table 2).

The performance turns out to be better in the DL07 case: The trend of vs. λmax is flatter, and the distribution is characterized by a smaller scatter at a given value of λmax,rest. For example, at a rest-frame λmax = 200 μm with DL07 modeling, it is always possible to constrain Mdust to a ~40% relative uncertainty, while the relative error on the MBB-based Mdust can reach values as large as 70%.

In the best case scenario of five to six bands available (covering up to 500 μm in the observed frame), the uncertainty on dust mass is typically on the order of 20% but can reach up to 3040%, depending on the S/N of the available photometry and on the adopted modeling.

7.2. Comparison of DL07 and MBB Mdust estimates

Figure 12 reports on the comparison of Mdust estimated with MBB and DL07 SED fitting. A median systematic offset of 50% between the two estimates is found. Studying objects with fully sampled SEDs (up to sub-mm wavelengths), Magdis et al. (2012, 2013) report a systematic offset of a factor ~2. Magnelli et al. (2012a) report a factor ~3 discrepancy. Comparing these different results is not as straightforward as it might seem. In fact, the underlying assumptions for the MBB modeling differ significantly. We have already mentioned in Sect. 5.2.1 there is a discrepancy between the β = 1.5 emissivity adopted by Magnelli et al. (2012b) and the β used by the DH02 models, which depends on the intensity of the radiation field. Finally, the κν (Li & Draine 2001; Draine 2003) adopted by most authors has a frequency dependency on the power of β = 2.08. Bianchi (2013) shows that when fitting local NGC galaxies and using β = 2.08, any offset should be washed away (see Sect. 3.2). Although we are using a consistent (β,κν) set for the two approaches, this incongruence with respect to Bianchi (2013) results might come from possible mismatched temperatures due to less well-sampled SEDs of high-z objects. It is therefore always important to report on the adopted setup when referring to this delicate comparison.

thumbnail Fig. 13

Comparison of three different Mgas estimates: the first is based on CO (rescaled to a common αCO,MW = 4.36, including a metallicity correction; Genzel et al. 2012, 2013, 2015); the second based on the scaling of depletion times, τdep; and the third derived from Mdust, using the δGDR-Z and M-Z relations in the PP04 metallicity scale. Left: datasets include GOODS-N (blue filled circles) and GOODS-S (red filled squares) 160 μm-detected sources; Herschel-detected SMGs (green filled pentagons; Magnelli et al. 2012a); CO-detected galaxies (different blue symbols, see below); and the Magnelli et al. (2014) stacked points (gray crosses). The right-hand panels only include CO-detected sources, namely: PHIBSS galaxies (blue triangles, Tacconi et al. 2013); BzK galaxies (blue squares, Daddi et al. 2010); other star-forming galaxies (blue filled circles, Magnelli et al. 2012b); lensed galaxies (empty blue circles, Saintonge et al. 2013); and SMGs (green upside down triangles, Bothwell et al. 2013). The dashed and dotted lines in the lower-right panel mark the ± 0.3 and ± 0.6 dex deviations from the 1:1 locus (solid line).

Open with DEXTER

Color coding the ratio of the two estimates by redshift, we see that the scatter in the distribution of points becomes very large above z ~ 1. At this redshift, one begins to lose PACS bands, on the blue side of the SED, because of k correction and sensitivity issues; at the same time, the SPIRE photometry also becomes poorer and more affected by confusion noise.

Similar results are obtained fitting the stacked photometry by Magnelli et al. (2014), which has the advantage of a more extensive SED coverage on the long-wavelength side (Fig. 12, right). In this case, it is seen that objects with longer λrest,max lie closer to the 1:1 locus. On the contrary, for poorly sampled SEDs, the discrepancy becomes larger (see also Sect. 6). We conclude that two concomitant effects contribute to the difference in MBB/DL07 mass ratios obtained by different authors: the underlying MBB model assumptions and the available spectral coverage.

7.3. Dust-based Mgas

Following Sect. 2.2, the DL07-based dust masses are converted into the molecular gas content of galaxies following δGDR-Z scaling of Magdis et al. (2012). Metallicities are computed from stellar masses with the parameterization of the M-Z relation by Genzel et al. (2015). Both of these relations are calibrated to the PP04 metallicity scale.

Figure 13 compares different Mgas estimates for the datasets in hand, all inclusive of the helium contribution. CO-based masses have been rescaled to a common Milky Way conversion factor αCO,MW = 4.36, and we have taken a metallicity correction into account as well. The adopted correction is the geometric mean of the Bolatto et al. (2013) and Genzel et al. (2012) dependencies of αCO on metallicity (see Genzel et al. 2015).

Gas masses based on the scaling of depletion times have been obtained adopting Eq. (3) in Sect. 2.1 (Genzel et al. 2015, and priv. comm.) and the Whitaker et al. (2014) definition of star formation main sequence.

In the case of CO-detected galaxies, the CO-based and τdep-based Mgas estimates are consistent overall within a factor of ~2 with only few exceptions. For the low-mass lensed galaxies, Saintonge et al. (2013) adopted the metallicity-dependent values of αCO derived using the Genzel et al. (2012) relation, nevertheless, the metallicities of most of their sources are out of the range where the Genzel et al. (2015)τdep scaling was calibrated (and holds). It is thus no surprise that these sources show an offset between CO- and τdep-based determinations of Mgas (bottom right panel of Fig. 13).

As far as CO- and dust-based estimates are compared (top-right panel of Fig. 13), for roughly one third of the sources the two Mgas estimates differ by more than 0.3 dex with dust systematically providing a higher Mgas than CO. The SEDs of these sources are missing some long-wavelength bands and, therefore, the Mdust estimate is affected by systematics, as seen in previous sections.

Dust-based Mgas of the stacked photometry by Magnelli et al. (2014) are in good agreement with τdep-based estimates with a relatively small scatter of ~0.23 dex, which is on the same order of the scatter in the adopted τdep(z,M,sSFR) scaling relation (Genzel et al. 2015). The situation for individual objects is more complex: at the low-mass end, the two estimates of Mgas are in fair agreement, but, at high masses, the dust-based estimate is significantly larger than τdep results. The mismatch between the two Mgas determinations becomes larger as the number of available bands decreases, mainly driven by limited wavelength coverage. This holds both for PACS-selected sources and for CO-detected objects.

For the general GOODS-N/S FIR population, λmax(rest) ≤ 200μm and can be as low as 100 μm even at intermediate-low redshift because of the large noise in SPIRE bands (see Table 2). According to our simulations, in this case dust masses can be overestimated (see Sect. 6.3), thus explaining the difference in Mgas in these cases. At higher redshift, λmax(rest) becomes even shorter if no sub-mm detections are available, and the overestimate of masses becomes more critical (see Sect. 5.2). Figure 15 limits the results on GOODS-N/S sources to λmax(rest) ≥ 160μm.

When an object benefits from sub-mm observations (e.g., as is the case for most CO-detected galaxies and SMGs), then the maximum rest-frame wavelength available is on average longer than for other PACS-selected galaxies. If a sufficient number of bands is available, and the SED coverage is fine enough, then τdep- and DL07-based masses turn out to be in good agreement. Outliers with Mgas(DL07) ≫ Mgas(τdep) suffer from as poor SED coverage as analogous GOODS-N/S cases.

It is worthwhile recalling that the adopted τdep scaling was calibrated for MS galaxies of nearly solar metallicity (Genzel et al. 2015). Therefore, we also expect a contribution to the Mgas(DL07) /Mgas(τdep) mismatch from the adopted τdep scaling. Figure 14 exemplifies the possible trends of this ratio as a function of distance from the MS of star formation in terms of Δ(log (sSFR))MS. The Whitaker et al. (2014) definition of MS has been adopted. The dependence of Mgas on the distance from the star-forming MS found by Genzel et al. (2015) does not play a role in this case because it is factorized out by taking the Mgas,DL07/Mgas,τdep ratio.

thumbnail Fig. 14

Comparison of Mgas as based on dust masses and on the scaling of τdep, as a function of distance from the main sequence of star formation. The Whitaker et al. (2014) definition of the MS has been adopted. Symbols are as in Fig. 13.

Open with DEXTER

thumbnail Fig. 15

Same as Figs. 13 and 14 for GOODS-S/N sources only, limiting the maximum available rest-frame wavelength to μm.

Open with DEXTER

In the case of the Magnelli et al. (2014) stacked photometry, the data show a trend of the τdep-based Mgas estimates of bins above the MS to be larger than those based on dust. Limiting to Δ(log (sSFR))MS ≤ ± 0.5 dex (vertical dashed lines), we restrict the analysis to the MS proper locus. The scatter of the Mgas(DL07) /Mgas(τdep) decreases to 0.19 dex for the stacked photometry. As a result of stacking, the SED is covered up to λmax(rest) ~ 160μm and up to z = 2. Our simulations (Sect. 6.3) prove that this is enough to avoid underestimates of Mdust and keep overestimates to less than a factor of 2 (details depending on the actual values of qpah, Umin, and γ).

As for individual sources, the majority of catastrophic failures lying at > 3σ (measured above using the result for stacking) above Mgas(DL07) /Mgas(τdep) = 1 turn out to be on the high-z tail of the Herschel sample. The SEDs of these galaxies are poorly sampled and are only covered up to λmax(rest) ≤ 100μm.

In summary, the analysis of the SEDs of real sources and stacked photometry by means of DL07 fitting confirms the expectations from the MC exercise discussed before. When the mid- and FIR rest-frame SEDs are covered up to at least 160 μm, the DL07 dust-based estimate of Mgas is reliable, not affected by systematics, and consistent with independent estimates based on τdep scalings or CO observations. If the SED wavelength coverage is poorer, i.e., limited to shorter wavelengths, dust masses can be overestimated. This might happen because the SED is redshifted and therefore the available bands sample shorter wavelengths and also because of a decreased S/N and increased confusion in the case of faint, high-z galaxies. If the maximum rest-frame wavelength available falls shorter than 100 μm, DL07-based dust masses are severely overestimated and should be not recommended anymore.

7.4. The dust-to-gas ratio

thumbnail Fig. 16

Gas-to-dust mass ratio, δGDR, of CO-detected galaxies with rest-frame λmax ≥ 200μm, computed by combining CO-based gas masses and DL07-based dust masses. Symbols are as in Fig. 13. Left: δGDR as a function of metallicity. The red lines represent the local relation by Magdis et al. (2012) and its 0.15 dex scatter. Right: δGDR as a function of other relevant parameters: redshift z; depletion times based on CO; M (Wuyts et al. 2011a); and .

Open with DEXTER

One assumes a value of the gas-to-dust mass ratio, δGDR to compute the molecular gas content of galaxies from their dust masses. As seen in Sect. 2.2, it is common to adopt the local scaling of δGDR with the metallicity derived by Leroy et al. (2011) or some a variant of it (e.g., Magdis et al. 2012). This procedure implicitly assumes that the relation holds regardless of redshift.

The sample of CO-detected sources in hand benefits from fully independent measurements of gas and dust masses. We combine CO-based Mgas with the DL07-based determinations of Mdust to derive δGDR. All gas masses have been renormalized to the Galactic αCO,MW = 4.36 and a metallicity dependence of αCO has been included as well (Genzel et al. 2015). The M-Z relation described in Eq. (6) is used with stellar masses from BC03 SED fitting (Wuyts et al. 2011a) to produce metallicities on the PP04 scale.

We only consider sources with μm hereto minimize the uncertainties on Mdust and avoid systematics (see Sects. 5 and 6). The left panel in Fig. 16 shows the resulting trend of δGDR vs. compared to the Magdis et al. (2012) local relation and its 0.15 dex scatter. A ±1.0 dex systematic uncertainty on its zero point is also reported by Magdis et al. (2012).

The data lie close to the local relation and are consistent with it within the uncertainties and possible systematics. At the high metallicity end, the scatter is very large and there exist cases with very small δGDR, which are significantly below the locus occupied by local galaxies. The right-hand panel of Fig. 16 seeks possible dependencies of δGDR on other derived quantities. No significant trend is found as a function of M, SFR, τdep, or redshift. The most critical outlier turns out to be a z ≃ 1.4 galaxy with a very high Mdust ~ 8 × 109 [M] and an SED with 10205 μm rest-frame coverage.

7.5. Synergies between Herschel and ALMA

The photometric determination of dust masses through SED fitting allows us, in principle, to obtain an estimate of gas masses for a very large number of objects via the δGDR scaling. This approach offers the advantage that photometric observations are still significantly faster than sub-mm spectroscopy, even with the last generation of FIR or sub-mm facilities, and especially for galaxies at z> 1.

Far-infrared photometry also has the advantage of sampling the dust emission of galaxies near the SED peak, providing a calorimetric measurement of their SFR (Elbaz et al. 2011; Nordon et al. 2010) for a relatively cheap time investment. Nevertheless, the analysis has shown the effects of the limited wavelength coverage: when the SED extends only up to λrest,max ≤ 200μm, the uncertainty on the derived dust mass can be very large and systematics might affect this measurement in a fraction of cases (see Sects. 7.1, 6).

On the sub-mm side, assuming the Rayleigh-Jeans (RJ) tail of SEDs is observed, Scoville et al. (2014) developed a strategy aimed at deriving gas (dust) masses on the basis of a one-band sub-mm continuum measurement. If confirmed reliable, this approach would be extremely competitive in terms of exposure time with respect to CO spectroscopy. Genzel et al. (2015) showed that this method can lead to incorrect (up to a factor >3) gas masses, even when applying corrections to take the RJ approximation into account. In fact, the proposed scaling of sub-mm luminosity holds for a specific dust temperature Tdust = 25 K, which is a condition that is only met in a few cases. Magnelli et al. (2014) have shown that Tdust varies across the MS of star formation, and increases as a function of Δ(log (sSFR))MS. A dust temperature measurement is hence necessary to reach a correct estimate of Mdust from the sub-mm emission, even in the simple MBB, RJ approximation.

If not using the available Tdust(z,M,sSFR) scaling (e.g., Magnelli et al. 2014), Genzel et al. (2015) suggest a strategy based on continuum observations in two distinct sub-mm bands. In this context, using ALMA bands 6 and 7 (centered at 1100 and 850 μm), the best relative uncertainty of Mdust that can be reached is on the order of 60%, for galaxies at z ~ 2 with 10σ detections in the two sub-mm bands. The uncertainty on dust mass becomes worse at lower redshift. In practice, the SED coverage is too small and too far off the SED peak to allow for a robust determination of Tdust and Mdust. Adopting bands 7 and 9 (at 850 and 450 μm) lowers to 3040% at best and using bands 6 and 9 finally brings it to 2030%, i.e., a better than 3σ estimate of Mdust. However, band 9 is particularly demanding in terms of atmospheric conditions. At z< 1, this approach is not necessarily faster than ALMA CO spectroscopy. We defer to Genzel et al. (2015) for further details.

It is interesting to study how FIR and sub-mm observations complement each other and how combining them improves Mdust measurements. We simulate the ALMA performance in conjunction with Herschel data in the following way. For each galaxy detected by PACS in the GOODS-S field, we assume that the best-fit DL07 model obtained in Sect. 7.1 represents the actual emission of the galaxy. This best-fit model is then convolved with a box ALMA passband, centered at ~1100 μm (i.e., in ALMA’s band 6). The artificial photometry thus obtained is added to the real Spitzer and Herschel data and the extended SED is fitted again following the usual procedure. We assume that a S/N of 5 holds in ALMA band 6.

Figure 17 compares the uncertainties on Mdust obtained with and without the artificial band 6 photometric point as a function of the maximum rest-frame wavelength covered by real data, λrest,max, and dust mass (see also Fig. 11). The availability of good quality (S/N = 5 in this case) photometry at λobs = 1100 μm reduces the uncertainty on Mdust to 33% for >85% of sources. For comparison, only ~20% of the sample have similarly good quality Mdust estimates when no ALMA data are available. Similar results are obtained for the Magnelli et al. (2014) stacked data, adding the ALMA band 6 filter to the simulation described in Sect. 5.

We then use the artificial ALMA band 6 flux density obtained above to derive a gas mass expectation applying the recipe by Scoville et al. (2015; see their Appendix A and also Scoville et al. 2014). As usual, DL07-based dust masses are transformed into Mgas adopting the Genzel et al. (2015)M-Z-z and the Magdis et al. (2012)δGDR-Z relations (see Sect. 2.2). The Mgas expectation based on band 6 flux only (Scoville et al. 2015) includes the RJ correction prescribed by the authors. Figure 18 compares the two.

thumbnail Fig. 17

Results of DL07 SED fitting to galaxies detected by Herschel in the GOODS fields, including artificial ALMA band 6 photometry. Color coding is based on the number of available photometric bands, ranging from only four bands (lighter colors) to eight bands (darker colors). This number reflects the maximum rest-frame wavelength available as well as redshift dependencies. Top: relative uncertainties on dust mass, obtained only with the real photometry (see also Fig. 11). Bottom: the same quantity obtained when adding artificial photometry in the ALMA band 6 (~1100 μm).

Open with DEXTER

thumbnail Fig. 18

Comparison of the DL07-based Mgas estimate obtained, including artificial ALMA band 6 photometry to expectations, based on the Scoville et al. (2015) recipe, and applied to the same artificial ALMA band 6 data. Color coding is based on the number of available photometric bands (see Fig. 17). Left: direct comparison: error bars only include statistic noise. Right: ratio of the two estimates as a function of redshift and DL07-based gas mass.

Open with DEXTER

The DL07-based Mgas estimate tends to be systematically larger than the expectation based on the recipe by Scoville et al. (2015). For a small fraction of the sample the opposite trend holds. The overestimate becomes larger at the high-mass end. The trend seen here is the combination of the systematic higher Mdust obtained with DL07 modeling with respect to MBB (see Sect. 7.2), the fact that dust temperature is varying as a function of the position in the M-SFR-z space (see Genzel et al. 2015; Magnelli et al. 2014), a slight difference in the β values underlying the two methods (β = 2.08 vs. 1.8), and possibly other systematics (e.g., related to metallicity dependencies). This is different from the problem encountered in Fig. 13 at large Mgas because now we are reasoning in relative terms (DL07 vs. Scoville methods), and the long wavelength side of the SED is constrained by the artificial ALMA photometry.

Following a similar path to Genzel et al. (2015), we now estimate the on source exposure time, texp, needed to reach S/N = 5 in band 6 for our Herschel-detected star-forming galaxies based on DL07 models.

To this aim, the ALMA sensitivity calculator7 is used with the standard Cycle 3 configuration via an array of 36 12 m antennas, requesting an angular resolution of 1 arcsec. The left-hand panels of Fig. 19 show texp for individual sources and the cumulative texp for sources below a given dust mass (solid histograms) along with the number distribution of objects (dotted histograms). The dispersion of the points in the upper panel reflects the spread in redshift of the sample. However, the sub-mm negative k correction implies that exposure times for a given mass are of the same order of magnitude at all redshifts.

Targeting all Herschel galaxies detected by the PEP/GOODS-H survey in the GOODS-S field down to Mdust = 109 [M] would require a few minutes of on-source exposure, i.e., without accounting for overheads. Reaching down to Mdust = 108.5 [M] requires roughly 2025 h on source with ALMA band 6 in the above mentioned configuration. Using ALMA band 7 would need a similar amount of time.

The right-hand panel of Fig. 19 includes the distribution of the expected ALMA band 6 fluxes as a function of dust mass. Instead of requesting a fixed S/N, a 0.38 mJy and 0.18 mJy depths are shown, similar to recent ALMA Cycle-3 approved surveys8. At these depths, ALMA detects the majority of Herschel sources above Mdust = 108.5 [M] with a fainter tail extending down to ~108 [M], and with the added value of possible undetected PACS objects, which were not included in the current analysis.

Assuming an average dust/gas ratio of 0.01, the limits mentioned above correspond to Mgas = 1010.5−11 [M]. Berta et al. (2013a) showed that the Schechter characteristic mass, of the molecular gas mass function lies between Mgas = 1010.3−10.9 [M] at z = 0.2−2.0. A measurement of the molecular gas mass function based on dust observations, modulo the δGDR scaling, and on the dust mass function itself down to and up to z ~ 2, is within reach of ALMA in 2025 h spent on source.

thumbnail Fig. 19

Further results of DL07 SED fitting including artificial ALMA band 6 photometry. Color coding is based on the number of available photometric bands (see Fig. 17). Left: exposure time estimate (without overheads) for ALMA band 6 observations of Herschel-detected galaxies. The texp estimate for individual sources (upper bigger panel, dots) are shown and summed up to compute the cumulative texp needed to observe all sources above a given value of Mdust (solid histograms). Dotted histograms represent the number distribution of sources. Right: distribution of ALMA band 6 expected fluxes, as a function of dust mass, for galaxies detected by Herschel in the GOODS-S field. The horizontal dashed lines correspond to 0.38 and 0.18 [mJy] flux limits (dark and light gray, respectively). The right side and bottom panel show the projected 1D distributions. The hatched dark and light gray histograms only include sources above the flux limits indicated.

Open with DEXTER

8. Conclusions

We have exploited the deepest FIR blank field maps available to date from the PEP, GOODS-H, and HerMES surveys to study the feasibility of deriving dust and gas masses via SED fitting for individually detected sources. In parallel, we built extensive Monte Carlo simulations to study the limitations of real data, and to understand how they influence the uncertainties and systematics on such dust mass determinations. We focused the analysis on two popular modeling approaches: SED fitting by means of a single-temperature modified blackbody and by means of the Draine & Li (2007) model. The main results of the analysis based on artificial sources and MC sampling are:

  • FIR SED fitting recovers dust mass consistently as long as the wavelength coverage offered by the data extends at least up to 160200 μm (rest frame). In this case, no systematics are expected on Mdust. In contrast, if the >3σ detections fall shorter than this wavelength limit, Mdust can be severely overestimated with the amplitude of the systematic effect depending on the details of the model.

  • The uncertainty of Mdust also strongly depends on wavelength coverage. As a rule of thumb, it is not possible to reach a ≥ 3σ determination of Mdust if the available photometry does not extend at least up to ~200 μm in the rest frame.

  • The determination of dust temperatures based on MBB fitting is rather stable, even if long-wavelength bands are missing. In addition, Tdust is always constrained to better than a 1020% uncertainty as long as the blue side of the FIR SED is constrained. Nevertheless, a small offset on Tdust can induce a large systematic error on Mdust.

  • The discreteness of the sampled wavelengths (photometric bands), combined with redshift, can cause significant residual systematics in the determination of Tdust, up to few K.

  • Reliability tests based on artificial catalogs tend to minimize systematics on Mdust if the SED shape adopted for SED fitting is consistent with the one used to produce the mock catalog. In other words, these kinds of simulations assume that the adopted SED library is consistent with real-world’ SEDs. On the other hand, small differences in SED shape (e.g., DH02 vs. DL07 models) can induce strong and not easily predictable systematic effects.

The GOODS fields benefit from an extensive multiwavelength coverage from the UV to FIR and sub-mm wavelengths. The SEDs of galaxies detected by Herschel can reach up to 200 μm rest frame and beyond. The SEDs of individual sources have been fitted with DL07 and MBB models. At the same time, stacked photometry of NIR selected galaxies has also undergone the same analysis, binned in z-M-SFR space. The main results of these pieces of analysis are:

  • At the depth of the deepest Herschel extragalactic surveys (GOODS-S as observed by PEP, GOODS-H, and HerMES), it is possible to retrieve dust masses with a S/N ≥ 3 for galaxies on the main sequence of star formation down to M ~ 1010 [M] up to z ~ 1. At higher redshift (z ≤ 2), the same goal is achieved for objects only lying increasingly above the MS at similar stellar masses or for galaxies at the tip of the MS (i.e., with higher M). At shallower depths (e.g., in the case of the COSMOS field), this reasoning shifts to even higher values of SFR. Dust temperatures (based on MBB fit) can be constrained within a 10% accuracy in most of the cases across the z-M-SFR space, modulo residual systematics due to the discreteness of SED sampling (see above).

  • As in the case of simulated data, Spitzer and Herschel data alone are not sufficient to produce an estimate of Mdust to better than 30% uncertainty if the maximum rest-frame wavelength covered by the data is shorter than ~160200 μm.

  • Comparing MBB- and DL07-based masses, the average offset between the two, regardless of redshift, is a factor ~1.5. At z> 1, the scatter of the MBB/DL07 mass ratio becomes very large, mainly because photometric points are progressively missed. We stress that to allow a direct and meaningful comparison, it is paramount to adopt a consistent set of parameters, taking special care in the value of β, i.e., the power of the dependence of dust emissivity on frequency.

  • Dust masses estimated with DL07 modeling are more robust than those based on MBB: relative errors are more mildly dependent on the maximum covered rest-frame wavelength and less scattered.

  • Dust mass estimates, based on DL07 SED modeling and on scaling of depletion times, τdep, are consistent with each other as long as the data guarantee a sufficient wavelength coverage. Applying the local dependence of δGDR on metallicity to transform Mdust into Mgas, these estimates are overall consistent with CO-based estimates for a small sample of star-forming galaxies. While comparing Mgas estimates based on different methods, it is important to adopt a consistent set of relations, calibrated to a common metallicity scale.

  • Using CO-based Mgas renormalized to αCO,MW = 4.36, Mdust obtained through DL07 modeling, and metallicities computed with the M-Z-z relation (Genzel et al. 2015), the δGDR of z> 1 galaxies depends on metallicity in a similar manner as for local galaxies within uncertainties and systematics.

CO-based Mgas estimates, which also represent an important anchor for validation of indirect methods despite αCO uncertainties, are still limited at z> 1 (see also Genzel et al. 2015). More CO observations of individual galaxies covering a wide range of parameters are highly desirable and are becoming more easily accessible with ALMA and NOEMA. In parallel, also the determination of gas-phase abundances and, therefore, of the redshift evolution of ISM physical properties (e.g., Kewley et al. 2013; Shapley et al. 2015), will undergo significant progress over the next few years thanks to NIR multiobject spectroscopy (e.g., with KMOS, MOSFIRE).

Finally, as Mdust estimates based on Spitzer and Herschel photometry are limited to cases with high quality SEDs available, we recalled the advantages and limitations of estimates including only sub-mm data (e.g., Genzel et al. 2015). For example, a scaling of sub-mm fluxes into Mdust can be affected by strong systematics if the characteristic dust temperature of the SED is not known. Continuum observations in two sub-mm continuum bands might eventually help, but they are time consuming and the results are still limited by large uncertainties because the wavelength range covered is relatively small.

We therefore explored a combined IR plus sub-mm approach, combining existing Herschel data to expected ALMA 850 or 1100 μm continuum fluxes. These single band observations allow one to reduce the uncertainties on Mdust down to <30% for virtually all Herschel-detected galaxies in the GOODS-S field. A direct measurement of the molecular gas mass function based on dust observations up to z ~ 2 will be soon within reach.


1

Using ex = 1 + x + x2/ 2 ! + x3/ 3 ! + ...

5

So that the first realization has all eight bands available (40850 μm); the second has seven (40500 μm); the third has six (40350 μm); and so on until the last with only two bands (4070 μm).

6

So that the first realization has all eight bands available (40850 μm); the second has seven (70850 μm); the third has six (100850 μm); and so on until the last with only two bands (500850 μm).

8

See list of ALMA Cycle-3 high priority projects: http://almascience.eso.org/observing/highest-priority-projects

Acknowledgments

The authors would like to thank the referee, Dr. Georgios Magdis, for his useful and effective comments that improved the quality of this manuscript, as well as for his positive and fast response. PACS has been developed by a consortium of institutes led by MPE (Germany) and including UVIE (Austria); KU Leuven, CSL, IMEC (Belgium); CEA, LAM (France); MPIA (Germany); INAF-IFSI/OAA/OAP/OAT, LENS, SISSA (Italy); IAC (Spain). This development has been supported by the funding agencies BMVIT (Austria), ESA-PRODEX (Belgium), CEA/CNES (France), DLR (Germany), ASI/INAF (Italy), and CICYT/MCYT (Spain). SPIRE has been developed by a consortium of institutes led by Cardiff University (UK) and including University of Lethbridge (Canada), NAOC (China), CEA, LAM (France), IFSI, University of Padua (Italy), IAC (Spain), Stockholm Observatory (Sweden), Imperial College London, RAL, UCL-MSSL, UKATC, University of Sussex (UK), Caltech, JPL, NHSC, University of Colorado (USA). This development has been supported by national funding agencies: CSA (Canada); NAOC (China); CEA, CNES, CNRS (France); ASI (Italy); MCINN (Spain); SNSB (Sweden); STFC, UKSA (UK) and NASA (USA).

References

Appendix A: Ancillary data for CO samples

thumbnail Fig. A.1

Comparison of input and output Mdust in DL07 simulations obtained with ten photometric bands. Left/right columns present catalog entries obtained removing bands from the long-wavelength and short-wavelength side of SEDs, respectively. Top/middel/bottom diagrams show results obtained with specific values of two parameters and sampling the 3rd vary within the ranges recommended by Draine et al. (2007, see also Sect. 3.3. When fixed, the adopted values are: qPAH = 0.47, Umin = 0.7, and γ = 0.11.

Open with DEXTER

The PHIBSS survey (Tacconi et al. 2013) performed CO spectroscopy of galaxies at z = 1−2, and included additional data from past work, for a total of 73 objects. The survey observed 38 sources in the Extended Groth Strip (EGS) field, also targeted by PEP and HerMES. Another 18 objects belong to the “Q-fields” (Steidel et al. 2004), Tacconi et al. (2013) also included six BzK sources by Daddi et al. (2010) and six sources by Magnelli et al. (2012b), plus a few other additional objects.

FIR photometry is sought in the EGS from the PEP and HerMES data (Lutz et al. 2011; Oliver et al. 2012), using a closest-neighbor algorithm and visual inspection of multiwavelength maps (IRAC, MIPS 24 μm, PACS, and SPIRE). The sources by Magnelli and Daddi, and few other isolated objects studied by Magnelli et al. (2012b) and Saintonge et al. (2013), and their photometry can be found in their works.

Out of the 38 PHIBSS sources in EGS, we have the following detection statistics: 37, 18, 15, 14, 15, and 13 are detected at 24, 100, 160, 250, 350, and 500 μm, respectively. Five out of the six BzK sources by Daddi et al. (2010) included in PHIBSS have 24 μm to 500 μm photometry by PEP and HerMES (Lutz et al. 2011; Oliver et al. 2012). The object BzK21000 benefits from millimeter photometry (1.3, 2.2, 3.3 mm) by Dannerbauer et al. (quoted as in prep. by Magdis et al. 2011), Carilli et al. (2010), and Daddi et al. (2009) as collected by Magdis et al. (2011). CB58 and the Cosmic Eye have photometry collected by Saintonge et al. (2013). For the other sources in the PHIBSS sample (mainly the Q-fields), no FIR photometry has been retrieved from the literature yet.

Magnelli et al. (2012b) included in their analysis six PEP sources, two HDF sources, and the well-studied object GN20. They actually included nine sources in their observations, but only six are used by Tacconi et al. (2013). For the six PEP objects, 24-to-500 μm photometry has been retrieved from PEP and HerMES catalogs (Lutz et al. 2011; Oliver et al. 2012). The two HDF169 and HDF242 objects also have 24-to-500 μm photometry Magnelli et al. (2012b). In addition to the 24-to-500 μm imaging, GN20 also has 850, 1100, 2200, 3300, 6600 μm photometry collected by Magdis et al. (2011) and obtained by Pope et al. (2006), Perera et al. (2008), Dannerbauer et al. (2009), Daddi et al. (2009), Carilli et al. (2011). For this source, a 1160 μm photometry is also available, obtained by Penner et al. (2011) from the combination of 1100 μm Atzec and 1200 μm MAMBO maps.

Saintonge et al. (2013) have observed 17 lensed galaxies, out of which ten have CO detections (and derived physical quantities). The available photometry includes 100-to-500 μm PACS and SPIRE data, and 1200 μm IRAM/MAMBO measurements for five of them.

Bothwell et al. (2013) observed 40 SMGs in CO and detected 32 of them.

They also included some additional data taken from the literature. Their sample is spread over several different sky areas: three sources in the Subaru deep field, or UDS; three sources in smaller fields; five sources in the Lockman Hole East (LH-East); 16 sources in the HDF; three sources in SSA-13; seven sources in ELAIS-N2 (seven sources); and three sources in SSA-22. The Subaru Deep field (SXDF or UDS) is part of HerMES; the LH-East is part of PEP and HerMES; the HDF is part of PEP and HerMES (GOODS-N); ELAIS-N2 is part of HerMES, but has not been released yet (at the time of writing, in DR2 and DR3). In the UDS field, only 250, 350, 500 μm bands from HerMES DR2, plus the 850 μm by Bothwell et al. (2013) are available.

In the HDF (GOODS-N) there are the 24, 100, 160, 250, 350, 500 μm data by PEP, GOODS-Herschel and HerMES, plus the 850 μm fluxes by Bothwell. In the LH-East there are the 24, 100, 160, 250, 350, 500 μm data by PEP + HerMES, plus the 850 μm fluxes by Bothwell. None of the Bothwell et al. (2013) sources is detected in the AzTEC 1.1 mm maps by Michałowski et al. (2012).

In synthesis, out of the 40 SMGs in the sample by Bothwell et al. (2013), 16 have a secure CO detection and six have a candidate detection. Out of these, 16 detected objects lie in UDF, LH, or HDF, and two candidates lie in HDF. Out of these 16+2, 13 have enough mid- and FIR photometric detections to ensure an estimate of dust mass based on DL07 SED fitting. Out of these 13, seven sources have enough data to be included in our analysis (i.e., Mgas,CO, Mdust, Mstars).

Appendix B: Detailed analysis of systematics

We dissect the fine details of possible systematic effects in the derivation of Mdust, as emerging from the MC runs with rest-frame photometry (see Sects. 6, 6.3, and 6.5).

thumbnail Fig. B.1

Comparison of input and output Tdust (top) and Mdust (bottom) in MBB simulations, split in bins of input dust temperature. Left/right panels again refer to cases with long- or short-wavelength bands removed, respectively. Color coding is the same as in Fig. 10

Open with DEXTER

Appendix B.1: DL07 systematics

The DL07 case is analyzed first. We only focus on dust mass estimates because the other parameters are much more poorly constrained (see Sect. 6.2).

Figure A.1 presents the comparison of Mdust,out and Mdust,in as a function of λmax,min, i.e., for cases obtained removing long-wavelength bands (left hand diagrams) and short-wavelength bands (right). Simulations including ten photometric bands are shown. Possible dependencies on qPAH,in, Umin,in, and γin are studied by fixing two parameters and splitting the analysis in bins of the third.

First of all, we focus on the top row of diagrams, showing cases at specific values of qPAH,in = 0.47 and Umin,in = 0.70 and sampling γin in the allowed range. When removing long-wavelength bands, there is a tendency to underestimate Mdust, if the bands coverage and the S/N are not adequate. This tendency becomes smaller for larger values of γin.

This is still a special case, and it is necessary to disentangle the effect of the other two parameters qPAH and Umin (see middle and bottom panels in Fig A.1). For example, moving to Umin,in = 2.0, the tendency is to overestimate Mdust, rather than to underestimate it. This difference is indeed mainly driven by Umin,in (see middle panels), but the overestimate becomes larger for larger γin. A similar effect is also shared by qPAH (bottom panels). A larger value of γin produces an increased chance to overestimate Mdust. Nevertheless, if there is a general tendency to underestimate Mdust (e.g., because Umin,in is small), then this underestimate becomes milder for the larger values of γin. If, instead, there is the tendency to overestimate Mdust (e.g., because of the large value of Umin,in), then this tendency becomes even worse with larger γin.

However, here we are speaking of large values of γin (>0.4), which are rarely seen in real galaxies (see, e.g., Tab. 4 in Draine et al. 2007). Keeping γin< 0.2, differences in trends are more difficult to appreciate.

When removing short-wavelength bands, Mdust is usually easily retrieved and there are not very significant trends related to γin in over- and underestimating Mdust. This holds also varying qPAH,in and Umin,in.

thumbnail Fig. B.2

Direct comparison of input/output Tdust vs. Mdust in MBB simulations, split in bins of input dust temperature. Left/right panels again refer to cases with long- or short-wavelength bands removed, respectively. Color coding is based on the number of available bands; from two (light pink) to eight (red).

Open with DEXTER

We now examine the middle panels in Fig. A.1, highlighting what happens when sampling Umin,in variations at qPAH,in = 0.47 and γin = 0.11. When removing long-wavelength bands, there are now big differences in the in/out trends as a function of λmax, depending on the actual value of Umin,in. At low Umin,in there is a tendency to underestimate Mdust if the band coverage is poor, while at high Umin,in there is a tendency to overestimate it.

This recalls what is happening with the MBB (see Sect. 6.5). Changing the input Umin,in, in fact, produces a change of the position of the FIR peak (see Fig. 1 and Sect. 3.3). Consequently, there is a tendency to overestimate Mdust when the combination of bands coverage and “equivalent temperature” is produces the effects described for the MBB case study (see Sect. 6.5). At lower Umin,in, the SED peaks at longer wavelengths and at some point it exits the range where overestimate of fluxes is “preferred”, causing underestimates of Mdust to dominate.

When removing short-wavelength bands the same effects are seen, although with a much smaller amplitude.

Finally, in the bottom panels of Fig. A.1, we focus on Umin,in = 0.7 and γin = 0.11, and follow the variation of qPAH,in. In this case, as in the case of variable γin, the trend to over- or underestimate Mdust changes as a function of the value of qPAH,in.

When removing long-wavelength bands, for the quoted choice of Umin,in and γin, dust masses can be underestimated for small values of qPAH,in, but this trend becomes much milder for large values. On the other hand, remember that for Umin,in = 2.0 and same range of γin, the tendency is to overestimate Mdust (see middle panels of Fig. A.1), and the tendency is to have a larger overestimate when qPAH,in becomes larger.

In summary, the larger qPAH,in, the more Mdust can be overestimated. If there is a tendency of underestimating Mdust (e.g., driven by a small value of Umin,in), then a large qPAH,in makes this underestimate milder. On the other hand, if Umin,in is large and thus there is already the tendency to overestimate Mdust, then a large qPAH,in makes things worse.

Using Eq. (20) to derive the mean radiation field U, a similar behavior of the ratio Uout/ ⟨ Uin is observed. Depending on the actual values of Umin, γ and qPAH, the chance of overstimate or underestimate U can be amplified or suppressed due to summing or compensating effects of the three parameters. Since U and Mdust is linked by Eq. (19), this explains the relative stability of the Mdust estimate on a more fundamental level.

When removing short-wavelength bands, as usual there are very small tendencies for under- or overestimates of Mdust. In this case, by varying qPAH,in these trends do not seem to change significantly, or at least not significantly enough to be appreciated by this analysis.

Simulations limited to a maximum wavelength of 250 μm, obviously behave equivalently when removing long-wavelength bands. On the other hand, when removing short-wave bands, things are different because the three longest wavelength bands are missing. Trends as a function of γin are still very small and are only slightly enhanced with respect to the ten bands case, but still much smaller than when removing long-wavelength bands. Trends as a function of Umin,in also start having a larger amplitude when removing short-wavelength bands, while trends with ten bands continued to be marginal. Finally qPAH,in trends are marginal, as in the case of ten bands.

Appendix B.2: MBB systematics

It is worthwhile to verify how reliably the MBB fit to the artificial photometry can retrieve the input values of dust temperature and mass. Figure B.1 shows the comparison of input and output quantities, as a function of λmax,min, split in bins of input dust temperature, Tdust,in.

As the average photometric uncertainty grows (lighter blue dots), there is the chance to systematically overestimate Mdust. The overestimate worsens when removing datapoints at the long-wavelength side, while the problem is relatively milder (but still significant) when removing bands at the short-wavelength side.

When removing datapoints at the long-wavelength side, Tdust is always well recovered, but when removing datapoints at the short-wavelength side it can be systematically underestimated or overestimated. This effect reflects the fact that the peak of the MBB shifts to longer wavelengths as dust temperature decreases, combined with the poor sampling of the SED obtained by removing datapoints.

Nevertheless, the incidence of catastrophic Mdust overestimates is small: only ~10% of cases turn out to have Mdust,out> 2 × Mdust,in. Finally, there exists a 0.1% of cases with the tendency to underestimate Mdust for the lowest Tdust,in considered, when removing datapoints at the short-wavelength.

Although affecting a limited fraction of cases, it is interesting to study in detail the nature of these systematic trends to understand what causes them. Figure B.2 directly compares the ratio of input/output Tdust to the same quantity for Mdust. Color coding is now based on the number of available bands.

As long as short-λ bands are progressively removed, it is possible to overestimate Tdust at low input temperatures (Tdust,in ≤ 20 K, for example). Consequently, the shape of the MBB SED changes, but the long-λ bands still give a good constraint on the model normalization, hence, on Mdust. Something different happens at higher input dust temperatures. At Tdust,in = 25−35 K, the output Tdust is more stable but there are still variations, which start driving an overestimate of dust mass. The peak of the MBB starts moving out of the covered wavelength range and the shape of the model in the covered range of wavelengths is roughly constant (Rayleigh-Jeans regime; RJ). Nevertheless, an underestimate of Tdust implies a lower emissivity, and thus a larger normalization (i.e., a larger dust mass) is needed to reproduce the “observed” fluxes. This causes an overestimate of dust mass.

When removing long-λ bands, the trend is again to overestimate Mdust, but the cause is more subtle. In fact, in this case, Tdust,out is rather consistent with Tdust,in. This is because the short-wavelength side of the SED is always constrained by the available photometry. There is still a possibility of slightly underestimating dust temperatures, however, because by anchoring the model at the short-wavelength side, it is still possible to reproduce the short-λ colors with a lower temperature, which allows for a different normalization (i.e., a different Mdust). The smaller the number of long-wave bands, the larger is the freedom in renormalization. However, this time a small effect on Tdust translates into a big effect on Mdust. A small variation of Tdust implies a significant change in the shape of the MBB SED blueward of the peak. Now the short-wavelength side of the SED is constrained by the data because we are removing long-λ bands. Therefore a small underestimate of Tdust induces a significant overestimate of Mdust.

In the previous case (i.e., removing short-λ bands), the opposite was happening: a big variation on Tdust implied a relatively smaller variation of Mdust because we were sampling the RJ side of the SED, where shape variations due to Tdust changes are less prominent than at the short-wavelength side. For example, focusing on input temperature in the range Tdust,in = 45−50 K one can note that:

  • When removing long-wavelength bands, a variation of 10% in Tdust turns out to produce a change of normalization of up to a factor of 10.

  • When removing short-wavelength bands, a variation of 30% in Tdust induces a change in normalization of up to a factor of 45 “only”.

All Tables

Table 1

Regions of the M-SFR-z space where good stacked SEDs are available (see Figs. 4 and 5 in Magnelli et al. 2014).

Table 2

Herschel PACS and SPIRE 3σ depths in the GOODS-S and COSMOS fields, adopted in Monte Carlo simulations.

All Figures

thumbnail Fig. 1

Effect of varying qpah (left), Umin (center), γ (right), in Draine & Li (2007) models, independently. On top of each panel, the values of the frozen parameters are quoted. The orange squares at the bottom indicate the position of the chosen rest-frame bands used to generate artificial catalogs.

Open with DEXTER
In the text
thumbnail Fig. 2

Relative uncertainty of dust masses (defined as σ(Mdust) /Mdust), based on DL07 SED modeling, as a function of position in the z-M-SFR space. The two diagrams belong to two independent simulations, obtained noise levels of individual detections in GOODS-S and COSMOS, applied to artificial photometry based on the SEDs by Magnelli et al. (2014). For reference, black lines denote the position of the main sequence (MS, solid line) of star formation (Elbaz et al. 2011) and ± 4, ± 10 MS levels (dashed).

Open with DEXTER
In the text
thumbnail Fig. 3

Relative uncertainty of dust masses (defined as σ(Mdust) /Mdust) at the GOODS-S depth, based on DL07 SED modeling, as a function of redshift and SFR. Color coding is as in Fig. 2. The dotted, solid, and dashed lines indicate the position of the MS of star formation at different values of stellar mass (Elbaz et al. 2011). These values were obtained from the dependence of the sSFRMS on redshift by adopting log (M/M) = 10.0, 10.5, and 11.0, respectively.

Open with DEXTER
In the text
thumbnail Fig. 4

Comparison of input and Monte Carlo (MC) output Mdust in DL07 simulations. Color coding is based on redshift, stellar mass, specific SFR, and Mdust relative uncertainty. As each dot belongs to an individual object, there is a general overlap of colors (green possibly hiding pink). The two diagrams belong to two independent simulations obtained with GOODS-S and COSMOS depths. Data points indicated in red have Mdust/σ(Mdust) < 3.

Open with DEXTER
In the text
thumbnail Fig. 5

Systematics on dust masses (defined as (OUT-IN)/OUT) as a function of position in the z-M-SFR space, in the simulation with GOODS-S depth. Black diagonal lines indicate the position of the MS (Elbaz et al. 2011) and ± 4, ± 10 MS levels (dashed).

Open with DEXTER
In the text
thumbnail Fig. 6

Same as Fig. 4, for GOODS-S only and for MC results obtained using DL07-based photometry instead of DH02-based photometry.

Open with DEXTER
In the text
thumbnail Fig. 7

Comparison of input and Monte Carlo (MC) output Tdust in MBB simulations at the COSMOS depth, with fixed β. Color coding is based on redshift, stellar mass, specific SFR, and Mdust relative uncertainty.

Open with DEXTER
In the text
thumbnail Fig. 8

Study of possible biases on Tdust in a MBB fit and dependence of DH02 templates (e.g., Magnelli et al. 2014). DH02 templates are characterized by their α parameter. Black filled squares and open circles represent the values of Tdust associated with each DH02 template by fitting it with a MBB model with β = 1.5 or 2.0, respectively. The light blue grid maps the variation of Tdust as a function of redshift for each given template. The inset includes a zoom on a single template, aimed at showing the details of the dependence of the Tdust derivation on redshift. Each redshift point was artificially shifted by a small amount along the x-axis to avoid overlapping. Black crosses refer to the result obtained by applying Wien’s displacements law to DH02 templates.

Open with DEXTER
In the text
thumbnail Fig. 9

Relative uncertainty on DL07 free parameters, as obtained with the rest-frame simulation. For simplicity, we show only the case with α = 2, Umax = 106, Umin = 0.7, qPAH = 0.47, and γ = 0.01−0.11. The left diagram is obtained using ten photometric bands from 8 to 850 μm; the right diagrams uses only six bands between 8 and 250 μm. Left/right columns belong to the case obtained by removing long- and short-wavelength bands. Color coding is based on the average photometric uncertainty as computed over all available bands. See Sect. 6.1 for more details.

Open with DEXTER
In the text
thumbnail Fig. 10

Relative error on dust temperature (top) and dust mass (bottom) obtained with MBB simulations, when erasing long-wavelength bands (left) and short-wavelength bands (right). Color coding is based on the average relative photometric uncertainty computed on the available bands, and is explained in the legend. See text for more details on how to read these diagrams.

Open with DEXTER
In the text
thumbnail Fig. 11

Relative uncertainty on Mdust for MBB (left) and DL07 (right) fits to the FIR SEDs of GOODS-N (blue) and GOODS-S (red) sources. Color coding is based on the maximum S/N of the available photometry (at λobs ≥ 100μm) for each object, ranging from a value of 3.0 to 1000.

Open with DEXTER
In the text
thumbnail Fig. 12

Comparison of MBB and DL07 dust estimates. The MBB fitting adopts β = 2.08 and a minimum rest-frame wavelength of 50 μm in the fit. The DL07 estimates are obtained with the whole 8500 μm photometric set for a maximum of nine bands, and DL07 parameters are limited as prescribed by Draine et al. (2007). Left: individual GOODS-N/S sources detected by Herschel. Right: stacked photometry by Magnelli et al. (2014) color coded on the basis of the maximum available rest-frame wavelength (see also inset).

Open with DEXTER
In the text
thumbnail Fig. 13

Comparison of three different Mgas estimates: the first is based on CO (rescaled to a common αCO,MW = 4.36, including a metallicity correction; Genzel et al. 2012, 2013, 2015); the second based on the scaling of depletion times, τdep; and the third derived from Mdust, using the δGDR-Z and M-Z relations in the PP04 metallicity scale. Left: datasets include GOODS-N (blue filled circles) and GOODS-S (red filled squares) 160 μm-detected sources; Herschel-detected SMGs (green filled pentagons; Magnelli et al. 2012a); CO-detected galaxies (different blue symbols, see below); and the Magnelli et al. (2014) stacked points (gray crosses). The right-hand panels only include CO-detected sources, namely: PHIBSS galaxies (blue triangles, Tacconi et al. 2013); BzK galaxies (blue squares, Daddi et al. 2010); other star-forming galaxies (blue filled circles, Magnelli et al. 2012b); lensed galaxies (empty blue circles, Saintonge et al. 2013); and SMGs (green upside down triangles, Bothwell et al. 2013). The dashed and dotted lines in the lower-right panel mark the ± 0.3 and ± 0.6 dex deviations from the 1:1 locus (solid line).

Open with DEXTER
In the text
thumbnail Fig. 14

Comparison of Mgas as based on dust masses and on the scaling of τdep, as a function of distance from the main sequence of star formation. The Whitaker et al. (2014) definition of the MS has been adopted. Symbols are as in Fig. 13.

Open with DEXTER
In the text
thumbnail Fig. 15

Same as Figs. 13 and 14 for GOODS-S/N sources only, limiting the maximum available rest-frame wavelength to μm.

Open with DEXTER
In the text
thumbnail Fig. 16

Gas-to-dust mass ratio, δGDR, of CO-detected galaxies with rest-frame λmax ≥ 200μm, computed by combining CO-based gas masses and DL07-based dust masses. Symbols are as in Fig. 13. Left: δGDR as a function of metallicity. The red lines represent the local relation by Magdis et al. (2012) and its 0.15 dex scatter. Right: δGDR as a function of other relevant parameters: redshift z; depletion times based on CO; M (Wuyts et al. 2011a); and .

Open with DEXTER
In the text
thumbnail Fig. 17

Results of DL07 SED fitting to galaxies detected by Herschel in the GOODS fields, including artificial ALMA band 6 photometry. Color coding is based on the number of available photometric bands, ranging from only four bands (lighter colors) to eight bands (darker colors). This number reflects the maximum rest-frame wavelength available as well as redshift dependencies. Top: relative uncertainties on dust mass, obtained only with the real photometry (see also Fig. 11). Bottom: the same quantity obtained when adding artificial photometry in the ALMA band 6 (~1100 μm).

Open with DEXTER
In the text
thumbnail Fig. 18

Comparison of the DL07-based Mgas estimate obtained, including artificial ALMA band 6 photometry to expectations, based on the Scoville et al. (2015) recipe, and applied to the same artificial ALMA band 6 data. Color coding is based on the number of available photometric bands (see Fig. 17). Left: direct comparison: error bars only include statistic noise. Right: ratio of the two estimates as a function of redshift and DL07-based gas mass.

Open with DEXTER
In the text
thumbnail Fig. 19

Further results of DL07 SED fitting including artificial ALMA band 6 photometry. Color coding is based on the number of available photometric bands (see Fig. 17). Left: exposure time estimate (without overheads) for ALMA band 6 observations of Herschel-detected galaxies. The texp estimate for individual sources (upper bigger panel, dots) are shown and summed up to compute the cumulative texp needed to observe all sources above a given value of Mdust (solid histograms). Dotted histograms represent the number distribution of sources. Right: distribution of ALMA band 6 expected fluxes, as a function of dust mass, for galaxies detected by Herschel in the GOODS-S field. The horizontal dashed lines correspond to 0.38 and 0.18 [mJy] flux limits (dark and light gray, respectively). The right side and bottom panel show the projected 1D distributions. The hatched dark and light gray histograms only include sources above the flux limits indicated.

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

Comparison of input and output Mdust in DL07 simulations obtained with ten photometric bands. Left/right columns present catalog entries obtained removing bands from the long-wavelength and short-wavelength side of SEDs, respectively. Top/middel/bottom diagrams show results obtained with specific values of two parameters and sampling the 3rd vary within the ranges recommended by Draine et al. (2007, see also Sect. 3.3. When fixed, the adopted values are: qPAH = 0.47, Umin = 0.7, and γ = 0.11.

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

Comparison of input and output Tdust (top) and Mdust (bottom) in MBB simulations, split in bins of input dust temperature. Left/right panels again refer to cases with long- or short-wavelength bands removed, respectively. Color coding is the same as in Fig. 10

Open with DEXTER
In the text
thumbnail Fig. B.2

Direct comparison of input/output Tdust vs. Mdust in MBB simulations, split in bins of input dust temperature. Left/right panels again refer to cases with long- or short-wavelength bands removed, respectively. Color coding is based on the number of available bands; from two (light pink) to eight (red).

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.