Benchmarking Dust Emission Models in M101
Jeremy Chastenet, Karin Sandstrom, I-Da Chiang, Brandon S. Hensley, Bruce T. Draine, Karl D. Gordon, Eric W. Koch, Adam K. Leroy, Dyas Utomo, Thomas G. Williams
DDraft version February 24, 2021
Typeset using L A TEX twocolumn style in AASTeX63
Benchmarking Dust Emission Models in M101
J´er´emy Chastenet , Karin Sandstrom , I-Da Chiang ( 江 宜 達 ) , Brandon S. Hensley , Bruce T. Draine , Karl D. Gordon ,
3, 4
Eric W. Koch , Adam K. Leroy , Dyas Utomo , andThomas G. Williams Center for Astrophysics and Space Sciences, Department of Physics, University of California, San Diego9500 Gilman Drive, La Jolla, CA 92093, USA Princeton University Observatory, Peyton Hall, Princeton, NJ 08544-1001, USA Space Telescope Science Institute, 3700 San Martin Drive, Baltimore, MD, 21218, USA Sterrenkundig Observatorium, Universiteit Gent, Gent, Belgium University of Alberta, Department of Physics, 4-183 CCIS, Edmonton AB T6G 2E1, Canada Department of Astronomy, The Ohio State University, 4055 McPherson Laboratory, 140 West 18th Ave, Columbus, OH 43210, USA Max Planck Institut f¨ur Astronomie, K¨onigstuhl 17, 69117 Heidelberg, Germany (Received; Revised; Accepted February 24, 2021)
Submitted to AJABSTRACTWe present a comparative study of four physical dust models and two single-temperature modifiedblackbody models by fitting them to the resolved WISE,
Spitzer , and
Herschel photometry of M101(NGC 5457). Using identical data and a grid-based fitting technique, we compare the resulting dust andradiation field properties derived from the models. We find that the dust mass yielded by the differentmodels can vary by up to factor of 3 (factor of 1.4 between physical models only), although the fits havesimilar quality. Despite differences in their definition of the carriers of the mid-IR aromatic features,all physical models show the same spatial variations for the abundance of that grain population. Usingthe well determined metallicity gradient in M101 and resolved gas maps, we calculate an approximateupper limit on the dust mass as a function of radius. All physical dust models are found to exceed thismaximum estimate over some range of galactocentric radii. We show that renormalizing the models tomatch the same Milky Way high latitude cirrus spectrum and abundance constraints can reduce thedust mass differences between models and bring the total dust mass below the maximum estimate atall radii.
Keywords:
Spectral energy distribution(2129); Interstellar dust(836); Dust continuum emission(412);Gas-to-dust ratio(638); Metallicity(1031); M101 (NGC 5457) INTRODUCTIONDust grains play key roles in processes that shape theinterstellar medium (ISM) and galaxy evolution. Theyrelease photo-electrons that participate in heating gas(e.g. Wolfire et al. 1995; Weingartner & Draine 2001a;Croxall et al. 2012), they shield dense molecular cloudsfrom stellar UV radiation and aid their collapse (e.g. Fu-magalli et al. 2010; Byrne et al. 2019), they catalyze a
Corresponding author: J´er´emy [email protected] number of chemical reactions, and offer surface area forthe production of H (e.g. Bron et al. 2014; Castellanoset al. 2018; Thi et al. 2020, see also a review by Wakelamet al. 2017). It is therefore critically important to un-derstand dust properties and abundance, and how dustaffects these processes.One main way to infer dust properties in externalgalaxies is to interpret infrared (IR) spectral energy dis-tributions (SEDs) with the aid of dust models. Thenear-to-mid-IR part of the spectrum is dominated bythe emission of stochastically heated grains that do notachieve a steady-state equilibrium with the incident ra-diation field. At longer wavelengths, the emission is al- a r X i v : . [ a s t r o - ph . GA ] F e b Chastenet, Sandstrom et al. most entirely due to grains in the thermal equilibrium,with a large enough radius to be constantly receivingand emitting photons. In this regime, the steady-stategrain temperature is set by the strength of the incidentradiation field.Modified blackbody models are a convenient paramet-ric representation of the emission from large grains inthermal equilibrium. As such, they provide good fitsto the far-IR SED and yield satisfactory dust massesif correctly calibrated (e.g. Bianchi 2013), since thesegrains contain most of the dust mass. Large grains inthermal equilibrium are reasonably well described by asingle temperature, in which case their emission can berepresented by a blackbody radiation, B ν , modified byan effective grain opacity, κ ν ( λ ), so that the grain emis-sion I ν ∝ κ ν ( λ ) B ν ( λ, T ). The opacity as a functionof wavelength is often described with a power-law withspectral index, β , such that κ ν = κ ( ν/ν ) β . Severalvariations to the modified blackbody model have beenused to fit the far-IR SED, for example multiple dustpopulations having different temperatures and β , or adifferent functional form for κ ν such as a broken power-law.Physical dust models aim to reproduce the IR emis-sion, extinction, and depletions, among other observa-tions, with a self-consistent description of dust prop-erties. Building these models requires specifying grainsizes, shapes, and chemical composition, which lead tothe optical properties and heat capacities. These grainpopulations are then illuminated by a radiation fieldwith a specified intensity and spectrum. Once the ra-diation field is modeled, one can compare the predictedand observed dust emission. Dust extinction measuredtowards specific lines-of-sight (not high A V ) helps con-strain the size distribution of grains, their composition,and total mass relative to H (e.g. Weingartner & Draine2001b). Depletion measurements provide important lim-its on the elemental abundance locked in dust grains(e.g. Jenkins 2009; Tchernyshyov et al. 2015). Experi-mentally, many studies rely on material thought to beISM dust analogs to provide laboratory constraints onoptical properties and heat capacities (e.g. Richey et al.2013; Demyk et al. 2017; Mutschke & Mohr 2019).Both physical and modified blackbody models are al-most always calibrated in the Milky Way (MW), wherethe the relevant observables of dust are well-constrained.This includes measurements of the diffuse emission andextinction per H column of the ISM at high Galacticlatitudes also referred to as the Milky Way cirrus (e.g.Compi`egne et al. 2011; Guillet et al. 2018). The cirrusis also a unique place where the interstellar radiationfield that heats dust grains, a necessary component to constrain models, can be estimated. The radiation fieldmeasured by Mathis et al. (1983) at the galactocentricdistance D G = 10 kpc is generally used to describe thestarlight heating dust grains.The stringency with which each model follows the el-emental abundances locked in dust grains from deple-tion studies varies. Some models use them as strict con-straints (e.g. Gordon et al. 2014). On the other hand,most physical models allow more flexibility in the massof metals locked up in dust grains, to more closely matchother, better constrained observables. For example theZubko et al. (2004) dust models follow the depletion con-straints strictly, while the Weingartner & Draine (2001c)can require up to ∼
30% (assuming Si/H = 36 ppm inthe model) more silicon than observed in the cold neutralmedium. However, the latter reproduce the observed ex-tinction to a better degree than the former.With the increasing number of observational con-straints on dust models, the complexity of physical dustmodels has grown. One of the earliest dust models de-scribed grains as a single mixture with a power-law sizedistribution (Mathis et al. 1977). Later on, very smallgrains known as the Polycyclic Aromatic Hydrocarbons(PAHs) were suggested as responsible for the mid-IRemission features (Leger & Puget 1984; Allamandolaet al. 1985), and included in dust models. While somedust models consider them to be defined by their size (inthe model description; e.g. Draine & Li 2007; Compi`egneet al. 2011; Galliano et al. 2011), other models identifythe mid-IR features carriers in the form of aromatic-richmantles onto amorphous grains (e.g. Jones et al. 2013,and their hydrogenated amorphous hydrocarbons com-ponent).The precise nature of large grains is also uncertain.The presence of amorphous silicate material in grainswas demonstrated early on by conspicuous absorptionfeatures (Mathis et al. 1977; Kemper et al. 2004) andincluded in dust models. But their exact compositionremains an active research topic (e.g. Zeegers et al.2019). While some models have a strong focus towardsreproducing the observations (e.g. Draine & Li 2007,“astrosilicates”), other models are closely tied to newlaboratory data (e.g. Jones et al. 2013, olivine and py-roxene). Finally, with the growing amount of far-IRpolarization data, new models emerge and take into ac-count this important grain property (e.g Guillet et al.2018; Hensley & Draine 2021). Future missions and in-struments are being developed to focus on polarization,and will bring new constraints to dust properties (e.g.SOFIA/HAWC+: Harper et al. 2018).There are now several physical dust models available,that are all different in some—sometimes small—ways. enchmarking Dust Emission Models ∼ . Table 1.
M101 (NGC 5457) properties.Property Value ReferenceRight ascension 14:03:12.6 Makarov et al. (2014) a Declination +54:20:57 Makarov et al. (2014) a Inclination 30 ◦ de Blok et al. (2008) b Position angle 38 ◦ Sofue et al. (1999)r . (cid:48) Makarov et al. (2014) a Distance 6.7 Mpc Tully et al. (2009) a from the HyperLEDA data base; http://leda.univ-lyon1.fr/ b Note the difference with the HyperLEDA database value(16 ◦ ). metallicity gradients of M101 also offer an independentroute to put an upper limit on its dust content. Thegalaxy has detailed measurements of metallicity fromauroral lines (e.g. Croxall et al. 2016; Berg et al. 2020),which puts good constraints on the gas phase metalabundance. It has also been extensively targeted fordeep H I and CO observations (Walter et al. 2008;Leroy et al. 2009), allowing us to account for all thegas (modulo any limitations in the ability of the COand 21 cm lines to trace the gas). Combining these,we can evaluate the impact of model choice on deriveddust properties across a wide range of environments withwell-understood metal and gas content.In Section 2 we present the photometry used to samplethe IR emission of M101. Section 3 describes the phys-ical and modified blackbody dust models and Section 4the technical aspects of the emission fitting technique.The results of these fits are presented in Section 5 withdiscussions of the differences in dust properties yieldedby the dust models. Finally in Section 6 we analyze thecalibration differences in the dust models themselves. DATAWe present our adopted distance and orientation pa-rameters for M101 in Table 1.In order to derive the dust properties in M101, weperform fits to its IR SED, comprised of measurementsat 16 different photometric bands: • µ m from the Wide-field In-frared Survey Explorer (WISE; Wright et al. 2010).We use the maps delivered by Leroy et al. (2019) inthe z = 0 Multiwavelength Galaxy Survey, alreadyconvolved to a 15 (cid:48)(cid:48) -FWHM Gaussian resolution; Chastenet, Sandstrom et al. • µ m from the Infrared ArrayCamera instrument (IRAC; Fazio et al. 2004), and24 and 70 µ m from the Multiband Imaging Pho-tometer for Spitzer instrument (MIPS; Rieke et al.2004) , both on-board the Spitzer Space Telescope (Werner et al. 2004). We used the product deliv-ery DR5 from the Local Volume Legacy surveyfor IRAC and MIPS maps (Dale et al. 2009); •
70, 100, and 160 µ m from the Photoconductor Ar-ray Camera and Spectrometer instrument (PACS;Poglitsch et al. 2010) and 250, 350, and 500 µ mfrom the Spectral and Photometric Imaging Re-ceiver instrument (SPIRE; Griffin et al. 2010),both on-board the
Herschel Space Observatory (Pilbratt et al. 2010). We downloaded the scansfrom the KINGFISH program (Kennicutt et al.2011) in the Herschel Science Archive, and pro-cessed them from L0 to L1 with
HIPE (v. 15; PACSCalibration v. 77; SPIRE Calibration v. 14.3; Ott2010) and
Scanamorphos (v. 25; Roussel 2013).2.1.
Data Processing
We correct the images for extended source calibrationappropriately and convert them to the same units inthe following way. We apply extended-source correctionfactors for the IRAC images. These are multiplicativecorrections of 0.91, 0.94, 0.695, and 0.74 at 3.6, 4.5, 5.8and 8.0 µ m, respectively, as suggested by the IRAC In-strument Handbook . For SPIRE, we adopt calibrationfactors of 90.646, 51.181, and 23.580 MJy/sr/(Jy/beam)at 250, 350, and 500 µ m, respectively, as suggested bythe SPIRE Handbook . These factors convert the datato MJy/sr from Jy/beam and also correct from pointsource to extended source calibration. The PACS datain units of Jy/pixel are converted to units of MJy/srusing the image pixel size (contained in the headers).We then process all the maps following the same steps,as described below. We remove a background in eachimage at their native resolution, by fitting a 2-D planeusing regions identified not to have significant galaxyemission. We find these regions with the following pro-cedure: 1) We use the r radius (r ∼ (cid:48) ) to first maskthe galaxy (this covers the visible SPIRE 500 emission We do not use MIPS 160 (e.g. as opposed to Aniano et al. 2020)to gain back some resolution (MIPS 160 PSF is ∼ (cid:48)(cid:48) ), withoutlosing wavelength coverage. This does not lead to major differ-ences. https://irsa.ipac.caltech.edu/data/SPITZER/LVL/LVL DR5 v5.pdf https://irsa.ipac.caltech.edu/data/SPITZER/docs/irac/iracinstrumenthandbook/29/ http://herschel.esac.esa.int/Docs/SPIRE/spire handbook.pdf completely); 2) We measure the median of the remain-ing pixels and the standard deviation, and clip all ofthose that are 3 standard deviations above that median;3) We iterate that clipping until the medians at itera-tions i and i + 1 differ by less than 1%. This clippingis only done for the purpose of measuring a backgroundlevel, and we do not keep this mask applied to the datafor the following steps. Table 2 lists the final standarddeviations of the background pixels in each bands.Each map is then convolved to the SPIRE 500 PSF(FWHM ∼ (cid:48)(cid:48) ), using convolution kernels from Ani-ano et al. (2011). Finally, all the maps are aligned andprojected onto the astrometric grid of the SPIRE 500image. In Figure 1, we show the 16 bands that we useto model the dust emission from 3.4 to 500 µ m.The final pixel size (9 (cid:48)(cid:48) ) oversamples the SPIRE 500beam size, to which all data are convolved. We take thisinto account by correcting by (cid:112) N pix /N beam when cal-culating uncertainties on quantities that use the valuesof multiple pixels.2.2. Ancillary data: neutral and molecular gas
In Section 6, we use gas surface density and metallic-ity measurements in additional to IR SEDs. We use thesame data as Chiang et al. (2018) to get the total gas sur-face density, Σ gas , combining gas surface density mapsfrom H I
21 cm and CO (2 →
1) emission converted toH . The H I
21 cm line data is from the THINGS collab-oration (The H I Nearby Galaxy Survey; Walter et al.2008). The map is converted from integrated intensityto surface density assuming optically thin 21 cm emis-sion, following Walter et al. (2008, Equ. 1 and 5, andmultiplying by the atomic mass of hydrogen).We use CO (2 →
1) emission from HERACLES (Leroyet al. 2009), and convert it to H surface density, as-suming a line ratio R = 0 .
7, and a α CO conversionfactor, in units of M (cid:12) /pc (K km/s) − , following twoprescriptions (both including helium): the first one isrepresentative of the MW CO-to-H conversion , α MWCO = 4 . , (1)and the second one follows Bolatto, Wolfire, & Leroy(2013), α BWL13CO = 2 . (0 . /Z (cid:48) ) (Σ ) − γ with (cid:40) γ = 0 . ≥ γ = 0 otherwise (2) A standard X CO = 2 × (K km/s) − is assumed for the col-umn density conversion factor. The mass of helium and heavierelements have been accounted for in α CO . enchmarking Dust Emission Models WISE 3.4 IRAC 3.6 IRAC 4.5 WISE 4.6IRAC 5.8 IRAC 8.0 WISE 12 WISE 22MIPS 24 MIPS 70 N G C PACS 70 PACS 100 h m s s m s s m s s D e c ( J ) PACS 160 SPIRE 250 SPIRE 350 SPIRE 500
Figure 1.
Emission of M101 (NGC 5457) from 3.4 to 500 µ m, in all the bands used in this study. The maps show the finalversion of each map: after extended source correction (IRAC and SPIRE bands), unit conversion (PACS and SPIRE bands),background removal, convolution to SPIRE 500 PSF ( ∼ (cid:48)(cid:48) ) and regridding. On the MIPS 70 panel, we also show: the locationof the pixel for the fit example (white cross; Figure 4), the region used for the IRS measurement (white rectangle; Section 5.6),and the contours for the 3 σ detection threshold. where Z (cid:48) is the metallicity relative to solar metallic-ity , Σ is the total surface density map in units of100 M (cid:12) /pc . The gas maps are convolved to a 41 (cid:48)(cid:48) Gaus-sian PSF (Chiang et al. 2020); since we plot the radialprofile of the gas maps by averaging in growing annuli,we find that the minor differences between the dust andgas maps resolutions are negligible. We measure an un-certainty of ∼ . →
1) map,and an uncertainty of ∼ . (cid:12) /pc for the H I
21 cmmap . We build a radial profile by averaging Σ gas ingrowing annuli from the center out to r .2.3. Ancillary data: metallicity Z (cid:48) = 10( (12+log(O / H)) − (12+log(O / H) (cid:12) ) ) with 12+log(O/H) (cid:12) =8 . The rms per channel is ∼ .
46 mJy/beam. Assuming σ z , gas =11 km/s (Leroy et al. 2008), the uncertainty is 0.96 M (cid:12) /pc . The5% calibration error (Walter et al. 2008) becomes significant inthe dense regions. Metallicity measurements 12 + log (O / H) are takenfrom the CHAOS survey (Berg et al. 2020), derived from72 H II regions of M101. In particular, we use their fittedmetallicity gradient12 + log (O / H) = 8 . ± . − . ± .
07 ( r/ r ) (3)to convert it to a radial profile of metal surface density(see Section 6.1). We adjusted the slope given by Berget al. (2020, 0.90) to account for the different r usedhere. Berg et al. (2020) find that the scatter in individ-ual region metallicities around the measured gradient is ∼ II regions, andfound minimal depletion of oxygen ( (cid:46) . Chastenet, Sandstrom et al.
Table 2.
Band-related details. References to the σ cal and σ sta coefficients: WISE: http://wise2.ipac.caltech.edu/docs/release/prelim/expsup/sec4 3g.html; IRAC: Reach et al.(2005) and https://irsa.ipac.caltech.edu/data/SPITZER/docs/irac/iracinstrumenthandbook/29/; MIPS: Engelbrachtet al. (2007) and Gordon et al. (2007); PACS: M¨uller et al.(2011) and Balog et al. (2013); SPIRE: Griffin et al. (2010)and Bendo et al. (2013).Band σ bkg a m cal b m sta c − MJy/sr % %WISE 3.4 0.170 2.4 10.0IRAC 3.6 0.149 9.0 1.5IRAC 4.5 0.107 6.0 1.5WISE 4.6 0.095 2.8 10.0IRAC 5.8 0.109 30 1.5IRAC 8.0 0.085 26 1.5WISE 12 0.114 4.5 10.0WISE 22 0.055 5.7 10.0MIPS 24 0.049 4.0 0.4MIPS 70 5.24 10.0 4.5PACS 70 10.7 10.0 2.0PACS 100 10.1 10.0 2.0PACS 160 7.95 10.0 2.0SPIRE 250 3.96 8.0 1.5SPIRE 350 2.73 8.0 1.5SPIRE 500 1.82 8.0 1.5 a The background standard deviation in thepixels considered to measure the back-ground covariance matrix, once the mapsare background-removed, convolved and pro-jected. b m cal is the error on the calibration used ineach instrument. The large errors of theIRAC bands are from the extended sourcecorrections, which we consider to be corre-lated calibration errors. c m sta measures the stability of an instrument,i.e. the scatter when measuring the same sig-nal. the case in Berg et al. 2020), but given the widespreaduse of O/H for extragalactic studies, we use the typical12 + log(O / H) tracer. DUST EMISSION MODELSThe goal of this study is to investigate the differencesin dust properties derived using physical and modifiedblackbody models. In this Section, we present the char- acteristics of each model. The parameters we use to fitthe IR SEDs are presented in Table 3. A summary tableof the following information is presented in Appendix A.3.1.
Modified blackbodies
We use single-temperature modified blackbodies mod-els. For these two models, we only fit photometry from100 to 500 µ m. At λ < µ m, stochastically heatedgrains contribute to the emission and are not well mod-eled by a modified blackbody. Assuming the opticallythin case, the dust emission I ν is described as I ν ( λ ) = κ ν ( λ ) Σ dust B ν ( λ, T dust ) , (4)where B ν ( λ, T dust ) is the Planck function at wavelength λ (in MJy/sr), T dust the dust temperature, Σ dust thedust surface density, and κ ν the opacity. We use thesimple power-law opacity (Equ. 5) and a broken power-law opacity (Equ. 6), which we normalize at λ = 160 µ m.The broken power-law model presented the best resultsin terms of quality of fits in the study by Chiang et al.(2018) and yielded physically reasonable Σ dust and T dust values in their study.We follow Gordon et al. (2014) and Chiang et al.(2018) and calibrate the opacity values for each model.We fit the modified blackbody model to the dust emis-sion per H column of the MW high-latitude cirrus de-scribed in Chiang et al. (2018) to derive the opacity. Theabundance constraints are based on a depletion strengthfactor typical for lines-of-sight with similar N H to theMW cirrus, e.g., F ∗ = 0 .
36 (Jenkins 2009). This sets theallowed dust mass per H atom. By fitting the tempera-ture for the MW cirrus we can then derive the opacityfor each modified blackbody model. More details on theopacity and comparison to the physical models can befound in Section 3.2.3.1.1.
Simple-Emissivity (SE)
In this case, the opacity is described as a single power-law: κ ν = κ ν (cid:18) νν (cid:19) β , (5)where β is the spectral index. For all blackbody models,we fix λ = 160 µ m, and ν = c/λ . The free parametersfor this model are then the dust surface density, Σ dust ,the dust temperature, T dust , and the dust spectral index, β . In this model the calibrated κ ν = 10 . ± .
42 cm /g(Chiang et al. 2018), from fitting the high latitude cirrusas described above.3.1.2. Broken-Emissivity (BE)
The broken-emissivity model stemmed from the iden-tification of the sub-millimeter excess (Gordon et al. enchmarking Dust Emission Models κ ν = κ ν (cid:16) νν (cid:17) β λ < λ c κ ν (cid:16) ν c ν (cid:17) β (cid:16) νν c (cid:17) β λ ≥ λ c (6)where λ c is the wavelength at which the opacity changes,and ν c its equivalent frequency. Following Chiang et al.(2018), we fix β = 2, and λ c = 300 µ m. The free pa-rameters are then the dust surface density, Σ dust , thedust temperature, T dust , and the second dust spectralindex, β . The calibrated opacity, κ ν , for this modelis 20 . ± .
97 cm /g (Chiang et al. 2018), as describedabove. 3.2. Opacity calibrations
The physical dust models used in our study have opac-ities that are set by each individual model’s calibra-tion procedure. All models use similar constraints fromthe MW high latitude cirrus (described in more detailin Appendix A). For ease of comparison to the opaci-ties we have derived for the modified blackbody models(10.1 cm /g for the simple-emissivity and 20.7 cm /g forthe broken-emissivity at 160 µ m), we list the opacity inthe physical models below all scaled to 160 microns forcomparison: • Draine & Li (2007): from the Weingartner &Draine (2001c) model updated in DL07 we report κ = 10 . /g (see also Bianchi 2013); • Compi`egne et al. (2011): from the work of Bianchi(2013) we report κ = 12 . /g; • Jones et al. (2017, THEMIS): from thework of Galliano et al. (2018) we report κ = 14 . /g; • Hensley & Draine (2021): the models uses κ ∼ .
95 cm /g (B. Hensley; priv. comm.).For all models, we use the listed opacity regardless ofthe specific environment in M101 we are studying. Fewof the currently available physical dust models have beencalibrated in any other environment than the MW cirrus(although the Draine & Li 2007 model does have Smalland Large Magellanic Cloud-like calibrations, thoughthey are not widely used even for these very galaxies;Sandstrom et al. 2010; Chastenet et al. 2019). Indeed,it is standard in current extragalactic applications to ap-ply MW cirrus R V = 3 . F ∗ = 0 .
36 value doesnot describe the depletion in the ISM over the full rangeof column densities probed in galaxies. However, for thepurposes of our comparative study of widely used dustmodels, we proceed by using the R V = 3 . Physical dust models
Physical dust models assume a composition, density,and shape for the dust grains and adopt heat capac-ities and optical properties from laboratory and the-oretical studies that are appropriate for such materi-als. For simplicity, most models assume spherical grainsor planar molecules for PAHs. In the case of PAHs,the grains/molecules are additionally described by anionization state, which changes the absorption cross-sections as a function of wavelength. The temperatureand emission of a grain of a given size, shape and com-position in a radiation field with a specified intensityand spectrum can then be calculated analytically (e.g.Draine & Lee 1984; Desert et al. 1986).The full dust population is represented by a grain sizedistribution and abundance relative to H for the speci-fied compositions. Physical dust models are calibratedby adjusting the grain size distributions and dust massper H to simultaneously match observations of extinc-tion, emission, and abundances (and more recently, po-larization) in a location where the underlying radiationfield that is heating the dust is well known. This hasgenerally been taken to be the high-latitude MW cirrus,where the radiation field intensity and spectrum are ap-proximately given by the Mathis et al. (1983) model forthe Solar neighborhood. The degree to which the mod-els must adhere to the somewhat uncertain abundanceconstraints varies model to model. For example, themodified blackbody models from Gordon et al. (2014)use the depletion measurements in the MW as a strictlimit. However, most physical models allow the final el-ement abundances to vary from depletions, using themonly as loose guide.In the following, we use 4 physical dust models:Draine & Li (2007), Compi`egne et al. (2011), THEMIS(Jones et al. 2017) and Hensley & Draine (2021). Herewe briefly describe these models, the key differencesbetween them. Details on their respective calibrationmethodologies can be found in Appendix A.3.3.1.
Draine & Li (2007)
Chastenet, Sandstrom et al.
In the Draine & Li (2007, hereafter DL07) model, dustis comprised of PAHs, graphite grains and amorphoussilicate grains. It stems from the original models pre-sented in Li & Draine (2001a). The carbonaceous dustoptical properties are adopted from Li & Draine (2001b)with updates to the PAH cross-sections and form of thegrain size distribution. A balance between ionized andneutral PAHs is assumed, following Li & Draine (2001b).The optical properties of silicate material are adoptedfrom the “astrosilicates” in the original model. We donot make use of the mass renormalization in Draine et al.(2014).The mass fraction of PAHs, q PAH , is described as themass of carbonaceous grains with less than 10 carbonatoms with respect to total dust mass. We effectivelyobtain q PAH in % by converting the part-per-million car-bon abundance b C to a PAH fraction, using the reference q PAH = 4 . ≡ b C = 55 ppm.The calibration of the model is described in severalpapers (Draine & Li 2001; Li & Draine 2001b; Wein-gartner & Draine 2001c). We use the DL07 Milky Waymodel, with R V = 3 .
1, which has a fixed ratio of silicateto carbonaceous grains. Details on the calibration canbe found in Appendix A.1.3.3.2.
Compi`egne et al. (2011)
The Compi`egne et al. (2011, hereafter MC11) modelis composed of PAHs, hydrogenated amorphous carbongrains, and amorphous silicate grains. The size distri-bution of the carbonaceous components includes PAHs,small amorphous carbon grains (SamC), and large amor-phous carbon grains (LamC). The PAH cross-sectionsand ionization as a function of size adopted in the modelare based on Draine & Li (2007) with slight modifi-cations to the cross sections of several bands. Amor-phous carbonaceous grains have optical properties fromZubko et al. (2004) and heat capacities from Draine &Li (2001). The amorphous silicates (aSil) have opticalproperties from Draine (2003) and heat capacities fromDraine & Li (2001). Details on the calibration can befound in Appendix A.2.The
DustEM tool allows both ionized and neu-tral PAHs to be fit independently. Because notall models allow that separation, we tie their emis-sion spectra together by summing them, hence keep-ing their ratio constant at roughly 60% neutral and40% ionized . Additionally, we fix the mass frac- DustEM The ratio of ionized and neutral PAHs is size-dependent. Thevalues are set in
DustEM from the
MIX PAH x .DAT files. tions of the large carbonaceous-to-silicate grains suchas ( M LamCdust /M H ) / ( M aSildust /M H ) = 0 . . Since theseshare a very similar far-IR spectral index, their respec-tive emission cannot be properly determined indepen-dently with our wavelength coverage, and would lead todegenerate abundances if fit separately.3.3.3. THEMIS
The Heterogeneous dust Evolution Model for In-terstellar Solids (THEMIS; Jones et al. 2013; K¨ohleret al. 2014; Ysard et al. 2015; Jones et al. 2017)is a core/mantle dust model, consisting of large sil-icate and hydrocarbons, both coated with aromatic-rich particles (Hydrogenated Amorphous Hydrocarbons,HACs). This model defines its dust components focusingstrongly on laboratory data, slightly adjusted to bettermatch observations. We use the diffuse ISM version ofthe model, described in Jones et al. (2013). The amor-phous carbon grain properties are size-dependent, andthere is no strictly independent population of grains re-sponsible for the mid-IR features, carried by aromaticclusters in the form of mantles and very small grains(sCM20). As they grow to larger grains (lCM20) insize, their core becomes aliphatic-rich, coated with anaromatic mantle. The silicate grains (aSilM5) in par-ticular are discussed in K¨ohler et al. (2014). They area mixture of olivine and pyroxene-type material, withnano-inclusions of Fe and FeS. Their mass ratio is keptconstant due to their extreme resemblance in emissionin the considered wavelength range. The dust evolutionmodels (from diffuse to denser medium) are discussedin Ysard et al. (2015), with the impact of aggregatesand thicker mantles on the final abundances. In thediffuse model we use, aSil have a 5 nm mantle, and hy-drocarbons have a 20 nm mantle. In our fitting, thelCM20-to-aSilM5 mass ratio is kept constant, such as( M lCM20dust /M H ) / ( M aSilM5dust /M H ) = 0 . . Details on thecalibration can be found in Appendix A.3.3.3.4. Hensley & Draine (2021)
Rather than employ separate amorphous silicate andcarbonaceous grain components, the Hensley & Draine(2021, hereafter HD21) model invokes a single homo-geneous composition, “astrodust” (Draine & Hensley2020), to model most of the interstellar grain mass. Inaddition to astrodust, the model incorporates PAHs us-ing the cross-sections from Draine & Li (2007) and asmall amount of graphite using the turbostratic graphitemodel presented in Draine (2016). The HD21 model was This value is that from the
DustEM file
GRAIN MC10.DAT . This value is that from the
DustEM file
GRAIN J13.DAT enchmarking Dust Emission Models U and q PAH . FITTING METHODOLOGY4.1.
Making matched model grids
Since each model was developed independently, it isnot always possible to create model grids that have ex-actly the same parameter sampling, because of the lackof parameter equivalence. However we attempt to do soas much as possible.
Radiation field —
We implement identical dust heat-ing in each of the physical models: a fraction γ of thedust mass in a pixel is heated by a power-law distribu-tion of radiation fields with U min < U ≤ U max (Daleet al. 2001), where U is the interstellar radiation field,expressed in units of the MW-diffuse radiation field at10 kpc from Mathis et al. (1983). The remaining frac-tion of dust (1 − γ ) is heated by the minimum radiationfield U min (Draine & Li 2007):1 M dust dM dust dU =(1 − γ ) U min + γ ( α − U − α min − U − α max U − α (7)We fix α = 2 for all models, U max = 10 for the DL07,THEMIS and MC11 models, and U max = 10 for HD21.This choice is constrained by the available parameterranges in the models: the U max values for the DL07 andHD21 models are fixed, and we do not have the freedomto change them; thanks to DustEM , we can adjust U max for MC11 and THEMIS . For THEMIS and the MC11model, which have multiple dust populations that can Using
DustEM and THEMIS, we checked the effect of usingU max = 10 or 10 , for a range of γ values. We find that most ofthe mid-IR bands used in this study are only minimally affected.The difference for IRAC 4.5 and WISE 4.6 can be significant at γ ≥ .
1, while the maximum γ values reach by the model fitsis 0.07 (for DL07). Additionally, these bands are dominated bystarlight, and mostly modeled by the 5,000 K blackbody. be independently heated, each population is heated bythe same radiation field.The minimum radiation field parameter grid is fixedby the values provided in the Draine & Li (2007) model.Thanks to the DustEM tool, we can use the exact samevalues with THEMIS and the MC11 model. These pa-rameter values in HD21 however are not exactly identi-cal to those in U DL07min . We therefore use the spectra withclosest U HD21min values. Although not strictly equal, theradiation field values in Hensley & Draine are within 5%of U DL07min .For the modified blackbody models, we use a singleradiation field intensity and translate it into a dust tem-perature. We use the relationship U ∝ T β +4dust with β = 2 (8)to convert U min to T dust and find the approximatelymatching sampling to use in the blackbody models. Weuse the normalization from Draine et al. (2014), i.e. U = 1 at T dust = 18 K, found using the same radia-tion field from Mathis et al. (1983). Mid-IR feature carriers —
The q PAH parameter iskept strictly identical between DL07 and HD21. Wechoose to use THEMIS and the MC11 models in asimilar fashion. Despite the definition of “PAHs” be-ing different in THEMIS, we parameterize the modelso that the fraction of small grains can vary, keepingthe large-carbonaceous-to-silicate grains ratio constant.The MC11 default model has 4 populations that canvary independently, and we choose to tie together SamC,LamC and aSil, to leave the amount of PAHs as a freeparameter. We explain these choices in more detail inSection 7.
Stellar surface brightness —
In addition to the dustmodel parameters, we use a scaling parameter of a5,000 K blackbody, Ω ∗ , to model the stellar surfacebrightness visible at the shortest wavelengths. This tem-perature is a good approximation as the shortest wave-lengths are nearly on the Rayleigh-Jeans tail. The freeparameter Ω ∗ scales the amplitude of the stellar black-body.In Table 3, we list the final free parameters we use foreach model. There are 5 free parameters for the phys-ical models, and 3 for the modified blackbody models.Figure 2 shows all the dust emission models used in thisstudy. The top-left panel shows the fiducial MW highgalactic latitude diffuse ISM models, labeled ‘GalacticSED’. The IR MW diffuse emission from Compi`egneet al. (2011) is also plotted. The bottom-left panel showsthe physical dust models at the same radiation field U min = 1. The other different panels detail the mod-els: THEMIS is divided in two grain population when0 Chastenet, Sandstrom et al.
Table 3.
Fitting parametersParameter Range Step UnitAll physical modelslog (Σ dust ) [-2.2, 0.1] 0.035 M (cid:12) /pc U min [0.1, 50] Irr a –log ( γ ) [-4, 0] 0.15 –log (Ω ∗ ) [-1, 2.5] 0.075 –Model-specificDL07, HD21: q PAH [0, 6.1] 0.25 %THEMIS: f sCM20 b [0, 0.5] 0.03 –MC11: f PAHs b [0, 0.5] 0.03 –Modified-blackbody modelslog (Σ dust ) [-2.2, 0.1] 0.035 M (cid:12) /pc T dust [12, 35] 0.3 KSE: β , BE: β [-1,4] 0.05 – aU min ∈ { } . b We use a “fraction” parameter so thatΣ X = f X × Σ d , where X = { sCM20 THEMIS ,PAHs
MC11 } . fitted to the SED of M101 (note that we tie the lCM20and aSilM5 population); for the MC11 model we tiethe SamC, LamC and aSil emissions together, and thusthere are two free parameters for the dust mass: f PAH or f sCM20 , and the total dust surface density. The twomodified blackbodies, are shown with their best fit val-ues to the MW diffuse SED: { T dust = 20 . β = 1 . } for the simple-emissivity model, and { T dust = 18 . β = 1 . } for the broken-emissivity model.4.2. Fitting tool
Bayesian fitting with
DustBFF
We use the
DustBFF tool from Gordon et al. (2014) tofit the data with the chosen dust models.
DustBFF is aBayesian fitting tool that uses flat priors for all param-eters. The probability that a model with parameters θ fits the measurements ( S obs ) is given by: P ( S obs | θ ) = 1 Q e − χ ( θ ) / (9)with Q = (2 π ) n det | C | χ ( θ ) = [ S obs − S mod ( θ )] T C − [ S obs − S mod ( θ )] (10) where S X the observed ( X = obs) or modeled ( X =mod) 16 band SEDs used here, and C is the covariancematrix that includes uncertainties from random noise,astronomical backgrounds and instrument calibration(described further below).To create S mod ( θ ) each model spectrum is convolvedwith the spectral responses for the photometric bandsused here. PACS and SPIRE band integrations are donein energy units, whereas all the others are done in pho-ton units, as necessitated by the instrument’s calibrationscheme. We also follow the reference spectra used for thecalibration of each instrument: the MIPS bands requirea reference shape in the form of the 10 K blackbodywhile the other bands use a reference shape as 1 /ν .4.2.2. Covariance matrices
The covariance matrix in the previous equations de-scribes the uncertainties on the measured flux, both dueto astronomical backgrounds, and instrumental noiseand the uncertainties in calibrating the instruments.This takes into account the correlation between the pho-tometric bands due to the calibration scheme and thecorrelated nature of astronomical background signals.We define a background matrix and an instrument ma-trix such that C = C bkg + C ins to propagate the cor-related errors of the background and noise, and of thecalibration uncertainties, respectively. Background covariance matrix C bkg — The “background”in our images encompasses astronomical signals that donot come from the IR emission of the target M101. It isdominated by different objects depending on the wave-length, and can therefore be correlated between bands.For instance, the background from 3.4 to 5.8 µ m ismostly from foreground stars, as well as zodiacal light;from 8 to 24 µ m, it is from evolved stars and back-ground galaxies; in the far-IR, it is dominated by Galac-tic cirrus emission and background galaxies. To includethis uncertainty, we measure this combination of signalsin “background pixels” using the processed data and amasking procedure described in Section 4.2.3.The elements of the background covariance matrix arecalculated as( C bkg ) i,j = (cid:80) N k ( S ki − (cid:104) S i (cid:105) )( S kj − (cid:104) S j (cid:105) )N − , (11)where S kX is the flux of pixel k , in band X , and (cid:104) S X (cid:105) is the average background emission in band X (close to0). The final number of 9 (cid:48)(cid:48) pixels used to measure thecovariance matrix is N = 18 , correlation matrix in Fig-ure 3 (not the elements’ absolute values but the Pearsoncorrelation coefficients; see Gordon et al. 2014). Three enchmarking Dust Emission Models I [ e r g / s / c m / s r / ( H / c m )] DL07; U min = 0.8HD21; U min = 1.6 THEMIS; U min = 1.0MC11; U min = 1.0 Galactic SED
Wavelength [ m]
Wavelength [ m]
SEBE1 10 100 1000
Wavelength [ m] I [ e r g / s / c m / s r / ( H / c m )] U min = 1 Figure 2.
The dust emission models used in this study: Draine & Li (2007), Jones et al. (2017, THEMIS), Compi`egne et al.(2011), Hensley & Draine (2021), simple-emissivity, and broken-emissivity.
Top left: all six models at the U min value that bestfits the calibration SED. Bottom left: the four physical models at U min = 1. The right side panels show each model at the U min value that best fits the calibration SED, and their break-down as they are used in this study. clear sets of positive correlations appear, like previouslyexplained: correlations are due to starlight in the near-IR, evolved stars and MW cirrus in the mid-IR, andbackground galaxies and MW cirrus in the far-IR. Ad-ditionally, some bands may be noise dominated: it isthe case for the PACS 70 band which is only weaklycorrelated with other bands. Instrument calibration matrix C ins — This matrix is calcu-lated as the quadratic sum of the correlated and uncor-related errors for each instrument. The correlated errors(matrix with diagonal and anti-diagonal terms) refer tothe instrument calibration itself ( m cal in Table 2, whilethe uncorrelated (matrix with diagonal terms only) ex-press the instrument stability, or repeatability ( m sta inTable 2). We use the same calibration errors reported inChastenet et al. (2017, see Table 2) for the Spitzer /MIPSand
Herschel bands. The correlated errors for the IRACbands were changed to take into account the uncertain-ties in the extended source correction factors, larger thanthe calibration error themselves. The errors due to re- peatability are unchanged. The errors for WISE bandswere taken from the WISE documentation .The elements of the calibration matrix C ins are calcu-lated “model-by-model” as( m ins ) i,j = S mod i ( θ ) S mod j ( θ ) ( m , i,j + m , i,j ) (12)with particular elements m cal , i,j = 0 if i, j belong to two different instruments; m sta , i,j = 0 if i (cid:54) = j. We fit all the pixels that are above 1 σ of the back-ground values, in all bands. We use these pixels to showparameter maps and radial profiles, while the galaxy-integrated values are calculated for pixels above 3 σ de-tection above the background (black contours in Appen-dices B and C).4.2.3. Stars in the background/foreground http://wise2.ipac.caltech.edu/docs/release/prelim/expsup/sec4 3g.html Chastenet, Sandstrom et al. W I S E . I R A C . I R A C . W I S E . I R A C . I R A C . W I S E W I S E M I P S M I P S P A C S P A C S P A C S S P I R E S P I R E S P I R E W I S E . I R A C . I R A C . W I S E . I R A C . I R A C . W I S E W I S E M I P S M I P S P A C S P A C S P A C S S P I R E S P I R E S P I R E Figure 3.
Background correlation matrix. The dark colorindicates a strong correlation between the bands. The matrixis symmetric and we show only half. The text indicates theastronomical signals that dominates the contoured bands,and that explains their strong correlation.
Here we describe the masking procedure to mea-sure the background covariance matrix in Section 4.2.2.To do so, we use the final images, i.e. background-subtracted, convolved and projected to the same pixelgrid.The covariance matrix elements are calculated withthe assumption of a Gaussian noise. While the assump-tion works well for faint and unresolved stars and thecosmic infrared background galaxies, it is no longer cor-rect if we include bright stars. Bright foreground starsonly must therefore be cut to measure this matrix. Thismasking has the effect of making the approximation ofGaussian noise for the remaining background more cor-rect. Note, however, than they are not masked for thefitting within the boundary of galaxy , but only for thepurpose of the covariance matrix measurement.In the process of masking, we first exclude the regionwithin r , to mask the galactic emission. We then maskthe brightest stars, using the star masks from Leroyet al. (2019) that leverage the known positions of stars We find no pixels showing a bad fit due to a foreground star.The Ω ∗ maps do not show conspicuous peaks, indicating that theforeground stars are not dominant, and the models successfullyfit the galaxy emission. from the GAIA and 2MASS catalogs . To match the fi-nal products, we convolve and regrid these star masks tothe SPIRE 500 resolution. After convolution, the maskvalues are no longer binary 0 and 1, but show inter-mediate values around the position of bright stars. Wemask pixels above 0.15 to exclude these regions wherebright foreground stars contaminate our measurementsfrom the covariance matrix calculation.One argument against the decision to mask brightstars to measure the covariance matrix would be thatthe emission from these stars is an astrophysical signalthat should be taken into account to propagate the noisefrom a band to another. This would require creating anew noise model and significant changes to the fittingmethodology. Rather than major changes to the fittingapproach, we decide to mask the stars to measure thebackground covariance matrix. RESULTSWe investigate the differences in some of the key pa-rameters from dust emission modeling, when the mod-els presented here are all used in an identical fittingframework. All residuals in this analysis are presentedas (Data-Model)/Model . For radial profiles, we use thepixel size as the annulus width and cover from the centerof M101 to r . Appendix B shows the resolved mapsof the fitted parameters for each model. Appendix Cshows the residual maps of each model.5.1. Quality of fits
Figure 4 shows an example of the fits for each model ina single pixel (marked by the cross in Figure 1). In eachpanel we plot the best fit spectrum (colored lines) to thedata SED (empty circles). The residuals are shown inthe colored symbols. Negative residuals mean that themodel overestimates the data. For example, the nega-tive (and decreasing) residuals from 250 to 500 µ m in thephysical models panels are representative of a system-atic overestimation of the data (present also in otherlocations; see Appendix Figures 19–22). In Figure 5,we show the fractional residuals (Data-Model)/Model ineach band, for the pixels above the 3 σ threshold. Ap-pendix D shows the reduced χ for all models.The bulk of the residuals in the short wavelengthbands are within the instrument uncertainties and cali-bration errors. For example, despite larger uncertaintydue to the extended source correction, the IRAC 8band shows residuals mostly within 10%. Below 4 µ m,THEMIS shows a clear offset that may be related to the https://irsa.ipac.caltech.edu/data/WISE/z0MGS/overview.html enchmarking Dust Emission Models I [ M J y / s r ] DL07
10 100-0.50.00.5 R e s i d u a l s HD21
10 1001.010.0 I [ M J y / s r ] THEMIS
10 100-0.50.00.5 R e s i d u a l s MC11
10 100
Wavelength [ m] SE
10 100 BE
10 100
Figure 4.
Example of the best fits in a pixel for the six models used in this study (color lines). The data is shown as theempty symbols, with 1 σ error-bars. The synthetic photometry from the model spectrum is shown with a cross symbol, due tothe band integration, this does not always sit exactly on the spectrum. The location of the measurement is marked by a crossin Figure 1. absence, or low amount, of an ionized component in theHAC population, which leads to enhance the mid-IR fea-tures. It is also worth noting these bands are dominatedby starlight, modeled by a 5,000 K blackbody, which isindependent from the dust model itself. The residualsat 12 µ m are systematically positively offset by less than10%, but all physical models show very narrow residualdistributions. This is in contrast with the broader resid-uals in IRAC 8 and WISE 22 bands, where we can seemore differences between models.All models show differences in the central values of theresidual distribution at 100 µ m. At 160 µ m, all resid-uals overlap and models perform fits of similar quality.At longer wavelengths, significant difference begin to ap-pear. At all far-IR wavelengths, the modified blackbodymodels reproduce the SEDs the best. This is likely be-cause of the additional parameter that can adjust thefar-IR slope ( β in the simple-emissivity model and β in the broken-emissivity model). The SPIRE 250 bandshows symmetrical residuals centered on 0 for the phys-ical models (except THEMIS), while the residuals getprogressively worse at 350 and 500 µ m for all physical models, showing that on average, the modeled far-IRslope of the SEDs is steeper than the data.In the SPIRE 350 and SPIRE 500 bands, the largenumber of pixels underestimated by the models showresiduals much larger than the uncertainties, ruling outstatistical noise and indicating that the models are notable to fit these wavelengths. In the SPIRE 500 band,some of the pixels show the so-called “sub-millimeterexcess” seen in other studies (e.g. Galametz et al. 2014;Gordon et al. 2014; Paradis et al. 2019).5.2. Total Dust Mass and Average Radiation Field
We compute several galaxy averaged quantities: thetotal dust mass, M dust , the dust mass-weighted aver-age radiation field (cid:104) U (cid:105) for the physical models, and themass-weighted average dust temperature (cid:104) T dust (cid:105) for theblackbody models. The average radiation field U is cal-culated for each pixel as U = (1 − γ ) U min + γ × U min ln(U max /U min )1 − U min / U max (13)4 Chastenet, Sandstrom et al. W I S E . DL07 HD21 THEMIS MC11 SE BE I R A C . I R A C . W I S E . I R A C . I R A C . W I S E W I S E M I P S M I P S P A C S P A C S -0.4 -0.2 0.0 0.2 0.4 P A C S -0.4 -0.2 0.0 0.2 0.4 S P I R E -0.4 -0.2 0.0 0.2 0.4 S P I R E -0.4 -0.2 0.0 0.2 0.4 S P I R E (Data - Model) / Model A r b i t r a r y un i t s Figure 5.
Fractional residuals (Data-Model)/Model for each model in each band. We plot the (Gaussian) kernel densityestimates of the the residual distributions. The WISE 12 band shows narrow, offset fits for all models while other mid-IR bandsshow clear over/under-estimations by some models. The physical models perform a good fit at 250 µ m that gets progressivelyworse towards longer wavelengths. Only the modified blackbody models show systematically good fits, within 10%, in all far-IRbands, likely because of the spectral index, β being a free parameter. since we fixed α = 2 (Aniano et al. 2020). The galaxy-integrated mass-weighted averages are calculated as (cid:104) U (cid:105) = (cid:80) j U j × Σ dust ,j (cid:80) j Σ dust ,j (cid:104) T dust (cid:105) = (cid:80) j T dust ,j × Σ dust ,j (cid:80) j Σ dust ,j . (14)The integrated values are calculated over the pixelsabove the 3 σ detection threshold. In Figure 6, weshow these measurements for each model as diagonalelements. We also provide the ratios between all modelsin M dust and (cid:104) U (cid:105) (or (cid:104) T dust (cid:105) ) to explicitly show their dif- ferences. These off-diagonal elements read as Y-model/ X-model (e.g. M HD21dust /M DL07dust = 0 . . × M (cid:12) (before their renor-malization) by fitting the DL07 at MIPS 160 resolu- enchmarking Dust Emission Models ∼ (cid:48)(cid:48) ) . Using the CIGALE
SED fitting tool(e.g. Boquien et al. 2019, and reference therein) andTHEMIS, Nersesian et al. (2019) found a total dust massof 4 . × M (cid:12) , which is similar to the 5 . × M (cid:12) from THEMIS in this study. It is worth noting thatNersesian et al. (2019) performed fits to the integratedSED, which could lead to a lower dust mass (Anianoet al. 2012; Utomo et al. 2019). However, low signal-to-noise pixels are included in the integrated SED, andexcluded in the resolved fits.It is interesting to note that despite the fairly goodagreement ( ∼ (cid:104) U (cid:105) and (cid:104) T dust (cid:105) . As expected, the HD21 modelrequires the highest radiation field, based on its colderdust. THEMIS and the MC11 model show similar valuesof (cid:104) U (cid:105) . Nersesian et al. (2019) found a dust tempera-ture of 21 . Dust Surface Density, Σ dust Figure 7 shows maps of the dust surface density, Σ dust ,for each model, as well as their ratios with each other.We can see that the physical dust models DL07, HD21,THEMIS and MC11 are all fairly close to each other(light colors), with variations in dust surface densitywithin a factor of 2. They all yield similar dust surfacedensity structures, and appear to vary from one anotherby a spatially smooth offset. THEMIS and the HD21model show the closest Σ dust values, but show an in-version of their ratio around 0.38 r/ r , where THEMISrequires less dust (see Figure 12). The HD21/DL07 and The difference between Aniano et al. (2020)’s dust mass and oursis due to the larger area used in the former for the total dust masscalculation. D L D L H D H D T H E M I S T H E M I S M C M C B E B E S E S E D L H D T H E M I S M C B E S E M dust [10 M ] D L H D T H E M I S M C B E S E D L H D T H E M I S M C B E S E UT [K] Figure 6.
Integrated values (over the 3 σ pixels). The di-agonal elements show the total dust mass (top panel), andradiation field (bottom panel, above the line) or tempera-ture (bottom panel, below the line). The off-diagonal ele-ments are ratios of the integrated values between differentmodels, and are read as “model Y-axis/model X-axis” (e.g. M HD21dust /M DL07dust = 0 . MC11/THEMIS ratio maps are particularly flat, withratios ∼ . . dust values in the center, but drops more6 Chastenet, Sandstrom et al. rapidly than any other model with increasing radius(Figure 12; also Chiang et al. 2018).A caveat of our approach is the difference in treatingthe heating of dust grains between physical dust modelsand modified blackbodies. In the latter, we only use asingle temperature, which is not equivalent to the en-semble of radiation fields used in the physical models.However, given the two extreme behaviors of the modi-fied blackbodies (the simple-emissivity model requiring ahigh dust mass, and the broken-emissivity model requir-ing the lowest), it is not obvious that the use of multipletemperatures (or radiation fields) drive the differences indust mass observed here. Rather, the calibration of themodified blackbodies, and their effective opacity, seemto be more important.5.4.
Average radiation field, U We perform the same ratio analysis with U , derivedfrom the fitted parameters in the physical dust models(Equ. 13 and 14). In Figure 8 we show the radial profilesfor U , as well as the parameter maps and the ratios ofeach model. In the radial profile, the thick lines stopwhere the selection effect due to fitting only bright pixelsbecomes important. Using the SPIRE 500 image, wefound that radial profile of IR emission for all pixelsand for fitted-pixels only differ significantly at ∼ (cid:48) (i.e.0.5 r ). The variations in U are reflective of those of U min , since the γ values are overall small, lending morepower to the delta-function than the power-law (Equ. 7).The overall variations of U appear to be rathersmooth, which is expected since it is dominated by thediffuse radiation field U min . However, we do find en-hanced values of U in H II regions. The ratio mapsdo not strongly show these peaks, which indicates thatall models behave similarly and require higher radia-tion field intensities in these regions. Like for the Σ dust parameter, the ratio maps for U do not display any con-spicuous spatial differences, but rather an offset betweeneach model. This is also visible in the upper panel ofFigure 8.The HD21 model shows the highest values of U . Thisis expected as the dust in the model is “colder” thanother models, due to very few large carbonaceous grains.This leads to a higher radiation field intensity requiredto reach the same luminosity.The spatial distributions of the γ parameter are sim-ilar in all physical models, and we chose not to showthem. Instead, we use the average radiation field, U ,that includes γ in its calculation. The HD21 modelshows the lowest values of γ , which, combined to thehighest values of U min , means it requires more power inthe delta function than the other physical models. 5.5. Fraction of PAHs
Three models have an explicitly defined PAH fraction.The DL07 and HD21 models define the parameter q PAH as the fraction of the dust mass contained in carbona-ceous grains with less than 10 carbon atoms, roughlyless than 1 nm in size. We use the f PAH parameter fromthe MC11 model to estimate a PAH fraction, i.e. massof grains with sizes from 0.35 to 1.2 nm, as defined bythe fiducial parameters. We refer to the mid-IR emis-sion features carriers in THEMIS as HACs. We use thedefinition in Lianou et al. (2019) to compute a fractionof HACs from THEMIS results, that can compare to thePAHs in other dust models: they found that this fractionof HACs corresponds to grains between 0.7 and 1.5 nmof the sCM20 component. It is important to remem-ber that the strict definition of the PAH/HAC fractionis different in each model, but its purpose—fitting themid-IR emission features—remains similar.We investigate the variations of the surface densityof the carriers, Σ
PAH and Σ
HAC , instead of their abun-dances. In the top panel of Figure 9, we show the radialprofiles of Σ { PAH; HAC } . Although the absolute valuesof the surface density of PAHs (HACs) differ by a factorup to ∼ . q PAH fit consistent with 0%. Tovisualize the variations between models, we normalizeeach parameter map to their mean value (as shown inthe color-bar labels). We are thus able to compare thespatial variations of the maps, and avoid the offsets dueto the definition differences of PAHs or HACs.The q PAH map in Aniano et al. (2020, using the DL07model) shows similar features to ours. A large por-tion of the disk of M101 has a rather constant distri-bution of q PAH , with conspicuous drops in H II regions.In their study using the Desert, Boulanger, & Puget(1990) dust model, Rela˜no et al. (2020) found a flat ra-dial profile of the small-to-large grain mass ratio, up to0.8 r ( ∼ . (cid:48) ). Our maps of the fraction of PAHs, orHACs, present a somewhat flat distribution (variationsless than 1%) out to ∼ . ( ∼ . (cid:48) ; see Appendix B),and a steep change further out.5.6. Reproducing the mid-IR emission features enchmarking Dust Emission Models D L H D T H E M I S M C B E DL07 S E HD21 THEMIS MC11 BE SE
Y axisdust / X axisdust log ( dust ) Figure 7.
Dust surface density maps (diagonal) and corresponding ratios with each model. The MC11 model requires thelargest dust mass, followed by the simple-emissivity and DL07 models. All physical dust models show very similar spatialvariations despite having different values of Σ dust . The HD21/DL07 and MC11/THEMIS ratio maps are particularly smoothacross the disk, but the MC11/DL07 and THEMIS/HD21 have the ratios closest to 1. The simple-emissivity model shows clearstructural differences with the other models, by requiring a lot of dust in the center, but rapidly dropping in the outskirts.
To investigate in more detail the ability of each phys-ical model to reproduce the PAH features, we performa fit on an integrated SED and compare the results tomeasurements from the Infrared Spectrograph (IRS; on-board
Spitzer ) in that same region (J. D. T. Smith, pri-vate communication). Figure 10 shows a zoom on themid-IR part of the models and the results of the fits tothe integrated SED. From that fit, it appears that allmodels are able to generally reproduce the mid-IR fea- tures, with their different parameters, but we can noticea few differences between models.From the residuals (bottom panel), we see that themodels perform similarly at 5.8, 8.0 and 12 µ m. At 22and 24 µ m, the offset between the measurements meansthe models tend to split the difference and sit betweenthe points. We note that all models appear to overes-timate the continuum around 7 µ m, in between the 6.2and 7.7 µ m PAH features.8 Chastenet, Sandstrom et al. r / r l o g ( U ) DL07HD21 THEMISMC11 D L H D T H E M I S DL07 M C HD21 THEMIS MC11 U Y axis / U X axis log ( U ) Figure 8.
Top:
Radial profile of (cid:104) U (cid:105) ( mass averaged radi-ation field, see Equ. 14) for each physical model. The thicklines stop where the radial profile is affected by the selectioneffect due to fitting only bright pixels (Section 4.2.2). Bot-tom:
Average radiation field, U , maps (diagonal) and cor-responding ratios with each model. The HD21 model showsthe highest values of U in all the pixels. Despite differentvalues, all models show very similar spatial distribution of U , including in H II regions. r / r { P A H ; H A C } [ M / p c ] DL07HD21 THEMISMC11 D L q PAH = 3.4% H D q PAH = 1.4% T H E M I S q PAH = 3.8%
DL07 M C HD21 THEMIS MC11 q PAH = 2.1% ( q PAH / q PAH ) Y axis / ( q PAH / q PAH ) X axis P q PAH P q PAH [%]
Figure 9.
Top:
Radial profiles for Σ
PAH , or Σ
HAC , inM (cid:12) /pc . The models yield very similar spatial variations.This indicates that the definition of the features carriers mat-ters in terms of their contribution to the total dust mass, butthey all reproduce the mid-IR features similarly. The thicklines stop where the radial profile is affected by the selec-tion effect due to fitting only bright pixels (Section 4.2.2). Bottom:
Fraction of PAHs (or HACs) (diagonal) centeredon their respective mean value, with boundaries at 5 and 95percentiles ( P , P ). The normalized ratios (normalized tothe mean; off-diagonal) show the spatial variations betweenmodels. enchmarking Dust Emission Models µ m than theother three models, despite a similar continuum at20 µ m. This rules out the higher U min found in theHD21 model as the reason for the higher flux at 10 µ m.Rather, in this model, the emission at 10 µ m is stronglydependent on the amount of nano-silicates, used to ac-count for the lack of correlation between PAH emissionand anomalous microwave emission (Hensley & Draine2017). The ratio between the flux at 20 and 10 µ m is ∼ . µ m, while the other mod-els do (although it has no impact on this particular fit).Around 7 µ m, all models show a higher flux than theone seen in the IRS spectrum.It is notable that, using only photometric bands fromWISE and Spitzer /IRAC in the fit, all models reason-ably well reproduce the mid-IR emission features, de-spite having different values of the PAH (or HAC) frac-tion and different definitions of the carriers. However,the comparison to spectroscopic measurements showthat there are still differences between models.5.7.
SPIRE 500 and Σ dust The monochromatic dust emission in far-IR wave-lengths has often been used as a mass tracer of the ISM(e.g. Eales et al. 2012; Berta et al. 2016; Scoville et al.2017; Groves et al. 2015; Aniano et al. 2020; Baes et al.2020). In Figure 11, we plot the emission of M101 at500 µ m as a function of the fitted Σ dust for each model(pixels above the 3 σ detection threshold), color-codedby the minimum radiation field U min , or dust tempera-ture T dust .In all cases, we can see two distinct relations as theSPIRE 500 emission increases. The majority of the fit-ted pixels show a linear scaling between the emission at500 µ m and Σ dust , while in some specific regions of thegalaxy, all models prefer a higher radiation field (or tem-perature) and a lower dust surface density. We providethe scaling relations between the emission at 500 µ m inMJy/sr, and the dust surface density in M (cid:12) /pc , foreach model. We measure the 5 th and 95 th percentilesof the data points above the 3 σ detection threshold (tokeep the bulk of the distribution only). We fit a linear slope to these points :log (Σ d ) DL07 = 1 . ± . × log ( I µmν ) − . ± . (Σ d ) MC11 = 1 . ± . × log ( I µmν ) − . ± . (Σ d ) THEMIS = 1 . ± . × log ( I µmν ) − . ± . (Σ d ) HD21 = 1 . ± . × log ( I µmν ) − . ± . (Σ d ) SE = 1 . ± . × log ( I µmν ) − . ± . (Σ d ) BE = 1 . ± . × log ( I µmν ) − . ± . . (15)To identify the pixels that “branch out” from the bulkwe select any pixel that falls below one standard devia-tion from the fit (dashed-lines). In each panel, we showthe spatial location of these pixels. It becomes clear thatthe regions that need a higher U min (or T dust ) are H II regions and in the outskirts of the galaxy. These pixelscan account between 4 and 11% of the pixels above the3 σ detection. The branching-out from the main rela-tion is likely the consequence of the fact that the dustin these H II regions is significantly hotter than average(as shown by the enhanced U in all models in Figure 8).Lower dust-to-gas ratios in the galaxy outskirts may alsocontribute to this trend.In Aniano et al. (2020), the authors found that thisrelation is well represented by a power-law scaling inthe KINGFISH sample, with slope of 0.942, which islower than our measured value of 1.2 for M101 with theDL07 model. A linear slope would be expected if thedust temperature and optical properties were uniformthroughout the galaxy, leading to a constant I ν / Σ dust ratio. The spatial distribution of temperatures through-out each individual galaxy leads to distinct slopes, andfor M101 we find that in regions of higher dust surfacedensity, the dust is also warmer, leading to more 500 µ memission (this can be seen in change in the color tableon Figure 11 at the highest Σ dust ). The branch of H II region points we see represents only a small fraction ofthe total data and may not be evident on the Anianoet al. (2020) plot. MODEL PERFORMANCE GIVEN ABUNDANCECONSTRAINTS ON DUST MASS The fit coefficients and uncertainties were measured using the numpy.polyfit procedure. Chastenet, Sandstrom et al. I [ M J y / s r ] DL07 - q PAH = 5.6 %HD21 - q PAH = 2.3 % THEMIS - q HAC = 5.4 %MC11 - f PAH = 3.3 % Photometry SEDIRS SpectrumIRS Integrated SED
Wavelengths [ m] -0.250.00.25 R e s i d u a l s Figure 10.
Fits to the integrated SED from the broad-band photometry (open circles) within a rectangle box (drawn inFigure 1). The
Spitzer /IRS spectrum and its corresponding SED (convolved in the 16 bands used here) are shown in gray (lineand filled circles). All models perform a good fit to the mid-IR part of the SED, but show different fractions of PAHs (or HACs).Differences can be noticed between models: a higher 10 µ m emission in the HD21 model, due to nano-silicates, or the lack ofan emission feature at 17 µ m in THEMIS. -2.25-2.0-1.75-1.5-1.25-1.0-0.75-0.5-0.25 l o g ( d [ M / p c ]) : . DL07 : . HD21 -0.2 0.0 0.2 0.4 0.6 0.8 1.0 log (SPIRE 500 [MJy/sr]) -2.25-2.0-1.75-1.5-1.25-1.0-0.75-0.5-0.25 l o g ( d [ M / p c ]) : . THEMIS -0.2 0.0 0.2 0.4 0.6 0.8 1.0 log (SPIRE 500 [MJy/sr]) : . MC11 -2.25-2.0-1.75-1.5-1.25-1.0-0.75-0.5-0.25 l o g ( d [ M / p c ]) : . SE -0.2 0.0 0.2 0.4 0.6 0.8 1.0 log (SPIRE 500 [MJy/sr]) -2.25-2.0-1.75-1.5-1.25-1.0-0.75-0.5-0.25 l o g ( d [ M / p c ]) : . BE l o g ( U m i n ) T d u s t [ K ] Figure 11.
SPIRE 500 emission as a function of the fitted Σ dust for each model, color-coded by the radiation field U min or thedust temperature T dust . We identify a separation in the linear scaling of Σ dust with the dust emission at 500 µ m, marked bythe dashed lines. The pixels below the lines are plotted in color in the maps, and are located in H II regions and surroundings.In these specific locations, the luminosity does not follow the same relation with Σ dust than the rest of the galaxy. enchmarking Dust Emission Models Maximum Dust Surface Density
The calibration of dust models involves a constraint onelements locked in grains (see Section 3). This step relieson depletion measurements, which characterize the dis-tribution of heavy elements between the gas and solidphases. The final amount of elements allowed in dustgrains varies between different physical dust models.The final dust masses derived by each model vary aswell, as discussed in Sections 5.2 and 5.3.A way to assess the performance of dust models is toverify that the required dust mass does not exceed theavailable heavy element mass, as constrained by metal-licity measurements (e.g. Gordon et al. 2014; Chianget al. 2018). We perform this test in M101 since itsmetallicity gradient has been thoroughly characterized(e.g. Zaritsky et al. 1994; Moustakas et al. 2010; Croxallet al. 2016; Berg et al. 2020; Skillman et al. 2020).We estimate the dust mass surface density upper-limitby assuming all available metals are in dust and calculat-ing the metal mass surface density from the metallicitygradient and observed gas mass surface density: M maxdust M gas = M Z M gas , (16)where M gas and M Z are determined as follows.The gas surface density is the sum of H I and H sur-face densities including a correction for the mass of He(we neglect the ionized gas contribution). The latter isbuilt from CO emission, assuming two prescriptions forthe CO-to-H conversion (see Section 2). We include theMW α CO prescription as it is widely used (Equ. 1). Wealso choose the α CO prescription from Bolatto, Wolfire,& Leroy (2013), which takes into account environmentalvariations of α CO with metallicity and surface density(Equ. 2). We emphasize that the Σ dust upper-limit inthis section is dependent on the choice of α CO to derive agas surface density. This is particularly true in the cen-tral region of M101, where H dominates (e.g. Schrubaet al. 2011; V´ılchez et al. 2019) and where the two α CO differ the most. Note also that this result differs withthat of Chiang et al. (2018). This is expected as the α CO conversion factor in their study (from Sandstromet al. 2013) is lower than the ones used in this study,which leads to a lower upper-limit.We use the 12 + log (O / H) radial profile from Berget al. (2020) and convert it to metallicity through: M Z M gas = m O m H (cid:0) (O / H) (cid:1) − . M O M Z , (17) with m O , m H the atomic masses of oxygen and hydro-gen, respectively; and the oxygen-to-metals mass ratio M O /M Z = 0 .
445 (Asplund et al. 2009).The top panels in Figure 12 show the radial profile ofΣ dust for each model and the Σ maxdust upper-limits yieldedby the two α CO prescriptions used here (black lines, dot-ted and dash-dotted). We see on the top left panel thatthe DL07 and MC11 models are above both Σ maxdust upper-limits in almost all the significant pixels. THEMIS andthe HD21 models are fairly in line with the Bolatto,Wolfire, & Leroy (2013) upper-limit, with the conser-vative assumption that 100% metals are in dust grains.These behaviors suggest that the dust emissivity in allphysical models is too low, which leads to requiring toomuch dust. We note that the large dust masses are likelydue to the opacity calibration, rather than a wrong fitin the far-IR bands: in Figure 5, all physical modelsshow a reasonable fit at 160 µ m, much closer to the IRpeak than the 500 µ m band. We are confident that theIR peak is correctly recovered, and that the high dustmasses are not due to the sub-millimeter excess. Basedon the reasonable quality of the fits, we believe that theexcessive mass is likely due to the opacity calibration.6.2. (Re-)Normalization Despite sharing a common calibration approach, thedetails of the opacity calibration in the dust models usedin this study vary in small but significant ways. Whileall models were calibrated to MW diffuse emission, theydid not use exactly the same high-latitude cirrus spec-trum. In addition, the M dust /M H adopted for the MWdiffuse ISM by the physical dust models varies. Addi-tionally, the radiation field that best reproduces the MWdiffuse emission, U MW , differs slightly from one modelto the next. Because of the relationship between dusttemperature and radiation field ( U ∝ T β dust ) and dusttemperature and luminosity, even a slight difference inthe assumed radiation field may lead to a significantchange in the model’s calibrated dust opacity.To investigate calibration discrepancies, we re-normalize each of the dust models via a fit to a commonMW diffuse emission spectrum using the same abun-dance constraints. We use the MW SED described inGordon et al. (2014), which we previously used to cali-brate κ ν of the modified blackbody models (Chiang et al.2018). This SED is the same as that used in Compi`egneet al. (2011), a combination of DIRBE and FIRAS mea-surements (e.g. Boulanger et al. 1996; Arendt et al.2 Chastenet, Sandstrom et al. d u s t [ M / p c ] DL07 maxdust with
BWL13CO
HD21 maxdust with
MWCO
THEMIS maxdust / gas MC11 BE SE
Re-normalized d r / r d u s t / g a s r / r Figure 12.
Top:
Radial profiles for Σ dust (colored lines), and dust surface density upper-limits (black dotted and dash-dotted lines). Upper-limits are estimated from gas and metallicity measurements, assuming all metals are locked in dust grains(Section 6.1). We have projected the dust surface density maps to the gas maps pixel grid, and masked the gas maps wherethere is no dust data (to ensure we are selecting identical pixels in building radial profiles). The upper-limits are invariantbetween the left and right panels.
Left: radial profiles from the fits. All physical models and the simple-emissivity model areeither above or similar to the upper-limits, out to 0.4 r /r . Right: renormalized dust surface densities for the physical models(Section 6.2). The renormalization forces physical models to the same abundance constraints ( M dust /M H = 1 / maxdust lines). Bottom:
Dust-to-gas ratios for each model, assuming the gas surfacedensity derived with α CO from Bolatto, Wolfire, & Leroy (2013). The gray line represents the upper-limit from Berg et al. (2020)and assuming a dust-to-metal ratio of 1 (Equ. 17). The thick lines stop where the radial profile is affected by the selection effectdue to fitting only bright pixels (Section 4.2.2), creating the conspicuous features in the H II region locations. The main bumpat 0.6 r/ r corresponds to the two H II regions NGC 5447 and NGC 5450. The less visible bump at 0.5 r/ r corresponds tothe H II region NGC 5462. . We do not use the ionized gas correction be-cause depletions measurements do not correct for it, and Although more recent measurements from
Planck are availablein the far-IR, we emphasize here that the important aspect isabout uniformity. We choose the DIRBE+FIRAS SED as it isconveniently the one used to calibrate the modified blackbodymodels. Additionally, the significant input brought by
Planck measurements are past the wavelength range used in our study(sub-milimeter and millimeter range). instead use a correction factor of 0.97 for molecular gasonly (Compi`egne et al. 2011). We integrate the SEDin the PACS 100, PACS 160, SPIRE 250, SPIRE 350and SPIRE 500 bands, so all models use the same wave-length coverage. We use a 2.5% uncorrelated and 5%correlated errors to account for FIRAS and DIRBE un-certainties (Gordon et al. 2014). The M dust /M H ratioset in the normalization is 1/150, as suggested by de- enchmarking Dust Emission Models F ∗ = 0 .
36 from Jenkins 2009, see alsoGordon et al. 2014).To perform the re-normalization using the MW SED,we use the same fitting technique as previously describedwith the following choices: 1) We do not use the com-bination of radiation fields, nor the stellar component(i.e. Ω ∗ = γ = 0). 2) We allow the minimum radia-tion field U min and the total dust surface density Σ dust to vary in each physical model. 3) We keep the rela-tive ratios between grain populations fixed and do notvary them independently (e.g. for each model, we usethe total spectra solid lines ‘DL07’, ‘HD21’, ‘MC11’ and‘THEMIS’ in Figure 2).The fits yield renormalization factors that correct allphysical models from their respective assumptions to adust-to-H ratio of 1/150. These corrections range from1.5 for THEMIS, to 3 for the DL07 model (see Ap-pendix A). With this normalization, we are able to meetthe metallicity constraints. The top right panel of Fig-ure 12 shows the radial profiles with the correction fac-tors applied to the surface densities of the physical dustmodels. The renormalization brings the models to lowerdust surface densities that agree with the upper-limitbased on the metal content. It is interesting to notethat the renormalized models now show three distinctbehaviors in the dust mass surface density radial pro-file: DL07, HD21 and the broken power-law emissivitymodified blackbody yield very similar results; THEMISand MC11 are very similar to each other and offset by afactor of ∼ U MWmin to fit the MW SED, and theirinitial spectrum used for calibration.In the bottom panels of Figure 12, we show the radialprofile of the dust-to-gas ratios (DGR) for each model,using the α BWL13CO conversion factor. The bottom leftshows the DGR with the derived dust surface densities,and the bottom right panels is after re-normalization.The thick gray line shows the upper-limit of the DGRusing the metallicity gradient from Berg et al. (2020) andEqu. 17, with α BWL13CO . We find the same abrupt changein the DGR as V´ılchez et al. (2019) around 0.5 r , al-though slightly lower values, consistent with the higherdust surface densities in this study. DISCUSSIONHere we discuss some next steps that should be un-dertaken to improve dust emission modeling. In the pre- vious sections, we investigated some systematic effectslinked to the modeling choices.We showed that physical dust models are likely torequire too much dust mass, exceeding what is avail-able based on metallicity measurements (for reasonablechoices of CO-to-H conversion factors, though we notethat, in the central region, this assertion is strongly de-pendent on this choice). This excess is linked to thecalibration of these models, in particular the elemen-tal abundances prescribed, and their assumed radiationfield. Combined with the growing number of metallic-ity measurements in nearby galaxies (e.g. Kreckel et al.2019; Berg et al. 2020), additional constraints for exter-nal environments (beyond the MW) may help to performbetter fits of dust emission.Several aspects of galaxy evolution studies rely ongrain properties. Dust evolution models heavily relyon observational constraints to find the parameters thatbest match the observed properties of dust. The bal-ance between destruction and formation processes indust evolution models are adjusted by observed dustmasses, that need to be accurately measured (e.g. byemission fitting). Similarly, dust-to-gas ratio evolutionwith metallicity is often derived using dust masses fromemission measurements (e.g. R´emy-Ruyer et al. 2014;De Vis et al. 2019; Nersesian et al. 2019; De Looze et al.2020; Nanni et al. 2020), and are subject to the system-atic biases found in this work.Our study was designed for a rigorous comparison be-tween models fit to mid- through far-IR SEDs. Severalchoices made in this study are justified by our imple-mentation of each model in the most similar frameworkpossible. This also requires using a uniform radiationfield description, and the parameters that go into de-scribing the physical dust models in their fiducial formare not always adapted to the choice of the radiationfield model of this work (Equ. 7).Because of the limited wavelength coverage and SEDsampling of this study, we “tie” together different grainpopulations, based on the similarity of their respectiveemission spectra. In the MC11 the large carbon (LamC)and silicate (aSil) grains have very similar slopes in thefar-IR, which would make these fit parameters stronglydegenerate if both were allowed to vary. In THEMIS, thekey difference between the large carbon grains (lCM20)and the large silicate grains (aSilM5) is their slopes inthe far-IR: the lCM20 grains have a flatter SED thanaSilM5. However, the spectral coverage used in thisstudy is too limited to properly constrain the emissionfrom the two grain populations. These choices have im-plications for the evolution of dust composition in the4 Chastenet, Sandstrom et al.
ISM, since the ratios of carbon-to-silicate in large grainsis assumed to be constant for each model of this study.Additionally, our further tests show that the γ param-eter is degenerate with the emission of some of the grainpopulation in MC11 and THEMIS. The abundance ofsmall amorphous carbon (SamC) in MC11 helps adjustthe slope between 24 and 70 µ m. In the radiation fieldparameterization chosen for this study, the γ parameterhas a similar impact on the shape of the dust emission.Keeping both the small amorphous carbon grains abun-dance and γ as free parameters introduces a degeneracyin the fitting. For this reason we choose to keep the fidu-cial relative abundances of SamC with respect to thatof big grains (LamC+aSil) fixed. In THEMIS, whenallowing both lCM20 and aSilM5 populations to vary(e.g. with f aSil , the fraction of large grains in the formof aSilM5), we also introduce a degeneracy with the γ parameter. A varying ratio of carbon-to-silicate grainsperforms a similar change of the SED shape and bothparameters, γ and f aSil become slightly degenerate. Fu-ture studies with more wavelength coverage and moredetailed constraints on individual elemental abundancesmay be able to allow for more free parameters in thefits. CONCLUSIONSIn this study, we compared the dust properties ofM101 derived from six dust models: four physical dustmodels and two blackbody models. We used the modelsfrom Draine & Li (2007), Compi`egne et al. (2011), Joneset al. (2017, THEMIS) and Hensley & Draine (2021),as well as a simple-emissivity and a broken-emissivitymodified blackbody models to assess the differences invarious dust properties yielded by fitting the mid- tofar-IR emission from
WISE , the
Spitzer Space Telescope and the
Herschel Space Observatory photometry. Ourmain conclusions are: • There are a few notable trends in the fitting residu-als (described as (Data-Model)/Model ; Figure 5).All physical models reproduce the mid-IR bandswithin 10%, with very similar residual distribu-tions in the WISE 12 band. All models performfits of similar quality at 160 µ m. While the mod-ified blackbody models can reproduce the data inall far-IR bands (residuals centered on 0), the fitsfrom physical models have large residuals at longwavelengths. This suggests that the flexibility toadjust the long wavelength slope of the opacity isimportant to reproduce the observed SEDs. • All physical models reproduce the mid-IR emis-sion features but yield different values of the massfraction of their carriers (Figure 9). Models that attribute the mid-IR emission features to PAHs orHACs do similarly well in reproducing the mid-IRspectrum. • We provide scaling relation of Σ dust = f ( I µmν ),and identified a diverging relation in H II regions,where hot dust changes the relationship betweendust emission and mass (Figure 11).Examining the fitting results of total dust masses anddust surface density distributions, we find: • Models yield different total dust masses, up to afactor of 1.4 between physical models, and up to 3including modified blackbodies (Figure 6), but allshow similar spatial distributions of dust surfacedensity (note the fairly low discrepancy betweendust masses from physical models, compared tomodified blackbody models). The MC11 modelrequires the highest dust mass, and the broken-emissivity model the lowest. • We use metallicity and gas measurements to cal-culate a dust surface density upper-limit (assum-ing all metals in dust) and show that all physicaldust models require too much dust over some ra-dial ranges in M101. Only the broken-emissivitymodified blackbody model is below the upper limitof Σ maxdust (Figure 12). This finding is dependent onthe chosen prescription for the CO-to-H conver-sion factor. • To investigate the differences between dust massesand their relationship to the available heavy ele-ments, we renormalized the models via fits to thesame SED of the MW diffuse emission, assuming astrict abundance constraint of M dust /M H = 1 / maxdust (Figure 12). We find that the choices made to cal-ibrate dust models have a non-negligible impacton the derived dust masses.To provide the strictest comparison, we do not alwaysuse dust models in their fiducial aspect, sometimes as-suming a fixed ratio between two dust grain populations.The observational constraints brought by IR emissionfitting are used to validate evolution models or derivescaling relations like the dust-to-gas ratio. Our resultsshow that these derived dust properties have systematicuncertainties that should be taken into account. Al-though there are still systematic uncertainties inherentin H II region metallicity measurements, resolved metal-licity gradients in nearby galaxies can be helpful for test-ing the opacity calibrations in dust models. enchmarking Dust Emission Models Herschel
Space Observatory.
Herschel is an ESA space observa-tory with science instruments provided by European-led Principal Investigator consortia and with impor-tant participation from NASA. This publication makesuse of data products from the Wide-field Infrared Sur- vey Explorer, which is a joint project of the Univer-sity of California, Los Angeles, and the Jet PropulsionLaboratory/California Institute of Technology, fundedby the National Aeronautics and Space Administra-tion. This work is based in part on observations madewith the Spitzer Space Telescope, which was operatedby the Jet Propulsion Laboratory, California Instituteof Technology under a contract with NASA. This re-search made use of matplotlib , a Python library forpublication quality graphics (Hunter 2007). This re-search made use of Astropy, a community-developedcore Python package for Astronomy (Astropy Collab-oration et al. 2018, 2013). This research made useof NumPy (Van Der Walt et al. 2011). This researchmade use of SciPy (Virtanen et al. 2020). This researchmade use of APLpy, an open-source plotting package forPython (Robitaille & Bressert 2012; Robitaille 2019).We acknowledge the usage of the HyperLeda database(http://leda.univ-lyon1.fr).REFERENCES
Allamandola, L. J., Tielens, A. G. G. M., & Barker, J. R.1985, ApJL, 290, L25, doi: 10.1086/184435Aniano, G., Draine, B. T., Gordon, K. D., & Sandstrom, K.2011, PASP, 123, 1218, doi: 10.1086/662219Aniano, G., Draine, B. T., Calzetti, D., et al. 2012, ApJ,756, 138, doi: 10.1088/0004-637X/756/2/138Aniano, G., Draine, B. T., Hunt, L. K., et al. 2020, ApJ,889, 150, doi: 10.3847/1538-4357/ab5fdbArendt, R. G., Odegard, N., Weiland, J. L., et al. 1998,ApJ, 508, 74, doi: 10.1086/306381Asplund, M., Grevesse, N., Sauval, A. J., & Scott, P. 2009,ARA&A, 47, 481,doi: 10.1146/annurev.astro.46.060407.145222Astropy Collaboration, Robitaille, T. P., Tollerud, E. J.,et al. 2013, A&A, 558, A33,doi: 10.1051/0004-6361/201322068Astropy Collaboration, Price-Whelan, A. M., Sip˝ocz, B. M.,et al. 2018, AJ, 156, 123, doi: 10.3847/1538-3881/aabc4fBaes, M., Trˇcka, A., Camps, P., et al. 2020, MNRAS, 494,2912, doi: 10.1093/mnras/staa990Balog, Z., M¨uller, T., Nielbock, M., et al. 2013,Experimental Astronomy,doi: 10.1007/s10686-013-9352-3Bedell, M., Bean, J. L., Mel´endez, J., et al. 2018, ApJ, 865,68, doi: 10.3847/1538-4357/aad908Bendo, G. J., Griffin, M. J., Bock, J. J., et al. 2013,MNRAS, 433, 3062, doi: 10.1093/mnras/stt948 Berg, D. A., Pogge, R. W., Skillman, E. D., et al. 2020,ApJ, 893, 96, doi: 10.3847/1538-4357/ab7eabBerta, S., Lutz, D., Genzel, R., F¨orster-Schreiber, N. M., &Tacconi, L. J. 2016, A&A, 587, A73,doi: 10.1051/0004-6361/201527746Bianchi, S. 2013, A&A, 552, A89,doi: 10.1051/0004-6361/201220866Bolatto, A. D., Wolfire, M., & Leroy, A. K. 2013, ARA&A,51, 207, doi: 10.1146/annurev-astro-082812-140944Boquien, M., Burgarella, D., Roehlly, Y., et al. 2019, A&A,622, A103, doi: 10.1051/0004-6361/201834156Boulanger, F., Abergel, A., Bernard, J. P., et al. 1996,A&A, 312, 256Bron, E., Le Bourlot, J., & Le Petit, F. 2014, A&A, 569,A100, doi: 10.1051/0004-6361/201322101Byrne, L., Christensen, C., Tsekitsidis, M., Brooks, A., &Quinn, T. 2019, ApJ, 871, 213,doi: 10.3847/1538-4357/aaf9aaCastellanos, P., Candian, A., Andrews, H., & Tielens,A. G. G. M. 2018, A&A, 616, A167,doi: 10.1051/0004-6361/201833221Chastenet, J., Bot, C., Gordon, K. D., et al. 2017, A&A,601, A55, doi: 10.1051/0004-6361/201629133Chastenet, J., Sandstrom, K., Chiang, I. D., et al. 2019,ApJ, 876, 62, doi: 10.3847/1538-4357/ab16cfChiang, I. D., Sandstrom, K. M., Chastenet, J., et al. 2018,ApJ, 865, 117, doi: 10.3847/1538-4357/aadc5f Chastenet, Sandstrom et al.
Chiang, I.-D., Sandstrom, K. M., Chastenet, J., et al. 2020,arXiv e-prints, arXiv:2011.10561.https://arxiv.org/abs/2011.10561Compi`egne, M., Verstraete, L., Jones, A., et al. 2011, A&A,525, A103, doi: 10.1051/0004-6361/201015292Croxall, K. V., Pogge, R. W., Berg, D. A., Skillman, E. D.,& Moustakas, J. 2016, ApJ, 830, 4,doi: 10.3847/0004-637X/830/1/4Croxall, K. V., Smith, J. D., Wolfire, M. G., et al. 2012,ApJ, 747, 81, doi: 10.1088/0004-637X/747/1/81Dale, D. A., Helou, G., Contursi, A., Silbermann, N. A., &Kolhatkar, S. 2001, ApJ, 549, 215, doi: 10.1086/319077Dale, D. A., Cohen, S. A., Johnson, L. C., et al. 2009, ApJ,703, 517, doi: 10.1088/0004-637X/703/1/517Davies, J. I., Baes, M., Bianchi, S., et al. 2017, PASP, 129,044102, doi: 10.1088/1538-3873/129/974/044102de Blok, W. J. G., Walter, F., Brinks, E., et al. 2008, AJ,136, 2648, doi: 10.1088/0004-6256/136/6/2648De Looze, I., Lamperti, I., Saintonge, A., et al. 2020,MNRAS, doi: 10.1093/mnras/staa1496De Vis, P., Jones, A., Viaene, S., et al. 2019, A&A, 623,A5, doi: 10.1051/0004-6361/201834444Demyk, K., Meny, C., Lu, X. H., et al. 2017, A&A, 600,A123, doi: 10.1051/0004-6361/201629711Desert, F. X., Boulanger, F., & Puget, J. L. 1990, A&A,500, 313Desert, F. X., Boulanger, F., & Shore, S. N. 1986, A&A,160, 295Draine, B. T. 2003, ApJ, 598, 1026, doi: 10.1086/379123—. 2011, Physics of the Interstellar and IntergalacticMedium—. 2016, ApJ, 831, 109, doi: 10.3847/0004-637X/831/1/109Draine, B. T., & Hensley, B. S. 2020, arXiv e-prints,arXiv:2009.11314. https://arxiv.org/abs/2009.11314Draine, B. T., & Lee, H. M. 1984, ApJ, 285, 89,doi: 10.1086/162480Draine, B. T., & Li, A. 2001, ApJ, 551, 807,doi: 10.1086/320227—. 2007, ApJ, 657, 810, doi: 10.1086/511055Draine, B. T., Aniano, G., Krause, O., et al. 2014, ApJ,780, 172, doi: 10.1088/0004-637X/780/2/172Eales, S., Smith, M. W. L., Auld, R., et al. 2012, ApJ, 761,168, doi: 10.1088/0004-637X/761/2/168Engelbracht, C. W., Blaylock, M., Su, K. Y. L., et al. 2007,PASP, 119, 994, doi: 10.1086/521881Fanciullo, L., Kemper, F., Scicluna, P., Dharmawardena,T. E., & Srinivasan, S. 2020, MNRAS, 499, 4666,doi: 10.1093/mnras/staa2911Fazio, G. G., Hora, J. L., Allen, L. E., et al. 2004, ApJS,154, 10, doi: 10.1086/422843 Finkbeiner, D. P., Davis, M., & Schlegel, D. J. 1999, ApJ,524, 867, doi: 10.1086/307852Fitzpatrick, E. L. 1999, PASP, 111, 63, doi: 10.1086/316293Fitzpatrick, E. L., Massa, D., Gordon, K. D., Bohlin, R., &Clayton, G. C. 2019, ApJ, 886, 108,doi: 10.3847/1538-4357/ab4c3aFumagalli, M., Krumholz, M. R., & Hunt, L. K. 2010, ApJ,722, 919, doi: 10.1088/0004-637X/722/1/919Galametz, M., Albrecht, M., Kennicutt, R., et al. 2014,MNRAS, 439, 2542, doi: 10.1093/mnras/stu113Galliano, F., Galametz, M., & Jones, A. P. 2018, ARA&A,56, 673, doi: 10.1146/annurev-astro-081817-051900Galliano, F., Hony, S., Bernard, J. P., et al. 2011, A&A,536, A88, doi: 10.1051/0004-6361/201117952Gordon, K. D., Engelbracht, C. W., Fadda, D., et al. 2007,PASP, 119, 1019, doi: 10.1086/522675Gordon, K. D., Roman-Duval, J., Bot, C., et al. 2014, ApJ,797, 85, doi: 10.1088/0004-637X/797/2/85Grevesse, N., & Sauval, A. J. 1998, SSRv, 85, 161,doi: 10.1023/A:1005161325181Griffin, M. J., Abergel, A., Abreu, A., et al. 2010, A&A,518, L3, doi: 10.1051/0004-6361/201014519Groves, B. A., Schinnerer, E., Leroy, A., et al. 2015, ApJ,799, 96, doi: 10.1088/0004-637X/799/1/96Guillet, V., Fanciullo, L., Verstraete, L., et al. 2018, A&A,610, A16, doi: 10.1051/0004-6361/201630271Harper, D. A., Runyan, M. C., Dowell, C. D., et al. 2018,Journal of Astronomical Instrumentation, 7, 1840008,doi: 10.1142/S2251171718400081Hensley, B. S., & Draine, B. T. 2017, ApJ, 836, 179,doi: 10.3847/1538-4357/aa5c37—. 2020a, arXiv e-prints, arXiv:2009.00018.https://arxiv.org/abs/2009.00018—. 2020b, ApJ, 895, 38, doi: 10.3847/1538-4357/ab8cc3—. 2021, arXiv e-prints. https://arxiv.org/abs/200x.xxxxxHunter, J. D. 2007, Computing In Science & Engineering,9, 90Jenkins, E. B. 2009, ApJ, 700, 1299,doi: 10.1088/0004-637X/700/2/1299Jones, A. P., Fanciullo, L., K¨ohler, M., et al. 2013, A&A,558, A62, doi: 10.1051/0004-6361/201321686Jones, A. P., K¨ohler, M., Ysard, N., Bocchio, M., &Verstraete, L. 2017, A&A, 602, A46,doi: 10.1051/0004-6361/201630225Kemper, F., Vriend, W. J., & Tielens, A. G. G. M. 2004,ApJ, 609, 826, doi: 10.1086/421339Kennicutt, R. C., Calzetti, D., Aniano, G., et al. 2011,PASP, 123, 1347, doi: 10.1086/663818K¨ohler, M., Jones, A., & Ysard, N. 2014, A&A, 565, L9,doi: 10.1051/0004-6361/201423985 enchmarking Dust Emission Models Kreckel, K., Ho, I. T., Blanc, G. A., et al. 2019, ApJ, 887,80, doi: 10.3847/1538-4357/ab5115Leger, A., & Puget, J. L. 1984, A&A, 500, 279Lenz, D., Hensley, B. S., & Dor´e, O. 2017, ApJ, 846, 38,doi: 10.3847/1538-4357/aa84afLeroy, A. K., Walter, F., Brinks, E., et al. 2008, AJ, 136,2782, doi: 10.1088/0004-6256/136/6/2782Leroy, A. K., Walter, F., Bigiel, F., et al. 2009, AJ, 137,4670, doi: 10.1088/0004-6256/137/6/4670Leroy, A. K., Sandstrom, K. M., Lang, D., et al. 2019,ApJS, 244, 24, doi: 10.3847/1538-4365/ab3925Li, A., & Draine, B. T. 2001a, ApJ, 554, 778,doi: 10.1086/323147—. 2001b, ApJ, 554, 778, doi: 10.1086/323147Lianou, S., Barmby, P., Mosenkov, A. A., Lehnert, M., &Karczewski, O. 2019, A&A, 631, A38,doi: 10.1051/0004-6361/201834553Makarov, D., Prugniel, P., Terekhova, N., Courtois, H., &Vauglin, I. 2014, A&A, 570, A13,doi: 10.1051/0004-6361/201423496Mathis, J. S. 1990, ARA&A, 28, 37,doi: 10.1146/annurev.aa.28.090190.000345Mathis, J. S., Mezger, P. G., & Panagia, N. 1983, A&A,500, 259Mathis, J. S., Rumpl, W., & Nordsieck, K. H. 1977, ApJ,217, 425, doi: 10.1086/155591Moustakas, J., Kennicutt, Robert C., J., Tremonti, C. A.,et al. 2010, ApJS, 190, 233,doi: 10.1088/0067-0049/190/2/233M¨uller, T., Nielbock, M., Balog, Z., Klaas, U., & Vilenius,E. 2011, PACS Photometer - Point-Source FluxCalibration, Tech. Rep. PICC-ME-TN-037, HerschelMutschke, H., & Mohr, P. 2019, A&A, 625, A61,doi: 10.1051/0004-6361/201834805Nanni, A., Burgarella, D., Theul´e, P., Cˆot´e, B., &Hirashita, H. 2020, arXiv e-prints, arXiv:2006.15146.https://arxiv.org/abs/2006.15146Nersesian, A., Xilouris, E. M., Bianchi, S., et al. 2019,A&A, 624, A80, doi: 10.1051/0004-6361/201935118Onaka, T., Yamamura, I., Tanabe, T., Roellig, T. L., &Yuen, L. 1996, PASJ, 48, L59, doi: 10.1093/pasj/48.5.L59Ott, S. 2010, in Astronomical Society of the PacificConference Series, Vol. 434, Astronomical Data AnalysisSoftware and Systems XIX, ed. Y. Mizumoto, K.-I.Morita, & M. Ohishi, 139.https://arxiv.org/abs/1011.1209Paradis, D., M´eny, C., Juvela, M., Noriega-Crespo, A., &Ristorcelli, I. 2019, A&A, 627, A15,doi: 10.1051/0004-6361/201935158 Peimbert, A., & Peimbert, M. 2010, ApJ, 724, 791,doi: 10.1088/0004-637X/724/1/791Pilbratt, G. L., Riedinger, J. R., Passvogel, T., et al. 2010,A&A, 518, L1, doi: 10.1051/0004-6361/201014759Planck Collaboration, Abergel, A., Ade, P. A. R., et al.2011, A&A, 536, A24, doi: 10.1051/0004-6361/201116485Planck Collaboration Int. XVII. 2014, A&A, 566, A55,doi: 10.1051/0004-6361/201323270Planck Collaboration Int. XXII. 2015, A&A, 576, A107,doi: 10.1051/0004-6361/201424088Planck Collaboration XI. 2018, arXiv e-prints,arXiv:1801.04945. https://arxiv.org/abs/1801.04945Poglitsch, A., Waelkens, C., Geis, N., et al. 2010, A&A,518, L2, doi: 10.1051/0004-6361/201014535Reach, W. T., Megeath, S. T., Cohen, M., et al. 2005,PASP, 117, 978, doi: 10.1086/432670Rela˜no, M., Lisenfeld, U., Hou, K. C., et al. 2020, A&A,636, A18, doi: 10.1051/0004-6361/201937087R´emy-Ruyer, A., Madden, S. C., Galliano, F., et al. 2014,A&A, 563, A31, doi: 10.1051/0004-6361/201322803Richey, C. R., Kinzer, R. E., Cataldo, G., et al. 2013, ApJ,770, 46, doi: 10.1088/0004-637X/770/1/46Rieke, G. H., Young, E. T., Engelbracht, C. W., et al. 2004,ApJS, 154, 25, doi: 10.1086/422717Robitaille, T. 2019, APLpy v2.0: The AstronomicalPlotting Library in Python, doi: 10.5281/zenodo.2567476Robitaille, T., & Bressert, E. 2012, APLpy: AstronomicalPlotting Library in Python, Astrophysics Source CodeLibrary. http://ascl.net/1208.017Roussel, H. 2013, PASP, 125, 1126, doi: 10.1086/673310Sandstrom, K. M., Bolatto, A. D., Draine, B. T., Bot, C.,& Stanimirovi´c, S. 2010, ApJ, 715, 701,doi: 10.1088/0004-637X/715/2/701Sandstrom, K. M., Leroy, A. K., Walter, F., et al. 2013,ApJ, 777, 5, doi: 10.1088/0004-637X/777/1/5Schlafly, E. F., Meisner, A. M., Stutz, A. M., et al. 2016,ApJ, 821, 78, doi: 10.3847/0004-637X/821/2/78Schruba, A., Leroy, A. K., Walter, F., et al. 2011, AJ, 142,37, doi: 10.1088/0004-6256/142/2/37Scott, P., Asplund, M., Grevesse, N., Bergemann, M., &Sauval, A. J. 2015a, A&A, 573, A26,doi: 10.1051/0004-6361/201424110Scott, P., Grevesse, N., Asplund, M., et al. 2015b, A&A,573, A25, doi: 10.1051/0004-6361/201424109Scoville, N., Lee, N., Vanden Bout, P., et al. 2017, ApJ,837, 150, doi: 10.3847/1538-4357/aa61a0Skillman, E. D., Berg, D. A., Pogge, R. W., et al. 2020,ApJ, 894, 138, doi: 10.3847/1538-4357/ab86aeSofue, Y., Tutui, Y., Honma, M., et al. 1999, ApJ, 523, 136,doi: 10.1086/307731 Chastenet, Sandstrom et al.
Tanaka, M., Matsumoto, T., Murakami, H., et al. 1996,PASJ, 48, L53, doi: 10.1093/pasj/48.5.L53Tchernyshyov, K., Meixner, M., Seale, J., et al. 2015, ApJ,811, 78, doi: 10.1088/0004-637X/811/2/78Thi, W. F., Hocuk, S., Kamp, I., et al. 2020, A&A, 634,A42, doi: 10.1051/0004-6361/201731746Tully, R. B., Rizzi, L., Shaya, E. J., et al. 2009, AJ, 138,323, doi: 10.1088/0004-6256/138/2/323Utomo, D., Chiang, I. D., Leroy, A. K., Sand strom, K. M.,& Chastenet, J. 2019, ApJ, 874, 141,doi: 10.3847/1538-4357/ab05d3Van Der Walt, S., Colbert, S. C., & Varoquaux, G. 2011,Computing in Science & Engineering, 13, 22V´ılchez, J. M., Rela˜no, M., Kennicutt, R., et al. 2019,MNRAS, 483, 4968, doi: 10.1093/mnras/sty3455Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020,Nature Methods, 17, 261,doi: https://doi.org/10.1038/s41592-019-0686-2Wakelam, V., Bron, E., Cazaux, S., et al. 2017, MolecularAstrophysics, 9, 1, doi: 10.1016/j.molap.2017.11.001Walter, F., Brinks, E., de Blok, W. J. G., et al. 2008, AJ,136, 2563, doi: 10.1088/0004-6256/136/6/2563 Weingartner, J. C., & Draine, B. T. 2001a, ApJS, 134, 263,doi: 10.1086/320852—. 2001b, ApJ, 548, 296, doi: 10.1086/318651—. 2001c, ApJ, 548, 296, doi: 10.1086/318651Werner, M. W., Roellig, T. L., Low, F. J., et al. 2004,ApJS, 154, 1, doi: 10.1086/422992Wolfire, M. G., Hollenbach, D., McKee, C. F., Tielens,A. G. G. M., & Bakes, E. L. O. 1995, ApJ, 443, 152,doi: 10.1086/175510Wright, E. L., Eisenhardt, P. R. M., Mainzer, A. K., et al.2010, AJ, 140, 1868, doi: 10.1088/0004-6256/140/6/1868Ysard, N., K¨ohler, M., Jones, A., et al. 2015, A&A, 577,A110, doi: 10.1051/0004-6361/201425523Zaritsky, D., Kennicutt, Robert C., J., & Huchra, J. P.1994, ApJ, 420, 87, doi: 10.1086/173544Zeegers, S. T., Costantini, E., Rogantini, D., et al. 2019,A&A, 627, A16, doi: 10.1051/0004-6361/201935050Zubko, V., Dwek, E., & Arendt, R. G. 2004, ApJS, 152,211, doi: 10.1086/382351 enchmarking Dust Emission Models A. CALIBRATION DETAILSWe present here the details of the calibration methodology used in each physical model, and a summary of thecalibration constraints. A.1.
Draine & Li (2007)
DL07 was calibrated using the following constraints. The extinction is described in Weingartner & Draine (2001c)and uses the Fitzpatrick (1999) extinction curves with a normalization of N H /E ( B − V ) = 5 . × H/cm or A V /N H = 5 . × − cm . The high-latitude cirrus emission per H observed by DIRBE (Diffuse Infrared BackgroundExperiment) and FIRAS (Far Infrared Absolute Spectrophotometer; Arendt et al. 1998; Finkbeiner et al. 1999) is usedas a reference for the far-IR emission, complemented by mid- and near-IR emission from IRTS (Infrared Telescope inSpace; Onaka et al. 1996; Tanaka et al. 1996). Weingartner & Draine (2001b) adopt solar abundances from Grevesse& Sauval (1998), assuming 30% of carbon is in the gas phase. They assume all silicon is depleted and has abundanceequal to the Solar value. DL07 uses M dust /M H = 1 . × − . The radiation field used in the Draine & Li (2007) modelis based on Mathis et al. (1983). A.2. Compi`egne et al. (2011)
MC11 was calibrated using the following constraints. Extinction constraints were taken from Mathis (1990) andFitzpatrick (1999), including the R V = 3 . N H /E ( B − V ) =5 . × H/cm . At λ > µ m, MC11 use the MW cirrus emission per H observed by COBE-DIRBE and WMAP(integrated in the Herschel and Planck/HFI bands; see MC11). At λ ≤ µ m, a compilation of mid-IR observations ofhigh latitude MW cirrus are used (combining measurements from AROME, DIRBE, and ISOCAM; we refer to readerto the Compi`egne et al. paper for details). They scale the emission SED by 0.77 to account for ionized and moleculargas not accounted in the H column. The allowed dust-phase abundances for C, O, and other dust components comefrom the difference between Solar (or F/G star) abundances and the observed gas phase abundances. In total, the M dust / M H = 1 . × − . MC11 assumes the Mathis et al. (1983, D G = 10 kpc) solar neighborhood radiation field toheat the dust grains. A.3. THEMIS
THEMIS was calibrated using the same constraints as Compi`egne et al. (2011) presented in the previous Section,with the addition of the far-IR-to-extinction relation τ /E ( B − V ) = 5 . × − (Planck Collaboration et al. 2011).In THEMIS, M dust /M H = 7 . × − . A.4. Hensley & Draine (2021)
The full set of observational constraints used to develop the model are described in Hensley & Draine (2020a).In brief, the extinction curve is primarily a synthesis of those of Fitzpatrick et al. (2019) in the UV and optical,Schlafly et al. (2016) in the optical and near-infrared, and Hensley & Draine (2020b) in the mid-infrared. Thenormalization N H /E ( B − V ) = 8 . × cm − mag − is used to normalize extinction to the hydrogen column (Lenzet al. 2017). The infrared emission in both total intensity and polarization are based on the analyses presented inPlanck Collaboration Int. XVII (2014), Planck Collaboration Int. XXII (2015), and Planck Collaboration XI (2018),including the normalization to N H . The solid phase interstellar abundances are re-determined in Hensley & Draine(2020a) using a set of Solar abundances (Asplund et al. 2009; Scott et al. 2015a,b), a measurement of Galactic chemicalenrichment from Solar twin studies (Bedell et al. 2018), and determination of the gas phase abundances from absorptionspectroscopy (Jenkins 2009). M dust /M H = 1 . × − in this model. HD21 assumes the same radiation field as theDL07 model to heat the dust grains, with updates from Draine (2011).0 Chastenet, Sandstrom et al. T a b l e . P h y s i c a l M o d e l s C a li b r a t i o nSu mm a r y D r a i n e & L i ( ) C o m p i ` e g n ee t a l. ( ) T H E M I S H e n s l e y & D r a i n e ( ) E x t i n c t i o n c u r v e F i t z p a t r i c k ( ) M a t h i s ( ) M a t h i s ( ) F i t z p a t r i c k e t a l. ( ) S c h l a fl y e t a l. ( ) H e n s l e y & D r a i n e ( b ) N H / E ( B − V ) . × H / c m . × H / c m . × H / c m a . × c m − m ag − E m i ss i o n s p e c t r u m O n a k a e t a l. ( ) c o m p il e d i n C o m p i ` e g n ee t a l. ( ) c o m p il e d i n C o m p i ` e g n ee t a l. ( ) P l a n c k C o ll a b o r a t i o n I n t . XV II ( ) T a n a k a e t a l. ( ) P l a n c k C o ll a b o r a t i o n I n t . XX II ( ) A r e nd t e t a l. ( ) P l a n c k C o ll a b o r a t i o n X I ( ) F i n k b e i n e r e t a l. ( ) M d / M H . × − . × − . × − . × − R a d i a t i o n F i e l d M a t h i s e t a l. ( ) M a t h i s e t a l. ( ) M a t h i s e t a l. ( ) D r a i n e ( ) R e n o r m a li z a t i o n : e m i ss i o n c o n s t r a i n t o n l y f r o m C o m p i ` e g n ee t a l. ( ) b a nd f o r c i n g M d / M H = . × − U M W m i n . . . . N o r m a li z a t i o n f a c t o r . . . . a a dd i t i o n a l c o n s t r a i n t : τ / E ( B − V ) = . × − ( P l a n c k C o ll a b o r a t i o n e t a l. ) b C o rr e c t e d f o r m o l e c u l a r ga s o n l y , n o t i o n i z e d ga s . N o t e — A ll m o d e l ss h a r e R V = . . enchmarking Dust Emission Models B. FITTED PARAMETER MAPSWe present the spatial variations of the fitted param-eters for all six models, in Figures 13 – 17. The contoursmark the 3 σ detection threshold. We use the same scalefor identical parameters, when possible. C. RESIDUAL MAPSWe present the spatial variations of the fractionalresiduals (Data-Model)/Model for all six models, in Fig-ures 18 – 22. The contours mark the 3 σ detection thresh-old. The so-called “sub-millimeter” excess is visible inmost maps at SPIRE 500 (blue shade). D. FITS QUALITYWe show the relative quality of the fits between eachmodel. The value displayed in the maximum likelihood,in arbitrary units. The simple-emissivity and broken-emissivity models show the least dynamic range butnever reach the highest values of the physical models.For the physical models, we can clearly see the H II re-gions showing fits with low confidence, likely related tothe issues mentioned in Section 5.7.2 Chastenet, Sandstrom et al.
15 20 25 30 35 T dust [K] log ( d [M pc ]) T dust [K] log ( d [M pc ]) Figure 13.
Maps of fitted parameters.
Top: simple-emissivity model, dust temperature ( T dust ), total dust surface density(Σ dust ) and spectral index ( β ). Bottom: broken-emissivity model, dust temperature ( T dust ), total dust surface density (Σ dust )and second spectral index ( β ); the breaking wavelength is fixed ( λ c = 300 µ m) as well as the first spectral index ( β = 2). enchmarking Dust Emission Models log ( U min ) log ( d [M pc ]) log ( ) q PAH [%] log ( * ) Figure 14.
Maps of the fitted parameters for the Draine & Li (2007) model: minimum radiation field ( U min ), total dust surfacedensity (Σ dust ), fraction of dust mass heated by a power-law distribution of radiation field ( γ ), PAH fraction (mass in grainswith less than 10 C atoms, q PAH ), and scaling parameter of surface brightness (5,000 K blackbody, Ω ∗ ). Chastenet, Sandstrom et al. log ( U min ) log ( d [M pc ]) log ( ) f PAH log ( * ) Figure 15.
Maps of the fitted parameters for the Compi`egne et al. (2011) model: minimum radiation field ( U min ), total dustsurface density (Σ dust ), fraction of dust mass heated by a power-law distribution of radiation field ( γ ), PAH fraction (withrespect to total dust mass, f PAH ), and scaling parameter of surface brightness (5,000 K blackbody, Ω ∗ ). enchmarking Dust Emission Models log ( U min ) log ( d [M pc ]) log ( ) f sCM20 log ( * ) Figure 16.
Maps of the fitted parameters for the Jones et al. (2013) model: minimum radiation field ( U min ), total dust surfacedensity (Σ dust ), fraction of dust mass heated by a power-law distribution of radiation field ( γ ), sCM20 fraction (small carbongrains, f sCM20 ), and scaling parameter of surface brightness (5,000 K blackbody, Ω ∗ ). Chastenet, Sandstrom et al. log ( U min ) log ( d [M pc ]) log ( ) q PAH [%] log ( * ) Figure 17.
Maps of the fitted parameters for the Hensley & Draine (2021) model: minimum radiation field ( U min ), total dustsurface density (Σ dust ), fraction of dust mass heated by a power-law distribution of radiation field ( γ ), PAH fraction (mass ingrains with less than 10 C atoms, q PAH ), and scaling parameter of surface brightness (5,000 K blackbody, Ω ∗ ). enchmarking Dust Emission Models P A C S P A C S S P I R E S P I R E S P I R E P A C S P A C S S P I R E S P I R E S P I R E . . . . . (Dat-Model)/Model - Blackbody models S i n g l e - e m i ss i v i t y B r o k e n - e m i ss i v i t y Figure 18.
Maps of the fractional residuals for the simple-emissivity and broken-emissivity models. Chastenet, Sandstrom et al. W I S E . I R A C . I R A C . W I S E . I R A C . I R A C . W I S E W I S E M I P S M I P S P A C S P A C S P A C S S P I R E S P I R E S P I R E . . . . . . . . . (Dat-Model)/Model - DL07 Figure 19.
Maps of the fractional residuals for the Draine & Li (2007) model. enchmarking Dust Emission Models W I S E . I R A C . I R A C . W I S E . I R A C . I R A C . W I S E W I S E M I P S M I P S P A C S P A C S P A C S S P I R E S P I R E S P I R E . . . . . . . . . (Dat-Model)/Model - MC11 Figure 20.
Maps of the fractional residuals for the Compi`egne et al. (2011) model. Chastenet, Sandstrom et al. W I S E . I R A C . I R A C . W I S E . I R A C . I R A C . W I S E W I S E M I P S M I P S P A C S P A C S P A C S S P I R E S P I R E S P I R E . . . . . . . . . (Dat-Model)/Model - THEMIS Figure 21.
Maps of the fractional residuals for the Jones et al. (2013) model. enchmarking Dust Emission Models W I S E . I R A C . I R A C . W I S E . I R A C . I R A C . W I S E W I S E M I P S M I P S P A C S P A C S P A C S S P I R E S P I R E S P I R E . . . . . . . . . (Dat-Model)/Model - HD21 Figure 22.
Maps of the fractional residuals for the Hensley & Draine (2021) model. Chastenet, Sandstrom et al.
SE DL07 THEMISBE MC11 HD21 R e d u c e d Figure 23.
Maps of the reduced χ in each pixel. The H II regions show the lowest confidence for the physical models, whilethe modified blackbodies show good fits on the entire disk. The contour marks the 3- σσ