A centimetre-wave excess over free–free emission in planetary nebulae

Abstract

We report a centimetre-wave (cm-wave, 5–31 GHz) excess over free–free emission in planetary nebulae (PNe). Accurate 31- and 250-GHz measurements show that the 31-GHz flux densities in our sample are systematically higher than the level of optically thin free–free continuum extrapolated from 250-GHz. The 31-GHz excess is observed, within one standard deviation, in all 18 PNe with reliable 31- and 250-GHz data, and is significant in nine PNe. The only exception is the peculiar object M2-9, whose radio spectrum is that of an optically thick stellar wind. On average the fraction of non-free–free emission represents 51 per cent of the total flux density at 31 GHz, with a scatter of 11 per cent. The average 31–250 GHz spectral index of our sample is 〈α²⁵⁰₃₁〉=−0.43 ± 0.03 (in flux density, with a scatter of 0.14). The 31–250 GHz drop is reminiscent of the anomalous foreground observed in the diffuse interstellar medium (ISM) by cosmic microwave background (CMB) anisotropy experiments. The 5–31 GHz spectral indices are consistent with both flat spectra and spinning dust emissivities, given the 10 per cent calibration uncertainty of the comparison 5-GHz data. However, a detailed study of the objects with the largest cm-excess, including the low-frequency data available in the literature, shows that present spinning dust models alone cannot explain the cm-wave excess in PNe. Although we have no definitive interpretation of our data, the least implausible explanation involves a synchrotron component absorbed by a cold nebular screen. We give flux densities for 37 objects at 31 GHz, and for 26 objects at 250 GHz.

radiation mechanisms: general, planetary nebulae: general, radio continuum: ISM, submillimetre

1 INTRODUCTION

An increasing amount of evidence supports the existence of a new continuum emission mechanism in the diffuse interstellar medium (ISM) at 10–30 GHz, other than free–free, synchrotron or a hypothetical Rayleigh–Jeans tail of cold dust grains (Leitch et al. 1997; Finkbeiner 2004; Watson et al. 2005). As proposed by Draine & Lazarian (1998), a promising candidate emission mechanism is electric dipole radiation from spinning very small grains (VSGs), or spinning dust. Cosmic Background Imager (CBI) observations of the dark cloud LDN 1622 linked, on morphological grounds, the cm-wave emitters to the VSGs (Casassus et al. 2006).

The SED of the Helix nebula [NGC 7293, an evolved planetary nebula (PN)] also exhibits a 30-GHz excess over the free–free continuum expected from 250-GHz measurements (Casassus et al. 2004). Resolved 31-GHz images of the Helix, obtained with the CBI, provided morphological evidence for a new radio continuum component.

PNe have traditionally served as test laboratories for nebular astrophysics. Their spectra are the product of the photoionization by a single exciting star of its surrounding point-symmetric nebula, in a characteristic radial and homologous expansion. PNe are archetypical free–free sources; the interpretation of their <14-GHz continuum in terms of free–free emission is paradigmatic (e.g. Condon & Kaplan 1998; Siódmiak & Tylenda 2001). Yet inconsistencies with the free–free paradigm exist in the literature, although not recognized as such. Hoare, Roche & Clegg (1992), in an investigation of the cm- to mm-wave SED of PNe, ascribed deviations from free–free spectra to instrumental uncertainties. Garay, Gathier & Rodriguez (1989) compared continuum and radio-recombination line maps, finding systematic differences interpreted as non-local thermodynamic equilibrium (LTE) amplification or large T_e variations.

Here we report on a survey of the cm- and mm-wave continuum in PNe acquired with the CBI (Padin et al. 2002) and with the SEST IMaging Bolometer Array (SIMBA; Nyman et al. 2001), at the Swedish–ESO Submillimetre Telescope (SEST). The PN sample is mostly selected from Casassus et al. (2001), and represents the population of compact and mid-IR-bright PNe (it is not a complete sample). Section 2 discusses PN electron temperatures and relevant free–free spectral indices. Sections 3 and 4 present the CBI and SIMBA data, respectively. Section 5 compares the two data sets, which are in turn compared to previous measurements in Section 6. Section 7 interprets our findings in terms of possible emission mechanisms. Section 8 analyses the SEDs of individual PNe. Section 9 concludes. All the data analysis was carried out using Perl Data Language (http://www.pdl.org), unless otherwise stated. The error bars on all figures are one standard deviation (±1σ; the error bars have a total length of 2σ).

2 ADOPTED PLANETARY NEBULA ELECTRON TEMPERATURES AND FREE–FREE SPECTRAL INDICES

The radio-frequency free–free spectral index is roughly constant at α∼−0.15 in flux density, F(ν) =F(ν_○) (ν/ν_○)^α. Yet α has a residual dependence on frequency, electron temperature T_e and helium contribution. Since our aim is to perform a precision analysis of the PN SEDs, we calculated accurate free–free indices for the frequencies relevant to this work. The values we obtained are exact for non-relativistic electrons (i.e. T_e < 10⁸ K). Full details are given in Appendix A.

The electron temperature in PNe has been the matter of significant debate. Peimbert (1971) first observed that optical recombination lines (ORLs) and collisionally excited lines (CELs) give broadly different T_e values. Since then, several new T_e diagnostics have been reported, all pointing towards the conclusion that ORLs and continuum diagnostics give systematically lower T_e values than [O iii]λ4363/λ4959, by a factor of ∼2 (and up to 10).

The current view, as summarized by Liu (2006), is that of a bi-abundance model in which PNe (and to a lesser extent H ii regions) contain dense, hydrogen deficient (with metallicity Z∼1) and cool inclusions accounting for most of the ORL emission. In the bi-abundance model the bulk of the nebula, by mass, is best represented by the CEL diagnostics, with T_e∼10 000–15 000 K. In the case of NGC 6153, the free–free luminosity of the ORL emitting gas is ∼1/10 that of the CEL emitting gas, and the hydrogen mass in the ORL gas is only ∼1–2 per cent that of the CEL gas (Mike Barlow, private communication; Liu et al. 2000, their model IH3). Assuming the ORL T_e is ∼1/2 the CEL T_e, an average T_e, weighted by free–free luminosity gives, in the case of NGC 6153, 〈T_e〉= 0.95 T_e(CEL).

Our strategy was to adopt a single conservative value for the whole PN sample. We thus adopted T_e= 7000 K, a value that approaches a canonical T_e for H ii regions (e.g. Dickinson, Davies & Davis 2003), but is unrealistically cold for PNe. This choice minimizes free–free deviations in the data set, as cooler T_e give steeper 31–250 GHz spectral indices. Since the emissivity of the He⁺ continuum gives somewhat steeper indices than for H, we also included a contribution from the He⁺ continuum, assuming a large but conservative 10 per cent abundance of He⁺⁺ relative to H⁺. Table 1 summarizes the free–free spectral indices relevant to our analysis (see Appendix A, Section A3 for details).

Table 1

Observed and free–free spectral indices at T_e= 7000 K.

α315	α3626	α25031
Free-free indices
−0.125a	−0.137	−0.152
PN sample
−0.06 ± 0.01(0.22)b,c	−0.18 ± 0.03(0.21)d	−0.43 ± 0.03(0.14)e

^aThe free–free indices bear a 1σ (one standard deviation) uncertainty of 10⁻³, stemming from the uncertain contribution of the He⁺ free–free continuum. ^bThe root-mean-square (rms) scatter is indicated in parenthesis. ^cFrom all 31 objects with 5-GHz flux densities and χ²_CBI < 2. ^d From all 32 objects with χ²_CBI < 2. ^e From all 17 objects with 250-GHz flux densities and χ²_CBI < 2.

3 CBI OBSERVATIONS

The CBI is a planar interferometer array with 13 antennas, each 0.9 m in diameter, mounted on a 6 m tracking platform, which rotates in parallactic angle to provide uniform uv coverage. Its uniform-weight synthesized beam is ∼6 arcmin. The CBI receivers operate in 10 frequency channels, with 1-GHz bandwidth each, giving a total bandwidth of 26–36 GHz. It is located in Llano de Chajnantor, Atacama, Chile. In this section we describe the analysis of total intensity and polarization observations of PNe.

3.1 Total intensity

We acquired CBI total intensity observations of a sample of 36 PNe, whose flux densities are listed in Table 2. The visibility data were reduced and edited using a special-purpose package (cbical, developed by T. J. Pearson). Flux calibration was obtained by reference to Mars or Jupiter, and tied to the WMAP (Wilkinson Microwave Anisotropy Probe) temperature for Jupiter, T_Jup= 147.3 ± 1.8 K (Readhead et al. 2004, and references therein). The fractional uncertainty on T_Jup, σ(T_Jup) = 1.2 per cent, affects the overall calibration of the CBI.

Table 2

Summary of the CBI and SIMBA measurements. Reduced χ² values, χ²_CBI, refer to parametrized models of the CBI visibilities. All flux densities are reported in mJy.

Name	PNG	5a	14a	31 GHz	α3626	χ2CBI	250 GHz	Δb	Δ/σ(Δ)	δc	α25031
SwSt1	001.5 − 06.7	130	240	190 ± 6	−0.33 ± 0.30	1.18	14 ± 102	−125 ± 102	(−1.2)	90 ± 73	−1.3 ± 3.5
NGC6369	002.4 + 05.8	2002	1718	1452 ± 30	−0.20 ± 0.06	1.20	479 ± 54	−583 ± 58	(− 10.0)	55 ± 5	−0.5 ± 0.1
Hb4	003.1 + 02.9	170	148	138 ± 5	0.07 ± 0.39	1.11
Hb6	007.2 + 01.8	243	241	195 ± 6	−0.40 ± 0.25	1.20
NGC6309	009.6 + 14.8	102	146	95 ± 7	−1.37 ± 0.87	1.08	16 ± 16	−54 ± 17	(− 3.2)	77 ± 23	−0.9 ± 0.5
M2-9	010.8 + 18.0	36	44	55 ± 4	0.05 ± 0.84	1.11	273 ± 48	233 ± 48	(4.8)	−582 ± 129	0.8 ± 0.1
NGC6537	010.1 + 00.7	610	557	477 ± 12	−0.63 ± 0.15	2.13	331 ± 111	−17 ± 111	(−0.2)	5 ± 32	−0.2 ± 0.2
NGC6818	025.8 − 17.9	281	270	270 ± 18	−1.04 ± 0.69	1.47	<36
NGC6741	033.8 − 02.6	220	183	181 ± 13	−1.42 ± 0.79	1.21
NGC6572	034.6 + 11.8	1260	-	1073 ± 22	−0.10 ± 0.06	1.16	442 ± 54	−343 ± 56	(− 6.1)	44 ± 7	−0.4 ± 0.1
NGC6790	037.8 − 06.3	240	256	267 ± 8	0.05 ± 0.24	1.20	118 ± 60	−77 ± 60	(−1.3)	40 ± 31	−0.4 ± 0.2
NGC7009	037.7 − 34.5	750	649	512 ± 11	−0.23 ± 0.08	1.26	157 ± 27	−217 ± 28	(− 7.7)	58 ± 7	−0.6 ± 0.1
CN3-1	038.2 + 12.0	65	62	55 ± 5	−0.39 ± 0.99	1.15
Vy2-2	045.4 − 02.7	50	-	234 ± 6	−0.20 ± 0.21	1.18	144 ± 104	−27 ± 104	(−0.3)	16 ± 61	−0.2 ± 0.3
M1-71	055.5 − 00.5	204	-	178 ± 7	0.37 ± 0.35	1.17
IC4997	058.3 − 10.9	100	127	81 ± 4	−0.20 ± 0.49	1.14	64 ± 64	4 ± 64	(0.1)	−7 ± 108	−0.1 ± 0.5
NGC6886	060.1 − 07.7	105	102	58 ± 5	−0.82 ± 1.00	1.15
NGC246	118.8 − 74.7	247	248	104 ± 12	−1.20 ± 1.30	0.32
IC418	215.2 − 24.2	1613	1529	1334 ± 30	−0.39 ± 0.11	1.22	704 ± 83	−271 ± 86	(− 3.2)	28 ± 9	−0.3 ± 0.1
NGC1360	220.3 − 53.9	-	-	132 ± 16	1.16 ± 1.35	1.19
IC2165	221.3 − 12.3	188	186	147 ± 7	0.18 ± 0.53	1.12
NGC2440	234.8 + 02.4	370	325	267 ± 13	0.50 ± 0.47	1.29	251 ± 64	56 ± 65	(0.9)	−29 ± 33	−0.0 ± 0.1
NGC3242	261.0 + 32.0	896	739	544 ± 12	−0.33 ± 0.09	1.25	212 ± 34	−186 ± 35	(− 5.3)	47 ± 9	−0.5 ± 0.1
NGC3132	272.1 + 12.3	235	198	184 ± 20	1.81 ± 1.19	1.34	161 ± 66	27 ± 68	(0.4)	−20 ± 51	−0.1 ± 0.2
NGC2867	278.1 − 05.9	252	265	194 ± 20	0.73 ± 1.08	1.31	21 ± 56	−121 ± 58	(−2.1)	85 ± 39	−1.1 ± 1.3
IC2501	281.0 − 05.6	261	236	167 ± 6	−0.03 ± 0.33	1.11	69 ± 29	−53 ± 29	(−1.8)	43 ± 24	−0.4 ± 0.2
Hen2-47	285.6 − 02.7	170	187	136 ± 9	0.59 ± 0.73	1.17
IC2621	291.6 − 04.8	195	175	165 ± 5	0.23 ± 0.27	1.10	22 ± 30	−99 ± 30	(− 3.3)	82 ± 25	−1.0 ± 0.7
NGC3918	294.6 + 04.7	859	765	626 ± 19	−0.44 ± 0.25	1.36	167 ± 40	−291 ± 42	(− 6.9)	63 ± 9	−0.6 ± 0.1
NGC5315	309.1 − 04.3	480	366	443 ± 22	0.51 ± 0.52	1.24	141 ± 38	−183 ± 41	(− 4.4)	56 ± 12	−0.5 ± 0.1
Hen2-113	321.0 + 03.9	115	160	133 ± 5	0.05 ± 0.34	1.15
Hen2-142	327.1 − 02.2	65	68	65 ± 4	2.12 ± 0.70	1.76
Pe1-7	337.4 + 01.6	117	119	111 ± 6	0.70 ± 0.54	1.36
NGC6153	341.8 + 05.4	477	559	475 ± 20	0.71 ± 0.43	2.40	208 ± 38	−139 ± 41	(− 3.4)	40 ± 11	−0.4 ± 0.1
NGC6072	342.1 + 10.8	-	-	64 ± 34	−1.03 ± 5.12	1.74	<87
NGC6302	349.5 + 01.0	3100	3034	2578 ± 53	0.19 ± 0.05	2.98	1963 ± 199	78 ± 203	(0.4)	−4 ± 11	−0.1 ± 0.0
H1-12	352.6 + 00.1	719	-	-	-	-	383 ± 65	−17 ± 76
M1-26	358.9 − 00.7	400	420	-	-	-	346 ± 81	124 ± 84
19W32	359.2 + 01.2	21	-	18 ± 3	−10.27 ± 1.13	4.34
Hb5	359.3 − 00.9	548	551	-	-	-	506 ± 71	201 ± 77

^a5- and 14-GHz measurements from the compilation of Acker et al. (1992), taken mostly from the Parkes 64-m data of Milne & Aller (1982), and bearing a ∼10 per cent calibration uncertainty. ^bΔ is the difference between the observed 250-GHz flux density and optically thin free–free emission extrapolated from 31 GHz, if available, or 5 GHz otherwise. ^cδ= 100 ×Δ/(F(250 GHz) −Δ) is the per cent fraction of 31-GHz emission not due to free–free emission. ^dThe 5-GHz flux for NGC 1360 is contaminated by a point source. ^eUpper limits are 3σ.

Approximate 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. In general we used an on-source integration time of 8 min. We found that for point sources the acquisition of a reference field is not mandatory, since restricting to the longer baselines avoids any ground contamination (e.g. Casassus et al. 2006). We thus used the differenced data set when available, and discarded baselines shorter than 200λ otherwise. The baseline-restricted and differenced data sets have similar signal-to-noise ratios on point source measurements.

We extracted integrated flux densities by fitting parametrized models to the visibilities. The uncertainties on the inferred flux densities correspond to the Δχ²= 1 contour. Point-source models gave better reduced χ² values than elliptical Gaussians or core-halo models, except for NGC 1360, which we fit with a 4.9 × 2.3 arcmin² elliptical Gaussian [full width at half-maximum (FWHM)], with its major axis at a position angle of 22°.9, positive east of north.

The use of point-source models to extract flux densities assumes there is no contaminating emission within the CBI primary beam, which has a FWHM of 45.2 arcmin. Additional sources within the primary beam result in high reduced χ² values (⁠ hereafter). We checked by inspection that extended emission is present for greater than ∼2.0, which is the cut-off we adopted to select reliable point source flux densities. We attempted to extract flux densities from images of the fields with extended emission, and obtained the same value as listed in Table 2. The extended emission that fills the CBI primary beam is featureless and difficult to separate from the PN. In addition to the extended emission, there are also bright

point sources near NGC 6537. 19W32 is altogether drowned by the diffuse background. Thus the sky-plane measurements have the same degree of unreliability as the visibility fits.

A measurement is considered reliable if its uncertainty, as estimated from the point-source models, indeed corresponds to a 1σ deviation. Several integrations on different nights for a subsample of 22 objects confirmed that, when χ²_CBI < 2.0, all nightly measurements are consistent within 3σ with the flux densities derived from model fitting the combined visibility data set. We also checked for consistency between the differenced and baseline-restricted data sets.

In addition to the statistical uncertainties due to thermal noise, the final error bars also include uncertainties in the overall calibration of the CBI for the PNe observations. The systematics stem in part from uncertainties in our primary calibrator (Jupiter), and in part from gain variations during a given night. We thus compared PN flux densities obtained with two different calibrators, during the same night. We measured the fractional flux density differences, Δ_gain= | F_cal1−F_cal2 |/F_cal1, when calibrating against calibrator 1 or 2. On the night of 2003 January 15, NGC 3918 was observed along with Jupiter and Tau A; on 2001 June 3 NGC 3242 was observed along with J1230+123 and Mars. Using the same reduction (cbical) scripts and analysis tools as those used to produce the flux density reported in Table 2, we obtain Δ_gain= 0.68 per cent for NGC 3918, and 0.80 per cent for NGC 3242. For the PNe observations, the uncertainty on the CBI response is thus of order σ_gain= 1 per cent.

Including the uncertainty on T_Jup, we take as a conservative estimate of fractional systematic uncertainties σ_systematic= 2 per cent. The final uncertainty on the flux densities F_CBI reported in Table 2 is σ²_total= (σ_systematic F_CBI)²+σ², where σ derives from model-fitting the combined visibility data sets.

The model-fit spectral indices over 26–36 GHz are all consistent with free–free emission (see Table 1). The accuracy of the CBI indices reaches ∼0.1 for NGC 6369, 6572 and 7009, IC 418 and NGC 3242.

The CBI 31-GHz and comparison 5-GHz flux densities are also consistent with optically thin free–free, within the 10 per cent calibration uncertainty of the literature data (Table 1). The only exceptions are Pe 1-7, M2-9, SwSt 1 and Vy 2-2, all with excesses at 31 GHz and point-source models with < 2.

In comparison with the 14-GHz data set, we obtain 〈α¹⁴₅〉= 0.00 ± 0.03 with an rms scatter of 0.19 and 〈α³¹₁₄〉=−0.23 ± 0.03 with a scatter of 0.37. is steeper than the free–free value by 3.3σ, if the uncertainties on the 14-GHz are indeed smaller than 10 per cent rms.

3.2 CBI polarization upper limits

Four PNe were observed by the CBI in polarization during 2002 November and 2005 June: NGC 7009, 1360, 246 and 7293 (the Helix nebula). The polarization capability of the CBI has been demonstrated with deep observations of cosmic microwave background (CMB; Readhead et al. 2004) and also in H ii regions (Cartwright et al. 2005; Dickinson et al. 2006). Details of the polarization calibration and data reduction procedures can be found in the aforementioned papers.

Natural-weighted maps of Stokes Q and U were made using the difmap package and combined to make polarization intensity maps, ⁠. All maps were consistent with noise – no significant polarization was detected in any of the four PNe. Table 3 summarizes the results and includes the 3σ upper limit on the polarization fraction, p=P/I, derived from the data.

Table 3

CBI 31 GHz polarization measurements for four PNe. Upper limits are 99 per cent CL.

PNe	Polarization fraction
NGC 7009	<8.5 per cent
NGC 1360	<8.5 per cent
NGC 246	<17 per cent
NGC 7293	<2.4 per cent

The 99 per cent confidence level (CL) upper limit on the polarization fraction was calculated assuming the Q and U maps have normal distribution with zero mean and a standard deviation of σ. The polarized intensity, ⁠, follows the Rayleigh distribution f(P). The CL F(P_u) associated to an upper limit P_u on the polarized intensity is thus ⁠, and a CL of 99 per cent corresponds to P_u≈ 3σ. Finally, the upper limit on the fraction of polarized intensity is P_u/I_max, where I_max is the peak total intensity (which in units of Jy beam⁻¹ is equal to the source flux density for unresolved sources).

Strong radio polarization is not expected for PNe since the emission is typically dominated by free–free emission, which is intrinsically unpolarized. However, the excess anomalous component could perhaps be significantly polarized, depending on the physical mechanism producing the emission. The most constraining and interesting result here is for the Helix nebula where significant excess emission was observed at 31 GHz with the CBI (Casassus et al. 2004). Since free–free radiation is not polarized, we can limit the polarization fraction on the excess emission in the Helix to <3.8–12 per cent, for 36–80 per cent free–free emission at 31 GHz.

4 SIMBA OBSERVATIONS

SIMBA, at SEST,1 is a 37-channel bolometer array, operating at 1.2 mm (250 GHz). The half-power beamwidth of a single element is 24 arcsec. We observed each PN with a scanning speed of 80 arcsec s⁻¹, in 51 steps of 8 arcsec orthogonal to the scan direction, thus obtaining 900 × 408 arcsec² maps. The SIMBA scans are reduced with the standard recipe for point sources using the mopsi package written by Robert Zylka (IRAM, Grenoble). Flux calibration is carried out by comparison with Uranus maps. We extracted the calibration flux densities of Uranus by integrating the Uranus specific intensity map in a circular aperture 1 arcmin in radius about the centroid of an elliptical Gaussian fit. A residual sky background was estimated from the median intensity in a ring 1 arcmin wide and immediately surrounding the flux extraction aperture.

We observed 35 PNe during three different observing runs in 2001 and 2002. Flux densities were measured by integrating each PN specific intensity map in a circular aperture 1 arcmin in radius about the centroid of an elliptical Gaussian fit to the PN image, and subtracting a residual sky background, as in the case of Uranus.

The uncertainty assigned to each individual integration is the rms noise in a circular annulus surrounding our extraction aperture (and taking into account correlated pixels). Up to three integrations from different days on individual objects were combined to reduce the uncertainties.

A 10 per cent calibration uncertainty affects SIMBA flux densities. The systematic error comes from the uncertainty in the brightness temperature estimate of Uranus and the uncertainty in the determination of telescope parameters, such as the aperture efficiency. A fractional systematic uncertainty with a standard deviation of 10 per cent was therefore added in quadrature to the combined statistical uncertainties derived from the aperture photometry.

Fig. 1 and illustrates our SIMBA images. All the SIMBA images can be approximated to elliptical Gaussians. The only case where the nebular solid angle Ω_N, as inferred from the Gaussian fits, is larger than 1.5 times the SIMBA beam, Ω_S, is NGC 6369, with Ω_N= 1.6 Ω_S. In the case of NGC 6302 (shown on Fig. 1) low-level emission can be resolved by SIMBA. The north-east–south-west (NE–SW) low-level extension in Fig. 1 follows the bipolar axis of NGC 6302. The stripes and negatives in the SIMBA maps are proportional to the object flux densities. They are artefact inherent to the technique of fast scanning, which involves a high-pass filter to cancel sky emission. These artefacts are also present on the Uranus calibration images, and are taken into account by the flux calibration procedure.

Figure 1

SIMBA map of NGC 6302. The grey-scale is in units of MJy sr⁻¹, and the peak intensity is 138 MJy sr⁻¹ (or 2.1 Jy beam⁻¹, for SEST's circular beam with a FWHM of 24 arcsec at 250 GHz). x-axis and y-axis show offset J2000 RA and Dec. in arcsec. Contour levels at −4.9, −3.6, 0.2, 7.5, 19.5, 37.4, 62.4, 95.7, 138.4 MJy sr⁻¹. The SIMBA beam is indicated on the lower-left corner of the map.

5 CBI–SIMBA COMPARISON

As summarised in Table 2, we have obtained 31- and 250-GHz flux densities for a total of 21 objects, of which three have ⁠. Optically thin free–free emission at 31 GHz can be extrapolated to 250 GHz with a spectral index of α=−0.15. We find that 20 out of 21 objects have 250-GHz data consistent, within 1σ, with a deficit over the free–free level extrapolated from 31 GHz. Only four objects have < 1σ excesses at 250 GHz rather than a deficit (i.e. Δ > 0, in the notation of Table 2). The 250-GHz deficit is significant in nine PNe, at 3σ. The only significant detection of an excess at 250 GHz is in M 2-9, with a 5σ excess at 250 GHz.

The CBI–SIMBA spectral index, averaged over the 17 measurements of α²⁵⁰₃₁ listed in Table 2 that have < 2, and excluding M2-9, is 〈α²⁵⁰₃₁〉=−0.43 ± 0.03, with a weighted rms scatter of 0.14. Thus, the observed 〈α²⁵⁰₃₁〉 is significantly different from the free–free value.

We searched for a correlation between the 31-GHz excess and the mid-IR continuum due to VSGs, as observed in the anomalous CMB foreground. Unfortunately the IRAS 12 μm band is contaminated by strong ionic lines as well as free–free emission, while ground-based 10-μm spectra available from the literature are extracted from varying apertures. We can still compare the 31-GHz excess to the IRAS 25-, 60- and 100-μm flux densities. The IRAS 60-μm band turned out to give the best correlation with the 250-GHz deficit. Fig. 4 plots the difference between the expected level of free–free emission at 250 GHz, extrapolated from 31 GHz, against IRAS 60 μm. The 10-μm dust emission features compiled from the literature are indicated by ‘C’ for polycyclic aromatic hydrocarbons (PAHs), ‘c’ for SiC, ‘O’ for silicates, ‘+’ for weak and featureless continuum and ‘*’ when no data are available (as in Casassus et al. 2001). A straight line fit to the data, required to cross the origin, gives a dimensionless slope of (−3.05 ± 0.18) × 10⁻³. The cross-correlation with IRAS 60 μm could simply reflect an overall scaling of the flux densities with nebular luminosity. The cm-excess does not seem related to the 10-μm dust emission features.

Figure 4

SED of NGC 6369. We plot flux density in Jy (y-axis) as a function of frequency in GHz (x-axis). The CBI and SIMBA datapoints are indicated by circles. (a) Synchrotron: the solid line includes a synchrotron component with a break at 31 GHz (dash–double-dotted), a core free–free component at 10⁴ K (dash–dotted), a sub-mm modified blackbody (dotted), and the emission expected from the free–free screen, chosen here to be a 20 K (dashed). (b) Spinning dust: the dashed line is a free–free component at 7400 K, with ⁠. A spinning dust component for the ‘WNM’ case of Draine & Lazarian (1998), is shown in thick dotted line and contributes to the total flux density shown in solid line. The other ISM environments considered by Draine & Lazarian (1998) are shown in dotted line, scaled by the same factor as the ‘WNM’ case. The low-frequency data cannot be reconciled with the spinning dust models. (c) Needles: the solid line includes extinction from very cold metallic needles, the dashed line is the free–free component, and the dotted line is the sub-mm tail of a modified blackbody, attenuated by grey needle extinction. Since needle absorption does not extend appreciably below 100 GHz, the dashed line indicates nebular emission had it been pure free–free.

We also calculated the 250-GHz deficit when extrapolating the free–free level from the 5-GHz data. The resulting deficits are still significant, albeit for a reduced number of objects. The reduced significance may be due in part to the 10 per cent calibration uncertainty of the 5-GHz data, but also to a 5–31 GHz spectral index greater than the optically thin value (see Section 3). In this case a straight line fit to the 250-GHz deficit as a function of IRAS 60-μm band flux density, for the same objects as for the 31-GHz cross-correlation of Fig. 2, gives a dimensionless slope of (−2.65 ± 0.27) × 10⁻³.

Figure 2

Deficit at 250 GHz from the expected level of free–free emission extrapolated from the CBI measurements, excluding M2-9 for clarity. The y-axis shows the difference between the SIMBA flux density and the free–free level extrapolated from 31 GHz, and the x-axis shows IRAS 60-μm flux density. The dotted line is a fit to the data.

Relative to the total 31-GHz flux density, the excess emission over the level of free–free emission extrapolated from 250 GHz represents ∼0.3–0.7. In order to avoid induced correlations from scaling effects, we plot in Fig. 3 the ratio of the 31-GHz excess to the total 31-GHz flux density against the ratio of the IRAS 60-μm flux density over the 5-GHz flux density. We have excluded from Fig. 3 all objects with large 1σ uncertainties on the fractional excess (>0.3). The remaining 10 objects, not including M2-9, have an average fractional excess of 0.51 ± 0.03 with an rms scatter of 0.11. Interestingly, the fractional excess seems to correlate with the ratio of far-IR and radio flux, as expected for a dust-related emission mechanism at 31 GHz.

Figure 3

The ratio of the 31-GHz excess (extrapolated from the 250-GHz measurement, assuming it is entirely due to free–free emission) over the total 31-GHz flux density, against the ratio of the IRAS 60-μm flux density over the 5-GHz flux density. The dashed line is a fit to the data that crosses the origin, with a slope of 4.54 ± 1.5 and reduced χ² of 0.5.

We caution that the IR-radio slopes derived in this section cannot be directly compared with the slopes observed in the ISM at large. PNe have much hotter dust temperatures (100–200 K, e.g. Casassus et al. 2001) than dark or cirrus clouds (15–25 K).

6 COMPARISON WITH PREVIOUS WORK

6.1 NGC 7027

An important PN which is absent from our sample is northern NGC 7027, the brightest at radio wavelengths. Terzian et al. (1974) and Hoare et al. (1992) constructed the SED of NGC 7027, finding its radio spectrum follows free–free emission. In this respect NGC 7027 is similar to NGC 6302 (discussed in Section 8.2.2), another high-excitation PN with a massive progenitor. Both are surrounded by molecular envelopes with copious amounts of dust. Although any cm-excess appears to be negligible in NGC 7027 and 6302, the mm-wave data could be affected by either molecular lines or cold dust, and definitive conclusions in these objects require further analysis. NGC 7027 and 6302 are reminders that PNe form a heterogenous class.

6.2 Previous mm-wave flux densities from the literature

Deviations from standard free–free spectra in the radio continua of PNe were already present in the SEDs obtained by Hoare et al. (1992), although not considered significant. It can be seen in fig. 1 from Hoare et al. (1992) that NGC 6572, 6302, 6537 and 6543 show deficits at 200–300 GHz. We confirm that our 1.2-mm flux densities are consistent with those measured at 1.1 mm by Hoare et al. (1992), within 2σ, for the four objects we have in common: NGC 6572, 6302, M2-9 and NGC 6537. These objects are studied in detail below.

6.3 Preliminary mm-wave flux densities from SEST

The 250-GHz deficit over the expected level of optically thin free–free emission was originally found by one of us (L-ÅN) with heterodyne and single pixel bolometer observations at SEST. However, in pointed observations the measurement of flux densities is hampered by pointing uncertainties. Until the development of bolometer arrays such as SIMBA, designed to produce images of the sky, it was difficult to ascertain whether the 250-GHz deficits were not due to pointing errors (which may also have been the problem in Hoare et al. 1992).

The SEST heterodyne observations were performed between 1995 and 1996 with the 15-m SEST. The data were obtained using simultaneously two single-channel SIS (superconductor-insulator-superconductor) receivers at 99 and 147 GHz with typical system temperatures of 150 K (SSB) above the atmosphere. The telescope beamwidth is 50 arcsec at 99 GHz and 34 arcsec at 147 GHz (FWHM). Two acousto-optical spectrometers (AOS) were used simultaneously with bandwidths of about 1 GHz, channel separations of about 0.7 MHz and resolutions of about 1.4 MHz. The telescope was used in a dual beam switch mode with the source alternately placed in each of the two beams, a method that yields very flat baselines, and allows the determination of continuum levels as the offset of the baseline from the zero level. The beam separation was about 11.5 arcmin. The flux densities were determined from the chopper-wheel corrected antenna temperature, ⁠, using aperture efficiencies (η_A), determined from observations of planets, of 0.60 and 0.52 at 99 and 147 GHz, respectively, and are estimated to be accurate to within 10 per cent (1σ). The pointing was checked during the observations using SiO maser sources, and the pointing accuracy is estimated to be better than ± 3 arcsec.

The SEST single pixel bolometer observations were performed between 1993 and 1995, using a single channel bolometer developed at MPIfR (Max-Planck-Institut für Radioastronomie), operating at a wavelength of 1200 μm (i.e. for SIMBA). Beam switching was done with a focal-plane chopper having a horizontal beam throw of 70 arcsec. The observations were done in the ON–OFF mode, alternating the source between the two beams after 10 s of integration time. Pointing was checked on nearby quasars and was found to be better that 3 arcsec rms. The atmospheric opacity was determined frequently through sky dips. Uranus was used as a calibrator source, and the estimated uncertainty in intensities are about 10 per cent (1σ). For single-pixel observations we report point-source flux densities, with uncertainties derived from the rms dispersion of all ON–OFF pairs. The maps were obtained by scanning the telescope in azimuth at a rate of 8 arcsec s⁻¹, with adjacent scans separated by 8 arcsec in elevation. The maps were reduced using the nod2 software (Haslam 1974) adapted for SEST. Flux densities were extracted with Gaussian fits.

We summarize in Table 4 the measurements from the preliminary SEST observations, both from the single pixel bolometer and the heterodyne backends. In general the previous pointed observations at 230 GHz are within 3σ with the 250 GHz values from Table 2. The only exceptions are Hb 5, NGC 6537 and 3132, where the pointed data give lower flux densities (consistent with a pointing error, or with nebular extensions larger than the beam), and NGC 7009, where the difference with the 230 GHz map is at 3.1σ. The concordance of the whole data set on NGC 6369 is particularly comforting, and will be commented further below.

Table 4

Previous flux densities from SEST. Values are reported in mJy. Systematic uncertainties are not included, and are assumed to be of ∼30 per cent.

PNe	Single pixel		Heterodyne
	On–off	Map	147 GHz	99 GHz
NGC 3132	30.1 ± 15.0
NGC 3242	180.2 ± 9.0
NGC 3918	149.4 ± 7.9		497.7	342.6
NGC 5315	136.4 ± 9.9
NGC 6072	16.6 ± 9.4
M2-9	276.8 ± 11.2		107.4	149.7
NGC 6302	2091.9 ± 13.7	1890	1890	2054
NGC 6369	396.5 ± 12.9	485	783	1094
M1-26	185.7 ± 15.6		368.1	487.7
Hb 5	89.7 ± 8.6		250.2	337.4
NGC 6537	80.1 ± 11.1
NGC 6572	670.8 ± 12.1	434	843	882.1
NGC 6790	154.8 ± 9.4		120.6	176.0
Vy 2-2	94.2 ± 15.2
NGC 7009	178.2 ± 16.5	279	231	379.6

6.4 Absence of radio PN haloes

Neither the CBI nor the SIMBA data sets show evidence for radio haloes. When compared with the point-source models, none of the core-halo models that we tried to fit the CBI visibilities gave improved χ²_CBI values. Extended haloes, about 4 arcmin in diameter and larger, have indeed been observed in a few PNe (Corradi et al. 2003). However, the optical Hα+[N ii] haloes are extremely faint. The halo intensities are of order 1/1000 the values in the central regions that dominate our radio measurements. The integrated halo line fluxes are about 1/100 the fluxes of the compact regions.

6.5 Is the PN cm-wave excess seen in H ii regions?

As natural comparison objects to PNe, H ii regions are also photoionized nebulae with thermal continua. Does their radio continua also bear the cm-excess? H ii regions are extended and complex objects; the measurement of precise radio SEDs is difficult. Single-dish maps at frequencies higher than ∼14 GHz require cancelling the sky emission with either fast scanning (as for SIMBA) or differencing, while interferometric data are affected by flux loss and a wavelength-dependent (u, v) coverage. Another difficulty with mm-wave continuum data of H ii regions is the possibility of molecular line emission from surrounding or coincident molecular material (Sutton et al. 1984).

Notwithstanding these difficulties, Dickinson et al. (2007) reported that the radio SEDs of six H ii regions constructed from CBI observations and literature data are compatible with free–free emission. They constrain the possibility of spinning dust, and note a 3σ 31-GHz excess in G284.3−0.3 similar to spinning dust in its narrow profile.

However, the SEDs of H ii regions at mm-wavelengths are not established. A close look at the Boomerang and SEST data reported by Coble et al. (2003) reveals that H ii regions may be affected by the same 250-GHz deficit seen in PNe. The SED of RCW 38 can be constructed with the Boomerang data (Coble et al. 2003) and the CBI data point.2 We fit the 31–400 GHz SED of RCW 38 with two power laws. The best-fitting index at high frequencies is α=+3.58 ± 0.69, which is consistent with greybody emission (with an emissivity index of β= 1.58 ± 0.69). The low-frequency data are fitted by α=−0.22^−0.08_−0.38± 0.14 (where we have indicated the Δχ²= 1 limits). Given the uncertainties such an index over 31–90 GHz is consistent with the indices seen in PNe.

The H ii region SEDs published by Malkamäki et al. (1979) also cover the mm-wave regime. The 10–100 GHz index is α∼−0.25 in W51, and is α∼−0.36 in W3. Uncertainties are difficult to estimate from the figures in Malkamäki et al. (1979), but our point here is that given the existing data, H ii regions may well show the same cm-wave excess as PNe.

7 CANDIDATE EMISSION MECHANISMS

In this section we speculate on the emission mechanisms that could explain the observed PN SEDs. We take the case of NGC 6369 as an example, since it is the object with the most significant 250-GHz deficit. We assume a distance for NGC 6369 of 330 pc, and a nebular diameter of 38 arcsec (from optical images, as compiled by Acker et al. 1992). Refer to Fig. 4(c) and Section 7.3 for a comparison of the data and a pure free–free spectrum. The free–free continua were calculated with the exact Gaunt factors for hydrogen.

7.1 A synchrotron component?

Assuming the 250 GHz continuum is due only to free–free emission and possibly to the thermal dust Rayleigh–Jeans tail, any synchrotron component should be very steep above 31 GHz, with a spectral index steeper than ⁠, or about −0.4 (Table 1). This is in contrast with 〈α³¹₅〉, which is essentially flat.

Synchrotron self-absorption is not expected in PNe. However, the flat index at low frequencies could be explained by free–free absorption with a turnover frequency at 20–30 GHz. For a typical emission measure of 10⁶ cm⁶ pc⁻¹, a turnover frequency of 40 GHz requires a free–free screen with an electron temperature T_e < 600 K for the bulk of the nebula.

7.1.1 Background

Synchrotron radiation was sought in PNe during the 1960s as a means to detect the magnetic field required by Gurzadyan (1969) in the theory of magnetic shaping (see also e.g. Pascoli 1985; Chevalier & Luo 1994). The evidence put forward by Gurzadyan (1969) was based on an UV excess in photospheric spectra, and by the mention of non-thermal emission being observed in a few objects. However, only the case of NGC 3242 appears to be backed by a reference. Although Menon & Terzian (1965) report a non-thermal spectrum for NGC 3242, their data are not confirmed by the flux densities given in Calabretta (1982), Milne & Aller (1982) and Condon & Kaplan (1998).

The possibility of non-thermal emission from the interaction of a magnetized and fast wind with dense knots has been investigated by Jones, Kang & Tregellis (1994) and Jones, Ryu & Tregellis (1996), in diverse astrophysical contexts. They expect synchrotron emission to arise at the shocked interfaces between knots and fast wind. Synchrotron should be strongest along the stretched magnetic field lines around the perimeter of the clumps, and increase at the time of knot disruption.

An interesting alternative has been proposed by Dgani & Soker (1998). In the specific case of PNe, Dgani & Soker (1998) argue that simultaneous cosmic ray acceleration and magnetic field enhancement could occur in the turbulent mixing layers at the sides of the knots and their turbulent wakes. They find that such synchrotron emission could be detectable, at a level of 1 mJy at ν∼1 GHz, in the case of hydrogen deficient nebulae, where the level of free–free emission is reduced. Dgani & Soker (1998) claim the detection in PN A30 of faint (0.1 mJy), compact (3 arcsec), non-thermal and circumstellar emission at 8 GHz. However, the data have not been published.

There is thus no data available to back the synchrotron hypothesis for the cm-excess. Recently, Cohen et al. (2006) have reported the detection of two non-thermal hotspots in the faint PN surrounding the OH/IR star V1018 Sco (with spectral indices of −0.85 ± 0.01 and −0.95 ± 0.11 over 5 and 10 GHz). Cohen et al. (2006) attribute the ionization of this object to the interaction between a fast asymptotic giant branch (AGB) wind that has been recently turned on and is overrunning a precursor slow wind. It is very different from the photoionized nebulae which make the focus of this work.

Magnetic fields have however been detected by Miranda et al. (2001) in the young PN K3-35, through the detection of a circularly polarized OH maser line at 1.665 GHz, implying field intensities of order ∼1 mG. Fields of ∼1 mG are in the range of the theoretical predictions of Dgani & Soker (1998), based on the equipartition of the magnetic field and fast wind energies.

Vlemmings et al. (2006) found water masers in the proto-PN W43A, whose polarization traces the precession of a strong collimating magnetic fields, >1 mG on 1000-au scales. Bains et al. (2003, 2004) provide further detections of magnetic fields in proto-PNe. Sabin et al. (2007) observed the nebular sub-mm dust polarization concomitant to a prevasive ∼1-mG field.

In fact, the energetic requirements of synchrotron radiation are easily met in PNe. Given a magnetic field and a synchrotron index α observed over a frequency range from ν₁ to ν₂, the total energy U_e in cosmic rays can be estimated from (Moffet 1975)

1

2

3

4

The observed synchrotron luminosity in NGC 6369, between 0.1 MHz and 31 GHz, is L_obs= 8.7 × 10³⁰ erg s⁻¹, as obtained by integrating a power-law spectrum with F_ν(31 GHz) = 1.5 Jy and α=−0.3. We find that, for a 1-mG field permeating the whole nebula, the total energy in cosmic rays required to account for L_obs is ∼5 × 10⁴² erg. Since the lifetime of a typical nebula is 10⁴ yr, such an amount of cosmic-ray energy requires a stellar luminosity of ∼4 × 10⁻³ L_⊙, much less than the observed tip-of-the-AGB luminosity of 10⁴ L_⊙.3

It is interesting to compare the above energetic requirement and the calculation in Dgani & Soker (1998). In their idea synchrotron should stem from the whole ensemble of globules and their overlapping turbulent wakes, encompassing a linear size R_ens. Dgani & Soker (1998) assume equipartition between magnetic field, cosmic ray electrons and fast wind. Their equation (9) gives a synchrotron luminosity of only 4.9 × 10²⁸(R_ens/5 × 10¹⁵ cm)³ erg s⁻¹. Accordingly Dgani & Soker (1998) predict synchrotron flux densities 1/100 of the values required for the cm-excess. However, a value of R_ens= 5 × 10¹⁵ cm is actually very small for a typical PN. As seen in the Helix nebula, the ensemble of knots covers at least 1/3 of the nebular solid angle, and the nebulae considered in this work are typically ∼0.1 pc. Given that the Dgani & Soker (1998) synchrotron luminosity scales as the cube of R_ens, a mere increase to 2.5 × 10¹⁶ cm, still much less than 0.1 pc, approximates the cm-excess.

7.1.2 NGC 6369

Given the detection of magnetic fields in PNe, we did the exercise of modelling the SED of NGC 6369 including a synchrotron component (Fig. 4a). The core free–free emission, which represents the bulk of the 250-GHz flux density, is approximated as a uniform slab at 10⁴ K. The synchrotron component is modulated by a cold screen, with a power-law distribution of opacities (as in the case of M2-9, see Section 8.2.1), which is required for NGC 6369 to fit the curvature of the <30-GHz data.4 The cold screen opacity includes both free–free and bound-free contributions, although the latter is negligible at frequencies where the screen opacity is ∼1.

We could only obtain a fit to the data by including a synchrotron break at 31 GHz. For a single power law to account for the whole SED, it would need to join the CBI and SIMBA data points, so that free–free emission would be altogether negligible compared to synchrotron. Thus, the best fit we obtained with a single power law has an index of −0.57, and a free–free screen at 15 K, with a power-law distribution of opacities (see Section 8.2.1). This model fits only marginally the low-frequency data. Note the need for a high-frequency break is alleviated by the inclusion of a spinning dust component (see Section 7.2).

Any free–free screen with a temperature below ∼1000 K can fit the data, but the electron density in the screen is an increasing function of temperature. At a fixed screen temperature of T_screen= 500 K, 〈N_e(screen)〉= 1221 ± 291 cm⁻³ and N_e(core) = 1446 ± 507 cm⁻³, while at T_screen= 10 K all of the 250-GHz flux density is due to the core free–free emission, N_e(core) = 2404 cm⁻³ and 〈N_e(screen)〉= 33 ± 20 cm⁻³. Perhaps an advantage of the lower temperatures is that the low-frequency synchrotron index can be somewhat steeper; at 1000 K the index is fixed at ∼−0.1 by the 5–31 GHz data.

For integrated flux densities the two-dimensional (2D) opacity fields can be modelled following Calabretta (1991). Section 8.2.1 on M2-9 gives more detail on the procedure, as well as an eloquent example of an opacity profile. Fig. 4(a) shows the fit to the SED of NGC 6369 with a fixed low-frequency synchrotron index of α_LF=−0.4, a fixed high-frequency synchrotron index of α_LF=−1.0, a screen at 20 K, which has an average density of 574 ± 88 cm⁻³ and a power-law opacity index γ= 6.0 (see Section 8.2.1). A modified blackbody with an emissivity index of β=+2, required to pass through the 100- and 60-μm IRAS flux densities, corresponds to T_dust= 52 K and is negligible at <250 GHz.

7.1.3 Predictions of the synchrotron interpretation

The main difficulty for the synchrotron interpretation is thus the need for a cold free–free screen. Screen temperatures of 1000 K are reminiscent of the ORL temperatures, but in the bi-abundance model the ORL gas only accounts for a very small fraction of the nebular plasma, while the 1000-K synchrotron screen represents about half the total nebular mass. In contrast, the average electron density of the very cold screen is a factor of ∼1/100–1/1000 less than the bulk nebular density. Such a 10–20 K screen could perhaps be found in the residual ionization of dense and largely neutral globules exposed to the central star wind, where Dgani & Soker (1998) expect cosmic-ray acceleration.

The cold free–free opacity profile, as a function of a 1D nebular solid angle parameter ω (Section 8.2.1), must be a power law to attenuate the synchrotron spectrum and mimic a free–free spectrum in emission below ∼31 GHz. Such opacity profiles occur naturally in PNe with a large system of globules, which seems to be the case in every object observed at sufficient resolution. A prototypical example is the Helix nebula, where knots are known to be dusty (Meaburn et al. 1992) and molecular (Huggins et al. 2002), with physical conditions approaching the predictions of Dyson et al. (1989): temperatures of T∼10 K, proton densities of up to ∼10⁶ cm⁻³. If the bulk of carbon in the knots is photoionized by the stellar UV photons with energies <16 eV, then the electron density in globules is about 10² cm⁻³. The C i free–free continuum could thus provide the cold free–free screen required to attenuate in situ the synchrotron background. If so, we expect concomitant C i radio recombination lines at the sites of synchrotron emission.

Synchrotron emission should be clumpy as well. The cometary globules providing the cold free–free screens must be coincident with the synchrotron emitting regions. The overlap between cold free–free and synchrotron regions is natural in models where synchrotron emission arises in the interaction between the fast central star wind and the clumps (Jones et al. 1996; Dgani & Soker 1998). In the clumpy model the free–free opacities should also include a filling factor (f, the fraction of nebular volume occupied by the clumps). Thus the opacity for a nebula with depth L is τ=f N_eN_iκ_ffL, where the free–free absorption coefficient is κ_ff. Values of f∼1/10 are allowed by the above synchrotron models with a factor of 2–3 increase in the electron density N_e. Smaller values of f∼1/100 can be obtained with shallower synchrotron backgrounds, with low-frequency indices of ∼−0.2.

If the cm-excess stems from the fast wind, then a relationship exists between the cm-excess, the fast wind pressure ⁠, and X-ray emission from the shocked cavities predicted by the interacting winds model (e.g. Kwok, Purton & Fitzgerald 1978). However, the sample of PNe where diffuse X-rays have been detected is too small to allow a cross-correlation with our sample. Of the PNe listed in Table 2 only NGC 7009 and 1360 are among the 13 PNe detected by ROSAT PSPC (Guerrero, Chu & Gruendl 2000). XMM and Chandra have allowed to resolve the diffuse X-rays and isolate the central star emission (Guerrero, Chu & Gruendl 2005), adding only NGC 3242 to the intersection with Table 2 (an up-to-date list of X-ray PNe can be found at http://www.iaa.es/xpn/).

The data on central star fast winds is also rather scarce, and no trends are easily found. For example, NGC 6572, with a large cm-excess, has a high mass-loss rate of in a rather slow wind with v_∞= 1190 km s⁻¹ (Modigliani, Patriarchi & Perinotto 1993). NGC 7009, also with a significant excess, has a faster wind 2770 km s⁻¹ but much smaller (Cerruti-Sola & Perinotto 1989), which gives a ram pressure ∼250 times weaker than in NGC 6572. IC 418, with a rather small cm-excess, shares a similar fast wind as NGC 7009 (v_∞= 940 km s⁻¹ and ⁠; Cerruti-Sola & Perinotto 1989).

A prediction of the synchrotron interpretation for the cm-wave continua of PNe is that the nebular continuum should be polarized at ∼60 GHz. In order to calculate the expected level of synchrotron polarization we assume the free–free screen is also embedded in a pervasive nebular magnetic field of 1 mG, and make use of the formulae given in Pacholczyk (1970, his equation 3.78). For a synchrotron spectral index of −0.5, the fraction of intrinsic (unabsorbed) polarization is Π= 69 per cent, in the case of a uniform magnetic field. Including a cold nebular screen gives Π= 69 per cent at 250 GHz, decreasing to Π∼3–15 per cent at CBI frequencies. The total nebular emission includes an unpolarized thermal component, which should account for most of the 250 GHz flux density to reconcile the continuum and recombination line data in the literature. Thus, for a uniform magnetic field, the total fraction of polarization Π_T should reach a maximum of Π_T∼ 20 per cent between 50 and 100 GHz.

However, PNe are often very point symmetric, so that the polarization fraction on flux densities integrated across the whole nebulae should vanish for perfect point symmetry. Moreover, in the interpretation of Dgani & Soker (1998), the magnetic fields required for synchrotron emission are amplified in turbulent wakes. This is a case of tangled magnetic field, without net polarization. Thus, very high angular resolution observations, at 50–100 GHz, are required to detect the polarization predicted by the synchrotron interpretation.

In summary the synchrotron interpretation predicts that the radio C i recombination lines should be coincident, on angular scales comparable to the cometary globules, with peaks in the continuum polarization between 50 and 100 GHz.

7.2 Is it spinning dust?

Spinning dust peaks at ∼30 GHz and is 30–70 times weaker at 5 GHz (depending on environment; Draine & Lazarian 1998). If spinning dust emission accounts for 30 per cent of the 31-GHz flux density, then α³¹₅= 0.045, which is consistent with 〈α³¹₅〉 within 2σ (see Table 1). Spinning dust drops with a Boltzmann cut-off above 30 GHz, and is >300 times weaker at 100 GHz. Could the 250-GHz deficit be explained by a 31-GHz excess due to spinning dust? VSGs are found aplenty in PNe, and the proximity of a hot star warrants such VSGs will be charged by the photoelectric effect.

However, taking all the low-frequency measurements into account, we find that the existing models for spinning dust cannot explain the SED of NGC 6369. Fig. 4(b) illustrates that the low-frequency data have too flat a spectral index to be explained by spinning dust. We used the spinning dust emissivity of Draine & Lazarian (1998, http://www.astro.princeton.edu/~draine/), for the warm neutral medium (WNM) case, which peaks at the lowest frequency. The spinning dust emissivity, being electric dipolar in origin, stems from the Larmor power radiated by individual grains, and has a very steep rise with frequency before peaking sharply at the Boltzmann cut-off. Magnetic dust, proposed by Draine & Lazarian (1999), also with a steep low-frequency rise, faces the same problem of not fitting the low-frequency data.

One problem with spinning dust as a candidate emission mechanism is that the current models give spectra that are too narrow in frequency. Arbitrarily changing the frequency peak does not help, since fitting the 5–10 GHz data would then miss 31 GHz. The spinning dust emissivities of other ISM environments considered by Draine & Lazarian (1998) all peak at slightly higher frequencies. Another problem is that the predicted flux densities, calculated from the highest diffuse ISM emissivities (the ‘RN’ case), barely reach peaks of 1 mJy if scaled to nebular densities of 10⁴ cm⁻³.

The PN environment may none the less lead to higher spinning dust emissivities, especially in the dense cometary knots, and the low-frequency levels could be explained by a synchrotron component. If so, spinning dust, in conjunction with a synchrotron component, could explain the cm-excess and alleviate the need for a synchrotron break.

7.3 Extinction due to mm-sized metallic needles?

Alternatively the 250-GHz deficit could be interpreted as extinction due to dust. For standard dust grains such millimetre extinction would require dust-to-gas mass ratios of order 1. However, metallic needles bypass the mass problem for extinction, and explain the 250-GHz deficit with grey extinction at wavelengths shorter than 1 mm, and an opacity decreasing as λ⁻² towards 1 cm. Using the formulae given in Dwek (2004a,b), a long-wavelength cut-off for the grey extinction of λ_○= 1 mm requires a needle aspect ratio l/a of ∼8000 for a resistivity of ρ_R= 10⁻⁶Ω cm. Imposing a unit needle opacity at 1 mm gives a total needle mass of 2 × 10⁻⁶ M_⊙ for a grain material density of ρ_m= 8.15 g cm⁻³, i.e. much less than typical nebular masses of ∼0.1 M_⊙. Yet, to be seen in absorption at 1 mm the needles would have to be extremely cold, about 1 K. Thermal balance gives lower temperatures for longer needles (Eli Dwek, private communication), and although needles are colder than normal interstellar grains (the needle temperature would be 8 K in Cas A; Dwek 2004a), temperatures < 1 K seem unphysical. A solution may lie in modelling the needle extinction in a distribution of needle aspect ratios, such that there exist a very small population of ultracold grains.

We attempted to fit the 31–250 GHz drop in NGC 6369 with needle extinction. The needle extinction is applied as a screen to the emergent flux. For simplicity we added the needle emission as if it were a foreground to the free–free nebula. A fit to the data is obtained with λ_○= 1.2 mm, maximum opacity 0.82 and T= 1.2 K. The free–free component corresponds to a uniform slab with N_e= 3373 ± 85 cm⁻³, T_e= 6580 ± 1031 K.

The SEST 99 and 147 GHz heterodyne data shown in Fig. 4(c) are in very good agreement with needle extinction. They were originally left out from the fit, and added at a later stage for comparison purposes. The SEST heterodyne data are subject to numerous uncertainties. The possibility of needle extinction will not be considered further in this work because it requires extreme temperatures.

8 COMMENTS ON INDIVIDUAL OBJECTS

In this section we fit the SEDs of selected PNe with analytic radiative transfer solutions. All fits (except the single-component free–free spectra) include modified blackbodies constrained to the 60- and 100-μm IRAS flux densities, with a dust emissivity index of β= 2.0. This choice for β may seem rather high compared to the canonical ISM value of β= 1.7. Yet higher emissivities minimize the contribution of the Rayleigh–Jeans tail to the 250-GHz flux density, thereby minimizing the 250-GHz deficit relative to the free–free level expected from 31 GHz. The PNe listed in Table 2 were selected on their being IR-bright. Their flux densities peak in the IRAS 25–IRAS 60-μm bands, and decrease from 60 to 100 μm.

Unless otherwise stated the low frequency measurements are extracted from Milne & Aller (1982, at 5 and 14 GHz), Calabretta (1982, at 408 MHz) and Condon & Kaplan (1998, at 1.4 GHz). We generally preferred single-dish measurements, when available, in the cases where the interferometric data are resolved. PN distances are notoriously difficult to measure, and are usually assigned 100 per cent uncertainties.

8.1 The largest 250-GHz deficits

8.1.1 NGC 7009

Few Galactic nebulae can be found near NGC 7009, at 21 h. Its CBI and SIMBA flux densities are the best measurements reported in this work. The 250-GHz value is averaged from three different observing runs. The 31-GHz value is averaged over seven different nights and spread over three years.

The SED of NGC 7009 in Fig. 5 is very similar to NGC 6369. As in Section 7.1, a fit to the data can be obtained with a broken synchrotron spectrum, with a low frequency index of α_LF=−0.4, and a high frequency index α=−1.0 absorbed by a free–free screen with 〈N_e(screen)〉= 345 ± 54 cm⁻³, T_e(screen) = 50 K and a power-law opacity distribution with γ= 6.0. The core free–free component would be N_e= 1456 ± 86 cm⁻³, T_e= 10⁴ K. We used a distance of 420 pc, and a nebular diameter of 28.5 arcsec (as compiled by Acker et al. 1992).

Figure 5

SEDs of the PNe with the largest 250-GHz deficits. We plot flux density in Jy (y-axis) as a function of frequency in GHz (x-axis). The CBI and SIMBA data points are indicated by circles. The solid line is a single uniform-slab free–free component, drawn to emphasized deviations from free–free spectra, and fitted to the ≤31-GHz data with the following parameters: NGC 7009, ⁠, ⁠; NGC 6572, ⁠, ⁠; NGC 3242, ⁠, ⁠; NGC 3918, ⁠, ⁠.

8.1.2 NGC 6572

The SED of NGC 6572 in Fig. 5 is again very similar to NGC 6369. We adopt the expansion parallax distance of 1.49 ± 0.62 kpc (Hajian, Terzian & Bignell 1995), and a nebular diameter of 8 arcsec (as observed with the Very Large Array by Zijlstra, Pottasch & Bignelli 1989).

We could fit the data with an attenuated and broken synchrotron spectrum, with a low-frequency spectral index of −0.4, and a high-frequency index of −1.0. The bulk free–free component has N_e= 10056 ± 162 cm⁻³ and T_e= 10⁴ K. We used a free–free screen at 50 K, with 〈N_e(screen)〉= 998 ± 248, and power-law opacity distribution with index γ= 5.2. The sub-mm dust has a temperature of 76 K.

The data from Hoare et al. (1992) confirm our calibration of the SIMBA flux densities. It can be seen in Fig. 5 that the SIMBA data point, at 250 GHz, is low but agrees with the Hoare et al. (1992) data within the uncertainties. The rising sequence with frequency could perhaps be due to the Rayleigh–Jeans tail of sub-mm dust extending to 250 GHz. We have also plotted on Fig. 5 data points from the preliminary SEST measurements (assigned exaggerated error bars).

8.1.3 NGC 3242

NGC 3242 has expansion parallax distance of 420 ± 160 pc (Hajian et al. 1995), and a nebular diameter of about 25 arcsec. We could fit the SED of NGC 3242 with the same synchrotron component as for NGC 6369, except the core free–free emission is constrained to the 250-GHz flux density, has N_e= 2401 cm⁻³, and is assigned a temperature of 20 000 K. A free–free screen at 20 K has an average density of 200 cm⁻³, and a power-law distribution of opacities with γ= 10. A single power-law fit at 5–31 GHz gives an exponent of {−}0.290^{−0.237}_{−0.340{/pm 0.051, which is different from the free–free α^{31}₅ value by more than 3σ.

8.1.4 NGC 3918

The high excitation PN NGC 3918 has been studied in depth by Clegg et al. (1987). We adopt their preferred distance of 1.5 kpc, and a nebular angular size of 15 arcsec. The SED of NGC 3918 shown in Fig. 4 can be fitted by a broken synchrotron spectrum, with spectral indices of −0.4 and −1.0. The free–free screen is at 50 K, and has an average density of 〈N_e(screen)〉= 371 ± 91. The bulk of the 250-GHz flux density is due to the traditional free–free emission at 10⁴ K and N_e= 2342 ± 306 cm⁻³.

8.2 Outliers

In this section we discuss the objects that we have considered as outliers from the correlations described in Section 5. They differ in physical characteristics and radio spectra from the rest of the sample.

8.2.1 M2-9

M2-9 is a peculiar object. Its ionizing central star is undergoing sustained mass loss, as revealed by broad Hα wings (e.g. Swings & Andrillat 1979; Arrieta & Torres-Peimbert 2003). As reported by Kwok et al. (1985), the radio SED of M2-9 is that of an optically thick stellar wind, and is in marked contrast with the other objects studied here.

In Fig. 6 we model the radio spectrum of M2-9 with a power-law distribution of opacities, following Calabretta (1991), except the opacity law we used is τ(ω) =δω^−γ. As explained in Calabretta (1991), the detailed spatial distribution of opacities is irrelevant to the observed free–free flux density. Only the distribution of opacities as a function of a single solid angle variable is required. In general the opacity will be a function of some 1D parameterization of solid angle, ω. As an example definition of ω, Calabretta (1991) pixelates the sky emission of the nebulae and orders pixels by their opacities. A natural choice for ω in spherically symmetric nebulae is the nebular radius. In the case of the optically thick stellar wind in M 2-9, the bulk of the emission at cm to mm wavelengths comes from the wind, with a power-law profile: ⁠. A clumpy free–free screen can also be modelled with a similar opacity profile as M 2-9, except the clumps are scattered across the nebula. Larger values of γ correspond to smaller fractions of the nebular pixels containing most of the opacity. Power-law profiles allow for a smoother transition to the optically thick frequency regime than uniform slabs (a proof is given in Calabretta 1991).

Figure 6

SEDs of outliers. We plot flux density in Jy (y-axis) as a function of frequency in GHz (x-axis). The CBI and SIMBA data points are indicated by circles. M 2-9: the components that sum up the model spectrum (in solid line) are a modified blackbody (dotted line), an optically thick stellar wind (dashed line) and a uniform slab (dash–dotted). NGC 6302: we use a core–halo model to fit the observed spectrum. The core component is shown by the dashed line (⁠⁠, ⁠), and the halo component is shown by the dashed-dotted line (⁠⁠, ⁠, with a diameter of 20 arcsec). This model is not unique and only serves the purpose of indicating a physical reference. Vy 2-2: line styles follow from the NGC 6302 spectrum. The flux densities are extracted from Purton et al. (1982) and Altenhoff et al. (1994).

The emergent flux density is

5

where Ω_N is the nebular solid angle, ω=Ω/Ω_N is a solid angle variable (for instance angular distance from the central star), ω_u≪ 1 is fixed by requiring and

6

We compiled the SED of M2-9 with the data in Hoare et al. (1992), our own CBI and SEST measurements, the 5- and 14-GHz data from Milne & Aller (1982) and the 2.7-GHz flux density from Calabretta (1982). We obtain a fit to the radio-IR spectrum of M2-9 with a superposition of a uniform slab nebula to account for the flattening at low frequencies, a power-law distribution with δ= 1.6 × 10⁻⁹, γ= 2.06 to model the optically thick core, an average core density of 〈N_e〉= 2093 cm⁻³, and a temperature T_e= 7700 K. The best-fitting average emission measure entering in the average opacity 〈τ_ν〉 is 〈EM〉= 1.3 × 10⁶ pc cm⁻⁶, with a total nebular diameter of 1 arcmin. We take the agreement between the SIMBA flux density of M2-9 with the data from Hoare et al. (1992) and Altenhoff, Thum & Wendker (1994), as well as with the model SED, as a confirmation of our data calibration and reduction procedure.

8.2.2 NGC 6302

NGC 6302 is the highest excitation PN known, and is probably derived from the upper mass range of PN progenitors (e.g. Casassus, Roche & Barlow 2000). As such it has been the focus of much study. The ionized nebula is surrounded by copious amounts of neutral and molecular material (Rodriguez & Moran 1982).

In Fig. 6 we summarize the available data set on NGC 6302. The flux densities have been extracted from Terzian, Balick & Bignell (1974), Milne & Webster (1979), Calabretta (1982), Milne & Aller (1982), Rodriguez et al. (1985), Zijlstra et al. (1989), Hoare et al. (1992) and Condon & Kaplan (1998). We adopt the expansion parallax distance of 1.6 ± 0.6 kpc (Gomez, Rodriguez & Moran 1993), and a core nebular diameter of 10 arcsec (Zijlstra et al. 1989), surrounded by a halo extending up to ∼1 arcmin. As in the case of M2-9, we find that the SIMBA data point agrees, within the error bars, with the data from Hoare et al. (1992).

We fit the SED of NGC 6302 with a core-halo model, as in Rodriguez et al. (1985), including the dust temperature as a variable. The IRAS 100 and 60 μm data alone give T_dust= 46 K, if the dust emissivity is β=+2. Including the United Kingdom Infrared Telescope (UKIRT) and James Clerk Maxwell Telescope (JCMT) data from Hoare et al. (1992), as well as our SIMBA point, gives a dust temperature of 42 K. We believe there is room for a colder dust component that could be present at ∼200 GHz.

The same situation occurs in NGC 6537, often considered as a slightly lower excitation twin of NGC 6302. The Rayleigh–Jeans tail of the sub-mm dust could contribute to the 250 GHz flux.

NGC 6302 is an important cross-check on the SIMBA calibration scale. Raising the SIMBA flux densities by 50 per cent would eliminate the cm-excess, and result in a mm-excess in some objects. However, the SIMBA data on NGC 6302 is consistent with the literature. Raising the NGC 6302 250-GHz flux densities by 50 per cent could be inconsistent with Hoare et al. (1992) by 5σ.

8.2.3 Vy 2-2

Vy 2-2 is a very young low excitation PN (e.g. Liu et al. 2001), with strong 10 μm silicate emission at ∼365 K (Aitken & Roche 1982). Faint, compact and opaque radio emission at the nucleus of Vy 2-2 is resolved into a 0.4-arcsec diameter ring at 22 GHz (Christianto & Seaquist 1998), and is surrounded by a faint 25 arcsec Hα halo.

Its distance is greater than 1 kpc (Christianto & Seaquist 1998), and traditional indicators place it at 7.9 kpc (Zhang 1995). Since Vy 2-2 is a young object, we adopt the distance of 3.9 kpc, given by the IRAS fluxes and the tip-of-the-AGB luminosity function (Casassus et al. 2001).

Vy 2-2 is very IR bright, while optically thick and faint at radio frequencies. The large uncertainty on its SIMBA flux density left it out of Fig. 3. Its IRAS 60-μm to 5-GHz ratio is 0.86 × 10⁻³, three–four times higher than even the largest values in Fig. 3 or in NGC 6302.

The radio SED of Vy 2-2 shown on Fig. 6 can be fitted with a simple free–free model consisting of a uniform slab, with T_e= 8522 ± 660, N_e= 3 × 10⁵± 4 × 10³ cm⁻³ and a depth of 0.01 pc. It can be appreciated this SED is very different from the rest of the PNe considered in this work.

8.2.4 NGC 2440

NGC 2440 is more evolved than the rest of our sample, and its mid-IR flux is due to ionic lines. Dust in this object radiates only in the far-IR, and is presumably located further from the central star than in the other PNe. The SIMBA data point is unfortunately very uncertain, and so is the <31 GHz data. We mention it here because the low-frequency data seems to peak at 5 GHz and drop with frequency steeper than allowed by free–free emission at low frequencies, in a fashion that is reminiscent of the 250 GHz deficit evident in NGC 6369.

9 CONCLUSIONS

We have found an inconsistency with the free–free paradigm for the radio continua of prototypical PNe, in the form of an excess over free–free emission at 31 GHz when compared to the 250-GHz flux densities. The fraction of the 31 GHz continuum not due to free–free emission, averaged over the 10 objects with less than 30 per cent uncertainties, is 51 ± 3 per cent, with a weighted scatter of 11 per cent. The cm-excess is present in all objects that are optically thin at 31 GHz, with the exception of NGC 6302 and 6537, whose sub-mm dust tails reach 250 GHz. The most significant detections are in NGC 6369, 7009, 3242 and 6572.

The evidence for the cm-excess stems mostly from the 250-GHz data acquired with SIMBA. Raising the SIMBA calibration scale by ∼50 per cent would bring some of the 250-GHz levels in agreement with free–free, but would also render the data on NGC 6302 and M 2-9 incompatible with the literature. Additional evidence can be found in the 5–31 GHz spectrum of NGC 3242, and in the average ⁠, which are steeper than free–free.

The 31-GHz excess follows a loose correlation with the IRAS flux densities. Yet the examination of individual objects shows that spinning dust alone cannot explain the low-frequency data, at least with the models currently available. Needle extinction at 250 GHz would require extreme temperatures. A synchrotron component is a possibility, if it is attenuated at low frequencies by a cold free–free screen (<500 K). The cm-wave–IR correlation could be due to synchrotron radiation originating in compact and dusty inclusions exposed to the fast stellar winds observed in PNe.

1

The SEST is situated on La Silla in Chile.

2

We assigned 10 per cent uncertainties to the CBI data point to allow for differences in the photometric extraction procedures: the CBI flux densities are obtained by Gaussian fitting, while Coble et al. (2003) used a circular photometric aperture encompassing most of the emission.

3

Synchrotron radiation from the nebular knots taps energy from the fast stellar wind, which in these hot stars is accelerated by radiation pressure coupled to the gas through ionic lines (e.g. Pauldrach et al. 1988). The radiation–fast wind coupling must be very effective since observations give P_rad∼P_wind.

4

This cold free–free component cannot, of its own accord explain the cm-excess, because the lowest possible free–free index is α=−0.175 at 1000 K (see Fig. A4), and because the nebulae, as cold blackbody radiators, would have to be much larger to reach the observed <5 GHz levels.

5

The LFL was originally derived by Elwert, as explained by Oster (1961)

6

We have not included the contribution of radio recombination lines (RRLs) to broadband data. RRLs could be important at ∼300 GHz for ⁠.

We are grateful to Pat Roche for useful discussions and for reminding SC about helium, to Eli Dwek for an e-mail exchange providing extensive information on his needle models, and to Wouter Vlemmings for feedback on the nebular magnetic fields. The Swedish-ESO was operated jointly by ESO and the Swedish National Facility for Radioastronomy, Onsala Space Observatory at Chalmers University of Technology. SC acknowledges support from FONDECYT grant 1060827, and from the Chilean Centre for Astrophysics FONDAP 15010003. We gratefully acknowledge the generous support of Maxine and Ronald Linde, Cecil and Sally Drinkward, Barbara and Stanely Rawn, Jr, Fred Kavli and Rochus Vogt. Part of the research described in this paper was carried out at the Jet Propulsion Laboratory, California Institute of Technology, under a contract with the National Aeronautics and Space Administration. CD thanks B. and S. Rawn, Jr, for funding a fellowship at the California Institute of Technology for part of this work.

NOTES ADDED IN PROOF

Gayet (1970) gives a hierarchy of approximations to calculate the hydrogenic free–free Gaunt factors.

The frequency of the SEST single-pixel bolometer data described in Section 7 is the same as for SIMBA, or 250 GHz, and varies by about 10 GHz as a function of atmospheric opacity. The frequency of the SEST heterodyne data is 230 GHz.

Appendix

APPENDIX A: THERMAL SPECTRAL INDICES

The free–free emissivity ε for electron–ion encounters for a dilute plasma (at frequencies much larger than the plasma frequency; Scheuer 1960) can be written as (e.g. Beckert, Duschl & Mezger 2000, per unit solid angle)

(A1)

where N_i is the ion density and g_ff(ν, T_e, Z_i) is a frequency- and temperature-dependent Gaunt factor.

A1 Exact Gaunt factors for free–free emission

We searched the literature for approximations to the exact g_ff(Z_i, ν), and found errors comparable to the level of precision in our spectral analysis (of ∼5–10 per cent). We thus chose to calculate exact Gaunt factors, using the formulae and procedures described in Grant (1958), which have also been reproduced in Beckert et al. (2000) and Muñoz (2006).

The form we used for the non-relativistic Gaunt factor in encounters between electrons and bare nuclei is that given by Menzel & Pekeris (1935). Velocity-dependent Gaunt factors for individual encounters involve the evaluation of complex hypergeometric functions ₂F₁(a, b, c, z). We used the hypergeometric series when possible (for |z| < min (0.1, |0.1 c/(ab)|, as prescribed by Beckert et al. 2000), and integrated the hypergeometric equation otherwise. We used the GNU Scientific Library routines for ODEs, finding that the Bulirsch–Stoer method gave best results. We used the prescription in Beckert et al. (2000) as starting point for the ODE integration. For |z| > 1 we used the ‘linear transformation formulae’ from section 15.3 of Ambramowitz & Stegun (1964) to keep the complex argument z within the unit circle.

Our hypergeometric routines are robust for any frequency for relatively fast electrons, when the parameter η_i=−αZ (c/v_i) < 1, in the notation of Grant (1958). However, the ODE integration crashed for large a and b, corresponding to slow incident electrons, or when the parameter η_i > 100. Therefore we tried the hypergeometric routines in Press et al. (1996, NR), edited to double precision and with the starting point of Beckert et al. (2000), finding the same numerical limitation as for GSL.

The domain of slow incident electrons corresponds to the so-called ‘classical’ limit, when the de Broglie wavelength is small compared to Z e²/(m v²_i), a measure of size of the Coulomb field (e.g. Brussaard & van de Hulst 1962). In what follows we will refer to the ‘classical’ limit as the WKB limit, and use the term classical for the classical mechanics calculation leading to the Gaunt factors used by Scheuer (1960) and Oster (1961). The WKB limit involves Hankel functions with imaginary arguments, and are rather difficult to evaluate.

Fortunately, a good and simple approximation for η_i > 100 can be obtained in the following way. When the product photon has a frequency close to the recombination (bound–free) limit, the WKB Gaunt factor reduces to a very simple expression given by equation (14) from Grant (1958). In the limit of low frequencies, the incident electron follows a classical trajectory unperturbed by the emission of the radio photon, and the Oster/Scheuer Gaunt factors are applicable. At frequencies between 0.1 and 1000 GHz η_i varies between 0.3 and 2500 when T_e= 1–50 000 K. Fig. A1 shows the free–free Gaunt factor for individual encounters involving incident electrons with η_i of −1, −10 and −100, as a function of frequency. Such values of v_i correspond to the rms velocities at temperatures of ∼10⁵, 10³ and 10 K, respectively. The maximum frequency corresponds to the recombination limit ν_g. We compare the exact calculation, in asterisks, with the low-frequency limit (hereafter LFL, as given in Beckert et al. 2000, their equation 23, but replacing the absolute value by a minus sign),5 the WKB and the classical Gaunt factor. Oster (1961) shows that the LFL converges to the classical limit for large η_i. It can be seen that as ν approaches ν_g the exact calculation scatters about the expected WKB value, especially for large η_i. This is because the ODE integration becomes numerically unstable, and its results cannot be trusted, as ν approaches ν_g. The ODE integration altogether crashes for even larger values of η_i. However, there is always a regime where an approximation can be substituted for the exact calculation, and this value is the largest between the classical and WKB Gaunt factors. Thus, we used the WKB Gaunt factor when ν > 0.01ν_g and η_i < −1, and the maximum between the classical and WKB Gaunt factors when η_i < −100. This is different from the choice of Beckert et al. (2000), who always used the LFL instead of the WKB approximation. The relative difference between the use of LFL or WKB is important only at very low temperatures, and is within ∼10⁻³ when T_e > 1000 K.

Figure A1

The free–free Gaunt factor for individual electron–proton encounters with incident velocities corresponding to η_i=−1, −10 and −100. We plot the exact calculation (with ‘x’, ‘+’‘*’ symbols, respectively), the LFL (dotted line), the classical (dashed line) and WKB Gaunt factors (solid line). The recombination limit Gaunt factor (equation 6 from Brussaard & van de Hulst 1962) is indicted by a circle.

The Maxwellian integration for the total rate coefficients (a.k.a. the average coefficients), equation (24) in Grant (1958), was carried out using the PDL interface to the GNU Scientific Library. In order to confirm our calculation we did the exercise of using the LFL of Beckert et al. (2000) for large |η_i | (i.e. not the WKB limit), and found that we obtain the same total Gaunt factors within a relative error of ∼10⁻⁵ as the calculations in Muñoz (2006), who used the mathematica package to directly evaluate the Menzel & Pekeris (1935) Gaunt factors (and following Beckert et al. 2000).

The use of common approximations to the radio-frequency Gaunt factors involve errors ∼5 per cent in the 1–300 GHz frequency range for traditional nebular temperatures. In Fig. A2 we have compared the exact Gaunt factors with the approximations in Scheuer (1960), Oster (1961), and the cloudy photoionization package (version c06.02c, last described by Ferland et al. 1998), which uses free–free Gaunt factors based on Hummer (1988). Note that the Altenhoff et al. (1960) free–free opacities, with the a(ν, T) factor correction from Mezger & Henderson (1967), are equal to free–free opacities calculated directly from the Oster (1961) Gaunt factors.

Figure A2

Comparison between various Gaunt factor formulae at ⁠. We plot the perceptual difference between the corresponding Gaunt factors g for Z=2 (dotted), the Oster (1961) (dashed) and cloudy (solid) approximations, and the exact value with Z=1 from Menzel & Pekeris (1935), g_MP. For reference we obtain g_MP(ν= 250 GHz, T_e= 7000 K, Z= 1) = 2.7954. The curve for Z=2 has been divided by 2.

For practical applications, in order to speed up the optimizations involved in fitting the PN SEDs, we gridded the exact free–free Gaunt factors in the [log(T_e), log(ν)] space. Fig. A3 shows the corresponding Gaunt factor images. The limit for a cold plasma, at high frequencies, is 1, as expected. Note that the Gaunt factors given by Oster (1970) tend to 0 for cold plasmas because they include the exp[−hν/(kT)] term, which is usually associated with the emissivity (equation A1) rather than with the Gaunt factor.

Figure A3

Exact (non-relativistic) Gaunt factors for Z= 1, as a function of log (ν/GHz) in x-axis and log(T_e/K) in y-axis.

A2 Free–bound emission

Since we are considering thermal emission from cold plasma, the possibility of radio free–bound emission must be taken into account.6 We use the formulae compiled in Brussaard & van de Hulst (1962), taking free–bound Gaunt factors of 1. Fig. 8 of Brussaard & van de Hulst (1962) shows that, in the radio–sub-mm regime, the free–bound Gaunt factors are essentially 1 for recombinations to any principal quantum numbers.

A3 Contribution from He and representative PN spectral indices

The free–free and free–bound emissivities depend on temperature. In Fig. A4 we summarize the variations of α³¹₅(H⁺), α³¹₅(He⁺⁺), α²⁵⁰₃₁(H⁺) and α²⁵⁰₃₁(He⁺⁺), in the case of optically thin emission and for a range of electron temperatures T_e.

Figure A4

Exact (non-relativistic) optically thin thermal spectral indices for encounters between electrons and bare H and He nuclei, evaluated for a range of temperatures, and two pairs of frequencies: 31–250 and 5–31 GHz.

The Gaunt factors for the He free–free continuum should be almost identical to the H values in the domain where the wavelength of the electron incident on He⁺ has a wavelength much larger than the He⁺ Bohr radius. Thus, we can neglect differences with He up to a temperature T_e= (2h/a_○)²/(3k m_e) ≈ 17 × 10⁶ K. As explained above, the Gaunt factors for the He⁺ free–free continuum are slightly different from H. Considering that the most He rich nebulae (of Peimbert type I) reach He abundance by number of 18 per cent relative to H, a fraction of which is doubly ionized, we can conservatively assume a contribution from the He⁺ continuum of up to 20 per cent, and a minimum of 0. This translates into the spectral index uncertainties given in Table 1, where we have assumed an average abundance of He⁺⁺ of 10 per cent relative to H⁺.