A New Class of Roche Lobe–filling Hot Subdwarf Binaries

Subdwarf B stars (sdBs) are stars of spectral type B with luminosities below the main sequence. The formation mechanism and evolution of sdBs are still debated, although most sdBs are likely helium (He)–burning stars with masses ≈0.5 M_⊙ and thin hydrogen envelopes (Heber 1986, 2009, 2016). A large fraction are found in binary systems (Maxted et al. 2001; Napiwotzki et al. 2004) and the most compact ones have periods ≲1 hr (e.g., Vennes et al. 2012; Geier et al. 2013; Kupfer et al. 2017b, 2017a, 2020). Systems with orbital periods ≲2 h at the exit of the last common envelope will overflow their Roche lobe while the hot subdwarf is still burning helium (e.g., Savonije et al. 1986; Tutukov & Fedorova 1989; Tutukov & Yungelson 1990; Iben & Tutukov 1991; Yungelson 2008; Piersanti et al. 2014; Brooks et al. 2015).

So far only a few hot subdwarf binaries have been identified that will start mass transfer while the hot subdwarf is still burning helium in its core. The first discovered system is CD–30°11223, which has an orbital period of P_orb = 70.5 minutes and will start mass transfer in ≈40 Myr (Vennes et al. 2012; Geier et al. 2013).

Most recently Kupfer et al. (2020) discovered ZTF J2130+4420, the first system where the sdOB fills its Roche lobe and has started mass transfer to its white dwarf (WD) companion. In this Letter, we report the discovery of ZTF J205515.98+465106.5 (hereafter ZTF J2055+4651), the second confirmed member of Roche lobe–filling hot subdwarf stars in compact binaries. ZTF J2055+4651 was discovered in a period search using the conditional entropy period-finding algorithm (Graham et al. 2013) of all Pan-STARRS objects with colors of g − r < 0.2 and r − i < 0.2, and more than 50 ZTF measurements at the time, as described in Burdge et al. (2019) and Coughlin et al. (2020). Conditional entropy is an information-theory-based algorithm that seeks to find the period that produces the most ordered (lowest entropy) arrangement of data points in the phase-folded light curve while also accounting for the phase space coverage of the data. Independently, on 2019 July 19, Gabriel Murawski reported the object to the International Variable Star Index (VSX).¹⁸ Also, Sandoval et al. (2020) report the absence of X-ray emission from 3.7 ks Swift observations and found a limit for the X-ray flux of <1 × 10⁻¹³ erg s⁻¹ cm⁻² leading to a limit of <6 × 10³¹ erg s⁻¹ at a distance of 2170 pc (Bailer-Jones et al. 2018). We present the system properties, compare them to ZTF J2130+4420, and discuss possible implications for this new class of Roche lobe–filling hot subdwarf binaries.

As part of the Zwicky Transient Facility (ZTF), the Palomar 48-inch (P48) telescope images the sky every clear night. ZTF J2055+4651 was discovered as part of a dedicated high-cadence ZTF survey at low Galactic latitudes (Bellm et al. 2019a; Graham et al. 2019), which either observed one field or alternated between two adjacent fields continuously for ≈1.5–3 hours on two to three consecutive nights in the ZTF-r band (Bellm et al. 2019b). The high-cadence data were complemented by data from the ZTF public survey that was made available after data release 2 on 2019 December 11. Image processing of ZTF data is described in full detail in Masci et al. (2019). The ZTF light curve of ZTF J2055+4651 consists of 224 observations in ZTF-g and 496 observations in ZTF-r.

Additionally on 2019 September 5, ZTF J2055+4651 was observed with the Gran Telescopio Canarias (GTC) using HiPERCAM, which is a five-beam imager equipped with frame-transfer CCDs allowing the simultaneous acquisition of u_s, g_s, r_s, i_s, and z_s¹⁹ band images at a rate of up to 1000 frames per second (Dhillon et al. 2016, 2018). ZTF J2055+4651 was observed at a 3 s cadence with a dead time of 10 ms for 1130 frames with HiPERCAM. The run lasted 1 hr, covering a little more than one 56 minute binary orbit. The cadence in u_s and z_s alone was a factor of 3 slower than the other bands to optimize signal to noise. The data were reduced using the dedicated HiPERCAM pipeline,²⁰ including debiasing and flat-fielding. Differential photometry was also performed. Due to technical problems the i_s-band light curve could not be used for the analysis.

On 2019 July 30 a sequence of five optical spectra was obtained with the Shane telescope and the KAST spectrograph at Lick observatory using a low-resolution mode (R ≈ 1200). The spectra were taken consecutively and had an exposure time of 12 minutes. Data reduction was performed with the standard IRAF packages. Although 12 minutes covers 1/5 of the orbit we detected large radial velocity shifts between the spectra, confirming the binary nature.

Additionally, on 2019 September 24 and 2019 September 25 phase-resolved spectroscopy of ZTF J2055+4651 was obtained using the Gemini-North Telescope, using Gemini Multi-Object Spectrograph (GMOS) in long-slit mode with the B600 grating (R ≈ 1400). We obtained a total of 64 spectra with an exposure time of 3 minutes per spectrum. CuAr arc exposures were taken at the position of the target before and after each observing sequence to account for telescope flexure. The data were reduced using the IRAF package for GMOS.

As we followed the same procedure as described in Kupfer et al. (2020), we only briefly describe the major steps here. To derive the orbital period of the system, we use the ZTF light curve with its multi-month baseline in combination with the HiPERCAM light curve. The Gatspy (Vanderplas 2015; VanderPlas & Ivezić 2015) module that utilizes the Lomb–Scargle periodogram (Lomb 1976; Scargle 1982) was used for period finding, yielding P_orb = 56.34785(26) minutes. Using the parameter from the best light-curve fit (see Section 4) and only fitting the zero-point of the ephemeris (T_o) to the HiPERCAM light curves, we find an ephemeris of

$$\begin{eqnarray}&&{T}_{o}(\mathrm{BMJD})=58731.9639(90)+0.03913045(18)E,\end{eqnarray} \tag{ 1 }$$

where E corresponds to the epoch.

We find spectral features of hydrogen and neutral (He i) and ionized (He ii) helium lines typical for sdOB stars. By fitting the rest-wavelength corrected average Gemini spectrum with metal-free NLTE model spectra (Stroeer et al. 2007), we determine T_eff = 33,700 ± 1000 K, log(g) =5.54 ± 0.11, $\mathrm{log}y=\mathrm{log}[n(\mathrm{He})/n({\rm{H}})]=0.06\pm 0.05$, and v_rot sin i = 201 ± 30 km s⁻¹. Additionally, we measured atmospheric parameters using NLTE models computed with TLUSTY (Hubeny & Lanz 1995) that include typical hot subdwarf abundances of nitrogen, carbon, and oxygen (Naslim et al. 2013). We find fully consistent results with our metal-free NLTE models. Figure 1 shows the fit to the Gemini spectrum. The occurrence of both He i and He ii as well as the increased helium abundance $\mathrm{log}\,y\gt 0$ classifies the hot subdwarf as an He-sdOB star (see Heber 2016).

Zoom In Zoom Out Reset image size

We folded the individual Gemini spectra on the ephemeris shown in Equation (1) into 20 phase-bins and coadded individual spectra observed at the same binary phase. This leads to a signal to noise per phase-bin of ≈30. We used the FITSB2 routine (Napiwotzki et al. 2004; Geier et al. 2011) to then measure radial velocities. Assuming circular orbits, a sine curve was fitted to the folded radial velocity (RV) data points (Figure 2) excluding the data points around Phase 0.8–0.2. We find a velocity semi-amplitude K = 404.0 ± 11.0 km s⁻¹. Similar to ZTF J2130, the velocity deviates significantly from a pure sine curve around the phase when the sdOB is farthest away from the observer, which can be explained with the Rossiter–McLaughlin effect. The red curve in Figure 2 show the residuals predicted from the Rossiter–McLaughlin effect calculated from our best-fitting light-curve model.

Zoom In Zoom Out Reset image size

The LCURVE code was used to perform the light-curve analysis of the HiPERCAM u_s, g_s, r_s, and z_s light curves (Copperwheat et al. 2010). We assume a Roche lobe–filling sdOB star, a hot WD, and an irradiated disk around the WD. The passband specific beaming parameter B (${F}_{\lambda }\,={F}_{0,\lambda }[1-B\tfrac{{v}_{r}}{c}]$; see Bloemen et al. 2011) was calculated following the approximation from Loeb & Gaudi (2003). We use B = 1.7, 1.5, 1.37, and 1.25 for u_s, g_s, r_s, and z_s, respectively. The passband specific gravity-darkening (β) and limb-darkening (a₁, a₂, a₃, a₄) were fixed to the theoretical values from Claret & Bloemen (2011) for T_eff = 35,000 K and log(g) = 5.00 and solar metallicity. We used β = 0.54, a₁ = 0.93, a₂ = −0.98, a₃ = 0.85, and a₄ = −0.29 for u_s; β = 0.42, a₁ = 0.96, a₂ = −1.05, a₃ = 0.87, and a₄ = −0.28 for g_s; β = 0.41, a₁ = 0.90, a₂ = −1.11, a₃ = 0.94, and a₄ = −0.31 for r_s; and β = 0.39, a₁ = 0.83, a₂ = −1.20, a₃ = 1.09, and a₄ = −0.38 for z_s. We also tried to set the limb- and gravity-darkening parameters according to the tables of Claret et al. (2020) for hydrogen-rich atmospheres and found no significant difference in the final results. In addition, to account for any residual airmass effect we added a first-order polynomial. The best value of χ² for this model was around 3500 for 3012 points, which also includes a weak eclipse of the hot WD. Although the eclipse is weak (≤1%), the χ² for the non-eclipsing solution is 3900, which is statistically significantly worse compared to the solution with the weak eclipse. For the final model we normalized the errors in the data to account for the small additional residuals and obtain a reduced χ² ≈ 1. To determine the uncertainties in the parameters we combine LCURVE with emcee (Foreman-Mackey et al. 2013), an implementation of a Markov Chain Monte Carlo sampler that uses a number of parallel chains to explore the solution space. We use 256 chains and let them run until the chains stabilized to a solution, which took approximately 1500 generations. Figure 3 shows the best fit for each HiPERCAM band.

Zoom In Zoom Out Reset image size

Although ZTF J2055+4651 is a single-lined binary, we can still constrain the masses and radii of the two stars by combining the light-curve modeling with the spectroscopic fitting. Parameters derived in this way by a simultaneous fit to the HiPERCAM u_s, g_s, r_s, and z_s light curves are summarized in Table 1.

Table 1.  Measured and Derived Parameters for ZTF J2055+4651

Right ascensiona	R.A. (hr)	20:55:15.98
Declinationa	decl. (°)	+46:51:06.46
Magnitudeb	g (mag)	17.72 ± 0.04
Parallaxa	ϖ (mas)	0.432 ± 0.098
Distance	d (kpc)	${2.3}_{-0.4}^{+0.7}$
reddeningc	E(g − r) (mag)	0.46 ± 0.03
Absolute magnitude	Mg (mag)	4.4 ± 0.6
(reddening corrected)d
Proper motion (R.A.)a	μαcos(δ) (mas yr−1)	−3.33 ± 0.17
Proper motion (decl.)a	μδ (mas yr−1)	−5.160 ± 0.167
Atmospheric and Orbital Parameters of the sdOB from Spectroscopic Fits
Effective temperature	Teff (K)	33,700 ± 1000
Surface gravity	log(g)	5.54 ± 0.11
Helium abundance	log y	0.06 ± 0.05
Projected rotational velocity	vrot sin i (km s−1)	201 ± 30
RV semi-amplitude (sdOB)	K (km s−1)	404.0 ± 11.0
System velocity	γ (km s−1)	−20.0 ± 6.0
Binary mass function	fm (M⊙)	0.267 ± 0.022
Derived Parameters from Light-curve Fit Combined with Spectroscopic Results
Ephemeris zero-point	T0 (MBJD)	58731.9639(90)
Orbital period	Porb (min)	56.34785(26)
Mass ratio	$q=\tfrac{{M}_{\mathrm{sdOB}}}{{M}_{\mathrm{WD}}}$	0.60 ± 0.03
sdOB mass	MsdOB (M⊙)	0.41 ± 0.04
sdOB radius	RsdOB (R⊙)	0.17 ± 0.01
WD mass	MWD (M⊙)	0.68 ± 0.05
WD radius	RWD (M⊙)	0.0148 ± 0.0020
Orbital inclination	i (°)	83.4 ± 1.0
Separation	a (R⊙)	0.50 ± 0.01
Derived Parameters from Spectrophotometry and Parallaxe
Color excess	E(44–55)	0.559 ± 0.021 mag
Angular diameter	$\mathrm{log}({\rm{\Theta }}\,(\mathrm{rad}))$	−11.479 ± 0.011
sdOB radius	R = Θ/(2ϖ) (R⊙)	${0.17}_{-0.04}^{+0.05}$
sdOB mass	M = gR2/G (M⊙)	${0.37}_{-0.14}^{+0.29}$
Luminosity	L/L⊙ = (R/R⊙)2(Teff/ Teff,⊙)4	${34}_{-12}^{+24}$

Notes.

^aTaken from Gaia DR2, epoch 2015.5 (Gaia Collaboration et al. 2016, 2018). ^bTaken from PanSTARRS DR1 (Chambers et al. 2016). ^cTaken from 3D extinction maps (Green et al. 2019). ^dReddening corrected with 3D extinction maps from Green et al. (2019). ^eThe given uncertainties are single-parameter 1σ confidence intervals based on χ² statistics. A generic excess noise of 0.026 mag has been added in quadrature to all photometric measurements to achieve a reduced χ² of unity at the best fit.

Download table as:  ASCIITypeset image

Our solution converges on a mass ratio of q = M₁/M₂ = 0.60 ± 0.03, with individual masses of M_He−sdOB = 0.41 ± 0.04 M_⊙ and M_WD = 0.68 ± 0.05 M_⊙. The inclination is i = 83.4 ± 10 and the He-sdOB is Roche lobe–filling and has a volumetric-corrected radius of R_sdOB = 0.17 ± 0.01 R_⊙ (Table 1). We find that the He-sdOB would have a projected rotational velocity v_rot sin i =215 ± 6 km s⁻¹ if synchronized to the orbit, consistent with the observed v_rot sin i = 201 ± 30 km s⁻¹. We detect a weak eclipse of the hot WD that constrains its properties (Figure 3). We derive a blackbody temperature of 63,000 ± 10,000 K for the WD, with a radius of 0.0148 ± 0.0020 R_⊙, which is 10%–50% larger than a fully degenerate WD. This is consistent with predictions for hot WDs as shown in Romero et al. (2019), who predict 25% increased radius for a carbon–oxygen WD with T_eff = 60,000 K and a hydrogen layer of M = 10⁻⁴ M_⊙.

The spectral energy distribution (SED) in combination with the Gaia parallax allows for an independent estimate of the stellar parameters of the hot subdwarf. To this end, we first construct the observed SED from Gaia (Evans et al. 2018) and PanSTARRS photometry (Chambers et al. 2016).²¹ Using the atmospheric parameters derived by spectroscopy, a synthetic SED is computed with an updated version of Atlas12 (Kurucz 1996; Irrgang et al. 2018), which is then matched to the observed SED by fitting the angular diameter Θ as a distance scaling factor and E(44–55) as an interstellar reddening parameter (see Heber et al. 2018 for details). The interstellar extinction law by Fitzpatrick et al. (2019) is used. The results are summarized in Table 1 and visualized in Figure 4. As expected for a star in the Galactic disk, reddening is large (E(44–55) = 0.559 ± 0.021 mag).

Zoom In Zoom Out Reset image size

In combination with the Gaia parallax ϖ and the atmospheric parameters (T_eff and log(g)) the angular diameter allows us to determine the stellar parameters (R_sdOB, M_sdOB, and luminosity L_sdOB). Results are listed Table 1. The mass derived for the hot subdwarf agrees well with that derived from the light curve, although the uncertainty of the former is large because of the ≈25% parallax error and 0.11 dex $\mathrm{log}(g)$ error.

Using 3D extinction maps provided by Green et al. (2019)²² toward the direction of ZTF J2055+4651, we find a reddening of E(g − r) = 0.46 ± 0.03 at a distance of ${2.3}_{-0.4}^{+0.7}$ kpc, which is somewhat lower than the reddening derived by the SED analysis. Nonetheless, ZTF J2055+4651 is highly reddened and has a PanSTARRS color of g − r = 0.07 mag (Chambers et al. 2016). Such a color puts it in the regime of A- to F-type stars with T_eff ≈ 10,000 K, which is significantly cooler than typical hot subdwarfs. With the 3D map extinction maps corrected magnitude, we find an absolute magnitude of M_g = 4.4 ± 0.6 mag, consistent with a hot subdwarf star (Geier et al. 2019).

To put constraints on the population origin of ZTF J2055+4651 we calculated its kinematics. The proper motions of the system is taken from the Gaia data release 2 catalog (Gaia Collaboration et al. 2018; μ_αcos(δ) = −3.326 ±0.167 mas yr⁻¹, μ_δ = −5.160 ± 0.167 mas yr⁻¹). We followed the same approach as described in Kupfer et al. (2020) following Odenkirchen & Brosche (1992) and Pauli et al. (2006). We use the Galactic potential of Allen & Santillan (1991) as revised by Irrgang et al. (2013). The orbit was integrated from the present to 3 Gyr into the past. ZTF J2055+4651 moves between a distance of 8.6 to 11.6 kpc from the Galactic center with an eccentricity of 0.07 within a height of 50 pc of the Galactic equator. Exactly as for the other Roche lobe–filling hot subdwarf binary, ZTF J2130+4420, ZTF J2055+4651 is a member of the Galactic thin-disk population.

5.1. Evolution

Following the same prescription as described in full detail in Section 6.2 in Kupfer et al. (2020) we employed MESA version 12115 (Paxton et al. 2011, 2013, 2015, 2018, 2019) to construct evolutionary models for ZTF J2055+4651. To match the observed properties of ZTF J2055+4651, we start with a MESA model of a progenitor star with a mass of 2.9 M_⊙, which has a main-sequence lifetime of 350 Myr. Once the main-sequence star starts to evolve onto the red giant branch (RGB) He core burning starts under nondegenerate conditions. We remove most of the hydrogen-rich envelope, leaving a He core–burning sdB with T_eff ≈ 25,000 K, a core mass of 0.36 M_⊙, and a thin H/He envelope of 0.03 M_⊙. The envelope is defined by where the H mass fraction is above 0.01 M_⊙. This subdwarf model is then placed in a binary with a 0.65 M_⊙ WD with an orbital period of ≈140 min. Gravitational-wave radiation removes angular momentum from the binary and shortens the period. The period does not shrink far enough during He core burning, which lasts ≈350 Myr, that the sdB fills its Roche lobe (see the left panel in Figure 5).

Zoom In Zoom Out Reset image size

Once He core burning finishes the core contracts and hydrogen burning starts, the radius of the hot subdwarf expands beyond its Roche radius, and mass transfer starts at an orbital period close to the observed orbital period. Mass transfer will continue for ≈1 Myr at a rate exceeding 10⁻⁹ M_⊙ yr⁻¹ until hydrogen shell burning is finished and the star contracts to become a carbon–oxygen WD with a thick helium layer and a small residual layer of hydrogen. The high accretion rate will heat the accreting WD significantly (Townsley & Bildsten 2003). Models predict a T_eff ≈ 50,000 K for accretion rates of 10⁻⁹ M_⊙ yr⁻¹ (Burdge et al. 2019) consistent with the high blackbody temperature of the accretor observed in ZTF J2055+4651. Accretion onto the WD companion at this rate will cause unstable hydrogen ignition after ≈10⁻⁴ M_⊙ accumulates, leading to classical novae eruptions (Nomoto 1982; Nomoto et al. 2007; Wolf et al. 2013). This accretion rate predicts a recurrence time of order 10⁵ yr for a total of approximately 10 novae. Our binary model suggests that this system is within the first ≈10% of the 1 Myr accretion phase. At the current state the orbit will shrink with $\dot{P}\approx -6\times {10}^{-13}$ s s⁻¹, which will be detectable after a few years of monitoring. The right panel of Figure 5 shows the evolution of the donor through this mass transfer phase. After accretion has ceased, the orbit of the system will continue to shrink due to the radiation of gravitational waves and the system will merge in ≈30 Myr. Our models predict that there is a substantial He layer of ≈0.05 M_⊙ left in the former hot subdwarf and the total mass of the system is relatively high (M_total ≈ 1.1 M_⊙). Recent models predict that such a system explodes as a subluminous thermonuclear supernova (Perets et al. 2019; Zenati et al. 2019). If the system avoids a thermonuclear supernova it will merge and could evolve into a rapidly rotating single high-mass carbon–oxygen WD (Saio 2008; Clayton 2012; Schwab 2019).

5.2. The Accretion Disk

5.3. Implications for a Population of Roche Lobe–filling Hot Subdwarfs

ZTF J2055+4651 was discovered as a P_orb = 56.34785(26) minute variable with a light-curve shape similar to ZTF J2130+4420. Follow-up observations show that ZTF J2055+4651 is the second Roche lobe–filling hot subdwarf binary known today. High signal-to-noise ratio photometry obtained with HiPERCAM in combination with Gemini spectroscopy constrain the system parameters to a mass ratio q = M_sdOB/M_WD = 0.60 ± 0.03, a mass for the sdOB M_sdOB = 0.41 ± 0.04 M_⊙, and a WD companion mass M_WD = 0.68 ± 0.05 M_⊙. Additionally, the HiPERCAM light curve revealed a weak (≈1%) eclipse of the hot WD that allowed us to constrain its blackbody temperature to 63,000 ± 10,000 K.

We use the stellar evolution code MESA to compute the evolution of the system and find that the hot subdwarf formed in a common envelope phase at a period of ≈140 minutes and filled its Roche lobe close to the observed period when the hot subdwarf expanded during residual hydrogen shell burning. We predict that the binary will transfer mass for ≈1 Myr and we currently observe the binary in its first 10% of this phase. Our models predict a mass transfer rate of more than 10⁻⁹ M_⊙ yr⁻¹. Due to its high total mass the system will merge in ≈30 Myr and either explode as subluminous type Ia supernova (Perets et al. 2019; Zenati et al. 2019) or merge and form a massive single WD (Webbink 1984; Saio 2008; Clayton 2012; Schwab 2019).

ZTF J2055+4651 and ZTF J2130+4420 are the only known binaries with a Roche lobe–filling hot subdwarf donor. Both systems have similar properties and we see mass transfer during the short-lived hydrogen shell–burning phase compared to the He core–burning phase, which lasts significantly longer. The most likely explanation is a selection criteria bias against more typical hot subdwarfs with T_eff ≈ 25,000 K due to severe reddening at low Galactic latitudes. Future studies require reddening-corrected samples. We predict that a substantial population of He core–burning Roche lobe–filling hot subdwarf binaries can still be discovered at low Galactic latitudes in young stellar populations.

Based on observations obtained with the Samuel Oschin Telescope 48-inch at the Palomar Observatory as part of the Zwicky Transient Facility project. ZTF is supported by the National Science Foundation under grant No. AST-1440341 and a collaboration including Caltech, IPAC, the Weizmann Institute for Science, the Oskar Klein Center at Stockholm University, the University of Maryland, the University of Washington, Deutsches Elektronen-Synchrotron and Humboldt University, Los Alamos National Laboratories, the TANGO Consortium of Taiwan, the University of Wisconsin at Milwaukee, and Lawrence Berkeley National Laboratories. Operations are conducted by COO, IPAC, and UW.

Based on observations obtained at the international Gemini Observatory, proposal ID GN-2019B-FT-102, a program of NSF's OIR Lab, which is managed by the Association of Universities for Research in Astronomy (AURA) under a cooperative agreement with the National Science Foundation. on behalf of the Gemini Observatory partnership: the National Science Foundation (United States), National Research Council (Canada), Agencia Nacional de Investigación y Desarrollo (Chile), Ministerio de Ciencia, Tecnología e Innovación (Argentina), Ministério da Ciência, Tecnologia, Inovações e Comunicações (Brazil), and Korea Astronomy and Space Science Institute (Republic of Korea). The authors thank the staff at the Gemini-North observatory for performing the observations in service mode.

Some results presented in this Letter are based on observations made with the Shane telescope. Research at Lick Observatory is partially supported by a generous gift from Google.

Based on observations made with the Gran Telescopio Canarias (GTC), installed at the Spanish Observatorio del Roque de los Muchachos of the Instituto de Astrofísica de Canarias, in the island of La Palma.

This research was supported in part by the National Science Foundation through grant ACI-1663688, and at the KITP by grant PHY-1748958. This research benefited from interactions that were funded by the Gordon and Betty Moore Foundation through grant GBMF5076.

We acknowledge the use of the Center for Scientific Computing supported by the California NanoSystems Institute and the Materials Research Science and Engineering Center (MRSEC) at UC Santa Barbara through NSF DMR 1720256 and NSF CNS 1725797.

HiPERCAM and VSD are funded by the European Research Council under the European Union's Seventh Framework Programme (FP/2007–2013) under ERC-2013-ADG grant agreement No. 340040 (HiPERCAM).

T.R.M. was supported by a grant from the United Kingdom's Science and Technology Facilities Council. P.S. acknowledges support from NSF grant AST-1514737. E.S.P.'s research was funded in part by the Gordon and Betty Moore Foundation through grant GBMF5076. D.S. acknowledges support by the Deutsche Forschungsgemeinschaft through grant HE1356/70-1.

This work has made use of data from the European Space Agency (ESA) mission Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and Analysis Consortium (DPAC, https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the Gaia Multilateral Agreement.

This work benefited from a workshop held at DARK in 2019 July that was funded by the Danish National Research Foundation (DNRF132). We thank Josiah Schwab for his efforts in organizing this.

Facilities: PO:1.2 m (ZTF) - , Gemini:Gillett (GMOS) - , Shane (KAST) - , GTC (HiPERCAM). -

Software: Lpipe (Perley 2019), Gatspy (VanderPlas & Ivezić 2015; Vanderplas 2015), FITSB2 (Napiwotzki et al. 2004), LCURVE (Copperwheat et al. 2010), emcee (Foreman-Mackey et al. 2013), MESA (Paxton et al. 2011, 2013, 2015, 2018, 2019), Matplotlib (Hunter 2007), Astropy (Astropy Collaboration et al. 2013, 2018), Numpy (Oliphant 2015), ISIS (Houck & Denicola 2000), TLUSTY (Hubeny & Lanz 1995).