Measurements of the Intensity and Polarization of the Anomalous Microwave Emission in the Perseus molecular complex with QUIJOTE

The new data presented in this article were acquired with the QUIJOTE experiment. QUIJOTE is a collaborative project that consists of two telescopes and three polarimeter instruments covering respectively the frequencies 10-20, 30 and 40 GHz, and located at the Teide Observatory (2400 m a.s.l.) in Tenerife (Spain). The main science driver of this experiment is to perform observations of the CMB polarization to constrain the B-mode signal down to $r=0.05$, and to characterize the polarization of low-frequency foregrounds, mainly synchrotron emission and the AME, so that this signal can be removed from the primordial maps to a level that will permit reaching the previous level of $r$. The two QUIJOTE telescopes are based on an offset crossed-Dragone design, with projected apertures of $2.25$ m and $1.89$ m for the primary and secondary mirrors, and provide highly symmetric beams (ellipticity $>0.98$) with very low sidelobes ($\leq-40$ dB) and polarization leakage ($\leq-25$ dB). The first instrument to be fielded on the QUIJOTE first telescope is a multi-frequency instrument (MFI) with 4 horns covering the frequency range $10-20$ GHz, and with angular resolutions close to one degree. These detectors are fitted with MMIC low-noise amplifiers (noise temperature better than 10 K), and use stepped polar modulators to measure the polarization of the incoming radiation, providing instantaneous sensitivities of $\approx 650~\mu$K$~s^{1/2}$ in four individual frequencies with nominal figures: 11, 13, 17 and 19 GHz. The median integrated PWV above the observatory is $\approx 4$ mm, giving a zenith atmosphere temperature of $\approx 2$ K at 11 GHz and $\approx 5$ K at 19 GHz. The MFI saw first light on November 2012 and since then it is performing routine observations of different Galactic and cosmological regions. The second instrument consists of 31 polarimeters at 30 GHz (TGI, thirty-gigahertz instrument), and is based on the same design of the MFI except that the polarization modulation is achieved electronically through phase switches. This instrument will be commissioned during 2015. Finally, the third instrument is planned to have 31 polarimeters at 40 GHz (FGI, forty-gigahertz instrument). Using the TGI and the FGI, which will provide instantaneous sensitivities of $50~\mu$K s^1/2, we plan to survey an area of $3000$ deg² down to a projected sensitivity of $\leq 1~\mu$K/beam. A more detailed description of the technical and scientific aspects of this project can be found in Rubiño-Martín et al. (2012a) or in Rebolo et al. (in preparation).

All the polarization data that will be used in this article comes from the QUIJOTE experiment. However, in order to obtain the full spectral energy distribution (SED) of G159.6-18.5, from which the residual AME fluxes will be inferred, we use ancillary data from other experiments. In the low-frequency range we use the Haslam et al. (1982) map¹¹1We use the map supplied by Platania et al. (2003). at 0.408 GHz, the Berkhuijsen (1972) map²²2We projected the map downloaded from http://www.mpifr-bonn.mpg.de/survey.html into HEALPix pixelization. at 0.820 GHz and the Reich & Reich (1986) map at 1.4 GHz.

At $10.9$, $12.7$, $14.7$ and $16.3$ GHz, similar frequencies to QUIJOTE, we use data from the COSMOSOMAS experiment (Watson et al., 2005). In order to minimize the $1/f$ noise, the data from this experiment was filtered by the suppression of the first seven harmonics in the FFT of the circular scans, which results in a flux loss on large angular scales. For this reason, the comparison with other experiments that preserve all the angular scales is not straightforward. It is necessary to account for the flux lost, as it was done in Planck Collaboration et al. (2011). The fluxes presented in that paper are already corrected, so we directly take those fluxes.

We also use data from the 9-year release of the WMAP satellite³³3Downloaded from the LAMBDA database, http://lambda.gsfc.nasa.gov/. (Bennett et al., 2013) to provide flux estimates at frequencies 23, 33, 41, 61 and 94 GHz. Recent data from the first release of the Planck mission⁴⁴4Downloaded from the Planck Legacy Archive, http://www.sciops.esa.int/index.php?project=planck&page=
Planck_Legacy_Archive. (Planck Collaboration et al., 2014c) cover the frequencies 28, 44, 70, 100, 143, 217, 353, 545 and 857 GHz. We also download the released Type 1 CO maps (Planck Collaboration et al., 2014e), which are used to correct the 100, 217 and 353 GHz frequency maps from the contamination introduced by the CO rotational transition lines (1-0), (2-1) and (3-2), respectively. Finally, in the far-infrared spectral range we use Zodi-Subtracted Mission Average (ZSMA) COBE-DIRBE maps (Hauser et al., 1998) at 240 $\mu$m (1249 GHz), 140 $\mu$m (2141 GHz) and 100 $\mu$m (2998 GHz).

Intensity and polarization fluxes will be calculated by applying an aperture photometry integration on the maps. This is a well known and widely used technique in this context (López-Caraballo et al., 2011; Dickinson et al., 2011; Génova-Santos et al., 2011), consisting of integrating temperatures of all pixels within a given aperture, and subtracting a background level calculated in an external ring. Instead of the mean, following Planck Collaboration et al. (2011), we chose to use the median of all the pixels in the external ring as the estimate of the background level. The median is a better proxy for the real level in cases of strongly variable backgrounds with many outlier pixels. The flux is then given by

where $n_{1}$ is the number of pixels in the aperture, and $T_{i}$ and $T_{j}$ represent respectively the pixel thermodynamic temperatures in the aperture and in the background annulus. The median is calculated over the $n_{2}$ pixels in this annulus. The function $a(\nu)$ gives the conversion factor from temperature to flux,

where $h$ and $k_{\rm b}$ are the Planck and Boltzmann constants, $c$ the speed of light, $T_{\rm cmb}=2.725$ K (Mather et al., 1999) the CMB temperature and $\Omega_{\rm pix}$ the solid angle subtended by each pixel.

The determination of the error associated with the previous estimate is crucial for the results of this paper. In a hypothetical case of perfect white uncorrelated noise, it could easily be estimated through:

where $\sigma$ represents the error of each pixel.

However, in QUIJOTE the instrument noise is correlated due to the presence of $1/f$ residuals, and also the galactic background fluctuations introduce, mainly in intensity, an important contribution to the error which is correlated on the order of the beam size. Ideally we should then use the covariance matrices of the instrument and background noises. The former can be extracted through a characterization of the $1/f$ noise spectrum, however the latter is difficult to determine. Instead, in the previous equation we can introduce in the denominator the number of independent pixels in the aperture and in the ring, which we will denote respectively as $n_{1}^{\prime}$ and $n_{2}^{\prime}$. The pixel variance will be calculated from the pixel-to-pixel standard deviation of all the pixel temperatures $T_{j}$ in the background, $\sigma(T_{j})$. Obviously, the standard deviation of the pixels in the aperture would be biased by the presence of the source. On the contrary, the standard deviation of the pixels in the ring gives a reasonable estimate of the contributions of the background and of the instrumental noise to the true error. Therefore, the final equation that we will use to estimate errors in this paper reads as:

In the case of the error being completely dominated by the background, then we could use for $n_{1}^{\prime}$ and $n_{2}^{\prime}$ the number of beams in the aperture and in the background, $n_{1}^{\rm b}$ and $n_{2}^{\rm b}$. However, while being particularly strong in intensity, the background fluctuations from the Galactic emission are not so important in polarization. For this reason, in this case the relative contribution from the $1/f$ residuals to the uncertainty is significant. To quantify this, we selected 20 random positions around our source, G159.6-18.5, and performed flux integration on those positions using the same aperture and ring sizes. The standard deviation of those fluxes gives a reasonable estimation of the true noise of our flux estimate. From this analysis we determined that for intensity $n_{1,2}^{\prime}=n_{1,2}^{\rm b}$, while for polarization $n_{1,2}^{\prime}=2n_{1,2}^{\rm b}$. This is what we will use in our estimation of the flux errors.

In cases of low signal-to-noise fluxes, or when placing upper limits on the polarized flux $P$, as it will be our case, it is necessary to de-bias the fluxes derived from the aperture photometry integration. This requirement comes from the fact that the posterior distribution of the polarized intensity $P$ does not follow a normal (Gaussian) distribution. Furthermore $P$ is a quantity that must always be positive, and this introduces a bias into any estimate. For any true $P_{0}$ we would expect to measure on average a polarization $P>P_{0}$. In order to get the de-biased fluxes, $P_{\rm db}$, from the measured ones, $P$, we choose the Bayesian approach described in Vaillancourt (2006) and in Rubiño-Martín et al. (2012b), consisting of integrating the analytical posterior probability density function over the parameter space of the true polarization. The same posterior can not be used for the polarization fraction, $\Pi=P/I\times 100$, as it follows a different distribution. As, to our knowledge, there is not in the literature any analytical solution for the posterior distribution of $\Pi$, we numerically evaluate this function by applying Monte-Carlo simulations. This approach has already been carried out in López-Caraballo et al. (2011) and in Dickinson et al. (2011).

In Fig. 1 we show the intensity map at 11 GHz resulting from combining 149 h of observations, where emissions from G159.6-18.5, the California nebula (NGC1499) and the 3C84 quasar are clearly visible. More detailed $I$, $Q$ and $U$ maps at our four frequencies around the position of G159.6-18.5 are shown in Fig. 2. The $Q$ and $U$ maps on this source are consistent with zero polarization, and therefore upper limits on the polarized intensity will be extracted in section 3.3. Some stripping is clearly visible in the maps, which is produced by the presence of regions with a higher noise due to a lower integration time per pixel and, to a lesser extent, to $1/f$ residuals. The inhomogenities in the coverage (integration time per pixel) maps are caused by the separation of the horns in the focal plane, which leads the sky coverage to be different when we observe the field before or after crossing the local meridian. In the $Q$ and $U$ maps at 11 and 13 GHz of Fig. 2, the two orthogonal stripes with clearly higher noise correspond to regions with integration of $\sim 3-7$ s/pixel (pixel size 6.9 arcmin). By comparison, in the central region inside the circle where we perform the aperture photometry, the integration time is $\sim 30-35$ s/pixel, resulting in a lower pixel-to-pixel dispersion of the data.

An important consistency test, that may reveal the presence of systematics and other spurious effects, is obtained through jackknife maps. We have uniformly split our full dataset in two halves in such a way that the maps of number of hits associated to these two halves are as similar as possible. The differences of the two halves divided by two, for the intensity and polarization maps at our four frequencies, are shown in Fig. 3. As expected, the intensity emission coming from G159.6-18.5 is consistently cancelled out in these maps. While a similar striping pattern to the maps in Fig. 2 is seen, these maps are dominated by instrumental noise. This is confirmed by the noise values shown in Table 1, where we compare the pixel-to-pixel RMS calculated in two different regions: the external ring that we will use for background subtraction when calculating the intensity and polarization fluxes, and a region of very low sky emission enclosed by the dashed circle represented in Fig. 3. The RMS values in $Q$ and $U$ are similar in the original maps and in the sum and difference of the two halves. They are typically $\sim 250~\mu$K/pixel (pixel size: 6.9 arcmin) or $\sim 25~\mu$K/beam (beam size: $1\sim$ degree). In the case of the $I$ maps, the RMS in the background ring are higher in the sum than in the difference because of the emission of the source. In the circle with low sky signal the values in the sum and in the difference maps are very similar.

The last column of Table 1 show the average RMS levels in the $Q$ and $U$ maps normalized by the integration time per pixel. As the number of hits per pixel is very inhomogeneous, to calculate these numbers we have made a realization of Gaussian noise in which we assign to each pixel a noise proportional to $t_{\rm pix}^{-1/2}$, where $t_{\rm pix}$ is the integration time per pixel, and then calculate the pixel-to-pixel RMS. The amplitudes of the white-noise in the spectra of the time-ordered-data range between 898 $\mu$K s^1/2 (at 16.7 GHz) and 1228 $\mu$K s^1/2 (at 11.2 GHz)⁵⁵5Note that in this article we are using data from only two out of the four horns of QUIJOTE. Then, the global sensitivities of the experiment are a factor $\sqrt{2}$ better, i.e. $\approx 650~\mu$K s^1/2, the number quoted in Section 2.1. $1/f$ residuals make the noises calculated on the maps only slightly higher, typically by a factor $\sim 15\%$, confirming our previous statement that these maps are dominated by white (Gaussian) noise.

Another important consistency test for the presence of systematics, and in particular for the $I$ to $Q/U$ polarization leakage, is to confirm that our polarization maps are consistent with noise in the position of unpolarized sources. This verification is provided by the nearby California HII region, which is also covered by our observations. As any standard HII region, it is dominated by free-free emission at the QUIJOTE frequencies, which is known to be practically unpolarized. In Fig. 4 we show intensity and polarization maps towards the California region, showing that the $Q$ and $U$ maps are consistent with zero polarization.

As was indicated in section 2.2, we take the COSMOSOMAS fluxes for G159.6-18.5 from Planck Collaboration et al. (2011), which have already been corrected for the flux loss caused by the filtering of COSMOSOMAS data. While here we will use aperture photometry to derive our fluxes, in Planck Collaboration et al. (2014a) they were obtained by fitting the amplitude of an elliptical Gaussian with a fixed size of $1.6^{\circ}\times 1.0^{\circ}$ (FWHM). In a first-order approximation, the fluxes obtained through Gaussian fitting will be equivalent to those obtained from aperture photometry using a given aperture size. Therefore, in order to get a reliable intensity SED we choose a size for the aperture that gives the most similar fluxes to those presented in Planck Collaboration et al. (2011) for the Haslam et al. (1982), Berkhuijsen (1972), Reich & Reich (1986), WMAP , Planck and DIRBE maps (it must be noted that the WMAP and Planck maps used in this work correspond to a different release to that used in Planck Collaboration et al. (2011), but this will have a negligible effect). After trying different values, we found that a radius of $1.7^{\circ}$ gives the best agreement, with a very low reduced chi-squared of $\chi^{2}_{\rm red}=0.098$. The median of the background is computed in an external ring between $1.7^{\circ}$ and $1.7\sqrt{2}^{\circ}$, which has the same area as the aperture. The derived fluxes in the QUIJOTE maps and in the other ancillary maps are listed in Table 2.

The final SED is depicted in Fig. 5, where the presence of AME clearly shows up at intermediate frequencies as an excess of emission over the other components. The intensities derived from these QUIJOTE observations trace, for the first time after the original measurements of the COSMOSOMAS experiment (Watson et al., 2005), the downturn of the spectrum at frequencies below $\sim 20$ GHz, as it is predicted by spinning dust models. In total, 13 data points are dominated by AME: the four QUIJOTE points, the four COSMOSOMAS points, the WMAP $22.8$, $33.0$ and $40.7$ GHz frequencies and the Planck $28.4$ and $44.1$ GHz frequencies. We perform a joint multi-component fit to all the data, consisting of free-free emission, which dominates in the low-frequency tail, spinning dust, which dominates the intermediate frequencies, a CMB component, and thermal dust, which dominates the high-frequency end. As it was done in Planck Collaboration et al. (2011), to avoid possible CO residuals we exclude from the fit the 100 GHz and the 217 GHz values. We fix the spectrum of the free-free using the standard formulae shown in Planck Collaboration et al. (2011), with a value for the electron temperature typical of the solar neighbourhood, $T_{\rm e}=8000$ K, and fit for its amplitude, which is parameterized through the emission measure ($EM$).

Following Planck Collaboration et al. (2011), we consider a high-density molecular phase and a low-density atomic phase which produce spinning dust emission, and fit their respective amplitudes, which are given by the hydrogen column densities $N_{\rm H}^{\rm mol}$ and $N_{\rm H}^{\rm at}$. The CMB amplitude is denoted by $\Delta T_{\rm cmb}$, and would correspond to the average of the primordial CMB fluctuations in the aperture. Finally, the thermal dust is modelled as a single-component modified blackbody curve, $\tau_{\rm 250}(\nu/1200~{\rm GHz})^{\beta_{\rm d}}B_{\nu}(T_{\rm d})$, which depends on three parameters: the optical depth at $250~\mu$m ($\tau_{\rm 250}$), the emissivity spectral index ($\beta_{\rm d}$) and the dust temperature ($T_{\rm d}$). Therefore, we jointly fit 7 parameters to all the data points: $EM$, $N_{\rm H}^{\rm mol}$, $N_{\rm H}^{\rm at}$, $\Delta T_{\rm cmb}$, $\tau_{\rm 250}$, $\beta_{\rm d}$ and $T_{\rm d}$.

To define the spinning dust spectra of the molecular and atomic phases we use the spdust.2 code⁶⁶6http://www.tapir.caltech.edu/$\sim$yacine/spdust/spdust.html (Ali-Haïmoud et al., 2009; Silsbee et al., 2011). Initially we use the same values as in Planck Collaboration et al. (2011) for the different parameters on which the spinning dust emissivity depends. This involves a modification of the code, as by default it uses $a_{0,1}=0.35$ nm and $a_{0,2}=3.0$ nm for the centroids of the two lognormal functions defining the dust grains size distributions, while the carbon abundance per hydrogen nucleus, $b_{\rm C}$, is selected from any of the values in table 1 of Weingartner & Draine (2001). Following Planck Collaboration et al. (2011), we use instead $a_{0,1}=0.58$ nm and $a_{0,1}=0.53$ nm, respectively for the molecular and atomic phases, and $b_{\rm C}=68$ ppm. Using these two models, we get a good fit to our full dataset, with $\chi^{2}_{\rm red}=0.99$, where the spinning dust component is clearly dominated by the molecular phase, as was found by Planck Collaboration et al. (2011).

To analyze the possibility of the existence of slightly different spinning dust models that could provide a better fit to the data, we take the AME residual fluxes from the previous fit, and produce a grid of models varying some of the parameters of the molecular phase component. In particular, we vary: i) the hydrogen number density $n_{\rm H}$ between 10 and 500 cm^-3 with a step of 5 cm^-3; ii) the kinetic gas temperature $T_{\rm g}$ between 5 and 200 K in steps of 5 K; iii) the intensity of the radiation field relative to the average interstellar radiation field, for which we consider only the values $G_{0}=1$ and 2; iv) and the hydrogen ionization fraction, for which we consider the values $x_{\rm H}=10$, $112$, $1000$ and $10000$ ppm. The best fit is obtained for $G_{0}=1$ and $x_{\rm H}=112$ ppm, the same values used in Planck Collaboration et al. (2011). As the fit is very degenerate, to define the most-likely values for $n_{\rm H}$ and $T_{\rm g}$ we set Gaussian priors on four different parameters. First, we put soft priors on $n_{\rm H}$ and $T_{\rm g}$ centred on the same values used in Planck Collaboration et al. (2011), $(n_{\rm H})_{0}=250$ cm^-3 and $(T_{\rm g})_{0}=40$ K, and with standard deviations $\sigma(n_{\rm H})=80$ cm^-3 and $\sigma(T_{\rm g})=60$ K. An additional prior on $N_{\rm H}^{\rm mol}$ can be derived from the canonical relation $2.13\times 10^{24}$ H cm${}^{-2}=1\tau_{100}$ (Finkbeiner et al., 2004). Using $\tau_{250}$ and $\beta_{\rm d}$ from our best-fit model of the thermal dust component, we extrapolate the optical depth to $100~\mu$m, and find $N_{\rm H}^{\rm mol}=2.909\times 10^{21}$ H cm^-2. We use this value to define the centre of the Gaussian prior, and $\sigma(N_{\rm H}^{\rm mol})=2\times 10^{21}$ H cm^-2. Finally, the ratio between the hydrogen column density and the hydrogen volume density, $z=N_{\rm H}^{\rm mol}/n_{\rm H}$, gives an estimate of the length along the line of sight of the spinning-dust-emitting region. We assume that this length might be of the order of the transverse angular size of the source. The source subtends an angle of around $2^{\circ}$, which at the distance of the Perseus complex, 260 pc (Cernicharo et al., 1985), corresponds to 9.08 pc. We therefore set a fourth prior on this quantity defined by $z_{0}=9.08$ pc and $\sigma(z)=4$ pc.

In Fig. 6 we show the marginalized likelihoods over $n_{\rm H}$ and $T_{\rm g}$. We define the best-values for these parameters from the 50% integrals of the probability distribution, and the confidence intervals from the region encompassing the 68% of the area around those values. We get $n_{\rm H}=223.2^{+69.5}_{-62.8}$ and $T_{\rm g}=60.2^{+45.3}_{-32.6}$ K. As it was said before, the values of the intensity of the radiation field and of the hydrogen ionization fraction that maximize the likelihood are $G_{0}=1$ and $x_{\rm H}=112$ ppm, respectively. We fix the other parameters of the molecular phase, and all the parameters corresponding to the atomic phase, to the same values that were used in Planck Collaboration et al. (2011). All these values are shown in Table 3. In this table, $x_{\rm C}$ represents the ionized carbon fractional abundance, $y$ the molecular hydrogen fractional abundance, and $\beta$ the average dipole moment per atom. The meaning of the other parameters have been explained before in the text. We then obtain the corresponding spinning dust spectra for the molecular and atomic phases using these parameters as inputs for spdust.2. Fixing these spectra, we perform a joint fit of the five aforementioned components, obtaining the best-fit values for the 7 parameters defining these models, which are also shown in Table 3. We get $\chi^{2}_{\rm red}=0.99$, the same value as before, so we do not manage to improve the quality of the global fit after improving the spinning dust models. This highlights the difficulty of constraining the parameters on which the spinning dust emission depends due to the the strong degeneracies between them.

The total hydrogen column density is $(4.40\pm 0.83)\times 10^{21}$ H cm^-2. This is a bit higher than the expected value of $2.89\times 10^{21}$ H cm^-2, which has been derived from the aforementioned $\tau_{100}-N_{\rm H}$ canonical relation, and extrapolating $\tau_{100}$ to $\tau_{250}$ using the fitted spectrum for the thermal dust emission. The inferred lengths of the two spinning dust phases along the line of sight, $z=N_{\rm H}/n_{\rm H}$, are respectively $z^{\rm mol}=5.63\pm 0.76$ pc and $z^{\rm at}<12.6$ pc (a 68.3% C.L. upper bound is used here, as the error bar is higher than the estimate). These values are of the order, or compatible, with the transverse size of the region, $9.08$ pc. Our fitted values for $\beta_{\rm d}$ and $T_{\rm d}$ are consistent within 1-sigma with those derived in Planck Collaboration et al. (2011). On the other hand, we get lower values for $EM$, $N_{\rm H}$ and $\tau_{250}$, but this is because these depend on the solid angle subtended by the region. This value is a factor $\approx 5$ higher in our case, because we are using aperture photometry instead of Gaussian fitting. When this factor is taken into account, then our values are brought into a better agreement with those of Planck Collaboration et al. (2011). Finally, we find a positive value for $\Delta T_{\rm cmb}$, whereas in Planck Collaboration et al. (2011) it was negative. However, we have checked that our value agrees with the average level of the CMB anisotropies within the aperture, which has been found to be $23.6~\mu$K in the Planck-DR1 CMB map resulting from the SMICA component separation method⁷⁷7Downloaded from the Planck Legacy Archive, http://www.sciops.esa.int/index.php?project=planck&page=
Planck_Legacy_Archive. (Planck Collaboration et al., 2014d).

In Fig. 7 we show the residual spinning dust spectrum, obtained after subtracting the best-fit free-free, CMB and thermal dust components.

As no clear emission is seen in the $Q$ and $U$ maps of Fig. 2, we derive upper limits on the polarization fraction of G159.6-18.5, following the procedure explained in section 2.3. In section 3.2 we defined the sizes of the aperture and of the background annulus so that we reproduced the fluxes in Planck Collaboration et al. (2011), which were obtained through Gaussian fitting. In order to minimize the error associated with the background subtraction, here we extend the size of the background annulus, and use a circular aperture with radius $1.5^{\circ}$, and a background ring between $1.5^{\circ}$ and $2.5^{\circ}$, with an extension towards the north-west as is shown in the maps of Fig. 2. At each frequency we calculate the RMS levels in this background annulus, to define the quantities $\sigma(T_{j})$ that are introduced in equation 5 to calculate the errors on the Stokes parameters $Q$ and $U$. These errors, together with the fluxes resulting from the aperture photometry integration, are quoted in Table 4. Note that here the intensity fluxes in the QUIJOTE frequencies are slightly different to those presented in Table 2 owing to the different apertures. In order to get the AME residual fluxes shown in Table 4, using the new fluxes, we repeat the same fit that was performed in section 3.2, considering the same spinning dust spectra for the molecular and atomic phases.

It is important to point out that, as it became clear when we discussed the results of the jackknife tests in section 3.1, the uncertainties on the $I$ fluxes are here biased high because some extended emission of the source leaks into the background annulus. This does not have significant implications in our analysis as the uncertainty in the polarization fraction is driven by the errors in $Q$ and $U$.

The $Q$ and $U$ fluxes shown in Table 4 are consistent with zero, and therefore we obtain de-biased (in section 2.3 we explained how this de-biasing is applied) upper limits at the 95% confidence level on the polarized intensity, $P_{\rm db}$. We also show in Table 4 upper limits at the 95% confidence level on the polarization fraction, taking as reference both the total intensity ($\Pi_{\rm db}$) and the residual AME intensity ($\Pi_{\rm AME,db}$). The values below the horizontal line in this table correspond to constraints obtained in maps that have been built by combining the two frequency bands of each horn. The most stringent upper limit we get on the polarization fraction is $\Pi_{\rm db}<2.85\%$, and is obtained after combining the maps at 16.7 and 18.7 GHz. Due to the decrease of the intensity flux at lower frequencies, the constraints at 11.2 and 12.9 GHz are less stringent.

Note that, under the reliable assumption that the free-free emission is unpolarized (Rybicki & Lightman, 1979), any possible detection of polarization below $\nu\sim 30$ GHz, where the thermal dust is clearly sub-dominant, should in principle be ascribed to AME. One caveat to this hypothesis is the possible presence of a Faraday Screen (FS) hosting a strong and regular magnetic field, which could induce a rotation of the background polarized emission. This idea was proposed by Reich & Reich (2009) to explain the high degree of polarization they detected towards G159.6-18.5 in observations from the Effelsberg telescope at 2.7 GHz. They suggest that this same mechanism could indeed be the responsible for the tentative polarized emission seen by Battistelli et al. (2006) at 11 GHz. Using as face value the rotation measure obtained by Reich & Reich (2009), RM$=190$ rad m^-2, López-Caraballo et al. (2011) estimated a polarization fraction of $\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{$\sim$}}\hss}\hbox{$<$}}}0.2$%, well below the upper limit at this frequency. Using this RM we estimate polarization fractions from the FS of $\approx 3.5\%$ and $\approx 0.4\%$, respectively at 11 and 19 GHz. These values are well below our upper limits, while the value at 11 GHz is compatible with the measurement of Battistelli et al. (2006).

In section 3.1 we mentioned that one important consistency test for our data processing is the verification that it shows no polarization at the position of un-polarized sources like the California HII region. In Table 5 we show the corresponding $I$, $Q$ and $U$ fluxes, and derived upper limits on $P$ and $\Pi$, which have been obtained using the apertures shown in Fig. 4. The limits on the polarization fraction stand at the level of $2$ to $4\%$ (95% C.L.), depending on the frequency band. This constrains the possible existence of a leakage from intensity to polarization to values below this level. In fact, as it was mentioned in section 2.1, we expect the leakage to be below 0.3% ($\leq-25$ dB). This number has been verified using observations of Cas A, where we recover a leakage pattern with a quadrupolar shape in the $Q$ and $U$ maps, with a level of around $0.2\%$ (further details will be presented in separate technical papers).

We plot our constraints, and previous results in the literature at different frequencies, in Fig. 8. It can be seen here that our observations fill the gap between previous results at frequencies below 11 GHz and above 20 GHz. The point at 9.65 GHz is an upper limit coming from GBT observations on LDN1622 (Mason et al., 2009). The value at 11 GHz represents the aforementioned tentative detection towards G159.6-18.5 at the $1.8\sigma$ level ($\Pi=3.4_{-1.9}^{+1.5}\%$) of AME polarization from the COSMOSOMAS experiment (Battistelli et al., 2006), which is still consistent with our upper bound at $12.0$ GHz. The points at $\nu>20$ GHz come from WMAP 7-year data on the Perseus (López-Caraballo et al., 2011) and $\rho$ Ophiuchi (Dickinson et al. 2011, they also present results on Perseus) molecular clouds and on the HII region [LPH96] 201.663+1.643 (Rubiño-Martín et al. 2012b, they also show less-stringent constraints from the Pleiades reflection nebula and from the dark nebula LDN1622). Currently, the most stringent constraint is that obtained by López-Caraballo et al. (2011): $\Pi<0.98\%$ at 23 GHz. This limit benefits from the fact that it is obtained at a frequency close to the peak ($\approx 28$ GHz) of the AME in G159.6-18.5. The constraints on the polarization fraction from QUIJOTE come from a spectral region where the intensity flux density drops. However, according to the model of Hoang et al. (2013), the polarization fraction of the spinning dust emission peaks at $\approx 5$ GHz, and therefore observations at lower frequencies are more appropriate to constrain this model. All measurements, including ours, are consistent with this model, except maybe the previous limit from López-Caraballo et al. (2011) and that from Mason et al. (2009), which are slightly below.

Our upper limits on $\Pi$ are not only affected by the decrease of AME flux at our frequencies, but also by the sparse sky coverage of our observations. We decided to survey an area large enough to cover the California HII region and also to get the source simultaneously in the four horns. Currently we are undertaking observations on a smaller area, centred in two individual horns, with the goal to increase by a factor $8$ the integration time per unit area. This would improve our map sensitivity by a factor $2.8$, which would help push our current upper limits down to $\approx 1.7\%$ at $18.7$ GHz, and $\approx 1.0\%$ after combining the $16.7$ and $18.7$ GHz frequency bands, which is well below the polarization degree predicted by the Hoang et al. (2013) model at these frequencies. Furthermore, future observations at 30 GHz with the upcoming TGI experiment, which is expected to be 13 times more sensitive than the current MFI, will push these upper limits down by a significant amount.

However, it is important to note that the curve of the polarization degree obtained by Hoang et al. (2013) represents an upper limit. They inferred the alignment efficiency of interstellar dust grains from observations of the UV polarization excess towards two stars to predict the polarization degree of the spinning dust emission. In so doing they assumed that the UV polarization bump is produced exclusively by polycyclic aromatic hydrocarbon (PAHs) molecules, those which are thought to be responsible for spinning dust emission. However, if the graphite grains were also aligned, contributing to the UV polarization bump, then the alignment efficiency of PAHs could actually be lower and so would be the inferred degree of spinning dust polarization.

The previous model from Lazarian & Draine (2000) for spinning dust polarization is also plotted in Fig. 8, for the case of a cold neutral medium (CNM), together with different models from Draine & Lazarian (1999) predicting the degree of polarization of the magnetic dipole emission for different grain shapes and compositions. In particular, we show the cases of grains made of metallic Fe, and of the hypothetical material X4 defined in Draine & Lazarian (1999). Concerning the grain geometries, different axial ratios $a1$:$a2$:$a3$ are shown, as indicated in the figure. All these models predict polarization fractions typically above $10\%$, much higher than that of the spinning dust emission, and are ruled out by our and previous observations. However, it must be noted that these models consider the dust grains to have a perfect ordering of the atomic magnetic moments in a single domain. A more realistic case considering randomly-oriented single-domain magnetic inclusions has been recently studied by Draine & Hensley (2013), resulting in much lower polarization degrees, of the order of $5\%$ in the range $10-20$ GHz, or below.