# Spatially Resolved [C ii] Emission in SPT0346-52: A Hyper-starburst Galaxy Merger at z ~ 5.7

Katrina C. Litke, Daniel P. Marrone, Justin S. Spilker, Manuel Aravena, Matthieu Bethermin, Scott Chapman, Chian-Chou Chen, Carlos de Breuck, Chenxing Dong, Anthony Gonzalez, et al.

► To cite this version:

Katrina C. Litke, Daniel P. Marrone, Justin S. Spilker, Manuel Aravena, Matthieu Bethermin, et al. Spatially Resolved [C ii] Emission in SPT0346-52: A Hyper-starburst Galaxy Merger at z ~ 5.7. The Astrophysical Journal Letters, 2019, 870 (2), pp.80. 10.3847/1538-4357/aaf057. hal-02015894

HAL Id: hal-02015894

https://hal.science/hal-02015894v1

Submitted on 13 Aug 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Spatially Resolved [C II] Emission in SPT0346-52: A Hyper-starburst Galaxy Merger at $z \sim 5.7$

Katrina C. Litke<fnref n="1" />, Daniel P. Marrone<fnref n="1" />, Justin S. Spilker<fnref n="2" />, Manuel Aravena<fnref n="3" />, Matthieu Béthermin<fnref n="4" />, Scott Chapman<fnref n="5" />, Chian-Chou Chen<fnref n="6" />, Carlos de Breuck<fnref n="6" />, Chenxing Dong<fnref n="7" />, Anthony Gonzalez<fnref n="7" />, Thomas R. Greve<fnref n="8" />,<fnref n="9" />, Christopher C. Hayward<fnref n="10" />, Yashar Hezaveh<fnref n="11" />, Sreevani Jarugula<fnref n="12" />, Jingzhe Ma<fnref n="13" />, Warren Morningstar<fnref n="11" />, Desika Narayanan<fnref n="7" />,<fnref n="8" />,<fnref n="14" />, Kedar Phadke<fnref n="12" />, Cassie Reuter<fnref n="12" />, Joaquin Vieira<fnref n="12" />, and Axel Weiss<fnref n="15" />

Received 2018 September 20; accepted 2018 November 9; published 2019 January 10

## Abstract

SPT0346-52 is one of the most most luminous and intensely star-forming galaxies in the universe, with $L_{FIR} > 10^{13} L_{\odot}$ and $\Sigma_{SFR} \approx 4200 M_{\odot} \text{yr}^{-1} \text{kpc}^{-2}$. In this paper, we present $\sim 0.1'15$ ALMA observations of the [C II] 158 $\mu$m emission line in this $z = 5.7$ dusty star-forming galaxy. We use a pixellated lensing reconstruction code to spatially and kinematically resolve the source-plane [C II] and rest-frame 158 $\mu$m dust continuum structure at $\sim 700$ pc ($\sim 0.1'12$) resolution. We discuss the [C II] deficit with a pixellated study of the $L_{[C II]}/L_{FIR}$ ratio in the source plane. We find that individual pixels within the galaxy follow the same trend found using unresolved observations of other galaxies, indicating that the deficit arises on scales $\lesssim 700$ pc. The lensing reconstruction reveals two spatially and kinematically separated components ($\sim 1$ kpc and $\sim 500$ km s$^{-1}$ apart) connected by a bridge of gas. Both components are found to be globally unstable, with Toomre Q instability parameters $\ll 1$ everywhere. We argue that SPT0346-52 is undergoing a major merger, which is likely driving the intense and compact star formation.

**Key words:** galaxies: high-redshift – galaxies: starburst

## 1. Introduction

Dusty star-forming galaxies (DSFGs) are among the most infrared-luminous ($L_{FIR} > 10^{12} L_{\odot}$, where $L_{FIR}$ is the luminosity integrated from 42.5 to 122.5 $\mu$m, Helou et al. 1988) and intensely star-forming ($\Sigma_{SFR} \sim 1000 M_{\odot} \text{yr}^{-1} \text{kpc}^{-2}$) galaxies in the universe (Greve et al. 2012; Casey et al. 2014; Ma et al. 2015). The origin of DSFGs is heavily debated (Sanders et al. 1988; Engel et al. 2010; Narayanan et al. 2010; Hayward et al. 2012, 2013; Chen et al. 2015; Oteo et al. 2016). It has been theorized that strong starbursts like DSFGs will eventually form the massive elliptical galaxies seen in the centers of galaxy clusters at $z < 1.5$ (Thomas et al. 2005, 2010; Kodama et al. 2007; Kriek et al. 2008; Zirm et al. 2008; Gabor & Davé 2012; Hartley et al. 2013; Toft et al. 2014).

Recent surveys with the 2500 deg$^2$ South Pole Telescope (SPT; Vieira et al. 2010; Carlstrom et al. 2011; Mocanu et al. 2013) have greatly expanded the number of known, bright, strongly lensed DSFGs up to $z \sim 7$ (Strandet et al. 2017; Marrone et al. 2018). One of the most extreme DSFGs discovered by SPT, with the highest $L_{FIR}$ and $\Sigma_{SFR}$ in the SPT sample, is SPT-S J034640-5204.9 (hereafter SPT0346-52). SPT0346-52 has been studied at

radio, infrared, optical, and X-ray wavelengths (Hezaveh et al. 2013; Vieira et al. 2013; Weiß et al. 2013; Gullberg et al. 2015; Ma et al. 2015, 2016; Spilker et al. 2015; Aravena et al. 2016; Spilker et al. 2016; Strandet et al. 2016).

SPT0346-52 is a gravitationally lensed galaxy at $z = 5.6559$ (Weiß et al. 2013), with lensing magnification $\mu = 5.6 \pm 0.1$ (Spilker et al. 2016). It has an apparent $L_{FIR} = 1.23 \times 10^{14} L_{\odot}$ (Spilker et al. 2015) and specific star formation rate sSFR $> 15.7 \text{ Gyr}^{-1}$ (Ma et al. 2015). This galaxy's star formation rate density, $\Sigma_{SFR}$, is $4200 M_{\odot} \text{yr}^{-1} \text{kpc}^{-2}$, one of the highest of any known galaxy (Hezaveh et al. 2013; Ma et al. 2015, 2016; Spilker et al. 2015).

Hezaveh et al. (2013) and Spilker et al. (2016) performed gravitational lensing reconstructions of the 860 $\mu$m continuum emission in SPT0346-52. This work was continued in Spilker et al. (2015), which reconstructed the CO(2–1) line in 200 km s$^{-1}$ channels. This lensing reconstruction showed that gas with velocities blueward of $-100$ km s$^{-1}$ was spatially offset from the rest of the emission, but it was unable to distinguish between a merging system of galaxies or a rotation-dominated system due to insufficient spatial resolution ($\gtrsim 0.1'5$).

<footnote n="1">¹Steward Observatory, University of Arizona, 933 North Cherry Avenue, Tucson, AZ 85721, USA; klcitke@email.arizona.edu</footnote>

<footnote n="2">²Department of Astronomy, University of Texas at Austin, 2515 Speedway Stop C1400, Austin, TX 78712, USA</footnote>

<footnote n="3">³Núcleo de Astronomía, Facultad de Ingeniería, Universidad Diego Portales, Av. Ejército 441, Santiago, Chile</footnote>

<footnote n="4">⁴Aix Marseille Univ., Centre National de la Recherche Scientifique, Laboratoire d'Astrophysique de Marseille, Marseille, France</footnote>

<footnote n="5">⁵Dalhousie University, Halifax, Nova Scotia, Canada</footnote>

<footnote n="6">⁶European Southern Observatory, Karl Schwarzschild Straße 2, D-85748 Garching, Germany</footnote>

<footnote n="7">⁷Department of Astronomy, University of Florida, Gainesville, FL 32611, USA</footnote>

<footnote n="8">⁸Cosmic Dawn Center (DAWN), DTU-Space, Technical University of Denmark, DK-2800 Kgs. Lyngby, Denmark; Niels Bohr Institute, University of Copenhagen, Juliane Maries vej 30, DK-2100 Copenhagen, Denmark</footnote>

<footnote n="9">⁹Department of Physics and Astronomy, University College London, Gower Street, London WC1E 6BT, UK</footnote>

<footnote n="10">¹⁰Center for Computational Astrophysics, Flatiron Institute, 162 Fifth Avenue, New York, NY 10010, USA</footnote>

<footnote n="11">¹¹Kavli Institute for Particle Astrophysics and Cosmology, Stanford University, Stanford, CA 94305, USA</footnote>

<footnote n="12">¹²Department of Astronomy and Department of Physics, University of Illinois, 1002 West Green Street, Urbana, IL 61801, USA</footnote>

<footnote n="13">¹³Department of Physics and Astronomy, University of California, Irvine, CA 92697, USA</footnote>

<footnote n="14">¹⁴University of Florida Informatics Institute, 432 Newell Drive, CISE Bldg. E251, Gainesville, FL 32611, USA</footnote>

<footnote n="15">¹⁵Max-Planck-Institut für Radioastronomie, Auf dem Hügel 69 D-53121 Bonn, Germany</footnote>
Ma et al. (2016) explored the origin of the high-luminosity surface density and star formation rate and found that the infrared luminosity is star-formation-dominated, with negligible contributions from a central active galactic nucleus (AGN).

In this paper, we use ALMA Band 7 observations of [C II] 158 µm (hereafter [C II]), a fine-structure line of singly ionized carbon, combined with an interferometric lensing reconstruction tool developed by Hezaveh et al. (2016), to study the structure of SPT0346-52. [C II] 158 µm is an ideal tracer of the gas in the interstellar medium (ISM) of galaxies. At 11.26 eV, neutral carbon has a lower ionization potential than neutral hydrogen, so [C II] can be found in many different phases of the ISM and trace regions inaccessible to observations of ionized hydrogen emission. [C II] is the dominant cooling line in far-UV heated gas (Hollenbach et al. 1991), making it an ideal line with which to study the structure of SPT0346-52. By studying the structure of this galaxy in [C II] we can begin to understand what drives the intense star formation rates observed.

The ratio of the [C II] line luminosity to the far-infrared (FIR) continuum luminosity has been observed to decrease as the FIR luminosity increases (e.g., Malhotra et al. 1997; Luhman et al. 1998, 2003; Díaz-Santos et al. 2013; Gullberg et al. 2015), forming the so-called “[C II] deficit”. Several different mechanisms to produce the observed [C II] deficit with respect to $L_{FIR}$ have been proposed. These include charged dust grains in high-UV radiation fields, self-absorption of [C II] or optically thick [C II], saturated [C II] emission in very high density photodissociation regions (PDRs), dust-bounded photoionization regions, and star formation rates driven by gas surface densities (Malhotra et al. 1997; Luhman et al. 1998, 2003; Muñoz & Oh 2016; Narayanan & Krumholz 2017). Pinpointing the origin of this deficit has been difficult, especially because the deficit is not always observed in DSFGs (e.g., Wagg et al. 2010; De Breuck et al. 2014). Spatially resolved studies of the [C II] deficit have recently become possible at high redshift (e.g., Rawle et al. 2014; Oteo et al. 2016) and should be able to provide a more comprehensive look at the gas conditions associated with the deficit. The analysis in this paper allows for a spatially resolved study of the [C II] deficit in SPT0346-52.

The [C II] observations from ALMA are described in Section 2. Section 3 describes the lensing reconstruction code used and the source-plane reconstruction of SPT0346-52. The [C II] deficit and a kinematic analysis of the results are described in Section 4 and further discussed in Section 5. A summary and conclusions are provided in Section 6. Throughout this paper, we adopt the cosmology from Planck Collaboration et al. (2016), with $\Omega_m = 0.309$, $\Omega_\Lambda = 0.691$, and $H_0 = 67.7$ km s$^{-1}$. At $z = 5.6559$, $1'' = 6.035$ kpc (Wright 2006).

## 2. ALMA Observations

ALMA Band 7 observations of SPT0346-52 were carried out on 2014 September 2 and 2015 June 28 (project ID: 2013.1.01231, PI: Marrone). SPT0346-52 was observed twice at different resolutions, with ~5 minutes on-source in both observations. The 2014 September 2 observation used 34 antennae with baselines up to 1.1km. For the 2015 June 28 observation, 41 antennae were used with baselines up to 1.6km, yielding higher resolution. The reference frequency (first local oscillator frequency) was 291.53 GHz. J0334-4008 was used as both the flux calibrator and the phase calibrator on both days. More information about the observations is available in Table 1.

<figure><img src="image_3.png" /><figcaption>**Figure 1.** High-resolution ALMA observation of SPT0346-52. The image shows the velocity-integrated [C II] line. The continuum is shown overlaid as white contours. The synthesized beam ($0'' \times 0''$ 19) is illustrated as the blue ellipse in the lower left corner.</figcaption></figure>

<table><caption>**Table 1**
ALMA Band 7 Observations of SPT0346-52</caption><thead><tr><th>Date</th><th># of Ant.</th><th>Resolution<br>(arcsec)</th><th>PWV<sup>a</sup><br>(mm)</th><th>t<sub>int</sub><sup>b</sup><br>(min)</th><th>Noise Level<sup>c</sup><br>(mJy/beam)</th></tr></thead><tbody><tr><td>2014 Sep 2</td><td>34</td><td>0.26 × 0.23</td><td>0.944</td><td>5.3</td><td>0.44</td></tr><tr><td>2015 Jun 28</td><td>41</td><td>0.19 × 0.17</td><td>1.315</td><td>5.2</td><td>0.54</td></tr></tbody></table>

**Notes.**

<sup>a</sup> Precipitable water vapor at zenith.

<sup>b</sup> On-source integration time.

<sup>c</sup> Root-mean-square noise level in continuum image.

The data were processed using the Common Astronomy Software Applications package (CASA; McMullin et al. 2007) pipeline version 4.2.2. Some additional flagging was carried out before processing the data with the pipeline. Images were made using the CLEAN algorithm within CASA, with Briggs weighting (robust = 0.5). The continuum was subtracted from the line cube using UVCONTSUB (fitorder = 1). The [C II] data were binned to 50 km s$^{-1}$ channels. The observed 158 µm dust emission and integrated [C II] emission are shown in Figure 1.

The observed [C II] spectrum is shown in Figure 2. It shows a profile with two peaks, one redshifted and one blueshifted relative to the [C II] rest frequency. This spectrum is obtained using the observed complex visibilities rather than CLEANed images of the [C II] line. The method used to create the spectrum is described further in the Appendix.

## 3. Lensing Reconstruction

Gravitational lensing is a useful phenomenon for observing faint emission. Lensing conserves the surface brightness of lensed background sources but it increases their apparent sizes, resulting in greater observed flux. However, strong lensing produces multiple distorted images of background sources and studying the intrinsic properties of lensed objects requires correcting for the lensing distortion.
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Spectrum of observed [C II] emission in SPT0346-52. The spectrum was obtained using the observed visibilities and the model continuum visibilities as described in the Appendix.</figcaption></figure>

### 3.1. Pixellated Lensing Reconstruction

To determine the structure of SPT0346-52, we use the pixellated lensing reconstruction code RIPPLES. This code is described in detail in Hezaveh et al. (2016), with the general framework of using pixellated sources described in Warren & Dye (2003) and Suyu et al. (2006). RIPPLES models the interferometric observations of lensed sources. It models the mass distribution in the lensing galaxy and the background source emission while accounting for observational effects such as those due to the primary beam.

Using a pixellated source reconstruction is advantageous, as it does not assume a specific source structure (i.e., the source is not constrained to follow, for example, a Gaussian or Sérsic surface brightness profile). Instead, it has the flexibility to model more complex source structures, especially when high-resolution data are available, because of the large parameter space of the source pixels and the less constraining priors. Using an inherently interferometric code such as RIPPLES also allows us to use all of the data available from an observing session with an interferometer like ALMA.

The model visibilities can be written as a linear matrix equation,

$$V = F(BLS). \qquad (1)$$

Light from the background source, S, is first lensed by the foreground galaxy. Pixels in the image (lensed) plane are mapped back to the source (de-lensed) plane for a given set of lens parameters. The lensing operation is a matrix represented by L in Equation (1), and depends on the mass distribution of the lensing galaxy. The lensed emission is then modified by the primary beam of the telescope (represented by the matrix B). Finally, we take a Fourier transform of the sky emission, F, to obtain complex visibilities of interferometric observations, V. The model visibilities are compared to the observed visibilities via a $\chi^2$ goodness-of-fit test. A Markov Chain Monte Carlo (MCMC) method is used to solve for the lens galaxy mass distribution parameters.

In addition to the lens parameters, there is a regularization term, $\lambda$. The regularization term acts to smooth the source and minimize large gradients between adjacent pixels in the source

plane. This prevents over-fitting of the data, or fitting to the noise in the source-plane image. $\lambda$ is determined by

$$N_s - \lambda \mathrm{Tr}([FL + \hat{C}_s]^{-1}\hat{C}_s^{-1}) - \lambda S^T \hat{C}_s S = 0, \qquad (2)$$

where $N_s$ is the number of source pixels and $C_s$ is the source covariance matrix. $\lambda$ scales an arbitrarily normalized source covariance matrix, $\hat{C}_s$. It is determined for a fixed lens model, rather than being simultaneously fit for with the lens parameters. We fit for $\lambda$, then run the MCMC with RIPPLES. These two steps are repeated until the chains have converged around a most-likely set of parameters.

After modeling the mass distribution of the lensing galaxy, we obtain a pixellated map of the source-plane emission, a model image, and model complex visibilities.

### 3.2. Reconstruction of SPT0346-52

We model the lensing galaxy as a singular isothermal ellipsoid at $z = 0.9$ with an external shear component. The initial parameters are taken from previous lensing reconstructions of SPT0346-52 by Hezaveh et al. (2013) and Spilker et al. (2015). The best-fit model was determined by fitting the 158 µm continuum data because the continuum has a much higher signal-to-noise ratio than the individual line channels. Figure 3 shows a probability density plot of the lens parameters with the results of the MCMC. The determined lens parameters are given in Table 2.

The best-fit model was then applied to the [C II] line in each 50 km s⁻¹ channel. This channel width was chosen to be wide enough to have high signal-to-noise to be able to reconstruct the source in each channel, while being narrow enough to study kinematic features in SPT0346-52. The source regularization, $\lambda$, for the [C II] line was fixed for all 50 km s⁻¹ channels. The original image, model image, and model source are shown for each channel in Figure 4.

Hezaveh et al. (2013) reconstructed the 860 µm continuum of SPT0346-52 only using short baselines, assuming a symmetric Gaussian source profile. Spilker et al. (2015) reconstructed the CO(2–1) line emission using the code VISILENS (Spilker et al. 2016). They used four 200 km s⁻¹ channels and assumed a symmetric Gaussian source-plane structure for each channel. Channels blueward of −100 km s⁻¹ were spatially offset from redder emission, with the same orientation and velocity ranges obtained with the reconstruction of the [C II] line. The −400 km s⁻¹ and +200 km s⁻¹ channels from the parametric reconstruction of CO(2–1) in Spilker et al. (2016) show disks with similar size and orientation as the two components in the [C II] pixellated reconstruction from this work (see Figure 4).

We combined the reconstructed channels from Figure 4 into a map of the [C II] emission (Figure 5). The top image shows the reconstructed continuum emission, while the bottom image shows the velocity-integrated reconstructed [C II] line. There are two lobes in the [C II] emission, with the lower left component much brighter than the upper right component. The continuum emission, arising mostly from dust, is more regular and is roughly elliptical. It is located near the center of the [C II] emission, between the two components, and is less extended than the [C II] emission. The dust continuum emission has been found to be more compact than the [C II] emission in other high-z galaxies. Gullberg et al. (2018) measured [C II]-emitting regions that are 1.6× more extended than the regions with dust
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Triangle plot with the model lens parameters from reconstructing SPT0346-52. *M* is measured in *M*☉ and is the mass enclosed within 10 kpc. *e*<sub>x</sub> and *e*<sub>y</sub> are the x- and y-components of the lens galaxy's ellipticity. *x* and *y* are the offset of the lens center in arcseconds. *γ*<sub>x</sub> and *γ*<sub>y</sub> are the x- and y-components of shear.</figcaption></figure>

continuum emission in four $z \sim 4.5$ DSFGs. Wang et al. (2013a) found dust continuum in 1.2–2.3 kpc regions and [C II] emission in 1.7–3.5 kpc regions in $z \sim 6$ quasars with vigorous star formation in the central region of the quasar host galaxies. Oteo et al. (2016) found similar sizes in the dust continuum and [C II] emission regions of SGP38326. An offset between the brightest [C II] emission and the center of the dust

emission, as is seen in the source-plane reconstruction of SPT0346-52, was also observed in the Seyfert 2 galaxy NGC1068 (Herrera-Camus et al. 2018a). We also find the dust continuum emission to be smoother than the [C II] emission. Smooth dust continuum emission with clumpy [C II] emission was also seen by Oteo et al. (2016) in SGP38326, a pair of interacting dusty starbursts.
<table><caption>**Table 2**
SPT0346-52 Lens Parameters</caption><thead><tr><th>Parameter</th><th>Value</th></tr></thead><tbody><tr><td>logMass[M<sub>☉</sub>]</td><td>11.43 ± 0.02</td></tr><tr><td>Ellipticity x-Component, e<sub>x</sub></td><td>-0.16 ± 0.02</td></tr><tr><td>Ellipticity y-Component, e<sub>y</sub></td><td>0.49 ± 0.06</td></tr><tr><td>Ellipticity, e<sup>a,c</sup></td><td>0.52 ± 0.06</td></tr><tr><td>Position Angle, φ<sub>e</sub> (E of N)<sup>b,c</sup></td><td>72 ± 6°</td></tr><tr><td>Lens x Position, x</td><td>0''076 ± 0''014</td></tr><tr><td>Lens y Position, y</td><td>0''31 ± 0''01</td></tr><tr><td>Shear x-Component, γ<sub>x</sub></td><td>0.035 ± 0.022</td></tr><tr><td>Shear y-Component, γ<sub>y</sub></td><td>-0.093 ± 0.015</td></tr><tr><td>Shear Amplitude, γ<sup>a,c</sup></td><td>0.10 ± 0.02</td></tr><tr><td>Shear Position Angle, φ<sub>γ</sub> (E of N)<sup>b,c</sup></td><td>111 ± 12°</td></tr></tbody></table>

**Notes.**

<sup>a</sup> α = $\sqrt{\alpha_x^2 + \alpha_y^2}$, where α = e or α = γ.

<sup>b</sup> φ<sub>α</sub> = arctan(-α<sub>γ</sub>/α<sub>x</sub>), where α = e or α = γ.

<sup>c</sup> Derived from best-fit parameters.

Figure 6 shows a spectrum of the reconstructed [C II] emission. The spectrum was obtained by summing the flux from the pixels in the source-plane reconstruction, while excluding pixels near the edges of the source plane. There are two clear peaks in the spectrum with similar maximum fluxes. A two-component Gaussian was fit to the spectrum; the results are overlaid on Figure 6. The blue component is centered at -287 ± 22 km s⁻¹ and has an FWHM of 337 ± 22 km s⁻¹. The red component is centered at +158 ± 22 km s⁻¹ with an FWHM of 319 ± 15 km s⁻¹. Similar velocity structure to that seen in the source-plane spectrum was also observed by Spilker et al. (2015), Aravena et al. (2016), Dong et al. (2018), and Y. Apostolovski et al. 2018, (in preparation).

### 3.3. Source-plane Resolution

In order to determine the resolution in the source plane, we created a set of mock visibilities for a point source in the source plane that was lensed by the best-fit lens model. RIPPLES was then applied to the mock visibilities to reconstruct the source and a 2D Gaussian was fit to the reconstructed source. This Gaussian is the effective resolution. This process was repeated with the location of the point source varying throughout the source plane to understand the variation in resolution across the source. We find that in regions away from the central diamond caustic the effective resolution is 0''13 × 0''15, while closer to the diamond caustic the effective resolution decreases to 0''12 × 0''12. Example effective resolution ellipses are shown in blue in Figure 5.

## 4. Analysis

### 4.1. [C II] Deficit

[C II] is usually the brightest coolant line of the ISM. While it can be emitted in a variety of ISM conditions, it is primarily produced in warm, diffuse gas at the edges of PDRs being heated by an external FUV radiation field, such as a star-forming region or AGN (Hollenbach et al. 1991; Malhotra et al. 1997; Luhman et al. 1998; Pineda et al. 2010). Pavesi et al. (2018) calculated that ~85% of [C II] emission in a z = 5.7 DSFG comes from PDRs. One of the more interesting aspects of the [C II] line is the so-called "[C II] deficit," in which the

$L_{[C II]}/L_{FIR}$ ratio has been found to decrease at high $L_{FIR}$, though this is not always the case. The deficit is often associated with AGN activity, though not all AGNs have a [C II] deficit (Sargsyan et al. 2012). Farrah et al. (2013) also showed that the deficit is stronger in merging systems, with no clear dependence on the presence of an AGN. The deficit was found to be strongest in the AGNs with the highest central starlight intensities, rather than those with the highest X-ray luminosities at low redshift (Smith et al. 2017). This is further supported by Lagache et al. (2018), who found that the [C II] deficit is correlated with the interstellar radiation field in their simulations. In resolved [C II] studies of the Orion Nebula in our galaxy and other DSFGs, the [C II] deficit has been shown to be strongest in regions with higher star formation rates (Goicoechea et al. 2015; Oteo et al. 2016).

Several mechanisms have been suggested as a cause for the [C II] deficit. The [C II] line may be optically thick or self-absorbed by foreground gas. Enhanced IR emission, from intense star formation or an AGN, can also lead to a deficit (Malhotra et al. 1997; Luhman et al. 1998). More recently, Narayanan & Krumholz (2017) proposed that increased surface densities in clouds and increased star formation rates cause a rise in the fraction of gas that is CO-dominated, rather than [C II]-dominated, leading to a [C II] deficit. In addition, Díaz-Santos et al. (2017) found a correlation between the UV-flux-to-gas-density ratio, $G/n_H$, and $\Sigma_{IR}$. They found a critical surface density, $\Sigma_{IR}^* \approx 5 \times 10^{10} L_{\odot} \text{ kpc}^{-2}$, below which $G/n_H$ remains constant. Above $\Sigma_{IR}^*$, they found that $G/n_H$ increases. They argued that the relation between $G/n_H$ and $\Sigma_{IR}$ links kiloparsec-scale galaxy properties to those of individual PDRs. Herrera-Camus et al. (2018a) also found a critical surface density, $\Sigma_{FIR} \gtrsim 10^{11} L_{\odot}$, above which the $L_{[C II]}/L_{FIR}$ ratio decreases, but with increased scatter. The [C II] deficit has also been found to correlate directly with $G/n_H$ (e.g., Malhotra et al. 1997).

Using the pixelated lensing reconstruction, we have resolved maps of the source-plane continuum and [C II] emission. This allows us to obtain a resolved map of the $L_{[C II]}/L_{FIR}$ ratio and probe the [C II] deficit to smaller scales than has previously been possible at high redshift.

In order to study the [C II] deficit, we assume that the 158 µm continuum flux density, $F_{\text{cont}}$, traces $L_{\text{FIR}}$. The measured continuum flux in each pixel, $F_{\text{cont},i}$, is scaled proportional to the total $L_{\text{FIR}}$, using $L_{\text{FIR}}$ from Gullberg et al. (2015) and corrected for lensing using the magnification from Spilker et al. (2016) such that

$$L_{\text{FIR},i} = F_{\text{cont},i} \left( \frac{L_{\text{FIR}}}{\Sigma F_{\text{cont},i}} \right). \quad (3)$$

This method assumes a constant flux density-to-luminosity ratio, and thus a constant dust temperature. We tested the effect of this assumption by determining the total FIR luminosity for a range of dust temperatures measured in DSFGs from the SPT sample (22 < T<sub>d</sub> < 57). $L_{\text{FIR}}$ was calculated by integrating the spectral energy distribution (SED; modeled by a modified blackbody; see Greve et al. 2012) from 42.5 to 122.5 µm and scaling the SED to go through the flux of SPT0346-52 at 158 µm. Resulting $L_{\text{FIR}}$ values were within a factor of ~2 of the luminosity measured by Gullberg et al. (2015), so the variation caused by variable dust temperatures in the galaxy is within a factor of ~2.
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Reconstruction of [C II] 158 µm emission in 50 km s⁻¹ channels. Left three-panel column: blue velocities. Right three-panel column: red velocities. Within each column: left: observed [C II] from ALMA; center: model sky emission; right: reconstructed source with lensing caustic. Lensed images are 5'' on a side. The source-plane images are 8.5 kpc a side. The light blue ellipses in the corners of the source-plane panels show the effective resolution where the emission is brightest in each channel. The data and model sky emission have one color scale, and images of the reconstructed source have another color scale. These color scales are the same across all channels.</figcaption></figure>

A map of the $L_{[C II]}/L_{FIR}$ ratio is shown in Figure 7. Typical values of the $L_{[C II]}/L_{FIR}$ ratio in the center of SPT0346-52 are around ~1 × 10⁻⁴. This value is consistent with other

ultra-luminous infrared galaxies (ULIRGs) and DSFGs that have the [C II] deficit (e.g., Maiolino et al. 2005; Iono et al. 2006; Oteo et al. 2016; Mazzucchelli et al. 2017; Decarli et al. 2017). The
<figure><img src="image_7.png" /><figcaption>**Figure 5.** Result of pixelated lensing reconstruction of SPT0346-52. Top: map of continuum emission at 158 µm. Bottom: map of integrated [C II] emission. The resolution varies across the source depending on the location relative to the caustics (white lines). The blue, green, and purple ellipses in the corner indicate the 1σ size of the 2D Gaussian fit to determine the effective resolution in the reconstruction at the locations of the dots in the maps of the same color. To see more of how the resolution varies across the source plane, see Figure 4. The color bar units are mJy/pixel.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Figure 6.** Spectrum of reconstructed [C II] emission in SPT0346-52. The solid black line shows the two-component Gaussian fit to the spectrum, with the individual components shown as red and blue dashed lines. The two components are centered at -287 km s⁻¹ and +158 km s⁻¹ relative to z = 5.6559. The purple dashed line shows the spectrum obtained from the (lensed) visibilities in Figure 2.</figcaption></figure>

higher values of $L_{[C II]}/L_{FIR}$ at the edges of the galaxy are due to the low amounts of continuum emission in those regions. Oteo et al. (2016) found a similar mapped distribution in a pair of

<figure><img src="image_9.png" /><figcaption>**Figure 7.** [C II] to FIR luminosity ratio in SPT0346-52. The center region, where the continuum emission is strongest, shows relatively uniform values for $L_{[C II]}/L_{FIR}$ in the center, ~1 × 10⁻⁴. The white contour traces the strongest region of continuum emission. The larger values at the edges of the galaxy are due to the falling continuum emission.</figcaption></figure>

interacting DSFGs and suggested it was due to the different morphology of the [C II] emission compared to the dust continuum emission. As with SPT0346-52, the sources studied by Oteo et al. (2016) do not show evidence for AGN activity (Oteo et al. 2016). The uniformity of the $L_{[C II]}/L_{FIR}$ ratio is similar to the merging system observed by Neri et al. (2014).

Figure 8 shows the $L_{[C II]}/L_{FIR}$ versus $\Sigma_{FIR}$ relation for SPT0346-52, pixels from the lensing reconstruction, other high-z sources, and ULIRGs from the GOALS survey (Díaz-Santos et al. 2013). As noted by Spilker et al. (2016), the $L_{[C II]}/L_{FIR}$ versus $\Sigma_{FIR}$ relation continues to higher values of $\Sigma_{FIR}$ for high-z sources. The tight relation continues to hold true at smaller physical scales (the purple diamonds and green shaded region in Figure 8 are individual pixels in resolved [C II] observations from this work and Oteo et al. 2016), with a similar scatter as previous, galaxy-averaged studies at high redshift.

The spatially resolved [C II] deficit was recently explored in nearby galaxies. Smith et al. (2017) measured $L_{C II}/L_{TIR}$ at 0.2–1.6 kpc scales in the KINGFISH sample and found that the $L_{[C II]}/L_{TIR}$ versus $\Sigma_{SFR}$ relation continues at lower values of $\Sigma_{SFR}$. The relation between $L_{[C II]}/L_{TIR}$ and $\Sigma_{SFR}$ found by Smith et al. (2017) is in good agreement with the star formation rate density and $L_{[C II]}/L_{FIR}$ ratio in SPT0346-52. This relation, and the similar $L_{[C II]}/L_{FIR}$ versus $\Sigma_{FIR}$ trend explored in this work, spans many orders of magnitude. It holds true both for spatially resolved regions and galaxy-averaged values at high and low redshift. The [C II] deficit appears to come from local conditions in the ISM because it continues to hold over smaller physical areas. Gullberg et al. (2018) reached the same conclusion in their resolved study of [C II] emission in four $z \sim 4.5$ DSFGs.

Muñoz & Oh (2016) argued that the [C II] deficit is the result of thermal saturation of the [C II] emission line. The relation of $L_{[C II]}/L_{FIR}$ versus $\Sigma_{IR}$ from Muñoz & Oh (2016) is plotted in Figure 8 for $f_{[C II]} \approx 0.13$ as a black dashed-dotted line. $f_{[C II]}$ is proportional to $L_{[C II]}/L_{FIR}$, and $f_{[C II]} \approx 0.13$ is the fraction of the total gas in a galaxy traced by [C II] for a typical DSFG. However, if we calculate $f_{[C II]}$ (Equation (5) from Muñoz & Oh 2016) for SPT0346-52 using $\alpha_{CO} = 2.2$, using $L'_{CO(2-1)}$
<figure><img src="image_10.png" /></figure>

**Figure 8.** $L_{[C\,II]}/L_{\text{FIR}}$ vs. $\Sigma_{\text{FIR}}$. The magenta star shows the galaxy-summed value for SPT0346-52 from this work. The purple diamonds show the individual pixels from the lensing reconstruction. The size of the diamonds is weighted by the continuum signal-to-noise. The right panel zooms in on the region from the left panel with the pixels from the SPT0346-52 reconstruction (light blue box). The green squares are the resolved DSFGs (SGP38326) from Oteo et al. (2016), while the green shaded region outlines the parameter space occupied by their individual pixels. The turquoise shaded squares represent spaxels from a spatially resolved survey of FIR lines in the $z \sim 0$ SHINING sample (Herrera-Camus et al. 2018a). The pink shaded region represents binned pixels from a spatially resolved study of the [C II] deficit in the KINGFISH galaxy sample (Smith et al. 2017). The blue triangles are SPT DSFGs (Gullberg et al. 2015). The black points are low-redshift star-forming galaxies and ULIRGs from the GOALS survey (Díaz-Santos et al. 2013). The light blue points represent global values from nearby galaxies in the SHINING survey (Herrera-Camus et al. 2018a). The steel blue points are additional nearby galaxies (Brauher et al. 2008; Farrah et al. 2013). The dark blue points are $z > 6$ quasars from Decarli et al. (2018) and Izumi et al. (2018). The gray points are additional high-redshift objects from the literature (Walter et al. 2009; Carniani et al. 2013; Riechers 2013; Wang et al. 2013a; De Breuck et al. 2014; Neri et al. 2014; Riechers et al. 2014; Yun et al. 2015; Díaz-Santos et al. 2016; Pavesi et al. 2018; Hashimoto et al. 2018). The gray circles represent the $L_{[C\,II]}/L_{\text{FIR}}$ ratio in the cores (measured within the continuum spectrum) and annuli (measured between the continuum and [C II] apertures) for three DSFGs at $z \sim 4.5$ in Gullberg et al. (2018). The solid gray lines connect the core of a given DSFG to the annulus of that DSFG. The core values have higher $\Sigma_{\text{FIR}}$ and lower $L_{[C\,II]}/L_{\text{FIR}}$ than their annulus counterparts. The gray solid line represents the relation found by Díaz-Santos et al. (2013) and extended by Spilker et al. (2016). The black dashed-dotted line represents the relation from Muñoz & Oh (2016) for saturation of the [C II] line for $f_{[C\,II]} = 0.13$, the value for typical DSFGs. The gray dashed-dotted line represents the relation from Muñoz & Oh (2016) for $f_{[C\,II]} = 0.2$, the value calculated for SPT0346-52. The solid line is the relation found by Díaz-Santos et al. (2017). The dotted line is the relation fit by Smith et al. (2017). $\Sigma_{\text{SFR}}$ values were converted to $\Sigma_{\text{FIR}}$ following Murphy et al. (2011) and Greve et al. (2012). The dashed line is the fit from Lutz et al. (2016).

from Spilker et al. (2015) and assuming $L'_{\text{CO(1-0)}} = L'_{\text{CO(2-1)}}$, we find that $f_{[C\,II]} = 0.21$. This moves the relation from Muñoz & Oh (2016) above the majority of the pixels in the reconstruction of SPT0346-52 in Figure 8 (gray dashed-dotted line). It should be noted that the other lines shown in Figure 8 are empirical fits to the data.

Herrera-Camus et al. (2018a) looked at the $L_{[C\,II]}/L_{\text{FIR}}$ ratio in the SHINING sample of nearby galaxies, with spatially resolved information for 25 of their galaxies. In Herrera-Camus et al. (2018b), they use a pair of toy models to explore the origin of the [C II] deficit, one with the ISM modeled as having OB stars and molecular gas clouds closely related, and the other with OB associations and neutral gas clouds randomly distributed throughout the ISM. In the former case, the [C II] intensity only weakly depends on $G_0$ (because the ionization parameter reaches a limit, $U \approx 0.01$) and $n_H$ (because the density of the neutral gas exceeds the critical density for collisional excitation of [C II]), and in the latter case, the [C II] intensity is nearly independent of $G_0$ (because photoelectric heating efficiency decreases), while the FIR intensity is proportional to $G_0$ in both scenarios. They conclude that the combination of both scenarios best replicates the observed [C II] deficit, including a critical luminosity surface density of $\Sigma_{\text{FIR}} \approx 10^{10} L_{\odot} \, \text{kpc}^{-2}$ above which the $L_{[C\,II]}/L_{\text{FIR}}$ ratio begins to decline.

## 4.2. Kinematic Analysis

In addition to the [C II] emission map shown in Figure 5, we calculate moment 1 (intensity-weighted average velocity, shown in the top left panel of Figure 9) and moment 2 (intensity-weighted velocity dispersion, shown in the top right panel of Figure 9).

Figure 9) of the reconstructed line. The velocity dispersions in the center of the system reach very high values ($\sigma > 200 \, \text{km s}^{-1}$). Extracting the velocities along the major axis of SPT0346-52 (dashed line in Figure 9, top right panel) reveals two spatially distinct velocity components. This is shown in the position-velocity diagram in the middle panel of Figure 9.

In order to separate the two velocity components seen in Figure 9, we fit the spectrum of each pixel with two Gaussian components. Each Gaussian is assigned to the appropriate galaxy component based on its velocity. The shape of the velocity-integrated [C II] emission in each of the two spatial components is fitted with an elliptical Gaussian distribution. The centers and elliptical FWHM shapes of these components are shown in Figure 10. The spatial distribution of the gas bridge, outlined in purple in Figure 10, is determined by selecting pixels with emission at velocities intermediate to the two main galaxy components.

# 5. Discussion

## 5.1. Merging Galaxies

Estimates of the fraction of DSFGs that have multiple components or are merging are varied. For example, from continuum emission only Spilker et al. (2016) found that only 13% of lensed DSFGs from the SPT sample showed strong evidence of having multiple components, while Bussmann et al. (2015) found 69% of DSFGs with multiple components. If only the intensity-weighted velocity (moment 1) map (see Figure 9) were considered when studying the kinematics of SPT0346-52, this system could appear to be a symmetric, rotating disk. However, only ~40%–80% of merging systems
<figure><img src="image_11.png" /><figcaption>**Figure 9.** Top: moments 1 and 2 of [CII] in SPT0346-52. Middle: position–velocity diagram of the major axis of SPT0346-52. The velocities were extracted along the dashed line shown in the moment 2 map. There are two spatially and kinematically components. Bottom: position–velocity diagram of the bridge connecting the two components of SPT0346-52. The velocities were extracted along the dotted line shown in the moment 2 map, with positions 0 kpc marked by the black x's. The units of the color bars in the position–velocity diagrams are Jy/pixel.</figcaption></figure>

show asymmetric kinematics in their star-forming gas (Hung et al. 2016), so symmetric gas kinematics is not a definitive way to determine that a galaxy is not a merging system.

The lensing reconstruction of SPT0346-52 reveals two separated components (see Figure 9). The centers of these components are separated by ~1 kpc and ~500 km s⁻¹. There is a significant decrease in emission at velocities between the center velocities of the two components, as shown in Figure 6. Both of these components are larger than the effective resolution in this region of the reconstructed source plane, and they are separated by ~2–3 resolution elements. Therefore, these two components are more likely to represent two separate structures, rather than a barely resolved rotating disk.

The two components overlap in the middle of the system, near the peak of the continuum emission. This region of

<figure><img src="image_12.png" /><figcaption>**Figure 10.** Separation of the two components in SPT0346-52. The blue and red ellipses outline the FWHM size of the blue and red components, and the centers are marked with an x. The amplitudes of the red and blue components of the Gaussian fits are shown as the red and blue images. The darker region in the center shows where the components overlap. The purple line outlines the pixels containing the “bridge” connecting the two components. Contours of the continuum emission are overlaid in white for reference.</figcaption></figure>

overlap has a more complex velocity structure and higher velocity dispersions (Figure 9, top panel). Because the overlap region is where there is the most dust continuum emission, the star formation is likely occurring most intensely in that region, as has been observed in other merging systems such as the Antennae Galaxies (Mirabel et al. 1998; Karl et al. 2010). Teyssier et al. (2010) also found that merger-induced star formation is relatively concentrated near the center of merging systems in their hydrodynamic simulations.

In addition to the red and blue components, a position–velocity slice through the more complex velocity structure reveals a bridge of gas connecting the two components. The extraction line is indicated by the dotted line in the moment 2 map in Figure 9, and the position–velocity diagram through this slice is shown in the bottom panel of Figure 9. The location of this bridge feature is also indicated by the purple contour in Figure 10. This structure resembles simulated tidal tails and observed tidal tails, such as in Arp105 (Bournaud et al. 2004), as well as the south tail in the Antennae (NGC 4038/9; Gordon et al. 2001) in position–velocity diagrams. Decarli et al. (2017) also found [C II] emission connecting a quasar host galaxy, PJ308-21, and a companion galaxy, though on much larger scales (25 kpc and 1000 km s⁻¹) than what is observed in SPT0346-52.

Both components of SPT0346-52 have large velocity dispersions (>200 km s⁻¹, determined by the Gaussian fits to the spectra in each pixel). These large turbulent motions can help stabilize disks against gravitational fragmentation (e.g., Westmoquette et al. 2012; Rangwala et al. 2015), see Section 5.2.

Several other merging DSFGs have been observed. For example, Neri et al. (2014) observed [C II] emission in HDF850.1 (z = 5.185) and found two components, one redshifted and one blueshifted, and separated by 2 kpc with radii ~1 kpc. These components are similar in size to the components observed in SPT0346-52. Neri et al. (2014) explored the idea that HDF850.1 was a rotating disk, but concluded that they observed a merger-driven starburst. Rawle et al. (2014) also observed a late-stage merging DSFG at z > 5 (HLS0918) with up to four components separated by <4 kpc. Engel et al. (2010) concluded
using CO observations that most bright DSFGs with $L_{IR} > 5 \times 10^{12} L_{\odot}$ are major mergers. This is consistent with the conclusion drawn from studies of stellar structures (e.g., Chen et al. 2015). At $z \sim 7$, Hashimoto et al. (2018) concluded that the Lyman-break galaxy B14-65666 was a merger-induced starburst galaxy based on the velocity gradient in the [C II] line and a two-component spectrum, whose spatial positions are consistent with two [C II] knots and UV emission peaks. Cosmological hydrodynamic galaxy formation simulations by Narayanan et al. (2015) have shown that many DSFGs have multiple components, though the intense star formation may be driven by stellar feedback rather than major mergers.

Mergers can trigger intense star formation activity without producing an obvious AGN in DSFGs (Wang et al. 2013b). Though many ULIRGs, which have similar $L_{IR}$ as DSFGs and enhanced star formation, have AGN activity that heats the dust and causes their high values of $L_{FIR}$, Younger et al. (2009) found that star formation alone can produce warm IR colors and produce UV radiation that is reradiated by hot dust. About 63% of luminous infrared galaxies (LIRGs) have multiple components (Engel et al. 2011; Haan et al. 2011). Merging ULIRGs, which have higher FIR luminosities than LIRGs and have FIR luminosities more similar to that of SPT0346-52, have small nuclear separations (average 1.2 kpc) and are in later merging systems (Haan et al. 2011). Similarly, SPT0346-52 could be a late-stage merger.

Pavesi et al. (2018) recently observed [C II] in a DSFG similar to SPT0346-52 at $z = 5.667$, COSMOS (FIR-)Red Line Emitter (CRLE), with SFR = 1500 $M_{\odot}$ yr$^{-1}$ and a diameter of 2.7 kpc. They determined that CRLE is an intermediate-stage merger. CRLE has a gas depletion timescale of 45 Myr. For SPT0346-52, we calculate a gas depletion timescale of $31 \pm 10$ Myr by dividing the gas mass from Spilker et al. (2015) by the star formation rate from Ma et al. (2015), similar to the depletion timescale calculated by Aravena et al. (2016) for this system.

An alternative explanation for the kinematic morphology in SPT0346-52 is that it is a rotating galaxy with a clumpy gas disk. Clumpy, rotating disks have been observed in DSFGs (e.g., Hodge et al. 2012; Iono et al. 2016; Dannerbauer et al. 2017; Tadaki et al. 2018). However, Hodge et al. (2016) searched for $\sim 1$ kpc clumps (comparable to the sizes of the clumps in GN20 and the components in SPT0346-52) in luminous DSFGs and found no significant evidence for clumping in most cases. Gullberg et al. (2018) looked at [C II] in four $z \sim 4.5$ DSFGs. They found three that showed a smooth morphology, while the fourth could be a clumpy disk, though they cannot rule out the possibility of it being a smooth disk. The data explored by Gullberg et al. (2018) and Hodge et al. (2016) did not have enough signal-to-noise to definitively show that the observed clumps were real, rather than noise fluctuations. The data presented in this work have a higher signal-to-noise ratio, allowing a more confident classification of this system as a merger rather than a clumpy, rotating disk.

## 5.2. Stability of Components

The Toomre $Q$ parameter describes the stability of rotating disk against gravitational collapse. It is given by

$$Q = c_s \kappa \frac{1}{\pi G \Sigma_{\text{gas}}}, \qquad (4)$$

where $c_s$ is the sound speed, $\kappa$ is the epicyclic frequency, and $\Sigma_{\text{gas}}$ is the gas surface density (Toomre 1964). In a system

dominated by turbulent pressure, rather than thermal pressure, this becomes

$$Q = \sqrt{c_s^2 + \sigma_T^2} \kappa \frac{1}{\pi G \Sigma_{\text{gas}}}, \qquad (5)$$

where $\sigma_T$ is the turbulent velocity dispersion. In the limit of high turbulence, $\sqrt{c_s^2 + \sigma_T^2} \sim \sigma_T$ (Hayward & Hopkins 2017). The gas is stable against gravitational collapse if $Q > 1$ and unstable for $Q < 1$, though observations of galaxies and simulations of thick disks place this threshold at $Q \sim 0.7$ (Kennicutt 1989; Kim & Ostriker 2007).

The sound speed and turbulent line width are difficult to measure directly, so the gas velocity dispersion, $\sigma_r$, is often used instead. In cosmological simulations from the FIRE (Feedback In Realistic Environments) suite, Su et al. (2017) found that stellar feedback, which would be an important factor in a rapidly star-forming system like SPT0346-52, increases the turbulent velocity dispersion by a factor of 2–3. Using the velocity dispersion instead of the sound speed and true turbulent velocity dispersion likely provides an upper limit to $Q$ (Prieto & Escala 2016). Because the observed velocity dispersion ($\sigma_r$) can include ordered motion such as rotation or outflows, it tends to overestimate $\sigma_T$. (Su et al. 2017).

The epicyclic frequency, $\kappa$, is $a\Omega$, where $\Omega$ is the rotational frequency and $a$ is a constant. For a flat rotation curve, $a = 1$. In general, $1 < a < 2$. Swinbank et al. (2015) and Oteo et al. (2016) used an intermediary value of $a = \sqrt{3}$; we use the same substitution here. The rotational frequency can be described as $\Omega = v_r/r$. Then, $\kappa = a\Omega \approx \sqrt{3} v_r/r$.

With the above substitutions, we calculate the Toomre $Q$ stability parameter using

$$Q \approx \sigma_r \frac{\sqrt{3} v_r}{r} \frac{1}{\pi G \Sigma_{\text{gas}}}. \qquad (6)$$

While we do not assume that the components of SPT0346-52 are disks, past spatially resolved calculations of $Q$ have found $Q < 1$ locally where there are star-forming regions and giant molecular clouds in other systems, even when the global disk has $Q > 1$ (i.e., Fisher et al. 2017; Genzel et al. 2011; Martig et al. 2009). The $Q$ parameter can therefore be used to find local instabilities independent of the global stability/instability of a system.

To calculate the gas surface density, $\Sigma_{\text{gas}}$, we assume the [C II] emission traces the gas. The total gas mass, $M_g = 1.5 \times 10^{11} M_{\odot}$, taken from Spilker et al. (2015), is divided among the pixels according to their [C II] luminosity. To convert to surface density, the gas mass in each pixel is divided by the area of the pixel.

The surface density is then given by

$$\Sigma_{\text{gas},i} \approx M_g \frac{S_i}{\sum_i S_i A_i}, \qquad (7)$$

where $S_i$ is the integral of the Gaussian component in each pixel for each component from Section 4.2, $\sum_i S_i$ is the total [C II] flux density, and $A_i$ is the area of a pixel.

The value of $\sigma_r$ used in Equation (6) is the standard deviation determined in the Gaussian line fitting described in Section 4.2.
<figure><img src="image_13.png" /><figcaption>**Figure 11.** Maps of the Toomre *Q* disk stability parameter. Left: blue component. The blue component pixels have [C II] intensity-weighted mean $\bar{Q} = 0.03$ and maximum $Q_{max} = 0.13$. Right: red component. The red component pixels have $\bar{Q} = 0.02$ and $Q_{max} = 0.06$. $Q \ll 1$ throughout the galaxy indicates that it is unstable against gravitational collapse, consistent with rapid star formation throughout both disks.</figcaption></figure>

To calculate $v_r$, we first created velocity fields for the two spatial/velocity components using the mean velocity determined by the two-Gaussian line-fitting described in Section 4.2. These velocity fields were then fit using the 2D tilted ring modeling in 3D-Barolo (Di Teodoro & Fraternali 2015). These model velocity fields are used as the values of $v_r$ throughout both components. The position, *r*, is defined relative to the center of each component, indicated by red and blue crosses in Figure 10.

A map of the Toomre *Q* stability parameter is shown in Figure 11. The individual pixels in the blue component have a [C II] intensity-weighted mean of $\bar{Q} = 0.03$ and a maximum value of $Q_{max} = 0.13$. The individual pixels in the red component have $\bar{Q} = 0.02$ and $Q_{max} = 0.06$. All values of *Q* are far less than one, indicating that the the system (separated into individual components) is unstable to gravitational collapse. As mentioned above, using $\sigma_r$ instead of $\sigma_T$ or $c_s$ gives the upper limit for *Q*. Thus, the result that $Q \ll 1$ everywhere and the disks are gravitationally unstable, does not depend on this substitution.

We also calculated values of *Q* that would be measured for SPT0346-52 if we could not spatially resolve its structure, as is typical of high-redshift galaxies observed to date. We also consider the unresolved estimates for the red and blue components alone. These values, along with the maximum rotational velocity, $V_{max}$, the mean velocity dispersion, $\bar{\sigma}$, and the radius of each component, *R*, are given in Table 3. The values of *Q* are low compared to previous studies of DSFGs. For example, Oteo et al. (2016) calculated $Q \sim 0.22$ and $Q \sim 0.35$ for an interacting pair of DSFGs at $z = 4.425$, and Swinbank et al. (2015) found $Q \sim 0.3$ in SDP.81. De Breuck et al. (2014) found a higher average in ALESS 73.1, a $z = 4.76$ DSFG, with average $Q = 0.58$, though with $Q < 1$ at all radii. In Arp220, $Q < 1$ only in the inner part of the disk, where the most intense star formation is occurring (Scoville et al. 1997). The values of *Q* in SPT0346-52 are consistent with studies of star-forming galaxies, where giant star-forming clumps and local overdensities were found to be unstable against fragmentation (Martig et al. 2009; Genzel et al. 2011; Westmoquette et al. 2012). Where $Q > 1$ in Seyfert galaxies and ULIRGs, the disks cannot fragment and form stars (Tacconi et al. 1999; Sani et al. 2012). The low values of *Q* throughout SPT0346-52 indicate that the components are very unstable against collapse, which is fully consistent with the observed high star formation rate.

<caption>**Table 3**</caption>

Kinematic Values of the Components of SPT0346-52

<table><thead><tr><th>Source</th><th>V<sub>c</sub><br>(km s<sup>-1</sup>)</th><th>V<sub>max</sub><br>(km s<sup>-1</sup>)</th><th>σ&#773;<br>(km s<sup>-1</sup>)</th><th>R<br>(kpc)</th><th>Q</th></tr></thead><tbody><tr><td>Both</td><td>...</td><td>263</td><td>307</td><td>1.6</td><td>0.30</td></tr><tr><td>Blue</td><td>-309</td><td>95</td><td>301</td><td>0.95</td><td>0.11</td></tr><tr><td>Red</td><td>158</td><td>30</td><td>313</td><td>0.94</td><td>0.05</td></tr></tbody></table>

### 5.2.1. The Future of SPT0346-52

While mergers can trigger the onset of an AGN (e.g., Wang et al. 2013b), SPT0346-52 has negligible AGN activity (Ma et al. 2016). However, many DSFGs and merging systems do have AGNs (e.g., Younger et al. 2009, 2008; Engel et al. 2011; Westmoquette et al. 2012; Carniani et al. 2013; Rawle et al. 2014). At $z = 3.351$, Marsan et al. (2015) found an AGN in an ultra-massive and compact galaxy at $z = 3.35$ whose stars formed in an intense starburst 300–500 Myr prior. It is possible SPT0346-52 currently has an AGN that is so heavily obscured that X-ray emission is not visible. SPT0346-52 may also host an AGN in the future.

DSFGs are thought to evolve to form the red sequence by $z = 2$. The stars in this red sequence would form in an intense, short, dissipative burst of star formation at $z > 4$ within a compact, $r_e \approx 1$ kpc, region (Kriek et al. 2008). This effective radius is similar to that of SPT0346-52. The models by Narayanan et al. (2015) suggest that by $z \sim 0$ DSFGs (like SPT0346-52) will reside in massive dark matter (DM) halos with $M_{DM} \approx 10^{14} M_{\odot}$, though not all of the intense star formation is driven by major mergers. These studies are in agreement with that of Cattaneo et al. (2013), who found that most $2 < z < 4$ star-forming galaxies with $M_* > 10^{10} M_{\odot}$ evolve into the most massive galaxies on the red sequence and had a phase of intense star formation at $z > 2$. Similarly, angular clustering analyses of $z > 2$ blank-field DSFGs have suggested that DSFGs evolve into present-day halos with masses of $10^{13}-10^{14} M_{\odot}$ (e.g., Chen et al. 2016; Wilkinson et al. 2017). Oteo et al. (2016) observed a pair of interacting DSFGs at $z = 4.425$. The system observed by Oteo et al. (2016) is at an earlier merger stage than SPT0346-52. They concluded that this system is likely the progenitor of a massive, red, elliptical galaxy. At $z = 2.3$, Fu et al. (2013) studied two interacting massive starburst galaxies, separated by 19 kpc and connected by a tidal tail or bridge. They similarly concluded that this system will deplete its gas reservoir in 200 Myr and merge to form an elliptical galaxy with $M_* \sim 4 \times 10^{11} M_{\odot}$.

SPT0346-52 is currently undergoing a phase of intense star formation. It may deplete its gas reservoir in $\sim 30$ Myr (Fu et al. 2013; Narayanan et al. 2015; Spilker et al. 2015; Aravena et al. 2016) and evolve into a red sequence galaxy.

## 6. Summary and Conclusions

In this paper, we presented a pixellated lensing reconstruction of high-resolution [C II] emission observed with ALMA toward the $z = 5.7$ dusty star-forming galaxy SPT0346-52. With this reconstruction, we mapped the integrated [C II] emission and dust continuum at rest-frame 158 µm in the (unlensed) source plane. We spatially resolved the $L_{[C II]}/L_{FIR}$ ratio in SPT0346-52 and showed that the $L_{[C II]}/L_{FIR}$ versus $\Sigma_{FIR}$ relation continues at smaller spatial scales.
We also obtained source-plane velocity information on SPT0346-52, including a demagnified spectrum and moment maps. The reconstruction revealed two spatially and kinematically separated components, one redshifted and one blue-shifted relative to the [C II] rest frequency. These components are connected by a bridge of gas. Each individual component is extremely unstable, with the Toomre Q stability parameter Q ≪ 1 throughout both components.

These components are in the process of merging. This merger is likely driving the intense star formation observed in SPT0346-52. SPT0346-52 may have an AGN in its future and evolve into a massive red sequence galaxy.

K.C.L. and D.P.M. acknowledge support from the U.S. National Science Foundation (NSF) under awards AST-1715213 and AST-1312950 and K.C.L. through award SOSPA4-007 from the National Radio Astronomy Observatory (NRAO). J.V. acknowledges support from NSF award AST-1716127. This material has made use of the El Gato (supported by the NSF award MRI-1228509) high-performance computer. This paper makes use of the following ALMA data: ADS/JAO.ALMA#2013.1.01231.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 SPT is supported by the National Science Foundation through grant PLR-1248097, with partial support through PHY1125897, the Kavli Foundation and the Gordon and Betty Moore Foundation grant GBMF 947. The Flatiron Institute is supported by the Simons Foundation. This research has made use of NASA Astrophysics Data System.

## Appendix
Using Visibility Data to Create the [C II] Spectrum

SPT0346-52 is an extended source with an irregular structure due to gravitational lensing. Therefore, there is no optimal aperture to contain all of the emission. When imaging these data, one has to make assumptions about the structure of the source. Different weightings of the visibilities emphasize different aspects of the galaxy's structure (i.e., faint emission or small structures) and suppress some of the information inherently available from the visibilities. In contrast, the observed complex visibilities contain all of the spectral line information.

We therefore obtain a spectrum of the observed [C II] emission from the observed complex visibilities. The flux density in a given channel, $F_v$, is determined by

$$F_v = \frac{\sum_i \frac{\tilde{v}_{v,i}}{\tilde{m}_i} |\tilde{m}_i|^2}{\sum_i |\tilde{m}_i|^2}, \quad (8)$$

where $\tilde{v}_{v,i}$ is the complex line data visibility and $\tilde{m}_i$ is the complex model visibility for the integrated [C II] line. Dividing data by the model in the numerator removes the spatial structure, transforming the observed visibilities to a point source. The data visibilities are then weighted by the amplitude of the model visibilities in the sum, which emphasizes the visibilities where line emission is expected, and minimizes the weight of the visibilities with little-to-no sensitivity to the

emission structure. We use a model of the complex visibilities for the integrated [C II] line created using a pixellated gravitational lensing reconstruction described further in Section 3.1.

The error on the flux density in each channel, $\sigma$, is the quadrature sum of the weights such that

$$\sigma^2 = \left( \sum_i |\tilde{m}_i|^2 \right)^{-1}. \quad (9)$$

The resulting spectrum is shown in Figure 2.

## ORCID iDs

Katrina C. Litke https://orcid.org/0000-0002-4208-3532
Daniel P. Marrone https://orcid.org/0000-0002-2367-1080
Justin S. Spilker https://orcid.org/0000-0003-3256-5615
Manuel Aravena https://orcid.org/0000-0002-6290-3198
Chian-Chou Chen https://orcid.org/0000-0002-3805-0789
Carlos de Breuck https://orcid.org/0000-0002-6637-3315
Chenxing Dong https://orcid.org/0000-0002-5823-0349
Anthony Gonzalez https://orcid.org/0000-0002-0933-8601
Thomas R. Greve https://orcid.org/0000-0002-2554-1837
Christopher C. Hayward https://orcid.org/0000-0003-4073-3236
Yashar Hezaveh https://orcid.org/0000-0002-8669-5733
Jingzhe Ma https://orcid.org/0000-0003-4178-0800
Warren Morningstar https://orcid.org/0000-0002-5153-0920
Desika Narayanan https://orcid.org/0000-0002-7064-4309
Joaquin Vieira https://orcid.org/0000-0001-7192-3871
Axel Weiss https://orcid.org/0000-0003-4678-3939

## References

Aravena, M., Spilker, J. S., Bethermin, M., et al. 2016, MNRAS, 457, 4406
Bournaud, F., Duc, P.-A., Amram, P., Combes, F., & Gach, J.-L. 2004, A&A, 425, 813
Brauher, J. R., Dale, D. A., & Helou, G. 2008, ApJS, 178, 280
Bussmann, R. S., Riechers, D., Fialkov, A., et al. 2015, ApJ, 812, 43
Carlstrom, J. E., Ade, P. A. R., Aird, K. A., et al. 2011, PASP, 123, 568
Carniani, S., Marconi, A., Biggs, A., et al. 2013, A&A, 559, A29
Casey, C. M., Narayanan, D., & Cooray, A. 2014, PhR, 541, 45
Cattaneo, A., Woo, J., Dekel, A., & Faber, S. M. 2013, MNRAS, 430, 686
Chen, C.-C., Smail, I., Swinbank, A. M., et al. 2015, ApJ, 799, 194
Chen, C.-C., Smail, I., Swinbank, A. M., et al. 2016, ApJ, 831, 91
Dannerbauer, H., Lehnert, M. D., Emonts, B., et al. 2017, A&A, 608, A48
De Breuck, C., Williams, R. J., Swinbank, M., et al. 2014, A&A, 565, A59
Decarli, R., Walter, F., Venemans, B. P., et al. 2017, Natur, 545, 457
Decarli, R., Walter, F., Venemans, B. P., et al. 2018, ApJ, 854, 97
Di Teodoro, E. M., & Fraternali, F. 2015, MNRAS, 451, 3021
Díaz-Santos, T., Armus, L., Charmandaris, V., et al. 2013, ApJ, 774, 68
Díaz-Santos, T., Armus, L., Charmandaris, V., et al. 2017, ApJ, 846, 32
Díaz-Santos, T., Assef, R. J., Blain, A. W., et al. 2016, ApJL, 816, L6
Dong, C., Spilker, J. S., Gonzalez, A. H., et al. 2018, ApJ, submitted
Engel, H., Davies, R. I., Genzel, R., et al. 2011, ApJ, 729, 58
Engel, H., Tacconi, L. J., Davies, R. I., et al. 2010, ApJ, 724, 233
Farrah, D., Lebouteiller, V., Spoon, H. W. W., et al. 2013, ApJ, 776, 38
Fisher, D. B., Glazebrook, K., Abraham, R. G., et al. 2017, ApJL, 839, L5
Fu, H., Cooray, A., Feruglio, C., et al. 2013, Natur, 498, 338
Gabor, J. M., & Davé, R. 2012, MNRAS, 427, 1816
Genzel, R., Newman, S., Jones, T., et al. 2011, ApJ, 733, 101
Goicoechea, J. R., Teyssier, D., Etxaluze, M., et al. 2015, ApJ, 812, 75
Gordon, S., Koribalski, B., & Jones, K. 2001, MNRAS, 326, 578
Greve, T. R., Vieira, J. D., Weiß, A., et al. 2012, ApJ, 756, 101
Gullberg, B., De Breuck, C., Vieira, J. D., et al. 2015, MNRAS, 449, 2883
Gullberg, B., Swinbank, A. M., Smail, I., et al. 2018, ApJ, 859, 12
Haan, S., Surace, J. A., Armus, L., et al. 2011, AJ, 141, 100
Hartley, W. G., Almaini, O., Mortlock, A., et al. 2013, MNRAS, 431, 3045
Hashimoto, T., Inoue, A. K., Mawatari, K., et al. 2018, arXiv:1806.00486
Hayward, C. C., Behroozi, P. S., Somerville, R. S., et al. 2013, MNRAS, 434, 2572

Hayward, C. C., & Hopkins, P. F. 2017, MNRAS, 465, 1682

Hayward, C. C., Jonsson, P., Kereš, D., et al. 2012, MNRAS, 424, 951

Helou, G., Khan, I. R., Malek, L., & Boehmer, L. 1988, ApJS, 68, 151

Herrera-Camus, R., Sturm, E., Graciá-Carpio, J., et al. 2018a, ApJ, 861, 94

Herrera-Camus, R., Sturm, E., Graciá-Carpio, J., et al. 2018b, ApJ, 861, 95

Hezaveh, Y. D., Dalal, N., Marrone, D. P., et al. 2016, ApJ, 823, 37

Hezaveh, Y. D., Marrone, D. P., Fassnacht, C. D., et al. 2013, ApJ, 767, 132

Hodge, J. A., Carilli, C. L., Walter, F., et al. 2012, ApJ, 760, 11

Hodge, J. A., Swinbank, A. M., Simpson, J. M., et al. 2016, ApJ, 833, 103

Hollenbach, D. J., Takahashi, T., & Tielens, A. G. G. M. 1991, ApJ, 377, 192

Hung, C.-L., Hayward, C. C., Smith, H. A., et al. 2016, ApJ, 816, 99

Iono, D., Yun, M. S., Aretxaga, I., et al. 2016, ApJL, 829, L10

Iono, D., Yun, M. S., Elvis, M., et al. 2006, ApJL, 645, L97

Izumi, T., Onoue, M., Shirakata, H., et al. 2018, PASJ, 70, 36

Karl, S. J., Naab, T., Johansson, P. H., et al. 2010, ApJL, 715, L88

Kennicutt, R. C., Jr. 1989, ApJ, 344, 685

Kim, W.-T., & Ostriker, E. C. 2007, ApJ, 660, 1232

Kodama, T., Tanaka, I., Kajisawa, M., et al. 2007, MNRAS, 377, 1717

Kriek, M., van der Wel, A., van Dokkum, P. G., Franx, M., & Illingworth, G. D. 2008, ApJ, 682, 896

Lagache, G., Cousin, M., & Chatzikos, M. 2018, A&A, 609, A130

Luhman, M. L., Satyapal, S., Fischer, J., et al. 1998, ApJL, 504, L11

Luhman, M. L., Satyapal, S., Fischer, J., et al. 2003, ApJ, 594, 758

Lutz, D., Berta, S., Contursi, A., et al. 2016, A&A, 591, A136

Ma, J., Gonzalez, A. H., Spilker, J. S., et al. 2015, ApJ, 812, 88

Ma, J., Gonzalez, A. H., Vieira, J. D., et al. 2016, ApJ, 832, 114

Maiolino, R., Cox, P., Caselli, P., et al. 2005, A&A, 440, L51

Malhotra, S., Helou, G., Stacey, G., et al. 1997, ApJL, 491, L27

Marrone, D. P., Spilker, J. S., Hayward, C. C., et al. 2018, Natur, 553, 51

Marsan, Z. C., Marchesini, D., Brammer, G. B., et al. 2015, ApJ, 801, 133

Martig, M., Bournaud, F., Teyssier, R., & Dekel, A. 2009, ApJ, 707, 250

Mazzucchelli, C., Bañados, E., Venemans, B. P., et al. 2017, ApJ, 849, 91

McMullin, J. P., Waters, B., Schiebel, D., Young, W., & Golap, K. 2007, in ASP Conf. Ser. 376, Astronomical Data Analysis Software and Systems XVI, ed. R. A. Shaw, F. Hill, & D. J. Bell (San Francisco, CA: ASP), 127

Mirabel, I. F., Vigroux, L., Charmandaris, V., et al. 1998, A&A, 333, L1

Mocanu, L. M., Crawford, T. M., Vieira, J. D., et al. 2013, ApJ, 779, 61

Muñoz, J. A., & Oh, S. P. 2016, MNRAS, 463, 2085

Murphy, E. J., Condon, J. J., Schinnerer, E., et al. 2011, ApJ, 737, 67

Narayanan, D., Hayward, C. C., Cox, T. J., et al. 2010, MNRAS, 401, 1613

Narayanan, D., & Krumholz, M. R. 2017, MNRAS, 467, 50

Narayanan, D., Turk, M., Feldmann, R., et al. 2015, Natur, 525, 496

Neri, R., Downes, D., Cox, P., & Walter, F. 2014, A&A, 562, A35

Oteo, I., Ivison, R. J., Dunne, L., et al. 2016, ApJ, 827, 34

Pavesi, R., Riechers, D. A., Sharon, C. E., et al. 2018, ApJ, 861, 43

Pineda, J. L., Velusamy, T., Langer, W. D., et al. 2010, A&A, 521, L19

Planck Collaboration, Ade, P. A. R., Aghanim, N., et al. 2016, A&A, 594, A13

Prieto, J., & Escala, A. 2016, MNRAS, 460, 4018

Rangwala, N., Maloney, P. R., Wilson, C. D., et al. 2015, ApJ, 806, 17

Rawle, T. D., Egami, E., Bussmann, R. S., et al. 2014, ApJ, 783, 59

Riechers, D. A. 2013, ApJL, 765, L31

Riechers, D. A., Carilli, C. L., Capak, P. L., et al. 2014, ApJ, 796, 84

Sanders, D. B., Soifer, B. T., Elias, J. H., et al. 1988, ApJ, 325, 74

Sani, E., Davies, R. I., Sternberg, A., et al. 2012, MNRAS, 424, 1963

Sargsyan, L., Lebouteiller, V., Weedman, D., et al. 2012, ApJ, 755, 171

Scoville, N. Z., Yun, M. S., & Bryant, P. M. 1997, ApJ, 484, 702

Smith, J. D. T., Croxall, K., Draine, B., et al. 2017, ApJ, 834, 5

Spilker, J. S., Aravena, M., Marrone, D. P., et al. 2015, ApJ, 811, 124

Spilker, J. S., Marrone, D. P., Aravena, M., et al. 2016, ApJ, 826, L12

Strandet, M. L., Weiss, A., De Breuck, C., et al. 2017, ApJL, 842, L15

Strandet, M. L., Weiss, A., Vieira, J. D., et al. 2016, ApJ, 822, 80

Su, K.-Y., Hopkins, P. F., Hayward, C. C., et al. 2017, MNRAS, 471, 144

Suyu, S. H., Marshall, P. J., Hobson, M. P., & Blandford, R. D. 2006, MNRAS, 371, 983

Swinbank, A. M., Dye, S., Nightingale, J. W., et al. 2015, ApJL, 806, L17

Tacconi, L. J., Genzel, R., Tecza, M., et al. 1999, ApJ, 524, 732

Tadaki, K., Iono, D., Yun, M. S., et al. 2018, Natur, 560, 613

Teyssier, R., Chapon, D., & Bournaud, F. 2010, ApJL, 720, L149

Thomas, D., Maraston, C., Bender, R., & Mendes de Oliveira, C. 2005, ApJ, 621, 673

Thomas, D., Maraston, C., Schawinski, K., Sarzi, M., & Silk, J. 2010, MNRAS, 404, 1775

Toft, S., Smolčič, V., Magnelli, B., et al. 2014, ApJ, 782, 68

Toomre, A. 1964, ApJ, 139, 1217

Vieira, J. D., Crawford, T. M., Switzer, E. R., et al. 2010, ApJ, 719, 763

Vieira, J. D., Marrone, D. P., Chapman, S. C., et al. 2013, Natur, 495, 344

Wagg, J., Carilli, C. L., Wilner, D. J., et al. 2010, A&A, 519, L1

Walter, F., Riechers, D., Cox, P., et al. 2009, Natur, 457, 699

Wang, R., Wagg, J., Carilli, C. L., et al. 2013a, ApJ, 773, 44

Wang, S. X., Brandt, W. N., Luo, B., et al. 2013b, ApJ, 778, 179

Warren, S. J., & Dye, S. 2003, ApJ, 590, 673

Weiß, A., De Breuck, C., Marrone, D. P., et al. 2013, ApJ, 767, 88

Westmoquette, M. S., Clements, D. L., Bendo, G. J., & Khan, S. A. 2012, MNRAS, 424, 416

Wilkinson, A., Almaini, O., Chen, C.-C., et al. 2017, MNRAS, 464, 1380

Wright, E. L. 2006, PASP, 118, 1711

Younger, J. D., Dunlop, J. S., Peck, A. B., et al. 2008, MNRAS, 387, 707

Younger, J. D., Hayward, C. C., Narayanan, D., et al. 2009, MNRAS, 396, L66

Yun, M. S., Aretxaga, I., Gurwell, M. A., et al. 2015, MNRAS, 454, 3485

Zirm, A. W., Stanford, S. A., Postman, M., et al. 2008, ApJ, 680, 224