Two Candidate KH 15D–like Systems from the Zwicky Transient Facility

Kearns-Herbst 15D (KH 15D, Kearns & Herbst 1998) represents a rare class of circumbinary disk system. Photometric observations that go all the way back to the 1950s show a complex light-curve behavior (Hamilton et al. 2005; Johnson et al. 2005; Maffei et al. 2005; Capelo et al. 2012; Aronow et al. 2018; García Soto et al. 2020), characterized by periodic (with a period of 48 days) dippings and decades-long dimming and rebrightening. Together with the spectroscopic observations of the central object (Johnson et al. 2004), studies have shown that the circumbinary disk is largely tilted relative to the central, highly eccentric binary (Chiang & Murray-Clay 2004; Winn et al. 2006; Silvia & Agol 2008; García Soto et al. 2020; Poon et al. 2021). Observations and more detailed studies of systems like KH 15D can provide useful constraints and insights into the evolution and dynamics of circumbinary disks (see Poon et al. 2021 and references therein).

Motivated by this, one of our co-authors, K. Bernhard, as an amateur astronomer, performed his search to identify similar systems in ongoing all-sky variability surveys. Specifically, his search was focused on the variable star catalog of Chen et al. (2020), which was based on data collected by the Zwicky Transient Facility (ZTF; Bellm et al. 2019; Masci et al. 2019). Several other authors of the present work were notified by Bernhard about the potentially KH 15D–like candidates later. This eventually led to further analysis and observations of the identified systems, as will be presented in the rest of this work.

The primary feature of the photometric light curve of KH 15D is its deep (∼4 mag in I) and long (∼50% of the binary period) occultation event on an otherwise photometrically relatively quiet star. ¹⁰ In addition, the light curve of KH 15D shows gradual changes in the baseline flux over a timescale of decades, due to the precession of the warped circumbinary disk. With only a relatively short (i.e., a few years) time baseline, this latter feature is not expected to be seen in the ZTF data.

The search starts from the ZTF variable star catalog of Chen et al. (2020), which contains about 780,000 periodic variables with classifications and another about 1,000,000 suspected variables with no classifications. Our prototype, KH 15D, belongs to the second catalog, probably because it is a rare type of variable and could not be classified as any of the variable types of Chen et al. (2020). Motivated by this, we focused the search for more KH 15D–like objects on the suspected variable catalog of Chen et al. (2020). ¹¹

Our final search follows closely the original procedure of K. Bernhard. First, variables with small photometric amplitudes (defined as <1.5 mag in both g and r bands) or short periods (i.e., <10 days in either g or r) are excluded. This substantially reduces the sample to 1041. Next, we perform the box least-squares (BLS; Kovács et al. 2002; Hartman & Bakos 2016) analysis on all survival variables to search for transit-like signals that closely resemble the occultation light curve of KH 15D. In the BLS analysis, we restrict the duration to be between 1 day and 90% of the searched period. The lower bound is effective in excluding the majority of eclipsing binaries that escaped from the classification of Chen et al. (2020), whereas the upper bound is useful in excluding variables with periodic outbursts (e.g., U Gem–type variables). Finally, all BLS analysis results are visually inspected to identify the probable candidates. The inspection primarily checks the unfolded and folded light curves and the level of variations outside of the best-fit "transit duration" window.

The systematic search reveals no more promising candidates other than KH 15D and the three candidates that were originally identified by K. Bernhard. After a closer look into the light curves, we accept the two most promising candidates for further analysis. These are assigned the names Bernhard-1 and Bernhard-2. The third candidate, ZTF J070412.91–112403.2, is rejected because the photometric scatterings in its out-of-occultation light curve are comparable (∼1 mag versus 1.5 mag in r) to the depth of the presumed occultation event.

We provide in Table 1 the relevant information of Bernhard-1 and Bernhard-2. The significance of the parallax measurement from the Gaia mission is at the level of 2σ–3σ, leading to unreliable distance inferences for both objects. Nevertheless, we find that Bernhard-1 may be associated with the Cygnus OB associations (Quintana & Wright 2021), whereas Bernhard-2 is not associated with any known star-forming groups or young star associations.

Table 1. Candidate Information

	Bernhard-1	Bernhard-2
ZTF identifier	J202055.22+381323.1	J071445.39-090152.1
R.A.J2000	20h20m5522	07h14m4539
Decl.J2000	$+38^\circ {13}^{{\prime} }23\buildrel{\prime\prime}\over{.} 1$	$-9^\circ {01}^{{\prime} }52\buildrel{\prime\prime}\over{.} 1$
Parallax a (mas)	0.59 ± 0.16	0.34 ± 0.16
P (days)	192.10 ± 0.02	63.358 ± 0.003
tin (MJD)	58227.81 ± 0.07	59155.24 ± 0.04
tout (MJD)	58339.29 ± 0.04	59181.18 ± 0.03
vin (R⋆/day)	0.134 ± 0.003	0.416 ± 0.010
vout (R⋆/day)	0.150 ± 0.002	0.291 ± 0.003

Note.

^a Taken from Gaia Early Data Release 3 (Gaia Collaboration et al. 2021). These values are not changed in Gaia DR3.

Download table as:  ASCIITypeset image

With the given coordinates we have also retrieved archival photometric observations. In particular, we found optical observations (g, r, i, z, and y) from the Panoramic Survey Telescope And Rapid Response System (Pan-STARRS, Chambers et al. 2016; Flewelling et al. 2020) between MJD = 55,300–56,900, which provided useful constraints on the occultation model (see Section 3). Additionally, we found near-infrared observations from 2MASS (J, H, and K_s ; Skrutskie et al. 2003, 2006) and WISE (W1–4; Wright et al. 2010, 2019) taken at multiple epochs. These observations extend the spectral energy distribution (SED) and reveal the existence of the circumstellar (or circumbinary) disks (see Section 3).

3.1. Occultation Model

In both candidate systems, the occultations last for ∼50% of the total period. This cannot be explained by a circumstellar disk occulting a binary star at an exterior orbit, as seen in EE Cephei, Aurigae, and a few other similar systems (e.g., Mikolajewski & Graczyk 1999; Dong et al. 2014; Zhou et al. 2018). Additionally, the occultation periods, 192 and 63 days, are too short to be explained by a precessing disk occulting one single central object. Therefore, we conclude that the circumbinary disk occulting a central binary is the most plausible explanation for both systems. They are therefore KH 15D–like.

Given the scarce photometric observations, we do not apply the detailed precessing disk model, as was developed for the case of KH 15D (e.g., Chiang & Murray-Clay 2004; Winn et al. 2006; Poon et al. 2021). Instead, a simplified one-sided, fully opaque screen model, as illustrated in Figure 1, is used to describe the time evolution of the occultation event. This model involves five primary parameters: t_in (t_out) and v_in (v_out) are the epoch and projected perpendicular velocity of the star at the middle of ingress (egress), respectively, and P is the period of the occultation (i.e., the inner binary).

Zoom In Zoom Out Reset image size

Specifically, we model the stellar brightness profile as

$$\begin{eqnarray}&&F(t)={F}_{1}\cdot f(t)+{F}_{2}.\end{eqnarray} \tag{ 1 }$$

Here, F₁ + F₂ and F₂ are the flux values of the system outside and inside the occultation, respectively. These linear parameters are derived analytically for any given set of model parameters (t_in, t_out, v_in, v_out, P) through the maximum likelihood method. The normalized light curve f(t) is given by

$$\begin{eqnarray}f(t)=\left\{\begin{array}{ll}1, & x\leqslant -1\\ \tfrac{1}{2}-\tfrac{1}{\pi }\left[x\sqrt{1-{x}^{2}}+\arcsin x\right] & ,-1\lt x\lt 1\\ 0, & x\geqslant 1\end{array}\right.,\end{eqnarray} \tag{ 2 }$$

where

$$\begin{eqnarray}x=\left\{\begin{array}{ll}{v}_{\mathrm{in}}(t-{t}_{\mathrm{in}}), & \mathrm{ingress}\\ {v}_{\mathrm{out}}({t}_{\mathrm{out}}-t), & \mathrm{egress}\end{array}\right..\end{eqnarray} \tag{ 3 }$$

Here we have assumed no limb-darkening effect for the stellar surface.

This one-sided screen model is applied to the ZTF data of both Bernhard-1 and Bernhard-2 objects. The emcee sampler from Foreman-Mackey et al. (2013) is used to perform the Markov Chain Monte Carlo (MCMC) analysis and derive the uncertainties on the model parameters. The results from this analysis are given in Table 1. We then apply the best-fit models to the archival Pan-STARRS data, which went back as long as ∼11 yr. For each target, we show the joint photometric light curves from ZTF and Pan-STARRS as well as the period-folded light curves of individual data sets. These are illustrated in Figures 2 and 3. The same models are also used to construct the SEDs when the targets are inside and outside of the disk occultation. The resulting SEDs are shown in Figure 4.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

3.2. Bernhard-1

3.3. Bernhard-2

Apart from a shorter period (63 days), Bernhard-2 shares several similar features to Bernhard-1: a flat light curve outside of occultation, asymmetric ingress/egress regions, and the lack of rebrightening at the middle of occultation. The ingress and egress regions in the period-folded ZTF light curve have large scatterings. Together with the fact that the best-fit model cannot explain the archival Pan-STARRS observations, it again may suggest the precession of the disk and/or the failure of the one-sided screen model.

We attempted to obtain follow-up, high-cadence photometric observations of Bernhard-2 when it was visible from the ground. A few nights of observations were obtained from the 32 inch telescope at the Post Observatory, the 16 inch telescope at the Remote Observatory Atacama Desert (ROAD, Hambsch 2012), and the 1 m telescopes of the Las Cumbres Observatory Global Telescope (LCOGT, Brown et al. 2013) in 2021 December, when the system was exiting the egress. In early 2022, more systematic observations were obtained on the 1 m LCOGT telescopes and captured a broader range of the egress region. The LCOGT observations are reduced by AstroImageJ (Collins et al. 2017) and shown in Figure 5, indicating a rather smooth and gradual transition from inside to outside of the occultation. Our follow-up photometric observations are available as data behind the figure.

Zoom In Zoom Out Reset image size

We observed Bernhard-2 with the High Resolution Echelle Spectrometer on the 10 m Keck I telescope (Keck/HIRES) at UT 07:42 on 2022 January 8, when the object was outside the occultation. At a seeing of 14, we integrated for 20 minutes without the iodine cell and adopted the master wavelength solution for that night. We achieved a signal-to-noise ratio (S/N) of ∼18 per resolution limit. To retrieve the spectroscopic parameters, we selected five spectrum segments in ranges of 5350–5500, 5500–5750, 5910–6180, 6320–6410, and 6670–6780 Å, respectively, and performed matches to the ELODIE library (Prugniel et al. 2007). This leads to a surface effective temperature T_eff = 4865 ± 82 K and a surface gravity of $\mathrm{log}g=4.37\pm 0.04$, and the object is found to be a K1.3 ± 0.5 type pre-main-sequence star (Pecaut & Mamajek 2013). This spectroscopic classification is consistent with the broadband colors of Bernhard-2 outside the occultation (see Table 2), assuming a small amount of extinction (Fang et al. 2017). Unlike the spectrum of KH 15D (e.g., Fang et al. 2019), the spectrum of Bernhard-2 contains no emission lines. This confirms the result from photometric observations that the object has no accretion activity and thus is probably an older star than KH 15D.

Table 2. SED Information of Both Candidate Objects

Filter	Bernhard-1		Bernhard-2
	In	Out	In	Out
g	21.28(6)	19.09(6)	20.46(8)	18.01(2)
r	19.42(2)	17.30(2)	19.14(4)	17.02(1)
i	18.301(7)	16.770(7)	18.36(2)	16.514(4)
z	17.67(1)	16.20(1)	17.80(2)	16.261(6)
y	17.27(2)	15.806(6)	17.658(3)	16.09(1)
J	15.64(6)	⋯	⋯	14.92(4)
H	14.80(7)	⋯	⋯	14.27(4)
Ks	14.5(1)	⋯	⋯	14.07(7)
W1	⋯	13.21(5)	15.10(12)	14.15(6)
W2	⋯	12.85(4)	14.21(15)	13.60(10)
W3	⋯	9.29(5)	⋯	10.56(25)
W4	⋯	7.4(1)	⋯	7.8(5)

Note. For each object, magnitudes inside and outside the occultation are given

Download table as:  ASCIITypeset image

The WISE W3- and W4-band observations had relatively low S/N values, as shown in the right panel of Figure 4. Nevertheless, the combined SED indicates the existence of a cold disk component, which is consistent with the expectation that the system is surrounded by a circumbinary disk.

This paper presents the discovery and analysis of two systems that show disk-occultation signatures. These systems, named Bernhard-1 and Bernhard-2, show periodic dimmings that have large amplitudes (>1.5 mag in both g and r) and relatively long durations (≳50% of the identified periods of 192 and 63 days, respectively). Both inside and outside of the occultations, the light curves appear to be fairly flat with no clear signs of active accretion. These features are best explained by the tilted disk occulting a central binary, similar to the famous case of KH 15D (e.g., Poon et al. 2021). SEDs of both systems indeed suggest the existence of cold disk components.

A one-sided, fully opaque screen model is used to fit the photometric observations from ZTF and Pan-STARRS, which have a time span of up to ∼11 yr. Although this model cannot explain the details of the photometric light curve, especially the portions of ingress and egress, it can reasonably reproduce the key features of the observed data. According to this model, the circumbinary disks in both systems may have been gradually precessing on a timescale of ∼10 yr. The current observations are too sparse, and the orbital properties of the central binaries remain unknown, preventing us from applying a more complicated occultation model. Both photometric and spectroscopic observations are encouraged in order to reveal the true nature of these systems.

Our findings represent a valuable addition to the rare class of KH 15D–like systems. Prior to this work, the only known objects in this class other than KH 15D were WL4 and YLW 16A in the ρ Oph star-forming region, both of which were identified in the near-infrared and highly extincted in the optical (Plavchan et al. 2008, 2013). Taking the number of candidate young stellar objects from data release 3 of Gaia (Marton et al. 2022), which is comparable to ZTF in terms of the magnitude limit, we find ∼3/80,000 for the rate of KH 15D–like objects in the optical survey down to ∼20 mag. This very rough estimate only applies to stellar binaries with relatively long (>10 day) periods, and it is possible that young binaries with shorter periods may also show disk-occultation signatures. We leave the more complete search for such systems to future works.

We thank Richard Post for contributing to the photometric observations of Bernhard-2 and Greg Herczeg for useful discussion. We also thank the anonymous referee for comments and suggestions on the manuscript. W. Zhu, W. Zang, and T.G. acknowledge support from the National Natural Science Foundation of China (grant No. 12133005). W. Zhu and S.D. are supported by the science research grant from the China Manned Space Project No. CMS-CSST-2021-B12. W. Zhu is also supported by the National Science Foundation of China (grant No. 12173021). J.Z. was supported by the Natural Sciences and Engineering Research Council of Canada (NSERC) under the funding reference # CITA 490888-16. This research uses data obtained through the Telescope Access Program (TAP), which has been funded by the TAP member institutes. Based on observations obtained with the Samuel Oschin Telescope 48-inch and the 60-inch Telescope at the Palomar Observatory as part of the Zwicky Transient Facility project. ZTF is supported by the National Science Foundation under grant Nos. AST-1440341 and AST-2034437 and a collaboration including current partners Caltech, IPAC, the Weizmann Institute for Science, the Oskar Klein Center at Stockholm University, the University of Maryland, Deutsches Elektronen-Synchrotron and Humboldt University, the TANGO Consortium of Taiwan, the University of Wisconsin at Milwaukee, Trinity College Dublin, Lawrence Livermore National Laboratories, IN2P3, University of Warwick, Ruhr University Bochum, Northwestern University and former partners the University of Washington, Los Alamos National Laboratories, and Lawrence Berkeley National Laboratories. Operations are conducted by COO, IPAC, and UW. This publication makes use of data products from the Two Micron All Sky Survey, which is a joint project of the University of Massachusetts and the Infrared Processing and Analysis Center/California Institute of Technology, funded by the National Aeronautics and Space Administration and the National Science Foundation. This publication makes use of data products from the Wide-field Infrared Survey Explorer, which is a joint project of the University of California, Los Angeles, and the Jet Propulsion Laboratory/California Institute of Technology, funded by the National Aeronautics and Space Administration. Some of the data presented herein were obtained at the W. M. Keck Observatory, which is operated as a scientific partnership among the California Institute of Technology, the University of California, and the National Aeronautics and Space Administration. The Observatory was made possible by the generous financial support of the W. M. Keck Foundation.