Towards understanding the Planck thermal dust models

TL;DR

This paper compares three Planck thermal dust emission models to identify key differences in their parameters, which is crucial for improving cosmic microwave background measurements, especially for detecting primordial B-modes.

Contribution

It provides a detailed comparison of three Planck dust models, highlighting differences caused by priors and background contributions, aiding future CMB analysis.

Findings

Significant differences in spectral index, optical depth, and dust temperature among models.

Differences stem from priors and cosmic infrared background contributions.

Insights help refine dust modeling for CMB experiments.

Abstract

Understanding the properties of dust emission in the microwave domain is an important premise for the next generation of cosmic microwave background (CMB) experiments, devoted to the measurement of the primordial $B$-modes of polarization. In this paper, we compare three solutions to thermal dust emission by the Planck Collaboration \cite{PlanckDust03,planck_16,planck_com} to point out significant differences between their respective parameters (the spectral index $\beta$, the optical depth $\tau$ and the dust temperature $T_d$). These differences originate from e.g. the priors on the parameters or the contribution of the Cosmic infrared background (CIB).

Figures

Tables

Table 1. Table 1: Cross-correlation coefficients between T d subscript 𝑇 𝑑 T_{d} and β 𝛽 \beta of the P13, P16 and C15 solutions for the full sky.

	$T_{d,P13}$	$T_{d,P15}$	$T_{d,P16}$	${\beta }_{P13}$	${\beta }_{C15}$	${\beta }_{P16}$
$T_{d,P13}$	1.00	0.71	0.67	-0.60	-0.16	-0.06
$T_{d,C15}$		1.00	0.19	-0.11	-0.04	0.52
$T_{d,P16}$			1.00	-0.30	-0.22	-0.51
${\beta }_{P13}$				1.00	0.27	0.51
${\beta }_{C15}$					1.00	0.21
${\beta }_{P16}$						1.00

Table 2. Table 2: Cross-correlation coefficients between T d subscript 𝑇 𝑑 T_{d} and τ 353 subscript 𝜏 353 \tau_{353} for the P13, P16 and C15 solutions. Note that here we use a ring mask to exclude ± 10 ∘ plus-or-minus superscript 10 \pm 10^{\circ} around the Galactic plane.

	$T_{d,P13}$	$T_{d,C15}$	$T_{d,P16}$	${\tau }_{353,P13}$	${\tau }_{353,C15}$	${\tau }_{353,P16}$
$T_{d,P13}$	1.00	0.64	0.66	-0.57	-0.55	-0.59
$T_{d,C15}$		1.00	0.07	-0.62	-0.63	-0.64
$T_{d,P16}$			1.00	-0.20	-0.18	-0.21
${\tau }_{353,P13}$				1.00	1.00	1.00
${\tau }_{353,C15}$					1.00	0.99
${\tau }_{353,P16}$						1.00

Table 3. Table 3: List of the monopole correction for the 357-3000 GHz maps by P13 and P16, in MJy/sr.

$\nu$ [GHz]	353	545	857	3000
P13	0.085	0.095	0.093	-0.174
P16	0.125	0.336	0.556	0.113

Equations

$$I_{\nu }(\mathbf{n})={\tau }_0(\mathbf{n}){\left(\frac{\nu }{{\nu }_0}\right)}^{\beta (\mathbf{n})}{B}_{\nu }(T_d(\mathbf{n}))$$

$$B_{\nu }(T_d(\mathbf{n}))=\frac{2h{\nu }^3}{c^2}{\left(e^{\frac{h\nu }{k{T}_d(\mathbf{n})}}-1\right)}^{-1}$$

$$\displaystyle\left\{\begin{array}[]{ccc}\tau(\textbf{n})\\ \beta(\textbf{n})\\ T_{d}(\textbf{n})\end{array}\right\}=\sum_{l=0}^{\infty}\sum_{m=-l}^{l}\left\{\begin{array}[]{ccc}\tau_{lm}\\ \beta_{lm}\\ T_{d,lm}\end{array}\right\}Y_{lm}(\textbf{n})$$

$$D^X(l)=\frac{l(l+1)}{2\pi (2l+1)}\sum _{-m}^m|X_{lm}{|}^2$$

$$I_{\nu }(\text{n})\simeq \overline{\tau }{\left(\frac{\nu }{{\nu }_o}\right)}^{\overline{\beta }}{B}_{\nu }\left({\overline{T}}_d\right)\times$$

$$\times \left\{1+\frac{\Delta \tau (\text{n})}{\overline{\tau }}+\Delta \beta (\text{n})\cdot \ln \left(\frac{\nu }{{\nu }_0}\right)+{\frac{d\ln B_{\nu }}{d\ln T_d}|}_{{\overline{T}}_d}\cdot \frac{\Delta T_d(\text{n})}{{\overline{T}}_d}\right\}.$$

$$I_{\nu }(\text{n})\propto \tau (\text{n})T_d(\text{n}){\left(\frac{\nu }{{\nu }_0}\right)}^{2+\beta (\text{n})}$$

Full text

I Introduction

The next generation of the Cosmic Microwave Background (CMB) experiments, such as LiteBird Hazumi et al. (2012); Matsumura et al. (2014, 2016), the CMB-S4 missions Abazajian et al. (2016), CORE (CORE Collaboration et al., 2016; Di Valentino et al., 2016), DeepSpace (The DeepSpace Collaboration, 2016), PIXIE (Kogut et al., 2011), and Polarbear The Polarbear Collaboration: P. A. R. Ade et al. (2014); Kermish et al. (2012) will attempt even more precise measurements of the CMB than available so far, targeting a detection of primordial $B$-modes. In combination with the ongoing LIGO experiments, these missions seek to complement gravitational wave astronomy with cosmological gravitational waves from the very beginning of the evolution of the cosmic plasma, thereby opening up new possibilities to understand the properties of gravity and matter in extreme conditions.

The current status of the search for cosmological gravitational waves can be summarized by an upper bound on the tensor-to-scalar ratio, $r$, measured to be $r<0.12$ at 95$\%$ confidence at a pivot scale of $0.05\,\text{Mpc}^{-1}$ Planck Collaboration et al. (2016). To allow for a tensor-to-scalar ratio as low as $r\sim 10^{-3}-10^{-4}$, predicted by e.g. Higgs driven inflation Bezrukov & Shaposhnikov (2008) or Starobinsky inflation Starobinskii (1979), drastic technological advancement in the sensitivity of polarization detectors is required. Furthermore, it must not be forgotten, that simultaneous progress in the separation of CMB and foregrounds is prerequisite for a sensible interpretation of the resulting $B$-modes.

Future $B$-mode CMB experiments are designed to cover mainly the frequency range around 100 GHz111Different experiments measure or are planned to measure in different frequency ranges, e.g. LiteBird in the range 50-320 GHz, the CMB-S4 missions most likely in the range 35-250 GHz, CORE in 45-795 GHz, DeepSpace currently measures at 10-15 GHz but will be extended to include frequencies up to 143 GHz, PIXIE is planned to measure from 30 GHz up to 6 THz, and Polarbear currently measures at 150 GHz., where the relative intensity of the CMB is known to be highest. One of the most dominant foreground components in this range is thermal dust emission. Naturally, to be sensitive to values of $r$ as quoted above, the contribution of synchrotron emission and at presumably lower polarization fractions free-free emission, CO-line emission, and Anomalous Microwave Emission must not be neglected either. Furthermore, one should bear in mind alternatives to the single-component thermal dust model (Finkbeiner et al., 1999; Meisner & Finkbeiner, 2015), as well as extensions to given parameterizations for the inclusion of averaging effects, spatially as well as along the line-of-sight (Chluba et al., 2017).

In this paper we revisit three state-of-the-art solutions to thermal dust emission, namely Planck’s 2013 full sky model of the first public release Planck Collaboration et al. (2014) (hereafter P13), Planck’s Commander solution from the second data release Planck Collaboration et al. (2016) (hereafter C15), and a CIB-free solution derived in the so called GNILC framework (Remazeilles et al., 2011) also from Planck’s second data release Planck Collaboration et al. (2016) (hereafter P16). Differences of the methods used for obtaining these solutions are briefly described below.

In this note, we restrict ourselves to intensity measurements only, an analysis of thermal dust polarization products must await the advent of the future data. However, observations in intensity already reveal a lot about the model’s parameters. We focus on highlighting the differences of the parameters derived in the three different approaches and their relation to each other. After a first comparison of the power spectra of the parameters of P13 and P16 which both are provided at high resolution, we justify a smoothing of the maps in section II. In section III we then proceed to compare all parameter maps, including those of C15, at a common resolution of $2^{\circ}$, regarding their visual appearance, their power spectra, their distributions and their cross-correlations. We identify those areas on the maps with $|b|\gtrsim 50^{\circ}$ to be the main cause for the discrepancies, precisely the regions in which the CIB begins to dominate. In addition, in section IV, we contrast two zones on the sky, namely the BICEP2 zone, which has been on everyone’s lips in the past, and a region along the north celestial pole (NCP) representative for many ground-based experiments, such as QUIJOTE Hoyland et al. (2012); Gomez et al. (2010); Gómez-Reñasco et al. (2012); Rubiño-Martín et al. (2012); Génova-Santos et al. (2015); Rubiño-Martín et al. (2010) or DeepSpaceThe DeepSpace Collaboration (2016). We summarize in section V.

II Smoothing of the Dust Model parameter maps

In the optically thin limit, the spectral energy density of thermal dust emission at frequency $\nu$ and in direction ${\bf n}$ is modeled as a so called modified black-body (MBB) Cambrésy et al. (2001); Draine (2003); Draine & Li (2007); Hildebrand (1983); Planck Collaboration et al. (2014):

[TABLE]

where $\tau_{0}({\bf n})$ is the dust optical depth at a reference frequency $\nu_{0}$, and $\beta({\bf n})$ is the dust spectral index.

[TABLE]

is the Planck black-body function with dust temperature $T_{d}({\bf n})$.

The $\beta({\bf n})$, $T_{d}({\bf n})$ and $\tau_{0}({\bf n})$ maps of the three solutions to be compared here have different angular resolution and angular smoothing. In specific, the C15 maps have been smoothed by a $1^{\circ}$ Gaussian beam, significantly larger than the $5^{\prime}$ those of P13 and P16 have. In order to make a comparison between the three solutions possible, their angular smoothing must agree.

Due to the strong non-linearity of the thermal dust model, eq. 1, this transition must be done with care. The most straightforward would be to degrade the angular resolution of all intensity maps used in the respective thermal dust fits down to $\Theta=2^{\circ}$ by applying a Gaussian filter $G\left[(\textbf{n}-\textbf{n'})^{2}/2\Theta^{2}\right]$ and subsequently repeat the parameter fits. Unfortunately, this procedure is extremely complex and should not be done here. Instead, we argue, that even resulting, higher resolution maps of $\tau(\textbf{n}),\beta(\textbf{n})$ and $T_{d}(\textbf{n})$ can be smoothed in regions where their variation is small compared to the mean value or monopole. For this, consider the following decomposition of the respective parameter maps into spherical harmonics $Y_{lm}(\textbf{n})$:

[TABLE]

where $\tau_{lm},\beta_{lm},T_{d,lm}$ are the corresponding harmonic coefficients.

We compute the power spectrum $D(\ell)$ for each of the parameters as:

[TABLE]

where $X$ stands for $\tau,\beta$ or $T_{d}$. The monopoles for the respective parameters are computed to be $\bar{\tau}^{P13}\simeq 1.97\cdot 10^{-5}$, $\bar{\beta}^{P13}\simeq 1.61$ and $\bar{T}^{P13}_{d}\simeq 19.7\,K$ for P13, and $\bar{\tau}^{P16}\simeq 1.89\cdot 10^{-5}$, $\bar{\beta}^{P16}\simeq 1.60$ and $\bar{T}^{P16}_{d}\simeq 19.4\,K$ for P16. The standard deviations of $\beta$ and $T_{d}$ lie around $6\%$ of their mean values for both P13 and P16. In Fig. 1 we show the full sky power spectra for all three parameters, from both P13 and P16. Note, that the power spectra $D^{\beta}(\ell)$ and $D^{\tau}(\ell)$ were amplified for the sole purpose of visualization. Already here we want to highlight that, while the respective monopoles hardly differ, the multipole dependence of $D^{X}(\ell)$ reveals a significant mismatch, mainly concentrated at multipoles $\ell\gtrsim 50$, especially distinct for the spectral index. These differences arise from the increasing contribution of the CIB at the North and South polar caps, $|b|\gtrsim 50^{\circ}$, as we shall see in the next section. For now, we point to the fact that the amplitudes of the power spectra $D^{X}(\ell)$ of all maps are smaller than the corresponding monopoles222$\tau$ is an exception here because its variation can be orders of magnitude higher than the average value; $\tau$ is linear in the dust model, thus the separation in mean value and variation always works regardless of its amplitude.. For this reason, the intensity of dust emission eq. 1 can be represented as a first order approximation,

[TABLE]

As Gaussian smoothing is applied to a map via a linear operator, we now see from eq. 11, that it is equivalent to smooth the intensity map as suggested before and to smooth the individual parameter maps $\tau(\textbf{n}),\beta(\textbf{n})$ and $T_{d}(\textbf{n})$ afterwards. Note that the problem of extracting physically meaningful parameters in averaged regions on the sky (via instrumental beams, or data processing) has also been discussed in more details in (Chluba et al., 2017). While we present above considerations mainly to justify the smoothing of the parameter maps, the authors of (Chluba et al., 2017) go further and not only compute general expansions (beyond first order) of spectral emission distributions but also present solution approaches to beam and line-of-sight averaging effects.

III Full-sky comparison of the solutions

As justified in the last section, we now compare the parameter maps of P13, P16 and C15 after having smoothed them with a $\Theta=2^{\circ}$ Gaussian beam, see in fig. 2. Even with the naked eye one notices that while the maps for the optical depth are very similar, both the spectral index and the temperature maps have strong differences, despite their similar mean values. We begin to analyze these differences with the power spectra of the maps, using the above definition, eq. 10, see fig. 3.

In comparison to P13 and P16, the C15 solution is characterized by a lower and flatter power spectrum for the spectral index $\beta$ at low multipoles, $\ell\lesssim 70$, where $D(\ell)\propto const$. This seems to be compensated by more power in the temperature map on all scales. The optical depth experiences a loss of power towards higher $\ell$ at around $\ell\gtrsim 40$. Needless to say, all the maps presented in fig. 2 reveal significant statistical anisotropy, clearly seen in the even-parity asymmetry of the power spectra (the amplitude of $D(\ell)$ for $\ell=$“even” is systematically greater then for $\ell=$“odd”), most pronounced in the map/the power spectrum of the optical depth, due to the concentration of the signal along the Galactic plane – a source of “even”-power dominance Kim & Naselsky (2010, 2010, 2011), which is most pronounced in the Commander maps.

We proceed to compare number densities, $P(X)$, of all the parameter maps(at a HEALPix resolution of $N_{\text{Side}}=128$) in the intervals $[X-\Delta X,X+\Delta X]$ with $\Delta X=10^{-2}X_{max}$, where $X_{max}$ is the maximum of the corresponding parameter, see fig. 4. The distributions are normalized to the total number of pixels in the maps. It becomes obvious that the C15 distributions reflect precisely the choice of priors. The Gaussian priors in C15 were chosen to be $\beta^{C15}=1.55\pm 0.1$ and $T^{C15}=23\pm 3$; $\tau$ is calculated from a combination of the previous parameters and the dust amplitude, a parameter only constrained to be positive. Note that in P13 the allowed ranges of variation are $1.0\leq\beta^{P13}\leq 2.5$ and $10\leq T^{P13}\leq 60$ K.

In spite of different treatment regarding the CIB, the temperature and optical depth distributions of P13 and P16 look very similar. Only the distributions of the spectral index exhibit a discrepancy at $\beta\simeq 1.7$; it stands to reason that this difference is likely caused by the CIB.

Finally, we compute full-sky cross-correlations between $T_{d}$ and $\beta$ as well as $T_{d}$ and $\tau$, for each of the solutions, see fig. 5. The negative correlation between $T_{d}$ and $\beta$ is well known from P13 Planck Collaboration et al. (2014) and is also reproduced by P16. Also from the P13 and P16 maps in fig. 5 we find that for most regions of the sky lower temperature regions are accompanied by higher spectral index regions and vice versa. While the C15 solution at low temperatures shows negative correlation as well – yet less pronounced – it turns to become strongly positive at temperatures $T_{d}\gtrsim 21\,K$. Again referring to fig. 2 we find that the regions in which dust temperatures this high occur, happen to lie at higher Galactic latitudes $|b|\gtrsim 50^{\circ}$, those regions where the CIB begins to dominate. The negative correlation between $T_{d}$ and $\tau$ present for the P13 solution, is well reproduced by C15, however marginally for P16.

IV Local comparison of the solutions

In this section we would like to underline that the variability of the dust emission over the sky and the difference of the solutions forbids a flat-rate use of one particular choice of dust parameters to describe different regions on the sky properly. For this purpose we pick two regions on the sky as examples for areas under investigation in ground-based CMB experiments: the BICEP2 region in the Southern hemisphere, and a second region around the North celestial pole (NCP), observed by, e.g. QUIJOTE Hoyland et al. (2012); Gomez et al. (2010); Gómez-Reñasco et al. (2012); Rubiño-Martín et al. (2012); Génova-Santos et al. (2015); Rubiño-Martín et al. (2010) or DeepSpaceThe DeepSpace Collaboration (2016), see fig. 6, where we show the ratio of the dust emission intensity in P16 to that in P13. While in the BICEP2 zone the intensities differ by a comparably constant offset, the NCP region shows a strong gradient from the Galactic plane, were the dust models agree fairly well, decreasing to higher latitudes, however, with many small features. Again, we investigate the distributions of the parameters of the different solutions, now separately for the two regions, see fig. 7. The results differ considerably. Especially the dust temperatures and the spectral indices in the BICEP2 zone disagree strongly. As the latter are used for the extrapolation of dust intensities from higher to lower, CMB-friendly frequencies, this difference will immediately correspond to a difference in the obtained dust contamination.

For future evaluation of data by planned CMB experiments, which due to their frequency coverage cannot perform independent dust emission estimates, it is crucial to find concordance in which dust parameters to use. For frequencies below $\nu_{0}=353\,$GHz and characteristic dust temperatures, eq. 1 reduces to

[TABLE]

The parameters in the BICEP2 zone follow $\beta^{C15}<\beta^{P16}$, $T_{d}^{C15}>T_{d}^{P16}$, and $\tau^{C15}\simeq\tau^{P16}$. According to eq. 12 the intensity of dust emission at $\nu_{0}<353\,$GHz would be systematically greater in the C15 solution than in P16.

In the NCP region the distribution of the thermal dust intensity is more complicated, as there is much more variation present than in the BICEP2 zone. From the number densities only, this is most apparent in the broad distribution of dust temperatures peaked at 17 and 23 K, corresponding to the low and the high Galactic latitude parts of the region.

The comparison of these two zones clearly illustrates that the evaluation of the dust emission intensity at different frequencies cannot be done independently of the region investigated via e.g. the use of average values for the dust parameters. The 5-20% variations of these parameters will generate artificial signals, interfering with any measurement of the $B$-mode of polarization.

For completeness, we also show the $T_{d}$-to-$\beta$ and $T_{d}$-to-$\tau$ correlations in fig. 8 for both zones in question.

V Discussion and Conclusion

In this paper, we compared the three full-sky solutions to thermal dust emission, P13, C15, and P16 from the Planck Collaboration and found significant variations in their respective parameters, the spectral index $\beta$, the optical depth $\tau$ and the dust temperature $T_{d}$. Thereby P16 is the only solution which was cleaned of the CIB. At high resolution, the difference to P13 is distinctly visible in their power spectra at $\ell\gtrsim 50$.

After having smoothed the maps, comparison with those of C15 shows a strong deviation of the distribution of parameters to those of both former ones, obviously caused by the priors imposed in the C15 method. Hints towards an influence of the CIB also in the C15 maps was found by the systematically higher dust temperature in regions of higher Galactic latitude $|b|\gtrsim 50^{\circ}$, which is positively correlated with the spectral index, opposite to the case at lower temperatures/lower Galactic latitudes. However, no upturn was found for the $T_{d}-\beta$ correlation of P13. We summarize all cross-correlations between $T_{d}$ and $\beta$, and $T_{d}$ and $\tau$ in Table I and Table II, where for the sake of completeness we also include the correlations across methods. Note, that the $T_{d}-\beta$ anti-correlations are significant only for the P13 and P16 models.

Also note, that the parameters of the P13 and P16 models critically depend on the monopole corrections, presented in Table 3. Potentially these corrections could be one of the major sources of uncertainties for future CMB $B$-mode experiments.

We have also investigated the distributions of the spectral index $\beta$, the optical depth $\tau$ and the dust temperature $T_{d}$ maps in all three models in two regions of the sky: the BICEP2 zone and one around the North Celestial Pole. We have shown that in the BICEP2 zone the spectral index for P16 model is shifted to values around $\beta\simeq 1.8$ compared to the C15 spectral index at $\beta\simeq 1.55$. In the NCP zone the behavior of the spectral indices and the dust temperatures are very different than those in the BICEP2 zone. Here C15 predict the temperature distribution function with two maxima at 17K and 23K. Most likely, the highest values of the dust temperature in NCP zone are related to the CIB contamination (mainly due to its monopole) of the C15 MBB parameters.

In conclusion, all methods compared here agree to about $5-20\%$. However, in the light of current ambitions of planned $B$-mode experiments, much higher accuracy is needed. It is curious to note, that the different solutions agree best along the Galactic plane at $|b|\lesssim 50^{\circ}$–those regions which are usually excluded from analyses due to the large dust emission amplitudes. The regions with $|b|\gtrsim 50^{\circ}$ show less emission by the Galaxy, however the increased relative contribution of CIB, instrumental noise possibly other components as anomalous microwave emission, in the context of the Commander model investigated in (von Hausegger & Liu, 2015), causes disagreement among the different dust emission solutions. Due to rapid changes of the dust parameters over the sky, we warn against the use of average values of the model’s parameters. The residuals introduced by this approach would themselves act as a ”new” foreground, leading to biases in the polarized component separation (Such biases are nicely shown in form of altered spectral shapes presented in (Chluba et al., 2017) for intensity measurements). Assuming that one knows the correct model(s) describing the foregrounds, supplemented by observations at high enough frequencies above 300 GHz, we strongly promote local (pixel- and multipole-based) fitting methods ensuring proper treatment of both thermal dust emission and CIB, despite its higher computational cost. Only then can one ensure the reliable extrapolation of foreground emission to the frequency range of interest, around $100\,$GHz, where ultimately also low frequency foregrounds must be treated to the same level of precision. In general, the next generation of $B$-mode missions needs more sophisticated methods of foreground removal, which carefully treat the averaging effects described here, as well as those along the line-of-sight, possibly based on further development of the Commander and P16 algorithms including the CIB and less strong priors on the parameters.

Acknowledgements.

We thank Carlo Burigana, Per Rex Christensen, Jacques Delabrouille, Hans Kristian Eriksen, Andrew Jackson, and Philipp Mertsch for comments on the draft and useful discussions. We also thank the anonymous referee for useful comments on improving the manuscript. This work was partially funded by the Danish National Research Foundation (DNRF) through establishment of the Discovery Center and the Villum Fonden through the Deep Space project. Hao Liu is supported by the National Natural Science Foundation for Young Scientists of China (Grant No. 11203024) and the Youth Innovation Promotion Association, CAS.