The TW Hya Rosetta Stone Project IV: A hydrocarbon rich disk atmosphere
L. Ilsedore Cleeves, Ryan A. Loomis, Richard Teague, Edwin A. Bergin, David J. Wilner, Jennifer B. Bergner, Geoffrey A. Blake, Jenny K. Calahan, Paolo Cazzoletti, Ewine F. van Dishoeck, Viviana V. Guzman, Michiel R. Hogerheijde, Jane Huang, Mihkel Kama, Karin I. Oberg, Chunhua Qi, Jeroen Terwisscha van Scheltinga, Catherine Walsh
DDraft version February 22, 2021
Typeset using L A TEX twocolumn style in AASTeX62
The TW Hya Rosetta Stone Project IV: A hydrocarbon rich disk atmosphere
L. Ilsedore Cleeves, Ryan A. Loomis, Richard Teague, Edwin A. Bergin, David J. Wilner, Jennifer B. Bergner, ∗ Geoffrey A. Blake,
7, 8
Jenny K. Calahan, Paolo Cazzoletti, Ewine F. van Dishoeck,
10, 11
Viviana V. Guzm´an, Michiel R. Hogerheijde,
10, 13
Jane Huang,
5, 9, ∗ Mihkel Kama,
14, 15
Karin I. ¨Oberg, Chunhua Qi, Jeroen Terwisscha van Scheltinga,
16, 10 andCatherine Walsh University of Virginia, Charlottesville, VA National Radio Astronomy Observatory, Charlottesville, VA Center for Astrophysics | Harvard & Smithsonian, Cambridge, MA University of Michigan, Ann Arbor, MI Center for Astrophysics | Harvard & Smithsonian, 60 Garden Street, Cambridge, MA 02138, USA University of Chicago, Department of the Geophysical Sciences, Chicago, IL 60637, USA Division of Chemistry & Chemical Engineering, California Institute of Technology, Pasadena CA 91125, USA Division of Geological & Planetary Sciences, California Institute of Technology, Pasadena CA 91125, USA Department of Astronomy, University of Michigan, 1085 South University Avenue, Ann Arbor, MI 48109, USA Leiden Observatory, Leiden University, PO Box 9513, 2300 RA Leiden, The Netherlands Max-Planck-Institut f¨ur Extraterrestrische Physik, Giessenbachstraße 1, D-85748 Garching bei M¨unchen, Germany Instituto de Astrof´ısica, Ponticia Universidad Cat´olica de Chile, Av. Vicu˜na Mackenna 4860, 7820436 Macul, Santiago, Chile Anton Pannekoek Institute for Astronomy, University of Amsterdam, Science Park 904, 1098 XH, Amsterdam, The Netherlands Department of Physics and Astronomy, University College London, Gower Street, London, WC1E 6BT, UK Tartu Observatory, University of Tartu, 61602 T˜oravere, Estonia Laboratory for Astrophysics, Leiden Observatory, Leiden University, PO Box 9513, 2300 RA Leiden, The Netherlands School of Physics and Astronomy, University of Leeds, Leeds LS2 9JT, UK
ABSTRACTConnecting the composition of planet-forming disks with that of gas giant exoplanet atmospheres,in particular through C/O ratios, is one of the key goals of disk chemistry. Small hydrocarbons likeC H and C H have been identified as tracers of C/O, as they form abundantly under high C/Oconditions. We present resolved c − C H observations from the TW Hya Rosetta Stone Project, aprogram designed to map the chemistry of common molecules at 15 −
20 au resolution in the TW Hyadisk. Augmented by archival data, these observations comprise the most extensive multi-line set fordisks of both ortho and para spin isomers spanning a wide range of energies, E u = 29 −
97 K. Wefind the ortho-to-para ratio of c − C H is consistent with 3 throughout extent of the emission, andthe total abundance of both c − C H isomers is (7 . − × − per H atom, or 1 −
10% of thepreviously published C H abundance in the same source. We find c − C H comes from a layer near thesurface that extends no deeper than z/r = 0 .
25. Our observations are consistent with substantial radialvariation in gas-phase C/O in TW Hya, with a sharp increase outside ∼
30 au. Even if we are notdirectly tracing the midplane, if planets accrete from the surface via, e.g., meridonial flows, then such achange should be imprinted on forming planets. Perhaps interestingly, the HR 8799 planetary systemalso shows an increasing gradient in its giant planets’ atmospheric C/O ratios. While these stars arequite different, hydrocarbon rings in disks are common, and therefore our results are consistent withthe young planets of HR 8799 still bearing the imprint of their parent disk’s volatile chemistry.
Keywords:
Protoplanetary disks – Astrochemistry – Exoplanet atmospheric composition
Corresponding author: L. Ilsedore [email protected] ∗ NHFP Sagan Fellow a r X i v : . [ a s t r o - ph . E P ] F e b Cleeves et al. INTRODUCTIONWe are entering an era where measurements of thecompositions of giant exoplanet atmospheres are be-coming increasingly common. A diversity of chemicalproperties (via carbon-to-oxygen ratios, or C/O) hasbeen seen (e.g., Madhusudhan et al. 2011; Madhusudhan2012; Kreidberg et al. 2015), including within a singleplanetary system (HR8799; Bonnefoy et al. 2016; Lavieet al. 2017; Lacour et al. 2019; Molli`ere et al. 2020). Tounderstand the origins of this diversity, we must studyplanets’ formation environments: gas rich protoplane-tary disks. The chemical properties of these disks aredriven by at least two factors, i. the make-up of themolecular cloud out of which the star and disk formed(e.g., Visser et al. 2009, 2011; Drozdovskaya et al. 2019),and ii. the disk physical properties (irradiation level,temperature, density, etc.) that can drive an activelyevolving chemistry prior to and during planet formation(e.g., Cleeves et al. 2014).The relative contribution of these factors, i.e., therole of inheritance versus later chemical reprocessingin the disk itself, changes with both radial and verti-cal location. For example, near the highly irradiateddisk surface and/or close to the central star, the chem-istry is effectively “scrambled” leaving little memory ofthe molecular composition of the cloud. Near the mid-plane, especially in the outer disk (beyond ∼
10 au),the high extinction levels provide a safer haven for someof the material originating in the molecular cloud to bepreserved, with further processing requiring moderateto high external irradiation or cosmic ray fluxes (e.g.,Bergin et al. 2014; Cleeves et al. 2014; Yu et al. 2017;Eistrup et al. 2018); however, it is unclear whether disksare sufficiently ionized to facilitate an active midplanechemistry (Cleeves et al. 2015).There is a growing body of evidence that the observ-able chemistry in planet forming disks, at least thatwithin the “warm molecular layer” (Aikawa et al. 2002),deviates from “typical” molecular cloud chemistry. Forexample, observations of CO emission and CO isotopo-logues are faint when compared to expectations based ondust masses from millimeter emission, a gas-to-dust cor-rection factor, and an interstellar CO abundance (e.g.,Favre et al. 2015; Schwarz et al. 2016; Ansdell et al.2016). Water was surprisingly challenging to detect inprotoplanetary disks with
Herschel , and when detected,fluxes were more than an order of magnitude lower thananticipated based on astrochemical modeling with UVirradiated water ice at interstellar abundances (Berginet al. 2010; Hogerheijde et al. 2011; Du et al. 2017).The abundant interstellar complex organic methanol,first detected in disks by the Atacama Large Millimeter Array (ALMA), was similarly quite faint (Walsh et al.2014), yet observations of CH CN appeared relativelybright, with nearly a 1:1 abundance ratio inferred be-tween them for the TW Hya disk (Loomis et al. 2018).Therefore we need better constraints on the chemi-cal composition of disk gas, especially at an elementallevel, to understand the chemical reservoir from whichplanets may accrete and what, e.g., C/O or N/O ratiothey might inherit at least initially. Chemical modelshave demonstrated that the abundances of simple hy-drocarbons like C H and c − C H are very sensitive tothe C/O ratio of the gas (e.g., Bergin et al. 2016; Kamaet al. 2016; Cleeves et al. 2018; Miotello et al. 2019;Fedele & Favre 2020). Pure freeze-out of solids from gascan vary C/O in the gas or ice from ∼ . Hwas extremely sensitive to the bulk C/O in volatiles un-til C/O > ∼ .
5, which overlaps with ob-served exoplanet atmospheric values. As a result, diskhydrocarbon studies open an exciting potential avenueto connect the composition of the gas in disks to thatmeasured in planets’ atmospheres. In addition, observa-tions of hydrocarbon emission in disks have been foundto be generally quite bright at sub-mm wavelengths (Qiet al. 2013a; Kastner et al. 2014, 2015; Bergin et al. 2016;Cleeves et al. 2018; Bergner et al. 2019), enabling smallsurveys ( (cid:46)
14 sources) of C H emission in disks (Guil-loteau et al. 2016; Miotello et al. 2019; Bergner et al.2019).The next step is localizing the distribution of hydro-carbons in the disk, to understand the range of pos-sible C/O values a planet could inherit from a singledisk environment. The first resolved image of C H wasmade of the TW Hya disk by Kastner et al. (2015) withthe Submillimeter Array, where it was found to have aring-like morphology. Its disk-averaged physical naturewas constrained using multiple transitions and its hy-perfine structure. However, due to the face-on natureof this disk 5 − ◦ (Qi et al. 2006, 2008; Huang et al.2018), degenerate excitation solutions were found to fitthe data, either a cold, dense solution or alternatively arelatively warmer, but low density solution. The formerwould suggest an enhanced (greater than solar) C/O ra-tio fairly deep into the disk gas, near the planet-formingregion, while the latter would suggest that the C/O en-hancement primarily is closer to the surface and perhapsless connected to the planet forming midplane. − C H in TW Hya H’s ring-like geometry, as well asobservations of c − C H toward TW Hya. The spatialdistribution of C H and c − C H were found to matchidentically in radial distribution, suggestive that this ra-dial region of the TW Hya disk supports a generally richhydrocarbon chemistry, consistent with an elevated C/Oratio in this region.The present paper uses multi-line observations of c − C H to better spatially constrain the nature of thehydrocarbon layer in the TW Hya protoplanetary disk,with the goal of improving our interpretation of the C/Oratio(s) of this disk. The observations were conducted aspart of the TW Hya as a Chemical Rosetta Stone Project(PI: Cleeves), and have been augmented with ALMAarchival data. The set of lines covers energies spanning29.1 K to 96.5 K, and with the relatively highly criticaldensity and range of opacities probed, these lines arewell suited to constraining the nature, and crucially thelocation, of small hydrocarbon chemistry in TW Hya. OBSERVATIONS2.1.
ALMA Observations
The observations of c − C H utilized in the presentstudy were carried out with ALMA as part of threedifferent observational programs. We present new ob-servations of c − C H taken as part of the TW Hyaas a Chemical Rosetta Stone Program (PI: Cleeves,2016.1.00311.S), augmented by archival observationsfrom 2013.1.00198.S (PI: Bergin) and 2013.1.00114.S(PI: ¨Oberg). The complete set of c − C H transitionsfrom these programs used in our analysis are listed inTable 1. The observations from 2016.1.00311.S werecarried out with 40 antennas on April 8, 2017 (C40-2; 15m – 390m baselines) and May 21, 2017 (C40-5;15m – 1124m baselines). The April 8 observation usedJ1037-2934 for bandpass and phase calibration, andJ1058+0133 for flux calibration. The May 21 obser-vation used J1037-2934 for bandpass, flux, and phasecalibration.All data were calibrated with the ALMA pipeline withCASA version 4.7.0. Prior to imaging, we phase self-calibrated the data using the line-free portions of thecontinuum, adopting a solution interval of 30s and av-eraging polarizations. In addition, spws were self cali-brated independently, and we set a minimum signal tonoise ratio of 3 and a minimum number of baselines perantenna of 6. The calibration process for the observa-tions from program 2013.1.00198.S and 2013.1.00114.Sare reported in Bergin et al. (2016) and ¨Oberg et al.(2017) and not repeated here. The projects 2013.1.00198.S and 2016.1.00311.S hadoverlapping coverage of the bright 10 – 9 and 9 – 8blends (see Table 1). The integrated line flux mea-sured between these two programs differed by about 10%(well within the quoted ALMA flux uncertainty), how-ever we adjusted the flux of 2013.1.00198.S to matchthe flux measured in 2016.1.00311.S flux. As a result,6 of the 7 lines have internally consistent fluxes, andtherefore RMS uncertainties reported here do not in-clude flux calibration uncertainty, as that will either netincrease/decrease measured flux but will not impact thespatially resolved shape of the line images.2.2. Images and Radial Integrated Flux Profiles
The continuum was subtracted in the uv-plane usingthe CASA task uvcontsub assuming a linear fit to thecontinuum shape. Imaging was carried out using the tclean task in a semi-automated fashion described here.The source velocity was assumed to be 2.84 km s − .For each line, using all available data for a given tran-sition, we begin by creating a dirty image to determinethe standard deviation per beam in line free channelsand the beam size for each transition. From these data,we create a mask based on expectations from Keplerianmotion using the code presented as part of Pegues et al.(2020), see also code reference jpegues (2020). We as-sume an inclination of 5 ◦ , a position angle of 152 ◦ , anda stellar mass of 0.8 M (cid:12) to create the Keplerian maskfollowing Huang et al. (2018). We assume a distanceof 59.5 pc (Gaia Collaboration et al. 2016). The maskis convolved with the respective beam for each data set(typically between 0 . (cid:48)(cid:48) − . (cid:48)(cid:48) . (cid:48)(cid:48)
5, we calculate the uv-taper that is nec-essary to create a beam of just below 0 . (cid:48)(cid:48) σ stan-dard deviation, with the line specific uv-taper applied.Finally, we take the cleaned image and apply the CASAtask imsmooth, which operates in the image plane, tocreate a final image with an exactly 0 . (cid:48)(cid:48) GoFish package(Teague 2019) to improve the signal to noise of the ra-dial profile for each transition.
GoFish works in the
Cleeves et al.
Table 1. c − C H ObservationsLine o or p ν A ij a E ua Channel width rms b Int. flux c Program d (GHz) (s − ) (K) (km/s) (mJy/beam) (mJy km/s)4 , − , o 227.16913 3.113E-04 29.07 0.16 3.2 80 ± , − , o 351.52327 1.237E-03 77.24 0.21 3.4 237 ±
13 28 , − , o 352.19364 1.734E-03 86.93 0.21 4.2 220 ±
30 35 , − , p 338.20399 1.598E-03 48.78 0.21 5.3 140 ±
16 38 , − , p 352.18551 1.735E-03 86.93 0.21 4.3 96 ±
15 39 , − , o (bl) 351.96593 2.117E-03 93.34 0.21 2.9 398 ±
14 2,39 , − , p (bl) 351.96594 2.117E-03 93.3410 , − , o (bl) 351.78157 2.439E-03 96.50 0.21 2.9 330 ±
15 2,310 , − , p (bl) 351.78157 2.439E-03 96.50 a Line catalogue data from CDMS (M¨uller et al. 2005) b In a 0 . (cid:48)(cid:48) c Measured within a Keplerian mask (see Section 3, and Figure 8). d
1) 2013.1.00114.S, 2) 2016.1.00311.S, 3) 2013.1.00198.S -2.5 0.0 2.5-2.5 0.0 2.5 ,29 K
0 2 4 6 -2.5 0.0 2.5-2.5 0.0 2.5 ,77 K
0 5 10 15 -2.5 0.0 2.5-2.5 0.0 2.5 ,87 K
0 4 8 12 m J y b e a m k m s a) ortho -2.5 0.0 2.5-2.5 0.0 2.5 ,49 K
0 4 8 12 -2.5 0.0 2.5-2.5 0.0 2.5 ,87 K
0 3 6 9 m J y b e a m k m s b) para -2.5 0.0 2.5 -2.5 0.0 2.5 ,93 K
0 6 12 18 -2.5 0.0 2.5 -2.5 0.0 2.5 ,97 K
0 6 12 18 m J y b e a m k m s c) blend Figure 1.
Moment zero maps for the five isolated transitions of a) ortho, b) para, and c) two blended transitions of c − C H that are analyzed in this work. Line ID and upper state energy are labelled in the upper right corner of each panel. − C H in TW Hya GoFish deprojection.The deprojection is done on the cleaned image cubesproduced with the method described above, and the out-put is a spectrum at each radius where the width of theemission is some combination of the thermal width andany non-thermal broadening or beam convolution broad-ening, since Keplerian motion has been accounted for.After shifting and stacking the azimuthal emission, weintegrate over a conservative velocity range of 2.5 km s − to 3.18 km s − and obtain the profiles shown in Figure 2.Note that even though the para transitions are faint inthe moment-0 maps in Figure 1, the ring-like structurebecomes clearer in the GoFish extraction and the peakof the resolved profile is detected at ≥ − σ . METHODSFrom this multi-line data set of c − C H rotationaltransitions spanning upper state energies from 29.1 Kto 96.5 K, we aim to use excitation to constrain the lo-cation of the hydrocarbon layer in the TW Hya disk.To fit the data, we use a simple slab model method(see e.g., Qi et al. 2008; ¨Oberg et al. 2017) and em-ploy a non-LTE line radiative transfer code based uponRatran (Hogerheijde & van der Tak 2000, AstrophysicsSource Code Library, record 0008.002) designed to han-dle blended transitions and be computationally efficient(https://github.com/ryanaloomis/nLTErt1d). The in-puts to this code are height ( z ), gas volumetric density,dust density, gas temperature, dust temperature, non-thermal line width, and abundance relative to hydrogenof the molecule of interest.We approximate the emission from TW Hya as a seriesof radial annular regions, where the physical conditionsas a function of height at each radial location are takenfrom the Cleeves et al. (2015) model of TW Hya. Thegas temperature, dust temperature, and gas density aretaken from this model, where the total mass of this TWHya disk model is 0.04 M (cid:12) .The “slab” c − C H distribution is bounded by an in-ner and outer radius, and a vertically computed upperand lower column density of H, N H . Note that for theISM, N H = 1 . × cm − = 1 A V . But since ourmodel has 6.7 × less small grains in the surface of thedisk due to the formation and subsequent settling oflarger grains, then 1 A V occurs at N H = 1 . × cm − . The general reasoning for this vertical parame-terization rather than a simple z/r cut is that we ex- Table 2.
Slab Model ParametersParameter Values UnitMin H Column 10 cm − Max H Column ( × ) 1, 2.5, 5, 7.5, cm −
10, 12.5, 18.75, 25,50, 100, 187.5 χ ( c − C H ) ( × − ) 1, 1.875, 3.75, 5, per H7.5, 10, 25, 50, 75,100, 250, 500, 10 OPR 1, 3 n/a R inner c − C H
20, 25, 30, 35, 40 au R outer c − C H
80, 90, 100, 110, 120 auTotal number of models 7150 pect the c − C H distribution to be driven by the stel-lar radiation field, largely UV (Du et al. 2015; Berginet al. 2016). The abundance is set to a constant valuewithin this “slab,” where we vary the value of this con-stant. The underlying physical structure and the calcu-lated N H contours bounding the c − C H distributionare shown in Figure 3. As can be seen in the right-handpanel, the slabs closer to the disk surface are dominatedby warmer temperatures, while deeper in (higher N H )the gas temperatures decrease.In addition to the four parameters described above,we also fit for the ortho to para ratio of c − C H , whereas can be seen in Table 1 our data set contains both spinisomers and two blended pairs of ortho and para tran-sitions (the brightest lines). We adopt the collisionalrate coefficients hosted on the Leiden LAMDA database(Sch¨oier et al. 2005) originally computed in Chandra& Kegel (2000) for ortho and para separately. For theblended transitions, the components of the blends hadidentical molecular parameters (e.g., Einstein A coeffi-cients, near identical frequencies), and therefore we sim-ulated emission for the blends by first estimating thelevel populations for each of the ortho and para compo-nents independently (while we did not impose it on thesolution, the emission was seen to be largely in LTE). Wethen took a weighted average of the ortho and para pop-ulations in the upper and lower states respectively (thepopulation fractions were the same for both spin iso-mers). This allowed us to simulate the blended c − C H emission as that from a “single” molecule. We took thisapproach since for some of the models the blended tran-sitions became optically thick, so we could not simplysum the line fluxes from the independent componentswithout potentially overestimating the blended line flux.To create the synthetic integrated radial intensity pro-file for each model, we reconstruct a 2D sky image of thevelocity integrated line brightness and convolve the im- Cleeves et al.
0 1 2 3 4 r ( ) -2 0 2 4 6 8 10 m J y b e a m k m s a. ortho
0 1 2 3 4 r ( ) -2 0 2 4 6 8 b. para
0 1 2 3 4 r ( ) -2 0 2 4 6 8 10 12 14 16 18 c. blends
10 9 (blend) (blend)
Figure 2.
Profiles extracted by the pixel stacking technique within the GoFish code (Teague 2019). Transitions indicated inthe legend. Shaded region indicates the error on the radial profile. Z ( a u ) N H = 1 . × . × . × . × . × . × . × . × . × . × . × l og ( n H / c m − ) l og ( T ga s / [ K ] ) Figure 3.
Slab c − C H model definitions. The vertical extent of c − C H is described by the vertically integrated columndensity derived from the Cleeves et al. (2015) TW Hya model shown in grey shaded contours. The upper limit for all modelsis N H = 10 cm − (dashed line). The bottom of the c − C H slab is described by the color line contours as labelled in theleft panel. The inner and outer radius parameters truncate these boundaries vertically, and the abundance is assumed to beconstant inside of the layer. age with a Gaussian 0 . (cid:48)(cid:48) GoFish intensity profile.Goodness of fit is assessed based on a reduced χ be-tween the beam-convolved model spectrally integratedintensity profile and the data. To estimate χ , the fluxin a given annulus for the data and model are measuredper beam, and we sum the square of the difference di-vided by the observed RMS uncertainty. The number ofannuli super samples the beam but is divided out by thereduced χ . All lines are treated equally in our overallassessment of fit for a given set of slab model parameters,even though some lines are brighter than others. Wemade this decision since the brighter lines are opticallythick blended transitions, and the lower SNR profiles provide important constraints on our fit. Note, we onlyconsider the RMS uncertainty in this estimate and donot include flux calibration error since eight out of ninelines have self-consistent fluxes, such that flux uncer-tainty will globally shift all radial line profiles upwardsor downwards rather than expanding the error bars uni-formly. The only exception is the ortho c − C H line4 , − , . We therefore have also examined our bestfit models from the profiles as extracted directly from GoFish , and then with a 10% uniform increase and de-crease on just the 4 , − , , with the rest of the linesfixed. We do not find a change in the best fit models, inpart since we have two other ortho lines. Furthermore,the models that fit the rest of the data give reasonablefits to the native 4 , − , observations without fluxscaling, and therefore we do not vary the flux for theremainder of the analysis described here. − C H in TW Hya N H = 1 × cm − to 1 × cm − , we decided tofix the upper boundary for the full grid of simulationsto the former value of 1 × cm − , which for the re-duced dust atmosphere corresponds to an A V of 0.0008,i.e., very high up in the atmosphere, where the H gasdensity is also very tenuous ( ∼ cm − ). Thus verylittle c − C H in this layer contributes to the emission.The full set of parameters and their values considered inthe simulation grid are provided in Table 2. Essentially,each model varies the inner and outer extent, lower layerextent via the column density, the abundance within thelayer, and the ortho-to-para ratio of c − C H . The rangeof values was explored is based on a combination of pre-vious modeling results (e.g., Bergin et al. 2016) alongwith making sure we explored a wide space around re-gions that gave better fits to assess degeneracies in theparameter space. RESULTS4.1.
Fiducial Model Results
A comparison between all of the 7,150 models in ourgrid and all observed transitions is illustrated in Fig-ure 4. Each point represents χ for all observedlines simultaneously. Due to the large dynamic rangeacross the model grid, we have split each panel into atop and bottom panel where the top row has log scale χ and the bottom row shows the same data butwith a linear scale for clarity. The best fit models (thosewith the lowest χ ) are those that appear on thebottom row. From Figure 4, a few key features becomeclear, and are enumerated here:1. The models that agree the best with the data havethe outer edge of the c − C H emission at > <
100 au.3. An inner edge of 25 and 30 au tends to be a betterfit for the c − C H distribution. We note that thisedge is close to the half-beam of the data (0 . (cid:48)(cid:48)
25 =15 au), and so cannot constrain this any furtherfrom existing data. We also note that the c − C H “empty” inner disk model is consistent with thenon-zero flux near the star seen in our radial pro-file due to beam convolution effects. Finally, wenote that both the value for the inner and outeredge only hold if a slab model is an adequate repre-sentation of the data, and certainly more complex structural representations could be explored withhigher SNR data.4. The flatness of the χ for the best modelssuggests some degeneracy between the maximumcolumn density, i.e., the depth of the layer, and theabundance of c − C H in the layer, which is to beexpected especially if the emission is thermalized.5. While the observations suggest the c − C H origi-nates in a region of sufficiently high gas density tobe thermalized, we can further constrain the verti-cal location of the c − C H layer by taking advan-tage of the wide span of upper state energies in ourdata set. We can rule out c − C H emitting fromthe upper surface (i.e., N H ≤ cm − ) or fromthe midplane, below N H ≥ cm − . These col-umn densities translate to A V ≤ .
05 and A V ≥ c − C H of 1. Ortho-to-para of 3is strongly favored over 1. We do not try to fitbetween these values given the signal to noise ofthe data.To provide an even clearer picture of our best fit mod-els, we extract the 30 models with the lowest χ fromour grid and plot a histogram of the model parameters(see Figure 5). We also plot the radial intensity profilesof these models in Figure 6. Based on this sub-selectionof models, an inner radius of 25 or 30 au is favored by2/3 of the 30 models. Nearly all of the best fit modelshave an outer radius of c − C H at >
100 au, with 110au favored slightly over 120 au. While there is somedegeneracy between the depth of the slab measured by N max and the abundance per H of c − C H , it is clearthat models with N max of (1 − × cm − are fa-vored by most. The abundance of c − C H in the slabis between (3 − × − per H atom. To test howmuch of the spread is due to the degeneracy between thethickness of the layer and the abundance, we also showin Figure 5 a histogram of the product of χ c − C H and N max and find that this distribution is much tighter,and corresponds to an approximate column density of c − C H in the slab of (1 − × cm − .The 30 best fit models are shown in Figure 6 in darkgrey, and we have additionally highlighted the four bestfit models in color and provided their physical param-eters in the key. Note all 30 models have an ortho topara ratio of 3, so that parameter is not listed in thefigure. The best fit models have similar characteristics,and abundance of (7 . − × − per H atom, and Cleeves et al. R e d u c e d Triangle: o/p = 1Square: o/p = 3Marker size: Inner edge r in = 20, 25, 30, 40 auColor: Abundance1.00e-111.87e-113.75e-115.00e-117.50e-111.00e-102.50e-107.50e-105.00e-101.00e-092.50e-095.00e-091.00e-08 R out : 80 au N / cm R e d u c e d R out : 90 au N / cm R out : 100 au N / cm R out : 110 au N / cm R out : 120 au N / cm Figure 4.
Reduced χ between the observed visibilities and modelled visibilities. Each point represents the model fit for allobserved lines as a summed reduced χ . The model parameters are described by the point color (abundance), size (inner radiusof c − C H ), and shape (ortho to para ratio). The columns indicate different outer radii for the c − C H slab, for the full rangeof H column densities considered on each x-axis. The poor fits are shown on a logarithmic scale on the top row, while thebetter fits are shown as a linear scale on the bottom row.
0 5 10
Inner Radius (au)
0 5 10 15
Outer Radius (au) . . . . .
0 5 10 15 20 log( N max , H / cm ) - . - . - . - . - .
0 5 10 15 log( C H per H) . . . . . .
0 10 20 log( C H × N max / cm ) Figure 5.
Histogram of the parameters of the 30 best fit models. Note 1 A V Corresponds to N H = 1 . × cm − giventhe reduced dust in the upper layers of TW Hya. extend down to an H column density of approximately2 ± . × cm − . This depth corresponds to a verti-cally integrated A V of about 1 due to the reduced smalldust mass in the disk upper layers in the Cleeves et al.(2015) TW Hya model. DISCUSSIONWe have carried out a resolved, multi-line study ofhydrocarbon emission in the TW Hya disk as traced by c − C H . The observations span a wide range of upperstate energy from 29 K to 96.5 K, and include both orthoand para forms of c − C H . Hydrocarbons like c − C H and the more widely studied C H have been suggestedto be good tracers of C/O ratios in disks, which is aparameter that we are now beginning to measure inthe atmospheres of gas-rich exoplanets, which should bestrongly influenced by (if not set by) the gas in their par- ent disk (e.g., ¨Oberg et al. 2011; ¨Oberg & Bergin 2016,and more). Therefore, it is crucial to understand the na-ture of the gas that these species trace, and specifically where in the disk they are located to understand howclosely it links to what forming planets might accrete.The present study augments previous studies of thechemically related species, C H, also observed towardthe TW Hya disk (Kastner et al. 2015; Bergin et al.2016). We confirm the earlier finding of Bergin et al.(2016) that the radial extent of c − C H (qualitativelydetermined by stacking) appears very similar to thebrighter C H emission. If we assume C H and c − C H originate from a similar layer, e.g., in z/r space or col-umn density space, the present study provides furthercontext on the vertical distribution of the hydrocarbonlayer in TW Hya, and what that means for eventuallyusing this species to measure C/O in this disk (see, e.g., − C H in TW Hya
0 1 2 3 4 -2 0 2 4 m J y b e a m k m s
3, 2
2, 1 (o)
0 1 2 3 4 0 2 4 6 8 10
3, 4
4, 3 (o)
0 1 2 3 4 -2 0 2 4 6 8 C H = 7.5 × 10 , N H = 1.875 × 10 , (20-110au): 2.42 C H = 7.5 × 10 , N H = 1.875 × 10 , (25-110au): 2.54 C H = 7.5 × 10 , N H = 2.500 × 10 , (25-110au): 2.52
3, 6
2, 5 (o)
0 1 2 3 4 -2 0 2 4 6 8 m J y b e a m k m s
5, 1
4, 0 (p)
0 1 2 3 4 -2 0 2 4
2, 6
3, 5 (p)
0 1 2 3 4 r ( )
0 5 10 15 m J y b e a m k m s
1, 10
0, 9 and
0, 10
1, 9 (b)
0 1 2 3 4 r ( )
0 5 10 15
1, 8
2, 7 and
2, 8
1, 7 (b)
Figure 6.
Radial profiles of the integrated emission from each transition, including error bars (shaded gray). Profiles extractedusing GoFish. The thirty best fit models are show as gray lines, while the top three models are shown in color, and haveparameters indicated in the legend in the lower right-hand corner.
Cleeves et al. 2018). In the following sections we discussour main findings.5.1. c − C H is located in and above the warmmolecular layer in a ring We find that all models that provide reasonable fitsto the full data set do not extend much further belowan H column density of a few × cm − . Given thatthe small dust mass in this model is only 15% of theinterstellar value, this layer corresponds to a maximum A V of 1 – 2. We do not have strong constraints on theupper extent of the layer, and where additional obser-vations for lines with higher upper state energies (andsmaller Einstein A coefficients) would be helpful. Whilethis lower boundary does not directly correspond to anormalized height value of z/r due to the nature of theassumed density profile (a power-law with an exponen-tial taper), we can say that within the radial boundariesof the emission, that the c − C H arises above a z/r value of approximately 0.25 based upon Figure 3, andno higher than 0.4. This range is consistent with the lo-cation of C H in this source as modeled in Bergin et al.(2016). We also constrain the radial distribution of the c − C H layer. The inner edge of the c − C H dis-tribution (25 – 30 au) is consistent with the CO snowline (Qi et al. 2013b; Zhang et al. 2017) and/or the dipin CO emission (Huang et al. 2018) and scattered light(Debes et al. 2013; van Boekel et al. 2016). The innerdeficit of c − C H might be signalling a large chemicalshift, perhaps a large reduction of C/O at this radius.Models have predicted an enhancement of CO interiorto the CO snow line due to dust evolutionary processes(Krijt et al. 2018); however, CO adds equal parts C andO and will not tend to reduce C/O especially if it issomewhat elevated (above solar, 0.54) but still below1 already. An alternative explanation is that the sameprocess that sequestered water ice into larger grains hasbegun in a delayed way on the carbon-bearing molecules(e.g., Cleeves et al. 2018). This process would be fastestin the inner disk and would move outward, which mightexplain why different disks have a wide variety of C Hring patterns (e.g., Bergin et al. 2016). Seeing whetherthese rings have a predictable time evolution for a sta-tistically significant sample of disks would help elucidate0
Cleeves et al. the cause (see also discussion in Bergner et al. 2019).It would likewise be interesting to compare with otherspecies that also show inner chemical deficits, like CN(Cazzoletti et al. 2018).The outer edge of c − C H is located at 110 – 120au. It is not clear why this radius is remarkable, besidesthat it corresponds also to the edge of C H (Bergin et al.2016). The millimeter pebble disk extends to approxi-mately 60 au, however there is a flatter weakly emitting“shoulder” in the ALMA observed 852 micron flux outto around 100 au (Andrews et al. 2012, 2016; van Boekelet al. 2016; Huang et al. 2018). The CO disk extends outbeyond 200 au (Andrews et al. 2012; Schwarz et al. 2016;Huang et al. 2018), though there is a sharp drop in COat around 100 au (Zhang et al. 2017). While the scat-tered light extends out to 200 au like the CO gas disk,it has some potentially interesting features around 100au. For example, there is a brightness peak in the ob-served NIR scattered light at about 100 au (van Boekelet al. 2016). Debes et al. (2013) attributes this changeat 100 au to be potentially related to sculpting by aNeptune-mass planet. In addition, there appears to bea shift in the “color” of the scattering at 100 au, whereit is neutral interior to 100 au and becomes more blueoutside of 100 au (Debes et al. 2013; van Boekel et al.2016). Therefore this location seems to be an importanttransition; however, its nature remains unclear. Per-haps it marks a significant change in gas density (i.e.,H ) or alternatively a chemical change where the diskbecomes so tenuous that the external UV field becomestoo harsh and densities too low to support chemistrymore advanced than CO.5.2. The rotational emission of c − C H appears to bethermalized We find that our best fit models have some spread in N max and χ c − C H that would be expected given thatwe are more fundamentally tracing a column density of c − C H with our observations, especially due to the faceon nature of this disk. Within these models, most of theemission is coming from the bottom of the layer, andthe emission appears well represented by LTE (see alsoTeague et al. 2018, for a non-LTE analysis of CS). Toprovide context, some models with the same c − C H column density placed very high up in the disk, wellabove A V = 1 with a very high c − C H abundance anda low value of N max do not provide good fits to the ob-servations. These models are also no longer thermalized.We conclude that the c − C H appears to be reasonablywell approximated by LTE. The critical density for thetransitions observed is typically around 10 cm − , andour results suggest there remains a substantial amount of H gas, > cm − , at high ( z/r > .
25) altitudesand moderately far out radii (25 – 110 au), consistentwith disk mass estimates provided by HD for this source(Bergin et al. 2013; Favre et al. 2013; Cleeves et al. 2015;Trapman et al. 2019; Calahan et al. 2020). How this rel-atively old disk maintains such a large amount of gas,however, is still unclear.5.3.
The abundance of c − C H relative to C H isconsistent with gas-phase chemical models
Bergin et al. (2016) modeled the related species C Hin TW Hya using full thermo-chemical models. Consis-tent with the results of Du et al. (2015), in this workit appeared clear that to reproduce the observed bright-ness of C H the carbon relative to oxygen ratio of thegas must be very high, greater than the solar ratio of0.54. To obtain a high C/O ratio, either there must bea source of carbon enhancement, such as carbon grainor PAH destruction (Anderson et al. 2017), and/or adeficit of oxygen, possibly due to grain growth and set-tling (Hogerheijde et al. 2011; Salinas et al. 2016; Cleeveset al. 2018). The same will be true for c − C H . Takingour best fit model abundances of c − C H and compar-ing it to the results of Bergin et al. (2016) for C H inTW Hya, we find c − C H ’s abundance is ∼ − H. If instead we use the column density, forBergin et al. (2016)’s C/O = 1 model, the c − C H toC H ratio is approximately 10%. While detailed chem-ical modeling is beyond the scope of the present paper,we can compare this percentage to a published modelof a different disk, IM Lup, where C/O was also varied(Cleeves et al. 2018) purely by removal of volatile oxy-gen (via water sequestration into large grains, presum-ably settled into the midplane as well as O-removal fromCO for the more extreme C/O ratios). In the z/r layerof 0.25 – 0.5, that model finds a c − C H to C H per-centage of 1 – 10%, broadly consistent with our findingshere, without the need for additional carbon sources. Incomparing these percentages with other star-forming en-vironments, we see that this ratio of c − C H to C H isalso broadly consistent with what is observed in PDRs(e.g., ∼
3% toward the Orion Bar Cuadrado et al. 2015).While the chemistry is consistent with gas phaseroutes, we cannot formally rule out carbon grain / PAHdestruction as a source of additional carbon enhance-ment (e.g., Kastner et al. 2015; Bergin et al. 2016; An-derson et al. 2017; Bosman et al. 2021). We also notethat the abundance ratio is consistent with a scenarioof a “dry” – i.e., H O ice/gas poor – surface of the TWHya disk, which also agrees with the
Herschel results forthe disk’s cold water vapor abundance Hogerheijde et al. − C H in TW Hya
30 40 50 60 70 80 E u / Kelvin20.521.021.522.022.5 l n ( N u / g u / c m ) o/p = 3.02 ± 0.78 T rot,ortho = 55 ± 13 K T rot,para = 43 ± 14 K orthopara Figure 7.
Disk-averaged rotational diagram of the un-blended lines. We find similar ortho and para rotationaltemperatures and a disk integrated ortho-to-para ratio con-sistent with 3. (2011) and scattered light constraints (e.g., Weinbergeret al. 2002; Debes et al. 2013).5.4.
The ortho-to-para ratio of c − C H is consistentwith a constant value of 3 Via the non-LTE fitting we find that the ortho-to-para ratio across the slab can be well reproduced witha single value of 3, however we only consider values of1 or 3 in our grid. To confirm this result, we have donean additional disk-averaged rotational diagram analysisusing the unblended ortho and para transitions, whichare also optically thin based on our model comparisons.Figure 7 shows the results of the fitting, and we findthat the ortho and para lines have consistent rotationaltemperatures, and the derived ortho-to-para ratio fromthe column densities is 3 . ± . c − C H mainly form viagas phase reactions (e.g., Vrtilek et al. 1987; Park et al.2006), in particular via C H + reacting with H via ra-diative association. Morisawa et al. (2006) demonstratethat the ortho-to-para ratio of c − C H in steady statecan be approximated by:[ o/p ] C H = (5 + 3 φ ) × [ o/p ] H + 3 φ + 3(1 + φ ) × [ o/p ] H + φ + 3 , (1)where φ essentially reflects the ratio of the rate of forma-tion of c − C H by electron dissociative recombinationof C H +3 compared to reactions of C H +3 with atoms(however, see also the discussion in Park et al. 2006,for additional effects). The parameter φ is equal to k DR [ n e ] /k ion − Y [ n Y ], where k DR is the rate of dissocia-tive recombination, n e is the number density of elec-trons, k ion − Y is the reaction rate of C H +3 with speciesY and [ n Y ] is the number density of said species. If the ortho-to-para ratio of H goes to 3, the dependence on φ drops out of Equation 1. If the ortho-to-para ratio of H is small due to a cold formation mechanism of H , thenthe ortho-to-para ratio of c − C H depends on φ , goingto 1 if φ is small and 3 if φ is large. Therefore, if thisequation holds, obtaining the observed ortho to para ra-tio of 3 in c − C H means either H has an ortho-to-paraof 3, or φ is large.We can estimate what φ regime c − C H may be inusing the Cleeves et al. (2018) model of a solar nebulalike disk. Taking approximate values from a locationof r = 50 au and z/r of 0.3 au, the electron density isapproximately 10 − per H, and is largely governed byphotoionization of atomic carbon. The recombinationrate with electrons is approximately 10 − cm s − at 50K (Loison et al. 2017). For the denominator, if we takereactions with C as representative, the rate is ∼ − cm s − , and the abundance is also 10 − per H. There-fore, φ comes down to the ratio of the rate coefficients,which is >>
1. Therefore it is not clear whether ourmeasured c − C H ortho to para ratio is shedding lighton the ortho-to-para ratio of H or if it is related to therelative rate of dissociative recombination.5.5. C/O implications for forming planets
While there are exciting prospects to use hydrocar-bons to constrain C/O ratios in disks to connect withplanets, we must understand what region our observa-tions fundamentally probe. Our results find that the c − C H emission largely comes from the UV irradiatedlayer, above an A V of 1 ( z/r > . Cleeves et al. gaps. These results imply the surface gas in the disk isfeeding potential planet forming regions, and that withour observations we could be tracing the same gas reser-voir that might feed proto-giant planets’ atmospheres(see also Cridland et al. 2020). Specifically in TW Hya,Teague et al. (2019) found evidence of vertical motionaround 90 au, near the van Boekel et al. (2016) scatteredlight gap. If such motions are common, this would beexciting as the layer traced by c − C H would be moredirectly related to the chemistry that sets planets’ at-mospheric compositions.Along these lines, the 20 Myr-old HR 8799 system’sfour outer planets have been independently chemicallycharacterized using VLT SPHERE, and their C/O ra-tios radially varies (Lavie et al. 2017; Bonnefoy et al.2016). The inner two planets ( d at 24 au and e at 15 au;Marois et al. 2008, 2010) have low C/O ratios – con-sistent with zero due to low C/H, while the outer twoplanets ( b at 38 au and c at 68 au) have far larger C/Oratios of 0.8 – 0.9 (Lavie et al. 2017). More recently,Molli`ere et al. (2020) reevaluated the atmospheric C/Oratios for planet e and found a C/O value of 0.6, so lessextreme of a jump than previously estimated. Theseresults emphasize the need for careful treatment of non-equilibrium atmospheric effects and care in interpretingindividual C/O estimates.Interestingly, HR 8799’s tentative radial increase fromsub-solar/solar up to super solar values in C/O is spa-tially consistent with TW Hya’s increase in hydrocar-bon emission. In fact, the radial separation of the innerplanets falls interior to the start of TW Hya’s hydro-carbon ring, while the outer planets’ orbits are withinthe radial bounds of the hydrocarbon ring. Of coursethere are many differences between TW Hya and theHR 8799 system, including that HR 8799’s host star hasa mass of 1.5 M (cid:12) compared to TW Hya’s 0.8 M (cid:12) , sosuch a connection is purely speculative. However, an in-teresting feature worth pointing out is that the youngerdisk IM Lup does not show the same large inner ring inhydrocarbons (Cleeves et al. 2018), and instead has amore centrally peaked/flattened distribution. Rings doappear to be common among intermediate aged disks(Bergin et al. 2016; Bergner et al. 2019), including someof the (cid:38) Acknowledgements:
This paper makes use of the fol-lowing ALMA data: ADS/JAO.ALMA
Aikawa, Y., van Zadelhoff, G. J., van Dishoeck, E. F., &Herbst, E. 2002, A&A, 386, 622 Anderson, D. E., Bergin, E. A., Blake, G. A., et al. 2017,ApJ, 845, 13 − C H in TW Hya Andrews, S. M., Wilner, D. J., Hughes, A. M., et al. 2012,ApJ, 744, 162Andrews, S. M., Wilner, D. J., Zhu, Z., et al. 2016, ApJL,820, L40Ansdell, M., Williams, J. P., van der Marel, N., et al. 2016,ApJ, 828, 46Bergin, E. A., Cleeves, L. I., Crockett, N., & Blake, G. A.2014, Faraday Discussions, 168, 61Bergin, E. A., Du, F., Cleeves, L. I., et al. 2016, ArXive-prints, arXiv:1609.06337Bergin, E. A., Hogerheijde, M. R., Brinch, C., et al. 2010,A&A, 521, L33+Bergin, E. A., Cleeves, L. I., Gorti, U., et al. 2013, Nature,493, 644Bergner, J. B., ¨Oberg, K. I., Bergin, E. A., et al. 2019, ApJ,876, 25Bertout, C., Siess, L., & Cabrit, S. 2007, A&A, 473, L21Bonnefoy, M., Zurlo, A., Baudino, J. L., et al. 2016, A&A,587, A58Bosman, A. D., Alarcon, F., Zhang, K., & Bergin, E. A.2021, arXiv e-prints, arXiv:2101.12502Calahan, J., Bergin, E. A., & al, e. 2020, ApJ, 0Cazzoletti, P., van Dishoeck, E. F., Visser, R., Facchini, S.,& Bruderer, S. 2018, A&A, 609, A93Chandra, & Kegel. 2000, Astron. Astrophys. Suppl. Ser.,142, 113. https://doi.org/10.1051/aas:2000142Cleeves, L. I., Bergin, E. A., Alexander, C. M. O. ., et al.2014, Science, 345, 1590Cleeves, L. I., Bergin, E. A., Qi, C., Adams, F. C., &¨Oberg, K. I. 2015, ApJ, 799, 204Cleeves, L. I., ¨Oberg, K. I., Wilner, D. J., et al. 2018, ApJ,865, 155Cridland, A. J., Bosman, A. D., & van Dishoeck, E. F.2020, A&A, 635, A68Cuadrado, S., Goicoechea, J. R., Pilleri, P., et al. 2015,A&A, 575, A82Debes, J. H., Jang-Condell, H., Weinberger, A. J., Roberge,A., & Schneider, G. 2013, ApJ, 771, 45Drozdovskaya, M. N., van Dishoeck, E. F., Rubin, M.,Jørgensen, J. K., & Altwegg, K. 2019, MNRAS, 490, 50Du, F., Bergin, E. A., & Hogerheijde, M. R. 2015, ApJL,807, L32Du, F., Bergin, E. A., Hogerheijde, M., et al. 2017, ApJ,842, 98Eistrup, C., Walsh, C., & van Dishoeck, E. F. 2018, A&A,613, A14Favre, C., Bergin, E. A., Cleeves, L. I., et al. 2015, ArXive-prints, arXiv:1503.02659Favre, C., Cleeves, L. I., Bergin, E. A., Qi, C., & Blake,G. A. 2013, ApJL, 776, L38 Fedele, D., & Favre, C. 2020, arXiv e-prints,arXiv:2005.03891Gaia Collaboration, Brown, A. G. A., Vallenari, A., et al.2016, ArXiv e-prints, arXiv:1609.04172Guilloteau, S., Reboussin, L., Dutrey, A., et al. 2016, A&A,592, A124Hogerheijde, M. R., & van der Tak, F. F. S. 2000, A&A,362, 697Hogerheijde, M. R., Bergin, E. A., Brinch, C., et al. 2011,Science, 334, 338Huang, J., Andrews, S. M., Cleeves, L. I., et al. 2018, ApJ,852, 122jpegues. 2020, jpegues/kepmask: Python Package:pykepmask v2.0.0, vv2.0.0, Zenodo,doi:10.5281/zenodo.3672187.https://doi.org/10.5281/zenodo.3672187Kama, M., Bruderer, S., van Dishoeck, E. F., et al. 2016,A&A, 592, A83Kastner, J. H., Hily-Blant, P., Rodriguez, D. R., Punzi, K.,& Forveille, T. 2014, ApJ, 793, 55Kastner, J. H., Qi, C., Gorti, U., et al. 2015, ApJ, 806, 75Kreidberg, L., Line, M. R., Bean, J. L., et al. 2015, ApJ,814, 66Krijt, S., Ciesla, F. J., & Bergin, E. A. 2016, ApJ, 833, 285Krijt, S., Schwarz, K. R., Bergin, E. A., & Ciesla, F. J.2018, ApJ, 864, 78Lacour, S., Nowak, M., Wang, J., et al. 2019, VizieR OnlineData Catalog, J/A+A/623/L11Lavie, B., Mendon¸ca, J. M., Mordasini, C., et al. 2017, AJ,154, 91Loison, J.-C., Ag´undez, M., Wakelam, V., et al. 2017,MNRAS, 470, 4075Loomis, R. A., Cleeves, L. I., ¨Oberg, K. I., et al. 2018, ApJ,859, 131Madhusudhan, N. 2012, ApJ, 758, 36Madhusudhan, N., Harrington, J., Stevenson, K. B., et al.2011, Nature, 469, 64Marois, C., Macintosh, B., Barman, T., et al. 2008, Science,322, 1348Marois, C., Zuckerman, B., Konopacky, Q. M., Macintosh,B., & Barman, T. 2010, Nature, 468, 1080Miotello, A., Facchini, S., van Dishoeck, E. F., et al. 2019,A&A, 631, A69Molli`ere, P., Stolker, T., Lacour, S., et al. 2020, A&A, 640,A131Morisawa, Y., Fushitani, M., Kato, Y., et al. 2006, ApJ,642, 954M¨uller, H. S. P., Schl¨oder, F., Stutzki, J., & Winnewisser,G. 2005, Journal of Molecular Structure, 742, 215¨Oberg, K. I., & Bergin, E. A. 2016, ApJL, 831, L19 Cleeves et al. ¨Oberg, K. I., Murray-Clay, R., & Bergin, E. A. 2011, ApJL,743, L16¨Oberg, K. I., Guzm´an, V. V., Merchantz, C. J., et al. 2017,ApJ, 839, 43Park, I. H., Wakelam, V., & Herbst, E. 2006, A&A, 449, 631Pegues, J., ¨Oberg, K. I., Bergner, J. B., et al. 2020, ApJ,890, 142Qi, C., ¨Oberg, K. I., Wilner, D. J., & Rosenfeld, K. A.2013a, ApJL, 765, L14Qi, C., Wilner, D. J., Aikawa, Y., Blake, G. A., &Hogerheijde, M. R. 2008, ApJ, 681, 1396Qi, C., Wilner, D. J., Calvet, N., et al. 2006, ApJL, 636,L157Qi, C., ¨Oberg, K. I., Wilner, D. J., et al. 2013b, Science,341, 630Salinas, V. N., Hogerheijde, M. R., Bergin, E. A., et al.2016, A&A, 591, A122Sch¨oier, F. L., van der Tak, F. F. S., van Dishoeck, E. F., &Black, J. H. 2005, A&A, 432, 369Schwarz, K. R., Bergin, E. A., Cleeves, L. I., et al. 2016,ApJ, 823, 91 Teague, R. 2019, The Journal of Open Source Software, 4,1632Teague, R., Bae, J., & Bergin, E. A. 2019, Nature, 574, 378Teague, R., Henning, T., Guilloteau, S., et al. 2018, ApJ,864, 133Trapman, L., Facchini, S., Hogerheijde, M. R., vanDishoeck, E. F., & Bruderer, S. 2019, A&A, 629, A79van Boekel, R., Henning, T., Menu, J., et al. 2016, ArXive-prints, arXiv:1610.08939Visser, R., Doty, S. D., & van Dishoeck, E. F. 2011, A&A,534, A132Visser, R., van Dishoeck, E. F., Doty, S. D., & Dullemond,C. P. 2009, A&A, 495, 881Vrtilek, J. M., Gottlieb, C. A., & Thaddeus, P. 1987, ApJ,314, 716Walsh, C., Millar, T. J., Nomura, H., et al. 2014, A&A,563, A33Weinberger, A. J., Becklin, E. E., Schneider, G., et al. 2002,ApJ, 566, 409Yu, M., Evans, Neal J., I., Dodson-Robinson, S. E.,Willacy, K., & Turner, N. J. 2017, ApJ, 850, 169Zhang, K., Bergin, E. A., Blake, G. A., Cleeves, L. I., &Schwarz, K. R. 2017, Nature Astronomy, 1, 0130 − C H in TW Hya c − C H transitions toward TW Hya andthe masks used to clean the data and generate the moment-0 maps. -2 0 2
3, 2
2, 1 (o)
0 4 8 12 16 m J y b e a m -2 0 2
3, 4
4, 3 (o)
0 6 12 18 24 m J y b e a m -2 0 2
3, 6
2, 5 (o)
0 7 14 21 28 m J y b e a m -2 0 2
5, 1
4, 0 (p)
0 7 14 21 28 m J y b e a m -2 0 2
2, 6
3, 5 (p)
0 4 8 12 16 m J y b e a m -2 0 2
1, 8
2, 7 and
2, 8
1, 7 (b)
0 8 16 24 32 m J y b e a m -2 0 2 -2 0 2
1, 10
0, 9 and
0, 10
1, 9 (b)
0 8 16 24 32 m J y b e a m Figure 8.
Channel maps for the transitions of c − C H2