Star Formation and ISM Properties in Host Galaxies of Three Far-infrared Luminous Quasars at $z\sim6$

TL;DR

This study uses VLA, Herschel, and ALMA observations to analyze the molecular gas, dust, and star formation in three luminous quasars at z~6, revealing ISM properties akin to local starburst galaxies.

Contribution

First detailed multi-wavelength analysis of ISM properties in z~6 quasars, combining CO line, dust, and star formation measurements.

Findings

Molecular gas masses of (0.8-4.3) x 10^10 M_sun in the quasars

Star formation rates of 2400-2700 M_sun/yr

ISM properties similar to local starburst galaxies

Abstract

We present Karl G. Jansky Very Large Array (VLA) observations of the CO (2$-$1) line emission towards three far-infrared luminous quasars at $z\sim6$: SDSS J231038.88$+$185519.7 and SDSS J012958.51$-$003539.7 with $\sim0\farcs6$ resolution and SDSS J205406.42$-$000514.8 with $\sim2\farcs1$ resolution. All three sources are detected in the CO (2$-$1) line emission -- one source is marginally resolved, and the other two appear as point sources. Measurements of the CO (2$-$1) line emission allow us to calculate the molecular gas mass even without a CO excitation model. The inferred molecular gas masses are (0.8$-$4.3) $\times$ 10$^{10}$ $M_{\odot}$. The widths and redshifts derived from the CO (2$-$1) line are consistent with previous CO (6$-$5) and [\ion{C}{2}] measurements. We also report continuum measurements using the Herschel for SDSS J231038.88$+$185519.7 and SDSS J012958.51$-$003539.7, and for SDSS J231038.88+185519.7, data obtained at $\sim140$ and $\sim300$ GHz using the Atacama Large Millimeter/submillimeter Array (ALMA). In the case of SDSS J231038.88+185519.7, we present a detailed analysis of the spectral energy distribution and derive the dust temperature ($\sim40$ K), the dust mass ($\sim10^{9}$ $M_{\odot}$), the far-infrared luminosity (8$-$1000 $\mu$m; $\sim10^{13}$ $ L_{\odot}$) and the star formation rate (2400$-$2700 $M_{\odot}$ yr$^{-1}$). Finally, an analysis of the photo-dissociation regions associated with the three high redshift quasars indicates that the interstellar medium in these sources has similar properties to local starburst galaxies.

Tables

Table 1. Table 1: Sample and observations

Source	Herschel						VLA			ALMA
	Herschel PACS			Herschel SPIRE
	OD	OBSIDs		OD	OBSIDs		$t_{\mathrm{obs}}$	Configuration		$t_{\mathrm{obs}}$
(1)	(2)	(3)		(4)	(5)		(6)	(7)		(8)
J2310$+$1855	1122	13422468221/1342246822		1314	1342257359		32 h	C		5.5 h
J0129$-$0035	-	-		1330	1342258372		15 h	C		-
J2054$-$0005	-	-		-	-		20 h	D		-

Table 2. Table 2: CO (2 − - 1) results compared with previous work

Species	Redshift	FWHM	Flux	Luminosity	$L_{\mathrm{CO}\left(1-0\right)}^{{}^{\prime }}$	$M_{\mathrm{gas}}$
		(km s-1)	(Jy km s-1)	(106 $L_{\odot }$)	(1010$K\mathrm{km}{s}^{-1}{\mathrm{pc}}^2$)	(1010$M_{\odot }$)
(1)	(2)	(3)	(4)	(5)	(6)	(7)
J2310+1855
CO (6$-$5)	6.0025 $\pm$ 0.0007	456 $\pm$ 64	1.52 $\pm$ 0.13	543.0 $\pm$ 46.4	6.5 $\pm$ 0.6	5.2
[C ii]	6.0031 $\pm$ 0.0002	393 $\pm$ 21	8.83 $\pm$ 0.44	8700 $\pm$ 1400	-	-
CO (2$-$1)	6.0029 $\pm$ 0.0005	484 $\pm$ 48	0.18 $\pm$ 0.02	21.3 $\pm$ 0.2	5.4 $\pm$ 0.5	4.3
J0129$-$0035
CO (6$-$5)	5.7794 $\pm$ 0.0008	283 $\pm$ 87	0.37 $\pm$ 0.07	125.0 $\pm$ 23.6	1.5 $\pm$ 0.3	1.2
[C ii]	5.7787 $\pm$ 0.0001	194 $\pm$ 12	1.99 $\pm$ 0.12	1800 $\pm$ 300	-	-
CO (2$-$1) peak	5.7783 $\pm$ 0.0004	195 $\pm$ 41	0.036 $\pm$ 0.005	4.0 $\pm$ 0.6	1.0 $\pm$ 0.1	0.8
J2054$-$0005
CO (6$-$5)	6.0379 $\pm$ 0.0022	360 $\pm$ 110	0.34 $\pm$ 0.07	122.5 $\pm$ 25.2	1.5 $\pm$ 0.3	1.2
[C ii]	6.0391 $\pm$ 0.0001	243 $\pm$ 10	3.37 $\pm$ 0.12	3300 $\pm$ 500	-	-
CO (2$-$1) peak	6.0394 $\pm$ 0.0004	270 $\pm$ 47	0.06 $\pm$ 0.01	6.6 $\pm$ 0.9	1.7 $\pm$ 0.2	1.3

Table 3. Table 3: Continuum photometry information

	Telescope	J2310$+$1855	J0129$-$0035	J2054$-$0005
$m_{1450\AA }$ (mag)	-	19.30g	22.28e	20.6c
$z$ (mag)	SDSS	19.31 $\pm$ 0.11a	22.16 $\pm$ 0.11e	20.72 $\pm$ 0.09c
$w_1$ (mag)	WISE	18.48 $\pm$ 0.05b	-	-
$w_2$ (mag)	WISE	18.73 $\pm$ 0.12b	-	-
$w_3$ (mag)	WISE	17.51 $\pm$ 0.44b,1	-	-
$F_{100\mu m}$ (mJy)	Herschel/PACS	6.5 $\pm$ 1.2	-	2.7h
$F_{160\mu m}$ (mJy)	Herschel/PACS	13.2 $\pm$ 2.8	-	9.8 $\pm$ 1.3h
$F_{250\mu m}$ (mJy)	Herschel/SPIRE	19.9 $\pm$ 6.0	12.2	15.2 $\pm$ 5.4h
$F_{350\mu m}$ (mJy)	Herschel/SPIRE	22.0 $\pm$ 6.9	11.4	12.0 $\pm$ 4.9h
$F_{500\mu m}$ (mJy)	Herschel/SPIRE	29.4	15.2	19.5h
$F_{496GHz}$ (mJy)	ALMA	24.89 $\pm$ 0.72j	-	-
$F_{350GHz}$ (mJy)	ALMA	14.54 $\pm$ 0.21	-	-
$F_{338GHz}$ (mJy)	ALMA	14.49 $\pm$ 0.21	-	-
$F_{295GHz}$ (mJy)	ALMA	11.56 $\pm$ 0.10	-	-
$F_{283GHz}$ (mJy)	ALMA	11.60 $\pm$ 0.10	-	-
$F_{263GHz}$ (mJy)	ALMA	8.91 $\pm$ 0.08g	2.57 $\pm$ 0.06g	2.98 $\pm$ 0.05g
$F_{250GHz}$ (mJy)	IRAM	8.29 $\pm$ 0.63g	2.37 $\pm$ 0.49f	2.38 $\pm$ 0.53d
$F_{147GHz}$ (mJy)	ALMA	1.59 $\pm$ 0.04	-	-
$F_{143GHz}$ (mJy)	ALMA	1.41 $\pm$ 0.04	-	-
$F_{135GHz}$ (mJy)	ALMA	1.24 $\pm$ 0.04	-	-
$F_{131GHz}$ (mJy)	ALMA	1.10 $\pm$ 0.05	-	-
$F_{91.5\mathrm{GHz}}$ (mJy)	ALMA	0.41 $\pm$ 0.03i	-	-

Table 4. Table 4: Physical properties derived from the SED fit of J2310 + + 1855

		P+CAT 3D$\_$n+MBB	P+CAT 3D$\_$w+MBB	P+CAT 3D$\_$n+DL07	P+CAT 3D$\_$w+DL07
${\alpha }_{\mathrm{UV}/\mathrm{opt}}$	(1)	$-$0.47${}_{}{}^{+0.08}{}_{-0.08}{}^{}$	$-$0.46${}_{}{}^{+0.09}{}_{-0.08}{}^{}$	$-$0.47${}_{}{}^{+0.08}{}_{-0.08}{}^{}$	$-$0.47${}_{}{}^{+0.08}{}_{-0.09}{}^{}$
$L_{\mathrm{UV}/\mathrm{opt}}$ (1046 erg s-1)	(2)	14.6${}_{}{}^{+1.0}{}_{-1.0}{}^{}$	14.4${}_{}{}^{+1.0}{}_{-1.0}{}^{}$	14.0${}_{}{}^{+1.0}{}_{-1.0}{}^{}$	13.8${}_{}{}^{+1.0}{}_{-1.0}{}^{}$
$a$	(3)	$-$1.25	$-$1.50	$-$1.75	$-$3.00
$h$	(4)	1.25	0.50	0.25	0.10
$N_0$	(5)	10.0	10.0	5.0	10.0
$a_w$	(6)	-	$-0.50$	-	$-$0.50
${\theta }_w$	(7)	-	30	-	30
${\sigma }_{\theta }$	(8)	-	7.50	-	10.00
$f_{\mathrm{wd}}$	(9)	-	2.00	-	0.15
$R_{\mathrm{out}}$	(10)	500.0	500.0	500.0	450.0
${\tau }_{\mathrm{cl}}$	(11)	50.0	50.0	50.0	50.0
Inclination (degree)	(12)	15.0	0.0	0.0	45.0
$\gamma$	(13)	-	-	0.17${}_{}{}^{+0.07}{}_{-0.07}{}^{}$	0.39${}_{}{}^{+0.13}{}_{-0.11}{}^{}$
$U_{\min }$	(14)	-	-	25.0	25.0
$U_{\max }$	(15)	-	-	106	106
$q_{\mathrm{PAH}}$ (%)	(16)	-	-	0.47	0.01
$T_{\mathrm{dust}}$ (K)	(17)	39${}_{}{}^{+3}{}_{-3}{}^{}$	40${}_{}{}^{+3}{}_{-2}{}^{}$	${\mathrm{𝟑𝟔}}_{-\mathrm{𝟐}}^{+\mathrm{𝟐}}$	${\mathrm{𝟒𝟎}}_{-\mathrm{𝟐}}^{+\mathrm{𝟐}}$
$L_{\mathrm{FIR}}$ (${10}^{13}{L}_{\odot }$)	(18)	${1.4}_{-0.3}^{+0.3}$	${1.6}_{-0.3}^{+0.3}$	${3.2}_{-0.7}^{+0.7}$	${3.6}_{-0.6}^{+0.6}$
SFR (${10}^3{M}_{\odot }{\mathrm{yr}}^{-1}$)	(19)	${2.4}_{-0.4}^{+0.5}$	${2.7}_{-0.5}^{+0.6}$	${5.5}_{-1.2}^{+1.1}$	${6.1}_{-1.0}^{+1.0}$
$M_{\mathrm{dust}}$ (${10}^9{M}_{\odot }$)	(20)	${1.7}_{-0.3}^{+0.4}$	${1.6}_{-0.3}^{+0.3}$	3.8${}_{}{}^{+0.2}{}_{-0.3}{}^{}$	3.9${}_{}{}^{+0.3}{}_{-0.3}{}^{}$

Equations

$$\nu F_{\nu ,\mathrm{model}}=\frac{M_{\mathrm{dust}}}{4\pi D_L^2}[(1-\gamma )\nu p_{\nu }^{(0)}(q_{\mathrm{PAH}},U_{\min })+\gamma \nu p_{\nu }(q_{\mathrm{PAH}},U_{\min },U_{\max },\alpha )]$$

$$\mathrm{Fraction}=\frac{F_{\nu ,\mathrm{obs}}}{F_{\nu ,\mathrm{intrinsic}}}=1-\frac{B_{\nu }[\lambda ,T_{\mathrm{CMB}}(z)]}{B_{\nu }[\lambda ,T_{\mathrm{dust}}(z)]}$$

$$U_{\mathrm{ave}}=(1-\gamma )\times U_{\min }+\frac{\gamma \times \ln (\frac{U_{\max }}{U_{\min }})}{\frac{1}{U_{\min }}-\frac{1}{U_{\max }}}$$

$$T_{DL07}\approx 17.0\times U_{\mathrm{ave}}^{\frac{1}{4+\beta }}\mathrm{K}(\mathrm{grain}\mathrm{size}>0.03\mu \mathrm{m})$$

$$M_{\mathrm{dust}}=\frac{F_{125\mu m}{D}_L^2}{{\kappa }_{125\mu \mathrm{m}}{B}_{\nu }(125\mu m,T_{\mathrm{FIR}})}$$

Full text

galaxies: evolution — galaxies: active — galaxies: high redshift — submillimeter: galaxies — quasars: general — radio line: galaxies

66affiliationtext: Max-Planck-Institut for Astronomie, Königstuhl 17, D-69117 Heidelberg, Germany1313affiliationtext: Astronomy Department, California Institute of Technology, MC105-24, Pasadena, CA 91125, USA\dagger$$\daggeraffiliationtext: E-mail: shaoyali0922@pku.edu.cn\bigstar$$\bigstaraffiliationtext: E-mail: rwangkiaa@pku.edu.cn

1 Introduction

More than 200 quasars have been discovered above redshift 5.7 (e.g., Fan et al. 2006; Mortlock et al. 2009, 2011; Willott et al. 2009, 2010a, 2010b; Venemans et al. 2013, 2015; Bañados et al. 2014, 2016, 2018; Jiang et al. 2015, 2016; Matsuoka et al. 2016, 2018; Mazzucchelli et al. 2017). These high redshift quasars are key to understand the co-evolution between supermassive black holes (SMBHs) and their host galaxies at the end of the reionization epoch. Observations of the dust, molecular and atomic gas content of these objects allow us to probe their star formation activity and derive their interstellar medium properties.

The rest frame far-infrared (FIR) continuum emission in these sources originates mainly from dust heated by ultraviolet (UV) radiation from young and massive stars in the host galaxies. Submillimeter and millimeter observations using, e.g., the instruments on Herschel, the Max Plank Millimeter Bolometer Array (MAMBO) on the IRAM-30 m telescope or the Submillimeter Common User Bolometer Array (SCUBA) on the James Clerk Maxwell Telescope (JCMT) (Bertoldi et al. 2003a; Petric et al. 2003; Robson et al. 2004; Beelen et al. 2006; Wang et al. 2007, 2008a, 2008b, 2011b; Leipski et al. 2013, 2014), and the ALMA (Decarli et al. 2018; Venemans et al. 2018), have detected dust continuum in the host galaxies of many $z\sim 6$ quasars, with FIR luminosities of $\sim 10^{11-13}$ $L_{\odot}$, and dust masses on the order of 107-9 $M_{\odot}$. The most luminous objects have FIR luminosities similar to those of the ultraluminous infrared galaxies (ULIRGs; $L_{\rm FIR}>10^{12}L_{\odot}$) and hyper-luminous infrared galaxies (HLIRGs; $L_{\rm FIR}>10^{13}L_{\odot}$) in the local universe (e.g., Rosenberg et al. 2015), indicating that they are forming stars with star formation rates (SFRs) of a few tens to thousands $M_{\odot}$ yr*-1*, co-eval with rapid SMBH accretion.

Most of these high-$z$ FIR-luminous quasars have been detected in carbon monoxide (CO). As principal molecular tracer, CO provides a tool to probe the physical conditions of star-forming gaseous reservoirs, through multiple transitions redshifted into the submillimeter/millimeter range.

Intermediate to high-$J$ transitions of CO from 3$-$2 to 9$-$8 have been detected with the VLA, the IRAM Plateau de Bure interferometer and its successor - the Northern Extended Millimeter array (NOEMA) and the ALMA in many high redshift quasars (Bertoldi et al. 2003a; Walter et al. 2003, 2004; Carilli et al. 2007; Riechers et al. 2009; Wang et al. 2010, 2011b, 2013, 2016; Venemans et al. 2017a, b; Feruglio et al. 2018; Li et al. in preparation). Walter et al. (2004) observed the CO (3$-$2) line emission from SDSS J114816.64+525150.3 (hereafter 1148+5251) at $z=6.42$ using the VLA, and measured a CO source size of 3.6 kpc $\times$ 1.4 kpc (full width at half maximum; FWHM). They also derived a dynamical mass of 4.5 $\times$ 1010 $M_{\odot}$, which is less than the stellar bulge mass (of order 1012 $M_{\odot}$) predicted by the present-day $M_{\rm BH}-M_{\rm bulge}$ relation (Kormendy & Ho, 2013), which may indicate a faster SMBHs evolution than their hosts.

Observations of the low-$J$ CO transitions ($J_{\rm upper}$ $\leqslant$ 2) in $z\sim 6$ quasars are very difficult due to the low flux density and limited telescope sensitivity. Only eight $z\sim 6$ quasars have been observed in CO (2$-$1) line emission ($\nu_{\rm rest}$ = 230.538 GHz) all using the VLA Ka band, and five of them have been detected (Wang et al. 2011a, 2016; Stefan et al. 2015; Venemans et al. 2017b). The VLA is the only instrument that can observe the CO (2$-$1) line with proper frequency coverage and high sensitivity for $z\sim 6$ objects. The CO (2$-$1) line emission allows us to measure the molecular gas mass directly with a molecular gas mass conversion factor. The CO (2$-$1) line is also a crucial tracer to probe the low-$J$ part of CO spectral line energy distributions (SLEDs), and the excitation, the spatial distribution and the surface density of the extended cold gas in star-forming quasar host galaxies. For example, by observing the $z\sim 5.7$ quasar SDSS J092721.82+200123.7, Wang et al. (2011a) constrained the CO (2$-$1) line excitation in the central $\sim 10$ kpc of the source, and estimated a molecular gas mass of order of 1010 $M_{\odot}$. The redshift, the line width and the gas mass derived from CO (2$-$1) in this quasar are consistent with those from CO (6$-$5) and CO (5$-$4) observations in Carilli et al. (2007). However, the CO (2$-$1) line has been detected only toward a few $\sim 6$ quasars. In this paper we increase the sample with CO (2$-$1) detections by a factor of 60$\%$. These data will be critical in investigations of the CO SLEDs with multiple CO transition observations in the future.

In this paper, we report VLA observations of CO (2$-$1) line emission in three FIR-luminous high redshift quasars: SDSS J231038.88$+$185519.7 (hereafter J2310$+$1855 at $z=$ 6.0), SDSS J012958.51$-$003539.7 (hereafter J0129$-$0035 at $z=$ 5.7), and SDSS J205406.42$-$000514.8 (hereafter J2054$-$0005 at $z=$ 6.0). We complement these data with continuum measurements from Herschel Photodetector Array Camera and Spectrometer instrument (PACS; Poglitsch et al. 2010) and Spectral and Photometric Imaging Receiver (SPIRE; Griffin et al. 2010) for J2310$+$1855 and J0129$-$0035, and in the case of J2310+1855, with ALMA. This paper is organized as follows. In Section 2 we describe the sample, the observations and the data reduction. In Section 3 we present our results. In Section 4, we first discuss the gas distribution and the gas mass of the three targets; then we analyze the dust temperature, the dust mass and the star formation rate of J2310$+$1855 through a fit to the continuum spectral energy distribution (SED); finally we discuss models for PDRs in the three targets. In Section 5, we present a short summary. Throughout this work we assume a $\Lambda$CDM cosmology with $H_{0}$ = 71 km s*-1* Mpc*-1*, $\Omega_{\rm M}$ = 0.27 and $\Omega_{\Lambda}$ = 0.73 (Spergel et al., 2007).

2 Observations and data reduction

We selected the three high redshift quasars - J2310$+$1855, J0129$-$0035 and J2054$-$0005 from our previous CO (6$-$5) surveys (Wang et al. 2010; Wang et al. 2011b). They are all FIR-luminous quasars and contain a lot of gas. We carried out VLA observations of the CO (2$-$1) line emission in the three objects. FIR continuum data for J2310$+$1855 and J0129$-$0035 were also obtained from Herschel PACS and SPIRE, and measurements of the continuum at 140 and 300 GHz were made using the ALMA in the case of J2310+1855. We list the details and the observations of the three sources in Table 1.

2.1 VLA

The CO (2$-$1) observations were performed using the Ka band receivers (centered at 32 GHz) of the VLA. The observations of J2310+1855 were carried out on 2012 January 28, and 2014 October 14 and 15 using the C-configuration. The total observing time was 15 hours and 8 hours on-source. For J0129$-$0035, the data were taken between 2013 July 8 and August 6 also in the C-configuration, with a total observing time of 32 hours and 17 hours on-source. Finally, J2054$-$005 was observed from 2014 June 29 to July 21 using the D-configuration, with a total observing time of 20 hours and 10.5 hours on-source. We used the WIDAR 8-bit samplers to maximize the line sensitivity. The total bandwidth is 2 GHz with 16 128 MHz-wide sub-bands in full polarization mode, and we centred the CO (2$-$1) line in one of the 128 MHz sub-bands. The reference redshifts of these sources were based on the previous CO (6$-$5) observations ($\nu_{\rm obs}$ = 32.922 GHz for J2310+1855; $\nu_{\rm obs}$ = 34.006 GHz for J0129$-$0035; $\nu_{\rm obs}$ = 32.757 GHz for J2054$-$005; Wang et al. 2010, 2011a). Flux calibrations were performed with the standard VLA calibrators: 3C286 and 3C48. The spatial resolutions achieved during these observations are $0\farcs 6$ for C-configuration and $2\farcs 1$ for D-configuration.

The data were reduced using the Common Astronomy Software Application (CASA111https://casa.nrao.edu/) standard pipeline. The final data cubes were reduced by the CLEAN task using a robust weighting factor of 0.5 to optimize the noise per frequency bin and the resolution of the final map.

2.2 Herschel

2.2.1 PACS

J2310$+$1855 was observed by Herschel PACS at 100 and 160 $\mu$m using the mini-scan map observing template. We executed the observations over two different scan angles with observing parameters as recommended in the mini-scan map Astronomical Observation Template (AOT) release note. Fourteen repetitions were employed for each scan direction. The final on-source integration time was 1792 seconds.

Data reduction was performed within the Herschel Interactive Processing Environment (HIPE; version 15.0.1; Ott 2010). We followed the standard pipeline for PACS mini-scan observations. Each scan direction was processed individually and mosaicked at the end of the procedure. Aperture photometry of the final mosaics was also performed using HIPE. We used aperture sizes of 6$\arcsec$ and 9$\arcsec$ radii in 100 and 160 $\mu$m bands, respectively. The residual sky emission was derived in a sky annulus between 20$\arcsec$ and 25$\arcsec$ (100 $\mu$m map) or 24$\arcsec$ and 28$\arcsec$ (160 $\mu$m map). Aperture corrections were determined from the encircled energy fraction of unresolved sources provided by calibrator observations. The photometric uncertainties cannot be measured directly from the pixel-to-pixel variations because the final PACS maps are heavily influenced by correlated noise. We followed Leipski et al. (2013, 2014) to determine the photometric uncertainties: first we masked the centre source with the same aperture size used in quasar photometry, then randomly placed $\sim 1000$ apertures on the images with the same diameter as used for photometry. In order to exclude the noisy edge of the maps, we restricted ourselves to the region of sky with 75$\%$ or more integration time compared with the position of the quasar. Then we fitted a Gaussian to the measured fluxes in 1000 apertures. We set the 1$\sigma$ photometric uncertainty to the sigma value of the Gaussian profile of the final map. The Herschel PACS photometric results for J2310+1855 are listed in Table 3.3.

2.2.2 SPIRE

Herschel SPIRE observations towards J2310$+$1855 and J0129$-$0035 were carried out at 250, 350, and 500 $\mu$m using small scan map mode with 9 repetitions for each object. The total on-source integration time per source was 370 seconds. SPIRE observations are dominated by the confusion noise (Nguyen et al. 2010; PSW: 5.8 mJy beam*-1*; PMW: 6.3 mJy beam*-1*; PLW: 6.8 mJy beam*-1*).

Data reduction was executed in HIPE (version 15.0.1) following the standard SPIRE small scan observation pipeline. We used the HIPE built-in source extractor called ‘sourceExtractorSussextractor’ (Savage & Oliver 2007) to determine the source location and the flux density. We estimated photometric uncertainties following the method in Leipski et al. (2013, 2014): first we created an artificial source image with all sources found by ‘sourceExtractorSussextractor’ in the calibrated map, then subtracted it from the observed map to get the residual map. On this residual map we determined the pixel-to-pixel $rms$ in a reasonable box with size of $\sim 8$ times the FWHM ($18\farcs 2$, $24\farcs 9$, and $36\farcs 3$ for default map pixel sizes of 6$\arcsec$, 10$\arcsec$, and 14$\arcsec$ at 250, 350, and 500 $\mu$m, individually) and centered at the position of the quasar in order to have enough sampling of quasar surrounding sky and avoid lower coverage regions of noisy edged sky. For the J2310$+$1855 SPIRE observations in particular, the sky is not uniform due to bright foreground sources, so we calculated the flux noise as follows: first we added 100 fake sources with flux densities equal to that of the target with a Gaussian scatter of 3 mJy in the raw data, and used standard procedures to reduce the data to generate the scientific image. The fake sources are in source-free areas with good coverage (e.g. $>$60$\%$) and are near targets (e.g. within $\sim 8$ times the FWHM with targets at the center). Then we used ‘sourceExtractorSussextractor’ to estimate source flux density. If the source was not detected, we measured the flux density at the input position. The quoted photometric uncertainty was the $rms$ difference between the input and output fluxes. The final Herschel SPIRE photometric results for J2310+1855 and J0129$-$0035 are presented in Table 3.3.

2.3 ALMA

We observed the dust continuum emission around 140 and 300 GHz towards J2310+1855 during ALMA Cycle 3 (program ID: 2015.1.01265.S; PI: Ran Wang). This program aims to observe CO (9$-$8), CO (8$-$7), [N II], and [O I] and their underlying continuum emission. We here only report on the continuum data and will report on the molecular and atomic line emission in a subsequent paper (Li et al. in preparation). The observations were carried out using Bands 4 and 6 with spatial resolutions of 0$\farcs$42 to 0$\farcs$65 with 36$-$44 12-m antennas. The on-source observation time was 5.5 hours and 3.4 hours, respectively. One of the 2 GHz spectral windows was tuned on the line and the other three on the continuum. The reference redshift of $z=$ 6.0031 is from previous ALMA [C ii] observations presented in Wang et al. (2013). The phase and flux calibrators were J2253+1608 and Pallas. The flux density scale calibration accuracy is better than 10$\%$. We reduced the data with the CASA pipeline. The noise levels of the dust maps are 0.02$-$0.04 mJy beam*-1*.

3 Results

3.1 J2310$+$1855

This source was discovered in the SDSS (Wang et al. 2013; Jiang et al. 2016). With $m_{\rm 1450\AA}$ = 19.30 mag, it is one of the brightest optical quasars among the known $z\sim 6$ quasars. The 250 GHz dust continuum flux density is 8.29 $\pm$ 0.63 mJy (Wang et al. 2013), which makes it one of the most FIR-luminous quasars known at $z\sim 6$. The top left panel of Figure 1 shows the CO (2$-$1) intensity map integrated from $-$333 to 323 km s*-1*. The white cross is the quasar position calibrated with nearby bright stars by Gaia astrometry. The uncertainties of the Gaia calibrated position are 96 mas in RA and 106 mas in Dec. Using a 2D Gaussian to fit the CO (2$-$1) velocity integrated map, the source is found to be marginally resolved with a source size of ($0\farcs 602$ $\pm$ $0\farcs 184$) $\times$ ($0\farcs 400$ $\pm$ $0\farcs 208$). This size is consistent within the errors with that measured by [C II] line emission (($0\farcs 56$ $\pm$ $0\farcs 03$) $\times$ ($0\farcs 39$ $\pm$ $0\farcs 04$); Wang et al. 2013). The corresponding physical source size is about (3.51 $\pm$ 1.07) kpc $\times$ (2.33 $\pm$ 1.21) kpc, which is comparable to other $z\sim 6$ quasars (Wang et al., 2013). The source position derived from the CO (2$-$1) line is consistent with these from [C II] and its underlying dust continuum emission in Wang et al. (2013) and the Gaia calibrated position. We fitted a Gaussian profile to the CO (2$-$1) spectrum (left panel of Figure 2), and measured a redshift of 6.0029 $\pm$ 0.0005, which is in agreement with the redshifts from other ISM tracers (6.0025 $\pm$ 0.0007 from CO (6$-$5); 6.0031 $\pm$ 0.0002 from [C II]; Wang et al. 2013). We fitted a line width (FWHM) of 484 $\pm$ 48 km s*-1*, which is consistent with the CO (6$-$5) line width (456 $\pm$ 64 km s*-1*) but somewhat larger (2$\sigma$) than the [C II] line width (393 $\pm$ 21 km s*-1*). These measurements are given in Table 3.3. The velocity map of CO (2$-$1) line emission (top right panel of Figure 1) shows a velocity gradient, which is consistent with that from [C II] (Wang et al., 2013).

J2310$+$1855 was detected in two Herschel PACS bands at $\sim 5\sigma$ and in two SPIRE bands at $\sim 3\sigma$. In the case of the SPIRE 500 $\mu$m band observations, we give 3$\sigma$ value as an upper limit. The photometric results are listed in Table 3.3.

The continuum flux densities using 2D Gaussian fits to the ALMA continuum maps at 140 GHz and 300 GHz are listed in Table 3.3. The continua are all spatially resolved at the observing frequencies, with consistent deconvolved source sizes of $\sim 0\farcs 2\times 0\farcs 2$. The derived source positions from the dust continua are consistent with that from CO (2$-$1) line. We present a dust continuum map at about 300 GHz as an example in Figure 3.

3.2 J0129$-$0035

This quasar was selected from the SDSS stripe 82 with $m_{\rm 1450\AA}$ = 22.28 mag (Jiang et al. 2009), making it the faintest source among our targets. The flux density at 250 GHz is $2.73\pm 0.49$ mJy (Wang et al., 2011b) yielding a FIR luminosity that is comparable to that of J2054$-$0005 but much smaller than that of J2310+1855. The bottom left panel of Figure 1 shows the CO (2$-$1) velocity-integrated map integrated from $-$126 to 85 km s*-1* with Gaia calibrated source position marked as white cross. The uncertainties of the Gaia calibrated position are 149 mas in RA and 114 mas in Dec (Wang et al. 2018 in preparation). A 2D Gaussian fit to the intensity map indicates that it is a point source. The inferred source position from the CO (2$-$1) line is marginally consistent with those from [C II] and its underlying dust continuum emission (Wang et al., 2013) and the Gaia calibrated position. A Gaussian fit to the CO (2$-$1) peak spectrum provides a redshift of 5.7783 $\pm$ 0.0004 and a line width of FWHM = 204 $\pm$ 45 km s*-1*, and these results are consistent with measurements derived from CO (6$-$5) and [C II] (Wang et al. 2011b; Wang et al. 2013).

We did not detect J0129$-$0035 in any of the Herschel SPIRE bands. The 3$\sigma$ values upper limits are listed in Table 3.3.

3.3 J2054$-$0005

This source was discovered by Jiang et al. (2008) from SDSS stripe 82 with $m_{\rm 1450\AA}$ = 20.60 mag. The bottom right panel of Figure 1 presents the CO (2$-$1) intensity map, integrated over the velocity range from $-$18 to 201 km s*-1*, where we also plotted the Gaia calibrated position as white cross. The uncertainties of the Gaia calibrated position are 25 mas in RA and 27 mas in Dec. We fitted a 2D Gaussian to the CO (2$-$1) velocity-integrated map and found that it is a point source with a peak value of 0.055 $\pm$ 0.008 Jy km s*-1* beam*-1*. The derived CO (2$-$1) source position is 20h54m06.5001s $\pm$ 284 mas of RA and $-$00d05m14.0622s $\pm$ 267 mas of Dec, where the errors are the sum of the fitting-type error and the position error caused by thermal noise (Reid & Honma, 2014). This position is $\sim 0\farcs 4$ away from those inferred from [C II] (position errors are $\sim$ 20 mas in both RA and Dec) and its underlying dust continuum emission (position errors are $\sim$ 10 mas in both RA and Dec) in Wang et al. (2013) and the Gaia calibrated position. High sensitivity observations of the molecular CO lines are needed to check this tentative offset. A Gaussian fit to the line spectrum yielded a redshift of 6.0394 $\pm$ 0.0004, consistent with CO (6$-$5) and [C II] measurements (6.0379 $\pm$ 0.0022 and 6.0391 $\pm$ 0.0001 respectively; Wang et al. 2010; Wang et al. 2013).

4 Analysis and Discussion

4.1 Gas distribution and gas mass in quasar host galaxies

J0129$-$0035 and J2054$-$0005 are point sources in the CO (2$-$1) line emission, and J2310$+$1855 is marginally resolved with a physical size of (3.51 $\pm$ 1.07) kpc $\times$ (2.33 $\pm$ 1.21) kpc. Figure 2 presents the CO (2$-$1) line spectra compared to CO (6$-$5) and [C II] emission lines of the three targets discussed in this paper. The line widths and redshifts are consistent with each other.

The CO (2$-$1) line fluxes of the three sources allow us to estimate molecular gas masses. Following Solomon et al. (1992) and assuming a conversion factor $\alpha_{\rm CO}\sim 0.8\ M_{\odot}\ (\rm K\ km\ {s}^{-1}\ {pc}^{2})^{-1}$ and $L^{{}^{\prime}}_{\rm CO(2-1)}\approx L^{{}^{\prime}}_{\rm CO(1-0)}$ (Carilli & Walter, 2013), we estimate the CO luminosities to be (21.3 $\pm$ 0.2) $\times$ 106 $L_{\odot}$, (4.0 $\pm$ 0.6) $\times$ 106 $L_{\odot}$ and (6.6 $\pm$ 0.9) $\times$ 106 $L_{\odot}$ and the gas masses to be (4.3 $\pm$ 0.4) $\times$ 1010 $M_{\odot}$, (0.8 $\pm$ 0.1) $\times$ 1010 $M_{\odot}$ and (1.3 $\pm$ 0.2) $\times$ 1010 $M_{\odot}$ for J2310$+$1855, J0129$-$0035, and J2054$-$0005, respectively. Considering the calibration errors, the gas masses are consistent with those derived from the CO (6$-$5) line emission in Wang et al. (2011a) (Table 3.3), where they assumed $L^{{}^{\prime}}_{\rm CO(6-5)}/L^{{}^{\prime}}_{\rm CO(1-0)}=0.784$ from a Large Velocity Gradient model of J1148+5251 at $z=6.42$ and $\alpha_{\rm CO}\sim 0.8\ M_{\odot}\ (\rm K\ km\ {s}^{-1}\ {pc}^{2})^{-1}$.

4.2 Dust temperature, dust mass and star formation rate of J2310$+$1855

In this section, we perform a SED fit of J2310$+$1855, for which we have photometry in multiple bands from both this paper and archival data which are all listed in Table 3.3. We should declare first that our measurements except those from ALMA are all from the unresolved images. We can assume that those measurements represent all the emission from the quasar-galaxy system. However, the ALMA measurements are from marginally or fully resolved data, which may resolve out the low surface density region and result in flux loss. By comparing two flux density measurements in two similar bands - ALMA 263 GHz and IRAM 250 GHz, we find consistent values, which may indicate that almost all the dust emission comes from the central part of the quasar host galaxy. In this case, we may ignore the possible flux loss here. We use the $emcee$ package (Foreman-Mackey et al., 2013) to do SED fit to the emission from both star formation and the central AGN based on the method in Leipski et al. (2013, 2014). In our fitting procedure, we follow Shangguan et al. (2018) to create a likelihood function for the flux upper limits.

We utilize a power law to present the UV/optical emission from the accretion disk and a Clumpy AGN Tori in a 3D geometry (CAT 3D) model (Hönig & Kishimoto, 2017) to represent the near-infrared (NIR) and middle-infrared (MIR) contributions from the AGN dust torus, which considers different sublimation radii of various particles. We use two CAT 3D models, one without and one with a polar wind. Unlike Leipski et al. (2013, 2014), we do not include a NIR hot black body component from the inner region of the AGN dust torus, but rather an interpolated CAT 3D model described above. We adopt two scenarios for the FIR dust continuum emission. The first scenario uses a modified black body (MBB) with fixed emissivity index $\beta$ of 1.6 which is a typical value from the high-redshift quasar sample of Beelen et al. (2006). The second scenario uses the Draine & Li (2007) dust model which is a linear contribution of PDR and diffuse component models that is described in Equation 1:

[TABLE]

$q_{\rm PAH}$ is the polycyclic aromatic hydrocarbon (PAH) mass fraction, $U_{\rm min}$ and $U_{\rm max}$ present the minimum and maximum starlight intensity relative to the local interstellar radiation field, $\alpha$ is the power law index of the starlight distribution, $\gamma$ is the dust fraction heated by starlight with $U_{\rm max}>U>U_{\rm min}$, $\nu p_{\nu}^{(0)}$ is the power produced by a single starlight intensity of $U_{\rm min}$ (diffuse component), and $\nu p_{\nu}$ is the power heated by starlight with intensity in the range from $U_{\rm min}$ to $U_{\rm max}$ (PDR component). The parameter details are described in Draine & Li (2007). As we lack data for the NIR range in which the PAH features appear, we cannot fully constrain the Draine & Li (2007) dust model. The results of our fits are given in Table 4.

As shown in da Cunha et al. (2013), the radiation of the cosmic microwave background (CMB) has two potential effects on our results. One is the heating on the galaxy dust which is important when the CMB temperature is comparable to the dust temperature. At $z$ = 6, the CMB temperature is 2.73 $\times$ (1 + $z$) = 19.11 K, well below our derived temperature of 40 K. So we can neglect the CMB heating. The other is an extra CMB background which will reduce the observed flux density where it is subtracted. The CMB as a background diminishes the observed flux by a fraction defined in Equation 2:

[TABLE]

$F_{\nu,\rm obs}$ and $F_{\nu,\rm intrinsic}$ are the observed and the intrinsic flux densities. $B_{\nu}[\lambda,\ T_{\rm CMB}(z)]$ and $B_{\nu}[\lambda,\ T_{\rm dust}(z)]$ are Planck functions with CMB and galaxy dust temperatures at a given wavelength. Figure 5 shows this effect at $z=6$. This figure shows that the longer the wavelength and the cooler the dust temperature, the higher the CMB influence is as a background. In our analysis, we correct for this effect and get a dust temperature $\sim 40$ K for the first scenario of the FIR component. And for the second scenario, we first correct the CMB effect as a background assuming the dust temperature of 40 K which is a typical value from the first scenario and then do the fit. We also estimate the dust temperature based on the average starlight intensity with Equation 3 and 4 and assuming the same $\beta$ = 1.6 with the first scenario in order to make a comparison.

[TABLE]

[TABLE]

The two scenarios both give a consistent dust temperature $\sim 40$ K.

The dust temperature ($\sim 40$ K) of J2310+1855 is at the low end of that found for quasars at $z$ $\geq$ 5 (Leipski et al. 2013, 2014). We calculate the FIR luminosity by integrating the FIR dust model (the MBB dust model or the Draine & Li (2007) dust model) from 8 to 1000 $\mu$m, and the results are shown in Table 4. The values from the MBB dust model are (1.4$-$1.6) $\times$ 1013 $L_{\odot}$, which would be smaller by a factor of 1.7 if we use a Chabrier initial mass function (IMF; Chabrier 2003). Assuming a Salpeter IMF and adopting Equation 4 in Kennicutt (1998), we derive SFRs of 2400$-$2700 $M_{\odot}$ yr*-1*. The FIR luminosities from the Draine & Li (2007) dust model are (3.2$-$3.6) $\times$ 1013 $L_{\odot}$, which are about two times higher than those from the MBB dust model. This may be due to the lack of data in the NIR and MIR parts, where we only have one WISE $w_{3}$ data point which is a 2.5$\sigma$ marginal detection and two PACS data points. These data can poorly constrain the CAT 3D dust torus model as well as the Draine & Li (2007) dust model which contains some PAH lines. More NIR and MIR data would allow us to distinguish these models.

We estimate dust masses (1.7 and 1.6 $\times$ 10${}^{9}M_{\odot}$ without and with a polar wind of a CAT 3D torus model) for the first scenario by Equation 5:

[TABLE]

where $F_{\rm 125\mu m}$ is the flux density at 125 $\mu$m in the rest frame, $D_{L}$ is the luminosity distance, $\kappa$125μm = 18.75 cm2 g*-1* is the dust absorption coefficient at 125 $\mu$m (Hildebrand, 1983), and $B_{\nu}$ (125 $\mu$m, $T_{\rm FIR}$) is the Planck value at 125 $\mu$m with FIR MBB temperature $T_{\rm FIR}$. In the second scenario, the dust masses are directly determined by the fitting process, giving 3.8 and 3.9 $\times$ 10${}^{9}M_{\odot}$, respectively. The diffuse component contributes more than $80\%$ of the dust mass and dominates the luminosity at longer wavelengths $>$100 $\mu$m, however the PDR component dominates the luminosity at shorter wavelengths $<$60 $\mu$m.

The gas-to-dust ratios (GDRs) are $26\pm 6$ with dust masses from the first scenario and $11\pm 1$ with dust masses from the second scenario. All these ratios are smaller than the widely adopted value - 100$-$150 for the Milky Way, which considers both warm ($>30$K) and cold dust and hydrogen in ionized, atomic and molecular phases. However, GDR values vary in a wide range in different galaxies. For example, Sandstrom et al. (2013) studied 26 nearby star-forming galaxies and proposed a ratio of 91.2. Devereux & Young (1990) calculated a gas to warm dust ratio of 1080 $\pm$ 70 for 58 spiral galaxies. Baes et al. (2014) derived a gas to cold dust ratio $<14.5$ for an early-type galaxy - NGC 5485. Magdis et al. (2011) calculated a ratio $\sim 75(35)$ assuming a metallicity of $Z$ = 8.8(9.2) of a starburst galaxy GN20 at $z$ = 4.05 by fitting to the local GDR$-Z$ relation. The large uncertainties in GDRs can arise from determinations of both the gas mass and dust mass. The gas is in multiple phases. A gas mass conversion factor $\alpha$ that just includes molecular gas results in an underestimation and the value of $\alpha$ is likely to be different for different types of galaxies. The dust mass is generally calculated from SED modeling, where small uncertainties in the dust temperature can give large errors. Based on Equation 5, if we increase the $T_{\rm FIR}$ from 40 K to 50 K, the dust mass will decrease by a factor of 2.

4.3 The Characteristics of the ISM

Kaufman et al. (1999, 2006) presented a series of PDR models by solving for the radiation transfer, chemistry equilibrium and thermal balance in a PDR layer. Each model can be described by a constant hydrogen nucleus density $n$ in units of cm*-3*, and the incident FUV intensity $G_{0}$ in units of the Habing field (= 1.6 $\times$ 10*-3* erg cm*-2* s*-1*; the local Galactic interstellar FUV field). These models can be recreated using the PDR Toolbox (PDRT222http://dustem.astro.umd.edu/pdrt/). In PDRT modeling, the FIR luminosity is integrated from 30 $\mu$m to infinity. For J2310$+$1855, we estimate the FIR luminosity based on the SED fit result. However, for the remaining two sources, we don’t have detailed SED information, so we assume a MBB with a dust temperature of 47 K (the average value of the high redshift quasar sample in Beelen et al. 2006), and an emissivity index of 1.6. And we also calculate the FIR luminosities of these two targets with a dust temperature of 40 K from our SED fit of J2310+1855. As recommended by Kaufman et al. (1999), we multiply the measured CO line flux by a factor of 2 considering line luminosity from both sides of each cloud (optically thick) when adopting these PDR models.

Figure 6 shows $L_{\textsc{[C\,{\sc ii}]\ }}$/$L_{\rm FIR}$ as a function of $L_{\rm CO(1-0)}$/$L_{\rm FIR}$ for our three targets and the PDR model grid from Kaufman et al. (1999, 2006). We get $L_{\rm CO(1-0)}$ based on an assumption of $L^{\arcmin}_{\rm CO(2-1)}\approx L^{\arcmin}_{\rm CO(1-0)}$ (Carilli & Walter, 2013). In this figure, luminous infrared galaxies (LIRGs; $L_{\rm FIR}>10^{11}L_{\odot}$) and ULIRGs, starburst nuclei and Galactic star forming regions are also plotted. Our targets and other high redshift quasars at $z>4$ in the literature (Maiolino et al. 2005; Iono et al. 2006; Wagg et al. 2012, 2014; Leipski et al. 2013; Stefan et al. 2015; Wang et al. 2016; Venemans et al. 2017b) have FIR luminosities ranging from 2.3 $\times$ 1012 $L_{\odot}$ to 4.1 $\times$ 1013 $L_{\odot}$ and bolometric luminosities in a wide range from 5.7 $\times$ 1012 $L_{\odot}$ to 4.3 $\times$ 1014 $L_{\odot}$ (Priddey & McMahon 2001; Willott et al. 2003; Mortlock et al. 2011; Wang et al. 2013; Wu et al. 2015). These high redshift quasars fall in the same part of the diagram as do the galactic star forming regions with mean and median observed $L_{\textsc{[C\,{\sc ii}]\ }}$/$L_{\rm CO(1-0)}$ values of $\sim 3700$ and $\sim 3500$, respectively ($\sim 4100$ is an empirical value for the local starburst galaxies). In comparing with the PDR models, we constrain the hydrogen nucleus density and the FUV intensity to be a few to ten times 105 cm*-3* and a few times 103 for our sample.

The flux ratios between CO (6$-$5) and CO (2$-$1) are $8.44\pm 1.18$, $10.28\pm 2.41$ and $5.67\pm 1.50$ for J2310+1855, J0129$-$0035, and J2054$-$0005, respectively. They are similar to that ($=7.08\pm 0.99$) in the well studied quasar J1148+5251 at $z=6.42$. This may indicate that the CO excitation is similar in these quasars at $z\sim 6$. With the objects in this paper, there are now eight quasars at $z\sim 6$ detected in CO (2$-$1) line emission, and seven of them have CO (6$-$5) detections (Bertoldi et al. 2003b; Carilli et al. 2007; Wang et al. 2010, 2011a, 2011b, 2016; Stefan et al. 2015). The mean and median values of the flux ratio between CO (6$-$5) and CO (2$-$1) of these quasars are $7.85\pm 0.98$ and 8.42, respectively. These values are larger than that ($=4.41\pm 1.48$) from the central starburst disk in M82 (Weiß et al., 2005) and the mean values ($2.1\pm 1.3$ and $3.8\pm 0.7$) from a (U)LIRG sample in Rosenberg et al. (2015) and a SMG sample in Bothwell et al. (2013), which shows that these $z\sim 6$ quasars are more highly excited starburst associated systems or that a central AGN may contribute additional heating. A detailed SLED analysis towards J2310+1855 including more high-$J$ CO lines is presented in Li et al. (in preparation).

5 Summary

We have reported VLA observations of the CO (2$-$1) line emission in three FIR-luminous $z\sim 6$ quasars. One target (J2310+1855) is marginally resolved at an angular resolution of $\sim 0\farcs 6$, and the other two (J0129$-$0035 and J2054$-$0005) are point sources at $\sim 0\farcs 6$ and $\sim 2\farcs 1$ resolutions respectively. CO (2$-$1) line emission is critical to trace the cool gas and to estimate the molecular gas mass directly. We have increased the number of CO (2$-$1) detected quasars from 5 to 8 at $z\sim 6$. The flux ratios between CO (6$-$5) and CO (2$-$1) of the three targets are consistent with that ($\sim 7$) of J1148+5251 at $z=6.42$. This may indicate that the CO excitation is similar in these $z\sim 6$ quasar host galaxies. The gas masses based on CO (2$-$1) of the three targets with typical masses of (1$-$4) $\times$ $10^{10}M_{\odot}$ are consistent with those derived by CO (6$-$5). Different ISM tracers (CO (2$-$1), CO (6$-$5), and [C II]) show similar line widths and redshifts of the three quasars. PDR analysis yields very intense FUV fields in all three sources, and indicates similar ISM properties as found in local starburst galaxies. We also presented a detailed analysis of the dust continuum of J2310+1855 based on Herschel and ALMA data. An SED fit yields a FIR dust temperature of $\sim 40$ K and a SFR of \sim 2400$$-2700 $M_{\odot}$ yr*-1*, and the derived dust mass is roughly 109 $M_{\odot}$. This is consistent with the strong star formation activity that has been seen in other quasar hosts at $z>5.7$. The IR data and an SED decomposition are critical to separate the emission from both the central AGN and the star formation activity in the host galaxy.

We acknowledge support from the National Key R$\&$D Program of China (2016YFA0400703). R.W. acknowledges supports from the National Science Foundation of China (NSFC) grants No. 11473004, 11533001, and the National Key Program for Science and Technology Research and Development (grant 2016YFA0400703). D.R. acknowledges support from the National Science Foundation under grant number AST-1614213. This work makes use of the ALMA data: ADS$/$JAO.ALMA$\#$2015.1.01265.S. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada), MOST and ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI$/$NRAO and NAOJ. The National Radio Astronomy Observatory is a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc. Facilities: Herschel, VLA, ALMA.