Optical and Radio Observations of the T Tauri Binary KH 15D (V582 Mon): Stellar Properties, Disk Mass Limit, and Discovery of a CO Outflow

For convenience, we have compiled in Table 1 all of the high-quality (CCD) photometry of KH 15D that has been obtained since 1995 (Hamilton et al. 2005; Kusakabe et al. 2005; Windemuth & Herbst 2014; Arulanantham et al. 2017). The previously unpublished data, spanning the interval 2015 September to 2017 March, were obtained with the ANDICAM instrument on the 1.3 m telescope at Cerro-Tololo Inter-American Observatory (CTIO), operated by the SMARTS consortium. KH 15D was observed nightly during each observing season in three optical filters (VRI) and three near-IR filters (JHK). These data are obtained within an interval of about 30 minutes and combined to yield a single nightly magnitude in each color. The observing and reduction process is identical to what has been used for many years and is described in more detail by Windemuth & Herbst (2014). The full data set is most complete in the I magnitude, which is plotted in Figure 1. The newest data, for the last three observing seasons, are shown in Figure 2 with enough time resolution that the individual orbital cycles can be seen.

Table 1.  VRIJHK Photometry of KH 15D

Julian Date	V	${\sigma }_{V}$	R	${\sigma }_{R}$	I	${\sigma }_{I}$	J	${\sigma }_{J}$	H	${\sigma }_{H}$	K	${\sigma }_{K}$	Obs.	Ref
2450017.84	⋯	⋯	⋯	⋯	14.45	0.011	⋯	⋯	⋯	⋯	⋯	⋯	VVO	H05
2450021.75	⋯	⋯	⋯	⋯	14.52	0.024	⋯	⋯	⋯	⋯	⋯	⋯	VVO	H05
⋯	⋯	⋯	⋯	⋯	⋯	⋯	⋯	⋯	⋯	⋯	⋯	⋯	⋯	⋯
2457834.54	15.47	0.019	14.74	0.019	14.16	0.019	13.25	0.041	12.54	0.036	12.52	0.095	CTIO	A17
2457835.54	15.45	0.019	14.70	0.019	14.13	0.019	13.19	0.042	12.46	0.037	12.41	0.095	CTIO	A17

Note.

^aAn excerpt from the full version (available in machine-readable version), which includes all VRIJHK photometry from 1995 to 2017 as well as the observatory from which the measurements were taken and the paper from which they were accessed. In the machine-readable version, nights with no available data are represented by magnitudes and errors of 99.99 and 9.999, respectively. The observatory and reference keys are also included in this version.

Only a portion of this table is shown here to demonstrate its form and content. A machine-readable version of the full table is available.

Download table as:  DataTypeset image

The SMA on Maunakea, Hawaii was used in a compact antenna configuration to observe KH 15D, primarily to search for possible 870 μm dust continuum emission. Table 2 provides a summary of the observational parameters for each of the two runs. The first run consisted of a short (2 hr) filler track on 2010 October 13. The setup used 128 channels (0.8125 MHz spacing) in each of 24 chunks that together spanned±(4–6) GHz from a local oscillator (LO) frequency of 342.05 GHz. The second run was a full (7.5 hr) track on 2011 October 25. This observation took advantage of expanded correlator capability that processed ±(4–8) GHz from an LO frequency of 342.05 GHz, which included the spectral line CO J = 3–2 (345.796 GHz) in the upper side band. The weather conditions were good for both of these observations, with stable atmospheric phase and low to modest opacity at 225 GHz, as measured by the tipping radiometer at the nearby Caltech Submillimeter Observatory. The basic observing sequence consisted of 2.5 or 3 minutes on each of the calibrators J0532+075 and J0750+125, and 10 minutes on KH 15D. Passband calibration was performed using observations of available strong sources, and the absolute flux scale was set using observations of Callisto and Uranus, with an estimated accuracy of 20%. The field of view is set by the 36'' FWHM primary beam size at the observing frequency. The synthesized beam sizes are ∼2''. All of the calibration was performed using the MIR software package. Continuum and spectral line imaging and CLEAN deconvolution were done in the miriad package.

Table 2.  SMA Observational Parameters

Parameter	2010 Oct 13	2011 Oct 25
No. of Antennas	7	7
${\tau }_{225\mathrm{GHz}}$	0.04	0.08
Pointing center (J2000)	R.A. ${06}^{{\rm{h}}}{41}^{{\rm{m}}}10\buildrel{\rm{s}}\over{.} 18$, Decl. $+{09}^{{\rm{d}}}{28}^{{\rm{m}}}35\buildrel{\rm{s}}\over{.} 5$
Min/Max baseline	15–76 m	11–76 m
Gain Calibrators	J0532+075 (0.64 Jy)	J0532+075 (0.72 Jy)
	J0750+125 (0.57 Jy)	J0750+125 (0.67 Jy)
Passband Calibrator	3C 454.3	Uranus, J0927+390
Flux Calibrator	Callisto	Uranus
Synthesized Beam FWHMa	$2\buildrel{\prime\prime}\over{.} 2\times 1\buildrel{\prime\prime}\over{.} 6$, p.a. 35°	$2\buildrel{\prime\prime}\over{.} 1\times 2\buildrel{\prime\prime}\over{.} 0$, p.a. $-16^\circ$
rms noise (continuum image)	1.3 mJy	0.9 mJy
rms noise (line images)	⋯	87 mJy

Note.

^aNatural weighting of the visibilities.

Download table as:  ASCIITypeset image

ALMA Band 4 observations of KH 15D were obtained on 2016 January 16 (ADS/JAO.ALMA#2015.1.01537.S.). The object was observed for a total of 11.2 minutes using forty-six 12 m antennas. The data set contains four spectral windows: spw0 (138.99–137.00 GHz), spw1 (140.92–138.94 GHz), spw2 (149.00–150.99 GHz), and spw3 (151.00–152.99 GHz), for a total usable continuum bandwidth of 7.5 GHz. Each spectral window contains 128 channels, each 15.625 MHz wide. Higher frequency (Band 7) observations were also approved but, unfortunately, not obtained due to higher priority demand.

The data were calibrated using the ALMA Science Pipeline. To create the Band 4 continuum image, all channels containing line emission were excised out from the calibrated visibilities and the remaining channels were then averaged together, using the CASA task split. The channels containing line emission were determined from sharp peaks in a plot of channel versus amplitude. To prevent the diffuse clouds in NGC 2264 from obscuring the emission from KH 15D, coverage from the shortest baselines (<40 kλ) was excluded from the data set during this step. The 40 kλ value was chosen to minimize the rms noise in the continuum image, providing the best balance between reducing the flux contribution from the extended cloud emission while not excising so much data that the noise increased due to missing visibilities. Note that a 40 kλ baseline corresponds to a cutoff of spatial scales greater than 52 (∼4000 au), and thus we do not anticipate any loss of sensitivity to extended KH 15D disk emission. This technique reduced the rms noise of the image from 220 to 46 μJy. Deconvolution was carried out using the CLEAN algorithm and applying natural weighting. The size of the synthesized beam is $1\buildrel{\prime\prime}\over{.} 4\times 1\buildrel{\prime\prime}\over{.} 8$.

In addition to the continuum image, channel maps were created for each of the two spectral features present in the visibilities (details on these spectral lines will be discussed in Section 3.3.2). The continuum emission was subtracted using the CASA task uvcontsub, and the CLEAN algorithm (natural weighting) was again used to deconvolve the channel maps from the dirty beam. The velocity resolution in both ALMA channel maps is coarse (68 km s⁻¹ in spectral window 0, 62 km s⁻¹ in spectral window 2) since the original science goal was to attempt to detect the source in the Band 4 continuum image.

In this new era for KH 15D, represented schematically in panel 4 of Figure 1, the entirety of Star B is visible above the disk at each apastron passage, while Star A remains occulted at all phases of the binary orbit. We, therefore, wished to measure for the first time the mean brightness and color of the fully unocculted Star B. To do this, we selected the flat portions of the light curve peaks in the 2016/2017 data, corresponding to light curve phases between 0.35 and 0.65 (see Figures 3 and 4). Zero phase corresponds to the midpoint of the eclipse. A similar sample for Star A was generated using data from 1998 July through 2005 May, when Star A was the lone visible component, selecting orbital phases between 0.31 and 0.70. The differences in phase between the two samples is due to the slightly different positions of the disk edges in relation to the projected orbit of the stars (leading edge for Star A, trailing edge for Star B). The resulting mean magnitude of each component in VRIJHK and the standard deviation of the mean is given in Table 3. These represent our best estimates of the mean unocculted apparent brightness of each component, uncorrected for any scattered light that may be present. The I value for Star B (14.20) is in excellent agreement with the prediction of Windemuth & Herbst (2014) of I = 14.19 ± 0.06 based on the brightness of Star A and archival photographic data (see below).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 3.  Brightness of KH 15D Components

	Apparent Mag			Extinction (mag)a	Absolute Magb
Band	Star A	Star B	Both Stars		Star A	Star B	Both Stars
V	16.039 ± 0.003	15.509 ± 0.009	⋯	0.233	6.756 ± 0.055	6.226 ± 0.056	⋯
R	15.257 ± 0.005	14.776 ± 0.008	⋯	0.175	5.916 ± 0.055	5.436 ± 0.055	⋯
I	14.489 ± 0.001	14.198 ± 0.009	13.57 ± 0.03c	0.111	5.085 ± 0.054	4.794 ± 0.055	4.17 ± 0.06
J	13.504 ± 0.017	13.285 ± 0.010	⋯	0.066	4.054 ± 0.057	3.835 ± 0.055	⋯
H	12.825 ± 0.015	12.570 ± 0.013	⋯	0.044	3.353 ± 0.056	3.099 ± 0.056	⋯
K	12.541 ± 0.014	12.421 ± 0.022	⋯	0.027	3.052 ± 0.056	2.932 ± 0.059	⋯

Notes.

^aR = 3.1, $E(B-V)=0.075$ ± 0.003, extinction law from Cardelli et al. (1989). ^bd = 800 ± 20 pc. ^cWinn et al. (2006).

Download table as:  ASCIITypeset image

Each star is likely to be variable at about the 0.1 mag level, as they are spotted, rotating T Tauri stars. In fact, the rotation period of Star A, 9.6 days, has been determined from its periodic variations seen during unocculted time intervals (Hamilton et al. 2005). This could account for most or all of the scatter seen in the magnitudes at maximum brightness of each star. We have searched for a rotation period for Star B, but the time interval that the star is definitely unocculted is currently still too short to be certain of any result. We do find a hint of a periodicity near 10 days. That would be interesting because, as Herbst & Moran (2005) have discussed, there is evidence that Star A has been pseudo-synchronized into its current rotation period through tidal interaction with Star B during periastron passages. If so, one might expect Star B to also be pseudo-synchronized with about the same period. No attempt has been made to identify the most "pristine" state of the star, i.e., when it is least spotted. We have simply averaged together all of the data during the unocculted phases, forming an "average" brightness of the photosphere, representing a typical spot coverage, whatever that might be.

Another factor influencing these results is the possible fine structure, or corrugation, of the disk edge (Herbst et al. 2002). This complicates a bit the task of identifying when, exactly, the star is unobscured. As may be seen from Figure 4, there has not been a completely steady progression of Star B toward its maximum brightness. In 2014/2015, the star appeared to reach full visibility during the last cycle of that season, but in 2015/2016, it had apparently become slightly occulted once again during the first cycle and progressively rose in brightness during that observing season, mimicking its 2014/2015 behavior. This could be due to the irregular nature of the disk edge (i.e., some hills and valleys) as discussed for Star A by Herbst et al. (2002). Or it could involve the spottedness of the star, especially if there is pseudo-synchronization into a 5:1 spin–orbital resonance, as suggested by Herbst & Moran (2005). Most likely, both effects contribute to the brightness pattern. The resolution of these issues must await the progression of the disk edge to the point where Star B is clearly unocculted for a sufficiently long period of time to study its rotational modulation.

In Figure 5, we plot the color of the system over the last two observing seasons against both its orbital phase and V magnitude. Near the system's bright phase (V mag < 15.6, phase between 0.35 and 0.65; see Figure 3), Star B's unobscured photosphere dominates the total magnitude, and the system's brightness and color is independent of phase. The new data are consistent with the overall picture described by Arulanantham et al. (2016), who included longer wavelength (Spitzer) data in their analysis. As discussed in that paper, the system first gets bluer as it fades (an effect most easily seen here in the V − I and V − J data), due to the wavelength dependence of the scattering phase function. When it is very faint, the system becomes extremely red, as is clearly seen in Figure 5. This indicates the likely presence of a third, low-temperature, light component in addition to Stars A and B. Windemuth & Herbst (2014) and Arulanantham et al. (2017) argue that this could be light from a brown dwarf or a luminous giant-planet component.

Zoom In Zoom Out Reset image size

The mean unocculted I magnitudes for Star A (14.49) and Star B (14.20) are plotted in Figure 1 and may be compared with the expectation based on the photographic photometry record that extends back to the time when both stars were visible together (Johnson & Winn 2004; Johnson et al. 2005; Maffei et al. 2005). Note that there may be a scattered light component to each of the stars' individual magnitudes, which could lead to overestimates of their luminosities. Until a more comprehensive scattering model is devised for the system, it is hard to know how to handle that potential contribution. Based on the individual stellar brightnesses determined here, ignoring any scattered light component, we would predict the brightness of the system to be I = 13.58 at times when both stars were unobscured. That value is shown in Figure 1. The last time both components were visible simultaneously and (photographic) photometry was obtained was during the period 1951–1991. Using the available archival photographic plate data from Johnson & Winn (2004), Johnson et al. (2005), and Maffei et al. (2005), Winn et al. (2006) determined the mean brightness of the combined starlight at this time to be I = 13.57 ± 0.03 mag, in excellent agreement with our result. There are severe difficulties in transforming from the photographic magnitude system to the CCD system in use today, and the authors of the archival work have noted that one might anticipate an uncertainty in the zero-point magnitude of ∼0.15 mag (Johnson & Winn 2004). Hence, the near perfect agreement between the summed CCD magnitudes and the historical photographic photometry may be somewhat fortuitous. However, it does indicate that the transformation to Cousins-I was successful and provides further support for our overall interpretation of the light curve of KH 15D.

To derive physical parameters for the stellar components, we need to adopt a distance and reddening for KH 15D. Table 4 provides a summary of all published distance measurements to NGC 2264, the cluster where KH 15D is a member, that we could find. In most of these studies, the distance was determined by spectrophotometric fitting of main-sequence stars; exceptions are mentioned in the table notes. Based on this set of studies, we adopt a distance of 800 ± 20 pc (the unweighted mean, omitting the discrepant Dzib et al. 2014 value). The mean foreground reddening to the cluster is only E(B − V) = 0.075 ± 0.003, as recently derived by Turner (2012). There is relatively little differential reddening within the cluster, and a comparison of the uncorrected colors and spectral types for Stars A and B suggests that, like its cluster neighbors, KH 15D experiences relatively little foreground extinction. We have adopted a reddening value of E(B − V) = 0.075 and a standard extinction law from Cardelli et al. (1989), which has a ratio of total to selective extinction of R = 3.1. The resulting extinction corrections and derived absolute VRIJHK magnitudes of both stars out of eclipse are presented in Table 3. The intrinsic V − I color was then used to estimate the effective temperature and bolometric correction, BC_V, of each star following Table 6 of Pecaut & Mamajek (2013). Stellar luminosities were calculated from the resulting bolometric magnitudes, and radii were determined from the effective temperatures and luminosities. The full set of derived properties for each star is given in Table 5. The errors listed in Table 5 result from the 1σ error on V − I, which corresponds to a spread in spectral type of approximately ±1 subclass. Masses and ages are, of course, model dependent, as will be discussed in what follows.

Table 4.  Distance and Reddening of NGC 2264

Reference	Distance (pc)	$E(B-V)$
Walker (1956)	800	0.082
Mendoza V. & Gomez (1980)	875	⋯
Kwon & Lee (1983)	700	0.07
Sagar & Joshi (1983)	800 ± 75	0.00–0.12
Perez et al. (1987)	950 ± 75	0.061 ± 0.021
Feldbrugge & van Genderen (1991)	700 ± 40	0.042 ± 0.009
Neri et al. (1993)	910 ± 50	⋯
Sung et al. (1997)	760 ± 85	0.071 ± 0.033
Park et al. (2000)	760 ± 35	0.066 ± 0.034
Baxter et al. (2009)a	913 ± 40	⋯
Sung & Bessell (2010)b	815 ± 95	⋯
Turner (2012)	777 ± 12	0.075 ± 0.003
Kamezaki et al. (2013)c	${738}_{-50}^{+57}$	⋯
Dzib et al. (2014)d	400	⋯
Gillen et al. (2014)e	756 ± 96	⋯

Notes.

^a $v\sin i$ modeling of NGC 2264 stars. Note that this technique has an additional systematic uncertainty of ±110 pc. ^bMean distance derived from SED fitting of cluster members. ^cParallax of water masers associated with NGC 2264. ^dCombination of the following three methods. (1) Fitting radial velocities and proper motions of sources in both NGC 2264 and (independently) HH 124 IRS, a radio cluster thought to be associated with NGC 2264, to an adopted galactic rotation curve. (2) Re-analyzing parallax measurements taken with Hipparcos of NGC 2264 members. (3) A convergent point analysis of NGC 2264 stars with known proper motions and radial velocities. ^eDerived from eclipsing binary CoRoT 223992193 in NGC 2264.

Download table as:  ASCIITypeset image

Table 5.  Derived Stellar Properties of KH 15D Components^a

	Star A	Star B
Spectral Type	K7b	K1c
(V − IC)0	1.671 ± 0.078	1.432 ± 0.078
Teff (K)	3970 ± 40	4140 ± 155
BCV	−1.14 ± 0.05	−0.95 ± 0.16
Luminosity (${L}_{\odot }$)	0.45 ± 0.03	0.61 ± 0.10
Radius (${R}_{\odot }$)	1.41 ± 0.05	1.52 ± 0.16
Mass (${M}_{\odot }$)	0.715 ± 0.005	${0.74}_{-0.04}^{+0.09}$

Notes.

^aResults for the adopted cluster parameters of d = 800 ± 20 pc and E(B − V) = 0.075 ± 0.003. Effective temperature and bolometric correction are derived from the intrinsic V − I_C color (Pecaut & Mamajek 2013). Solar parameters of ${M}_{\mathrm{bol},\odot }=4.74$ and ${T}_{\mathrm{eff},\odot }=5772.0$ were adopted for the luminosity and radius calculations (Prša et al. 2016). The model-dependent age is 2.5–5.0 Myr for both stars. ^bHamilton et al. (2001). ^cCapelo et al. (2012).

Download table as:  ASCIITypeset image

The main goal of the radio observations was to detect the inferred KH 15D accretion disk directly through its expected thermal continuum emission from millimeter-size grains present in the outer disk. At the outset, we knew this would be difficult because of the distance to the system and the confusing region, full of diffuse clouds, in which it is embedded. Scaling from disk surveys of nearby star-forming regions, we predicted a less than 50–50 chance that accretion disk emission would be detected for this WTTS. However, because of the importance of the object, discussed above, we felt it was worth the attempt.

3.2.1. SMA Continuum Results

3.2.2. ALMA Non-detection

ALMA observations in Band 4 (1959–2188 μm; 137–153 GHz) did not reveal a source at the location of KH 15D. We may use the non-detection to set an upper limit on the flux from the disk and a constraint on its mass. The final cleaned ALMA Band 4 image has an rms noise of $\sigma =46$ μJy, with no significant emission at the coordinates of KH 15D. We can therefore place a 3σ upper limit on the 2 mm flux density of the disk at 0.14 mJy. Note that we do not expect the disk to be spatially resolved, as the angular resolution of the data (15) corresponds to a linear scale of ∼1200 au at the distance of KH 15D. A non-detection is a little surprising given the 4σ detection with SMA, at lower sensitivity. Possible explanations for this apparent contradiction will be examined in Section 4.2.

If we assume that the emission is from an isothermal and optically thin portion of the disk, we can estimate an upper limit for the total mass of the disk using the following relation (Hildebrand 1983; Andrews & Williams 2005):

$$\begin{eqnarray}&&{M}_{d}=\displaystyle \frac{{d}^{2}{F}_{\nu }}{{\kappa }_{\nu }{B}_{\nu }({T}_{c})},\end{eqnarray} \tag{ 1 }$$

where d is the distance to the source, F_ν is the flux density at the observed frequency ν, ${\kappa }_{\nu }$ is the opacity, and ${B}_{\nu }({T}_{c})$ is the Planck function at a characteristic temperature T_c. The opacity was estimated to be ${\kappa }_{\nu }\approx 3(1\,\mathrm{mm}/\lambda )$ cm² g⁻¹. This approximation was derived by Beckwith & Sargent (1991) and assumes grains are between 0.5 and 20 μm in size. A characteristic temperature of 25 K (the approximate mean temperature in a survey of disks in Taurus–Auriga, conducted by Andrews & Williams 2005) was chosen, resulting in an upper mass limit of $0.0016\,{M}_{\odot }$, corresponding to 1.7 Jupiter masses. This mass estimate is sensitive to assumptions about the opacity and characteristic temperature. Further, we can only detect particles smaller in size than the wavelength of observation. Note, however, that this is still fairly large for a disk mass; the upper limit is larger than 55% of the known disk masses in the Andrews & Williams (2005) survey.

While the radio observations were designed to search for thermal emission from the accretion disk, which was not detected, they also permitted a search for molecular line emission, which was detected in both data sets. The SMA data revealed a CO J = 3–2 line source that is spatially coincident with the filamentary jet of shocked gas seen in the near-IR and very likely associated with the KH 15D system, if not the binary itself. The ALMA data revealed two spectral features that probably arise within the diffuse cloud emission in NGC 2264, but appear unrelated to KH 15D.

3.3.1. SMA Line Emission

The SMA velocity channel maps of the CO J = 3–2 (Figure 6) reveal a detailed image of the ambient cloud emission throughout the region. The clouds emit across kinematic local standard of rest (LSR) velocities of 10–21 km s⁻¹. Additionally, there is a linear structure present in the blueshifted half of the CO line, at 3–7 km s⁻¹. These maps suggest a spatial and kinematic association between the structure and KH 15D. They align perfectly with the filamentary jet of shocked H₂ emission discovered by Tokunaga et al. (2004) and span a velocity range close to the binary's systemic velocity (3.4 ± 1.3 km s⁻¹) inferred from the radial velocity measurements of Star A (Johnson et al. 2004) and a model of the system by Winn et al. (2006). Thus, we believe that this structure is associated with the KH 15D system, and that it could be a detection of the binary outflow first proposed by Hamilton et al. (2003) and discussed in more detail by Mundt et al. (2010).

Zoom In Zoom Out Reset image size

A close-up of the channels containing emission from the molecular outflow is shown in Figure 7. The northern half is slightly redshifted compared to the southern half. Together with the low radial velocities of the CO emission (jet speeds are typically 100–200 km s⁻¹; Mundt & Fried 1983), this implies that the jet is oriented nearly perpendicular to the line of sight, with the southern portion inclined slightly toward the direction of the observer. This reinforces the jet inclination angle of 84° proposed by Hamilton et al. (2003).

Zoom In Zoom Out Reset image size

In Figure 8, we overlap the image of the H₂ emission from Tokunaga et al. (2004) with the CO emission integrated across the channels exhibiting jet-like structure. The images are aligned using an average of the optical (6:41:10.31 +9:28:33.2; SIMBAD) and infrared (6:41:10.34 +9:28:33.5; CSI NGC 2264) positions of KH 15D (Cody et al. 2014). The shape of the emission aligns nearly perfectly along the axis of the outflow. Even the faint nebulosity visible in the southern border of the H₂ image, speculated by Tokunaga et al. (2004) to be a re-emergence of the hidden counterjet, aligns with a knot of CO emission. The similarities in the shape and knots of the emission between the two linear structures make the case that they are associated with one another and with KH 15D. There is, nonetheless, an intriguing small, but real, offset ($\sim 1^{\prime\prime}$ or 800 au) of the stellar image from the central line of the molecular emission, which we discuss further in Section 4.3.

Zoom In Zoom Out Reset image size

3.3.2. ALMA Line Emission

The rise of Star B sets the stage for an improvement to the system model (Winn et al. 2006) that should yield highly accurate stellar characteristics for T Tauri stars. There is only one other known low-mass, PMS binary in NGC 2264 with determined stellar parameters (CoRoT 223992193), which also shows evidence of an accretion disk (Gillen et al. 2014). Thus, the KH 15D system provides a valuable set of data points for improving PMS stellar models. As a prelude to the upcoming KH 15D model, we explore what the photometry alone predicts about the stellar components. The full model requires radial velocity information on Star B that has not yet been obtained, as well as a prescription for the scattered light component. Here we see what constraints the data in hand place on the stars independently of those considerations.

An output of the full model to come will be the stellar radii, which are highly constrained by the shape of the eclipse light curve. In combination with the stellar luminosities, these will yield effective temperatures of high precision. With the data in hand, we reverse this procedure and use derived effective temperatures to predict the stellar radii. For main-sequence stars, effective temperatures are usually estimated from their spectral type and/or color. Stars A and B have been reported to be K7 and K1 T Tauri stars, respectively, by Hamilton et al. (2001) and Capelo et al. (2012). However, these determinations are based on different spectral regions (optical and near-IR, respectively), and correlations between spectral type and effective temperatures for PMS stars are notoriously hard to quantify, given the complex, magnetically active photospheres of young stars (Cohen & Kuhi 1979; Hillenbrand 1997). Some authors have attempted to adapt the spectral type versus effective temperature relationship of main-sequence stars for young stars. For example, Pecaut & Mamajek (2013) proposed effective temperatures for 5–30 Myr stars of 4920 K and 3970 K for K1 and K7 stars, respectively. Other authors, such as Cohen & Kuhi (1979) and Hillenbrand (1997), simply adopt dwarf scales and apply larger uncertainties. Given the uncertainties inherent in these methods, as well as the lower precision data based on spectral types, we prefer to use the intrinsic color as our primary guide to effective temperatures. This choice is reflected in Table 5.

The stellar radii presented in Table 5 (1.41 ± 0.05 ${R}_{\odot }$ for Star A, 1.52 ± 0.16 ${R}_{\odot }$ for Star B) suggest that the two stars are about the same size. This appears sensible; the stars are presumably the same age, and the relatively small difference in brightness and color between the components, along with the Winn et al. (2006) model result, points to roughly equal masses (with Star B being slightly brighter, bluer, and therefore more massive). Using the Pecaut & Mamajek (2013) calibration, the V − I colors imply spectral types of K7 and K5, respectively, for Stars A and B. This result is consistent with the spectral type determination of K7 by Hamilton et al. (2001) based on optical spectra. It is marginally inconsistent with the spectral type of K1 determined by Capelo et al. (2012) based on near-IR spectra. This may simply reflect the difficulties of assigning consistent spectral types to T Tauri stars. If the Pecaut & Mamajek (2013) effective temperature and bolometric correction corresponding to spectral type K1 are used in place of those determined from the mean V − I color of Star B, the models give ages that differ by 36 Myr for the two stars in the system (using PMS evolutionary tracks from Marigo et al. 2017; see Figure 9). This reinforces our decision to use the stars' colors instead of their spectral types for the effective temperature estimates. It would be of interest to obtain additional optical spectra of Star B while it, alone, dominates the spectrum of KH 15D, for the purpose of obtaining a better spectral classification.

Zoom In Zoom Out Reset image size

Alternatively, the effective temperatures can be estimated by fitting a blackbody curve to the full VRIJHK spectral energy distribution of each star. As in the V − I estimate, we adopted the mean reddening to NGC 2264, E(B − V) = 0.075, from Turner (2012) and the extinction law from Cardelli et al. (1989), yielding ${A}_{{\rm{V}}}=0.23$. For each fit, the stellar angular diameter was allowed to vary as a function of the effective temperature and measured stellar magnitude, using the bolometric correction calibration from Hillenbrand (1997). This results in ${T}_{\mathrm{eff}}=4076\pm 18$ K and $R=1.19\pm 0.03\,{R}_{\odot }$ for Star A, and ${T}_{\mathrm{eff}}=4185\pm 63$ K and $R=1.38\pm 0.07\,{R}_{\odot }$ for Star B. The errors on the temperatures are representative of the errors on the fit and do not encompass the much larger systematic errors such as the assumption that each star is a perfect blackbody. The V − I color approach was favored because of the smaller systematic errors, which make it easier to properly quantify our (large) uncertainties. Additionally, the temperatures derived from the color provide a nicer fit to the Marigo et al. (2017) isochrones (see Figure 9). As noted previously, future dynamical models of the KH 15D system will derive the radii of the components to remarkable accuracy for PMS stars (within 10%). This presents the opportunity to work backwards to derive the effective temperatures, creating valuable data points to better characterize PMS stars in the future.

The effective temperatures of the binary components, derived using each of the three approaches described above, are plotted against their corresponding luminosities on an HR diagram in Figure 9. Also plotted are isochrones from Marigo et al. (2017), which build off the Padova-Trieste Stellar Evolution Code (PARSEC) previously edited by Bressan et al. (2012). A metallicity of Z = 0.0152, $Y=0.2485+1.78Z$ (the solar values adopted in Marigo et al. 2017) is used. The (unweighted) mean metallicity of NGC 2264 is $[\mathrm{Fe}/{\rm{H}}]\,=-0.08$ (Paunzen et al. 2010), and so we find it reasonable to adopt a solar metallicity. The age of the system is likely between 2.5 and 5.0 Myr, with 3.5 Myr providing the best simultaneous fit to both stars. We note that this representative age range does not reflect any systematic effects that may be present due to inaccuracies of the models or transformations from measured to physical variables. The Marigo et al. (2017) stellar models imply masses of $0.71\mbox{--}0.72\,{M}_{\odot }$ for Star A and $0.70\mbox{--}0.83\,{M}_{\odot }$ for Star B (Table 5). As described in more detail in the caption of Figure 9, the spectral type and blackbody SED fit methods do not produce stellar parameters that can be fit by a common age above the 3σ level.

Effective temperature scales are model dependent and have been shown to overestimate temperatures by about 200 K for PMS stars with masses less than $1\,{M}_{\odot }$ (Kraus et al. 2015; David et al. 2016). A temperature suppression of this magnitude can be explained by the exclusion of starspots and magnetic fields from the stellar models (Somers & Pinsonneault 2015; Feiden 2016). To address this concern, we shift the measured KH 15D stellar effective temperatures by 200 K toward cooler temperatures and again plot our results alongside the Marigo et al. (2017) isochrones in Figure 10. The shifted plot suggests a slightly younger age of 1.7–4.0 Myr for the system, as well as smaller masses ($0.57\mbox{--}0.65{M}_{\odot }$ for Star A and $0.55\mbox{--}0.73{M}_{\odot }$ for Star B). The discrepancy between the shifted and non-shifted data demonstrates the effect that systematic uncertainties in PMS star effective temperature scales have on our ability to accurately estimate physical parameters. The shifted values appear slightly more sensible in the context of the PARSEC models than their unshifted counterparts. Nonetheless, we proceeded in our analysis with the original measurements.

Zoom In Zoom Out Reset image size

Our values of radius and luminosity for Star A are entirely consistent with those determined by Hamilton et al. (2001) and adopted in the Winn et al. (2006) model. The Star B parameters also agree with those derived by Winn et al. (2006), although we note that the error bars in both studies are substantial. The effective temperature measurement continues to be of interest, since it is the most difficult parameter to determine just from photometric and spectroscopic data. The M_B/M_A ratio predicted by the Winn et al. model (M_B/M_A = 1.2 ± 0.1) is slightly larger than the mass ratio we determined here using the Pecaut & Mamajek (2013) effective temperature scale, but the discrepancy is only at the 2σ level. Use of the shifted temperature scale reduces the minor discrepancy to within 1σ. We expect the updated dynamical model to strongly constrain the radii and, thereby, the effective temperatures of these stars. With that anticipated precision, Stars A and B of the KH 15D system should be extremely useful for testing models of PMS evolution.

It was disappointing and slightly surprising not to detect thermal emission from the KH 15D accretion disk with ALMA, since the SMA data suggested a flux density that would be easily seen. We are not sure how to reconcile these facts but note the following. (1) The SMA detection did not reach the usual 5σ level that traditionally signifies a solid detection. (2) The frequencies of the SMA and ALMA observations were different, but it would require a physically unrealistic spectral index (around $\alpha =4.4$) to account for the ALMA non-detection on this basis. (3) Perhaps the SMA detection was of non-thermal emission associated with the source. Some WTTS are non-thermal radio emitters; synchrotron emission from events in their magnetospheres can cause detectable radio outbursts (Bieging & Cohen 1989). The stars in KH 15D probably pass close enough to one another at perihelion that their magnetospheres could interact, resulting in "pulsed" eruptions, such as those seen in DQ Tau (Tofflemire et al. 2017). To check on this possibility, we looked at the orbital phase distribution of the radio measurements. Unfortunately for this idea, it turns out that both the ALMA and 2011 SMA data were taken during an optically bright phase of the system near aphelion, whereas the 2010 SMA data were at an intermediate brightness phase. In short, we have no explanation for the apparent detection of KH 15D with SMA and clear lack of detection with ALMA. Obviously, it would be of interest to repeat the SMA measurements at some point in the future, with the benefit that it may soon be possible to detect changes in the filamentary CO outflow.

By far, the most exciting new result from the radio observations is the discovery in the SMA data of a CO outflow in the J = 3–2 transition that matches the near-IR molecular feature so nicely (Figure 8). This CO outflow extends 15'' to both the north and the south, corresponding to a full linear extent of approximately 24,000 au at the estimated distance of KH 15D. Although most CO outflows have wide opening angles, this structure is tightly collimated, with a collimation factor (length:width) of 13.6. Tight collimation is closely correlated with high jet velocities according to Bachiller (1996). Assuming a relatively modest constant jet speed of 100 km s⁻¹, we can estimate the timescale of the jet to be around 550 years. Such a short timescale suggests that we should be able to see changes in the structure of the outflow with repeated observations of the CO J = 3–2 line, perhaps over just the six-year period between the SMA observations and today. Spatially resolving these structural changes in the future would make it possible to estimate the translational velocity of the CO outflow, which combined with the radial velocities would give the true outflow speed.

The newly detected CO outflow is most certainly associated with the previously known H₂ feature. As can be seen in Figure 8, the morphology of the two filaments is precisely the same in almost every detail. The main difference is that the full CO filament can be seen in the radio, whereas the southern portion of the near-IR outflow is obscured by an intervening cloud. We have carefully measured the position angle of the well-delineated outflow and find it to be 176 ± 03.

Is this filament associated with KH 15D? The only reason to doubt that is the small, but real, offset between the optical/near-IR image of the star and the central line of the outflow, which is easily seen in Figure 8 and amounts to about 1'' or 800 au. By contrast, the following facts argue that the KH 15D system, if not the optically visible binary, is the source of the molecular filament. (1) The star is located almost at the central point along the CO filament, showing that the lack of a "counterjet" in the near-IR image is an extinction effect. (2) The radial velocity of the feature is consistent within the errors with the systemic velocity of KH 15D and may, in fact, prove to be an important element in ultimately establishing that systemic velocity for the full dynamical model to come. (3) The position angle of the CO outflow, 176 ± 03, agrees remarkably well with the position angle of the polarization of light from KH 15D during the eclipse measured by Agol et al. (2004), 17° ± 1°. Assuming the scattered light comes from forward scattering off of the CB ring (Silvia & Agol 2008), one would expect it to be polarized in a direction perpendicular to the ring/disk, i.e., in the same direction that one would expect an outflow to emerge. (4) The small range in velocities for what must be a high-speed outflow, given its tight collimation, indicates that the outflow must be nearly perpendicular to our line of sight, again the orientation one would expect given that the ring/disk must be nearly edge-on to cause the optical obscuration effects we observe.

In short, it is hard to believe that the spatial, kinematic, and orientation coincidences between KH 15D and the molecular filament are caused by some chance superposition of two unrelated physical entities. On the other hand, the small offset of the star from the central line of the CO element is definitely real and demands an explanation. Close inspection of Figure 8 or, more definitively, the original H₂ images displayed by Tokunaga et al. (2004) does show a wisp of nebulosity that emerges directly from the location of the binary star in a nearly exactly northward direction, then appears to bend and merge into the main molecular feature with a 176 position angle. This detailed feature of the outflow may be telling us something about the jet launching mechanism in this source (and in binary T Tauri stars in general). On the other hand, there is no hint of that feature in the radio map, so it could be an artifact. A speculative possibility is that the jet is actually associated with the putative third luminous object in the system mentioned in the introduction. In a study of spectroscopic binaries with periods between 12 and 30 days, Tokovinin et al. (2006) found 34% had a tertiary companion. However, we note that this occurrence rate was found to decline with increasing orbital period, and that KH 15D's 48.37 day period is over 50% longer than the 30 day upper limit for the sample. As discussed above, the timescales on which the outflow features are expected to evolve are sufficiently short that it is probably worth revisiting the source to obtain high spatial resolution modern images of the near-IR and radio molecular filaments.