© ESO, 2012

1. Introduction

Although it has been common to characterize the large-scale Galactic diffuse emission by treating each component of the Galaxy individually – the neutral gas distribution, HII region maps, the magnetic field structure – this emission is usually caused by a convolution of different combinations of these components. Far-infrared dust emission, for example, traces the dust distribution heated by stars; radio free-free emission depends specifically on the local density of main-sequence massive stars; radio synchrotron emission traces the local cosmic ray electron density (and energy spectrum) convolved with the magnetic field structure; and so on. Since our knowledge of the distribution of any one of these components is still rather incomplete, modeling emission that arises from the coupling of two or more of them is a daunting task, requiring numerous assumptions to develop a complete model that can be compared with observations.

The value of such an exercise lies in identifying the factors that have the largest influence on determining the properties of the observed emission, with the hope that eventually all of the different diffuse emission measurements can be combined to yield a more complete picture of the Galaxy. In this paper, we present models of the stellar and diffuse infrared emission in the Galactic plane. Although investigations of the diffuse infrared emission have been carried out in the past, to our knowledge this is the first work to focus on the mid-infrared (3.6–8.0 μm) diffuse emission from polycyclic aromatic hydrocarbons (PAHs), which is much more sensitive to the distribution of the young stars than the longer wavelength dust emission probed by IRAS and COBE/DIRBE.

The large-scale distribution of stars and diffuse dust emission in the Milky-Way has been the subject of numerous studies. On one hand, powerful three-dimensional stellar distribution models have been developed to simulate the inventory of stellar populations as a function of location in the Galaxy and position on the sky. The Besançon model (Robin & Creze 1986; Bienayme et al. 1987; Robin et al. 1996, 2003), the SKY model (Wainscoat et al. 1992; Cohen 1993, 1994, 1995), and the TRILEGAL model (Girardi et al. 2005) are three of the most widely known such models aimed at reproducing the stellar populations of the Galaxy. Such models typically include various components to describe the Galaxy, such as a bulge, halo, disk(s), spiral arms, and so on, and include stars at various masses, metallicities, and stages of evolution.

On the other hand, analytic models of the diffuse dust emission in the Galaxy have also been developed to reproduce the large-scale mid- and far-infrared emission. For example, Sodroski et al. (1997) used COBE/DIRBE data from 12−240 μm to determine the best-fitting abundances and temperature of large dust grains, the abundance of very small transiently heated dust grains and PAHs, and the energy density of the interstellar radiation field as a function of Galactocentric radius. Similarly, Davies et al. (1997) fit a model for the emission from cool (18–22 K) dust to the COBE/DIRBE far-infrared data. More recently, the Planck Collaboration (2011) fit the observations of the Galaxy by the Planck satellite, using a model similar to that used by Sodroski et al.

Drimmel (2000) carried out a qualitative comparison of COBE/DIRBE K-band and 240 μm observations, and found that while the near-infrared shows evidence only for two stellar arms, the far-infrared observations are consistent with four arms (in agreement with the influential HII region study of Georgelin & Georgelin 1970), leading to the conclusion that the main Galactic potential is two-armed, but that the structure of the dust and gas responding to the potential is more complex, and can be adequately characterized as four-armed. Drimmel & Spergel (2001) extended this by developing a quantitative model to fit the above data, which consisted of a parametric model of the dust distribution and associated far-infrared emission. This dust distribution was then used by Drimmel et al. (2003) to predict the three-dimensional K-band extinction.

What has been lacking in previous models is a self-consistent treatment of the source of the dust heating – the photons produced by stars of various temperatures and luminosities – and the target of the heating – grains of different sizes and PAHs – distributed throughout the disk of the Galaxy. In this paper we present a radiative transfer model of the Galaxy at infrared wavelengths, considering both stellar and dust emission. The model is self-consistent in the sense that the heating of the dust does not follow an analytical prescription, but is instead directly computed from the stellar populations. The main aims of this paper with this model are the following:

• To determine whether an existing Galactic model for stellarpopulations and dust – here the SKY model – used in conjunctionwith the dust properties determined by Draine &Li (2007), which include transiently heated verysmall grains and PAH molecules, adequately reproduces thestellar and diffuse emission seen at mid- and far-infraredwavelengths.

• To determine whether any modifications are needed to this model to reproduce the observations, and whether any new insights on Galactic structure can be gained. However, we stress that the aim is not to provide an exhaustive parameter space study, as given the number of stellar populations and components, the number of parameters describing the model is on the order of several hundred. Rather, the aim is to determine possible and realistic modifications that improve the model.

• To determine the relative contributions of resolved and unresolved stellar flux, emission from PAH molecules and dust grains, and scattering at infrared wavelengths.

• To determine the relative importance of various stellar populations in heating dust grains and exciting PAH molecules.

In Sect. 2 we describe the observational datasets used in this paper. In Sect. 3, we describe the initial model used, the dust properties, and the radiative transfer code. In Sect. 4, we show results for both the initial model (Sect. 4.1) as well as an improved model (Sect. 4.2). In Sect. 5 we examine the main contributors to the observed flux and to the heating of the dust (Sect. 5.1), we study how much IRAC stellar flux is likely to be unresolved (Sect. 5.2), and we generate images of the Galaxy from an external viewpoint (Sect. 5.3). Finally, in Sect. 6 we summarize our findings.

2. Observations

In this paper we make use of the Spitzer GLIMPSE (Benjamin et al. 2003; Churchwell et al. 2009) and MIPSGAL (Carey et al. 2009) surveys of the Galactic mid-plane, specifically GLIMPSE/MIPSGAL I, which cover 10° ≤ |ℓ| ≤ 65° and |b| ≤ 1°, and GLIMPSE/MIPSGAL II, which fill in the region for |ℓ| < 10°, with |b| ≤ 1° for |ℓ| > 5°, |b| ≤ 1.5° for 2° < |ℓ| ≤ 5°, and |b| ≤ 2° for |ℓ| ≤ 2°. The total area covered by these surveys is 274 deg².

The GLIMPSE surveys were carried out using the InfraRed Array Camera (IRAC; Fazio et al. 2004) at 3.6, 4.5, 5.8, and 8.0 μm. The MIPSGAL surveys were carried out using the Multiband Imaging Photometer for Spitzer (MIPS; Rieke et al. 2004) at 24 and 70 μm, though only the processed 24 μm data are available and are included in this study.

In addition to the Spitzer data, we include data from the Improved Reprocessing of the IRAS Survey (IRIS; Miville-Deschênes & Lagache 2005) at 60 and 100 μm to provide further constraints on the dust emission.

At the different wavelengths considered, the total emission in a given direction consists of stellar emission resolved into point sources, unresolved stellar emission, and diffuse interstellar emission. The relative contribution of these three types of emission change both as a function of wavelength (diffuse emission increases with wavelength), and position in the Galaxy (unresolved stellar emission decreases with increasing angle from Galactic center). Here, longitude and latitude profiles were constructed for the four IRAC bands, the MIPS 24 μm data, and the IRIS data, by resampling the full data into 3′ by 3′ bins, including the resolved stellar flux. The re-binned data were then collapsed into selected longitude and latitude profiles to compare the observations to the radiative transfer model presented in Sect. 4. In this paper, we only use the rectangular region for |ℓ| < 65° and |b| < 1° to facilitate computing longitude and latitude profiles.

3. Model

Our self-consistent radiative transfer model of the stellar and dust emission at infrared wavelengths consists of two key ingredients: a model for the distribution of the stellar populations, and a model for the characteristics and distribution of the dust. The initial distributions used are described in Sect. 3.1 and in Sect. 3.2, but these are then modified in Sect. 4.2 to provide a better fit to the observations. The dust properties assumed are described in Sect. 3.3. Finally, the radiative transfer model is described in Sect. 3.4.

3.1. Stellar populations

To model the distribution of various stellar populations, we implemented the SKY model from Wainscoat et al. (1992, hereafter W92) and Cohen (1993, 1994, 1995). The choice of the SKY model over more recent models such as the Besançon or TRILEGAL models is that the former provides a prescription for spiral arms, which are essential for the present work, since we expect the diffuse PAH emission at IRAC wavelengths to depend on the distribution of massive stars, which are concentrated in the spiral arms. In addition, the three-dimensional distribution of the stellar populations in the SKY model is separated into spectral classes, making it straightforward to represent in a radiative transfer model.

The SKY model includes five components to describe the Galaxy: an exponential disk, a bulge, a halo, spiral arms, and a ring. We do not include the halo: as it is defined in W92 by a projected rather than a spatial distribution, it is more difficult to include in the radiative transfer model. However, it has a very low normalization factor relative to the disk (1:1250; Cohen 1995) and would have no noticeable impact on the Spitzer or IRAS observations in the Galactic mid-plane. The remaining four main components are briefly described in the following sections.

The model includes 87 different types of Galactic sources, including pre-main-sequence stars, main-sequence stars, giants, asymptotic giant branch (AGB) stars, supergiants, and planetary nebulae. Each component includes a normalization factor (ρ_D, ρ_B, ρ_A, ρ_r), and in some cases other parameters that depend on the source type, denoted S. For example, O, B, and A stars are absent from the bulge (typically composed of more evolved stars), while lower-mass main sequence stars are absent from the spiral arms (typically composed of relatively young stars). In the following sections, all parameter values that are functions of the source type S are given in Table 2 of W92. The following sections give an overview of the distributions used; the reader is referred to W92 for more details. We note at the outset that the parameters characterizing these components, and in some cases even their existence (such as the molecular ring), are uncertain. In the discussion that follows, we will explore the effect of changing some of these components.

3.1.1. The exponential disk

The exponential disk has a density distribution given in cylindrical polar coordinates by where h_z(S) is the disk scale height, which depends on the spectral class S, R₀ = 8.5 kpc is the distance from the Sun to the Galactic center, and h = 3.5 kpc is the radial scale-length. The exponential disk is truncated at R_max = 15 kpc.

3.1.2. The bulge

The bulge has a density distribution given by where In the W92 model, the bulge axis ratio k₁ is set to 1.6, and the bulge radius R₁ is set to 2 kpc.

3.1.3. The spiral arms

The original W92 model included four main logarithmic spiral arms, and one spur in the Solar neighborhood. Cohen (1994) subsequently replaced the single local spur by two narrower spurs on either side of the Sun. The shape of the arms is parameterized as where α is the winding constant, R_min is the inner radius, and θ_min is the angle at the inner radius. The arms go out to R_max, but are truncated above a certain length. The parameters for the arms (compiled from Table 2 of W92 and Cohen 1994) are given in Table 1. The arms follow the same radial and vertical exponential distributions as the disk, and have a constant width. Thus, their density distribution is effectively given by inside the arms, and ρ(R,z,S) = 0 outside. As mentioned previously, the stellar population in this component is dominated by young, massive stars.

3.1.4. The stellar (molecular) ring

The SKY model contains a stellar ring with a density given by where R_r = 6.75 kpc is the ring radius, and σ_r = 0.96 kpc is the ring width. Since this component was motivated by an annulus of molecular gas in the inner galaxy, the stellar populations in this component are also weighted to young massive stars, as is the case for the spiral arms.

3.2. Dust distribution

The second main ingredient for the Galactic model is a distribution of dust, which is initially taken from the extinction distribution of W92. For the dust properties described in Sect. 3.3, the extinction distribution from W92 corresponds to (1)where h_d = 3.5 kpc, z_d = 100 pc, and ρ_d = 10^-25 g/cm³.

Given the clumpiness of the molecular gas of the Galaxy, the smooth exponential distribution of dust is clearly an oversimplification, but serves as a useful starting point for the model presented here.

3.3. Dust properties

The dust model adopted is from Draine & Li (2007). The model consists of a mixture of carbonaceous and amorphous silicate grains with a size distribution from Weingartner & Draine (2001) that reproduces interstellar extinction curves with R_V = 3.1. The model uses the renormalization relative to H from Draine (2003). The size distribution includes a population of small grains in two log-normal components in the size distribution, to represent transiently heated very small grains and PAH molecules. The model from Draine & Li used here is the one with a mass fraction of PAHs set to q_PAH = 4.6%.

Draine & Li computed emissivities for dust grains illuminated by the interstellar radiation field in the solar neighborhood from Mathis et al. (1983) scaled by a factor – denoted U – which varied between 0.5 and 10⁷. As described in Draine & Li (2007), the ionization fraction of the PAHs is a function of grain size. The emissivities and opacities used were identical to those publicly available¹ but were separated into three size regimes (B. Draine, priv. comm.): ultra small grains (a < 20 Å), very small grains (VSGs; 20 Å < a < 200 Å), and big grains (a > 200 Å). These three regimes contain 5.86%, 13.51%, and 80.63% of the mass respectively in the dust model. Since the ultra-small grains are dominated by PAH molecules, we will refer to them as PAHs in the remainder of this paper (but we note that they are not strictly equivalent).

3.4. Radiative transfer model

3.4.1. Model set-up

The radiative transfer models presented in this paper were computed with the Monte-Carlo three-dimensional radiative transfer code Hyperion² (Robitaille 2011). This code can read in any arbitrary dust geometry, with multiple dust components, and multiple sources, making it perfectly suited to this project. It also allows images to be produced relative to an observer inside the density grid, which we make use of in this paper. The code implicitly solves the full equation of radiative transfer, including optical depth effects, scattering, polarization, and emission under the assumption of radiative equilibrium.

The stellar and dust densities were discretized onto a cylindrical polar grid with 200 radial, 50 vertical, and 100 azimuthal cells. The radial cells were distributed linearly between 0 and 15 kpc. The vertical cells were logarithmically spaced from 1 pc to 3 kpc above and below the mid-plane, with the region inside  ± 1 pc divided into two cells. The azimuthal cells were uniformly distributed from 0 to 2π.

Rather than treating sources individually, which would be unfeasible computationally, “diffuse” sources of emission were used to represent populations of stars: for each spectral class, each cell in the grid was given a probability for emission, and any photon emitted from that source was given the spectrum corresponding to the spectral class. Each diffuse source was given a total luminosity derived from the total number of “real” sources.

By default, the calibrated spectra used for each spectral class were derived from the absolute magnitudes given in Table 2 of W92 for BVJHK and 12 μm and 25 μm. However, for main-sequence stars, giant, and supergiants, we used Allen (1973) to transform the spectral types into effective temperatures and surface gravities, which we subsequently used to interpolate stellar photosphere models from Castelli & Kurucz (2004). These model spectra were then scaled to the absolute magnitudes given in Table 2 of W92. The motivation for doing this is that the B magnitudes given in the SKY model are not sufficient to properly characterize the UV portion of the spectrum, which is crucial for the heating of the PAHs and VSGs. For the “Young OB” source type in W92, which is meant to represent OB associations, we fit the Castelli & Kurucz spectra to the magnitudes supplied by W92, and found that an effective temperature of 15 000 K provided the best fit to all wavelengths.

Each of the dust types (PAHs, VSGs, and big grains) was given the spatial distribution from Eq. (1), or the modified spatial distribution defined in Sect. 4.2.2, scaled by the mass fraction of the dust types (cf. Sect. 3.3).

3.4.2. Caveats

Fig. 1 Spectrum of the average surface brightness inside −65° < ℓ < 65° and −1° < b < 1° for the initial model (left) and the final modified model (right). The thick colored solid lines represent different components of the model: supergiants (purple), AGB stars (blue), giants (teal), main sequence stars (green), other stellar types (dark gray), dust emission from PAHs (red), VSGs (orange), and big dust grains (yellow), and scattered light (black). The major components are labeled in the left panel. The shaded gray area shows the total flux. The open black circles show the model fluxes after convolution with the transmission curves, and the filled black circles show the observations.

The Hyperion radiative transfer code uses pre-computed emissivities for the grains for a range of interstellar radiation field strengths, rather than computing the transient heating of very small grains and PAH molecules exactly. As described in Sect. 3.3, the dust emissivities were pre-computed for a fixed illuminating spectral shape – but variable intensity U – so for regions in the Galaxy with different radiation field shapes, the shape of the emissivity will likely be different. The primary effect of a different spectral shape is to change the amount of radiation absorbed by and therefore heating the PAHs. Since PAHs have an opacity that is highest in the UV, harder radiation fields with a fixed intensity U will lead to higher excitation levels. This can be accounted for by quantifying the radiation field not by U (which ignores the wavelength dependence of the emission) but by the power of the radiation field absorbed by the grains (per unit mass): where J_isrf is the mean intensity of the radiation field. Thus, two radiation fields with the same intensity U but different spectral shapes – which would excite the PAHs by different amounts – would have correspondingly different values of Ȧ. Therefore, selecting the emissivities based on Ȧ is a better approximation than simply using U and takes into account the effects of spectral shape to first order.

With this effect taken into account, the difference between the emissivities due only to the change in spectral shape for a constant Ȧ is quantitatively small for our purposes. As shown in Fig. 7b of Draine (2011), for two very different radiation fields with the same³Ȧ, the difference in the PAH emissivities is less than a factor of 1.5. In this case the emissivities are systematically offset, and the colors do not change noticeably. Since at least some of the ISM will in fact have an interstellar radiation field close to that used to compute the dust emissivities, this effect may be smaller in reality.

As shown in see Fig. 7a of Draine (2011), the difference in the emissivities between fully ionized and fully neutral PAHs is more important. In particular, different levels of ionization can lead to differences in emissivities of factors of two, as well as significant changes in colors.

We conclude that the assumption of the specific dust model used, which includes an assumption for the ionization level of the PAHs, is the main caveat in the radiative transfer models, while the method of choosing the emissivities based on Ȧ, while approximate, appears to be adequate for our purposes.

3.4.3. Model output

Model images were computed for the survey area and collapsed into longitude and latitude profiles in the same way as the data. The latitude resolution of the profiles is 3′, while the longitude resolution is 1°. The images were computed for 160 wavelength bins logarithmically spaced from 3 to 140 microns, and were subsequently convolved with the transmission curves for IRAC, MIPS, and IRAS (following Appendix A of Robitaille et al. 2007). The convolution is especially important because the PAH features in the emissivities cause the flux to vary rapidly with wavelength in the mid-infrared, and picking a single wavelength instead of taking into account the proper transmission curve could result in the model fluxes being wrong by a factor of two or more.

For all four models presented in Sect. 4 we have computed the variance between the data and the model in log₁₀ space to quantify the goodness of fit, and we list these in Table 2. We do not compute χ² values, since the model is clearly misspecified – that is, the model is only a rough approximation of reality because it does not take into account the small-scale features in the emission (which we are not attempting to reproduce) – so that absolute χ² values and associated likelihoods would not be meaningful. We compute the variances separately for longitude and latitude profiles, and provide both the variances for individual wavelengths, and the overall values.

4. Results

4.1. Original model

In this section, we show the results from the model as described in Sect. 3 with no modifications.

Figure 1 shows the average spectrum of the surface brightness inside −65° < ℓ < 65° and −1° < b < 1° for both the observations and the model. The model underestimates the surface brightness at IRAC wavelengths and overestimates it in the far infrared at IRAS wavelengths.

In addition to the averaged surface brightness spectrum across the survey area, we also examine the spatial distribution of the emission. Figure 2 shows the observed and model longitude and latitude profiles from 3.6 to 100 μm. The emission in the model is on the right order of magnitude, and reproduces a few important features of the observations, including the scale height of the emission (longwards of IRAC 5.8 μm) and the flattening of the IRAC 5.8 and 8.0 μm emission inside |ℓ| < 30°, but there are also a number of issues:

• the model surface brightness is too low by a factor of two or moreforℓ < −30° for IRAC 3.6 μm to MIPS 24 μm, and for ℓ > 30° for IRAC 5.8 μm to MIPS 24 μm. Since ℓ = 30° corresponds to the tangency of the ring component for our model, this suggests that the model for the disk and/or spiral structure outside this radius are not adequate.

• the model surface brightness is too high for |ℓ| < 30° for IRAS 60 μm and 100 μm, which suggests that our model contains too much dust inside R = 4 kpc. It has previously been suggested by numerous authors that a hole exists in the distribution of gas within a few kpc of the Galactic center (cf. Ortiz & Lepine 1993) – presumably due to the bar – which would explain the excess flux predicted by the model here.

• the IRAC 3.6 and 4.5 μm latitude surface brightness profiles are not as peaked as the data. This suggests that either the scale height of the stellar populations is not small enough, or that dust extinction is too large in the mid-plane, since dust extinction is included in the model. The shape of the latitude profiles is dominated by the inner Galaxy (|ℓ| < 30°) since that is where most of the emission originates, so changes in this range will have the largest impact on the latitude profiles.

Fig. 2 Longitude and latitude surface brightness profiles in the range −65° < ℓ < 65° and −1° < b < 1° for the observations and for the initial model (Sect. 4.1). The thick black lines show the observed surface brightness, while the colored lines and the gray shaded area show the model surface brightness, with the same colors as used in Fig. 1. The “traditional” spiral arm tangencies are indicated in the top panel.

Fig. 3 Longitude and latitude surface brightness profiles as for Fig. 2, but for the model including two major and two minor spiral arms instead of four main spiral arms (Sect. 4.2.1).

Fig. 4 Longitude and latitude surface brightness profiles as for Fig. 3, with a deficit in the dust distribution in the central few kpc of the Galaxy (Sect. 4.2.2).

Fig. 5 Longitude and latitude surface brightness profiles as for Fig. 4, with 50% more PAHs (Sect. 4.2.3).

4.2. Improved model

In this section, we examine how the model presented in Sect. 4.1 changes as we modify our initial assumptions in ways that have been motivated by previous investigations of Galactic structure. The modifications explored are changes to the parameters for the spiral arms (Sect. 4.2.1 and Fig. 3), modifications of the dust distribution interior to R = 4 kpc (Sect. 4.2.2 and Fig. 4), and modifications to the abundance of PAHs (Sect. 4.2.3 and Fig. 5). We also explored the effect of concentrating the dust in the arms, but found that the simulations most consistent with the observations were the ones in which the dust is substantially broader than the stellar arms.

4.2.1. Modified spiral structure

Our modifications to the spiral structure model are motivated by the low model surface brightnesses for |ℓ| > 30° at IRAC and MIPS wavelengths relative to the observations

Previous research has shown that when using stellar tracers in the near- or mid-infrared, there is only evidence for two of the spiral arms in the Galaxy (e.g. Drimmel 2000; Drimmel & Spergel 2001), specifically the Scutum-Centaurus arm and the Perseus arm. On the other hand, surveys in the far-infrared and radio show that models with four arms are a better match. This suggests that the main potential of the Galaxy is in fact two armed, but that the non-linear response of the gas and dust flowing through the stellar potential is more complex, and includes arms, or large spurs, between the two major spiral arms.

We therefore modify the model so that the Scutum-Centaurus and Perseus arms (2 and 2′ in Table 1) become the dominant arms by increasing the normalization factor of the O and B spectral class populations in these arms by a factor of two. Rather than completely eliminating the Norma and Sagittarius arms (1 and 1′ in Table 1), which does not provide a good fit to the dust-dominated bands, we eliminate all stars except those belonging to the O and B spectral classes, which we leave unchanged, under the assumption that stars are indeed forming in these secondary arms but do not stay there in the long term.

The shape of these four arms is also changed to Gaussian profiles in the radial direction, with σ = 550 pc. Finally, we modify the spiral arm parameters slightly to provide a better match to the shape of the longitude profiles. The updated spiral arm parameters are listed in Table 3. The local spurs are left unchanged. The resulting model longitude and latitude profiles, shown in Fig. 3 provide a better fit to the observations at IRAC and MIPS wavelengths for |ℓ| > 30°, although we note that the IRAC 5.8 and 8.0 μm fluxes are still slightly too low in that longitude range. This is reflected in the decrease in the variance between the model and the data listed in Table 2 for IRAC and MIPS wavelengths.

4.2.2. Modified dust distribution

In order to provide a better fit to the IRAS 60 and 100 μm longitude profiles inside |ℓ| < 30°, we also experimented with including a hole in the dust within a few kpc of the Galactic center by modifying the radial dependence of the dust distribution from a simple exponential to a function with the following form: where μ₀ and σ₀ are parameters, and f₀ and R_smooth are defined as the normalization constant and transition radii for which the two functions transition smoothly, meaning that the functions and their derivatives are equal. These values are: and We found that values of μ₀ = 4.5 kpc and σ₀ = 1 kpc provided a model that fit the IRAS profiles well, while maintaining the fit to the IRAC and MIPS profiles. A comparison of the new density profile with the original one is shown in Fig. 6. The surface brightness profiles are shown in Fig. 4, and the corresponding average spectrum is shown in Fig. 1. The MIPS 24 μm, IRAS 60 μm, and IRAC 100 μm now appear to fit the observed profiles well, while the IRAC 5.8 μm and 8.0 μm remain systematically low.

As shown in Fig. 4 and Table 2, the addition of the hole in the dust distribution also improves the fit of the IRAC 3.6 and 4.5 μm latitude profiles, which are now peaked in the same way as the observations, thanks to the reduced mid-plane extinction toward the Galactic center.

4.2.3. Modified PAH abundance

The fact that the systematic offset at IRAC wavelengths is seen at all longitudes suggests that it cannot be solved by changes in the geometry alone, but rather by a change in the large-scale properties of the dust or radiation field. For example, one solution might be to increase the abundance of PAHs relative to larger grains. To test this, we computed models with an increased dust density for the PAHs (strictly speaking the dust grains with a < 20 Å), and found that a model with 50% more PAHs (corresponding to q_PAH ~ 6.7%) provides a good fit to the observations⁴. The resulting integrated spectrum is shown in Fig. 1 and the profiles are shown in Fig. 5: the fit to the IRAC 5.8 μm and 8.0 μm profiles is significantly improved (as seen in the decrease in the variance for these bands in Table 2), and the model now provides a reasonable fit at all wavelengths. This result suggests that – if the model for the stellar populations is correct – the fraction of PAHs in the ISM may be higher than previously found, although we caution that due to the large number of parameters involved, and given the caveats of the model (Sect. 3.4.2), this result is only tentative.

We note that the original value of q_PAH = 4.6% found by Draine & Li (2007) was also determined by using emission from the Galactic plane, including data from the GLIMPSE survey. Draine & Li implicitly assume that the emission and colors of the GLIMPSE emission are due to optically thin emission, but our model indicates that in the mid-plane, extinction can affect the IRAC fluxes by up to 50%. We ran a model without increasing the PAH abundance, instead computing the profiles as if the Galaxy was optically thin, and found that the profiles and total flux were a better match to the observations, and fit almost as well as when increasing the PAH abundance. Therefore, it could be that our result is consistent with Draine and Li, and that the true PAH fraction is higher than q_PAH = 4.6%, but appears to be lower due to the extinction of the IRAC wavelengths relative to the far-infrared.

Fig. 6 Relative radial dependence of the dust density for the original model presented in Sect. 4.1 (gray) and the model with an inner dust hole from Sect. 4.2.2 (black).

We also note that there may be other solutions to the systematic offset at IRAC wavelengths. For example, the clumpy nature of the ISM may affect the average energy input of the various dust types in a systematic way, resulting in different relative contributions of emission from the different dust types at different wavelengths, but we defer such an investigation to a future study.

Fig. 7 Fractional contribution of various stellar populations and dust grain sizes to the total flux in the survey area as a function of wavelength. The colors are as in Fig. 1. Only the main contributing spectral types and dust types are labeled.

Fig. 8 Fractional contribution to dust heating as a function of individual spectral type. The colors are as in Fig. 1, and the darker vertical lines indicate separations between individual spectral types. Only the main contributing spectral types are labeled.

5. Analysis

In this section, we use the improved model as described in Sect. 4.2, but we note that the qualitative results in Sect. 5.1 and Sect. 5.2 are not significantly changed if we use the initial model instead.

5.1. Flux and energy breakdown

Fig. 9 Ratio of stellar emission measured in the GLIMPSE point source Archives to the total stellar flux predicted by the final model. This suggests that unresolved stellar flux contributes a large component of the total stellar flux inside |ℓ| < 20°.

Figures 1 to 5 show the breakdown of the surface brightness into contributions from the different components of the models, split up into main-sequence stars, giants, AGB stars, supergiants, other spectral classes, the three dust size ranges, and scattered light. We summarize this information into a more intuitive representation in Fig. 7. Approximately half of the IRAC 3.6 and 4.5 μm flux is from giants, with the second and third biggest contributors being PAHs and main-sequence stars respectively. On the other hand, the IRAC 5.8 and 8.0 μm bands are strongly dominated by emission from PAHs. MIPS 24 μm is dominated by emission two thirds from PAHs (in the form of continuum emission rather than features) and one third from VSGs, while IRAS 60 and 100 μm are dominated by thermal emission from the larger dust grains.

A similar analysis can be done to understand which stellar populations contribute the most to the energy injected into and reprocessed by the dust populations. Figure 8 shows the breakdown of this energy by individual spectral type (with the main contributing spectral types indicated). Heating of the dust grains is dominated by B-type stars and to a lesser extent, OB associations. Since O stars are rarer and have shorter lifetimes, the heating from the O stars is likely to be less uniform, since O stars are mostly found in clusters. In fact, the strong peaks in emission seen in the longitude profiles, which correspond to regions of massive star formation, are likely to be heated by a larger fraction of O-type stars. For the larger dust grains, which are less sensitive to UV radiative than the PAHs and VSGs, giants provide around a third of the heating.

5.2. Unresolved flux

The model derived in Sect. 4.2 can be used to understand how much of the diffuse emission in the IRAC bands is due to unresolved stellar flux. In order to study this, we first computed the mean surface brightness of point sources in the GLIMPSE Archive source lists in 1° bins in longitude, and compared this to the model stellar flux (excluding the dust emission). Figure 9 shows the ratio of stellar flux measured in the GLIMPSE Archives compared to the model flux. For |ℓ| > 20°, the agreement between the resolved stellar flux and the model is good, suggesting that there is little or no unresolved stellar flux at these longitudes. Inside |ℓ| < 20°, the observed fluxes are as low as half of the predicted flux in places, suggesting that in that longitude range, the unresolved stellar flux could be as important a component to the total stellar flux as the resolved stellar flux.

5.3. External viewpoint

Fig. 10 Color composite view of the Galaxy model from an external viewpoint, viewed along the North Galactic Pole, for the initial model (left) and the final modified model (right). The Sun is located at (x, y) = (0, −8.5 kpc). The colors show IRAC 8.0 μm (red), IRAC 4.5 μm (green), and IRAC 3.6 μm (blue), all on a square-root intensity scale.

As well as computing the appearance of the model from the position of the Sun, we can also compute the appearance of the Galaxy from an external viewpoint. In Fig. 10, we show what the model would look like at IRAC wavelengths viewed along the North Galactic pole. The spiral arms are prominent in PAH emission, while the bulge is bright in IRAC 3.6 μm and 4.5 μm. The molecular ring also features prominently, although as noted previously, there is debate as to whether this is in fact a real component of our Galaxy, or whether it is simply formed by an inward extension of the major spiral arms (see e.g. Dobbs & Burkert 2012).

6. Summary

A radiative transfer model of the Galaxy which uses the SKY model in conjunction with the dust properties from Draine & Li (2007) was developed to self-consistently calculate the heating of dust grains, and verify whether it is able to reproduce the observed surface brightness from 3.6 μm to 100 μm. The main findings presented in this paper are the following:

• 1.

  The initial model is able to roughly reproduce the order ofmagnitude of the observed surface brightnesses observed, butthere are disagreements between the model and observations,notably for |ℓ| > 30° for the IRAC and MIPS bands, and inside |ℓ| < 30° for IRAS 60 μm and 100 μm.

• 2.

  By modifying the model to incorporate two major stellar spiral arms and two secondary spiral arms with only young massive stars, as well as removing dust from the central few kpc of the Galaxy, we are able to significantly improve the quality of the fit. A slight systematic offset remains at IRAC 5.8 and 8.0 μm. By increasing the abundance of PAHs by 50%, we are able to eliminate this systematic offset, though we caution that other effects may explain this offset.

• 3.

  The flatness of the IRAC 5.8 and 8.0 μm and MIPS 24 μm emission is not directly due to a hole in the dust distribution, but is a consequence of the lack of strong UV sources within the inner few kpc of the Galactic center. The evidence for a deficit of dust is only apparent at longer wavelengths, in the IRAS 60 and 100 μm bands.

• 4.

  Since the stellar populations are not symmetrically distributed, the heating of the dust is not uniform, and therefore the spiral arms appear in diffuse emission at infrared wavelengths without the need for a non-symmetrical dust distribution (i.e. the model does not require additional dust in the spiral arms). This is not to say that the distribution of dust is not enhanced in the spiral arms, but that there is no evidence from the data that it is.

• 5.

  The overall flux at IRAC 3.6 and 4.5 μm is dominated by giant stars, PAHs, and main-sequence stars, while the flux from IRAC 5.8 μm to IRAS 100 μm probes mostly dust, with the size of the dust grains probed increasing with wavelength.

• 6.

  On large scales, transiently heated very small grains and PAH molecules are predominantly heated by B-type stars. The larger dust grains also have a significant component of heating from giant stars.

• 7.

  The Spitzer/IRAC bands may contain as much unresolved stellar flux as resolved stellar flux for |ℓ| < 20°.

The models presented here, while simple, help us quantify the main factors that determine the large-scale Galactic mid-plane emission. We plan to carry out significant improvement to the models in the future, including a more systematic exploration of parameter space, a more realistic clumpy distribution of dust, and an improved treatment of the inner Galaxy, specifically relating to the molecular ring and the Galactic bar. In addition, the modeling of observations can be extended to spatial and wavelength regions outside that covered here. For instance, one could include observations for the outer Galaxy (GLIMPSE 360) or at longer wavelengths (Herschel HiGal). Once the data from the Wide-field Infrared Survey Explorer (WISE) are fully released, it will even be possible to model the all-sky observations from 3.5 μm to 100 μm when combined with the IRAS all-sky data.

Online material

Scripts used to produce the results Access here

Also available via GitHub

Acknowledgments

We wish to thank the referee for a careful review and for comments that helped improve this paper. We thank Bruce Draine and Aigen Li for useful discussions relating to PAH emission. This work is based in part on observations made with the Spitzer Space Telescope, which is operated by the Jet Propulsion Laboratory, California Institute of Technology under a contract with NASA. T.R. was supported by NASA through the Spitzer Space Telescope Fellowship Program, and B.W., E.C., M.M., and B.B. gratefuly acknowledge support from NASA through awards issued by JPL/Caltech (1368699 and 1367334). The models in this paper were run using Hyperion 0.9.0, and the scripts used to produce the results presented here are available as supplementary online material.

References

All Tables

All Figures

Fig. 1 Spectrum of the average surface brightness inside −65° < ℓ < 65° and −1° < b < 1° for the initial model (left) and the final modified model (right). The thick colored solid lines represent different components of the model: supergiants (purple), AGB stars (blue), giants (teal), main sequence stars (green), other stellar types (dark gray), dust emission from PAHs (red), VSGs (orange), and big dust grains (yellow), and scattered light (black). The major components are labeled in the left panel. The shaded gray area shows the total flux. The open black circles show the model fluxes after convolution with the transmission curves, and the filled black circles show the observations.

In the text

Fig. 2 Longitude and latitude surface brightness profiles in the range −65° < ℓ < 65° and −1° < b < 1° for the observations and for the initial model (Sect. 4.1). The thick black lines show the observed surface brightness, while the colored lines and the gray shaded area show the model surface brightness, with the same colors as used in Fig. 1. The “traditional” spiral arm tangencies are indicated in the top panel.

In the text

	Fig. 3 Longitude and latitude surface brightness profiles as for Fig. 2, but for the model including two major and two minor spiral arms instead of four main spiral arms (Sect. 4.2.1).
In the text

	Fig. 4 Longitude and latitude surface brightness profiles as for Fig. 3, with a deficit in the dust distribution in the central few kpc of the Galaxy (Sect. 4.2.2).
In the text

	Fig. 5 Longitude and latitude surface brightness profiles as for Fig. 4, with 50% more PAHs (Sect. 4.2.3).
In the text

	Fig. 6 Relative radial dependence of the dust density for the original model presented in Sect. 4.1 (gray) and the model with an inner dust hole from Sect. 4.2.2 (black).
In the text

	Fig. 7 Fractional contribution of various stellar populations and dust grain sizes to the total flux in the survey area as a function of wavelength. The colors are as in Fig. 1. Only the main contributing spectral types and dust types are labeled.
In the text

	Fig. 8 Fractional contribution to dust heating as a function of individual spectral type. The colors are as in Fig. 1, and the darker vertical lines indicate separations between individual spectral types. Only the main contributing spectral types are labeled.
In the text

	Fig. 9 Ratio of stellar emission measured in the GLIMPSE point source Archives to the total stellar flux predicted by the final model. This suggests that unresolved stellar flux contributes a large component of the total stellar flux inside |ℓ| < 20°.
In the text

	Fig. 10 Color composite view of the Galaxy model from an external viewpoint, viewed along the North Galactic Pole, for the initial model (left) and the final modified model (right). The Sun is located at (x, y) = (0, −8.5 kpc). The colors show IRAC 8.0 μm (red), IRAC 4.5 μm (green), and IRAC 3.6 μm (blue), all on a square-root intensity scale.
In the text