ALMA RESOLVES THE TORUS OF NGC 1068: CONTINUUM AND MOLECULAR LINE EMISSION

NGC 1068 is the prototypical Seyfert 2 galaxy and a prime example for active galactic nuclei (AGNs) unifying schemes since the discovery of polarized optical continuum and broadline emission in this source (Miller & Antonucci 1983; Antonucci & Miller 1985). Its central engine is thought to be hidden behind a screen of obscuring material located in a dusty molecular torus of a few parsecs size. Given the distance to the galaxy (D ∼ 14 Mpc; Bland-Hawthorn et al. 1997), which implies a scale ∼70 pc/'', spatially resolving the molecular/dusty torus has been challenging, as angular resolutions <01 are required.

Single-dish observations in the near- and mid-infrared (NIR and MIR) have shown extended emission around the central engine (Bock et al. 2000; Tomono et al. 2001; Rouan et al. 2004; Galliano et al. 2005; Gratadour et al. 2006, 2015). Higher-spatial resolution interferometric observations in the NIR and MIR identified two components (Jaffe et al. 2004; Weigelt et al. 2004; Raban et al. 2009; Burtscher et al. 2013; López-Gonzaga et al. 2014). First, a compact 0.5–1.4 pc sized core at the AGN with a position angle PA ≃ 140° that contains hot (≃800 K) dust co-spatial with the H₂O megamaser disk (Greenhill et al. 1996, hereafter GR96; Greenhill & Gwinn 1997; Gallimore et al. 2001, hereafter GA01; Gallimore et al. 2004).

A second component of warm (≃300 K) dust extends over 3 to 10 pc along the north-south axis, and may correspond to polar dust emission in the ionization cone. Due to their limited (u, v) coverage, IR interferometric data have to be processed with reconstruction techniques to provide a best-fit model of the putative torus in the (u, v) plane. This limitation can be to a large extent circumvented in the submillimeter using interferometers like ALMA, whose large number of antennas and baselines assure a much more complete coverage of the (u, v)-plane.

We used ALMA in Cycle 0 to image the dust continuum and CO(J = 6 − 5) line emissions at 689 GHz in the circumnuclear disk (CND) of NGC 1068 with a spatial resolution of ∼20 pc (García-Burillo et al. 2014, hereafter GB14). The CND appeared as a 300 pc × 200 pc sized off-centered ellipsoidal ring with two prominent knots located east and west of the AGN. Although significant continuum and line emissions were also detected at the position of the AGN, the insufficient spatial resolution of these observations prevented us from isolating the torus. The new ALMA Cycle 2 observations of the CND of NGC 1068 presented in this Letter have a spatial resolution of ∼4 pc. This factor of ≥20 smaller beam area compared with GB14 allowed us to resolve the CND and image, and for the first time, the dust emission and the distribution and kinematics of molecular gas from the torus of NGC 1068 in the submillimeter wavelength range.

We observed the CO(J = 6 − 5) line and its underlying continuum emission (at 432 μm) in NGC 1068 with ALMA during Cycle 2 using the Band 9 receiver (project ID: #2013.1.00055.S). The projected baselines range from 41 m to 2251 m. The data were calibrated using the ALMA reduction package CASA (McMullin et al. 2007). The calibrated UV tables were exported to GILDAS ²¹ for mapping and CLEANing. One track was observed with 35 antennas during September 2015 using the extended configuration of the array, which allowed us to reach a spatial resolution of 007 × 005 using natural (NA) weighting, and a field of view of 9'' that covers the entire CND of the galaxy. Four spectral windows with a spectral bandwidth of 1.875 GHz were placed at sky frequencies of 687.970, 688.888, 691.701, and 706.221 GHz. The second sub-band was centered at the redshifted frequency of the CO line. To optimize the deconvolution and cleaning, we used polygons to restrict the regions where cleaning components are identified. The polygons were defined for the CO line image for each channel based on the detection of significant (≥5σ) emission identified in our Cycle 0 project (GB14). The single polygon used in the cleaning of the continuum image was based on a similar detection criterion using the Cycle 0 data. The sensitivity in the line data cube is 4 mJy beam⁻¹ in channels of 12.8 km s⁻¹ width. We obtained images of the continuum emission by averaging in each sub-band the channels free of possible line emission. The resulting maps were combined to produce an image of the continuum emission at 694 GHz. The corresponding point source sensitivity for the continuum is 0.5 mJy beam⁻¹. Given the (u, v)-coverage, we expect to filter a significant amount of emission on scales ≥2''.

The phase tracking center of the central field is the same as in GB14 (α₂₀₀₀, δ₂₀₀₀) = (02^h42^m40.771^s, −00°00'4784) and the reference velocity is v_o(HEL) = 1136 km s⁻¹.

(1) The CND and the torus. Figure 1 shows the continuum emission of the CND imaged by the new ALMA observations. While ∼90% of the total continuum flux is spread over a large number of faint emission clumps dispersed throughout the CND, the strongest emission peak corresponds, within 002 (∼1–2 pc), to the position of the AGN (see Figures 1(a), (b)). The extended emission connecting the AGN with the CND in the map of GB14 is likely resolved out in the present data. The AGN is identified as the source S1 in the radio-continuum maps of Gallimore et al. (1996) and Gallimore et al. (2004).²²

Zoom In Zoom Out Reset image size

We have used the GILDAS task UV_FIT to find the best fit to the visibilities inside a radius r ≤ 02 (14 pc) from the AGN. The fit indicates that half of the total flux in this region stems from a spatially resolved ellipsoidal disk, characterized by a deconvolved diameter of 7 ± 1 pc oriented along PA^major = 142° ± 23°. Hereafter, we refer to this fitted disk as the dust torus. The orientation of the dust torus is similar to that of the H₂O maser disk (PA^maser = 140° ± 5°; GR96, GA01). The minor-to-major axis ratio of the dust torus is ≃0.8 ± 0.1, an indication that the disk is spatially resolved along its minor axis. Moreover, the residuals of the fit show polar emission that extends south to west out to a distance of ∼10 pc relative to the torus major axis (Figure 1(b)). The enhanced spatial resolution of the uniformly weighted data set better illustrates the goodness of the fit (Figure 1(c)).

Figure 2(a) shows the CO(6–5) intensity map of the CND obtained by integrating the line emission in the velocity interval $v-{v}_{{\rm{o}}}=[-250,+250]$ km s⁻¹. The CO emission from the east and west knots of the CND is comparatively more prominent than in the continuum image shown in Figure 1, which likely reflects the higher percentage of flux recovered for the line on scales ≥2''.

Zoom In Zoom Out Reset image size

Figure 2(b) shows the mean-velocity field of the gas in the CND. The gas kinematics derived from the new observations are compatible with the picture drawn from GB14: gas motions in the CND show a superposition of outward radial flows on rotation. The molecular outflow is responsible for tilting the kinematic major axis by ≃70° from PA ≃ 260° in the outer galaxy disk to PA ≃ 330° in the CND.

Furthermore, ALMA detects the CO emission from a spatially resolved lopsided disk located at the AGN, as illustrated in Figure 3. The size (diameter ≃ 10 ± 1 pc, minor-to-major axis ratio ≃0.5±0.1) and orientation $({{PA}}^{\mathrm{major}}=112^\circ \pm 20^\circ )$ of the CO disk, hereafter referred to as the CO torus, are close to those derived above for the dust torus. The orientation of the CO torus is virtually identical to that of the radio-continuum source S1 resolved in the VLBA map of Gallimore et al. (2004). Figure 3 shows no significant CO counterpart for the polar emission of the dust torus.

Zoom In Zoom Out Reset image size

(2.) The mass of the torus. From the continuum flux of the dust torus at 432 μm, ≃13.8 ± 1 mJy, we estimate a dust mass of ${M}_{{\rm{dust}}}^{{\rm{torus}}}\sim 1600\;{M}_{\odot }$, assuming that a maximum of ∼18% of the flux at 694 GHz may come from other mechanisms different than thermal dust emission (GB14). To estimate ${M}_{{\rm{dust}}}^{{\rm{torus}}}$ we used a modified black-body function with a dust temperature T_dust ≃ T_gas = 150 K, an emissivity index β = 2, and a value for the dust emissivity κ_{694 GHz} = 0.34 m² kg⁻¹. The values of T_dust and κ_{694 GHz} are based on the analysis of the line excitation and SED fitting done in GB14 and Viti et al. (2014), and on the adoption of a value for κ_{352 GHz} = 0.0865 m² kg⁻¹ (Klaas et al. 2001). Assuming that the gas-to-dust ratio in the central 2 kpc of NGC 1068 is ≃60 (GB14), we derive a molecular gas mass for the torus ${M}_{{\rm{gas}}}^{{\rm{torus}}}=(1\pm 0.3)\times {10}^{5}\;{M}_{\odot }$.

We combined the ALMA Band 9 continuum thermal flux with the NIR/MIR continuum and MIR interferometry data obtained by several groups in apertures ≤005–03 to construct the nuclear SED of NGC 1068 from 1.65 μm to 432 μm (Marco & Alloin 2000; Tomono et al. 2001; Weigelt et al. 2004; López-Gonzaga et al. 2014; this work). We also included as an upper limit the ALMA Cycle 0 Band 7 continuum flux from GB14. The 8–13 μm spectrum is the sum of components 1 and 2 in López-Gonzaga et al. (2014). We note that this nuclear SED, shown in Figure 4, corresponds to the smallest spatial scales probed to date for NGC 1068.

Zoom In Zoom Out Reset image size

We fitted the nuclear SED with the CLUMPY torus models of Nenkova et al. (2008a, 2008b), adopting a Bayesian approach using the BayesClumpy tool (Asensio Ramos & Ramos Almeida 2009). This scheme was applied to model the SED of other Seyferts (Alonso-Herrero et al. 2011; Ramos Almeida et al. 2011). Using this approach we can specify a priori information about the six model parameters. We restricted the torus angular width σ to be 40° < σ < 60° to match the opening angle of the ionization cones (80°; Das et al. 2006). The torus radial extent, defined by the outer–to–inner radius ratio Y, is in the range 5 < Y < 30 based on the torus size measured by ALMA at 432 μm. In this first fit, we restricted the inclination angle of the torus to 60° < i < 90° based on the H₂O maser detection for this galaxy (GR96 ) We translate the results into the corresponding SED shown in Figure 4. The NIR fluxes are underpredicted, an indication that we need to include some contribution from starlight or a distinct component to reproduce the NIR bump (e.g., Kishimoto et al. 2012; Lira et al. 2013).

To compare our torus model with direct measurements from the ALMA observations presented here, we estimated from the fitted parameters the gas mass ${M}_{{\rm{gas}}}^{{\rm{torus}}}$, the outer radius R_torus, and the inclination i of the torus. Using the median values of the fitted parameters, we derive ${M}_{{\rm{gas}}}^{{\rm{torus}}}=0.9{\pm }_{0.4}^{0.2}\times {10}^{5}\;{M}_{\odot }$, R_torus = 4 ± 1 pc, and $i=66^\circ {\pm }_{4}^{9}$. These values are consistent with those derived in Section 3 from the Band 9 observations alone. Nevertheless, if we allow the viewing angle of the torus to have smaller values in the range 30° < i < 90° (see Section 5), the best-fit median values are ${M}_{{\rm{gas}}}^{{\rm{torus}}}=0.6{\pm }_{0.2}^{0.4}\times {10}^{5}\;{M}_{\odot }$, R_torus = 2 ± 1 pc, and $i=34^\circ {\pm }_{2}^{5}$.

Figure 5 shows that ALMA has spatially resolved the molecular gas kinematics in the torus. Although the mean velocities derived from CO, shown in Figure 5(a), span a moderate range ($v-{v}_{{\rm{o}}}=[-30,+30]$ km s⁻¹), we identify spatially resolved velocity gradients within the torus. Overall, the line emission appears redshifted to the south and blueshifted to the north. Therefore, the apparent kinematic major axis of the CO disk, at PA^kin ≃ 180°, is tilted by a significant angle (≃53° ± 15°) relative to the morphological major axis of the CO/dust torus, oriented at a weighted average PA^major ≃ 127° ± 15° (see Section 3). Taken at face value, the tilt between both axes indicates that strong non-circular motions are superposed to the slow rotation pattern of the CO torus.

Zoom In Zoom Out Reset image size

Figure 5(b) shows the position–velocity (p–v) diagram taken along the morphological major axis of the CO torus (PA^M = 112° ± 20°). The bulk of the CO emission in the torus spans a wide range of velocities: $v-{v}_{{\rm{o}}}\simeq [-80,80]$ km s⁻¹. Figure 5(b) compares the observed kinematics with the velocity range expected from the H₂O sub-Keplerian rotation curve of Gallimore et al. (1996) projected along the CO major axis (${v}_{{\rm{rot}}}\propto {v}^{-\alpha }$, with α = 0.31; see also Lodato & Bertin 2003). Even in the most favorable case of a purely Keplerian extrapolation (α = 0.5), the CO velocities are still a factor ≥2 lower compared with the H₂O-based kinematics for an edge-on disk (i = 90°). To reach a better agreement, we would need to incline the CO torus down to i = 30°–40°, i.e., close to the inclination of the galaxy disk (GB14). Although these small values for i can be discarded for the mostly edge-on maser disk (at r ≤ 1 pc), we cannot rule out that the CO torus is tilted toward a more face-on orientation (at r ≥ 2−4 pc). Furthermore, this scenario is still within the range of solutions found in Section 4 for the clumpy torus. However, as argued below, the intricate kinematics of the CO torus cannot be attributed alone to the changing orientation of a purely rotating disk.

Figures 5(a), (b) show that the CO emission is redshifted on the southeast side of the major axis at Δx ≃ −006. The emission is nevertheless blueshifted northwest at Δx ≃ +006. This velocity gradient, measured on scales Δx ≃ ±4 pc from the AGN, is inverted relative to that measured for the H₂O maser disk at Δx ≃ ±1 pc. The apparent counterrotation of the outer disk relative to the inner disk could be caused by an enhancement of the gas turbulence in the torus. A rough estimate of the ratio of the deprojected mean velocity (v_mean) to the velocity dispersion (σ) of the gas along the major axis yields very low values: v_mean/σ < 1. However, Figure 5(c) shows that a comparable spatially resolved velocity gradient is measured along the minor axis p–v diagram. The tilted disk scenario or the simple addition of turbulence to rotation cannot alone explain the observed pattern along the minor axis, an indication that non-circular motions in the torus need to be invoked.

Papaloizou & Pringle (1984) first predicted that non self-gravitating AGN tori should likely develop low-order non-axisymmetric instabilities on a dynamical timescale. This type of instability, known as the Papaloizou–Pringle instability (PPI), has been studied in isolated self-gravitating tori (Kiuchi et al. 2011; Korobkin et al. 2013) and also in tori perturbed by the accretion of matter coming from the outer boundary of the system (Dönmez 2014). These new simulations show the growth of a long-lasting m = 1 mode related to the onset of the PPI. The m = 1 mode generates a lopsided distribution, as well as enhanced turbulence and strong non-circular motions in the torus gas. Evidence of inflowing gas connecting the CND with the torus has been found in the 2.12 μm H₂ lines by Müller Sánchez et al. (2009). The remarkably lopsided distribution and complex kinematics of the ∼7–10 pc CO torus of NGC 1068 seen by ALMA could be the footprint of the PPI in this Seyfert. Alternatively, the velocity gradient observed along the minor axis could be explained by gas being entrained in the outflow, as similarly detected by Cecil et al. (2002) in high-ionization lines close to the AGN. This radial shift is noticeably reversed farther out at the CND where the CO outflow follows the pattern of low-ionization lines (GB14).

We mapped with ALMA the dust continuum at 432 μm and the CO(6–5) molecular line emission in the CND of NGC 1068 with a spatial resolution of ∼4 pc. These observations allowed us to resolve the CND and image the dust emission and also the distribution and kinematics of molecular gas from a 7–10 pc diameter disk. This is the first detection of a submillimeter counterpart of the putative torus of NGC 1068. As expected for the contribution of a comparatively cooler component, the CO/dust torus extends on spatial scales that are twice larger compared with the brightest MIR sources detected by the VLTI at the position of the central engine of NGC 1068 (Burtscher et al. 2013; López-Gonzaga et al. 2014).

The ALMA measurements have been combined with the NIR/MIR continuum and MIR interferometry data available from 1.65 μm to 432 μm. We fitted the nuclear SED with the CLUMPY torus models of Nenkova et al. (2008a, 2008b). The mass and the radius of the best-fit torus are consistent with the values derived from the ALMA data alone: ${M}_{{\rm{gas}}}^{{\rm{torus}}}=(1\pm 0.3)\times {10}^{5}$ M_☉ and R_torus = 3.5 ± 0.5 pc.

We spatially resolved the kinematics of molecular gas in the torus. The overall slow rotation pattern of the disk is perturbed by strong non-circular motions and enhanced turbulence. Furthermore, while the maser condition of the H₂O line suggests that the inner molecular disk should be oriented mostly edge-on at r ≤ 1 pc, we found evidence that the molecular torus could be tilted toward lower inclination angles (i = 34°–66°) at r ≥ 2−4 pc. The lopsided morphology and complex kinematics of the NGC 1068 torus could be the signature of the PPI instability predicted by Papaloizou & Pringle (1984). This instability was foreseen to rule the dynamical evolution of AGN tori.

The results of this work are relevant to caution against the use of CO lines as a direct tracer of the supermassive black hole mass (M_BH) in galaxies (e.g., Davis et al. 2013). In cases like NGC 1068 where even mid-J CO lines trace the outer boundary of a kinematically distorted tilted torus, careful modeling would be mandatory to reliably estimate M_BH.

S.G.B., A.U., L.C., I.M., and A.F. acknowledge support from Spanish grants AYA2012-32295 and AYA2013-42227-P. We thank the referee, Makoto Kishimoto, for his constructive review. This paper makes use of the following ALMA data: ADS/JAO.ALMA#2013.1.00055.S. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan); NRC (Canada), NSC, 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.