ALMA reveals extended cool gas and hot ionized outflows in a typical star-forming galaxy at $z=7.13$

We present spatially-resolved morphological properties of [CII] 158 $\mu$m, [OIII] 88 $\mu$m, dust, and rest-frame ultraviolet (UV) continuum emission for A1689-zD1, a strongly lensed, sub-L* galaxy at $z=7.13$, by utilizing deep Atacama Large Millimeter/submillimeter Array (ALMA) and Hubble Space Telescope (HST) observations. While the [OIII] line and UV continuum are compact, the [CII] line is extended up to a radius of $r \sim 12$ kpc. Using multi-band rest-frame far-infrared (FIR) continuum data ranging from 52-400 $\mu$m, we find an average dust temperature and emissivity index of $T_{\rm dust} = 41^{+17}_{-14}$ K and $\beta = 1.7^{+1.1}_{-0.7}$, respectively, across the galaxy. We find slight differences in the dust continuum profiles at different wavelengths, which may indicate that the dust temperature decreases with distance. We map the star-formation rate (SFR) via IR and UV luminosities and determine a total SFR of $37\pm 1~M_\odot~{\rm yr}^{-1}$ with an obscured fraction of $87\%$. While the [OIII] line is a good tracer of the SFR, the [CII] line shows deviation from the local $L_{\rm [CII]}$-SFR relations in the outskirts of the galaxy. Finally, we observe a clear difference in the line profile between [CII] and [OIII], with significant residuals ($\sim 5\sigma$) in the [OIII] line spectrum after subtracting a single Gaussian model. This suggests a possible origin of the extended [CII] structure from the cooling of hot ionized outflows. The extended [CII] and high-velocity [OIII] emission may both contribute in part to the high $L_{\rm [OIII]}$/$L_{\rm [CII]}$ ratios recently reported in $z>6$ galaxies.


INTRODUCTION
The formation and evolution of galaxies across cosmic time is shaped by the "baryon cycle," namely the complex interplay between gas accretion, star formation, stellar and active galactic nucleus (AGN) feedback, and gas outflows (e.g. Davé et al. 2012;Anglés-Alcázar et al. 2017;Tumlinson et al. 2017). These various processes serve to transfer gas between the interstellar medium (ISM) and the surrounding circumgalactic medium (CGM). Key to understanding the role of the baryon cycle in the formation of the first galaxies are studies investigating the morphology and structure of the ISM and CGM at high redshift, where robust observational constrains are rare.
At high redshift, far-infrared (FIR) fine-structure lines are typically used to trace the ISM and the galaxies' star-formation activity. In particular, the [C ii] 2 P 3/2 → 2 P 1/2 transition at 157.74 µm (hereafter [C ii]), with its low ionization potential of 11.2 eV (compared to 13.6 eV of Hydrogen) is a dominant coolant of the neutral ISM (Stacey et al. 1991;Wolfire et al. 2003), and dense photodissociation regions (PDRs; Hollenbach & Tielens 1999) associated with giant molecular clouds. The [C ii] line, arising from singly ionized carbon, has been used to probe ISM properties (e.g. Capak et al. 2015;Pentericci et al. 2016;Knudsen et al. 2016;Carniani et al. 2017;Bakx et al. 2020;Matthee et al. 2020) and kinematics (e.g. Jones et al. 2017;Smit et al. 2018;Rizzo et al. 2020Rizzo et al. , 2021 in dozens of high redshift galaxies. While the [C ii] line traces multiple phases of the ISM, it has been shown to primarily trace cool, neutral gas through PDRs (Cormier et al. 2019). By contrast, the [O iii] 3 P 1 → 3 P 0 transition at 88.4 µm (hereafter [O iii]) primarily traces hot, diffuse, ionized gas, as the ionization potential of O + (35.1 eV) is significantly higher than that of Hydrogen. Moreover, both the [C ii] and [O iii] lines are known to trace the star-formation rate (SFR; Kapala et al. 2014;De Looze et al. 2014;Herrera-Camus et al. 2015).
The advent of Atacama Large Millimeter/submillimeter Array (ALMA) opens our observational window to detect both the [C ii] and [O iii] lines at z > 4, and the number of detections risen rapidly in recent years (see Hodge & da Cunha 2020). Of particular relevance, recent ALMA stacking studies have discovered extended [C ii] structures up to r ∼ 10-15 kpc around high redshift galaxies (Fujimoto et al. 2019;Ginolfi et al. 2020a). The extended structure, called the "[C ii] halo" (Fujimoto et al. 2019), shows a ∼ 5 times larger effective radius than the central galaxy as seen in the UV and a potential association with the Lyα halo (e.g. Leclercq et al. 2017). Similar extended [C ii] structures have also been identified in observations of individual star-forming galaxies at z > 4 (Fujimoto et al. 2020;Herrera-Camus et al. 2021) and massive quasar host galaxies at z > 6 ( Maiolino et al. 2012;Cicone et al. 2015;cf. Meyer et al. 2022).
In several cases, observations of the extended [C ii] structure have additionally detected the "broad-wing" feature in the [C ii] spectrum Ginolfi et al. 2020a;Herrera-Camus et al. 2021), characteristic of ongoing high-velocity outflows, which have been hypothesized to drive carbon-enriched gas out into the CGM. The existence of extended [C ii] emission is not fully reproduced by the current galaxy formation models (Fujimoto et al. 2019;Arata et al. 2020;Pallottini et al. 2019;Katz et al. 2019aKatz et al. ,b, 2022 suggesting that there may be missing physics in simulations at the Epoch of Reionization (Pizzati et al. 2020). Fujimoto et al. (2019) discuss several possible scenarios for what physical processes may power the extended [C ii] emission, including outflows, satellite galaxies (in which [C ii] emission traces their SFR), and photoionization via radiation from the central galaxy. Ginolfi et al. (2020a) and Fujimoto et al. (2020) find correlations between galaxy SFR and the detection of the broad-wing feature and [C ii] extension, suggestive of the SF-driven outflow scenario. However, this result may be biased by the signal to noise, as objects with higher SFR are likely brighter, making faint components in the spectra easier to detect (Ginolfi et al. 2020a;Fujimoto et al. 2020). Ginolfi et al. (2020b) introduce an additional scenario, tidal disruption, which can produce [C ii] emission on CGM scales through turbulence and shocks through gravitational interactions. Moreover, while the existence of the extended [C ii] structure implies the presence of dust in the CGM, high redshift ALMA observations have generally found a compact dust continuum (Fujimoto et al. 2019;Ginolfi et al. 2020a). In contrast, observations of local galaxies suggests that dust entrained by galactic winds can reach ∼ 10 kpc in scale (Kaneda et al. 2009;Meléndez et al. 2015;McCormick et al. 2018;Yoon et al. 2021) and theory suggests that stellar feedback is efficient at driving dust-enriched, multiphase winds (Kannan et al. 2021). ALMA observations of quasar host galaxies at z 6 have found [C ii] and dust emission extended out to r ∼ 10 kpc, and interpret both as tracing the star-forming ISM of these massive galaxies (Novak et al. 2020). A proper determination of the physical processes producing the extended [C ii] structure in line with the outflow thus depends on our ability to measure, with sufficient depth, the distribution of neutral/ionized gas and dust continuum emission. Simultaneous ALMA observations of the [C ii] and [O iii] fine-structure lines, and their underlying continua, are well-suited to accomplish this task.
In this work, we utilize deep ALMA observations of the galaxy A1689-zD1, a strongly lensed, sub-L* galaxy at z = 7.13 (Bradley et al. 2008;Watson et al. 2015;Knudsen et al. 2017). The high magnification (µ ∼ 9.3) of this galaxy makes it an ideal target for a spatiallyresolved investigation of the ISM in the abundant population of sub-L* galaxies in the Epoch of Reionization, and indeed, recent high-resolution follow-up observations have successfully detected bright [C ii] 158 µm and [O iii] 88 µm emission lines (Knudsen et al. in prep). Moreover, rich FIR continuum data has also been taken in bands 3, 6, 7, 8, and 9, providing valuable constraints on the FIR luminosity and dust temperature (Bakx et al. 2021). The availability of deep, high resolution, multiwavelength ALMA data make A1689-zD1 the ideal target to address the aforementioned issues and examine the extended [C ii] emission around high redshift galaxies.
The paper is structured as follows: In Section 2 we describe the data reduction steps for both ALMA and HST data. We report the results for the spatial extent of the [C ii], [O iii], UV, and FIR emission in Sections 3.1 and 3.2. In Section 3.3 we spatially-resolve the dust temperature, which we use in Section 3.4 to derive the spatially-resolved L [C ii] -SFR and L [O iii] -SFR relations for A1689-zD1. We present the line spectra of [C ii] and [O iii] in Section 3.5. Finally, in Section 4 we discuss our results, including the interpretation of the extended emission as a "[C ii] halo" vs. an extended disk and the physical origin of the emission. Throughout this work, we assume a ΛCDM cosmology with H 0 = 67.4 km s −1 Mpc −1 , Ω M = 0.315, Ω Λ = 0.685 (Planck Collaboration et al. 2020) and a Kroupa (2001) initial mass function (IMF).

ALMA Data
We utilize ALMA observations of A1689-zD1 (R.A. 13 h 11 m 29.9 s , Dec. −1 • 19 18.7 ) in bands 3, 6, 7, 8, and 9. The observations in bands 6 and 8 (ALMA program IDs 2015.1.01406.S and 2017.1.00775.S, respectively) will be presented in Knudsen et al. (in prep). We use archival observations for bands 7 (Knudsen et al. 2017), and 9 (Bakx et al. 2021), and refer the reader to these papers for more details. While previous studies of A1689-zD1 have adopted the optical redshift of z = 7.5 determined from an apparent continuum break (interpreted as the Lyman-α; Watson et al. 2015 52 µm line will be presented in a separate, upcoming paper. Band 3 observations were tuned to cover the CO(7-6) and CO(6-5) lines.
The Common Astronomy Software Applications package (CASA; McMullin et al. 2007) was used for reduction, calibration, and imaging. The data were first processed via the standard ALMA pipeline reduction, which was sufficient for these observations. All images were produced with the tclean task in CASA using natural weighting to maximize sensitivity. Cleaning was done using the "auto-multithresh" automasking algorithm (Kepley et al. 2020) initially down to the 3σ level but expanded to a lownoisethreshold of 1σ. In order to image both compact and extended emission, we apply the multiscale deconvolver (Cornwell 2008) with scales of 0 (delta-function) 1, and 3 times the beam size (i.e., the psf). Continuum maps were produced using the multi-frequency synthesis (mfs; Conway et al. 1990) mode in CASA over a frequency range identified as lacking line emission. To produce moment 0 maps for the [C ii] and [O iii] lines, we first subtract the continuum using the CASA task uvcontsub, with order 0, and produce a map using the mfs mode. While this mode is typically used for continuum imaging, we employ it to produce an average line intensity map, as a pseudomoment 0 map, in order to maximize the signal-to-noise of the map prior to cleaning and ensure that we detect and deconvolve any faint and/or diffuse signals in the tclean process. We integrate channels between 233.4 and 234.0 GHz for the [C ii] map and between 416.6 and 417.8 GHz for the [O iii] map. Both ranges correspond to roughly ±400 km/s from the line center.
For [C ii] and [O iii] line intensity maps we obtain an rms of 18.8 µJy beam −1 and 63.0 µJy beam −1 . For bands 3, 6, 7, 8, and 9 continuum maps we obtain an rms of 10.4, 5.3, 44.9, 24.5, and 181.2 µJy beam −1 . All maps are primary-beam corrected, though the reported rms is derived prior to primary-beam correction. All ALMA and HST images are aligned astrometrically using the reproject package in Python.
With natural weighting, the synthesized beam sizes for band 6 and 8 observations are 0.30 × 0.28 and 0.46 ×0.41 , respectively. To perform a fair comparison based on the same spatial resolution, we apply a uvtaper of 0.3 for the band 6 data, which increases the beam size to 0.44 × 0.40. For the analysis of extended emission (Section 3.2) we apply a uv-taper of 0.5 (0.7 ) to observations in band 8 (band 6). Table 1 summarizes Note-Columns: (1) ALMA band, (2-3) Central wavelength and range, (4-5) Beam size, for "fiducial" maps (with no uv-taper) and "tapered" maps (with a uv-taper of 0.5 ), (6-7) RMS noise, (8) Reference for ALMA data. † For band 6, the "fiducial" and "tapered" maps have uv-tapers of 0.3 and 0.7 , respectively, to more closely match the other bands in resolution. ‡ Band 6 and band 8 data are also presented in Wong et al. (2022).
the beam size and depth of both the fiducial (uv-taper of 0.3 for band 6, untapered for other bands) and tapered (0.7 for band 6, 0.5 for other bands) continuum data in each band.

HST Data
In addition to deep ALMA observations of FIR line and continuum emission, we use archival Hubble Space Telescope (HST) data to probe the rest-frame UV continuum of A1689-zD1 (Watson et al. 2015). We obtain images taken with Hubble WFC3/IR in the J (F125W) and H (F160W) bands. At z = 7.13, these bands trace the UV continuum at 1350Å< λ < 2050Å. These wavelengths are redder than the Lyman limit and unaffected by the Lyα forest, making this a good probe of the restframe UV continuum and the unobscured SFR.
Recent works (e.g. Rujopakarn et al. 2016;Dunlop et al. 2017;Fujimoto et al. 2019;Tamura et al. 2019;Herrera-Camus et al. 2021) have shown that there sometimes exist offsets of 0.1-0.3 arcseconds between HST and ALMA astrometry, which can significantly impact results of comparisons between high-resolution observations. The HST data we use in this work has been corrected for this astrometric offset by comparison to astrometry from the Gaia survey (Gaia Collaboration 2018), which is accurate to ALMA astrometry to the milliarcsecond scale. The correction shifts the HST astrometry −0.11 in right ascension and +0.06 in declination.
As noted in previous analysis of A1689-zD1 (Watson et al. 2015;Knudsen et al. 2017), there is an additional source ∼ 1.5 westwards of the target that has been identified as a low-redshift (z ∼ 2) galaxy. Other sources nearby A1689-zD1 can be identified as low-z interlop-ers galaxies by their detection in shorter-wavelength ACS/F814W observations, where A1689-zD1 disappears due to the Lyman-break nature of the galaxy. We regard pixels that have S/N > 3 in F125W+F160W and are more than 1.1 from A1689-zD1, that are also detected in ACS/F814W, as interlopers and replace them with randomly drawn values from the image rms noise. Finally, in order to present a fair comparison with the same spatial resolution between HST and ALMA images of A1689-zD1, we convolve the HST map with a Gaussian kernel constructed to match the idealized (i.e., Gaussian model) ALMA beam for the band 8 data. For the fiducial (tapered) maps, this kernel has a size of 0.44 ×0.39 (0.70 × 0.67 ).

Source Plane Reconstruction
To analyze the intrinsic ISM morphology of A1689-zD1, we reconstruct the ALMA and HST images into the source plane, correcting for the lens magnification and deflection.
The magnification of A1689-zD1 is derived to be a factor of µ = 9.3 ± 0.5 based on the mass model originally published in Limousin et al. (2007), and refined to include spectroscopic information of multiple systems as described in Bina et al. (2016). This parametric strong lensing model is constructed using the lenstool ) software which provides a set of models sampling the posterior distribution of each parameter of the mass distribution, allowing us to derive the statistical error on the magnification. The uncertainty of 0.025 dex of the magnification at the location of the source was also confirmed independently as described in Watson et al. (2015), and the magnification is quite uniform across the source, with variations around 5% (see  in Appendix A). We show source plane reconstructed maps for each emission in Figure A.2; we employ these source plane maps when measuring distance scales (i.e. Sections 3.2.2 and 3.2.3). This is because the asymmetrical deflection of the source confuses a simple relation between angular and physical scale in the image plane-the full reconstruction is necessary to assess the size/shape of the source.
When quoting the luminosity or flux from the source, we employ the notation that, for example, L IR represents the observed IR luminosity and L IR µ −1 represents the intrinsic IR luminosity, corrected for magnification. These images indicate that the [C ii] emission is likely extended relative to the [O iii], FIR, and UV. [C ii] emission is detected at > 2σ out to r ∼ 2-3 from the galaxy center. Furthermore, the dust continuum at 90 µm (band 8) appears more compact than at 163 µm (band 6), which we will revisit in Section 3.3. z

uv-plane Visibility Profiles
We also analyze the spatial structures of the different emission in the uv-visibility plane using the python package uvplot (Tazzari 2017). In the bottom panel of Figure A.2, we show the real part of the visibility amplitude in bins of the uv-distance, which is inversely proportional to the observed angular scale θ.
We apply a single and double Gaussian fits to the [O iii], [C ii], and dust continuum profiles (shown in blue and red, respectively). Based on the reduced chi-squared (χ 2 ν ) statistic, the single Gaussian is a reasonably good model for the [O iii] line and the dust continuum in both bands, with consistent FWHMs of ∼ 1.0 ; on the other hand, for the [C ii] line, we find a significant offset from the best-fit single Gaussian at the shortest uv-distance regime much beyond the error scales. The double Gaussian shows a better fit for the [C ii] line, lowering the χ 2 ν . For the compact and extended [C ii] components, we derive a best-fit FWHM of 1.05 ± 0.03 and 4.0 ± 0.6 , respectively. The FWHM of the compact [C ii] component is almost the same as those of the [O iii] line and the dust continuum, while the extended [C ii] component is ∼ 4 times larger than the compact component. These results suggest that there exists an independent extended component surrounding the central galaxy, which is reminiscent of the Lyα halo universally identified around normal star-forming galaxies at z ∼ 3-6 (e.g., Momose et al. 2014;Leclercq et al. 2017). The size ratio of the extended and compact components is comparable to the previous ALMA deep [C ii] stacking results of ∼ 5-6 estimated from galaxies and quasars at z ∼ 6 (Fujimoto et al. 2019;Novak et al. 2020).
We note that the existence of the extended [C ii] emission, despite no extended component in the [O iii] emission, is likely not just a result of the data depth or maximum recoverable scale. The S/N at the peak pixel in the [C ii] intensity map is 42.9, almost comparable to [O iii] (37.7), only different by ∼ 10%. Moreover, as the histogram of uv-dist in Figure 1 indicates, the band 8 observations cover more short baselines (larger angular scales) than band 6, due to the different antenna configurations (C-4 and C-3 for bands 6 and 8, respectively).
Despite the larger error bars, we note that the profiles of the dust continuum in both bands show some excess emission at large angular scale (low uv-dist), similar to the [C ii] line. Although the lower S/N of this excess trend prevents us from giving a definitive conclusion for the existence of the extended dust component, we examine the potential extended dust structure also in the image-based analysis (Section 3.2.2).

Image-based Radial Profiles
We derive radial surface brightness profiles, based on images for both emission lines and the UV and FIR continuum. We conduct this analysis in the source plane, as described in Section 2.3. For the radial profiles, we use the midpoint between the two UV-bright peaks as the center. Figure 2 shows radial profiles in the source plane for [C ii], [O iii], and UV continuum (top) and for the dust continuum at 90 µm and 163 µm (bottom), normalized at their peaks. Points show the mean surface brightness in circular annuli whose width increases at increasing distances, from 0.4 kpc to 1 kpc. Error bars show the standard error of the mean. For ease of readability, we show these radial profiles out to the first point at which  158 µm, band 8 dust continuum, band 6 dust continuum, and the rest-frame UV continuum in A1689-zD1. In each panel, white dashed contours show the −3σ and −2σ levels, while colored contours show the 2σ and 3σ levels and then increase following the Fibonacci sequence. A slight uv-taper of 0.3 is applied to both band 6 maps to match the resolution of the band 8 maps, and the ellipses on the bottom left indicate the beam size. The rest-frame UV continuum is measured by stacking the HST H and J-bands, masking out nearby sources, and convolving with a Gaussian kernel constructed to match the ALMA beam in band 8. For clarity, we show the unconvolved (high resolution) HST map in the inset panel for the image plane. We pinpoint the two peaks of the UV continuum emission on each panel in the image plane with two blue crosses. All panels show the same 5.6 × 5.6 area on the sky. Bottom: Visibility profiles of each ALMA observation. The amplitude is extracted from the real part of the complex visibility in bins of uv-distance, which corresponds inversely to the observed angular scale θ. We fit each visibility profile to a single Gaussian, and find that this is a good fit for the [O iii] line and dust continuum. For [C ii], however, a double Gaussian model provides a better fit, indicating the presence of significant extended emission. The bottom right panel shows histograms of uv-dist for the two ALMA bands. The detection of the extended [C ii] structure is not due to an increased sensitivity to large-scale structure (i.e., small uv-dist) in band 6. the mean surface brightness drops below 0. At this point and beyond, we consider the profile dominated by noise fluctuation and instead show only 1σ upper limits.
Both [O iii] and UV continuum emission are compact (though extended relative to the PSF) dropping to below 0 at radii of ∼ 4.5 and ∼ 3 kpc, respectively. In contrast, [C ii] emission is detected out to r ∼ 12 kpc. Furthermore, the [C ii] profile shape is inconsistent with the [O iii] profile out to ∼ 6 kpc, and inconsistent with the UV profile out to ∼ 8 kpc.
Dust continuum emission appears to be extended with a similar radial profile as the [C ii] line. This is true at both 163 µm and 90 µm, though the latter wavelength drops to 0 at ∼ 8 kpc while the former extends out to ∼ 12 kpc. 1 Both continuum profiles are inconsistent with UV and [O iii] at radii of r ∼ 2-4 kpc, more in line with the [C ii] profile. As with the visibility-based profiles, the large error bars at larger scales prevents robust evidence of extended emission. Nevertheless, the consistent positive surface brightness in the dust continuum may suggest some real extended emission. The presence of such extended dust continuum emission would be in contrast to previous ALMA results for star-forming galaxies at z ∼ 4-7, and we discuss this further in Section 4. While radial profiles give a quantitative sense as to the extent of [C ii] emission, they may also be sensitive to the "clumpiness" of the emission. Recent analysis has shown that A1689-zD1 is a highly clumpy object, requiring multiple spatial/spectral components to fit its [C ii] and [O iii] datacubes (Wong et al. 2022, Knudsen et al. in prep). As such, the precise choice of center may impact the shape of the surface brightness profile. In the imaging analysis presented thus far, we adopt the midpoint between the two UV-bright peaks as the cen- ter; however, in order to test against this possible bias, we additionally compute the Gini coefficient for each emission. The Gini coefficient (G) is a metric borrowed from econometrics, and was originally created to summarize the degree of income inequality in a population. It was first used to assess galaxy morphology by Abraham et al. (2003), who showed that it can, to first order, be interpreted similarly to other morphological indicators such as concentration. The advantage of using G in this work is that it is non-parametric: it does not depend on the galaxy having a well-defined center, but rather only indicates whether the bulk of the total intensity is contained within a few pixels or spread evenly across many. A high value of G (∼ 0.75) indicates that the emission from the object is concentrated in one or more bright nuclei, while a low value of G (∼ 0.25) indicates that the emission is more uniform. We compute the Gini coefficient following Lotz et al. (2004) by first sorting all n pixel values X i into increasing order and then computing

Gini Coefficients
We compute the error on G in a Monte-Carlo fashion by varying the underlying map by a randomized "noise  . The far-infrared SED of A1689-zD1. We fit a modified blackbody spectrum to ALMA photometry from bands 9, 8, 7, and 6, derived from 2.8 × 1.8 elliptical apertures placed over the continuum maps shown in the right panels. We additionally constrain the SED fitting according to the 3σ upper-limit from the non-detection of continuum in band 3 (rest-frame ∼ 400 µm). The grey shaded region shows the 68% range of best-fit SEDs determined from varying the flux values according to their uncertainties in a Monte-Carlo fashion. The inset panel plots the 1σ contour of the probability distribution of T dust and β values derived from the Monte-Carlo method, and highlights the degeneracy between the two parameters. We find a best-fit dust temperature of T dust = 41 +17 −14 K and dust emissivity index of β = 1.7 +1.1 −0.7 , where the error represents the 68% confidence interval. The images on the right show the four continuum detections, with contours denoting the 2, 3, 5, 8, 12, 16, 20, 25 × σ levels and dashed yellow lines denoting the aperture.
map," convolved with the ALMA beam. Since the Gini coefficient it sensitive to the number of pixels included, we compute G based on all pixels within a circular aperture of radius r = 2 for each map. Moreover, while the Gini coefficient can be sensitive to the S/N of the map, this effect is negligible above a per-pixel S/N of ∼ 10 (Lisker 2008), and the sensitivity to the S/N is accounted for by the Monte-Carlo error computation. Figure 3 shows the Gini coefficient computed in the source plane for each map. Our Gini coefficient results reaffirm the conclusions from the radial profiles: [O iii] and UV emission are compact and concentrated in 1 or 2 bright nuclei, while [C ii] and dust continuum emission are spread more evenly. Interestingly, the Gini coefficient for the dust continuum suggests a difference between the different wavelengths: the dust continuum at shorter wavelengths is more concentrated than at higher wavelengths. This may suggest that the dust temperature may vary in space, with hotter dust (which peaks at shorter wavelengths) more concentrated than colder dust (which peaks at higher wavelengths). The spatially-resolved dust temperature distribution is addressed further in Section 3.3.

Spatially-Resolved Dust Temperature
The rich FIR continuum data for A1689-zD1 facilitates accurate determinations of the dust temperature. Importantly, band 9 (rest-frame 53 µm) data is taken, which helps significantly to constrain the FIR SED at shorter wavelengths. However, the low S/N of the band 9 observations prohibits using this lower wavelength data to constrain the dust temperature in a spatiallyresolved sense. We opt instead to use the ratio between bands 6 (163 µm) and 8 (90 µm) to spatially-resolve the dust temperature; however, this method depends on the assumed dust emissivity index, β, which usually ranges from 1-2.5 (Casey 2012). Therefore, we first determine the best-fit β value via a modified blackbody fit on the galaxy as a whole, using continuum data from bands 3, 6, 7, 8, and 9. We then derive the spatially-resolved dust temperature map, assuming that β remains constant in space. Figure 4 shows the FIR SED of A1689-zD1. Measurements of the dust continuum flux from bands 6 (163 µm), 7 (107 µm), 8 (90 µm), and 9 (53 µm) are shown, and a 3σ upper limit is shown for the continuum in band 3 (∼ 400 µm), which is not detected. All photometry is computed in the 2.8 ×1.8 elliptical aperture shown over the dust continuum images on the right. This aperture is adopted consistently across all bands and is determined by eye to correspond roughly to the [O iii]-emitting region. The aperture is large enough to cover the bulk of the dust continuum in bands 6-8, but not so large as to include outskirt noise fluctuation in band 9 and significantly reduce the S/N in that band. Since band 9 holds most of the constraining power in our MBB fit, it is critical to focus on this central region where the band 9 continuum is robustly detected. Error bars in Figure 4 show the combined uncertainty, including both statistical uncertainty measured from the image rms, and an additional systematic uncertainty of 10% (20% in band 9) arising from flux calibration, as specified in the ALMA Cycle 8 proposer's guide. 2 We fit a modified blackbody spectrum to the data points, correcting for the decreasing temperature contrast between the target and background due to the heating of dust by the CMB following da Cunha et al. (2013). We perform this fitting using a Monte-Carlo approach, allowing T dust and β to vary freely with flat priors. We determine a best-fit dust temperature of T dust = 41 +17 −14 K and dust emissivity index β = 1.7 +1.1 −0.7 . The inset panel shows the distribution of T dust and β values derived from the MCMC fitting, and highlights the inherent degeneracy between the two parameters. The larger uncertainty in our T dust and β measurements compared to Bakx et al. (2021) is due to our use of a 20% uncertainty on the band 9 measurement.
To extend the analysis in a spatially-resolved fashion, we then fix β = 1.7 as a fiducial value. We note that, with variations in the dust chemical composition and crystallinity, it would be possible for β to vary in T CMB (z) T dust = 41 K, χ 2 ν = 3.40 Figure 6. Radial profile of the dust temperature in A1689-zD1. We derive the dust temperature from the ratio of the continuum at 90 µm to the continuum at 163 µm, correcting for CMB effects. As in Section 3.2, we use uv-tapered maps.
Error bars show 1σ uncertainty on the inferred dust temperature derived in a Monte-Carlo fashion from the uncertainty on the continuum ratio. We show the furthest bin as an open diamond, as the S/N is too low to provide a stringent constraint on the dust temperature. Despite this uncertainty, the dust temperature appears to decrease to larger radii, approaching the CMB temperature.
space with T dust fixed. However, we have fit the full FIR SED in two spatially-resolved apertures and found < 1σ difference in β between the two. Therefore, we consider it reasonable to fix β and explore the resulting variation in T dust . We do so by deriving a relationship between the observed rest-frame 90 µm to 163 µm continuum ratio (corrected for CMB effects) and the dust temperature. Figure 5 shows the spatially resolved dust temperature map for pixels where the dust continuum is detected (at 2σ) at both 163 µm and 90 µm. We show in the right panel the fractional uncertainty on the inferred dust temperature. While we see three high dust temperature regions, the two on the far northeast and northwest sides of the galaxy are likely a product of noise fluctuation, as the uncertainty is > 60%. We overplot contours showing the stacked HST J +Hband image. The southeast UV-bright region corresponds to a peak in the dust temperature of ∼ 50 K, while the northwest region corresponds primarily to a noisy peak. However, the high dust temperature in the galaxy center stretches out towards both noisy peaks, with low errors. Similarly, the far southwest of the map shows a higher temperature than the minimum, with low errors. These regions likely correspond to distinct clumps in the galaxy ISM, which have been identified in the galaxy spectrum (Wong et al. 2022, Knudsen et al. in prep).
In addition to estimating T dust on a pixel-by-pixel basis, we also endeavor to estimate the spatially-resolved T dust on kpc scales. Figure 6 shows the inferred dust temperature as a function of distance, in arcsec, based on image-plane radial profiles of the continuum maps at both wavelengths. We use circular annuli of variable width out to 2 from the galaxy center, which corresponds to ∼ 7 kpc when accounting for the gravitational lensing. We use image-plane maps for this analysis to minimize any additional uncertainty brought about by the source-plane reconstruction. Error bars show 1σ uncertainty on the dust temperature. The radial profile suggests that the dust temperature decreases with increasing distance from the galaxy center, starting at ∼ 50 K within 0.25" and dropping to ∼ 30 K beyond 1". At further distances, it appears to increase again, though our signal-to-noise becomes too low to robustly measure the dust temperature.
We explore the feasibility of a T dust gradient by applying a constant model in which the dust temperature is fixed at the best-fit spatially-integrated value (T dust = 41 K) everywhere. This model results in a reduced chi-squared of χ 2 ν = 3.4. In contrast, a simple linear model for the decreasing trend agrees within the error in every annulus. These results indicate that a model with constant T dust is unlikely to be consistent with the data. While the errors in our results are large, this analysis suggests that deep ALMA observations, with the aid of gravitational lensing, may be able to observe the T dust gradient at z > 7.

The Σ [C ii] -Σ SFR and Σ [O iii] -Σ SFR Relations
Spatially-resolving the dust temperature facilitates robust estimates of the spatially-resolved total IR luminosity (8-1000 µm) and obscured SFR. For each pixel, we compute the total IR luminosity by integrating a modified blackbody with β = 1.7 and T dust set to the inferred dust temperature shown in Figure 5, or, if the temperature was not well determined, 3 with T dust = 30 K. We then compute the obscured SFR following the where L UV is the UV luminosity measured from the stacked HST J+H band image. This probes the restframe wavelength range 1350Å< λ < 2050Å, comparable to the GALEX FUV transmission curve used in the calibration of Murphy et al. (2011). The calibrations from Murphy et al. (2011) assume a Kroupa (2001) IMF and a constant star-formation history with an age of ∼ 100 Myr. Uncertainties on the derived SFRs include the uncertainties on the ALMA/HST maps as well as the uncertainty in the dust temperature (for SFR IR ). Summing over the elliptical aperture used for photometry in Section 3.3, and assuming β = 1.7, we derive an obscured SFR of SFR IR = 32 ± 1 M yr −1 and an unobscured SFR of SFR UV = 4.95 ± 0.03 M yr −1 , resulting in a total SFR of ∼ 37 ± 1 M yr −1 and an obscured fraction of ∼ 87%. The high obscured SFR is inline with recent estimates for A1689-zD1 from Bakx et al. (2021). We derive a total infrared luminosity of L IR µ −1 = (2.16 ± 0.07) × 10 11 L .
To compare the derived SFR to [C ii] and [O iii] emission, we place circular annuli of varying width across the tapered maps, and compute the SFR and line surface densities in each annulus. Figure 7 shows the spatiallyresolved Σ [O iii] -Σ SFR and Σ [C ii] -Σ SFR relations with the empirical, local universe relation from De Looze et al. (2014) and (for Σ [C ii] ) the results for z ∼ 2-6 dusty star-forming galaxies from Spilker et al. (2016). We see that the [O iii] line roughly traces the SFR following the relation for local, low-Z, dwarf galaxies, with a slight (∼ 0.3 dex) excess in regions with high Σ SFR (i.e., the galaxy center). By contrast, we see a ∼ 0.7 dex deficit of [C ii] in the galaxy center, relative to local dwarfs. This deficit is consistent with the [C ii]-FIR deficit in high redshift DSFGs observed by Spilker et al. (2016), which has also been observed in local IR luminous galaxies (Díaz-Santos et al. 2014). The deficit at high Σ SFR diminishes at lower Σ SFR , with the two largest annuli showing [C ii] detections despite non-detections of UV or IR continua. The [O iii] line is also not detected in these annuli.

Central [C ii] and [O iii] Line Profiles
We also investigate the [C ii] and [O iii] line profiles in a spatially-resolved manner. In these emission line spectra, high-velocity outflows in different gas phases may be identified via the "broad-wings" that extend past the extent of the best-fit Gaussian profile (see e.g. Veilleux et al. 2020). We produce line profiles from the [C ii] and [O iii] cubes in a small (0.25 radius) aperture centered at the the peak of the Σ SFR map, i.e. over the UV+FIR-bright region. This small radius is roughly the beam size of the fiducial ALMA maps. We focus   on this region for two reasons: 1) the clumpy nature of A1689-zD1 confuses identification of the broad-wing feature with distinct spectral components, and 2) in a SF-driven outflow scenario, the faint broad-wing feature would be most easily detected in regions of highest Σ SFR (Davies et al. 2019). We have verified that this aperture is small enough to probe only one of the five kinematic component components identified in Knudsen et al. (in prep Katz et al. 2019bKatz et al. , 2022, as the [O iii] tends to trace younger starforming regions which tend to have a lower velocity dispersion than the system as a whole. The larger line width for [C ii] may also be related to the spatial extent of the [C ii] emission. If the extended emission were spherically symmetric, then the larger volume of gas observed may capture more complicated kinematics along the line of sight and thus contribute to the larger line width. To test for the outflow signature, we fit a Gaussian to each line with the FWHM fixed at the measured value. We see that the [C ii] line is well fit by a single Gaussian, with a reduced chi-squared of χ 2 ν = 1.18. By contrast, the [O iii] line is not well fit by a single Gaussian, with a reduce chi-squared of χ 2 ν = 2.16 (p < 10 −7 ). This owes to significant excess flux in high-velocity regions, despite the narrow line width overall. We show in the bottom panels of Figure 8 the residuals on each Gaussian fit. For [O iii], we find significant (> 5σ) residuals when integrating across the high-velocity tails from ±100-400 km s −1 . As [C ii] is well-fit by a single Gaussian, we find no significant residuals in these high-velocity tails.

DISCUSSION
We have shown in Section 3 that the [C ii] emission (and possibly also the dust continuum emission) in A1689-zD1 is extended on ∼ 12 kpc scales relative to [O iii] and UV emission. The [C ii] line emission in uv space is well-fit by a double Gaussian model, with a extended/compact component ratio of ∼ 4. This [C ii] emission deviates slightly from the local Σ [C ii] -Σ SFR relation, with a "deficit" of [C ii] in the galaxy center and a slight excess in the furthest apertures. We additionally detect significant large-velocity residuals in the [O iii] line spectrum at the galactic core, indicating the possible presence of outflow dominated by hot, ionized gas. We now discuss in depth the physical interpretation of these observations; in particular, 1) the description of the extended [C ii] emission as a part of the central galaxy versus a "[C ii] halo," separate from a single galaxy disk 4 We determine the FWHM via linear interpolation between the data above and below 50%. We determine the error on the FWHM in a Monte-Carlo fashion, by varying each data point in the line spectrum according to the uncertainty. and reminiscent of the Lyα halo, and 2) the physical origin of the extended emission.

[C ii] halo vs. Extended disk
The results presented thus far are reminiscent of the "Lyman-α halo" universally identified around normal star-forming galaxies at z ∼ 3-6 (e.g., Momose et al. 2014;Leclercq et al. 2017). However, another possibility is that we are observing the outskirts of the central galaxy disk even out to r ∼ 10 kpc, and that starformation in this disk is the source of the extended [C ii] emission (e.g., Novak et al. 2020). Due to the detection of extended FIR emission (Fig. 2), we cannot immediately rule out the latter scenario.
While the rest-UV and [O iii] sizes are similarly compact, we see potential evidence of extended dust emission in both uv -and image-based analyses (e.g., Figures 1 and 2). The presence of extended dust continuum emission would be in contrast to previous ALMA results for star-forming galaxies at z ∼ 4-7 (Fujimoto et al. 2019;Ginolfi et al. 2020a;Herrera-Camus et al. 2021), and would, at first glance, suggest significant obscured star-formation in the outskirts of the galaxy.
One possible explanation for this inconsistency with previous observations is that this extended dust emission is common in the early universe, but too cold and faint to be detected, being embedded in the warm CMB (e.g. da Cunha et al. 2013;Vallini et al. 2015;Zhang et al. 2016;Pallottini et al. 2017;Lagache et al. 2018). The dust continuum profiles for A1689-zD1 are similar to the UV and [O iii] line up to r ∼ 3 kpc, and only diverge at larger distances. This faint extended structure may have been missed in previous stacking results, especially if the observed dust temperature gradient is common, as it is in local galaxies (Hunt et al. 2015) and as is predicted in cosmological zoom-in simulations for early galaxies (e.g., Behrens et al. 2018;Sommovigo et al. 2020). For this study, the deep integration times in bands 6 and 8, combined with the strong lensing effect (with magnification µ ∼ 10) yields an effective sensitivity ∼ 10 times deeper than previous stacking results (e.g. Fujimoto et al. 2019;Ginolfi et al. 2020a). Under this interpretation, many high-z star-forming galaxies may host extended dust continuum emission corresponding to their extended [C ii] emission. In fact, recent stacking results for z 6 massive quasar host galaxies find similar radial profiles for [C ii] and dust (Novak et al. 2020), although the physical mechanisms producing extended dust emission might be different between massive quasar host galaxies and sub-L* galaxies.
It is worth noting that the detection of extended dust continuum emission implies a larger dust mass than what has previously been estimated (e.g. M dust = 1.7 +1.3 −0.7 × 10 7 M from Bakx et al. 2021). Indeed, the elliptical aperture applied in Section 3.3 fails to capture ∼ 13% of the band 6 continuum flux observed out to r ∼ 1.5 . Based on our spatially-resolved T dust estimates, fitting this excess flux to a modified blackbody with T dust = 30 K and β = 1.7 yields a dust mass of M dust ∼ 7 × 10 6 M . This suggests that the dust mass estimate could be increased by ∼ 40% compared with the previous estimate of Bakx et al. (2021) by including the extended dust component. If the extended dust component ubiquitously exists around early galaxies, this additional, previously unaccounted for dust mass may have implications for dust evolution in the early universe (e.g. Mancini et al. 2015;Micha lowski et al. 2019;Dayal et al. 2022).
The extended [C ii] and dust continuum emission could be explained by localized star-formation, i.e., an extended SF disk rather than an independent halo component. Given that the dust temperature at the outskirts of the galaxy approaches the CMB temperature ( Figure 6), it is possible that we underestimate Σ IR in these regions despite correcting for the decreasing contrast following da Cunha et al. (2013). However, even in this scenario, the compact rest-UV size of the galaxy implies that the obscured fraction increases up to ∼ 100% towards the outskirts of the galaxy. We show in Figure 9 the UV and IR SFRs (top) and obscured fraction (bottom) as a function of radius in several circular annuli (as in Figure 7). While the obscured fraction decreases with increasing radius up to r ∼ 1 , the detections of IR SFR despite non-detections of UV SFR yield 3σ lower limits on the obscured fraction of ∼ 93%. Such a reversal of the decreasing trend is in contrast to the metallicity gradient expected from simulations, which would continue to decrease to ∼ 1% at the outskirts .
Therefore, a realistic picture might include some dust heating mechanisms other than localized SF activity in the outskirts of the galaxy, such as heating from the interstellar radiation field (ISRF) or by hot outflows. In this scenario, we would need to correct the SFR estimates by removing the extended FIR emission counted as SFR. Then, the Σ [C ii] -Σ SFR relation in the furthest annuli would deviate even more from the local relation, strengthening the interpretation that extended [C ii] emission is powered by other heating mechanisms than localized SF.
At the same time, the [C ii] emission in the outskirts of the galaxy can also be suppressed by the CMB, depending on the gas temperature and density (Vallini et al. 2015;Kohandel et al. 2019). Since in general the gas  Figure 9. ΣSFR and obscured fraction assuming that the extended FIR continuum emission is solely due to localized SF. Top: Surface density of star formation derived from the UV (blue) and IR (red), as a function of radius in several circular annuli. Bottom: Obscured fraction of star formation (SFRIR/SFRtot) as a function of radius. In both panels, arrows indicate 3σ upper/lower limits. Non-detections of the UV SFR in the outer annuli, where IR emission is detected, constrain the obscured fraction to 93%. Such a high obscured fraction in the outskirts is unexpected, and may suggest that the extended FIR emission is not entirely due to SF. temperature ( 10 2 K) is much higher than T CMB , the detection of the extended [C ii] emission may indicate i) the extended C + gas maintains a relatively high density to be not significantly affected by the CMB suppression, or ii) the majority of the [C ii] emission in the outskirt regions is actually suppressed and the intrinsic [C ii] luminosity is even greater. The latter case would yield even further deviation from the local Σ [C ii] -Σ SFR relation (right panel of Figure 7. As A1689-zD1 is a low-mass (∼ 10 9 M ), sub-L* galaxy (Watson et al. 2015), the existence of an extended disk of r ∼ 10 kpc would be striking. The effective radius for such sub-L* galaxies is generally estimated   Figure 10 shows these fits. For the [O iii] line, we find that including a second component is unnecessary, as the emission is well fit by a single Sérsic profile with R eff = 1.21 ± 0.01 kpc and n = 0.59 ± 0.02. This is generally consistent with the rest-frame optical and UV sizes measured by Shibuya et al. (2015), though perhaps a bit more extended. For the [C ii] line, the single Sérsic model fails to capture the extended emission; instead, we need a central component with R central eff = 1.36±0.04 kpc and n = 0.55 ± 0.03 and an exponential halo component with R halo eff = 3.5 ± 0.3 kpc. In addition to the visibilitybased analysis (Section 3.2.1), the fitting clearly confirms the existence of two separate components, which is reminiscent of the Lyα halo identified around high-z SFGs.
In summary, the extended dust continuum and [C ii] emission might be explained by the outskirt of an extended galaxy disk. But this interpretation suggests an increasing trend of the obscured fraction as a function of radius up to ∼ 100% at large radii. Instead, the extended [C ii] structure is more likely explained by an independent gaseous halo, which is evident from both 5 While the rest-frame UV is not necessarily a perfect tracer of the overall source size, Shibuya et al. (2015) find consistent sizes with rest-frame optical observations out to z ∼ 1. Similarly, Jiménez-Andrade et al. (2021) find consistent sizes between UV, optical, and radio wavelengths at z ∼ 3.
visibility-and image-based profile fitting for the [C ii] emission, and we may call it "[C ii] halo." The [C ii] halo and the extended dust emission would have formed via some physical mechanism(s) other than localized SF, which we will discuss in the next section.

The Physical Origin of the [C ii] Halo
Assuming that the extended [C ii] emission represents a distinct [C ii] halo, we now discuss the physical origin of this feature. Fujimoto et al. (2019) present five possible explanations for [C ii] emission on circumgalactic scales: A) star formation from satellite galaxies; B) circumgalactic-scale photodissociation region (CG-PDR); C) circumgalactic-scale H ii region (CG-H ii); D) cold streams; and E) outflow. For merging galaxies in the center of the protocluster, Ginolfi et al. (2020b) argue a sixth possible explanation: F) turbulence and shocks associated with tidal stripping.
In all cases, some physical process is required to enrich the CGM with carbon and produce the halo of C + that we observe. Assuming that cold streams are composed primarily of pristine gas, only scenarios A), E), and F) include some mechanism to enrich the CGM with carbon. Scenarios B) (CG-PDR) or C) (CG-H ii) may certainly ionize the carbon and power [C ii] emission, provided the CGM was pre-enriched by some other process. However, given that A1689-zD1 is a sub-L*, pre-mature galaxy at z = 7.13 whose formation history is relatively short, we assume that the metal-enrichment process is related to on-going activity and focus our discussion on the three scenarios of A) star-formation, E) outflow, and F) tidal disruption.

Star Formation from Satellite Galaxies
As discussed in Section 4.1, we find that the spatiallyresolved Σ [C ii] -Σ SFR relation shows a slight excess of [C ii] emission relative to the local relation at low SFR (at the outskirts of the galaxy). Such deviation from the local relation suggests that the extended [C ii] emission is unlikely to be explained by individual, lowmass, low-SFR satellite galaxies. Moreover, following the mass-metallicity relation (e.g. Mannucci et al. 2010), such faint, low-mass satellite galaxies generally have low metallicity. At this low-Z, the [C ii] line emissivity for a given stellar continuum decreases (Vallini et al. 2015;Olsen et al. 2017;Katz et al. 2022), which would in fact have the opposite effect as we observe.

Outflows
Another possible scenario for the origin of the [C ii] halo is that high-velocity outflows drive carbon out to r ∼ 10 kpc. These outflows are responsible for the metal enrichment of the CGM, and this interpretation has been supported by observations of ongoing outflows, as traced by the [C ii] line, in high-SFR galaxies Ginolfi et al. 2020a;Herrera-Camus et al. 2021). Since the majority of [C ii] is thought to be radiated from PDRs in the cold, neutral hydrogen gas clouds, the broad wing feature observed in [C ii] suggests that the outflowing gas is dominated by cool, neutral gas.
However, galactic outflows in general are known to be made up of multiple phases of gas, with the archetypal example of M82 showing a hot phase observed in X-ray emission (T 10 6 K), a warm ionized phase (T ∼ 10 4 K), cold and warm molecular phases and a cold atomic phase (Heckman & Thompson 2017). The multiphase nature of galactic outflows suggests that the [C ii] line may not reliably trace outflows in all galaxies. In fact, Spilker et al. (2020) find that molecular outflows are ubiquitous in z > 4 galaxies, but not always identifiable by [C ii] emission. We show in Section 3.5 that the [O iii] line shows significant ∼ 5σ residuals in the high-velocity region after fitting a single Gaussian model. By contrast, we find no significant residuals in [C ii]. This may suggest the presence of ongoing outflows dominated by hot, ionized gas, rather than cool, neutral gas in A1689-zD1.
The theoretical connection between these hot outflows and the observed, cool, [C ii] halo depends on the cooling time of the hot ionized gas (t cool ). Optically thin intergalactic gas at z ∼ 3.5 (with temperature T ∼ 10 6 K and an overdensity δ ∼ 1-10) is predicted to have a cooling time of t cool 4 Gyr, which is longer than the Hubble time at this redshift (Madau et al. 2001). How-ever, this depends on the specific properties of the gas, e.g. the density, metallicity, and multiphase structure of the outflows. In fact, recent "cooling outflow" models implemented into the analytical models of Pizzati et al. (2020) have successfully reproduced the [C ii] halo out to r ∼ 10 kpc. In their models, hot ionized outflows experience catastrophic cooling shortly after their initial launch. Depending on the timescale on which this hot outflowing gas cools, we may or may not expect to observe the broad-wings in [C ii] under this scenario. This might contribute to the diversity of the presence or the absence of the [C ii] broad-wing feature in recent ALMA high-z studies (e.g. Decarli et al. 2018;Gallerani et al. 2018;Fujimoto et al. 2019;Stanley et al. 2019;Ginolfi et al. 2020a;Bischetti et al. 2019;Novak et al. 2020;Herrera-Camus et al. 2021;Spilker et al. 2020;Izumi et al. 2021a,b). It is worth noting that in the case of A1689-zD1 the [C ii] line profile may also be composed of narrow+broad components, as with the [O iii] line, but degenerate with a single Gaussian profile. In this scenario, we are observing multi-phase outflows in an early galaxy including both cold neutral and hot ionized gases.
We also note that outflows may also be capable of forming the extended dust structure potentially observed around A1689-zD1 (Section 3.2). Though the outflow may also contribute to destruction of dust, the existence of dust from ∼ 20 kpc to a few Mpc from galaxies has been confirmed from the reddening measurements of background objects (e.g., Ménard et al. 2010;Peek et al. 2015;Tumlinson et al. 2017). Fujimoto et al. (2020) report that 1 out of 23 isolated ALPINE galaxies shows a dust continuum profile similar to the [C ii] emission, despite compact UV continuum. Interestingly, the object also shows a large Lyα EW even with the large amount of dust, and the authors argue a potential mechanism of outflows dominated by hot ionized gas, which may serve to eject dust out to the CGM, keep the dust hot, and allow for the Lyα line to escape from the system. If the hot-mode outflows last only for a short timescale, this extended dust continuum structure could be somewhat rare in the early universe.

Tidal stripping
As previously noted, A1689-zD1 is a clumpy object with multiple spatial/spectral components necessary to fit the [C ii] and [O iii] datacubes (Wong et al. 2022, Knudsen et al. in prep). This raises the possibility that A1689-zD1 is actually multiple merging systems, and that the [C ii] halo arises from the dissipation of mechanical energy via turbulence and shocks associated with tidal stripping. Such tidal effects would also be expected to eject material from the galaxy, including dust and old stars that would contribute to the ISRF in the CGM. These effects have been explored by Ginolfi et al. (2020b), who observe [C ii] emission from circumgalactic gas structure around a massive (∼ 10 10 M ), merging galaxy system at z = 4.57 in the center of a protocluster environment.
The close projected distances (∼ 0.5 ) of the different components in A1689-zD1 make difficult a robust distinction between multiple, merging systems and a single complex, clumpy object. Moreover, A1689-zD1 is a sub-L*, low-mass (∼ 10 9 M ) galaxy, and the galaxy overdensity is not identified around A1689-zD1. As such, the situation may be different from a system of massive merging galaxies in the center of a high redshift protocluster. Nevertheless, minor mergers have been predicted to drive significant kinematic disruption of gas in low-mass galaxies (see e.g. the simulations of Kohandel et al. 2020), which could power [C ii] emission. While further analysis of the kinematics of the object is needed to explore the merging system scenario (and will be explored in Knudsen et al. in prep), it should be noted that tidal disruption in a merging system may play a role in producing the [C ii] halo around A1689-zD1.

Summary: origin of the [C ii] halo
Our observational results suggest the presence of outflows, dominated by hot, ionized gas, in and around the star-forming core of A1689-zD1. In conjunction with the analytical models of Pizzati et al. (2020), our results imply that the [C ii] halo likely formed via some combination of a) hot, ionized outflows and subsequent catastrophic cooling, or b) past outflow activity including the cool, neutral gas. Our results do not rule out the possibility that the CGM was pre-enriched by past outflow/merger activity and the observed [C ii] halo is powered by strong ionizing photons from the galaxy core or extended H ii regions. In any case, our results support the connection between the [C ii] halo and the outflow history of the galaxy, and highlight the importance of considering the multiphase nature of galactic outflows.
However, it remains a possibility that the central and secondary components trace a kinematically complex bulge and diffuse disk of the galaxy, respectively, where the latter is missed in previous HST observations of the rest-frame UV continuum. Rest-frame optical observations with JWST, which have been scheduled for Cycle 1, will map out the entire stellar distribution and help to answer definitively whether the extended [C ii] emission is part of the galaxy or the CGM. Recent evidence from ALMA has shown that high redshift galaxies tend to have (∼ 3-10 times) higher L [O iii] /L [C ii] ratios than local galaxies (Harikane et al. 2020;Hashimoto et al. 2019;Carniani et al. 2020;Bakx et al. 2020). This may suggest different properties of the ISM at high-z, such as a higher ionization parameter, PDR deficit, higher burstiness parameter (Vallini et al. 2021), or lower C/O abundance ratio than solar (Becker et al. 2011;Katz et al. 2022).
Using the best sensitivity [O iii] and [C ii] maps so far obtained at this redshift, we have shown the existence of the [C ii] halo and [O iii] outflows in A1689-zD1, both of which may be additional factors to explain the high L [O iii] /L [C ii] ratios in z > 6 galaxies. Although our results rely on one galaxy, A1689-zD1 is just one of the abundant population of sub-L* galaxies whose properties may be a representative picture in the Epoch of Reionization. Therefore, we discuss here the contribution from the [C ii] halo and [O iii] outflows to the L [O iii] /L [C ii] ratio in A1689-zD1. A full analysis of the spatially-resolved L [O iii] /L [C ii] ratio in A1689-zD1 will be presented in Knudsen et al. (in prep) 4.3.1. Contribution from [C ii] halo One possible explanation for the high L [O iii] /L [C ii] ratios in z > 6 galaxies is the systematic existence of, and failure to account for, the extended [C ii] emission on CGM scales. This has been investigated by Carniani et al. (2020), who find that correcting for surface brightness dimming in high-resolution, low S/N observations can lower the L [O iii] /L [C ii] ratios by a factor of ∼ 2. Figure 11 shows the L [O iii] /L [C ii] ratio as a function of SFR for A1689-zD1, in multiple circular apertures of different sizes. These different apertures each encompass different amounts of the [C ii] halo. We additionally show the 9 galaxies studied in Harikane et al. (2020). By accounting for the full extent of the [C ii] halo (< 3.5 ), we find a lower L [O iii] /L [C ii] ratio of ∼ 2, compared to ∼ 4 for the central region of the galaxy (< 0.5 ). That is, accounting for extended [C ii] emission can help to reconcile the high L [O iii] /L [C ii] ratios with lower ratios found in local galaxies, but cannot bring these values completely in line with local relations.

Contribution from [O iii] outflows
Another possible explanation for the high L [O iii] /L [C ii] ratio, at least in the case of A1689-zD1, is the contribution of the "broad-wing" feature in the [O iii] spectrum to the total [O iii] luminosity. If this emission indeed corresponds to hot, ionized gas in outflows, rather than H ii regions around massive stars, we would expect to see a higher L [O iii] /L [C ii] ratio than local universe galaxies not undergoing this outflow.
To examine this, we analyze the [O iii] line profile in A1689-zD1. Since the single Gaussian fitting in [O iii] shows significant residuals (Section 3.5), we perform two component (narrow+broad) Gaussian fitting to the [O iii] line profile. Figure 12 shows the [O iii] line profile in the same central aperture as Figure 8, with a two-component Gaussian fit applied. We see that the broad component dominates the overall spectrum, with the narrow component responsible for only about 13% of the total luminosity.
The bottom panel of Figure 12 shows how the L [O iii] /L [C ii] ratio changes when calculated using the total [O iii] luminosity (narrow+broad) vs. the luminosity of only the narrow component. We use the same SFR for both points, derived from the UV+FIR emission as described in Section 3.4. We see that, while this central region shows an overall L   Figure 11, for the total [O iii] luminosity as well as only the narrow component. (Inoue et al. 2016). Importantly, the S/N of the [O iii] line in our data is at least 4-8 times larger than previous studies (e.g. Harikane et al. 2020), which could explain the non-detections of the [O iii] outflow signature.

SUMMARY
In this paper, we present deep ALMA observations of the strongly lensed, z = 7.13 galaxy A1689-zD1, focusing on the [C ii] 158 µm and [O iii] 88 µm fine-structure lines as well as the FIR continuum. Utilizing archival HST data to probe the UV continuum, we compare the morphology and spatial extent of each wavelength of emission to examine the multiphase nature of the ISM. Our major findings are summarized below.  Figure A.1 shows a map of the magnification µ across the field. While the magnification is highest in the southwest corner of the field, and lowest in the northeast, it is quite uniform across the source.
We show in Figure A.2 the 2D maps of each emission, in the image plane (as in Figure 1) as well as reconstructed into the source plane, with a uv-taper of 0.5 (0.7 for band 6) applied. Much of the analysis in this paper is conducted on these uv-tapered maps to optimize the S/N of the extended, diffuse emission. Moreover, we derive radial profiles, gini coefficients, and the dust temperature profile from the uv-tapered source plane maps, as they provide a direct measurement of the radial distance, corrected for the lens deflection.