Dust-correlated centimetre-wave radiation from the M78 reflection nebula

Abstract

An anomalous radio continuum component at cm-wavelengths has been observed in various sources, including dark clouds. This continuum component represents a new property of the interstellar medium. In this work, we focus on one particular dark cloud, the bright reflection nebula M78. The main goal of this work is to investigate the cm-wave continuum emission in a prominent molecular cloud, nearby and with complementary observational data. We acquired Cosmic Background Imager (CBI) visibility data of M78 at 31 GHz with an angular resolution of ∼ 5.8 arcmin, and CBI2 data at an angular resolution of ∼ 4.2 arcmin. A morphological analysis was undertaken to search for possible correlations with templates that trace different emission mechanisms. Using data from Wilkinson Microwave Anisotropy Probe and the Rhodes/HartRAO 2326-MHz survey, we constructed the spectral energy distribution (SED) of M78 in a 45-arcmin circular aperture. We used results from the literature to constrain the physical conditions and the stellar content. The 5–31 GHz spectral index in flux density (α= 1.89 ± 0.15) is significantly different from optically thin free–free values. We also find closer morphological agreement with IR dust tracers than with free–free sources. Dust-correlated cm-wave emission that is not due to free–free is significant at small scales (CBI resolutions). However, a free–free background dominates at cm-wavelengths on large scales (∼ 1°). We correct for this uniform background by differencing against a set of reference fields. The differenced SED of M78 shows excess emission at 10–70 GHz over free–free and a modified blackbody, at 3.4σ. The excess is matched by the spinning dust model from Draine & Lazarian.

radiation mechanisms: general, ISM: clouds, dust, extinction, infrared: ISM, radio continuum: ISM

1 INTRODUCTION

Since 1996, experiments designed to measure the CMB anisotropy have reported an anomalous diffuse foreground in the range of 10–60 GHz. This diffuse emission is correlated with thermal emission from dust grains at 100 μm. The spectral index (considering S_ν∝ν^α) of the radio–IR correlated signal is α_radio/IR∼ 0 in the range 15–30 GHz, as for optically thin free–free (Kogut et al. 1996). But α_radio/IR∼− 0.85 between 20 and 40 GHz, for high-latitude cirrus (Davies et al. 2006). Additionally the absence of Hα emission concomitant to radio free–free emission would require an electron temperature T_e≥ 10⁶ to quench H i recombination lines (Leitch et al 1997). Another emission mechanism was presented by Draine & Lazarian (1998), who calculated that spinning interstellar dust grains produce rotational dipole emission in the range from 10 to 100 GHz, at levels comparable to the excess diffuse foreground observed over the free–free component.

Observations of specific targets may shed light on the anomalous foreground, whose existence in the diffuse interstellar medium (ISM) has been inferred statistically. Anomalous cm-wavelength radiation has been found in about a dozen clouds (Finkbeiner 2004; Watson et al. 2005; Casassus et al. 2006; Casassus et al. 2008; Dickinson et al. 2009; Scaife et al. 2009a,b; Vidal et al. 2010).

The search by Finkbeiner et al (2002) for anomalous cm-wave emission in nine dark clouds, resulted in one detection: LDN 1622.¹ The radio spectrum of LDN 1622 matches that of spinning dust emission (Finkbeiner 2004). Cosmic Background Imager (CBI; Padin et al. 2002) observations of LDN 1622 linked, on morphological grounds, the cm-wave emitters to the smallest interstellar dust grains (i.e. the very small grains, VSGs) – as required for spinning dust (Casassus et al. 2006). Interestingly, LDN 1622 is notable in being exposed to UV radiation from the adjacent hot stars in the Ori OB 1b association.

The prototypical dark cloud LDN 1688, close to the ρ Oph star, was recently found to be bright at cm-wavelengths (Casassus et al. 2008). The ρ Oph cloud is undergoing intermediate-mass star formation. UV radiation from its hottest young stars heats and dissociates exposed layers, but does not ionize hydrogen. Only faint radiation from the Rayleigh–Jeans tail of ∼10–100 K dust is expected at wavelengths longwards of ∼3 mm. Spinning dust comfortably explains the radio spectrum of ρ Oph W, as in LDN 1622. However, spinning dust encounters difficulties in explaining the cm-wave morphology of the ρ Oph cloud. Diffuse mid-IR (10–20 μm) emission from ρ Oph is interpreted as stochastic heating of VSGs (Bernard Boulanger & Puget 1993). Dust emissivities are proportional to both the dust grain density and the local UV energy density. Spinning dust was introduced to explain the 100-μm–1-cm correlation in the diffuse ISM, but does this correlation extend to smaller scales in denser environments?Casassus et al. (2008) report that the mid-IR intensity peaks in ρ Oph, from the circumstellar nebulae about Oph S 1 and Oph SR 3, have no 31-GHz counterparts. They find that the predicted spinning dust intensities towards S 1 are in excess over the observed values by a factor of >40. On the basis of Spitzer infrared spectrograph (IRS) spectroscopy, this discrepancy could marginally be explained by VSG depletion near S 1 and SR 3.

As part of an effort to understand the emission mechanisms in environments giving rise to the 31-GHz continuum, we have studied a dark cloud with extensive complementary data: the reflection nebula M78, with good constraints on its physical environment. Several observations are available in the literature that constrain the physical conditions in M78, making it a promising testbed of candidate emission mechanisms. As a reflection nebula M78 has relatively low levels of free–free emission but also contains stars bright enough to excite circumstellar dust grains. Reflection nebulae are regions of star formation, and they contain large amounts of dust.

In this report we refer to M78 as a group of reflection nebulae, all of them part of LDN 1630. This is a large dark cloud in the constellation of Orion, north from the belt, which also contains NGC 2023. LDN 1630 corresponds to the Orion B molecular cloud, which has recently been studied in molecular lines by (Buckle et al. 2009). The area that hereafter will be called M78 contains the following nebulae: NGC 2068 (sometimes also referred to as M 78), NGC 2071, NGC 2064 and NGC 2067 (see Fig. 1).

Figure 1

M78 area taken from Digitized Sky Survey (DSS) blue band survey. The largest and brightest nebula is NGC 2068/M78, the arc-shaped nebula to the west is NGC 2067. North-east from NGC 2068 is NGC 2071 and the faintest nebula, south-west from NGC 2068 is NGC 2064. The crosses indicate the positions of the brightest stars in the field. The dashed lines point to these stars and the nebulae are pointed to with solid lines. The coordinates are in right ascension and declination offsets from (⁠⁠, J2000), in degrees of arc.

Section 2 refers to the physical conditions and different phases of gas and dust found in M78. In Section 3 we will describe the observations used for this work. Section 3 also includes a discussion about the CBI flux loss. In Section 4 we present our results, derived from a comparison between the CBI data and various tracers of dust and free–free emission. In Section 5 we summarize our conclusions, quantifying the 31-GHz radiation from M78.

2 ENVIRONMENT AND PHYSICAL CONDITIONS

L1630 is embedded in the Orion Complex, at a distance of about 400 pc (Anthony-Twarog 1982). It is part of Orion Molecular Cloud II, among other star-forming regions and recently formed clusters. We also find Barnard's Loop in this region; it passes at ∼ 40 arcmin from the northern end of NGC 2071. As a reflection nebula M78 has low levels of free–free emission but is exposed to stars bright enough to excite IR emission from dust. Reflection nebulae are regions of star formation, so they contain young stars and dense molecular cores. As we can see in Fig. 2, although the large-scale radio emission is dominated by a uniform free–free background, the individual features of M78 become more pronounced with increasing frequency and angular resolution.

Figure 2

Large-scale structure around M78 in various frequencies. The solid line circle is the aperture where we integrated to extract flux densities. The crosses indicate the centre of the reference fields used for differentiation. The bright source at south-west from M78, common to each image, is NGC 2023. (a) At 2326 MHz, from Rhodes/HartRAO survey, we see clearly the arc to the north-east which corresponds to part of the Barnard's Loop, but there is no structure at the position of M78. (b) At 93.4 GHz, corresponding to WMAP W band, we see a faint feature corresponding to M78. (c) Image at 100 μm from IRIS, we see that the structure of M78 corresponds to that observed in (b). x- and y-axis show right-ascension and declination offsets from (⁠ J2000), in degrees of arc in degrees of arc.

The average H-nucleus density value in L1630, derived from CO low resolution maps, is n∼ 10³ cm⁻³, but towards NGC 2068 and NGC 2071, higher densities are found (Strom et al. 1975). Lada, Evans & Falgarone (1997) obtained, for these clumps, densities of log (n/cm⁻³) = 5.3 for NGC 2068 and log (n/cm⁻³) = 5.9 for NGC 2071. These values were derived from an excitation analysis of CS line fluxes for three transitions J = 2 → 1, J = 3 → 2 and J = 5 → 4. Buckle et al. (2009) have recently reported on a large-scale CO(3–2) isotopologue survey of Orion B, including the M78 area, and focusing on the star formation activity. Depending on the isotopologue, the molecular mass of NGC 2071 is 0.4 to 3.7 10³ M_⊙. NGC 2071 dominates NGC 2068 by mass.

M78 is still undergoing star formation (Flaherty & Muzerolle 2008). M78 contains two embedded clusters, in both NGC 2071 and NGC 2068. These clusters are considered to have an age of about 2 ± 1 Myr, and for a significant fraction of stars in the cluster there is evidence of accretion discs.

The illuminating stars for NGC 2068 and NGC 2071 have been described by Strom et al. (1975). NGC 2068 hosts the star HD38563N (B2 iii) and embedded in NGC 2071 is HDE290861 (B2-B3). Other early-type stars in NGC 2068 are HD28563S (B3-B5) and HD38563C (A0 ii). These stars appear to be very young, with ages around 10⁵ yr. For NGC 2064, the illuminating star is SSCV 111 (Strom et al. 1975), and its spectral type is B3 v (Chini et al. 1984).

The photo-dissociation-region (PDR) surrounding HD38563N was studied by Owl et al. (2002), finding a dust temperature of 41 K, a gas density of 5000 cm⁻³ and a gas temperature of 250 K. They also give G_○ = 2800 for the UV radiation incident on NGC 2068, where G_○ = 1 corresponds to the average interstellar radiation field of 1.6 × 10⁻³ erg cm⁻² s⁻¹ (Mathis, Mezger & Panagia 1983).

M78 has an associated C ii region. According to Pankonin & Walmsley (1978), the C ii region extends ∼ 20 arcmin along a line between NGC 2071 and NGC 2068. The electron temperature ranges from 20 to 50 K and the electron density is 0.2–1.0 cm⁻³ (Brown et al. 1975; Pankonin & Walmsley 1978). Both densities and temperatures were derived from the ratio of the power in two carbon recombination lines. According to Silverglate & Peter (1983), this C ii region is unusual because it does not present a sulfur line in its spectrum. This sulfur line is found towards all other known C ii regions. Another datum is that the carbon radio recombination line is the narrowest yet observed for reflection nebulae and has fainter brightness temperatures compared to other C ii regions.

Besides the C ii regions, a compact H ii region (∼ 1 arcmin) was reported by Matsakis et al. (1976) around the brightest star in NGC 2068, HD38563N. Matsakis et al. (1976) measured 0.1 Jy at 2.4 GHz. By assuming a temperature of ∼7000 K, they found an electron density of 40 cm⁻³. This region is responsible for most of the free–free emission of the nebula on 1–10 arcmin scales, although on 1° scale it is weaker than the background diffuse H ii region.

M78 also hosts H₂O masers (Campbell 1978). These masers are considered tracers of protostellar discs, given that they coincide with IR sources with associated IR excesses (i.e. circumstellar discs) and collimated outflows (e.g. as in class-II protostellar objects).

3 OBSERVATIONS

3.1 Cosmic Background Imager data

The CBI is described in Padin et al. (2002). The CBI is a planar interferometer array with 13 antennas, each 0.9 m in diameter, mounted on a 6-m tracking platform. The CBI receivers operate in 10 frequency channels covering 26–36 GHz. It is located in Llano de Chajnantor, Atacama, Chile. Cancellation of ground and Moon contamination was obtained by differencing with a reference field at the same declination but offset in hour angle by the duration of the on-source integration. We used an on-source integration time of 8 min, with a trailing reference field.

We used the CBI to observe M78 on 2004 December 20. The primary beam is at full width at half-maximum (FWHM) of the best-fitting Gaussian at the reference frequency of 31 GHz, and the synthesized beam is 5.94 × 5.70 arcmin² with uniform weights. We also observed M78 with CBI2, an upgrade where the 0.9-m CBI dishes were changed for 1.4-m dishes. In the case of CBI2, since the primary beam FWHM is 28.2 arcmin (this beam diverges from a Gaussian outside ∼FWHM/2), we mosaiced two pointings on NGC 2071 and NGC 2068. The uv coverage was different in the two fields, thus they have different synthesized beams (using uniform weights); for NGC 2068 we have 4.89 × 3.97 arcmin², while for NGC 2071 the values are 4.37 × 4.00 arcmin². The visibility coverage in the uv plane is displayed in Fig. 3.

Figure 3

The uv coverage for (a) CBI when taking M78, (b) for the observation of NGC 2068 with CBI2 and (c) for CBI2 in NGC 2071.

The data were reduced and edited using a special-purpose package (cbical, developed by T.J. Pearson). Flux calibration was performed using Tau A, whose fluxes are in turn calibrated against Jupiter (with a temperature of 146 ± 0.75 K; Hill et al. 2009). The flux calibrator is also used as the reference for an initial phase calibration (restricting to baselines longer than 2 m to avoid ground spill over).

The electronic phase stability of the CBI system is better than 10° over the time-scales between primary phase calibrations. However the pointing accuracy of the mount is approximately 0.5-arcmin rms. Due to the co-mounted nature of the CBI antennas, this appears as correlated phase errors in the visibilities and translates to an equivalent pointing uncertainty in the map. The pointing error can be measured, and removed, by interleaved observations of a nearby phase calibrator. For the CBI observations on 2004 December 20, the calibrator J0541−051 was observed, and the resulting pointing error correction was 62 arcsec. For the CBI2 observations, no sufficiently bright calibrator was observed. We note that these pointing errors and corrections are in any case a small fraction of the synthesized beam size (6 arcmin for CBI and 4.5 arcmin for CBI2) and do not significantly affect the positional matching of features in the CBI maps with data at other wavelengths described below (see Casassus et al. 2006, for tests on the pointing accuracy of a similar data set).

We fitted the CBI visibility data with three elliptical Gaussians, corresponding to NGC 2068, NGC 2071 and NGC 2064. The restored image is shown in Fig. 4(a). In the case of CBI2, the components are coincident with the four nebulae: NGC 2068, NGC 2071, NGC 2067 and NGC 2064. The Gaussian model convolved with the synthesized beam plus residuals of the CBI and CBI2 data are shown in Figs 4(a) and 5(a), respectively. Hereafter we refer to these images as ‘model-fitted’ images. The CBI2 image shows clearly the improved resolution of CBI2, and the increased flux loss of CBI2 due to reduced coverage at short spacings.

Figure 4

Comparison of CBI contours and various templates. The units are Jy beam⁻¹. The dashed circle corresponds to the CBI primary beam, the small circle outlined is the CBI synthesized beam. (a) Restored CBI image. Contours are taken at 95, 80, 65, 52, 43 and 35 per cent of the peak emission. White crosses indicate the brightest stars: from north to south, HDE290861, HD38563N, HD38563S and HD38563C. (b)–(e) Restorations of CBI-simulated visibilities on the IRIS images at 12, 25, 60 and 100 μm, with contours from CBI. (f) PMN image at 4.85 GHz with CBI contours. Coordinates follow from Fig. 2.

Figure 5

Comparison of model-fit and MEM reconstructions of CBI2 data for M78. Contours were taken at 95, 80, 65, 52, 43 and 35 per cent of the peak emission for each image. The dashed line corresponds to a noise level of twice the minimum noise. (a) Model-fitted reconstruction with contours. The maximum in the image is 0.32 Jy beam⁻¹. (b) MEM model with contours. The maximum is of 0.23 MJy sr⁻¹. (c) Restored MEM image, the model was convolved with the clean beam and residuals were added, the peak intensity in this image corresponds to 0.32 Jy beam⁻¹. All images are corrected for the primary beam (including the MEM model). Coordinates follow from Fig. 2.

For CBI2 we also implemented the Maximum Entropy Method (MEM) algorithm (as in Casassus et al. 2006). The MEM reconstruction provides a model sky consistent with the data, prior to convolution with the synthetic beam. We show a comparison of the two models in Fig. 5.

The expected thermal noise for the CBI image is 12.78 mJy beam⁻¹; for CBI2 we have two values, one for M78 of 9.94 mJy beam⁻¹ and another for NGC2071, of 11.91 mJy beam⁻¹. The achieved noise can be measured by taking the rms dispersion of a region within the half-power contour of the mosaiced primary beams, before correcting for it. The noise in Fig. 4(a), for CBI1, is 15 mJy beam⁻¹, and the noise in Fig. 5(a), for CBI2, is 20–30 mJy beam⁻¹ (the combination of two fields results in a varying noise map, we take 30 mJy beam⁻¹ as a conservative value).

3.2 Auxiliary data

The images that are used to make a morphological comparison were taken from Improved Reprocessing of the IRAS Survey (IRIS) (Miville-Deschenes & Lagache 2005), a re-processing of the IRAS survey (Wheelock et al. 1991) at 12, 25, 60 and 100 μm. We use these IR data to test for correlation between the 31-GHz emission and the emission from interstellar dust. The mid-IR emission at 12 μm is thought to be produced by stochastic heating of VSGs (e.g. Draine & Li 2001). The heat capacity of these grains is small enough so that the absorption of a single UV photon will increase its temperature to above 100 K, but their small cross-section makes it unlikely that they encounter many simultaneously. Thus VSGs are not in equilibrium with the radiation field. At longer IR wavelengths, (such as in the 60- and 100-μm IRIS bands), the emission is due to classical grains at T∼ 20–50 K.

The resolution of the IRIS data is not fine enough to allow a comparison with CBI2.² To trace small-scale structure, that could be compared with the improved resolution of CBI2, we used archival data acquired with the infrared array camera (IRAC) instrument aboard the Spitzer Space Telescope. This image is taken at 8 μm and traces VSGs and polycyclic aromatic hydrocarbons (PAHs), like IRIS at 12 μm, but with a finer resolution of 4 arcsec that allows proper point-source subtraction and comparison with CBI2 data. The IRAC image is shown in Fig. 6(a).

Figure 6

(a) Comparison between IRAC 8 μm (in grey-scale) and the CBI2 MEM restored image (in contours). The intensities are in MJy sr⁻¹ (b) In grey-scale we have the CBI2-simulated IRAC 8-μm image, scaled by the correlation slope with CBI2, against the restored image from CBI2, both reconstructed using model-fitting. (c) As in (b), but for simulations of CBI2 observations on the ratio of IRAC 8 μm and the UV field G_○. Coordinates follow from Fig. 2.

Point sources from IRAC were subtracted using a median filter. For the remaining sources, they were subtracted individually by fitting a Gaussian to the source and then subtracting it. CBI observations were simulated on both the IRIS and IRAC images using the mockcbi package (Pearson 2000, private communication), in order to have equivalent visibility sampling when comparing with CBI and CBI2 images. This will remove any uniform background from IRIS and IRAC. When simulating IRAC data, we embedded the image in a background taken from IRIS 12 μm, to extend the field of the IRAC 8-μm mosaic to at least three times the size of the primary beam. When embedding IRAC 8 μm in IRIS 12 μm, we fitted a plane to the difference of emission from the border of the images in order to remove most of the discontinuities and to remove the contribution from zodiacal light in the 8-μm image (IRIS images have only residual contribution from zodiacal light). The simulation is also useful to acquire visibility–plane data, which we used to make linear correlations to obtain a quantitative comparison between dust and 31-GHz emission. The simulations obtained will hereafter be referred to as CBI-simulated or CBI2-simulated images.

We also carried out a morphological comparison of the CBI2 data with Parkes-MIT-NRAO (PMN) (at 4.85 GHz) (Condon, Broderick & Seielstad 1991, with a resolution of 4.1 arcmin, see Section 4.6), NARO VLA Sky Survey (NVSS) (1.4 GHz) (Condon et al. 1998, with a resolution of 45 arcsec) and Southern H-Alpha Sky Survey Atlas (SHASSA) (at 6563 Å, and with a resolution of 48 arcsec; Gaustad et al. 2001). These templates trace free–free emission, which is present, at some level, at 30 GHz. It is necessary to quantify the intensity and spatial distribution of the free–free component in order to test for the existence of an excess at 30 GHz. These comparisons are discussed in Section 4.6.

For measurements of flux density, we have used images taken from Wilkinson Microwave Anisotropy Probe (WMAP), which give information about the continuum radiation on larger scales. The WMAP channels are at 22.8, 33.0, 40.7, 60.7 and 93.5 GHz, with resolutions of and ⁠, respectively (Limon et al. 2009). We also used data from the Rhodes/HartRAO survey at 2.326 GHz with a resolution of 20 arcmin (Jonas, Baart & Nicolson 1998).

3.3 Flux loss correction and spectral energy distribution

The CBI, being an interferometer with no total power measurements, is insensitive to an uniform background. The lack of total power also removes part of the extended emission, thus high-pass filtering the sky intensity field.

In order to estimate the level of flux loss, due to missing short spacings, we used the IRIS images as templates, shown in Fig. 4 to correlate with the CBI signal. We calculated the integrated fluxes in a 45-arcmin aperture from both raw and CBI-simulated data at each of the IRIS wavelengths, and compared them. This showed an average flux loss of 34 ± 8 per cent among the IRIS bands. This means that most of the flux in the IRIS images is within the compact objects corresponding to NGC 2071, NGC 2068 and NGC 2064.

Table 1 summarizes the spectral energy distribution (SED) of M78, for flux densities extracted in a circular aperture 45 arcmin in diameter, corresponding to the CBI primary beam diameter, and centred on RA and Dec. ⁠. Our choice of aperture contains the bulk of the M78 region, and matches the signal seen by CBI. Using the same aperture diameter, we also calculate the flux density in eight reference fields which are marked as crosses in Fig. 2. Subtracting the flux density in each reference field from the flux density centred on M78, we remove the uniform background, and estimate the error in this background level with the rms dispersion of the eight different background-subtracted flux densities.

Table 1

Flux densities extracted in a circular aperture 45 arcmin in diameter and centred on M78.

νa	Fbν	Fνc
2.326	7.08 ± 0.21	0.89 ± 0.48
22.8	4.40 ± 0.31	2.34 ± 0.66
33.0	4.34 ± 0.47	2.79 ± 0.81
40.7	4.37 ± 0.59	3.18 ± 0.98
60.7	5.07 ± 0.69	4.68 ± 1.44
93.5	20.20 ± 1.47	19.55 ± 3.60
31	3.96 ± 0.30d	2.73 ± 0.28d
3000	32418 ± 3220	29700 ± 3310

^aFrequency of observation in GHz.

^bFlux density in Jy.

^cFlux density in Jy after differencing with the reference fields.

^dThis value is obtained from IRIS 100-μm–31-GHz correlated emission, using the slope found from the visibility correlation.

^eOrigin of the data, for each frequency (in GHz): 2.326, Jonas et al. (1998); 22.8, 33.0, 40.7, 60.7, 93.5, Limon et al. (2009); 31, this work; 3000, Miville-Deschenes & Lagache (2005).

The error bars in Table 1 give 1σ, with ⁠, where

• σ_diff is the rms dispersion of the results from differentiations against the different reference fields,

• ⁠, where ΔΩ is the solid angle of one pixel, I_rms is the nominal thermal noise in specific intensity, when available, or the rms dispersion of a ∼ square box devoid of signal, N_beam is the number of pixels in a synthesized beam for each image, N_apert is the number of pixels in the aperture, in this case a diameter of 45 arcmin.

An additional term σ_cal could have been added in quadrature to the 1σ errors on a flux density F, σ_cal = 0.1 F, as a representation of systematic errors in calibration. However σ_diff dominates the error budget, and σ_cal is negligible when added in quadrature.

4 DISCUSSION

To characterize the emission from M78 that may not be explained as free–free and/or vibrational dust, i.e. the Rayleigh–Jeans tail of submm dust emission, we will perform a morphological and spectral analysis. First, we will summarize a comparison in the sky plane with various templates that trace different phases of dust and gas, to find qualitatively which corresponds better with the emission found at 31 GHz. Secondly, we will correlate the CBI data with the CBI-simulated templates, both in the visibility and sky planes. The visibility cross-correlations give information on the average properties over the CBI primary beam. Thirdly, another diagnostic is fitting an SED. This will give us the spectrum of M 78, which allows us to estimate free–free and vibrational dust contributions at 31 GHz, and infer the level of any excess.

4.1 Qualitative CBI–IR comparisons

In Figs 4(b)–(e) we compare the CBI and IRIS images. There is a qualitative correlation between CBI and IRIS, with some variations. At 25 μm, NGC 2071 is more prominent than NGC 2068, which contrasts with IRIS 12, 60 and 100 μm. For the other IRIS images it is difficult to tell by visual inspection which correlates best with CBI.

In Fig. 6, we compare the CBI2 MEM model in contours with the IRAC 8-μm image in grey-scale. This shows us that the spatial correlation obtained between IRIS 12 μm and CBI still holds on smaller scales. We also see a radio peak with no IR counterpart towards NGC 2067. In IRAC this nebula appears as diffuse emission, while in CBI2 this peak is more intense than the one observed towards NGC 2064. In Fig. 6(b), we show the reconstructed CBI2-simulated image from IRAC obtained by model-fitting compared with the CBI2 restored MEM image (Fig. 5). We scaled the CBI2-simulated visibilities by the cross-correlation slopes obtained in Section 4.2 (see Table 3) to bring them to the same range as the CBI2 data, and averaged 90 reconstructions, corresponding to different realization of Gaussian noise as given by the CBI2 visibility weights. This noise-simulation process applied to the CBI2-simulated IRAC visibilities, and subsequent averaging, avoids differences in convergence with the models fit to CBI2, as could arise from the varying visibility weights.

Table 3

Linear correlations results for CBI2 versus IRAC 8 μm.

Template	χ2/νa	rvisb	Cc	rskyd
NGC 2068	3.44	0.254	0.473 ± 0.017	0.739 ± 0.020
NGC 2068e	3.92	0.129	0.409 ± 0.015	0.677 ± 0.019
NGC 2071	3.99	0.205	0.347 ± 0.028	0.530 ± 0.034
NGC 2071e	4.08	0.180	0.402 ± 0.031	0.516 ± 0.036

^aReduced χ² of the visibility correlation, where ν is the number of degrees of freedom.

^bLinear correlation coefficient, obtained from the visibility correlation.

^cNormalized cross-correlation, from the sky plane correlation.

^dLinear correlation coefficient, obtained from the sky plane correlations.

^eComparisons made with UV field corrected IRAC images.

We see in Fig. 6 that emission from NGC 2064 is at 24 ± 1 per cent of that found for NGC 2068 in the CBI2-simulated 8-μm image, while in CBI2 this percentage is 70 ± 7 per cent. We also note that the radio emission peaks are not coincident with illuminating stars, particularly in NGC 2064 where the radio nebula is offset to the south-east of the IR counterpart. If the emission is due to spinning dust, and if the spinning dust emissivities are fairly independent of the radiation field (e.g. Ali-Haımoud Hirata & Dickinson 2009; Ysard & Verstraete 2009), we expect the ratio of mid-IR and CBI images to be proportional to the averaged UV intensity. So we do not expect a direct correspondence between the radio and the IR. For instance, in the vicinity of the brightest stars the mid-IR intensities will be enhanced relative to the spinning dust intensities. The UV field will also affect the grain size distribution, directly impacting the spinning dust spectrum. Mid-IR specific intensities from dust grains are expected to be proportional to the UV field, as well as to VSG column densities. This may explain the fact that NGC 2067 is fainter at 8 μm than at 31 GHz because this nebula does not contain stars that are bright enough to have a strong UV field, but still has enough dust column to radiate conspicuously in the radio.

We find that the morphological differences between IRAC 8-μm and the 31-GHz continuum cannot be explained only in terms of the proportionality of IR emission and the UV field. Following the procedure described in Section 4.3, we have estimated the UV field in M78, parametrized as G_○. Fig. 6(c) presents the simulation of CBI2 observations on the ratio of IRAC 8 μm to G_○. It is apparent from Fig. 6(c) that the radio/IR correlation worsens in M78 after division by G_○. The peak in NGC 2068 is shifted to the south when compared to the observed from the UV field uncorrected IRAC image and the one observed with CBI2. We also observe that emission corresponding to NGC 2067 appears clearly in Fig. 6(c), while in Figs 6(a) and (b) it is missing.

4.2 CBI–IR correlations

For a quantitative analysis, we performed linear correlations between CBI and IRIS data and between CBI2 and IRAC at 8 μm, both in the visibility and sky planes. We used different diagnostics in each case. On the visibility plane, we used reduced χ² and the linear correlation coefficient r, which we calculated using:

1

where N is the number of visibilities, V₃₁ are the visibilities at 31 GHz, V_t are the visibilities from the simulated template and the sums extend over every visibility. This last coefficient will be referred to as r_vis.

The sky plane correlations were carried out within circular areas, with radii of 45 arcmin (i.e. with diameters twice the CBI primary beam FWHM), and centred on each region of interest (i.e. M78 for CBI1 and NGC 2068, NGC 2071 for CBI2). These sky plane correlations were made with model-fitted images. For these correlations we also use the r coefficient, changing now visibilities for intensities and N for the number of pixels. This coefficient will be referred to as r_sky. We also use the normalized cross-correlation coefficient,

2

where I₃₁ and I_t are the intensities at 31 GHz and from the template respectively and α is defined by I₃₁=αI_t and determined by linear regression (this definition of C is as used in Casassus et al. 2006, except there is a typo in their equation 1, where α is missing). Here again, the sums are taken over all pixels falling within a circular aperture 90 arcmin in diameter. This choice of aperture is meant to include all of the signal seen by CBI.

To obtain statistics on r_sky and C, we averaged 1000 values of r_sky and C taken on different reconstructions of the CBI data, each with the addition of different realizations of Gaussian noise to the visibilities (given by the visibility weights). The values shown in Tables 2 and 3 correspond to the mean value of a normal distribution fitted to the data so obtained, and the errors correspond to 1σ.

Table 2

Linear correlations results for CBI versus IRIS, both for visibility and sky planes.

Template	χ2/νa	rbvis	Cc	rdsky
12 μm	1.40	0.776	0.683 ± 0.024	0.786 ± 0.019
12 μme	1.39	0.779	0.774 ± 0.025	0.786 ± 0.018
25 μm	1.59	0.768	0.644 ± 0.023	0.739 ± 0.021
25 μme	2.15	0.719	0.627 ± 0.022	0.717 ± 0.017
60 μm	1.58	0.755	0.651 ± 0.020	0.787 ± 0.017
60 μme	1.28	0.789	0.750 ± 0.024	0.819 ± 0.018
100 μm	1.39	0.777	0.699 ± 0.021	0.802 ± 0.019
100 μme	1.58	0.754	0.761 ± 0.027	0.761 ± 0.018

^aReduced χ², taken in the visibility plane.

ν is the number of degrees of freedom, i.e. twice the number of visibilities (counting real and imaginary parts).

^bLinear correlation coefficient from the visibility correlations.

^cNormalized cross-correlation, from the sky plane correlation.

^dLinear correlation coefficient, obtained from the sky plane correlations.

^eComparisons made with UV field corrected IRIS images.

Table 2 summarizes numerical results from the correlations between CBI and IRIS. As we can see the best results independently of the diagnostic (reduced χ², C, r_vis or r_sky) correspond to 12 and 100 μm, which both, given the uncertainties, correlate equally well with the CBI data. The worst correlations are obtained at 25 μm, except for the r_vis diagnostic, where it gives a value closer to one than the obtained result at 60 μm.

For correlations with the finer resolution data of CBI2 we used IRAC 8 μm as a template corresponding roughly to IRIS 12 μm. However, the two channels do not probe exactly the same medium, we have different widths and line contributions. In the case of 12 μm, the channel ranges from λ = 7 to 15.5 μ m (Δλ = 8.5 μm) and includes PAHs features at 7.7 and 11.3 μm and H₂ lines at 8.026, 9.665 and 12.3 μm. For 8 μm, the range is from λ = 6.5 to 9.5 μ m with contributions from the PAH feature at 7.7 μm and H₂ lines at 6.91 and 8.026 μm.

The correlation diagnostics worsen when testing with CBI2-simulated data from IRAC 8 μm, as summarized in Table 3. The degradation of the correlation when comparing CBI2 and IRAC 8 μm might be due to the intrinsic morphological differences pointed out in Section 4.1, which stand out with the finer resolution of CBI2 and IRAC, compared to the CBI and IRIS tests.

4.3 Stellar UV field and VSG column density map

Under the spinning dust hypothesis, we expect both the mid-IR and radio emission to be proportional to the column density of VSGs. We assume that mid-IR emission is due to stochastically heated VSGs. Thus the mid-IR emission will also be proportional to the UV field. We also assume that spinning dust is proportional to the column of VSGs, and is fairly independent of the local UV field.

Then to obtain a more accurate template of the VSG column we must divide the mid-IR images by an estimate of the UV field. We estimated the dust temperature T_d with the IRIS images at 60 and 100 μm. The resulting T_d map shown on Fig. 7(a) was obtained by fitting a modified blackbody to each pixel using an emissivity index for big grains of β = 2 (representative of the emissivity index in the far-IR). The T_d map can be converted to an equivalent radiation field using G_○ = (T_d/17.5 K)^β+4 (as in Ysard, Miville-Deschenes & Verstraete 2009). The result is shown in Fig. 7(b).

Figure 7

Proxies for the temperature and UV field in M78. (a) Dust temperature map obtained from IRIS 60/100, (b) implied G_○ map with a dust emissivity β = 2. (c) IRAC 8-μm filtered to remove point sources (i.e. stars), and resampled at IRIS resolutions. (d) Expected G_○ map from the stellar content of M78, and in the absence of extinction. The grey-scale has been chosen to match (b); G_○ diverges as the inverse-squared projected distance to the UV sources. Coordinates follow from Fig. 2.

We also fitted a modified blackbody to the integrated flux at 100 and 60 μm obtaining an average temperature of 32 K and a spectral index of 1.83, consistent with big grains. The spinning dust emission should arise from the VSG population that is traced by mid-IR emission. Calculating the expected emission from this modified blackbody at 12 μm we found that it is negligible (∼ 10⁻⁵ Jy), while the value extracted from the 12-μm image is ∼ 1200 Jy. Thus we expect a significant population of VSGs. This is also deduced from the 8-μm image which shows significant emission from the 7.6-μm PAH band.

The UV luminosities of HD38563N, HD38563S, HD38563C, HDE290861 and SSCV 111 are extracted from the ATLAS12 model atmospheres (Castelli 2005). The values we used for the stellar parameters are listed on Table 4. In connection with the discussion on C i continuum in Section 4.5, we include the Strömgren radius of the C ii region around each star using values of n_e = 1 cm⁻³ and T_e = 50 K, and recombination coefficients from Nahar & Pradhan (1997).

Table 4

Parameters for the brightest stars embedded within M78.

Name	Spectral type	Temperaturea	Surface gravityb	Luminosityc	C ii Strömgren radiusd
HD38563N	B2 iii	20000	−1	1.6 104	6.06
HDE290861	B2 v	22000	−0.5	6.6 103	4.84
HD38563S	B5 v	15000	−0.5	6.9 102	1.39
HD38563C	A0 ii	10000	−1.5	1.3 103	0.11
SSCV 111	B3 v	19000	−0.5	2.7 103	3.16

^aIn K; ^bIn log (g/g_⊙); ^cIn L_⊙; ^dIn pc using ⁠.

Knowledge of the stellar content can be used to estimate the circumstellar UV field, by diluting the stellar UV luminosities according to the inverse-square of the projected distance. This rough estimate of the G_○ map, shown in Fig. 7(d), assumes that the layers of M78 exposed to the stellar UV field do not overshadow. In other words, Fig. 7(d) assumes that the mid-IR emission stems from regions directly exposed to the radiation from the illuminating stars, without intervening extinction. The similarity of the maps in Figs 7(b) and (d) helps as a sanity check on Fig. 7(b).

With an estimate of the UV field in M78, as parametrized by G○, we can divide the IRAC 8-μm image by the G○ map to produce a VSG column density map. This estimative column density map, in arbitrary units, is shown on Fig. 7(c). The spike near the origin of coordinates in Fig. 7(c) corresponds to HD38563N.

4.4 CBI and CBI2 correlations with the VSG column density map

In Fig. 8(a), we see the Gaussian distributions fitted to the histograms of the r_sky coefficient for the different IRIS bands before performing the UV field correction. The two Gaussians fitted to the histograms at 12 and 60 μm are almost coincident, while at 100 μm we find the best correlation, r_sky = 0.802 ± 0.019. Finally at 25 μm we find the worst correlation result, as was expected because this is the only template where NGC 2071 appears brighter than NGC 2068.

Figure 8

(a) Fitted distributions to Monte Carlo noise simulation of the correlation coefficient r_sky between CBI and IRIS for M78. The solid lines corresponds to the 12-μm templates, the dashed lines to 25 μm, the dot–dashed lines to 60 μm and the dotted lines to 100 μm. The Gaussian parameters are listed in Table 2. (b) Distributions obtained from UV-field corrected templates. The corresponding parameters are also listed in Table 2.

In Table 2, we give the results for IRIS–CBI correlations, corrected by our estimate of the UV field. We see that the best correlation now corresponds to 60 μm, which improves from the uncorrected version, as reflected in the 3σ increase in r_vis, C and r_sky. For 12 μm the improvement is slight. In general we see that the values found for 25 and 100 μm suffer a worsening when compared with the UV-field uncorrected results. In the case of 25 μm, we see again the worst correlation results, as expected for here the brightest nebula is NGC 2071 and not NGC 2068. The Gaussian fit to the distribution of the r_sky coefficient is shown in Fig. 8, where we can compare cross-correlations with templates corrected and uncorrected by the UV field.

The UV field corrected results for IRAC–CBI2 correlations are given in Table 3. In Fig. 9 we show the distribution of the correlation coefficient r_sky obtained for each nebula, with and without the UV field correction. In NGC 2071 the correlation indicators are generally slightly lower after division by G_○, but the decrement is within the noise. In NGC 2068, all correlation indicators decrease markedly, with significant decrements in C and r_sky over 3σ.

Figure 9

Distribution obtained from Monte Carlo noise simulations of the correlation coefficient r_sky between CBI2 and IRAC for NGC 2068 and NGC 2071. The solid line histogram corresponds to the data obtained without the UV-field correlation, while the dashed line one is the corrected distribution. We also show Gaussian fits for each case, whose parameters are listed in Table 3.

The estimates of the VSG column density maps produced by dividing the IR images by the UV-field should be taken with some caution. Apart from the uncertainties in the UV-field estimates, there is also the possibility that the grain population may be influenced by the UV field, in particular, in its size distribution. To make a more complete comparison is beyond the scope of this paper.

4.5 Expected C i continuum

Estimates of the C i continuum level can be made under the assumption of ionization balance for the stellar content of M78 with T_e ranging from 20 to 50 K (see Section 2). If the stellar UV³ radiation is absorbed by carbon in the M78 PDRs, and the kinetic temperature is T_e = 50 K, we expect a flux density of F_ν = 2.7 Jy at 31 GHz for the whole of M78. Details are given in Casassus et al. (2008; note that there is a typo in their equation 1, for a missing D² in the denominator). We use the recombination coefficients tabulated in Nahar & Pradhan (1997).

However, Casassus et al. (2008) overlooked the competition of dust and carbon for the carbon-ionising continuum, which may be estimated by comparing opacities of dust, τ_d, and carbon, τ_C, to C-ionising radiation. We take a constant photoionization cross-section for carbon of 1.6 × 10⁻¹⁷ cm² between 11.2 and 13.6 eV, and use the relation A_V=N_H/2 10²¹ (Bohlin, Savage & Drake 1978, in CGS units). We find that τ_d/τ_C≈A_λ/(A_V 3.2) ≈ 2.0 for λ≈ 1000 Å, for a standard R_V = 3.1, and using the extinction curves given in Draine (2003). About one-third of all carbon-ionising photons is absorbed by carbon, although this result is sensitive to the exact value for R_V.

The expected C i continuum levels at 31 GHz, in a dust-free nebula, are F_ν = 2.1, 2.7 and 3.3 Jy, for T_e = 20, 50 and 100 K. But taking into account dust absorption of the carbon-ionising radiation, these values are reduced by 1/3, to not more than 1 Jy. This, however, is not negligible compared to the observed level of anomalous emission calculated from the differenced data set, of 2.79 ± 0.81 Jy.

For a uniform slab nebula, seen face-on as a disc with a radius of ⁠, and at a distance of 400 pc, the corresponding electron density at T_e = 50 K is 15 cm⁻³, which is n_H = 1.5 × 10⁵ cm⁻³ if the carbon in the PDR gas is all singly ionized. Ionization balance indicates that the carbon continuum could reach observable levels, given the observed physical conditions in M78 (Section 2). However, the contribution of optically thin carbon continuum at 5 GHz is probably negligible since the 5–31 GHz index is α³¹₅ > 0 (we do not see the nebulae in the PMN maps). To reach unit opacity at 5 GHz, we require electron densities of order ∼ 100 cm⁻³, corresponding to H-nucleus densities of ∼ 10⁶ cm⁻³ if all of carbon is ionized, with a filling factor ∼ 10⁻⁴.

Alternatively, the lack of detectable optically thin C i continuum at 5 GHz could be used to place constraints on the dust extinction law in M78. This requires full-blown PDR models, and is beyond the scope of this work.

4.6 Free–free specific intensities

In Fig. 10 we compare the CBI contours and the free–free tracers on small scales. The 4.85-GHz PMN image traces H i free–free. No 4.85-GHz emission is evident at the location of NGC 2071. The 4.85-GHz intensity peak, located in NGC 2068, is 0.073 ± 0.011 Jy beam⁻¹. If we assume an electron temperature of T_e = 7000 K, then the spectral index of optically thin free–free emission between 5 and 31 GHz is α=− 0.12 (e.g. Dickinson, Davies & Davis 2003; Casassus et al. 2007; Wilson, Rohlfs & Httemeister 2009), and the peak free–free intensity should be ∼ 0.055 Jy beam⁻¹ at 31 GHz. Hα, also a tracer of H i free–free, bears similar properties as the PMN map – no Hα is seen in NGC 2071, albeit faint wings from Barnard's loop.

Figure 10

CBI contours over free–free tracer images. We find that the peaks coincide, but there is no further correlation between CBI and any of the templates. (a) CBI contours against NVSS image. The crosses indicate the positions of the exciting stars. (b) SHASSA image with CBI contours. In this image we can appreciate the extended Hα emission that surrounds M 78. (c) A zoomed image from PMN with CBI contours. Coordinates follow from Fig. 2.

In each template, i.e. NVSS, SHASSA and PMN, there is compact free–free emission coincident with the brightest star embedded in NGC 2068. It is the compact H ii region reported by Matsakis et al. (1976). We see that the peak is coincident with HD38563N for both CBI and PMN images.

In order to make a proper comparison between PMN and CBI/CBI2 specific intensities, we have to smooth them to a common resolution. Any high-pass filtering from PMN will not be significant in this comparison since, as explained in Condon et al. (1994), the PMN survey filters out scales above 20 arcmin, and since the shortest CBI2 baseline is on 20-arcmin scales, CBI2 and PMN have fairly well-matched beams. Being an interferometer, the CBI filters on all scales above its point spread function (∼8 arcmin), so we expect the comparison with PMN to give conservative values in terms of the detection of an excess over optically thin free–free.

The specific intensity units of the PMN data retrievable from SkyView⁴ are not defined precisely. It is stated in the header of the PMN images that 200 000 beams are about 1 sr (so a beam of 7.2 arcmin). To obtain the exact size of the PMN beam, and to cross-check on the Jy beam⁻¹ units, we fitted a Gaussian to 3C273, which is at a declination similar to that of M78 (RA Dec. ⁠) so both images were obtained using similar beams. Although an elliptical Gaussian fit to 3C273 has a minor axis of 3.8-arcmin FWHM, trials on neighbouring and fainter radio point sources led to values around 3.4 arcmin. We assume the original PMN map we are using has a resolution of 3.7 arcmin. Given the resolution of the CBI maps used here, of ∼ 5.8 arcmin, we smoothed the resolution of PMN to that of CBI by convolving with a Gaussian whose FWHM is FWHM ²_ker = FWHM²_CBI− FWHM²_PMN. We also use this value for the beam units in the PMN image.

After smoothing, and scaling Jy beam⁻¹ units to the CBI beam, we find that the peak intensity in the M78 PMN image is of 0.015 ± 0.004 Jy beam⁻¹, while the peak found in the CBI image is of 0.492 ± 0.031 Jy beam⁻¹. Using a spectral index α=− 0.12 and assuming pure free–free contribution at 4.85 GHz, we find an expected emission of 0.012 ± 0.003 Jy beam⁻¹ at 31 GHz. The difference between the actual emission detected at 31 GHz and the optically thin free–free level expected from 4.85 GHz is of 0.480 ± 0.031 Jy beam⁻¹. The spectral index between 4.85 and 31 GHz is α = 1.89 ± 0.15.

From archival fluxes from PMN and NVSS, we find that the peak specific intensity between 1.4 and 4.85 GHz has a spectral index of 0.00 ± 0.23. This is an upper limit since NVSS filters out large scales (NVSS has more flux loss than PMN). We can none the less conclude that the free–free emission seen at 1.4 and 4.85 GHz is optically thin.

4.7 Low-frequency spectrum of NGC 2071

Under the hypothesis that the absence of NGC2071 at 5 GHz is due to an optically thick spectrum, with the observed 31-GHz flux, of ∼0.1 Jy we expect 16 mJy beam⁻¹ at 5 GHz. This is just barely at 2σ in PMN, so could be missed. Also there is hot dust in NGC 2071, as seen in the 25 μm image, perhaps due to classical grains heated by a hot star.

For optically thick emission, the expected flux density at 5 GHz from NGC 2071, 16 mJy beam⁻¹, corresponds to a diameter of 0.05 arcsec. Such a compact and optically thick source would stand out in Very Large Array (VLA) observations at frequencies higher than 5 GHz. Unfortunately we found no VLA observations that cover the 31-GHz centroid of NGC 2071.

However, the bulk of the 31-GHz flux density seen towards NGC 2071 is probably not due to an UCHII region. The MEM model of the CBI2 visibilities finds structure in NGC 2071, so that its extension is probably closer to 4 arcmin. If opacity is 1 at 31 GHz (so that the 5–31 GHz index is +2), then a distance of 400 pc and an angular extension of 4 arcmin, or a linear depth of 0.5 pc, gives an electron density of 10⁵ cm⁻³ if T_e = 10⁴ K. Such densities are found only in the most compact UCHII regions, under 0.01 pc in diameter. In fact, no UCHII larger than 0.1 pc are tabulated in Wood & Churchwell (1989), in their table 16. More recently, Murphy et al. (2010) have listed UCHII region properties in their table 3, where an upper limit of 0.1 pc can be found on the linear sizes of regions with N_e≈ 10⁵ cm⁻³. The ionized mass for NGC 2071 would have to be around 125 M_⊙, while the largest UCHII regions reach only 0.5 M_⊙.

4.8 Radio stars

There are several protostars in this region, which contribute to the total emission at 31 GHz. In most cases, the clumps associated with stellar formation do not have a related radio continuum source (to a limit of 0.1 mJy at 8.4 GHz; Gibb 1999). Table 2 in Gibb (1999) gives information on the known 8.4-GHz radio continuum sources within M78. We also included a source found in the VLA archive not included in Gibb (1999); this sources is roughly coincident with NGC 2071 IRS 1.

The radio continuum sources will contribute to some extent to the emission found at 31 GHz. To derive this contribution we use the formula S_ν∝ν^α, with α being the spectral index. Collimated, ionized stellar winds (sources for radio emission from young stellar objects) have spectral indices constrained by −0.1 < α < 2, with the highest values corresponding to totally opaque sources (e.g. Reynolds 1986). Using this and considering the highest possible spectral index, we find that the combined contribution from these radio continuum sources is 237.98 ± 11.22 mJy, far below the observed emission at 31 GHz, which is of 2.73 ± 0.28 Jy. Thus we neglect any contribution from radio stars to the 31-GHz flux density in a 45-arcmin aperture.

4.9 Spectral energy distribution fits

In the previous section, we used specific intensity (surface brightness) to make the comparison between 31 GHz and free–free templates. Here we use the integrated flux densities to build an SED. For this we use the values shown in Table 1 for differenced fluxes, unless otherwise indicated (see Section 3.3). We will fit a model of three components of emission: free–free, vibrational dust and spinning dust.

In Fig. 11, we see a fit to the data with two emission mechanisms, not including spinning dust emission. One is free–free emission for which we used α=− 0.12, adjusted the datum at 2.326 GHz. We also fitted a modified blackbody to the IRIS data points at 100 and 60 μm. This fit gives us a warm dust temperature of 32 K and a spectral index β = 1.83[we take the spectral index from S_ν∝ν^βB_ν(T)]. To check this value, we also fitted a modified blackbody to the IRIS images that were CBI-simulated, and found consistent values of T = 32 K and β = 1.82. We also fitted a Rayleigh–Jeans tail to the points of WMAP which accounts for cold dust, assuming that there is no contribution from spinning dust. The characteristic values for the spectral index β (defined by S_ν∝ν^2+β) of ISM dust in the Rayleigh–Jeans regime range from 1 < β < 2 (Shetty et al. 2009). To measure the confidence level of our fit we used a reduced χ² test. For this fit, we used a value of β = 0.6, lower than characteristic values found in the ISM, but also the highest value for which the confidence level is ≥ 5 per cent. We see from this fit that a model which uses only Rayleigh–Jeans from cold dust is unable to explain the emission from M78 at 95 per cent confidence.

Figure 11

M78 SED and model components. The spectrum of the radio flux in a 45-arcmin aperture is modelled using two components made up from free–free emission (dotted line) and a Rayleigh–Jeans tail fitted to the WMAP data (dashed line). Diamonds correspond to flux densities measured by WMAP, the triangle is the flux density obtained by Rhode/HartRAO and the asterisk is the flux density from CBI–IRIS 100-μm correlated emission. The solid line is the sum of the three emission components considered.

We used the spinning dust model from Draine & Lazarian (1998) for dark clouds, which assumes n_H = 10⁴ cm⁻³ and T_d = 10 K. We fitted this model to the flux densities from differentiated WMAP data, except the W band (see Table 1). For a given spinning dust emissivity per H-nucleus (which depends on the physical environment), the spinning dust intensities are determined by the column density of hydrogen nuclei, N_H. We fit the WMAP SED by varying N_H, and obtain a value of (7.99 ± 2.34) 10²¹ cm⁻² for N_H=fn_HL, where f is the filling factor and L is the nebular depth. The uncertainties on N_H take into account uncertainties on the level of free–free and the level of residual Rayleigh–Jeans dust emission at 31 GHz. We see that the cm-wave excess, parametrized by the spinning dust model, is significant at 〈N_H〉 = 3.4 σ(N_H).

The value for f can be constrained using the densities assumed by the spinning dust model, and those given by Lada et al. (1997) and Strom et al. (1975). We assume that n_H ranges from 10³ to 10⁶ cm⁻³. We then consider the nebula as a sphere, with characteristic size of 27 arcmin which is the angular distance from the north-east extreme NGC 2071 to the south-west end of NGC 2064. When we consider a distance of 400 pc, we see that N_H would range from 10²² to 10²⁵ cm⁻², which constrains −3.25 < log (f) < 0.01. When using the extinction maps from Schlegel, Finkbeiner & Davis (1998), we find H-nucleus column densities within a smaller range, with N_H = (2.28 ± 1.40) × 10²³. For this we considered a constant ratio A_V/E(B−V) = 5.5 and the relation N_H = 10²²A_V/5.3 cm⁻² (Kim & Martin 1996). This new range in H-nucleus column densities gives us the value log (f) =− 1.46 ± 0.30.

Fig. 12 shows the SED calculated from WMAP, Rhodes/HartRAO and CBI, this last one extracted from IRIS 100 μm using the coefficient obtained from the visibility–plane correlation (see Section 3.3). The free–free emission was derived from the datum at 2.326 GHz using a spectral index α=− 0.12. We assume that the excess in WMAP at 94 GHz from the modified blackbody fitted using IRIS at 100 and 60 μm traces the submm cool dust grey-body. We fixed to this datum point the Rayleigh–Jeans tail of vibrational dust using a characteristic β = 1.8, which lies within the range established by Shetty et al. (2009).

Figure 12

M78 SED and model components. The spectrum of the radio flux in a 45-arcmin aperture is modelled using three components made up from the vibrational blackbody emission of classical dust grains, and that from free–free and spinning dust. Diamonds correspond to flux densities measured by WMAP, the triangle is the flux density obtained by Rhodes/HartRAO and the asterisk is the flux density from CBI–IRIS 100-μm correlated emission. The dotted line is the free–free spectrum derived from 2.326 GHz, the dashed line corresponds to the submm dust tail, traced at 94 GHz and the dot–dashed line represents a spinning dust model. The solid line is the sum of the three emission components considered.

We assume that spinning dust is traced by dust-correlated radio emission, using the IRIS 100-μm template, so that the 31-GHz emission is given by F_{31 GHz}=a F_{100 μm}, where a = (0.092 ± 0.001) × 10⁻³ is the slope of the correlation in the visibility plane. If all visibilities were weighted equally, this would be equivalent to correcting the CBI flux density for flux loss.

If we again use 2.326 GHz to trace free–free emission, but this time we differentiate to remove the diffuse background, we find that the free–free level is 0.65 ± 0.35 Jy at 31 GHz. Compared with the value from dust-correlated emission at 31 GHz, which is of 2.73 ± 0.28 Jy, we have an excess of 2.08 ± 0.45 Jy over the expected level of optically thin H i free–free. The emission at 31 GHz is not dominated by free–free, confirming the results from Section 4.5.

Part of the Rayleigh–Jeans tail of submm dust emission could reach 31 GHz. Fixing the level of the Rayleigh–Jeans flux density to the point at 93.5 GHz, with β = 1.6 ± 0.3 (which extends conservatively to shallow values), we expect 0.40 ± 0.19 Jy at 31 GHz. We see that the submm dust cannot account for the observed 100 μm correlated flux density at 31 GHz of 2.73 ± 0.28 Jy (close to the 33-GHz flux density from WMAP 33 GHz, of 2.79 ± 0.81). When we consider the addition of vibrational dust emission and free–free emission, as previously derived, we find a 31-GHz excess of F_{31 GHz} = 1.68 ± 0.49 Jy, significant at 3.4σ. Comparing this value with the expected 2.73 ± 0.28 Jy, we can say that at least 62 ± 19 per cent of the observed 31-GHz emission cannot be explained by H i free–free and/or vibrational dust emission.

The addition of a spinning dust component improves significantly the confidence level (as derived from a reduced χ² test) for the SED fit, when compared with models that only consider the Rayleigh–Jeans tail from submm dust to account for the spectra at 20–60 GHz. The only free parameter for a spinning dust component is fN(H), whose value of (7.99 ± 2.34) 10²¹ cm⁻² differs significantly from zero. From this we can say that spinning dust is a significant component in our SED for M78 at angular scales of ∼ 8 arcmin.

4.10 Comparison with other radio–IR proportionality coefficients

For comparison with previous work on the radio–IR correlations in other objects, we converted our correlation slopes to T_{31 GHz}/I_{100 μm}, i.e. μK (MJy sr⁻¹)⁻¹ instead of dimensionless units.⁵ These values are summarized in Table 5. We found that the highest slopes correspond to the fainter nebulae, NGC 2064 and NGC 2067, with values of 7.69 ± 0.96 and 7.21 ± 1.72 μK (MJy sr⁻¹)⁻¹, while for NGC 2071 and NGC 2068 we found 4.51 ± 0.14 and 4.29 ± 0.25 μK (MJy sr⁻¹)⁻¹. We see that all these values are within the range constrained by 3.3 ± 1.7 μK (MJy sr⁻¹)⁻¹, derived by Dickinson et al. (2007) from 6 H ii regions, and the value of 11.2 ± 5.0 μK (MJy sr⁻¹)⁻¹ for cool dust obtained from 15 regions by Davies et al. (2006). Regions with higher temperatures tend to have lower slopes. The radio-correlated dust emission stems from hotter dust in H ii regions than in M78. Also we see that the correlation slopes in NGC 2071 and NGC 2068 are lower than in NGC 2067 and NGC 2064, as expected since the exciting stars in NGC 2071 and NGC 2068 are hotter.

Table 5

Radio–IR correlation slopes for individual clouds in M78.

Source	μK (MJy sr−1)−1
NGC 2068	4.29 ± 0.25
NGC 2071	4.51 ± 0.14
NGC 2067	7.21 ± 1.72
NGC 2064	7.69 ± 0.96
6 H ii regionsa	3.3 ± 1.7
15 Cool dust regionsb	11.2 ± 5.0

^aDickinson et al. (2007); ^bDavies et al. (2006).

5 CONCLUSIONS

We have found evidence of dust-correlated emission towards M78 at a frequency of 31 GHz. Its spatial distribution is coincident with that from dust, especially as traced at 12 and 100 μm, for which we found similar linear correlation coefficient in the visibility plane. The morphology at 31 GHz is qualitatively inconsistent with pure optically thin free–free emission, as can be traced by PMN, NVSS and SHASSA, where we only find significant emission around HD38563N.

The radio–IR correlations worsen at finer angular resolutions. The CBI2 and IRAC 8-μm correlation tests worsen relative to those found from CBI and IRIS. The correlation tests further degrade when dividing the IR images by the UV field. Thus the radio/IR differences cannot be explained as being due to modulation by the UV field.

The specific intensity levels seen by CBI at 31 GHz are higher than those expected for thin free–free emission derived from PMN data (smoothed to the CBI resolution). This 31-GHz excess is 0.480 ± 0.031 Jy beam⁻¹. However, the flux densities measured in a 45-arcmin circular aperture are dominated by diffuse free–free emission, as suggested by flux measurements with WMAP and Rhodes/HartRAO. When we differentiate against a reference field within the Orion B complex, we find an excess flux density at 31 GHz of 2.08 ± 0.52 Jy over the expected free–free in a 45 arcmin circular aperture. We need high-resolution data to subtract more accurately the background free–free emission and confirm the values found in this work.

When fitting a Rayleigh–Jeans law to WMAP differentiated fluxes, we derive a spectral index β = 0.6. This value is lower than those typically found in the ISM, where 1 < β < 2, but even allowing for this value of β the confidence level of this model spectrum, given the observed SED, is only of 5 per cent. We can neglect the contribution of a Rayleigh–Jeans tail from submm dust at cm-wavelengths.

The addition of a spinning dust component to the SED fit results in a significant improvement in confidence level. The only free parameter in our fit, the column of H-nuclei which scales the model spinning dust emissivities per H-nucleus, differs significantly from zero, thus confirming that the higher confidence level achieved is due to a non-negligible component. This is also supported by the fact that the combined contributions of Rayleigh–Jeans and free–free at 31 GHz yield a total flux that is 1.68 ± 0.49 Jy below the emission detected by CBI. A spinning-dust-like component is needed in order to explain the observed emission at cm-wavelengths in M78.

1

And also another detection, which was later refuted by Dickinson et al. (2006).

2

As the IRIS resolution, of 4 arcmin, approaches that of CBI2, the simulation of CBI2 observations on the IRIS images will yield coarser resolutions than real CBI2 maps.

3

The total number of carbon ionizing photons is 6.33 10⁴⁷ s⁻¹.

4

http://skyview.gsfc.nasa.gov/

5

The radio–IR correlation slopes in units of μK (MJy sr⁻¹)⁻¹ are sometimes referred to as ‘emissivities’.

Thorough comments from an anonymous referee improved this manuscript. We thank Evelyne Roueff for interesting discussions on PDRs. SC acknowledges support from a Marie Curie International Incoming Fellowship (REA-236176), from FONDECYT grant 1100221, and from the Chilean Center for Astrophysics FONDAP 15010003. CD acknowledges an STFC Advanced Fellowship, and an ERC grant under the FP7. This work has been supported in part by the Strategic Alliance for the Implementation of New Technologies (SAINT – see http://www.astro.caltech.edu/chajnantor/saint/index.html). This work has been carried out within the framework of a NASA/ADP ROSES-2009 grant, no. 09-ADP09-0059. The CBI observations were made possible, thanks to the engineering team: C. Achermann, R. Bustos, C. Jara, N. Oyarce, R. Reeves, M. Shepherd and C. Verdugo. We acknowledge the use of the Legacy Archive for Microwave Background Data Analysis (LAMBDA). Support for LAMBDA is provided by the NASA Office of Space Science. We used data from the Southern H-Alpha Sky Survey Atlas (SHASSA), which is supported by the National Science Foundation.