© The Authors 2026

1. Introduction

The loss of orbital energy through gravitational wave emission by low-mass compact binaries, such as binary neutron stars, neutron star-black hole systems, or double white dwarf binaries, causes the orbits of these systems to gradually tighten, leading many of them to conclude their evolution in a merger event (Nelemans et al. 2001; Blanchet 2002; Shen 2015). Depending on their masses, the merger of two white dwarfs can result in different outcomes, including a Type Ia supernova (Webbink 1984), a collapse into a neutron star (Nomoto & Iben 1985), the formation of a hot subdwarf (Heber 2009), or the formation of a R Coronae Borealis star (Clayton et al. 2007), but the most common outcome is the formation of another more massive white dwarf (Toonen et al. 2012; Kilic et al. 2024).

For mergers involving a black hole or neutron star, several studies (see, e.g., Rosswog 2007) have shown that the late-time transient behavior of the merger products, namely, the long-lasting high-energy (X-ray) emission that follows the short gamma-ray bursts produced in the process of merging, is likely due to the fallback accretion of gravitationally bound post-merger ejecta. In the case of GW170817, representing the first detected gravitational counterpart of a binary neutron star merger event (Abbott et al. 2017), several studies (see, e.g., Metzger & Fernández 2021; Ishizaki et al. 2021) have shown that this process neatly explains the year-scale X-ray excess. Double white dwarf mergers have also been theorized to be capable of hosting fallback accretion, with ejecta masses of the order of ∼10⁻³ M_⊙ (see, e.g., Guerrero et al. 2004; Lorén-Aguilar et al. 2009; Dan et al. 2014) being gravitationally bound and eventually re-accreted onto the newly formed remnant. Similar to the case of binary neutron star or neutron star-black hole mergers (Musolino et al. 2024), the fallback of the ejecta on a double white dwarf merger remnant is expected to not only drive the emission of soft X-rays (Rueda et al. 2019) but also optical and infrared emission (Rueda et al. 2019; Sousa et al. 2023). However, these predictions remain largely unconstrained, as no non-explosive double white dwarf merger events have been detected in transient surveys.

Potential white dwarf merger remnants may still act as a valuable source of information regarding the merging process. Models of white dwarf merger events tend to produce massive and rapidly rotating remnants as a consequence of mass and angular momentum conservation during the merging process (Schwab 2021). Additionally, due to the strong dynamos that arise during the merger, the remnants are predicted to be highly magnetized (Tout et al. 2008; García-Berro et al. 2012; Briggs et al. 2015). Thus, the combination of high mass, short spin period, and high magnetic field in a white dwarf are plausible fingerprints of a double white dwarf merger product. Examples of this are systems such as EUVE J0317−855 (Ferrario et al. 1997; Burleigh et al. 1999; Vennes et al. 2003), ZTF J1901+1458 (Caiazzo et al. 2021; Desai et al. 2025), and SDSS J2211+1136 (Kilic et al. 2021). Variations in the magnetic field strength and structure over the surface of these white dwarfs can be revealed by their spectroscopic and photometric variability (see, e.g., Vennes et al. 2003); therefore, time-resolved observations on the spin period are key in the study of these systems.

In this paper, we present the discovery and characterization of the white dwarf ZTF J200832.79+444939.67, hereafter ZTF J2008+4449, as a double white dwarf merger remnant showing the presence of circumstellar material. The existence of such material is supported by a host of evidence: the emission of non-photospheric X-rays, a significant increase in the short rotation period of the white dwarf (indicating that material is being ejected from the system, possibly in a propeller mechanism, and carrying away angular momentum), and the intriguing presence of a half ring of ionized gas trapped in the magnetosphere of the white dwarf. We pay particular attention to the analysis and modeling of the latter due to the unique nature of such a structure of circumstellar gas. Moreover, the absence of a stellar or brown dwarf companion (constrained by infrared observations) raises the question of the origin of this circumstellar material, and we propose a handful of possible answers: the fallback of post-merger ejecta, the disruption of a planetary body, or a magnetically driven wind from the white dwarf itself. ZTF J2008+4449 has also been identified as a potential member of the RSG 5 open cluster (Bouma et al. 2022; Prišegen & Faltová 2023; Miller et al. 2026; Yan et al. 2025), although we argue that it is likely an interloper. Given the unusual characteristics of this system, we believe it may be an invaluable probe into both the evolution and phenomenology of double white dwarf merger remnants. Additionally, the recent analysis of the merger remnant ZTF J1901+1458 by Desai et al. (2025) shows striking similarities to the white dwarf analyzed in this paper, as the object is a fast-spinning, highly magnetized white dwarf with similar X-ray emission properties to ZTF J2008+4449. This suggests their likely membership to the same new class of merger remnants.

The paper is structured as follows. In Section 2, we catalog and present all of the relevant available data on ZTF J2008+4449 and preliminarily describe the distinguishing features in these data. Section 3 is dedicated to a detailed account of our scientific methods and analysis of the data characterizing the white dwarf. In Section 4, we present the significance of our results and provide an extensive discussion of their physical implications for this system. Finally, in Section 5, we offer a summary of our findings and an outlook on future studies needed to fully understand this fascinating object.

2. Discovery and observations

ZTF J2008+4449 (Gaia DR3 2082008971824158720) was discovered as part of a search for periodic variability among Gaia white dwarfs (Gentile Fusillo et al. 2019, 2021) within the Zwicky Transient Facility archive (ZTF; Bellm et al. 2019; Graham et al. 2019; Masci et al. 2019; Dekany et al. 2020), which has already yielded several results, such as finding a large number of variable magnetic white dwarfs and compact binaries (Burdge et al. 2019a,b, 2020a,b, 2022a,b; Caiazzo et al. 2021, 2023). Previously, ZTF J2008+4449 was also identified as a large-amplitude variable within the Gaia Data Release 2 by Mowlavi et al. (2021).

In the Gaia color-magnitude diagram (CMD), ZTF J2008+4449 appears as a hot (≈30 000 K) and massive (≈1 M_⊙) isolated white dwarf (see Fig. 1). In our search, the object stood out because of its strong photometric variability in the optical with a very short period of ≈6.6 minutes. To understand the origin of the variability, we acquired time-resolved spectroscopic and photometric data necessary for the accurate characterization of this system. The array of available data is described in the following sections.

Fig. 1. Gaia CMD for white dwarfs within 100 pc from Earth and within the SDSS footprint (Kilic et al. 2020), where the x-axis depicts the difference between the Gaia BP and RP bands, and the y-axis the absolute magnitude in the Gaia G filter. Solid lines show theoretical cooling tracks (interpolated over the white dwarf’s mass with the WD_models package, Cheng 2024) for white dwarfs with masses between 0.6 M⊙ (top) and 1.28 M⊙ (bottom), equally spaced in mass; the atmosphere is assumed to be hydrogen-dominated (Holberg & Bergeron 2006) and the interior composition to be carbon-oxygen (Fontaine et al. 2001) for M < 1.1 M⊙ and oxygen-neon (Camisassa et al. 2019) for M > 1.1 M⊙. ZTF J2008+4449 is shown as a red star, and its location in the CMD reveals its high mass and high temperature. Reddening corrections were applied only to ZTF J2008+4449; as the objects in the sample are close to the Solar System, reddening is expected to be small. 1σ error bars are smaller than the size of the colored marker, and are omitted for the black background dots for clarity.

2.1. Optical light curve

The ZTF archive currently holds optical photometric data on ZTF J2008+4449 in the g and r bands throughout a roughly 6-year time span. Additionally, we obtained high-speed photometric data using the Caltech HIgh-speed Multi-color camERA (CHIMERA; Harding et al. 2016) on the 200 inch Hale telescope at Palomar Observatory. We obtained these observations on eight separate nights in the g and r filters. The full list of observations and exposures is shown in Table 1. The absolute local time of each exposure was measured to better than millisecond accuracy using a GPS receiver. We applied bias and flatfield corrections to the images and to perform aperture photometry we used the ULTRACAM pipeline (Dhillon et al. 2007). We used one reference star to construct differential light curves. In all observations, including spectroscopic and photometric data, we converted times into modified Barycentric Julian Date in the Barycentric Dynamical Time scale using the python library astropy.time (Astropy Collaboration 2022).

The phase-binned CHIMERA light curve, phase-folded using the period and period derivative for each epoch that we derive in Section 3.1, is shown in Fig. 2. The shape of the light curve is quite unusual compared to typical magnetic rotating white dwarfs, which usually show roughly sinusoidal light curves. Furthermore, the data show some aperiodic variations in the bright flat portion of the light curve (between phases 0.1 and 0.7), which results in increased scatter in the phase-folded light curve (we observed that, although it is brighter, the flat portion of the light curve is noisier than in the rest). The un-binned light curves are shown in the Appendix in Fig. C.1.

Fig. 2. Binned CHIMERA light curve with 3-second exposures, phase-folded at the correct period for each epoch as derived in Section 3.1 and normalized to the mean of the light curve in the g′-band (upper) and in the r′-band (lower).

2.2. Optical spectra

We acquired time-resolved optical spectra with the Low Resolution Imaging Spectrometer (LRIS; Oke et al. 1995) on the Keck I telescope on Maunakea, on six separate nights over a period of ∼five years (see Table 1). We reduced all LRIS observations using the publicly available LPIPE automated data reduction pipeline (Perley 2019). The full optical spectrum averaged over all exposures is shown in the top-left panel of Fig. 3 and shows broad and shallow absorption features that are typical of highly magnetized white dwarfs. In particular, the absorption feature at ≈8500 Å, with its peculiar step-like shape, is often observed in hydrogen-dominated white dwarfs with magnetic fields in excess of a few hundred MG (see, e.g., Caiazzo et al. 2021), as it corresponds to a blend of two strong Zeeman components of Hα, as illustrated in the lower section of the same plot in Fig. 3. Lastly, the spectrum shows weak emission in the Balmer lines Hα and Hβ. The emission lines are double-peaked, indicating the presence of Doppler shifts as high as ≈2000 km s⁻¹, but no Zeeman splitting due to a strong magnetic field. It is therefore likely that the narrow emission comes from a region far from the surface of the white dwarf, where the local magnetic field is weak.

Fig. 3. Top left: Phase-averaged Keck/LRIS spectrum of ZTF J2008+4449 and (below) the predicted positions of Hα and Hβ absorption lines as a function of magnetic field. The broad and shallow absorption lines, in particular at ≈8500 Å marked in blue, indicate a surface magnetic field B > 400 MG, while the Balmer emission lines suggest ongoing accretion. The red and blue dashed lines mark the maximal Doppler shifted wavelengths of Hα and have the same significance in all panels in this plot. Bottom left: Enlarged view of the phase-averaged Hα emission line, binned (unbinned) according to the instrument resolution shown as the solid black (faded gray) line. Bottom middle: Phase evolution of the trailed spectrum around Hα over the 6.6-minute period. The heat map reflects the mean-subtracted flux. Right: Trailed spectra of the Hα line, with spin phase progressing from bottom to top, as denoted by the numbers on the right hand side of the plot.

Trailed optical spectra were obtained by binning the data in ten phase bins over the optical modulation period, corrected for the epoch considered using the period and period derivative that we obtain in Section 3.1 (Fig. 3). The final spectrum of each phase bin was obtained by averaging all spectra in that bin and propagating their errors in quadrature. Upon inspecting the behavior of the Balmer emission lines in these phase-resolved spectra, it is clear that their shape in the averaged spectrum is caused by the periodic transition of one emission peak between maximally blue-shifted and maximally red-shifted wavelengths, as shown in the right and bottom-middle panels of Fig. 3 for the case of Hα. Interestingly, the variation in wavelength does not appear to be an S-curve, which is typically seen for Doppler shifts caused by rotational or orbital velocities. Instead, the transition between these states looks almost abrupt, with each one lasting roughly half of the period cycle (see the middle low panel of Fig. 3). The maximally blue-shifted emission corresponds to the minimum brightness in the ZTF and CHIMERA light curves (phase 0 in the right panel of Fig. 2). As a note to the reader, throughout this paper, all figures use a common zero phase.

2.3. UV spectra

We obtained ultraviolet (UV) spectra of ZTF J2008+4449 with the Cosmic Origins Spectrograph (COS; Green et al. 2011) and the Space Telescope Imaging Spectrograph (STIS; Woodgate et al. 1998) on the Hubble Space Telescope (HST, GO program 17720; Caiazzo et al. 2024). We used the low-resolution grating G140L on COS, with central wavelength 800 Å, which provides a spectral coverage over the range ≃900−1900 Å at a resolving power of ≈2000. For STIS, we used the G230L NUV MAMA grating, covering the range ≃1600−3000 Å. The COS observations spanned 5 orbits, for a total of ∼10 000 seconds, and were acquired in TIME-TAG mode, thus recording the arrival time of each individual photon. This event stream allowed us to construct the phase-binned spectra for 10 phase bins over the 6.6-minute optical modulation cycle.

In Fig. 4 we show the COS data, while in Fig. 5 we show the full UV spectrum and the STIS data. The phase-averaged COS spectrum is depicted in the upper panel of Fig. 4. The two strong emission features at ≈1210 and ≈1300 Å that exceed the upper limit of the plot are due to geocoronal airglow. The narrow absorption lines indicated in cyan are likely due to either interstellar or circumstellar carbon and silicon absorption. The phase-averaged UV spectrum of ZTF J2008+4449 shows an array of deep and broad absorption features clustered in the wavelength range between 1000−1350 Å. In blue in the figure, we indicate a feature at 1344 Å that could correspond to a hydrogen absorption line (a Zeeman Lyman-α component) in a strong magnetic field (in agreement with the Zeeman Hα component detected in the optical at 8500 Å). In the lower part of the same panel, we show the Lyman transitions as a function of magnetic field: the stationary Zeeman component at ≈1344 Å coincides with the observed absorption line. We do not have a clear identification for the other absorption lines, which we highlight in red. As we discuss in Section 3.2, these lines are likely absorption lines of metals. In the top panel of Fig. 5 we show the full COS+STIS phase-averaged spectrum. In the STIS range, the spectrum does not show any additional strong absorption features.

Time-binning the UV spectra once again reveals variability over the 6.6 min period. As illustrated in the other panels of Fig. 4 for COS, two types of variability stand out in the phase-resolved COS spectra. Firstly, the absorption line at 1344 Å (the possible Lyman-α component, highlighted in blue) varies in strength with phase, completely disappearing for half of the period. Secondly, we observed a strong sinusoidal modulation (with a peak-to-peak amplitude of 30%) in the continuum blueward of 1300 Å, with the brightest continuum phase corresponding to the strongest absorption in the 1344 Å line. A slight shift in phase at different wavelengths can also be seen; the minimum in the 1100−1200 Å range coincides with the minimum in the optical light curve (that we define as phase 0, see Fig. 2), while the minimum in the 1250−1280 Å is shifted later in phase by approximately 0.1 (we compare all light curves in all wavelengths in Fig. 16). The variability in STIS is more subtle. Sinusoidal variability can be seen in the 1800−2250 Å range, with an amplitude peak-to-peak of about 10% (see the middle right panel of Fig. 5).

Fig. 4. Upper panel: Phase-averaged COS spectrum in Fλ (solid black line, errors in gray). Strong absorption lines are highlighted. A possible Lyα component is in blue, and non-identified metal lines are in red (≈1092, 1108, 1174, 1190, 1193, 1234, 1253, 1277, 1328 Å; see Fig. 8). Narrow carbon and silicon absorption lines (possibly from the interstellar medium or circumstellar material) are highlighted in cyan. In the lower part, the Lyα component is compared to predicted line positions of Lyα and Lyβ as a function of magnetic field, where dashed lines indicate forbidden components (Ruder et al. 1994). Middle left: COS phase-resolved spectrum (five phase bins) over one period, zoomed in on the region with absorption lines. Each phase is shifted vertically by the same amount. Most of the absorption lines do not change significantly with phase, with the exception of the Lyα line, which disappears for half of the period. Middle right: COS light curves in two wavelength ranges: 1100−1200 Å (blue) and 1250−1280 Å (orange). Bottom two panels: Phase resolved COS spectra plotted over two periods. In the upper panel, the color bar indicates the flux in units of 10−15 erg s−1 cm−2 Å−1. In the lower panel, the phase-resolved spectra have been divided by the phase-averaged spectrum, so the color bar shows the variation with respect to the mean. Two sources of variability can be seen: a strong variation of the continuum between 1100 and 1300 Å, and a variation in the strength of the line at 1344 Å. The absorption is most pronounced when the continuum is at maximum brightness.

2.4. Photometry

We acquired new photometric data for ZTF J2008+4449 and compiled it together with pre-existing photometry. The photometric data are presented in the top panel of Fig. 9 and cataloged in Table 2. In the optical wavelength range, we used the photometry available in the PanSTARRS catalog (Chambers et al. 2016). As the errors in the Pan-STARRS photometry are lower than the photometric variation, we used an error of 0.1 magnitudes instead to take into account the error induced by variability. We also employed the Wide-field InfraRed Camera (WIRC) (Wilson et al. 2003) on the Hale Telescope at Palomar Observatory to obtain infrared photometry on two separate nights, in WIRC filters; J, H and Ks (shown in light blue). The near ultraviolet (NUV) photometry (dark blue) was acquired with the Neil Gehrels Swift Observatory (Swift) (Gehrels et al. 2004) through the Target of Opportunity proposals with numbers 14960 and 17885, target ID 13925. Finally, we converted the previously described COS UV spectrum to AB magnitudes to fill in the far ultraviolet (FUV) end of the photometric data (this is shown in cyan in Fig. 9).

2.5. XMM-Newton

The program that provided the HST observations (Caiazzo et al. 2024) was a joint program with XMM-Newton (Jansen et al. 2001). We requested 50 ks of non-continuous observation, with simultaneous Fast-Mode photometry with the Optical Monitor (OM; Mason et al. 2001). The observation was carried out in two separate visits (Obs IDs 0952990101 and 0952990201), each with a duration of approximately 32 ks. The source was clearly detected in all three cameras of the European Photon Imaging Camera (EPIC) in both visits: PN (Strüder et al. 2001), MOS 1 and MOS 2 (M1 and M2, Turner et al. 2001). The EPIC data were reduced using v21.0.0 of the Science Analysis Software (SAS) designed for the reduction of XMM-Newton data (Gabriel et al. 2004). We filtered the event files using the routine evselect and good time intervals (GTI) that were defined with standard defaults and thresholds on the count rate chosen by visual inspection of the total light curve for each observation; specifically, for Obs 01 we used 0.3, 0.1, and 0.15 cts/s for PN, M1, and M2, respectively, whilst for Obs 02 we used 0.26, 0.15, and 0.15 cts/s. This procedure resulted in total effective exposure times of 17.9, 21.4, and 21.9 ks for PN, M1, and M2, respectively, in Obs 01, and 19.3, 28.6, 28.9 ks for Obs 02. We extracted spectra for the sources and adjacent background regions using evselect, with background regions placed on the same chip and pixel columns as the source. The source and background circular apertures had radii of 20 arcsec and 80 arcsec, respectively. A source aperture of 20 arcsec should enclose 75−80% of source photons for an on-axis pointing in all three EPIC cameras¹. We list the total source and background counts in Table B.2, as well as the significance of detection in each band. In the PN camera, the source was detected with high significance (> 3σ) in all standard science bands below 4.5 keV (0.2−0.5, 0.5−2.0, and 2.0−4.5 keV) and in the broad band (0.2−10.0 keV), while M1 and M2 only detected the source below 2.0 keV and in the broad band. The X-ray spectral and timing analysis is presented in Section 3.5.

Light curves were obtained simultaneously to the X-ray observations with the Optical/UV Monitor (XMM-OM) in Fast Mode in the U and the UVM2 bands, in Obs 01 and 02, respectively. The signal-to-noise in the UVM2 band is very low and the light curve does not show any significant variability. On the other hand, the phase-folded U-band light curve, shown in Fig. 6, presents a similar shape as the optical light curves, albeit with a smaller amplitude and a small shift in phase (the orange dashed line indicates the minimum in the optical light curve, phase 0 in all plots).

Fig. 5. Data of STIS and COS UV for ZTF J2008+4449. Upper panel: Full phase-averaged UV spectrum in Fν (black solid line, errors in gray). Middle left: Phase-resolved spectrum (five phase bins) in the STIS band over one period. Each phase is shifted vertically by the same amount. There is very little variability with phase in the STIS band. Middle right: Comparison between the COS light curve in the 1100−1200 Å range (blue, same as Fig. 4) and the STIS light curve in the 1800−2250 Å range (red), which is where the STIS spectrum is most variable. The STIS light curve is also roughly sinusoidal, with a lower amplitude than the COS one, and slightly out of phase. Note that the error bars on the STIS light curve are smaller than the symbols used to plot it here. Bottom panel: STIS phase resolved spectra plotted over two periods and divided by the phase-averaged spectrum (the color bar shows the variation with respect to the mean). A weak variability can be observed, mostly in the 1800−2250 Å range.

Fig. 6. Binned XMM/OM light curve, shown with fewer (more) bins in blue (gray), phase-folded at the correct period for each epoch as derived in Section 3.1 and normalized to the mean of the light curve in the U-band.

3. Analysis

3.1. Period and period derivative

The Lomb-Scargle periodogram of the ZTF light curves reveals a strong optical modulation period at 6.55577 ± 0.00002 minutes and no additional significant periodicities (see Fig. 7). When phase-folded at the ZTF period, the CHIMERA light curve shows a significant shift in phase between each epoch, indicating a possible period variation. By measuring the shift in phase of the minimum in the light curve, we can measure the period derivative. For a small and constant period derivative (Ṗ ≪ P), the shift in the time of arrival of the minimum in the light curve after a number of cycles n can be expressed as

Fig. 7. Left: ZTF periodogram for ZTF J2008+4449 in the r-band. The main period and first harmonic (half of the period) are highlighted by vertical red dashed lines. The inset shows the complex structure of the main peak, due to a combination of the light curve windowing function and the fact that the period is changing significantly over the ZTF baseline. Right: O–C diagram for ZTF J2008+4449 for two high cadence ZTF observations (in green) and the 8 CHIMERA observations (in blue), all in the r band. On the y-axis, we plot the difference in seconds between the observed time of the minimum in the light curve and the expected time calculated assuming a constant period, while on the x-axis, we plot the elapsed time since our reference time T0 = 59401.76. For the second ZTF observation, the error in the x direction indicates the total duration of the ZTF data used in the fit (70 nights); for the other ones the duration is too small to show (three nights for the first ZTF observation, less than one night for the CHIMERA ones). The quadratic behavior indicates the presence of a significant period derivative.

$$\begin{array}{c}\hfill T_n-T_0\approx n{P}_0+\dot{P}{P}_0\frac{n^2}{2},\end{array}$$(1)

where T_n is the time of the minimum at the cycle n, T₀ is the time of the minimum at cycle n = 0, P₀ is the period at the time T₀, and Ṗ is the period derivative, which is assumed to be constant.

We took the eight CHIMERA observations in the r band and measured the time of the first minimum in each observation. This time was determined by period-folding the entirety of the observation and fitting a skewed inverted Gaussian profile to the dip-like feature in the light curve. As the change in period over the course of one individual observation is negligible, this method does not bias our measurement while it increases its precision. Fortunately, ZTF observed the white dwarf in high-cadence mode in two instances in 2018 and in 2019, allowing us to determine the phase shift in those times and therefore extend our baseline by three years. In particular, in November 2018, ZTF observed the target in the r band in “deep-drilling mode” (248 epochs in three nights, compared to the normal cadence of about one epoch every two nights, Kupfer et al. 2021), while between August and November 2019, ZTF observed the target in a high-cadence mode (310 epochs in 70 nights, an average of 4.4 epochs per night).

In the right panel of Fig. 7 we show the O–C (observed minus calculated) diagram for different observations in the r band of ZTF J2008+4449. On the y-axis, we plot the phase shift, T_n − T₀, in seconds, while on the x-axis, we plot the elapsed time since our reference time, T₀ = 59401.76, in days (in Barycentric Dynamical Time). If there was no period derivative, the O–C should follow a constant or linear trend, whilst for a (constant) period derivative, the O–C should follow a quadratic relation. The quadratic fit to the data provides an excellent match for the O–C, yielding a period derivative of Ṗ = (1.80 ± 0.09) × 10⁻¹² s/s, and a refined period estimate of P₀ = 6.55576688 ± 0.00000008 minutes at the epoch T₀. The positive period derivative provides clear evidence that the period is increasing. The CHIMERA optical light curves, phase-folded with the correct period at each epoch, are shown in Fig. 2.

3.2. Physical properties of the white dwarf

The broad and shallow absorption features in the optical spectrum of ZTF J2008+4449 (see left panel of Fig. 3) are characteristic of highly magnetized white dwarfs (> 100 MG). Similarly, the broad and shallow absorption lines in the ultraviolet spectrum are likely high-field Zeeman components of hydrogen and other metals. As we show in Figs. 3 and 4, the absorption lines at 8500 Å and at 1344 Å are likely Zeeman components of Hα and Lyα, respectively. If this is the case, they indicate an average magnetic field on the surface of the white dwarf in the range 500−600 MG (see Fig. 4). To explore the origin of the other lines, we calculated the transition energies in the Paschen-Back regime for a few metals (carbon, silicon and nitrogen) for which absorption lines have been observed in the spectra of dwarf novae in quiescence (see, e.g., Godon & Sion 2022) and for which lines can be expected at the effective temperature of the white dwarf. We explain our calculations in more detail in Appendix A and we show the different transitions in Fig. 8. Although our calculations are approximate (the perturbation theory that we employ might break down at such high field) and a more careful approach is needed to obtain confident identification in the lines, we observed that the π Zeeman components (the line transitions with no change in the m quantum number, Δm = 0) at high magnetic field cluster in the same wavelength range where we detect the absorption lines. The other Zeeman components, on the other hand, vary rapidly with field and can be strongly affected by magnetic broadening, producing extremely shallow and broad features that would be hard to detect.

Fig. 8. Top: Phase-averaged COS spectrum in Fλ (in black), with errors in light gray. The vertical bands in both panels mark the position of the absorption lines, with the blue one corresponding to the absorption line at 1344 Å, as in the top panel of Fig. 4; Bottom: Positions of a suite of C II–IV (black), Si II–III (purple), and N V (red) ionic transition components as a function of magnetic field strengths. Only transitions with a relatively large oscillator strength at 500 MG have been plotted. It can be seen that the π components of the Zeeman-split lines do not deviate far from the wavelength range where most of the absorption lines in the UV spectrum are located.

Another indication that the star is highly magnetized is in the shape of its spectral energy distribution (SED). As can be seen in Fig. 9, the SED shows a break between the optical and the ultraviolet that is often seen in highly magnetized white dwarfs (Green & Liebert 1981; Schmidt et al. 1986; Gänsicke et al. 2001; Desai et al. 2025), but that is hard to explain for non-magnetic stars. As we explain in more detail in our companion paper (Desai et al. 2025), this flattening of the UV spectrum in high-field white dwarfs is due to the effect of the strong field on the continuum opacities. The 550 MG magnetic hydrogen-dominated spectral model shown in black in Fig. 9, described below in this section, captures nicely the SED of the white dwarf across all wavelengths (with the exception of the unidentified absorption in the FUV). For comparison, we also show the non-magnetic, pure-hydrogen model (dot-dashed line Tremblay et al. 2011) at the same physical parameters, demonstrating that a model that does not account for the strong magnetic field will either be too shallow in the optical, or too steep in the UV.

Fig. 9. Top: SED of ZTF J2008+4449 plotting the AB magnitude against wavelength. The light blue, orange, and deep blue data points respectively show the IR (WIRC), optical (PanSTARRS) and NUV (Swift) photometry. The cyan line shows HST STIS and COS concatenated spectra spanning the NUV and FUV, with the COS spectrum truncated at the lowest wavelengths. The solid black and dashed gray lines show the best fitting pure Hydrogen magnetic and non-magnetic atmosphere models respectively. Lastly, a 1700 K (800 K) brown dwarf model is shown in purple (pink), by itself (faded line) and added to the best fitting magnetic white dwarf model (solid line). Bottom: Irradiated face temperatures for a series of Roche-lobe-filling brown dwarf and planetary companions of varying ages against their orbital separation from the white dwarf primary. The upper y-axis shows the orbital period given by each respective orbital separation. The four different tracks correspond to four different brown dwarf or planet ages (0.01, 0.1, 1 and 10 Gyr), marked next to each track. Brown dwarf models (≥13 MJup) are plotted in colors corresponding to their respective masses; all planetary models (≤13 MJup) are plotted in black. The size of each data point is representative of (but not linearly proportional to) the model radius. The golden (black) horizontal dashed line marks the minimal irradiation temperature for a brown dwarf (planet) among the evolutionary models utilized, i.e., ≈1700 K (≈900 K). The irradiated temperatures and corresponding orbital separations were calculated using the Sonora Bobcat brown dwarf evolutionary models grid (Marley et al. 2021).

Lastly, the continuum variability in the UV spectrum of ZTF J2008+4449 manifests as two broad humps: a strong hump in the FUV, blueward of 1300 Å (see bottom panel of Fig. 4) and a weaker one between ≈1800−2250 Å (see bottom panel Fig. 5). Their appearance as broad humps and the stark contrast between their large-amplitude modulation on the rotation period of the white dwarf and the rest of the weakly variable UV continuum (see Fig. 16), suggest that the emission in these ranges may be caused by cyclotron emission close to the surface of the white dwarf. As we do not have a good explanation for the rest of the photometric and spectroscopic variability of this white dwarf (see Section 4.1), we remain agnostic; however, we explore the possibility that indeed these humps are caused by cyclotron emission. For a magnetic field of ≈400−600 MG, the location of the observed humps would correspond to the energies of the second and first cyclotron harmonics, as shown in Fig. 10. Alternatively, the two ranges in wavelength could also correspond to the second and third harmonic for a lower field of ≈300 MG (in this second case, the brightening in the optical g-band at a phase of ∼0.5 could be due to emission in the first harmonic). Assuming that there is no contribution from cyclotron emission to the minimum of the light curve, we can extract a lower bound on the cyclotron luminosity from the difference between the minimum and the maximum fluxes of the strongest of the two humps; we find L_cyc ≳ 3 × 10³⁰ erg s⁻¹ = 0.001 L_⊙. If the cyclotron emission is due to inflowing material in free-fall close to the surface of the white dwarf, this luminosity would correspond to an accretion rate of Ṁ_cyc ≈ 2 × 10⁻¹³ M_⊙/yr.

Fig. 10. Central wavelengths of the first nine cyclotron harmonics (labeled 1 to 9) as a function of magnetic field (depicted by the black lines). The faded red and orange regions mark the ranges covered by the variable humps seen in STIS and COS, respectively. The dashed lines mark the centroid of each hump. For magnetic fields ranging from ≈400 − 600 MG (faded blue region), the positions of the two UV humps are consistent with the first and second cyclotron harmonics.

As we have just shown, non-magnetic atmosphere models are insufficient in constraining the physical properties of ZTF J2008+4449 from its optical and UV SED. We therefore developed state-of-the-art atmosphere models that include continuum and line opacity sets (e.g., Lamb & Sutherland 1974; Jordan 1992; Merani et al. 1995; Schimeczek & Wunner 2014) under a variable magnetic field strength. For simplicity, we employ a pure hydrogen composition, as we do not have a clear identification of the other metal lines. We describe the models in detail in our companion paper (Desai et al. 2025); the models employed here assume a simple dipolar structure with an average magnetic field across the surface of 550 MG and a field at the magnetic pole of 884 MG (this is quite arbitrary as we do not have constraints on the structure of the field on the surface of the star). We employed these models to fit the photometric and HST data and constrain the radius and temperature of the white dwarf, as well as the amount of interstellar reddening toward the star. Since we are using hydrogen-dominated models, we chose the spectral phase in which the likely hydrogen absorption is stronger in the UV spectrum (around phase 0.6) and we excluded from the fit wavelengths bluer than 1300 Å, where absorption by other elements dominates the spectrum and where the continuum variability is stronger. We find the effective temperature, stellar radius and interstellar reddening to be T_eff = 35 500 ± 300 K, R_WD = 4800 ± 300 km and E(B − V) = 0.042 ± 0.001, respectively, employing the photogeometric distance estimate from Bailer-Jones et al. (2021, 350 ± 20 pc). The best-fit model is shown as a black solid line in the upper plot of Fig. 9, while we show the corner plot for the fitting in Fig. D.1 in the Appendix. The derived reddening agrees perfectly with the value reported by the Bayestar19 dust map (Green et al. 2019). From these estimates, we employed the carbon-oxygen core evolutionary models of Bédard et al. (2020) to derive a mass of M_WD = 1.12 ± 0.03 M_⊙ and a cooling age of 60 ± 10 Myr (the cooling age is defined as the time since the birth of the object as a white dwarf, and therefore in this case it is the time elapsed since the merger). We also attempted to do the same by making use of the oxygen-neon core evolutionary models of Camisassa et al. (2019). However, the best-estimate temperature and radius from the SED fitting lie juts outside of the model grid, which is only calculated for masses above 1.1 M_⊙, and therefore we cannot derive a consistent estimate for the oxygen-neon composition. However, the radius of the lowest-mass model in the grid at the temperature of ZTF J2008+4449 is very close to our estimates, and we do not find a significant difference in the other parameters such as the cooling age.

Since it is useful for our photoionization model (see Section 4.5), we also derived the number of ionizing photons emitted from the white dwarf for the best-fitting model:

$$\begin{array}{c}\hfill L_{\mathrm{ion}}=4\pi R_{\mathrm{WD}}^2{\displaystyle {\int }_{{\nu }_0}^{\infty }}\frac{F_{\nu }}{h\nu }d\nu =(2.4\pm 0.1)\times {10}^{42}{\mathrm{s}}^{-1},\end{array}$$(2)

where F_ν is the flux from the white dwarf model as a function of frequency and ν₀ is the photoionization threshold for hydrogen (13.6 eV). All the parameters derived from the SED fitting are listed in Table 3. The quoted errors on the parameters are purely statistical in nature, and do not reflect the likely much larger systematic uncertainties introduced by the fact that we do not know the magnetic field strength and structure, and we just employ a fixed simple magnetic field structure in our models, while we know that the SED is significantly affected by the magnetic field strength and geometry. In future work, we plan to obtain an accurate identification of the metal absorption lines, which will yield a much more precise constraint on the field strength and geometry, and therefore we will be able to provide more accurate physical parameters for the white dwarf.

3.3. No trace of a brown dwarf companion

Although we rule out the presence of a companion on a 6.6 minute orbital period (see Section 4.1), the Balmer emission, the strange aperiodic variability in the optical light curve, and the possible cyclotron emission in the FUV could be signatures of accretion or, more generally, of ionized material in the relative vicinity of the white dwarf. We thus consider the possibility of a Roche-lobe-filling, sub-stellar companion, such as a brown dwarf, to be the mass donor in this system. To verify the presence of a companion brown dwarf, we obtained infrared photometry of the system (see Section 2.4 and Fig. 9). The white dwarf is sufficiently hot that a brown dwarf close enough to be filling its Roche lobe would be heated by irradiation. Additionally, as the brown dwarf would likely be tidally locked, we would expect strong variability at infrared wavelengths on the orbital period as the irradiated face of the brown dwarf comes in and out of our line of sight. To estimate the minimum effective temperature for the irradiated face of the brown dwarf companion, we employed a grid of brown dwarf evolutionary models (the Sonora Bobcat grid; Marley et al. 2021), from which we obtained the radius and intrinsic effective temperature of the brown dwarfs as a function of their age. We then assumed that the brown dwarf is overfilling its Roche lobe and therefore that it is orbiting the white dwarf at the distance for which its radius is equal to the Roche lobe (assuming a circular orbit). As the Sonora models extend into planetary regimes, we took the lower brown dwarf mass limit to be 13 Jupiter masses (M_Jup) and considered all models below this limit as planets (see Fig. 9); nonetheless, our convention of substellar and planetary regimes does not impact the underlying physics behind the following calculation of the companion irradiated face temperatures.

The effective temperature of a brown dwarf can be calculated from the total flux, F, leaving the atmosphere as

$$\begin{array}{c}\hfill T={\left(\frac{F}{\sigma }\right)}^{\frac{1}{4}},\end{array}$$(3)

where σ is the Stefan–Boltzmann constant. If the brown dwarf is irradiated, we can separate the total flux, F, into F = F_e + F_i, where F_i is the flux from the brown dwarf interior and F_e is the fraction of the white dwarf flux that is absorbed and re-emitted uniformly by the brown dwarf (Jermyn et al. 2017). The former of the two contributions is obtained from the intrinsic, non-irradiated effective temperature, T_BD, of the brown dwarfs in the Sonora Bobcat grid as F_i = σT_BD⁴. The latter contribution is related to the white dwarf luminosity by

$$\begin{array}{c}\hfill F_e=\frac{L_{\mathrm{WD}}}{16\pi a_{\mathrm{orbit}}^2}(1-A_B)=\frac{\sigma R_{\mathrm{WD}}^2{T}_{\mathrm{WD}}^4}{4a_{\mathrm{orbit}}^2}(1-A_B)\end{array}$$(4)

(Jermyn et al. 2017), where $L_{\mathrm{WD}}=4\pi R_{\mathrm{WD}}^2\sigma T_{\mathrm{WD}}^4$ is the white dwarf luminosity and a_orbit is the orbital separation between the white dwarf and brown dwarf. We also added the (1 − A_B) factor, where A_B is the Bond albedo (see, e.g., Karttunen et al. 2017) of the brown dwarf, to account for the fraction of the incoming radiation from the white dwarf that is scattered into space and does not contribute to the re-emission by the brown dwarf. Finally, we used a conservatively large value of 0.5 for the Bond albedo (implying a lower absorbed flux fraction and hence a lower irradiated temperature), as per the study of Marley et al. (1999) and in agreement with the study by Hernández Santisteban et al. (2016).

The resulting effective temperatures of the brown dwarf and planetary irradiated faces are plotted in the lower panel of Fig. 9, against both the orbital separation at which the companions fill their Roche lobes and the corresponding orbital period. The different tracks show distinct model ages, the marker sizes are proportional to their radii and the color bar shows the brown dwarf masses. Models corresponding to the planetary regime (< 13 M_Jup) are shown in black. It can be seen that the minimum brown dwarf irradiation temperature attainable is ≈1700 K, as marked by the golden horizontal line, for an orbital period of ≈160 minutes. We also mark in black the lowest mass planet in the Sonora grid (≈0.5 M_Jup), which would fill its Roche lobe at a period of ≈1180 minutes with a temperature of ≈900 K. The purple solid line in the upper panel of Fig. 9 shows the expected combined spectrum for the white dwarf with a 1700 K companion brown dwarf. We observed that the infrared photometry is inconsistent with such a hot brown dwarf. The coldest brown dwarf or planet that would be consistent with the WIRC IR photometry within three standard deviations has a temperature of ≈800 K (solid pink line): this is significantly below the lower limits on irradiation temperatures for any brown dwarf and for planets above 0.5 M_Jup. The WIRC J and H data were obtained over a span of two consecutive hours; furthermore, we have an additional two-hour long J and K observation taken a few days later and J does not show any hint of variability between the two epochs. Given the possible brown dwarf orbital periods, it is highly unlikely that we have only observed the cold, non-irradiated face of a potential brown dwarf. Hence, we confidently rule out the presence of a brown dwarf companion and, implicitly, that of accretion due to ongoing mass transfer as part of a cataclysmic variable evolution scenario. It is however still possible for a cold planet to be the mass donor in the system, given the orbital period at which it would be filling its Roche lobe.

3.4. Doppler tomography of the Balmer emission

From Fig. 3 we observed that the Hα emission presents Doppler shifts as high as ≈±2000 km s⁻¹, and that there is a narrow peak that jumps from blue-shifted to red-shifted values at ≈1700 km s⁻¹. These velocity features as a function of phase are directly related to the spatial geometry of the emitting gas (see, e.g., Marsh & Horne 1988); in order to visualize this geometry, we employed the technique of Doppler tomography (Marsh 2001). Doppler tomography, first introduced by Marsh & Horne (1988), is a technique designed to reconstruct the spatial geometry of emitting sources (e.g., circumstellar material) from the period-resolved variations of spectral features associated with them. Numerically, Doppler tomography is usually performed via a maximal entropy algorithm: individual values representing emission densities are assigned to each pixel among a 2D spatial canvas via the maximization of an appropriately defined associated entropy. This is done under the constraint that the resulting distribution of emission sources is a solution to (i.e., can produce) the phase-resolved variation of its emission features. The need to define an entropy associated to the source distribution arises from the spatial degeneracy of the solution and ensures that the final result is minimally assumptive, i.e., it is reconstructed making use of no information other than that contained by the variations of the spectral features. However, as these variations only encode information regarding the emitting sources’ radial velocity with respect to the observer, Doppler tomography only recovers the relative distribution of emission density in velocity Vsin(i) space (where i is the inclination of the white dwarf rotation axis with respect to the line of sight). This velocity-space solution must then be appropriately interpreted as one of multiple, doubly degenerate, real-space distributions of emitting sources.

We used the Doppler tomography software DOPTOMOG 2.0 (Kotze et al. 2015, 2016) to map the emission density from the shape of the Hα line variability. Fig. 11 depicts the obtained solution in the inclination-adjusted velocity space. The values of velocities depicted in the tomogram are as seen by the observer, with the white dwarf lying at the origin of the [V_x − V_y]sin(i) coordinate system (it is stationary with respect to the observer). The tomographic image in Fig. 11 is shown at rotation Phase 0; the amplitudes and phase evolution of the velocity vectors in the observer’s frame can be recovered by rotating the coordinate system over the obtained distribution, uniformly and clockwise over one variation cycle. (For a more detailed explanation of tomography, see, e.g., Kotze et al. 2015, 2016). The tomogram shows that the material is distributed in a half-circular structure centered around a velocity amplitude of ≈1700 km s⁻¹, dominated by a slightly asymmetric over-dense region.

Fig. 11. Doppler tomogram of the Hα emission feature. The minimally assumptive Hα emission density distribution, extrapolated from the emission line variation on the 6.6-minute white dwarf spin period, is shown in the v = Vsin(i) velocity space (i = inclination of the white dwarf rotation axis). The white dwarf is located at the origin at the coordinate system (v, θ). The evolution of the obtained source distribution with rotation phase can be recovered by rotating the coordinate system uniformly and clockwise over one variation cycle. The concentric dotted lines mark regions of constant velocity separated by 1000 km/s from each other, while the colors illustrate the density of emission. The emission distribution resembles a half ring, with velocities ranging between ≈1500 and ≈2500 km/s and peaking at a velocity of ≈1700 km/s.

Accounting for inclination, we note that a circular-arc distribution of emitting material in velocity space (i.e., spanning a single velocity amplitude) may straightforwardly (but not necessarily) correspond to a circular-arc distribution in real-space; the velocities of both circular Keplerian movement and rigid corotation only depend on the radial coordinate in the orbital plane. As we explain in more detail below, the velocities are too high to be explained by Keplerian orbits (the material would be close enough to the star to be forced to corotate with the strong magnetic field and we would be able to detect Zeeman splitting in the emission line). We therefore suggest that the velocities reflect rigid corotation of the emitting material in the magnetosphere of the white dwarf; as corotation velocities increase with distance from the star, in this case, the tomogram also provides an accurate representation of emission density in real-space, in a frame corotating with the white dwarf.

3.5. X-ray analysis

We conduct an analysis of the X-ray spectral data using the Bayesian X-ray Analysis package (BXA; Buchner 2016), which integrates the nested sampling algorithm ULTRANEST (Buchner 2016) with an X-ray spectral fitting platform. For the spectral modeling, we use XSPEC (Arnaud 1996) accessed through its Python interface, PyXSPEC. Our procedure involves jointly fitting the background-subtracted data from all three XMM EPIC cameras (PN, MOS1, and MOS2). To ensure reliable results in the low-count regime, we employ the Cash (C-) statistic for fitting (Cash 1979). Before fitting, we combined the spectra using the SAS command epicspeccombine and binned according to the optimal binning algorithm described by Kaastra & Bleeker (2016). We combined the PN spectra for both epochs, and all MOS (M1 and M2) for both observations, resulting in two optimally binned spectra which we fit simultaneously. We tested two different spectral models: (i) a simple power-law and (ii) a two-temperature, optically thin plasma model (APEC) with Solar abundance. Both models incorporate an absorption component to represent Galactic n_H absorption, as well as any intrinsic absorption within the system. Rather than rely on dust maps to estimate the Galactic absorption, we used the value of extinction derived from our optical and UV SED fitting (see Section 3.2). We converted this to hydrogen column density according to Equation (1) of Güver & Özel (2009). This has the advantage of implicitly accounting for both Galactic absorption and absorption intrinsic to the system, assuming that the X-ray and UV emission originates from the same place. From our SED analysis, we determine an extinction of E(B − V) = 0.042 ± 0.001, from which we find a hydrogen column density of n_H = (2.74 ± 0.69)×10²⁰ cm⁻². We adopted a Gaussian prior on the absorption, centered on this value with a standard deviation equal to the uncertainty.

For the two components of the optically thin plasma model, we find temperatures of kT₁ = 3 ± 1 keV and $k{T}_2=0.{23}_{-0.03}^{+0.04}$ keV. For the power law, we recover a photon index of Γ = 2.2 ± 0.1. The Bayes factor for the comparison of the two models is given by K = Z_pow/Z_apec, where Z denotes the Bayesian evidence (see, e.g., Buchner 2016) for each model, from which we find a value of K ≈ 2 which is sufficiently small that there is no Bayesian evidence that either model is favored. We checked whether the two best-fitting models provide a reasonable fit to the data using the XSPECgoodness routine, which uses Monte Carlo simulations to estimate the fraction of synthetic datasets that would yield a better fit statistic than the observed data. For the two-temperature model, from 1000 samples we find that 58.4% of the simulated spectra have lower, or better, fit statistics compared to our best-fitting model, implying that the model is a reasonable description of the observed data. For the power law model, we find 66.6% of the simulated spectra to have favorable fit statistics compared to our best-fitting model. We thus conclude that both the power law and two-temperature APEC models provide an adequate fit for the measured X-ray spectroscopic data. In both cases, the posterior on the absorption, n_H, closely resembles the prior informed by UV/optical extinction. For both models, we fit the combined spectra, but also the spectra from the two epochs separately. The top and bottom panels of Fig. 12 show the best-fitting model, convolved with the instrumental response, for the APEC and power law models, respectively. The posterior distribution of fit parameters are shown in the Appendix in Figs. D.2 and D.3. In black in the Appendix plots, we show the posterior distributions resulting from the fit of the combined spectra, whilst in blue and orange we show the fits to each epoch separately. From this analysis, we conclude that there is no evidence of variability in the X-ray spectral properties over the two-day period separating the two epochs.

Fig. 12. XMM spectral fit. The X-ray spectra are shown as a photon count rate as a function of energy (black data points) for the XMM PN detector (left column) and the combined M1 and M2 detectors (right column). The best fitting two-temperature APEC and power law emission models, convolved with the instrumental response, are shown in blue on the top and bottom rows respectively, with the 3-σ contours shown as shaded regions. In the case of the two-temperature APEC model, the two individual temperature components are shown in green and purple, along with their 3-σ contours. The residuals of the fit models, in both detectors, are shown on the bottom of their respective plots.

We performed a period search in the XMM-Newton EPIC data by phase-folding and binning the event times at a range of trial periods from 0.1−1 hours, with five phase bins. We combined both epochs of XMM data, making use of the barycenter-corrected event times, selecting events with energies in the range 0.2−2.0 keV to maximize the signal-to-noise. We fit each trial lightcurve with a constant, flat line and assess the goodness-of-fit using the χ²-statistic, such that large values of χ² indicate a lightcurve less-well described as being constant. Fig. 13 shows the results of this period search for EPIC-pn and MOS (M1 and M2 combined). We indicate the 3 and 5 standard deviations in χ² value with the horizontal dashed line. We marginally recover the white dwarf rotation period (shown in the vertical dotted line) with a confidence of 2σ, which is inadequate to confirm the white dwarf rotation period signal from the current X-ray observations. Fig. 14 shows the phase-folded and phase-binned X-ray events in the 0.2−2.0 keV band, folded on the rotation period. Although a hint of variability on the rotation period (mostly in the more sensitive PN) can be seen, a longer and deeper observation is needed to confirm it.

Fig. 13. Period search in XMM-Newton EPIC data using a χ2 test at a range of trial periods to search for a periodic signal. We show the resulting periodagram for the EPIC-pn and MOS data on the left and right, respectively. The vertical dotted line shows the period of the white dwarf at the epoch of the observation.

Fig. 14. Phase-folded light curve for the XMM-Newton EPIC-pn and MOS data. Barycenter-corrected events have been folded on the photometric rotation period, with the ephemeris defined for the epoch of the XMM observations, such that phase zero is consistent with all other plots in the paper (the minimum in the optical light curve).

3.6. RSG 5 cluster membership

ZTF J2008+4449 is identified by Miller et al. (2026) as a high-probability white dwarf member of the open cluster RSG 5. The identification is part of a search for high-mass cluster white dwarfs using the Hunt & Reffert (2023) cluster catalog and the Gaia EDR3 white dwarf catalog from Gentile Fusillo et al. (2021), with the goal of refining the initial-final mass relation for massive white dwarfs. However, they do not include ZTF J2008+4449 in their cluster-based white dwarf initial-final mass relation for two reasons: because of its likely nature as a merger remnant and because the white dwarf appears too old to be part of the cluster according to the cluster age derived in Hunt & Reffert (2023), and therefore they classify it as a likely interloper. It is worth mentioning that ZTF J2008+4449 was also independently identified as a candidate member of RSG 5 by Bouma et al. (2022) and Yan et al. (2025).

Although it is not relevant to the rest of the analysis presented in this paper, the possibility of cluster membership for the white dwarf would be interesting because it would provide a constraint on the total age of the system and therefore on the evolution of the binary before the merger. In order to explore further this possibility, we re-derived the age of RSG 5 with PARSEC 2.0 isochrones to better characterize its uncertainty (Bressan et al. 2012; Tang et al. 2014; Chen et al. 2014; Fu et al. 2018; Nguyen et al. 2022). We used stars with at least 50% probability of being cluster members in Hunt & Reffert (2023), and we applied the heuristic isochrone-fitting technique of Miller et al. (2026), which uses a visual fit to the main-sequence turnoff (MSTO) to select the best age and then derives 1σ uncertainties via a weighted χ² calculation that also heavily emphasizes the MSTO. Our fit yields a solar metallicity cluster with an age of ${45}_{-11}^{+8}$ Myr, in strong agreement with the pre-main sequence age of ${46}_{-7}^{+9}$ Myr reported by Bouma et al. (2022). The cluster CMD and best-fit isochrone are shown in Fig. 15.

Fig. 15. De-reddened Gaia DR3 CMD for RSG 5 (black) compared with the CMDs of three reference clusters: Alpha Persei (81 Myr, orange), IC 2391 (43 Myr, cyan), and NGC 2451A (50 Myr, red). The members of each cluster are those with at least a 50% membership probability in the Hunt & Reffert (2023) Milky Way cluster catalog, with clusters de-reddened using the catalog’s mean Av. Overlaid are PARSEC isochrones (Bressan et al. 2012; Tang et al. 2014; Chen et al. 2014; Fu et al. 2018; Nguyen et al. 2022) for our best-fit cluster age of RSG 5 (45 Myr) as well as the kinematic ages from Heyl et al. (2021) for the comparison clusters, with the colors matching those of the CMDs.

Since the MSTO of RSG 5 is sparsely populated, we also overlaid its CMD onto those of IC 2391, NGC 2451A, and Alpha Persei, each of which have estimated kinematic ages, independent of isochrone fitting, from Heyl et al. (2021). By visually comparing the lower main sequences, we find that RSG 5 falls squarely between IC 2391 (43 Myr) and NGC 2451A (50 Myr), and is clearly offset from Alpha Persei (81 Myr), as illustrated in Fig. 15. Given this agreement among ages from the MSTO, lower main sequence, and pre-main sequence, we conclude that the age estimate is solid and that there is no evidence for the true age of RSG 5 to substantially exceed ≈50 Myr.

Our estimate for the cooling age of ZTF J2008+4449, which is the time since the star became a white dwarf and therefore the time since the merger, is t_cool = 60 ± 10 Myr (see Section 3.2). This is of the order of the age of the cluster, or slightly higher. However, the cooling age is not the total age of the system, as it only indicates the time since the white dwarf was born as a merger remnant, and does not include the main sequence lifetime of the progenitor stars nor the binary inspiral time due to gravitational wave emission. Supposing ZTF J2008+4449 was the product of single star evolution, its mass (1.12 ± 0.03 M_⊙) would imply a progenitor mass of ≈6 M_⊙ (Miller et al. 2026; Cunningham et al. 2024), corresponding to a main-sequence lifetime of ≈70 Myr from PARSEC isochrones at solar metallicity, which one would have to add to the cooling age to obtain the total age of the system. As we think the star is the product of a double white dwarf merger, we need to consider the evolution of the binary as well. In some binary evolution scenarios, the merger time delay can be very small (Toonen et al. 2012; Yan et al. 2025), especially if the progenitor stars are massive; however, even in these scenarios, the time needed to create the merger remnant is at least a few tens of million of years, and there would not be enough time afterward for the white dwarf to cool down to its current effective temperature. We therefore conclude that the total age of ZTF J2008+4449 is likely higher than 100 Myr and therefore exceeds the age of RSG 5 significantly, indicating that the object is not a cluster member.

4. Discussion

The rapid spin of the white dwarf, together with the high magnetic field threading its surface, strongly suggest that ZTF J2008+4449 was born in a double white dwarf merger, and the absence of infrared excess exclude the possibility of mass transfer from a stellar or brown dwarf companion. However, the detection of Balmer emission and X-rays, together with the measured Ṗ, suggest the presence of circumstellar material. In the next sections, we discuss the origin of the 6.6 minute variability and the presence of circumstellar material trapped and ejected by the spinning magnetosphere, the observational constraints on its properties, and possible scenarios for its origin and for the nature of its interaction with the magnetosphere.

4.1. Variability

ZTF J2008+4449 shows variability on the 6.6 minute period at all wavelengths from the optical to the FUV, and possibly the X-rays (see Figs. 14 and 16). The optical light curves present a broad, dip-like feature that differs in depth and in shape between filters: it is deeper in the r-band (≈20%) and it becomes shallower at bluer wavelengths (≈17% in g and ≈10% in U), and the “ingress” is slower in the r-band than in the g-band, while the minimum appears shifted by about 0.1 in phase in the U-band. Additionally, the g-band displays a “brightening” peak at 0.5 phase away from the “dip”. The light curve shows aperiodic variability in the bright flat portion of the light curve (between phases 0.1 and 0.7), which results in enhanced scattering in the phase-folded light curve (see also Fig. C.1). We currently do not have a model to explain the unusual shape of the light curve, which is morphologically different from the sinusoidal light curves usually observed in magnetic white dwarfs. One possibility would be that some material trapped in the magnetosphere is partially eclipsing the white dwarf; in fact, the optical light curves bear some similarity to the light curves of an M dwarf that shows evidence of material trapped in the magnetosphere, for which this possibility has been suggested (Bouma & Jardine 2025). However, the minimum of the dip-like feature corresponds to the maximum blueshift in the Hα emission, indicating that at least the glowing disk (see Section 4.3) cannot be responsible for the dip. Material trapped in the magnetosphere further out could be the cause of the dips.

Fig. 16. Optical and UV light curves for ZTF J2008+4449, phased at the correct period for each epoch as derived in Section 3.1. Phase 0 corresponds to the minimum in the r-band. While for the CHIMERA and XMM/UVOT data the filter is indicated in the upper-left corner of each plot, for the STIS and COS data we extracted light curves by averaging the flux in the wavelength range indicated in the bottom right of each plot.

Another issue with the dip scenario is that, as we go to bluer wavelengths, the amplitude of the variability decreases and the dip-like feature disappears. The spectra in the COS and STIS wavelength ranges show very little continuum variation (see Figs. 4 and 5), with the exception of two humps that show quasi-sinusoidal variations in the 1100−1280 Å range and in the 1800−2250 Å range. In Section 3.2 we point out that this variation could be due to cyclotron emission close to the white dwarf surface as those energies would correspond to the first and second cyclotron harmonics in a magnetic field of about ≈500 MG (see Fig. 10).

Finally, the white dwarf shows strong variability in the Balmer emission lines, which we analyze in detail in the next few sections, and in one of the absorption lines in the UV. The optical absorption lines are extremely weak and only detectable in the phase-averaged spectrum (see Fig. 3). The signal-to-noise is too low in the phase-resolved spectra to constrain variability. In the ultraviolet, on the other hand, we detect stronger absorption lines in the 1000−1300 Å range (see Fig. 4). As we explain in Section 3.2 and in Appendix A, the absorption features most likely correspond to Zeeman components of absorption lines of metals (such as silicon, carbon and nitrogen), shifted by the strong magnetic field. The only exception is the line at 1344 Å, which, as we show in Fig. 4, is likely a Zeeman component of Lyman-α. None of the candidate metal lines show strong variability (once one subtracts the continuum variation) either in shape or in wavelength, while the Lyman-α component varies strongly in equivalent width, completely disappearing at Phase 0 (see Fig. 4). Strong variation in hydrogen absorption has been observed in other rapidly spinning white dwarfs: the emerging new class of Janus-like double-faced white dwarfs (see, e.g., Caiazzo et al. 2023; Koen et al. 2024; Cheng et al. 2025; Moss et al. 2024, 2025). In these systems, hydrogen and helium on the surface vary in anti-correlation, with the most extreme cases showing one of the two elements disappearing when the other is at maximum. It has been suggested that the variation is due to the presence of a magnetic field threading the surface of the white dwarf and a possible asymmetry in hydrogen content due to a merger origin (Bédard & Tremblay 2025); the variation that we see in the Lyman-α component of ZTF J2008+4449 could be due to a similar mechanism.

Due to the presence of circumstellar material and the strange shape in the optical light curve, we considered the possibility that the periodic variability on the 6.6 minute is connected to the orbital period of a binary. However, we can exclude this possibility. For a white dwarf in a compact binary, an orbital period of 6.6 minutes would be too short to host a stellar companion (the period minimum for cataclysmic variable systems is ≈80 minutes). If the 6.6 minute corresponded to an orbital period, it would require the white dwarf companion to be a compact object, such as a black hole, a neutron star or another white dwarf. We can exclude a neutron star or a black hole companion since it would not explain the Balmer emission, as the white dwarf would not be mass transferring (it would not fill its Roche lobe) and because we would see strong radial velocity shifts in the UV absorption lines. A companion white dwarf could be the donor in a 6.6-minute period, such as in the direct-impact accretor HM Cancri (Israel et al. 2003). However, the donor would have to be a low-mass white dwarf to fill its Roche lobe, and it would outshine the high-mass white dwarf accretor. Additionally, we would expect much lower velocity shifts in the Balmer emission (Roelofs et al. 2010). Finally, the modulation in the Lyman-α absorption in phase with 6.6 minutes strongly hints to a rotational rather than orbital modulation. The 6.6 minute periodic variability is therefore most likely associated with the rotation of the white dwarf, rather than the orbital period of a possible binary system.

4.2. Origin of the Balmer emission

The Balmer emission indicates the presence of hot ionized gas surrounding the white dwarf. The absence of Zeeman splitting, the peculiar jumping behavior between maximally red-shifted and maximally blue-shifted velocities of the peak emission, and the peak velocities can provide hints on the location and geometry of the emitting material. A key parameter of this discussion is the extent of the white dwarf magnetosphere, i.e., the region around the white dwarf in which the magnetic field is strong enough to force material to move along magnetic field lines. The extent of this region is indicated by the magnetospheric radius, which is usually taken as a fraction, ξ, of the Alfvén radius (Pringle & Rees 1972; Ghosh & Lamb 1979; Papitto & de Martino 2022):

$$\begin{array}{cc}\hfill R_{\mathrm{m}}& \approx \xi R_{\mathrm{A}}\approx \xi {\left(\frac{B_p^2{R}_{\mathrm{WD}}^6}{4\dot{M}\sqrt{2G{M}_{\mathrm{WD}}}}\right)}^{\frac{2}{7}}\hfill \\ \hfill & \approx 200R_{\mathrm{WD}}\left(\frac{\xi }{0.5}\right){\left(\frac{B_p}{500\mathrm{MG}}\right)}^{\frac{4}{7}}{\left(\frac{R_{\mathrm{WD}}}{4800\mathrm{km}}\right)}^{\frac{12}{7}}\hfill \\ \hfill & \times {\left(\frac{M_{\mathrm{WD}}}{1.12{\mathrm{M}}_{\odot }}\right)}^{-\frac{1}{7}}{\left(\frac{\dot{M}}{{10}^{16}\mathrm{g}/\mathrm{s}}\right)}^{-\frac{2}{7}},\hfill \end{array}$$(5)

where B_p is the magnetic field at the surface and Ṁ is the mass flux rate at the magnetospheric radius. Although we do not have a good determination of the magnetic field strength, we are confident that it is of the order of a few hundreds of MG or more (see Section 3.2). Employing our best estimate for the field from the absorption lines (≈500 MG), and the white dwarf parameters we derived in Section 3.2, we find that the magnetospheric radius is hundreds of white dwarf radii, even for a very high mass accretion rate.

Two possible scenarios can explain the velocities detected in Hα: the material could be orbiting the white dwarf at Keplerian velocities, or it could be rigidly corotating with the white dwarf, trapped by the strong magnetic field. Employing the white dwarf mass and radius estimated in Section 3.2, we find that a Keplerian velocity ($v_k=\sqrt{G{M}_{\mathrm{WD}}/d}$) of 1700 km s⁻¹ would indicate, in the case of a perfectly edge-on system, a distance from the white dwarf of d ≈ 11 R_WD (see the blue vertical line in Fig. 17). For an inclined system, the material would have to be orbiting even closer to the star (blue shaded region highlighted by the blue arrow in Fig. 17). However, at this distance, we are deep into the magnetosphere of the white dwarf and the magnetic field would be strong enough to force the ionized material into rigid-body corotation with the white dwarf (and therefore move slower than the Keplerian velocity at this distance). Moreover, we would be able to detect the Zeeman splitting in the emission line as the local magnetic field would be higher than 0.1 MG. We therefore exclude that the line-emitting material is orbiting the white dwarf in a Keplerian disk. If some material is present within the corotation radius, it could fall toward the white dwarf following the field lines, and therefore the line-of-sight velocity that we measure could be a mixture of corotation velocity and infall velocity. The same would apply if the material was moving outward in a wind, following the field lines. However, it would be hard to explain the precise jumping behavior in the line-of-sight velocity that we detect in both Hα and Hβ. In the next section, we show that a half ring of ionized gas trapped in corotation with the white dwarf can instead well reproduce the observed Hα profile and phase variation.

Fig. 17. Diagram showing the possible location of the line emitting material. On the x-axis, we plot the distance from the white dwarf (d) in white dwarf radii (RWD ≈ 4800 km). The black solid line indicates the strength of the magnetic field at the distance d, assuming a field of dipole strength at surface of ≈500 MG. We highlight with a dash-dotted vertical line where the field becomes weaker than one MG. The top x-axis marks the corotation velocity corresponding to each radial coordinate. The ratio of centrifugal to gravitational accelerations is unity at the Keplerian corotation radius (vertical dashed line). The vertical blue line indicates the distance at which the Keplerian velocity equals the observed velocity in Hα; the blue shaded region shows where the emitting material would lie if it was moving in Keplerian orbit, depending on inclination (with the vertical blue line indicating the maximum distance for an inclination i = 90°). We think the material is instead corotating with the white dwarf, with the vertical orange line indicating the distance at which the corotation velocity equals the observed velocity in Hα. The orange shaded region shows where we expect the material to lie, depending on inclination (with the vertical orange line indicating the minimum distance for i = 90°).

The velocity of material forced into corotation (v_c = Ωd), increases with distance, and a velocity of 1700 km s⁻¹ in an edge-on system would indicate a distance of d ≈ 22 R_WD, while the distance would be larger for an inclined system. At this distance, the magnetic field is less than 50 kG, and it can be easily hidden in the low-resolution Hα spectrum. As can be seen in Fig. 17, at this distance we are outside the Keplerian corotation radius (r_k = GM_WD/Ω²), which is the radius at which the Keplerian and corotation velocities are equal and the centrifugal force and gravity are in balance at the rotational equator. Outside the Keplerian corotation radius, for material in rigid corotation with the white dwarf, the centrifugal force dominates over the gravitational attraction toward the star; however, the magnetic field strength in this region is high enough to trap the material on the closed field lines, preventing it from being ejected.

Several studies on magnetic OB stars and M dwarfs (Townsend et al. 2005a,b, 2007; Townsend 2008; Bouma & Jardine 2025) show that, for a sufficiently large magnetic field strength, ionized material can be trapped in a stellar magnetosphere at radii larger than the Keplerian corotation radius. Within the magnetospheric radius, where ionized gas is forced to move along field lines, gas that is also within the corotation radius, where gravity is stronger than the centrifugal force, will follow field lines toward the star and accrete. On the other hand, if the material is outside the corotation radius, the centrifugal force will push the material outward as much as the closed field line will allow. Inside the magnetosphere and outside the corotation radius, the maximum distance from the white dwarf at which the material can be outwardly driven by centrifugal forces is set by the shape and finite extent of the field line to which the material is bound. More precisely, spatial potential minima will naturally occur beyond r_k in the case of large magnetic fields, as a consequence of magnetic pressure acting on the ionized gas. These spatial minima align into accumulation surfaces or volumes into which circumstellar emitting material may settle. This steady accumulation leads to the development of a warped disk, which is supported by the centrifugal force and corotates rigidly with the star (Townsend 2008). Although this mechanism was proposed to explain circumstellar emissions of rotating magnetic OB stars, it is not specific only to this class of objects. As such, we consider this to also be the most probable explanation for the consistent radial location of the circumstellar emitting material in ZTF J2008+4449.

4.3. A half ring in corotation

The double-peak profile of the Balmer emission in the phase-averaged spectrum of ZTF J2008+4449 resembles closely the profile of a thin, photoionized disk (see, e.g., Rybicki & Hummer 1983; Horne & Marsh 1986). However, when we analyzed the phase-resolved spectra, we observed that the bulk of the emission shifts abruptly between a maximally red-shifted velocity and maximally blue-shifted one every half a period (see Fig. 3). This behavior can be explained if the ionized gas is not distributed in a full disk, but is confined in a fraction of approximately f_d ≈ 50% of disk that is thin in the radial direction (a half ring). This is also suggested by the Doppler tomogram that we derived in Section 3.4.

We here build a simple toy model to reproduce the observed phase-resolved variations of Hα and Hβ. We take a half flat ring in corotation with the white dwarf, emitting uniformly, and we indicate the inner and outer radius of the half ring as R_in and R_out. We consider two simple cases for the nature of the emission. If the material is optically thin to Hα photons, we expect the emission from each section of the ring to be isotropic. If, on the other hand, the material is optically thick, we employ the Sobolev approximation (Hummer & Rybicki 1985) and impose an intensity that depends on the velocity gradient (along the line-of-sight) within the ring (Rybicki & Hummer 1983). Under this assumption, and assuming an edge-on inclination, the intensity at any point in the disk is given by

$$\begin{array}{c}\hfill I(\theta )=I_{0}sin\theta cos\theta ,\end{array}$$(6)

where θ is the azimuthal angle along the ring defined with respect to the line-of-sight. In the optically thin case, the intensity is given simply by I(θ) = I₀, where I₀ is the local intensity which in the following we assumed has no radial dependence.

We then simulated the emission from the half ring as it rigidly corotates around the white dwarf. In all calculations, we assumed an inclination of 90° of the normal to the ring with respect to the line-of-sight; for a lower inclination i, the model remains the same, but the inner and outer radii are shifted outward by a factor 1/ sin i. Fig. 18 shows the trailed spectrogram for Hα and Hβ in ZTF J2008+4449 compared with the phase-resolved intensity profile of our model in the optically thin and the optically thick case for a choice of R_in = 20R_WD and R_out = 35R_WD. We observed that the jumping behavior of Hα is well reproduced in both models.

Fig. 18. Comparison of the trailed spectrogram from LRIS for the Hα and Hβ emission of ZTF J2008+4449 (upper-left and lower-left plots, respectively, in blue) with our model for the optically thin (central plots) and optically thick (right plots) emission from a half ring with Rin = 20RWD and Rout = 35RWD.

We employed our simple model to fit the phase-averaged profile of the Hα emission. In the fit, we leave as free parameters the normalization (which we define as the inverse of the intensity, I₀), and the inner and outer radius of the ring. Additionally, we also fit the system velocity, v_sys, allowing the central wavelength of the profile to shift from the rest wavelength of Hα. We performed a fit for the optically thin case, for the optically thick case and for a combination of both; in this last fit, we allowed the fraction of optically thin to optically thick emission to vary freely as an additional parameter. The results of the three fits are shown in Fig. C.2 in the Appendix. We find that the best-fitting models are the optically thin one or the mixed model, while the fully optically thick model underpredicts the central section of the profile. Independently of the model, we recover very similar parameters for R_in, R_out, and v_sys. In Fig. 20, we show the phase-averaged profile for the best-fit for the mixed model, while we show the posteriors of the fit in the Appendix Fig. D.4. The best-fit parameters are: R_in = 20.1 ± 0.2R_WD, R_out = 34.9 ± 0.3R_WD and $v_{\mathrm{sys}}=0.8_{-8.9}^{+9.6}$ km/s. From the profile, we also extract the total flux in Hα to be $F_{\mathrm{H}\alpha }=(5.{02}_{-0.05}^{+0.06})\times {10}^{-16}$ erg s⁻¹ cm⁻². Fig. 19 shows the phase-resolved LRIS spectra of the Hα line (left panel) next to the phase-resolved profile of the best-fit for the mixed model (middle panel) and their comparison (right panel). We observed that, although ours is a simple and likely imperfect model (as suggested by minor discrepancies, such as excess emission at positive and negative velocities in the data at phases 0 and 7 respectively, in Fig. 19), it performs extremely well in reproducing the behavior of the Hα profile.

Fig. 19. Left: Observed phase-resolved LRIS spectra for the Hα profile of ZTF J2008+4449, divided into 12 phases (similar to right panel of Fig. 3). Middle: Phase-resolved profile for the best-fit model shown in the upper panel of Fig. 20. Right: Comparison of the observed phase-resolved profile with the model.

The dominant emission feature detected in this system is the Hα profile. The only other detected emission feature is Hβ, which is detected at much lower signal-to-noise. Nonetheless, the phase-resolved behavior of the higher-order Balmer line closely follows that of Hα (see Fig. 18). We therefore also performed a fit to the phase-averaged Hβ profile with the best-fitting Hα model in order to estimate the Hβ flux, leaving only the normalization as a free parameter (middle panel of Fig. 20). The posterior distribution from the fit yields a best-fit flux of F_Hβ = (2.8 ± 0.1)×10⁻¹⁶ erg s⁻¹ cm⁻². Compared to the Hα flux, this yields an intensity ratio in photon flux of R = 2.4 ± 0.1. We corrected this value for extinction by scaling the measured line fluxes by the best-fit extinction curve found in our SED fitting procedure (see Section 3.2), finding a line ratio of R = 2.3 ± 0.1. The bottom panel of Fig. 20 shows the measured, extinction-corrected line ratio and associated 1-, 2-, and 3-σ limits, compared to the simulated line ratios of Storey & Hummer (1995). In Section 4.5, we present a photoionization model of the Balmer emission, deriving an estimated electron density of n_e ≈ 10¹² cm⁻³. The corresponding curve in Fig. 20 is indicated with the thick, orange, solid line. The measured line ratio favors temperatures in the T_e ≈ 10 000 − 15 000 K range, consistent with the expectation for circumstellar photoionized hydrogen (e.g., Spitzer 1978; Cunningham et al. 2015), with a larger range of 5000−20 000 K permitted at 3-σ. In Section 4.5, we also estimate the vertical thickness of the ring, from which we derive a temperature, assuming hydrostatic equilibrium and an aligned dipole for magnetic field, of ≈9000 K. An improved radiative transfer model, accounting for detailed geometry, will be necessary to better constrain the ring temperature.

Fig. 20. Top: Five-parameter fit to the Hα emission profile employing our half-ring model and a mix of optically thick and optically thin emission, same as right-middle and right-bottom panels of Fig. C.2. The upper panel shows the comparison of our model (in red) to the phase-resolved Hα profile observed with LRIS for ZTF J2008+4449 (in blue), while the lower panel shows the residuals. Middle: One-parameter fit to the Hβ emission profile for ZTF J2008+4449, fixing all parameters to those found for Hα (see Fig. D.4), except the normalization, which is the only fit parameter. Bottom: Line ratios of Hα/Hβ as presented by Storey & Hummer (1995) for a range of gas temperatures and electron densities, shown with points and connecting lines. The measured Hα/Hβ ratio from the emission line fitting (top two panels) is shown with the horizontal dotted blue line. The 1-, 2-, and 3-σ contours are shown in the filled horizontal bands. The simulation data (points and connected lines) most closely matching the inferred electron density from our photoionization ring model are plotted with greater thickness.

4.4. Hydrogen absorption

The phase-averaged spectroscopy reveals a narrow absorption feature near the rest wavelength of Hα. In the left panel of Fig. 21, we show a Gaussian fit to the narrow Hα absorption feature. We find an equivalent width of EW = 8.4 ± 2.5 km s⁻¹. In the same figure, we also show the curves of growth of the Hα line as a function of optical depth and column density in the middle and right panels, respectively (see, e.g., Eq. 9.27 of Draine 2011). The measured equivalent width implies that the absorbing medium is optically thin to Hα for thermal broadening v > 5 km s⁻¹ (τ_Hα < 2). This value, however, is only a lower limit to the optical depth, because the half ring could be inclined with respect to our line of sight and it is likely thin enough that it does not entirely occult the star. For comparison, we also indicate the optical depth derived from the photoionization model (in the next Section, we find an optical depth in emission of ≈2, vertical dashed line in the middle panel of Fig. 21).

Fig. 21. Left: Continuum-normalized spectrum near the rest wavelength of Hα. Before the continuum normalization, the fit emission feature (Fig. C.2) was subtracted. The dashed line shows a Gaussian fit to the data. The fitting procedure yields an equivalent width of EW = 8.4 ± 2.5 km s−1 and central velocity of v0 = −37 ± 13 km s−1. Middle and right: Curves of growth for Hα. Middle: Evolution of the equivalent width as a function of optical depth. Right: Similar to middle panel, but as a function of column density. The dotted horizontal line indicates the measured equivalent width from the left panel. The vertical dashed line indicates the optical depth and column density inferred from the photoionization model.

In the COS UV spectrum, we also detect narrow absorption lines of carbon and silicon at their zero-field rest wavelengths (cyan lines in Fig. 4), which cannot be photospheric. These could be due to the interstellar medium or also originate from the trapped circumstellar material, similar to the Hα absorption. However, the uncertainty on the continuum level and blending with photospheric absorption lines does not allow us to derive a confident equivalent width. Future higher resolution observations could help understand if the lines are indeed from the half ring, allowing to better constrain its composition.

4.5. Properties of the glowing ionized gas

In this section, we estimate the properties of the material glowing in Hα and Hβ by assuming that the gas is ionized by the hot white dwarf and by assuming a half-ring geometry. Also, for simplicity, we assume a pure-hydrogen composition. The Balmer emission is very weak, and we cannot exclude the presence of other elements. Helium, for example, can be easily hidden in emission because the number of ionizing photons from the white dwarf is significantly lower than for hydrogen (from the best-fitting magnetic model, the number of ionizing photons is 20 times smaller for hydrogen than helium). Other elements, such as carbon, might be present in absorption in the UV spectrum, but we do not have a strong constraint on their abundance. The magnetically broadened photospheric metal lines observed in absorption in the UV spectrum do not provide additional information on the composition of the accreted material, as they could be easily explained by accretion of material at different metallicities (even solar) or by radiative levitation of metals from the interior of the white dwarf. We leave a detailed photoionization model, which would take into account different compositions, for future studies, and we here focus on deriving a rough estimate of the physical properties of the gas.

For the geometry of the half ring, we imposed the inner and outer radius of the ring (R_in and R_out) to be the ones derived in Section 4.3, based exclusively on the velocity profile. For the volume of the ring, we used the simple expression

$$\begin{array}{c}\hfill V=2\pi f_{d}H(R_{\mathrm{out}}^2-R_{\mathrm{in}}^2),\end{array}$$(7)

where H is half of the vertical thickness of the ring and f_d is the fraction of the ring, which we assumed to be one half, based on the trailed spectrogram and Doppler tomogram. The actual geometry of the trapped material could be much more complex, as it is likely determined by the morphology of the magnetic field and by the inclination of the magnetic axis with respect to the rotation axis (see Townsend et al. 2005a); however, with our simple model, we can still obtain an estimate of the properties of the trapped gas, in particular its density and optical depth.

The inner section of the ring, where the bulk of the Hα is emitted, is likely almost fully ionized by the incoming radiation, while at the outer edge of the ring, we expect a thin transition region between mostly ionized material and mostly neutral material (the edge of a “Strömgren disk”). Additional material might be trapped in the magnetosphere at larger radii (because of the likely presence of some small ionization fraction which would couple the material to the field), but we do not detect it in Hα because the ionization fraction is too low. We can write the equilibrium between ionization and recombination in the volume, V, of the ring as

$$\begin{array}{c}\hfill Q_H=V{\alpha }_{\mathrm{B}}{n}_{\mathrm{p}}{n}_{\mathrm{e}}=V{\alpha }_{\mathrm{B}}{n}_e^2,\end{array}$$(8)

where Q_H is the rate of ionizations in the ring and ${\alpha }_{\mathrm{B}}(T)\approx 2.56\times {10}^{-13}{T}_4^{-0.83}$ cm³ s⁻¹ is the recombination rate coefficient in the so-called Case B (which assumes that all Lyman-line photons are absorbed in the medium and converted into lower-series photons, Osterbrock & Ferland 2006); the assumption of Case B also involves the “on-the-spot” approximation, in which recombinations to the ground state do not contribute to the recombination rate as they produce photons just above the hydrogen ionization edge, which are reabsorbed locally. n_p and n_e are the number densities of protons and electrons in the ring, assumed to be uniform, and, since we are assuming a pure hydrogen composition, they are equal to each other and to the total ionized density. As we expected all the ionizing photons to be absorbed within the ring, we could impose that the rate of ionization is given by the number of ionizing photons entering the ring. If we take the total number of ionizing photons emitted by the white dwarf that we derive from our best-fitting magnetic model, L_ion (see Section 3.2 and Table 3), the number of photons entering the ring is given by L_ion times the ratio of the areas of the inner rim of the ring and that of the full sphere at R_in. The area of the inner rim of the ring can be expressed as A_rim = 2πR_in ⋅ 2H ⋅ f_d. We obtained

$$\begin{array}{c}\hfill Q_H=\frac{L_{\mathrm{ion}}{A}_{\mathrm{rim}}}{4\pi R_{\mathrm{in}}^2}=L_{\mathrm{ion}}{f}_d\frac{H}{R_{\mathrm{in}}}·\end{array}$$(9)

From Eqs. (8) and (9), we derived the ionized density in the ring, which does not depend on the thickness or the angular extent of the ring, as

$$\begin{array}{c}\hfill n_e=\sqrt{\frac{L_{\mathrm{ion}}}{2\pi {\alpha }_{\mathrm{B}}{R}_{\mathrm{in}}(R_{\mathrm{out}}^2-R_{\mathrm{in}}^2)}}\approx 1\times {10}^{12}{\mathrm{cm}}^{-3}.\end{array}$$(10)

Here, the value was calculated by employing the R_in and R_out obtained in Section 4.3. From the ionized density in the ring, we derived the neutral density assuming local equilibrium between ionization and recombination as

$$\begin{array}{c}\hfill n_0=\frac{{\alpha }_{\mathrm{B}}{n}_e^2}{q_{\mathrm{ion}}}\approx 2\times {10}^7{\mathrm{cm}}^{-3},\end{array}$$(11)

where q_ion = (∫L_ionσ_pdν)/(4πR_in²) is the ionization rate along the photon path in the radial direction, ${\sigma }_{\mathrm{p}}={\sigma }_{0\mathrm{p}}{(\nu /{\nu }_0)}^{-3}$ is the photoionization cross section for hydrogen (σ_0p ≈ 6.4 × 10⁻¹⁸ cm²), and ν₀ is the photoionization threshold for hydrogen (13.6 eV).

To estimate if the ring is optically thick to Hα photons, we need to estimate the population in the second excited level of hydrogen. Since the only relevant transitions are captures of electrons into excited levels and downward radiative transitions to lower levels, we can write the ratio of densities in the second and ground level, adopting the low density, coronal approximation (in which collisional processes may be neglected), as

$$\begin{array}{c}\hfill \frac{n_2}{n_1}={\alpha }_{\mathrm{B}}{n}_e\frac{1}{A_{21}^{\mathrm{eff}}}={\alpha }_{\mathrm{B}}{n}_e\frac{{\tau }_{\mathrm{L}\alpha }}{A_{21}},\end{array}$$(12)

where α_Bn_e is the recombination rate and $A_{21}^{\mathrm{eff}}$ is the effective transition rate, which can be expressed as the transition rate (A₂₁ = 4.41 × 10⁷ s⁻¹) times the escape probability of Lyman-α photons (ϵ ≈ 1/τ_Lα). We note that τ_Lα is the optical depth to Lyman-α photons in the ring. If we define δr = R_out − R_in as the radial extent of our ring, we can find the optical depth to Lyman-α photons in the radial direction as

$$\begin{array}{c}\hfill {\tau }_{\mathrm{L}\alpha }=\delta r{n}_0{\sigma }_{\mathrm{L}\alpha },\end{array}$$(13)

where σ_Lα is the absorption cross section for Lyman-α photons.

The frequency dependence of the absorption cross section for Lyman alpha, σ_Lα is well described by the Voigt function, ϕ(x) (see, e.g., Dijkstra 2014), such that

$$\begin{array}{c}\hfill {\sigma }_{\mathrm{L}\alpha }(x)={\sigma }_{\mathrm{L}\alpha ,0}·\varphi (x),\end{array}$$(14)

where σ_Lα, 0 is the absorption cross section at rest frequency, and x is the dimensionless frequency defined as x ≡ (ν − ν_α)/Δν_D, where thermal broadening is given by

$$\begin{array}{c}\hfill \mathrm{\Delta }{\nu }_D={\nu }_0\sqrt{\frac{2kT}{m{c}^2}}\equiv \frac{{\nu }_0{v}_{\mathrm{th}}}{c}·\end{array}$$(15)

The absorption cross section at the rest frequency is given by

$$\begin{array}{cc}\hfill {\sigma }_{\mathrm{L}\alpha ,0}& =\frac{f_{12}}{\sqrt{\pi }\mathrm{\Delta }{\nu }_{\mathrm{D}}}\frac{\pi e^2}{m_{e}c}\hfill \\ \hfill & =5.88\times {10}^{-14}\left(\frac{f_{12}}{0.416}\right)\left(\frac{{\lambda }_0}{1215.16\mathrm{\AA }}\right)\hfill \\ \hfill & \times {\left(\frac{T}{{10}^4\mathrm{K}}\right)}^{-1/2}{\mathrm{cm}}^2\hfill \\ \hfill & =3.79\times {10}^{-14}\left(\frac{f_{12}}{0.416}\right)\left(\frac{{\lambda }_0}{1215.16\mathrm{\AA }}\right)\hfill \\ \hfill & \times {\left(\frac{v_{\mathrm{th}}}{20\mathrm{km}/\mathrm{s}}\right)}^{-1}{\mathrm{cm}}^2,\hfill \end{array}$$(16)

where f₁₂ = 0.416 is the Lyman-α oscillator strength and λ₀ is the transition rest wavelength. The Voigt profile core can be well approximated by ϕ(x)≈e^−x2 in the limit that the broadening is dominated by velocity broadening rather than broadening of the intrinsic line profile. We adopted a normalization such that ϕ(0) = 1, and the integrated profile becomes $\int \varphi (x)dx=\sqrt{\pi }$. We can thus express the Lyman-α absorption cross section as

$$\begin{array}{cc}\hfill {\sigma }_{\mathrm{L}\alpha }& =\sqrt{\pi }{\sigma }_{\mathrm{L}\alpha ,0}\hfill \end{array}$$(17)

$$\begin{array}{cc}& =6.7\times {10}^{-14}\left(\frac{f_{12}}{0.416}\right)\left(\frac{{\lambda }_0}{1215.16\mathrm{\AA }}\right)\hfill \\ \hfill & \times {\left(\frac{v_{\mathrm{th}}}{20\mathrm{km}/\mathrm{s}}\right)}^{-1}{\mathrm{cm}}^2.\hfill \end{array}$$(18)

From Eqs. (12), (13) and (18), we obtained

$$\begin{array}{c}\hfill \frac{n_2}{n_1}={\alpha }_{\mathrm{B}}{n}_e\frac{\delta r{n}_0{\sigma }_{\mathrm{L}\alpha }}{A_{21}}\approx 4\times {10}^{-5},n_2\approx 7\times {10}^2{\mathrm{cm}}^{-3}.\end{array}$$(19)

We could then derive the optical thickness of the ring in Hα as

$$\begin{array}{cc}\hfill {\tau }_{\mathrm{H}\alpha }& =\delta r{n}_0\left(\frac{n_2}{n_1}\right){\sigma }_{\mathrm{H}\alpha }\hfill \\ \hfill & ={(\delta r{n}_0)}^2{\alpha }_{\mathrm{B}}{n}_e\frac{{\sigma }_{\mathrm{L}\alpha }{\sigma }_{\mathrm{H}\alpha }}{A_{21}}\approx 2,\hfill \end{array}$$(20)

where the Hα absorption cross section is given in a manner similar to Eq. (18) as

$$\begin{array}{cc}\hfill {\sigma }_{\mathrm{H}\alpha }& =\sqrt{\pi }{\sigma }_{\mathrm{H}\alpha ,0}\hfill \end{array}$$(21)

$$\begin{array}{cc}& =5.5\times {10}^{-13}\left(\frac{f_{23}}{0.641}\right)\left(\frac{{\lambda }_0}{6562.8\mathrm{\AA }}\right)\hfill \\ \hfill & \times {\left(\frac{v_{\mathrm{th}}}{20\mathrm{km}/\mathrm{s}}\right)}^{-1}{\mathrm{cm}}^2.\hfill \end{array}$$(22)

The value of order unity for the optical depth of the Hα line is consistent with the fit to the Hα profile derived in Section 4.3, where we derived that the best-fitting model is a combination of optically thick and optically thin emission. On the other hand, the equivalent width of the narrow Hα absorption feature (see Section 4.4), implies an optical depth of τ_Hα ≈ 0.3, suggesting the absorbing medium is optically thin. This discrepancy could indicate that the absorbing material has a lower column density than the emitting material. Such a situation could arise, for instance, in a thin, uniform density (partial) ring at a modest inclination of ≈1 deg away from edge-on, since the absorption is likely absorption of the stellar continuum, for which the circumstellar gas must be between the star and the viewer. The Balmer emission, on the other hand, could be seen at any inclination. This discrepancy could also be explained by uncertainties arising from the simplistic approximate treatment of geometry and radiative transfer presented in this study. For example, the gap between the n_e = 10¹¹ cm⁻³ and n_e = 10¹² cm⁻³ Balmer line ratio curves in the bottom panel of Fig. 20, which shows significant overlap with the rough range of expected electron densities of the ionized gas, probably arises due to collisional (de-)excitation effects becoming important. In this analysis, however, we do not account for this effect. Such uncertainties arising from the simplicity of the employed models may be tested via future studies exploring detailed geometric and radiative transfer models, although this is beyond the scope of this paper.

In all calculations until now, we did not have to specify the vertical thickness, H, or the angular extent, f_d, of the ring. However, if we want to determine the total mass contained in the ring, we need to derive its volume. One way to estimate the thickness of the ring is to calculate how many ionizing photons need to be captured in the ring to produce the observed Hα emission. We can express the total Hα luminosity (in photons per second) as

$$\begin{array}{cc}\hfill L_{\mathrm{H}\alpha }& =4\pi j_{\mathrm{H}\alpha }V\hfill \\ \hfill & =\left(\frac{4\pi j_{\mathrm{H}\beta }}{n_e{n}_p}\right)\left(\frac{j_{\mathrm{H}\beta }}{j_{\mathrm{H}\alpha }}\right)n_e{n}_{p}V\hfill \\ \hfill & =\zeta (T)n_e^{2}V,\hfill \end{array}$$(23)

where j_Hα and j_Hβ are the emissivities of the two lines and V is the volume of the ring. We grouped the two terms in parentheses into the function ζ(T), highlighting that they mostly depend on the temperature of the gas, although weakly (they also have an extremely weak dependence on density, which we ignore). The values for these two terms for Case B are tabulated for different temperatures in Table 4.4 of Osterbrock & Ferland (2006).

By combining Eqs. (8), (9), and (23), we derived the following expression for H:

$$\begin{array}{c}\hfill H=\frac{{\alpha }_{\mathrm{B}}(T)}{\zeta (T)}\frac{L_{\mathrm{H}\alpha }}{L_{\mathrm{ion}}}\frac{r_i}{f_d}\approx 0.1R_{\mathrm{WD}},\end{array}$$(24)

where we used values for α_B and ζ at a temperature of 10 000 K; however, the two parameters have a similar dependence on temperature, so the ratio of the two is almost independent of temperature. This value of H is small compared to the inner radius of the ring (≈20 R_WD), indicating a very thin disk; this is expected because the centrifugal force, which pushes the material toward the apex of each magnetic field line, is strong and acts together with gravity to reduce the thickness of the disk. Townsend et al. (2005a) derive a modified expression for the hydrostatic scale height of the trapped disk in the case of a dipole aligned with the rotation axis (their Eq. 28), which includes the effects of the centrifugal acceleration. Using the same expression, we derived the temperature of the gas that would result in the scale height derived in Eq. (24):

$$\begin{array}{c}\hfill T_{\mathrm{ring}}=H^2\frac{\mu }{2k_{\mathrm{B}}}\frac{G{M}_{\mathrm{WD}}}{r_K^3}(3-\frac{2r_K^3}{R_{\mathrm{in}}^3})\approx 6000\mathrm{K},\end{array}$$(25)

where μ indicates the mean molecular weight (which we took as one for pure hydrogen) and k_B is the Boltzmann constant. Although we did not expect the magnetic field to be a perfectly aligned dipole, we found a consistent temperature to what is expected from circumstellar photoionized gas and with the temperature estimate that we derived in Section 4.4 from the Hα absorption.

Finally, we estimated the total mass in the ring to be

$$\begin{array}{c}\hfill M_{\mathrm{ring}}\sim n_{e}V{m}_p\approx 2\times {10}^{-17}{\mathrm{M}}_{\odot }.\end{array}$$(26)

We note that this mass estimate serves as a lower limit, since this corresponds to the mass of Balmer-emitting hydrogen. The true total mass could be considerably higher when accounting for the presence of non-emitting hydrogen (either due to low ionization fraction, or low temperatures), and the presence of heavier elements such as helium and metals.

4.6. Mass depletion timescale in the half ring

In the model presented in the previous Section, the ionized gas glowing in Hα is currently trapped in a stable configuration in corotation with the white dwarf because ionized particles strongly feel the effect of the magnetic field. Upon recombination, however, neutral hydrogen atoms can decouple from the magnetic field line to which they were previously bound and move outwardly according to their orbital velocity at the time of recombination (i.e., tangentially to the circular trajectory of the ionized arc) until they are reionized. This effect will create a drift between neutral and ionized particles and therefore an outward motion that will slowly deplete the reservoir of ionized gas. In reality, the neutral particle will feel a drag from the charged particles around it, so its motion will not be truly ballistic. The ion-neutral drag is dominated by charge transfer and elastic collisions. At densities below 10¹¹, there are no collisions within the time for photoionization, and the motion is expected to be ballistic. At densities above 10¹⁴, on the other hand, there are thousands of collisions in a photoionization time, and the motion is diffusive. The density of 10¹² that we just derived is in the more difficult intermediate range; calculating an accurate estimate of the timescale of outward motion of the gas would require an estimate of the strength of the drag between ionized and neutral particles in this intermediate regime, which is beyond the scope of this paper. However, we expect the timescale to be shorter than the cooling age of the white dwarf (≈60 Myr), and therefore the ionized gas trapped in the magnetosphere would have had to be replenished by material entering the magnetosphere from outside or from outflows from the white dwarf itself.

4.7. Period increase: A wind or a propeller mechanism

In the previous sections, we focused only on the ionized gas glowing in Hα. The X-ray emission of ZTF J2008+4449 and the detected increase in its spin period, however, highlight the presence of additional material interacting with the rapidly spinning white dwarf. In particular, the observed Ṗ = (1.80 ± 0.09) × 10⁻¹² s/s is much larger than the period increase expected from magnetic dipole radiation:

$$\begin{array}{c}\hfill {\dot{P}}_{\mathrm{dipole}}=\frac{8{\pi }^2}{3c^3}\frac{{sin}^2\alpha R_{\mathrm{WD}}^6{B}_p^2}{I_{\mathrm{WD}}P}\approx 5\times {10}^{-14}\mathrm{s}/\mathrm{s},\end{array}$$(27)

where I_WD = (βR_WD)²M_WD is the moment of inertia of the white dwarf, β is the gyration radius (Claret 2023), and α is the inclination of the magnetic field to the rotational axis. We have used a value of β ≈ 0.39, as per the models of Claret (2023) and for α we have taken the mean of sin²α = 0.5. This value of Ṗ_dipole implies that the change in spin period is most likely driven by interaction with the circumstellar material, and, in particular, material being ejected from the system, carrying away angular momentum. Furthermore, as we commented in the previous Section, the ionized gas trapped in the magnetosphere likely needs to be replenished over smaller timescales than the age of the white dwarf.

If ionized gas can be extracted from the surface of the white dwarf, as for example in a magnetically driven wind, some could get trapped in the magnetosphere, populating the half ring, and some could escape, carrying away angular momentum. Alternatively, if diffuse infalling material is present outside the magnetosphere of the white dwarf (in Section 4.9 we discuss the possible origins of this material), the loss of angular momentum necessary to produce the observed spin-down rate may be driven by the binding (and subsequent outward propelling) of this material into rigid body corotation with the white dwarf once it makes contact with the magnetosphere (i.e., at the magnetospheric radius, Eq. 5). A similar propelling mechanism has been suggested for the white dwarf propellers such as AE Aquarii (Wynn et al. 1997; Eracleous & Horne 1996) and the so-called transitional pulsars (see Campana et al. 2018; Papitto & de Martino 2022, and references therein). The amount of angular momentum carried away by the ejected material in this mechanism depends on how deeply it can plunge within the magnetosphere and on how efficiently it is accelerated before being ejected. Similarly, in the case of a wind, the loss of angular momentum depends on how far from the star the material detaches from the field.

One way to estimate the minimum infalling or outflowing mass flux, Ṁ, would be by assuming that the material is accelerated into rigid corotation and ejected at the magnetospheric radius. If we equate the angular momentum transfer necessary to bring the infalling magnetospheric radius to the angular momentum losses of the white dwarf, we can express Ṗ as

$$\begin{array}{c}\hfill \frac{dln(P)}{\mathit{dt}}=\frac{\dot{P}}{P}\approx \frac{\dot{M}{R}_A^2}{I_{\mathrm{WD}}}\approx \frac{\dot{M}}{{\beta }^2{M}_{\mathrm{WD}}}{\left(\frac{R_{\mathrm{m}}}{R_{\mathrm{WD}}}\right)}^2.\end{array}$$(28)

Since the magnetospheric radius depends on Ṁ, we could in principle substitute the equation for the magnetospheric radius (Eq. 5) in Eq. (28) and derive an Ṁ ≈ 1 × 10⁻¹⁷M_⊙ yr⁻¹. However, at such a small mass loss rate, the magnetospheric radius lays beyond the light cylinder, which is the distance from the white dwarf at which the corotation velocity becomes equal to the speed of light. A significant transfer of angular momentum between the star and the circumstellar material could only happen if the infalling material reaches within the light cylinder (Papitto & de Martino 2022). Similarly, a wind would have to detach from the field at or within the light cylinder. Therefore, the ejected mass has to be higher to explain the observed Ṗ. A better way to estimate the mass loss rate is to calculate how much material would have to be ejected at the speed of light at the light cylinder to carry away enough angular momentum from the star to explain the Ṗ. We can derive this value by using the expression for the light cylinder radius (R_LC = cP/(2π), where c is the speed of light) in Eq. (28) instead of R_m. We obtain Ṁ ≳ 2 × 10⁻¹⁵M_⊙ yr⁻¹. In reality, the material will likely plunge or become detached much deeper in the magnetosphere and will be ejected at much lower velocities, so this is only an absolute minimum to Ṁ. A third way of estimating Ṁ is to impose that the magnetospheric radius is smaller than the light cylinder radius, which yields a limit Ṁ ≳ 8 × 10⁻¹⁵M_⊙ yr⁻¹. This is likely the more realistic and stringent constraint. For both estimates, only a tenth or less of the spin-down luminosity of the white dwarf (that we estimate at $L_{\mathrm{sd}}=I_{\mathrm{WD}}\omega \dot{\omega }\approx 9\times {10}^{31}$ erg/s) would be converted into the kinetic energy carried away by the ejecta.

Comparing these Ṁ limits with the estimate of the mass trapped in the half ring that we derive in Section 4.5, we find that the material glowing in Hα is only a very small fraction of the material that is being ejected further out. So even in the case of a very efficient propeller, the ring could be replenished by material entering the magnetosphere from outside. Close to the rotation axis, the centrifugal acceleration for material in rigid corotation with the white dwarf is weak out to large radii; some material could therefore trickle down to field lines much closer to the white dwarf than the magnetospheric radius from above the rotation axis. Some of this material can be captured by field lines that extend right outside the corotation axis and therefore replenish the ring that we see glowing in Hα, while some material could reach field lines within the corotation radius, and accrete onto the white dwarf.

4.8. X-ray emission mechanism

The observed X-ray emission (L_X = (2.3 ± 0.4)×10²⁹ erg s⁻¹) could be powered by a fraction of the spin-down luminosity, by the gravitational energy of material infalling into the potential well of the white dwarf, or by a combination of both. The X-ray spectrum is equally well described by a power-law model and by a two-temperature optically thin plasma model (see Section 3.5), indicating that the emission could be either thermal or nonthermal. If the X-rays are produced by the interaction of the magnetic field with the infalling (or outflowing) material at the edge of the magnetosphere, a thermal X-ray emission could be produced by shock-heated plasma where the material is being accelerated, as it has been suggested for the soft X-ray component of the white dwarf propeller AE Aquarii (Oruru & Meintjes 2012; Kitaguchi et al. 2014). Alternatively, high-velocity material trapped in the field at the edge of the magnetosphere could deform the field lines, causing continuous local reconnections; these magnetic reconnections could heat the plasma to high temperatures or accelerate particles to relativistic speed and lead to synchrotron emission (as it has been suggested for the pulsed X-ray emission of the binary white dwarf pulsar AR Scorpii, Takata et al. 2018).

Another possible origin for the X-ray emission is the region close to the surface of the white dwarf. If some material can reach field lines within the corotation radius, it can accrete onto the surface of the white dwarf; the observed X-ray emission could therefore be due to shock-heated plasma close to the surface of the white dwarf or by hot spots, similar to what is observed in polars. The Ṁ needed to power this emission would be of the order of 10⁻¹⁴ M_⊙ yr⁻¹. This scenario can also explain what could be cyclotron humps in the UV spectrum. Alternatively, Alfvén waves excited by the interaction of the material with the field at the magnetospheric radius can carry energy down to the white dwarf surface, giving rise to hotspots or local magnetic reconnections.

Variability of the X-ray emission could provide some important clues on the origin and nature of the emission. As we show in Section 3.5, the XMM light curves present a hint of variability on the period of the white dwarf. This could indicate that part of the emission is eclipsed by the white dwarf, hinting at an emission region close to the white dwarf surface, or that the X-ray emission is beamed, indicating the presence of high velocities in the emitting material. Deeper X-ray observations could provide a higher signal-to-noise X-ray light curve, needed to confirm the variability and constrain its amplitude.

4.9. Origin of the circumstellar material

Despite the host of evidence for the presence of circumstellar material, in Section 3.3 we have shown that the presence of a stellar or brown dwarf companion filling its Roche lobe and therefore of ongoing mass transfer as the source of this material can be excluded. There is still the possibility, as we show in Section 3.3, for a Roche-lobe-filling gas giant planet to be the source of mass transfer; however, because of the origin of the star as a double white dwarf merger, we find it unlikely that a planet would survive the evolution of the star and find itself at just the right distance from the white dwarf. Also, accretion from the interstellar medium would not be enough to explain the observed spin-down and X-ray emission (Wesemael 1979). The trapped material glowing in Hα, as well as the infalling (or outflowing) material that is at the origin of the spin-down and of the X-rays must therefore have a different origin.

One possible origin for the circumstellar material is the tidal disruption of a planetary body (e.g., a comet, asteroid, or planet). Accretion of planetary material onto white dwarfs is a common occurrence, as over a quarter of all white dwarfs show the presence of metals in their atmospheres (Koester et al. 2014; Zuckerman et al. 2003, 2010), and recently, X-rays produced by the accretion of planetary material have been detected on the prototypical polluted white dwarf G 29–38 (Cunningham et al. 2022). Due to the high surface gravity of white dwarfs, heavy elements sink quickly out of the observable layers, so the presence of these metals indicate recent accretion of debris from planetary bodies that have been scattered toward the central stellar remnant and tidally disrupted (e.g., Zuckerman et al. 2007; Xu et al. 2017; Vanderburg et al. 2015; Gänsicke et al. 2019; Manser et al. 2019). The scattering of planetary bodies to high eccentricities, where they can be tidally disrupted by the white dwarf, is thought to be caused by the chaotic evolution of the orbits of the planetary system after its central star loses mass in the last stages of the asymptotic giant branch, on its way to becoming a white dwarf (for a review, see Veras 2021; Veras et al. 2024). In the case of ZTF J2008+4449, if indeed the white dwarf is the product of a double white dwarf merger, the system has undergone two planetary nebulae and a merger event. No studies have been done on the effects of these combined events on the dynamics of a planetary system; however, we can imagine that, if a reservoir of planetary bodies has survived the previous stages of evolution of the system without being ejected or disrupted, the merger event might have perturbed the orbits of the remaining objects in the system, leading some of them to high enough eccentricities to be tidally disrupted by the white dwarf.

One possible caveat of this scenario is the composition of the infalling material. In the half ring, we only see emission of hydrogen and, although we cannot exclude the presence of other elements (see discussion in Section 4.4), it is unlikely for the composition of the ionized gas to be hydrogen poor. Planetary debris observed around polluted white dwarfs usually present a rocky composition (see, e.g., Xu et al. 2019). In those instances in which a gaseous, ionized disk is detected around polluted white dwarfs, it usually shows emission of calcium and other metals, but not hydrogen (e.g., Gänsicke et al. 2006; Manser et al. 2020). There is one notable exception, WD J0914+1914, a white dwarf accreting from an evaporating ice giant; however, even in this case, the hydrogen emission is accompanied by emission in silicon, oxygen and iron (Gänsicke et al. 2019).

Since ZTF J2008+4449 is the remnant of a double white dwarf merger, another possible explanation for the detected circumstellar gas is the fallback of material ejected during the merger event. Studies such as Rosswog (2007) and Ishizaki et al. (2021) propose the fallback of gravitationally bound ejecta following the merger of compact binary systems, in order to explain the long-lasting X-ray activity of young binary neutron star mergers such as GW 170817. Furthermore, most simulations of double white dwarf mergers predict that right after the merger, the remnant white dwarf is surrounded by a thick disk and an extended tidal tail (Guerrero et al. 2004; Lorén-Aguilar et al. 2009; Dan et al. 2014). Most of the material in the disk will be accreted rather quickly, on a viscous timescale, while material at high eccentricity in the bound tails accretes at longer timescales (Rosswog 2007; Lorén-Aguilar et al. 2009). Most simulations, however, only explore the first few minutes or days after the merger (see for example Lu & Quataert 2023), and no studies have been performed to explore if a very low level of accretion could still be ongoing 60 Myr after the merger. Furthermore, as in the planetary scenario, we expect the composition of the fallback material to be dominated by carbon, oxygen and other metals, with very little hydrogen.

A third possibility would be a wind of ionized gas flowing from the surface of the white dwarf. Radiatively driven winds are observed in hot pre-white dwarfs at low surface gravities (Werner et al. 1995; Dreizler et al. 1995), and for one system it has been suggested that some of the wind might be trapped in the magnetosphere of the star (Reindl et al. 2019). In the case of ZTF J2008+4449, however, because of its high surface gravity (log g ≈ 8.8) and relatively low temperature, a sizable wind cannot be driven by radiation pressure alone, even if metals are present in the atmosphere (see, e.g., Vennes et al. 1988; Unglaub 2008). As the star is rapidly spinning and highly magnetized, it is possible that the strong magnetic field could power a wind, as plasma is expected to flow outward on open field lines; moreover, if infalling material is creating stresses on the magnetosphere, Alfvén waves could induce reconnections close to the surface, accelerating material into a wind (in which case we would still need an external source of circumstellar material). Our estimate for the mass trapped in the magnetosphere is very small (2 × 10⁻¹⁷ M_⊙), and, even assuming a fast depletion timescale, only a weak wind would be required to replenish the half ring; the mass loss required to explain the period increase, on the other hand, is orders of magnitude larger. Current studies of magnetically driven winds mainly focus on main sequence and evolved stars and mostly on how the magnetic field accelerates winds that are powered by other mechanisms, such as radiation or reconnections (Weber & Davis 1967; Belcher & MacGregor 1976; Goldreich & Julian 1970; Thirumalai & Heyl 2010; Vidotto et al. 2014). Further studies are needed to explore if a wind can be extracted from such an extreme system in terms of surface gravity, stellar rotation and magnetic field strength as ZTF J2008+4449.

4.10. Evolutionary history of the rotation period

From the period and period derivative measured in Section 3.1, we measured a spin-down age of $P/\dot{P}\sim 7$ Myr, which is an order of magnitude lower than the cooling age of the white dwarf. However, we do not necessarily expect the period derivative to have been constant in the past. White dwarfs polluted by planetary material are believed to undergo episodic accretion, in which a small planetary body is scattered within the Roche radius of the white dwarf, where it is tidally disrupted and ground to dust, forming a debris disk around the white dwarf (Veras et al. 2021; Cunningham et al. 2025). The lifetime of these disks is uncertain, but it is expected to be of the order of one Myr (Girven et al. 2012; Cunningham et al. 2021).

If, on the other hand, the circumstellar material originated from the fallback of post-merger ejecta, the spin evolution is expected to change depending on the accretion rate (Sousa et al. 2022). At early times, when the high accretion rate pushes the magnetospheric radius within the corotation radius, the accretion of material onto the surface of the white dwarf can transfer angular momentum onto the white dwarf, resulting in a negative period derivative (the white dwarf spins up). As the accretion rate diminishes, the period derivative increases until the magnetospheric radius becomes larger than the Keplerian corotation radius. At this point, the material is unable to accrete and is instead ejected, carrying away angular momentum and braking the spin of the white dwarf (Campana et al. 2018; Papitto & de Martino 2022; Sousa et al. 2022).

4.11. Comparison with ZTF J1901+1458

ZTF J2008+4449 is not the only merger remnant to display the signature of circumstellar material in the absence of a companion star or brown dwarf. In a companion paper (Desai et al. 2025), we show that another merger remnant, ZTF J1901+1458, presents very similar characteristics to ZTF J2008+4449 as well the presence of circumstellar material through similar X-ray emission properties.

Just as the white dwarf treated in this work, ZTF J1901+1458 is massive (nearing the Chandrasekhar mass at ≈1.3 M_⊙), with a radius of ≈2600 km, and rapidly rotating on a very similar period of ≈6.9 minutes. These characteristics, together with its strong magnetic field of ≈800 MG point toward the fact that ZTF J1901+1458 is also most likely a double white dwarf merger remnant (Caiazzo et al. 2021). Most significantly, both systems are sources of soft X-ray emissions, and their X-ray spectra are strikingly similar to each other: when fit by a two-temperature optically thin plasma model, we recover similar temperatures (the soft component is $kT=0.{23}_{-0.03}^{+0.04}$ for ZTF J2008+4449 and kT = 0.23 ± 0.03 for ZTF J1901+1458, while the hard component is $kT=3.6_{-0.7}^{+1.0}$ for ZTF J2008+4449 compared to $kT=2.9_{-0.8}^{+1.2}$ for ZTF J1901+1458). Finally, for both systems, we have strong constraints from infrared photometry showing the absence of a stellar or sub-stellar companion.

Although we do not have a strong constraint on the origin of the circumstellar material in the two systems, the striking similarities between ZTF J2008+4449 and ZTF J1901+1458, in terms of mass, magnetic field, spin period, and X-ray emission, hint to the two system being part of the same class of merger remnants interacting with circumstellar material. However, the two systems show some differences. First, the X-ray luminosity is ≈100 times larger in ZTF J2008+4449 (L_X = (2.3 ± 0.4)×10²⁹ erg/s) as compared to that in ZTF J1901+1458 ($L_{\mathrm{X}}=1.{28}_{-0.15}^{+0.18}\times {10}^{27}$ erg/s). If the X-ray emission in both systems is at least partly powered by accretion, this could imply a larger mass influx rate in the former. Also, ZTF J1901+1458 does not display the Hα emission seen in the optical spectrum of ZTF J2008+4449 and we do not detect a significant spin-down, although our observational constraints for the Ṗ in ZTF J1901+1458 are much worse than for ZTF J2008+4449, and we can only derive a limit of Ṗ < 10⁻¹¹ s/s. These differences could be due to the difference in age between the two systems: the cooling age of ZTF J1901+1458, estimated through SED fitting, is 490 ± 10 Myr, about ten times older than ZTF J2008+4449. If, for example, the X-ray activity is due to fallback material from the merger, it would be expected for the accretion rate to diminish with time (although it is still to be understood if fallback accretion can happen at such late times). In the case of a mass loss from the surface of the white dwarf, on the other hand, the higher surface gravity of ZTF J1901+1458 might be the reason for a weaker wind.

5. Summary and outlook

We have presented our discovery and characterization of the hot, massive, and highly magnetized variable white dwarf ZTF J2008+4449 and provided a description of the data, analysis, and an interpretation of the physical properties of the white dwarf and surrounding material. In Table 4, we list the parameters that we have derived for the white dwarf, and we summarize our results in the following:

• The exceptionally high magnetic field (≈400 − 600 MG), short rotation period (P ≈ 6.6 minutes), and large mass (1.12 ± 0.03 M_⊙) of the isolated white dwarf ZTF J2008+4449 classify it as a likely double white dwarf merger remnant.

• ZTF J2008+4449 shows spectral and photometric variability across a broad range of wavelengths, from the optical to the FUV (and possibly X-rays), on its 6.6 minute rotation period. Although it is difficult to find a definitive explanation for all the variable features, as an ensemble they point toward the possibility of cyclotron emission, complex magnetic field geometry, and possibly the variation of chemical abundances across the surface of the white dwarf.

• The significant spin-down detected from the optical light curve, Ṗ = (1.80 ± 0.09) × 10⁻¹² s/s, along with the detection of Balmer emission lines and X-rays, indicate the presence of circumstellar material in close proximity that interacts with the white dwarf and its magnetosphere.

• The Doppler variation seen in the Balmer emission lines is consistent with that of a half ring of ionized gas in corotation with the white dwarf at a distance of ≈20 R_WD and a radial extent of ≈15 R_WD, implying the trapping of circumstellar gas in the magnetosphere at a distance larger than the Keplerian corotation radius.

• The detected spin-down can be explained by some material being ejected from the system and carrying away angular momentum from the star. We find that to explain the spin-down of the white dwarf, the mass loss is likely higher than Ṁ> 3 × 10⁻¹⁴ M_⊙/yr. The observed X-rays could be associated with the ejection region, being produced by shocks or magnetic reconnection. Alternatively, if a fraction of the infalling material manages to reach the surface of the white dwarf, the X-rays could be due to accretion onto the magnetic poles of the white dwarf.

• Infrared photometry excludes the presence of a Roche-lobe-filling brown dwarf companion as the source of the material through mass transfer.

• We identify three possible sources for the circumstellar material: the tidal disruption of a planetary object (a comet, asteroid, or small planet) scattered at high eccentricity toward the white dwarf, fallback accretion of gravitationally bound ejecta from the merger event, or a wind from the surface of the white dwarf.

• Although Gaia astrometry of ZTF J2008+4449 suggests that the white dwarf might be a member of the young open cluster RSG 5 (Miller et al. 2026; Bouma et al. 2022), we exclude this possibility because our re-derived age (Section 3.6) shows that the cluster cannot be substantially older than 50 Myr, which is comparable to the cooling age of the white dwarf (60 ± 10 Myr) and therefore likely much less than its total age.

The discovery of circumstellar material in two white dwarf merger remnants (ZTF J2008+4449 and ZTF J1901+1458) potentially opens a novel avenue toward a deeper understanding of compact binary evolution and the physics of compact mergers. Additionally, their study could also provide significant insight into accretion phenomena onto highly magnetic objects, given the likely complex nature of the interaction between the intense magnetic fields of the white dwarfs and their circumstellar material. Lastly, the striking similarities between these two uniquely behaving systems raise the possibility of their mutual membership in a new class of compact merger remnants interacting with circumstellar material.

A further in-depth analysis of both systems is required to confirm their nature and fully determine their physical properties. As such, we intend to conduct a more detailed photoionization study of the half ring of material glowing in Hα in ZTF J2008+4449 in order to produce a physically accurate model of the spatial geometry of the emitting region of ionized gas and to constrain its composition. Understanding the geometry and properties of the emission region will be a crucial step in constraining the ongoing accretion processes and relating them to the geometry of the magnetic field, the spin decay, and the X-ray emission of ZTF J2008+4449. Another key aspect in understanding the nature of ZTF J2008+4449 will be a better determination of the magnetic field strength on the surface of the white dwarf. The calculations reported in this work for the transitions of metals in high magnetic fields show that the observed absorption lines in the UV are likely stationary Zeeman components of elements such as carbon, nitrogen, and silicon (see Fig. 8); however, more detailed calculations are needed to correctly identify the lines and determine the magnetic field strength. An in depth analysis of the phase-resolved UV spectra with bespoke magnetic models of different geometries would also allow the geometry of the magnetic field to be constrained, which may be much more complex than that of a simple dipole, as hinted by the complex photometric variability (see Section 4.1) and the asymmetric structure of the trapped ionized gas.

Additional observational campaigns will also be crucial in understanding this peculiar white dwarf. Deep observations in the near and mid-infrared would allow us to constrain the presence of a cold disk or cloud of material surrounding the white dwarf or of planetary objects feeding the accretion. Radio observations would reveal possible nonthermal emission as the pulsating synchrotron emission observed in white dwarf pulsars. At the other end of the spectrum, deeper X-ray observations could better constrain spectral properties of the white dwarf and reveal if the X-rays are variable on the 6.6 spin period, providing additional information on the emission region.

Finally, discovering new white dwarfs such as ZTF J2008+4449 and ZTF J1901+1458 will allow the origin of these fascinating objects to be constrained and enhance our understanding of how they evolve. It is likely that many more of these systems are waiting to be discovered in the archives of time-domain surveys such as ZTF, ATLAS, and TESS and that many more are bound to be discovered with the Vera C. Rubin Observatory.

Acknowledgments

We thank Lynne Hillenbrand and Soumyadeep Bhattacharjee for helpful discussions, and Kishalay De for his help with the WIRC reduction pipeline. IC was supported by NASA through grants from the Space Telescope Science Institute, under NASA contracts NASA.22K1813, NAS5-26555 and NAS5-03127. TC was supported by NASA through the NASA Hubble Fellowship grant HST-HF2-51527.001-A awarded by the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., for NASA, under contract NAS5-26555. This project has received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (Grant agreement No. 101020057). This work was 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 Grants No. AST-1440341, AST-2034437, and currently Award #2407588. ZTF receives additional funding from the ZTF partnership. Current members include Caltech, USA; Caltech/IPAC, USA; University of Maryland, USA; University of California, Berkeley, USA; University of Wisconsin at Milwaukee, USA; Cornell University, USA; Drexel University, USA; University of North Carolina at Chapel Hill, USA; Institute of Science and Technology, Austria; National Central University, Taiwan, and OKC, University of Stockholm, Sweden. Operations are conducted by Caltech’s Optical Observatory (COO), Caltech/IPAC, and the University of Washington at Seattle, USA. 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. The Pan-STARRS1 Surveys (PS1) and the PS1 public science archive have been made possible through contributions by the Institute for Astronomy, the University of Hawaii, the Pan-STARRS Project Office, the Max-Planck Society and its participating institutes, the Max Planck Institute for Astronomy, Heidelberg and the Max Planck Institute for Extraterrestrial Physics, Garching, The Johns Hopkins University, Durham University, the University of Edinburgh, the Queen’s University Belfast, the Harvard-Smithsonian Center for Astrophysics, the Las Cumbres Observatory Global Telescope Network Incorporated, the National Central University of Taiwan, the Space Telescope Science Institute, the National Aeronautics and Space Administration under Grant No. NNX08AR22G issued through the Planetary Science Division of the NASA Science Mission Directorate, the National Science Foundation Grant No. AST–1238877, the University of Maryland, Eotvos Lorand University (ELTE), the Los Alamos National Laboratory, and the Gordon and Betty Moore Foundation. This work made use of Astropy (http://www.astropy.org): a community-developed core Python package and an ecosystem of tools and resources for astronomy (Astropy Collaboration 2013, 2018, 2022).

References

Appendix A: Zeeman shifted metal lines in the UV spectrum

As described in Section 2.3, the UV spectrum shows a suite of lines which we do not identify with either Hydrogen or Helium transitions in high magnetic fields. These lines are shaded in red in the topmost panel of Fig. 4, which we show once more in Fig. 8. We thus consider the possibility that these absorption lines are produced by the presence of light metallic elements in the atmosphere of the white dwarf. In a first attempt to constrain the ionic species present that could contribute to the lines’ formation and subsequently place limits on the surface magnetic field strength, we have investigated the behavior of the most frequently present light element transitions in hot dwarf novae with temperatures similar to that of our white dwarf. For this, we have compiled a list of transitions seen in such hot (≈ 30000 K) dwarf novae (see, e.g., Godon & Sion 2022) in the approximate wavelength range that contains the aforementioned unidentified absorption lines in our spectrum. These transitions are compiled in Table A.1 along with their rest wavelengths and belong to carbon, nitrogen, and silicon.

We then proceeded to calculate the transition energy of the Zeeman components of our compiled lines as a function of magnetic field from perturbation theory, including linear and quadratic operators in the Paschen-Back regime (Kemic 1975; Vennes et al. 2024). The calculations assume < r²> expectation values, or mean square radii (Cowan 1981), involved in the calculation of the quadratic term for individual levels in a zero-field environment (see Vennes et al. 2024). Self consistent calculations of atomic wave functions and < r²> expectation values at large magnetic fields should be performed to improve these results. Our current results are shown in Fig. 8 for line components with relatively significant oscillator strengths at a field of 500 MG. However, we remark that these results may not be fully accurate at the highest of fields, as a result of the perturbative methods utilized, which may break down at high magnetic field values. Unfortunately, more accurate calculations at high fields (see Hollands et al. 2023) are not currently available for these line transitions.

Our calculations show that the π Zeeman components of each line (the transitions with no change in the m quantum number, Δm = 0) do not deviate far from their zero-field wavelength, and they all cluster in the UV wavelength range where we detect the broad absorption lines. The other Zeeman components, on the other hand, vary rapidly with field and, thus, could be much more strongly affected by magnetic broadening, making them harder to detect. As such, it is plausible that the observed features are produced by some of the central components in these transitions, indicating that light metals could indeed be present at the surface of the white dwarf, in particular carbon and silicon.

Appendix B: XMM-Newton data and data statistics tables

Appendix C: Additional plots

Fig. C.1. CHIMERA light curve for ZTF J2008+4449. In blue, we show the unbinned data for the g (upper) and r (lower) filters, phase-folded using the correct period and period derivative for each epoch as derived in Section 3.1 In black, we show a binned light curve.

Fig. C.2. Fit of our toy model for an half ring to the phase-averaged Hα profile for ZTF J2008+4449. The left plots are for optically thin emission, the middle are for optically thick and the right ones are for a mix of optically thin and optically thick emission. The upper panels show the phase-resolved profiles for each model, the middle panels show the best fit phase-averaged compared to the data and the lower panels show the residuals and the reduced χ2 for each fit. It can be seen that the mixed model is a marginally better fit to the data than the optically thin model, and both are much better fits than the fully optically thick model.

Appendix D: Corner plots

We here show the corner plots for different analysis in the paper. The first corner plot in Fig. D.1 shows the results of the SED fitting discussed in Section 3.2. Then we show the corner plots for the X-ray analysis presented in Section 3.5. The first is for the X-ray spectral fit with a two-temperature Apec model (Fig. D.2) and the second is for the power law model (Fig. D.3). Finally, the last corner plot (Fig. D.4) shows the posterior distributions for the fit of the half-ring model to the Hα emission profile presented in Section 4.3.

Fig. D.1. Corner plot for the SED fitting presented in Section 3.2 for ZTF J2008+4449. As both the distance and the radius of the white dwarf act as a normalization factor on the atmosphere models, we fit for the ratio R*/D, where R* is the radius of the star in km and D is the distance to the star in parsecs. We then convert the value obtained to a radius for the white dwarf by employing the parallax measured by Gaia DR3 and factoring in the error on the parallax. The other two parameters, effective temperature and reddening, are fit directly. We employed an uniform prior on all parameters.

Fig. D.2. Corner plot for the X-ray spectral fit presented in Section 3.5 for a two-temperature Apec model. The different colors indicate fits on the two epochs and on the combined spectra; the data show no evidence of variability.

Fig. D.3. Corner plot for the X-ray spectral fit presented in Section 3.5 for a power-law model. The different colors indicate fits on the two epochs and on the combined spectra; the data show no evidence of variability.

Fig. D.4. Posterior distributions for the MCMC fit to the phase-averaged Hα emission profile. The fitting procedure includes five free parameters; 1) inner ring radius, 2) outer ring radius, 3) fractional contribution from optically thin emission, 4) absolute normalization, and 5) system velocity. The posteriors on the inner and outer ring radii yield a constraint on the radial extent of the (partially emitting) ring of δr = 14.8 ± 0.4 RWD. The fitting procedure favored a contribution dominated by optically thin emission. The system velocity was found to be consistent with zero.

All Tables

All Figures

Fig. 1. Gaia CMD for white dwarfs within 100 pc from Earth and within the SDSS footprint (Kilic et al. 2020), where the x-axis depicts the difference between the Gaia BP and RP bands, and the y-axis the absolute magnitude in the Gaia G filter. Solid lines show theoretical cooling tracks (interpolated over the white dwarf’s mass with the WD_models package, Cheng 2024) for white dwarfs with masses between 0.6 M⊙ (top) and 1.28 M⊙ (bottom), equally spaced in mass; the atmosphere is assumed to be hydrogen-dominated (Holberg & Bergeron 2006) and the interior composition to be carbon-oxygen (Fontaine et al. 2001) for M < 1.1 M⊙ and oxygen-neon (Camisassa et al. 2019) for M > 1.1 M⊙. ZTF J2008+4449 is shown as a red star, and its location in the CMD reveals its high mass and high temperature. Reddening corrections were applied only to ZTF J2008+4449; as the objects in the sample are close to the Solar System, reddening is expected to be small. 1σ error bars are smaller than the size of the colored marker, and are omitted for the black background dots for clarity.

In the text

	Fig. 2. Binned CHIMERA light curve with 3-second exposures, phase-folded at the correct period for each epoch as derived in Section 3.1 and normalized to the mean of the light curve in the g′-band (upper) and in the r′-band (lower).
In the text

Fig. 3. Top left: Phase-averaged Keck/LRIS spectrum of ZTF J2008+4449 and (below) the predicted positions of Hα and Hβ absorption lines as a function of magnetic field. The broad and shallow absorption lines, in particular at ≈8500 Å marked in blue, indicate a surface magnetic field B > 400 MG, while the Balmer emission lines suggest ongoing accretion. The red and blue dashed lines mark the maximal Doppler shifted wavelengths of Hα and have the same significance in all panels in this plot. Bottom left: Enlarged view of the phase-averaged Hα emission line, binned (unbinned) according to the instrument resolution shown as the solid black (faded gray) line. Bottom middle: Phase evolution of the trailed spectrum around Hα over the 6.6-minute period. The heat map reflects the mean-subtracted flux. Right: Trailed spectra of the Hα line, with spin phase progressing from bottom to top, as denoted by the numbers on the right hand side of the plot.

In the text

Fig. 4. Upper panel: Phase-averaged COS spectrum in Fλ (solid black line, errors in gray). Strong absorption lines are highlighted. A possible Lyα component is in blue, and non-identified metal lines are in red (≈1092, 1108, 1174, 1190, 1193, 1234, 1253, 1277, 1328 Å; see Fig. 8). Narrow carbon and silicon absorption lines (possibly from the interstellar medium or circumstellar material) are highlighted in cyan. In the lower part, the Lyα component is compared to predicted line positions of Lyα and Lyβ as a function of magnetic field, where dashed lines indicate forbidden components (Ruder et al. 1994). Middle left: COS phase-resolved spectrum (five phase bins) over one period, zoomed in on the region with absorption lines. Each phase is shifted vertically by the same amount. Most of the absorption lines do not change significantly with phase, with the exception of the Lyα line, which disappears for half of the period. Middle right: COS light curves in two wavelength ranges: 1100−1200 Å (blue) and 1250−1280 Å (orange). Bottom two panels: Phase resolved COS spectra plotted over two periods. In the upper panel, the color bar indicates the flux in units of 10−15 erg s−1 cm−2 Å−1. In the lower panel, the phase-resolved spectra have been divided by the phase-averaged spectrum, so the color bar shows the variation with respect to the mean. Two sources of variability can be seen: a strong variation of the continuum between 1100 and 1300 Å, and a variation in the strength of the line at 1344 Å. The absorption is most pronounced when the continuum is at maximum brightness.

In the text

Fig. 5. Data of STIS and COS UV for ZTF J2008+4449. Upper panel: Full phase-averaged UV spectrum in Fν (black solid line, errors in gray). Middle left: Phase-resolved spectrum (five phase bins) in the STIS band over one period. Each phase is shifted vertically by the same amount. There is very little variability with phase in the STIS band. Middle right: Comparison between the COS light curve in the 1100−1200 Å range (blue, same as Fig. 4) and the STIS light curve in the 1800−2250 Å range (red), which is where the STIS spectrum is most variable. The STIS light curve is also roughly sinusoidal, with a lower amplitude than the COS one, and slightly out of phase. Note that the error bars on the STIS light curve are smaller than the symbols used to plot it here. Bottom panel: STIS phase resolved spectra plotted over two periods and divided by the phase-averaged spectrum (the color bar shows the variation with respect to the mean). A weak variability can be observed, mostly in the 1800−2250 Å range.

In the text

	Fig. 6. Binned XMM/OM light curve, shown with fewer (more) bins in blue (gray), phase-folded at the correct period for each epoch as derived in Section 3.1 and normalized to the mean of the light curve in the U-band.
In the text

Fig. 7. Left: ZTF periodogram for ZTF J2008+4449 in the r-band. The main period and first harmonic (half of the period) are highlighted by vertical red dashed lines. The inset shows the complex structure of the main peak, due to a combination of the light curve windowing function and the fact that the period is changing significantly over the ZTF baseline. Right: O–C diagram for ZTF J2008+4449 for two high cadence ZTF observations (in green) and the 8 CHIMERA observations (in blue), all in the r band. On the y-axis, we plot the difference in seconds between the observed time of the minimum in the light curve and the expected time calculated assuming a constant period, while on the x-axis, we plot the elapsed time since our reference time T0 = 59401.76. For the second ZTF observation, the error in the x direction indicates the total duration of the ZTF data used in the fit (70 nights); for the other ones the duration is too small to show (three nights for the first ZTF observation, less than one night for the CHIMERA ones). The quadratic behavior indicates the presence of a significant period derivative.

In the text

Fig. 8. Top: Phase-averaged COS spectrum in Fλ (in black), with errors in light gray. The vertical bands in both panels mark the position of the absorption lines, with the blue one corresponding to the absorption line at 1344 Å, as in the top panel of Fig. 4; Bottom: Positions of a suite of C II–IV (black), Si II–III (purple), and N V (red) ionic transition components as a function of magnetic field strengths. Only transitions with a relatively large oscillator strength at 500 MG have been plotted. It can be seen that the π components of the Zeeman-split lines do not deviate far from the wavelength range where most of the absorption lines in the UV spectrum are located.

In the text

Fig. 9. Top: SED of ZTF J2008+4449 plotting the AB magnitude against wavelength. The light blue, orange, and deep blue data points respectively show the IR (WIRC), optical (PanSTARRS) and NUV (Swift) photometry. The cyan line shows HST STIS and COS concatenated spectra spanning the NUV and FUV, with the COS spectrum truncated at the lowest wavelengths. The solid black and dashed gray lines show the best fitting pure Hydrogen magnetic and non-magnetic atmosphere models respectively. Lastly, a 1700 K (800 K) brown dwarf model is shown in purple (pink), by itself (faded line) and added to the best fitting magnetic white dwarf model (solid line). Bottom: Irradiated face temperatures for a series of Roche-lobe-filling brown dwarf and planetary companions of varying ages against their orbital separation from the white dwarf primary. The upper y-axis shows the orbital period given by each respective orbital separation. The four different tracks correspond to four different brown dwarf or planet ages (0.01, 0.1, 1 and 10 Gyr), marked next to each track. Brown dwarf models (≥13 MJup) are plotted in colors corresponding to their respective masses; all planetary models (≤13 MJup) are plotted in black. The size of each data point is representative of (but not linearly proportional to) the model radius. The golden (black) horizontal dashed line marks the minimal irradiation temperature for a brown dwarf (planet) among the evolutionary models utilized, i.e., ≈1700 K (≈900 K). The irradiated temperatures and corresponding orbital separations were calculated using the Sonora Bobcat brown dwarf evolutionary models grid (Marley et al. 2021).

In the text

Fig. 10. Central wavelengths of the first nine cyclotron harmonics (labeled 1 to 9) as a function of magnetic field (depicted by the black lines). The faded red and orange regions mark the ranges covered by the variable humps seen in STIS and COS, respectively. The dashed lines mark the centroid of each hump. For magnetic fields ranging from ≈400 − 600 MG (faded blue region), the positions of the two UV humps are consistent with the first and second cyclotron harmonics.

In the text

Fig. 11. Doppler tomogram of the Hα emission feature. The minimally assumptive Hα emission density distribution, extrapolated from the emission line variation on the 6.6-minute white dwarf spin period, is shown in the v = Vsin(i) velocity space (i = inclination of the white dwarf rotation axis). The white dwarf is located at the origin at the coordinate system (v, θ). The evolution of the obtained source distribution with rotation phase can be recovered by rotating the coordinate system uniformly and clockwise over one variation cycle. The concentric dotted lines mark regions of constant velocity separated by 1000 km/s from each other, while the colors illustrate the density of emission. The emission distribution resembles a half ring, with velocities ranging between ≈1500 and ≈2500 km/s and peaking at a velocity of ≈1700 km/s.

In the text

Fig. 12. XMM spectral fit. The X-ray spectra are shown as a photon count rate as a function of energy (black data points) for the XMM PN detector (left column) and the combined M1 and M2 detectors (right column). The best fitting two-temperature APEC and power law emission models, convolved with the instrumental response, are shown in blue on the top and bottom rows respectively, with the 3-σ contours shown as shaded regions. In the case of the two-temperature APEC model, the two individual temperature components are shown in green and purple, along with their 3-σ contours. The residuals of the fit models, in both detectors, are shown on the bottom of their respective plots.

In the text

	Fig. 13. Period search in XMM-Newton EPIC data using a χ2 test at a range of trial periods to search for a periodic signal. We show the resulting periodagram for the EPIC-pn and MOS data on the left and right, respectively. The vertical dotted line shows the period of the white dwarf at the epoch of the observation.
In the text

	Fig. 14. Phase-folded light curve for the XMM-Newton EPIC-pn and MOS data. Barycenter-corrected events have been folded on the photometric rotation period, with the ephemeris defined for the epoch of the XMM observations, such that phase zero is consistent with all other plots in the paper (the minimum in the optical light curve).
In the text

Fig. 15. De-reddened Gaia DR3 CMD for RSG 5 (black) compared with the CMDs of three reference clusters: Alpha Persei (81 Myr, orange), IC 2391 (43 Myr, cyan), and NGC 2451A (50 Myr, red). The members of each cluster are those with at least a 50% membership probability in the Hunt & Reffert (2023) Milky Way cluster catalog, with clusters de-reddened using the catalog’s mean Av. Overlaid are PARSEC isochrones (Bressan et al. 2012; Tang et al. 2014; Chen et al. 2014; Fu et al. 2018; Nguyen et al. 2022) for our best-fit cluster age of RSG 5 (45 Myr) as well as the kinematic ages from Heyl et al. (2021) for the comparison clusters, with the colors matching those of the CMDs.

In the text

Fig. 16. Optical and UV light curves for ZTF J2008+4449, phased at the correct period for each epoch as derived in Section 3.1. Phase 0 corresponds to the minimum in the r-band. While for the CHIMERA and XMM/UVOT data the filter is indicated in the upper-left corner of each plot, for the STIS and COS data we extracted light curves by averaging the flux in the wavelength range indicated in the bottom right of each plot.

In the text

Fig. 17. Diagram showing the possible location of the line emitting material. On the x-axis, we plot the distance from the white dwarf (d) in white dwarf radii (RWD ≈ 4800 km). The black solid line indicates the strength of the magnetic field at the distance d, assuming a field of dipole strength at surface of ≈500 MG. We highlight with a dash-dotted vertical line where the field becomes weaker than one MG. The top x-axis marks the corotation velocity corresponding to each radial coordinate. The ratio of centrifugal to gravitational accelerations is unity at the Keplerian corotation radius (vertical dashed line). The vertical blue line indicates the distance at which the Keplerian velocity equals the observed velocity in Hα; the blue shaded region shows where the emitting material would lie if it was moving in Keplerian orbit, depending on inclination (with the vertical blue line indicating the maximum distance for an inclination i = 90°). We think the material is instead corotating with the white dwarf, with the vertical orange line indicating the distance at which the corotation velocity equals the observed velocity in Hα. The orange shaded region shows where we expect the material to lie, depending on inclination (with the vertical orange line indicating the minimum distance for i = 90°).

In the text

	Fig. 18. Comparison of the trailed spectrogram from LRIS for the Hα and Hβ emission of ZTF J2008+4449 (upper-left and lower-left plots, respectively, in blue) with our model for the optically thin (central plots) and optically thick (right plots) emission from a half ring with Rin = 20RWD and Rout = 35RWD.
In the text

	Fig. 19. Left: Observed phase-resolved LRIS spectra for the Hα profile of ZTF J2008+4449, divided into 12 phases (similar to right panel of Fig. 3). Middle: Phase-resolved profile for the best-fit model shown in the upper panel of Fig. 20. Right: Comparison of the observed phase-resolved profile with the model.
In the text

Fig. 20. Top: Five-parameter fit to the Hα emission profile employing our half-ring model and a mix of optically thick and optically thin emission, same as right-middle and right-bottom panels of Fig. C.2. The upper panel shows the comparison of our model (in red) to the phase-resolved Hα profile observed with LRIS for ZTF J2008+4449 (in blue), while the lower panel shows the residuals. Middle: One-parameter fit to the Hβ emission profile for ZTF J2008+4449, fixing all parameters to those found for Hα (see Fig. D.4), except the normalization, which is the only fit parameter. Bottom: Line ratios of Hα/Hβ as presented by Storey & Hummer (1995) for a range of gas temperatures and electron densities, shown with points and connecting lines. The measured Hα/Hβ ratio from the emission line fitting (top two panels) is shown with the horizontal dotted blue line. The 1-, 2-, and 3-σ contours are shown in the filled horizontal bands. The simulation data (points and connected lines) most closely matching the inferred electron density from our photoionization ring model are plotted with greater thickness.

In the text

Fig. 21. Left: Continuum-normalized spectrum near the rest wavelength of Hα. Before the continuum normalization, the fit emission feature (Fig. C.2) was subtracted. The dashed line shows a Gaussian fit to the data. The fitting procedure yields an equivalent width of EW = 8.4 ± 2.5 km s−1 and central velocity of v0 = −37 ± 13 km s−1. Middle and right: Curves of growth for Hα. Middle: Evolution of the equivalent width as a function of optical depth. Right: Similar to middle panel, but as a function of column density. The dotted horizontal line indicates the measured equivalent width from the left panel. The vertical dashed line indicates the optical depth and column density inferred from the photoionization model.

In the text

	Fig. C.1. CHIMERA light curve for ZTF J2008+4449. In blue, we show the unbinned data for the g (upper) and r (lower) filters, phase-folded using the correct period and period derivative for each epoch as derived in Section 3.1 In black, we show a binned light curve.
In the text

Fig. C.2. Fit of our toy model for an half ring to the phase-averaged Hα profile for ZTF J2008+4449. The left plots are for optically thin emission, the middle are for optically thick and the right ones are for a mix of optically thin and optically thick emission. The upper panels show the phase-resolved profiles for each model, the middle panels show the best fit phase-averaged compared to the data and the lower panels show the residuals and the reduced χ2 for each fit. It can be seen that the mixed model is a marginally better fit to the data than the optically thin model, and both are much better fits than the fully optically thick model.

In the text

Fig. D.1. Corner plot for the SED fitting presented in Section 3.2 for ZTF J2008+4449. As both the distance and the radius of the white dwarf act as a normalization factor on the atmosphere models, we fit for the ratio R*/D, where R* is the radius of the star in km and D is the distance to the star in parsecs. We then convert the value obtained to a radius for the white dwarf by employing the parallax measured by Gaia DR3 and factoring in the error on the parallax. The other two parameters, effective temperature and reddening, are fit directly. We employed an uniform prior on all parameters.

In the text

	Fig. D.2. Corner plot for the X-ray spectral fit presented in Section 3.5 for a two-temperature Apec model. The different colors indicate fits on the two epochs and on the combined spectra; the data show no evidence of variability.
In the text

	Fig. D.3. Corner plot for the X-ray spectral fit presented in Section 3.5 for a power-law model. The different colors indicate fits on the two epochs and on the combined spectra; the data show no evidence of variability.
In the text

Fig. D.4. Posterior distributions for the MCMC fit to the phase-averaged Hα emission profile. The fitting procedure includes five free parameters; 1) inner ring radius, 2) outer ring radius, 3) fractional contribution from optically thin emission, 4) absolute normalization, and 5) system velocity. The posteriors on the inner and outer ring radii yield a constraint on the radial extent of the (partially emitting) ring of δr = 14.8 ± 0.4 RWD. The fitting procedure favored a contribution dominated by optically thin emission. The system velocity was found to be consistent with zero.

In the text