© The Authors 2024

1 Introduction

Close white dwarf binaries (CWDBs), consisting of at least one white dwarf in the system, are possible descendants of common envelope evolution (e.g., Gokhale et al. 2007; Anguiano et al. 2022). They are categorized as one type of post-common envelope binary (Gokhale et al. 2007). Due to the complex mass loss mechanism and short timescale, the progenitors of CWDBs – common envelope binaries (CE) – lack well-measured stellar parameters, and thus lead to insufficient theoretical interpretation. The observational constraints from post-common envelope binaries greatly illuminate the studies of CE (Han et al. 2020).

On the other hand, the detection of CWDBs may also identify the progenitors of many exotic objects, such as over-luminous type IA supernovae.

When both components of CWDBs are white dwarfs, they are named double white dwarfs (DWDs), which are the main source of galactic gravitational wave foreground radiation in the frequency ranging from 10⁻⁴ to 0.1 Hz. A number of different population models predict that ~10⁴ of close DWDs (0.1% of Galactic ultra-compact binaries) will be resolved by the Laser Interferometer Space Antenna (LISA, Georgousi et al. 2022; Amaro-Seoane et al. 2017). But until the year 2022, only 25 DWDs had well-measured orbital parameters and mass ratios (Kilic et al. 2020a), and only 150 DWDs had known orbital parameters (Korol et al. 2022).

Many CWDB candidates have been discovered by time-domain photometric surveys, such as the Kepler K2 mission (Howell et al. 2014; Hallakoun et al. 2016; Van Sluijs & Van Eylen 2018) and the Zwicky Transient Facility (ZTF; see Keller et al. 2021; Burdge et al. 2020b,a, 2019). Ren et al. (2023) used photometric data from Gaia and ZTF to identify a group of CWDB candidates. Their orbital periods and classifications were determined by fitting light curves. Their results provide valuable samples for further spectroscopic study. In the near future, the upcoming Chinese Space Station Telescope (CSST) will detect faint sources down to ~26–27 mag in multiple pass-bands (ultraviolet to optical), greatly expanding the samples of CWDB candidates.

Follow-up confirmation of CWDB candidates generally relies on high-resolution spectra (HRS, R ≳ 20 000). In this case, the double lines from both components are mostly resolved, and thus precise stellar parameters and chemical abundances can be measured. For example, Bédard et al. (2017) demonstrated that it is possible to derive atmospheric parameters of DA+DA binaries by combining spectroscopic (R ~ 18 500), photometric, and astrometric measurements¹. Their work broke the parameter degeneracies by fitting Balmer lines and spectral energy distributions (SEDs) simultaneously. They identified 15 double degenerate binary candidates out of 219 CWDBs. In a follow-up study, Kilic et al. (2020a) performed multi-epoch high resolution (R ~ 18 000 for the Gemini Multi-Object Spectrograph, R ~ 37 000 for the Keck I HIRES echelle spectrograph) spectral analysis on four double degenerate binary candidates. They first determined orbital parameters of each binary with radial velocity (RV) measurements, and then fit the Balmer lines plus the SED to obtain effective temperatures and surface gravities. Finally, the white dwarf (WD) masses were determined through evolutionary tracks.

However, the high-resolution spectrum is so expensive that the resulting sample size is generally small. To fully embrace the arrival of large surveys with a lower spectral resolution, such as the Sloan Digital Sky Survey V (SDSS-V, Kollmeier et al. 2019) and the Large Sky Area Multi-Object Fiber Spectroscopic Telescope (LAMOST, Zhao et al. 2012), a few attempts have been made. It is challenging to decompose the flux contributions from both components since the line-blending effect in low-resolution spectra (LRS) is severe. The parameter degeneracy is attenuated in CWDBs with different spectral type components, because of the significant SED discrepancy between the two components. For example, the white dwarf + dwarf M (WD+dM) binary can be decomposed by fitting the red and blue ends of the spectrum individually, as WDs dominate the flux at the blue end, while M dwarfs dominate the red end. For example, Ren et al. (2013) utilized principal component analysis and template-matching techniques to obtain stellar parameters for ten WD+dM binaries based on LAMOST LRS.

During the process of estimating the stellar parameters of binary systems, two strategies are usually adopted to generate the model synthetic spectra: template-matching (e.g., Ren et al. 2013) and real-time computing (e.g., iSPEC, Blanco-Cuaresma 2019). Matching spectra with a pre-computed template grid limits the precision of the template resolution, while resolving the radiation transfer equation in real time is time-consuming. A trade-off method is to apply a machine-learning method or polynomial interpolation to parameterize the synthetic spectra. With these techniques, synthetic spectra can be generated faster and their parameter coverage is more continuous. Artificial neural network (ANN) parameterization is a popular algorithm, which employs codes such as PAYNE (Ting et al. 2019; Xiang et al. 2022) and wdtools (Chandra et al. 2020). The ANN behaves better than 21-term polynomial interpolation (e.g., Sharma et al. 2020, TGM2 in Ulyss).

Compared to HRS, LRS generally have the advantages of a wider wavelength range, a higher signal-to-noise ratio (S/N), and a shorter exposure time. A short exposure time is particularly important when measuring the RV of close binaries with short periods (of the order of minutes to hours), given that long exposure times can smear the RV signal completely (Rebassa-Mansergas et al. 2019). The information loss due to the line-blending effect is compensated by a wider wavelength range coverage and higher S/N. If we have robust models and the ability to fit all stellar labels simultaneously, the uncertainty should be entirely independent of the resolution. This holds over spectral resolutions from R ~ 1000 to R ~ 100 000 (Ting et al. 2017). At the same time, the large number of unresolved binary candidates with LRS call for a robust pipeline to determine the parameters of both components.

This work aims to develop a pipeline to derive the stellar parameters of both components of CWDBs based on LRS. The inputs of our code include: flux-calibrated spectra and parallax. The outputs include: the atmospheric parameters (T_eff, log 𝑔 for WDs, T_eff, log 𝑔, and [Fe/H] for MS), mass, and radius of each component. Our code is innovative compared to previous studies; for example, (1) implementing the binary decomposition technique of Kilic et al. (2020a) in LRS for the first time, (2) using the flux-calibrated spectra, instead of SED plus spectra lines (the latter requires multi-epoch observations), and (3) applying an ANN in binary decomposition for the first time, significantly improving the efficiency and accuracy of generated spectra. We structure this paper as follows. In Sect. 2, we outline the procedures of our code and show how to generate model spectra using an ANN. In Sect. 3, we show how to determine stellar parameters by simultaneously fitting line profiles and continuum. After testing the credibility of our code, we apply it to the SDSS LRS of 14 WDB candidates in Sect. 4. In Sect. 5, we discuss our results and compare them to previous works. In Sect. 6, we briefly summarize this work and give a future perspective.

2 Generating model spectra

The procedures of our code are the following: (1) For a given CWDB, we first evaluated the possible stellar combinations according to its location in the color magnitude diagram (CMD) (see figures of Inight et al. 2021). (2) We generated synthetic spectra of both stellar components based on known spectral libraries. In this work, we used the ANN stellar spectrum generator. (3) We then fit the flux-calibrated spectra to derive the atmospheric parameters of both components. (4) We finally estimated the mass and radius by interpolating between evolutionary tracks. Given that our code assumes a combination of two stellar components, it is not suitable to fit binary systems with significant emission lines, which indicate prominent binary interaction. We first describe how to generate model spectra using the ANN stellar spectrum generator in this section².

Fig. 1 Conceptual illustration of the ANN architecture. X1 is the number of input parameters, while X2 is the number of output pixels.

2.1 Spectral grids

Main-sequence stars: The Coelho synthetic stellar library (Coelho 2014) is a high-resolution (constant wavelength sampling of ∆λ = 0.02 Å) theoretical spectral database for main-sequence stars. This library covers B-type to M-type main-sequence stars in the following parameterspace: $$\begin{array}{c}2500Å\le \lambda \le 9000Å\hfill \\ 3000\text{\hspace{0.17em}K}\le T_{\text{eff\hspace{0.17em}}}\le 25000\text{\hspace{0.17em}K}\hfill \\ 4.0\le \log g\le 5.5\hfill \\ -1.3\le [\text{Fe}/\text{H}]\le 0.2\hfill \\ [\alpha /\text{Fe}]=0.0,\mathrm{0.4.}\hfill \end{array}$$

White dwarfs: Among WDs, there are ~80% DA-type WDs (DAs)³ and ~20% DB/O-type WDs (DB/O)⁴ (Sion et al. 1983; Kong & Luo 2019). Thus, we have ignored other type of WDs in this work. Here, we used the grid of pure-H and pure-He atmospheres from Koester (2010).

The pure-H templates⁵ cover the following parameter space: $$\begin{array}{l}900Å\le \lambda \le 30000Å\hfill \\ 5000\text{\hspace{0.17em}K}\le T_{\text{eff\hspace{0.17em}}}\le 80000\text{\hspace{0.17em}K}\hfill \\ 6.5\le \log g\le \mathrm{9.5.}\hfill \end{array}$$

The pure-He templates cover the following parameter space: $$\begin{array}{l}3000Å\le \lambda \le 9000Å\hfill \\ 8000\text{\hspace{0.17em}K}\le T_{\text{eff\hspace{0.17em}}}\le 40000\text{\hspace{0.17em}K}\hfill \\ 7.0\le \log g\le \mathrm{9.5.}\hfill \end{array}$$

2.2 Artificial neural network stellar spectrum generator

Inspired by the work of Ting et al. (2019), we adopted an ANN architecture. It is fully connected with three hidden layers, each consisting of 256 neurons (Fig. 1). The normalized flux (f_λ), as a function of the stellar parameter label (l), can be written as follows: $$f_{\lambda }=w_{\lambda }^{i}A\left(w_{\lambda }^{j}A\left(w_{\lambda ij}^k{l}_k+b_{\lambda ij}\right)+b_{\lambda i}\right)+{\overline{f}}_{\lambda },$$(1)

where λ represents the index of the output wavelength pixels, and indices k, j, i represent the neuron orders of the first, second, and third ANN layers, respectively. A represents the activation function of a neuron, where a leaky ReLU activation function (Eq. (2)) is employed. $$A(x)=\{\begin{array}{cc}x\hfill & x>0\hfill \\ \alpha x\hfill & x\le 0,\alpha =0.1\hfill \end{array}.$$(2)

The optimizer seeks the coefficients $\left(w_{\lambda }^i,w_{\lambda i}^j,w_{\lambda ij}^k,{\overline{f}}_{\lambda },b_{\lambda i},b_{\lambda ij}\right)$ that best fit the training spectra as a function of labels.

Next, we randomly selected a set of parameters as inputs for the spectrum generator. The resulting flux-calibrated spectrum was compared with theoretical spectrum in the spectral library with the same parameters. Our spectrum generator performs smoothly across most of the parameter space with high accuracy (See Appendix A for details).

One of our major improvements is the capability to generate flux-calibrated spectra, which are essential for the following fitting procedures. Previous studies always applied continuum normalization to the synthetic spectra before training (e.g., Chandra et al. 2020; Ting et al. 2019), which fails to generate flux-calibrated spectra. The technique of continuum-normalization can be described as $$F_{\text{norm\hspace{0.17em}}}(\lambda )=\frac{F(\lambda )-F_{\min }}{F_{\max }-F_{\min }},$$(3)

where F_max and F_min represent the minimum and maximum flux values of a given spectrum, respectively. Without this normalization, the neural network might struggle to achieve accurate interpolation results for some stellar parameters. For example, surface gravity mainly affects the line widths, which have little contribution to the mean square error (MSE) loss function in the neural network. However, the generated spectra after such flux normalization are also scaled, which is not suitable for our code.

We used a new normalization formula instead of Eq. (3): $$F_{\text{norm,new\hspace{0.17em}}}(\lambda )=\frac{F(\lambda )-F_{\lambda ,\min }}{F_{\lambda ,\max }-F_{\lambda ,\min }},$$(4)

where F_{λ, mаin} (F_{λ, max}) represents the minimum (the maximum) flux value among the entire spectral template sets at wavelength λ.

The flux values of any single parameter set share the same scaling factor (F_{λ, min}, F_{λ, max}) at wavelength λ. This allows us to save the “minimum” and “maximum” scaling spectra, F_{λ, min} and F_{λ, max}. We can restore the flux-calibrated spectra by reversing the scaling.

3 Determining stellar parameters

3.1 Spectral preparation before fitting

In this work, we intend to use LRS as the input, but the blended line profiles are insufficient to decompose both components. Inspired by the work of Kilic et al. (2020a), we split the LRS into line segments and continuum segments, and prepared them separately before fitting. We note that the continuum segments are defined as the spectral regions outside of the line segments. In total, 234 frequently used feature lines (of 38 types of atoms or ions) were selected from the National Institute of Standards and Technology (NIST) Atomic Spectra Database⁶.

The RV, caused by the Doppler effect or the gravitational effect, should be determined prior to fitting stellar parameters, given that the RV changes the observed wavelengths of the lines. The most persistent feature lines, such as Na D for M stars, Fe I for F-K stars, Hα for A stars, and He for O-B stars (Hilditch 2001), were used to determine the RVs.

We prepared the line segments and the continuum segments separately. For the former, we used the normalized spectra to maximize the changes in the line profiles. For the latter, we used flux-calibrated spectra. The model spectra can be expressed as $$F_{\text{line\hspace{0.17em}}}(\lambda )=\frac{R_1^2{S}_1(T,\log g,\lambda )+R_2^2{S}_2(T,\log g,\lambda )}{D^{2}P(\lambda )},$$(5) $$F_{\text{cont\hspace{0.17em}}}(\lambda )={\left(\frac{R_1}{D}\right)}^2{S}_1(T,\log g,\lambda )+{\left(\frac{R_2}{D}\right)}^2{S}_2(T,\log g,\lambda ),$$(6)

where R₁ and R₂ are the radii determined by the evolution tracks of each stellar component (i.e., R = R(M(T, log 𝑔), log 𝑔), see Sect. 3.3 for details), D denotes the distance, S₁, S₂ refer to the synthetic spectra of both components, and P(λ) represents the pseudo-continuum. Given that the distance, D, can be derived from trigonometric parallax measurements, the right-hand sides of these equations only depend on the stellar parameters of both components. The observational spectra were also divided into the corresponding continuum and normalized line segments: $$F_{\text{line,observed\hspace{0.17em}}}(\lambda )=\frac{O(\lambda )}{P(\lambda )},F_{\text{cont,observed\hspace{0.17em}}}(\lambda )=O(\lambda ),$$(7)

where O(λ) are the observed spectra and P(λ) is the pseudo-continuum.

3.2 Atmospheric parameters

A first guess of the parameters is preferred before the fitting, since it makes the later optimization convergence faster. Zhang et al. (2019) statistically confirmed the hypothesis that the temperatures of binary spectra derived by Ulyss⁷ are approximately equal to those of the luminosity-dominant component. Otherwise, if their flux contributions are similar (i.e., not dominated by any components), they probably have a similar temperature. In addition, the metallicities of both components are thought to be equal if they come from the same cloud.

The atmospheric parameters were obtained by simultaneously fitting the line segments and the continuum segments. We used the “minimize” optimizer from the SciPy package in this work. Due to the discrepancy in pixel numbers and flux values between F_line(λ) and F_cont(λ), one cannot simply add up the residuals of the two segments. We solved this problem by adjusting the weights of the residuals of F_cont(λ) and F_line(λ). Generally, we decreased the weight of the F_cont(λ), so that the line segments and the continuum segments could be fitted equally well. The predicted spectrum is expressed as Y_pređict = XL, where L represents the input parameters and X are the spectral models. The matrix of the squared residuals can be written as $$R(L)=\frac{1}{2}{\left(Y_{\text{predict\hspace{0.17em}}}-Y\right)}^T\left(Y_{\text{predict\hspace{0.17em}}}-Y\right)=\frac{1}{2}{(XL-Y)}^T(XL-Y).$$(8)

where Y is the observed spectrum. To represent the weights of different segments, we introduced a transformation matrix, α, to the residual: $$(\mathit{X}\mathit{L}-\mathit{Y})\to \mathit{\alpha }(\mathit{X}\mathit{L}-\mathit{Y}).$$(9)

The new squared residual, $\tilde{R}(L)$, is expressed as $$\tilde{R}(L)=\frac{1}{2}{(XL-Y)}^T{\alpha }^T\alpha (XL-Y).$$(10)

The minimal squared residual is located where its gradient equals zero. The gradient is written as $$\begin{array}{ll}\frac{\partial }{\partial L}\tilde{R}(L)\hfill & =X^T{\alpha }^T\alpha (XL-Y)\hfill \\ & =X^{T}A(XL-Y),\hfill \end{array}$$(11)

where each element of the transformation matrix, A = α^T α, can be regarded as the weight of each wavelength pixel. The transformation matrix, α, is set to be $$\alpha =\left(\begin{array}{cccccc}\frac{1}{F_{\text{med\hspace{0.17em}}}}& 0& 0& \cdots & 0& 0\\ 0& \frac{1}{F_{\text{med\hspace{0.17em}}}}& 0& \cdots & 0& 0\\ ⋮& & \ddots & & & \\ 0& \cdots & 0& 1& \cdots & 0\\ 0& \cdots & 0& 0& \ddots & 0\\ 0& \cdots & \cdots & \cdots & 0& \frac{1}{F_{\text{med\hspace{0.17em}}}}\end{array}\right).$$(12)

This is a diagonal matrix: elements corresponding to the continuum segments were set to 1/F_med (F_med is the median flux of the continuum segments), while elements corresponding to the line segments were set to 1. This transformation reduces the weights of the continuum segment residuals, leading to more reasonable fitting results.

While looking for the minimal squared residual, Eq. (11) indicates that the model, X, plays an important role. The spectral model, X, affects the predicted feature line strength (part of spectrum Y_predict at a given parameter, L. For example, Fig. 1 of Xiang et al. (2019) shows the derivatives of stellar spectra with respect to 24 parameters. In other words, some feature lines may be particularly strong (or weak) in one type of stars, which helps to break the degeneracy. In that sense, binaries with significantly different atmospheric parameters are expected to be better decomposed (and vice versa; see Sect. 4 for fitting results).

3.3 Stellar mass and radius

Assuming that the evolutionary tracks do not overlap, we can determine stellar mass on a three-dimensional grid (log(T_eff), log 𝑔, [Fe/H]). Additionally, stellar mass is used to determine the radius of the star and then constrain the model flux through Eq. (5). Therefore, the stellar radius is not a free parameter in our code, reducing the chance of non-physical results from Markov chain Monte Carlo fitting. In this work, we utilized the linear N-dimensional interpolator package (LinearNDInterpolator)⁸. The white dwarf evolutionary tracks were obtained from the Python package WD_models (Cheng et al. 2020)⁹. The MS stellar evolutionary tracks were picked from MIST¹⁰, which covers the following parameter space: $$\begin{array}{l}0.1\le M_{\text{initial\hspace{0.17em}}}/M_{\odot }\le 10\hfill \\ -1.5\le [\text{Fe}/\text{H}]\le 0.5\hfill \\ [\alpha /\text{Fe}]=0\hfill \\ V/V_{\text{crit\hspace{0.17em}}}=0,0.4,\hfill \end{array}$$

where M_initial is the initial stellar mass and V_crit is the critical linear velocity at the surface (Choi et al. 2016).

We confirmed the reliability of the interpolation function by reproducing the mass-luminosity and mass-radius relations (Huang 2012). Next, we calculated the mass uncertainty numerically using the gradients of the evolutionary tracks, and the radius uncertainty was derived through uncertainty propagation.

4 Verification and results

Our code can be applied to LRS with various resolutions. In this work, we chose SEGUE LRS from the Sloan Digital Sky Survey (Gunn et al. 2006, SDSS). These spectra cover a wavelength range of 3800–9200 Å with a spectral resolution of R ~ 1800–2200.

4.1 Test sample

To validate our code, we selected two CWDBs with accurate stellar parameter estimation as the reference sample (Table 1). They show distinct double lines in HRS; moreover, time-domain measurements, such as light curves and RVs, are available for these sources. The test sample size is small as the aforementioned high-quality datasets are challenging to obtain.

The first CWDB, WD 1534+503, was identified as a binary system with DWDs. Kilic et al. (2021) obtained the orbital and atmospheric parameters of WD 1534+503 with high confidence (Table 1). In this work, our prediction successfully reproduce the flux-calibrated LRS, given that most of the residuals are within the 3σ level. In particular, the Balmer lines are distinguishable in the residual spectrum, indicating an accurate fitting (the bottom panel of Fig. 2). The fitted parameters of the primary star and the fitted effective temperature of the secondary star are consistent with the literature, while the surface gravity of the secondary star is not (see Fig. 2). The latter also leads to a significant mass-radius discrepancy (see Fig. 3).

The second CWDB, PG 1224+309, is a DA+dM4 binary according to Orosz et al. (1999). In their study, the masses of both components were obtained from the orbital parameters of this binary system. Given the low luminosity of the M dwarf, its temperature was fixed to 3100 K during the fitting (Table 1).

In this work, our prediction generally reproduces the flux-calibrated LRS, since most of the residuals are within the 3σ level (see Fig. 4). We also notice two features in the residual spectrum: (1) distinct features near 7000–8000 Å, probably due to the inaccurate modeling of the molecular lines in M dwarfs; (2) a strong Hα emission line, which is also visible in the observed spectrum. This is not surprising, given that PG 1224+309 is a short-period eclipsing binary. These defects may be the reason for the inconsistency in the stellar parameters of the secondary star (Table 1).

Through these two CWDBs with high-quality datasets, we prove that our code works reasonably well for LRS with R ~ 2000: (1) the stellar parameters of the primary stars are accurate; (2) given the lower luminosities of the secondary stars, their derived atmospheric parameters may be slightly less accurate, leading to larger errors in the estimated masses and radii.

4.2 Results of 14 white dwarf binary candidates

After estimating the reliability of our code on LRS, we plan to apply it to photometric identified CWDB candidates with LRS observations. As was mentioned in Sect. 1, Ren et al. (2023) have provided a catalog of 423 short-period eclipsing binary candidates. Among them, 28 sources have SDSS LRS. We further excluded: (1) spectra with S/N below 10; (2) emission-line dominated spectra; (3) previously identified subdwarf+M binaries (J162256.66+473051.1, J011339.09+225738.9, and J082053.53+000843.4, Dai et al. 2022; Devarapalli et al. 2022). We treated the remaining 14 sources as our final sample (Table B.1).

As was mentioned in Sect. 3.2, we first estimated the possible stellar combinations based on their CMD locations (Fig. 5). Our sample sources are located between the MS and WD cooling sequence. Based on their colors, the possible primary stars should be MSs or WDs earlier than G-type stars Zhang et al. (2019), while the secondary stars have no strong constraint on their spectral type. We further reduced the searching parameter space using the results of Inight et al. (2021), which show the possible CMD location of CVs, DWDs, WD+AFGK, and WD+M (Fig. 5).

For each possible stellar combination, we generated the spectra for both components, and added up their fluxes according to Eqs. (6) and (5). The most probable stellar combination was found by fitting the observed spectra (Eq. (3.2)). Though our sources show detectable luminosity variation, the secondary star may not be significant enough to be decomposed from single-epoch LRS (too faint or blocked by the primary star). In this work, we consider a source to be a single star if the fitted temperature of the secondary star shows a large uncertainty (σ_T > 30% T). The fitted parameters and spectra of our sample are shown in Table B.1 and Fig. B.1, respectively.

Fig. 2 Fitting results for WD 1534+503. The top panel shows the SDSS spectrum (black line), the convolved smoothed continuum (dashed orange line), the pseudo-continuum obtained by polynomial fitting (dashed red line), the best-fit spectrum of Component 1 (blue line), the best-fit spectrum of Component 2 (solid orange line), and the composite spectrum (green line). The bottom panel shows the residual spectrum. The ±3σ levels are indicated by the dashed lines.

Fig. 3 Kiel diagram of three DAs in the test sample. The evolutionary tracks provided by Cheng et al. (2020) are shown as small dots, where their masses are indicated by the color bar. The red symbols with error bars represent the expected values and uncertainties of our fitting results. Similarly, orange symbols and green symbols show the literature results of Kilic et al. (2021) and Orosz et al. (1999), respectively. Obviously, the error in log 𝑔 is propagated to the stellar mass, and thus the radius. As the surface gravity is expressed on a log scale, and the mass and radius are expressed on a linear scale, the mass-radius errors are magnified.

5 Discussion

5.1 Residual spectra

The residual spectra reveal possible features of interacting binaries (accretion disk, Roche lobe overflow, stellar wind, etc.). A few sources exhibit multiple metal or molecular lines around 8000-9000 Å (e.g., J034137.68+011027.6, J082145.28+455923.4, J110045.15+521043.8, and J093934.56+531751.5) in the residual spectra. These DAZ-like features suggest an accretion of tidally disrupted asteroids and planets (Bédard 2024).

5.2 Statistics

Though it has a limited sample size, it is interesting to discuss our fitting results in the context of other WD+MS studies.

Our analysis identifies one MS+DB, five DA+MS, one DA, two MSs, and five MS+MS. We note that only six out of 14 CWDB candidates have been confirmed, proving the necessity of LRS follow-up observations. Our DA mean mass is 0.509 M_⊙, which is consistent with the peak value of WD mass distribution in Ren et al. (2018), while our MS mean mass is 0.77 M_⊙. We further divided DA+MS into DA+dM and DA+AFGK systems, to compare them with literature studies. For the two DA+dM binaries in our samples, the mean WD temperature is 54 612 K, which is much higher than those reported in two previous studies (~ 15 000 K, Nayak et al. 2023; Ren et al. 2018). This discrepancy is caused by their different colors: the targets in our work are bluer (G_bp – G_rp < 1) compared to the aforementioned studies (G_bp – G_rp > 1). Binaries with hotter primary components tend to have bluer colors. As for the three DA+AFGK binaries in our sample, the mean temperature of DA and MS agrees with the peak value of the DA+AFGK temperature distribution given by (Ren et al. 2020).

Fig. 4 Fitting results for PG 1224+309. Symbols are the same as Fig. 2.

Fig. 5 Our final samples on the Gaia CMD (colored symbols). Different types of binary systems identified by Ren et al. (2023) are denoted by different colored symbols. The black dots denote Gaia samples within 100 pc with small errors. The region outlined by cadet dashed blue lines was set by Ren et al. (2023) to select their candidates. The color-shaded regions indicate the locations of gold samples (CVs-green, DWDs-blue, WD+AFGK-orange, WD+M-red) compiled by Inight et al. (2021).

5.3 Individual analysis

ZTF J110045.15+521043.80 is identified as a hot DA WD accompanied by an invisible star. The WD dominates almost all the luminosity. Residual analysis reveals the presence of metal and OH I emission lines in the spectral region of >8000 Å, suggesting a potential accretion of surrounding dust and cool dense molecular clouds. Our fitting results are consistent with those of Kepler et al. (2019), T = 23 083 ± 208K, log 𝑔 = 7.525 ± 0.028.

ZTF J140847.17+295044.87 is identified as a DA+dM system in this work. Our fitting WD temperature (T_DA = 31118K) is consistent with the literature (T_DA = 30 000K in Morgan et al. 2012 and T_WD = 29 050 ± 484 K in Parsons et al. 2013), while the fitting WD surface gravity and mass (log 𝑔_DA = 8.058, M_DA = 0.68 M_⊙) are slightly higher than these studies (log 𝑔_DA ~ 7.75, M_DA ~ 0.5 M_⊙). Its spectrum displays narrow Balmer emission lines, similar to PG 1224+309 (Sect. 4), indicating that the near side of the MS star is heated by the WD (Sing et al. 2004).

ZTF J082145.29+455923.46 is classified as a DA+MS system in this work. Our classification is consistent with Rebassa-Mansergas et al. (2010) and Silvestri et al. (2007). A single WD spectroscopic fitting suggests that T = 30 384 ± 152K and log 𝑔 = 7.980 ± 0.034 (Kepler et al. 2019). On the other hand, Parsons et al. (2013) classified this source as a WD+M2 binary system using light curve fitting. Their derived WD temperature is T_WD = 80 938 ± 4024K, with M_WD = 0.66 ± 0.05 M_⊙, which are consistent with our prediction. The discrepancy between our measurements and those of the former study is possibly caused by (1) the orbital phase variation or (2) the defect of single component fitting.

ZTF J034137.67+011027.88 is identified as an MS star with an invisible companion in this work. The spectrum lacks prominent Balmer lines, suggesting there is no strong WD contribution. The spectrum shows metal emission lines (e.g., Ca I, Na D, Ca II, and Fe I) at the red end, indicating a binary accretion.

ZTF J143257.08+491143.01 is classified as an MS+MS system in this work. Single MS spectroscopic fitting suggests that T_eff = 8419 ± 31 K, log 𝑔 = 4.2 ± 0.03 and [Fe/H] = −1.03 ± 0.05 (Chen et al. 2020), which is consistent with our prediction of the primary star.

ZTF J140118.80-081723.51 is classified as a DA+MS system in this work. We note that this source was identified as a different type of system in other studies: an extremely low mass (ELM) WD (Brown et al. 2020) or an sdA (Kepler et al. 2019). The former study derived T = 8813 ± 90 K, log 𝑔 = 5.731 ± 0.048 by fitting its spectrum with a grid of pure hydrogen-atmosphere models. The masses of both components were determined by RV measurements: M_WD = 0.216 ± 0.042 M_⊙ is M_{2, min}¹¹ = 0.79 M_⊙. Our fitted results do not agree with them. The residual spectrum (Fig. B.1) shows non-negligible Balmer absorption lines, which may be attributed to the lack of the ELM or sub-dwarf models in our work.

The remaining sources have undergone no previous spectro-scopic analysis, but some of them are listed in several catalogs (Heinze et al. 2018; Marsh et al. 2017; Drake et al. 2014; Kepler et al. 2019). Readers are referred to these works for more details.

To summarize, we find that our code performs reasonably well when comparing to other studies statistically or comparing fitting results individually. Our fitted residual spectra are also promising to study mild binary activities. We again recognize the disadvantage of our code: the lack of stellar models out of binary evolution; for example, involving subdwarfs and ELM WDs.

6 Conclusion

Close white dwarf binaries play a significant role in understanding the evolution of stellar post-common-envelope phases, and they are also progenitors of several gravitational-wave-emitting objects. The follow-up identification of CWDBs often relies on HRS. However, HRS are too expensive and insufficient to satisfy the growing demands. On the other hand, LRS, like those from SDSS and LAMOST, are more accessible. The LRS also have the advantages of a wider wavelength range, higher S/N, and shorter exposure time. In particular, upcoming photometric surveys (e.g., CSST, LSST, Euclid) will provide a large sample of faint CWDB candidates beyond the magnitude limit of HRS, even with a 10-meter telescope (e.g., Keck, GTC). The most feasible option is degrading the spectral resolution to achieve enough S/N. Therefore, decoupling CWDB components from LRS is a promising task.

In this work, we develop a new code for fitting LRS (R ~ 2000) of CWDB candidates. Our code estimates the atmospheric parameters and masses of both components based on flux-calibrated spectra and parallax. We first generated the spectra of WDs and MS stars using an ANN. We then obtained the stellar parameters by fitting the feature lines and continuum simultaneously. Our code has several advantages: (1) We applied the binary decomposition method of Kilic et al. (2020b) in LRS for the first time. (2) We used flux-calibrated spectra, rather than line profile plus SED data for binary decomposition. This not only avoids the phase variation between spectra and SED observations, but also increases the accessibility of the required data. (3) The ANN spectrum generator has the advantages of the template-matching technique and the synthetic spectra computing technique, which can generate the model spectra both efficiently and accurately.

By comparing to two CWDBs with well-measured components, we have proven that our code performs well when estimating the stellar parameters of the primary stars. The estimation is less accurate for the secondary stars, given their lower luminosities.

We then applied our code to decompose the stellar parameters of 14 white dwarf binary candidates suggested by Ren et al. (2023). We find that only six candidates are CWDBs, proving the necessity of LRS spectroscopic decomposition. The fitting results from our code agree with others from the literature statistically. The individual fitting results are generally consistent with other studies. Two caveats are noted with regard to the implementation of our code: (1) our code is not suitable to fit binary systems with significant emission lines; (2) our code lacks stellar models outside of binary evolution, e.g., subdwarf and ELM ones. These should be the goals of future improvements.

To sum up, our code offers a fully automated approach to identify the most probable stellar combinations and atmospheric parameters of both components. It is a valuable tool for conducting pilot spectroscopic study on CWDB candidates rapidly.

Given the large number of LRS available in current and future spectroscopic studies (e.g., LAMOST-II, DESI), our code will provide a large sample of CWDBs for studying stellar evolution and monitoring gravitational waves.

Acknowledgements

We thank Jingkun Zhao for helpful discussions. G.L., B.T., and C.L. gratefully acknowledge support from the China Manned Space Project under grant NOs. CMS-CSST-2021-B03, CMS-CSST-2021-A08, etc., the National Natural Science Foundation of China through grants No. 12233013, and the Natural Science Foundation of Guangdong Province under grant No. 2022A1515010732. J.N.F. acknowledges the support from the National Natural Science Foundation of China (NSFC) through the grants 12090040 and 12090042, and the science research grants from the China Manned Space Project. This work utilizes data from the Gaia satellite, a space mission of the European Space Agency (ESA). The Gaia data is provided by the Gaia Data Processing and Analysis Consortium (DPAC). The funding for DPAC comes from national agencies, particularly institutions involved in the Gaia Multilateral Agreement. This work also utilizes filter curve and theoretical spectral model data from various surveys collected and organized by the Spanish Virtual Observatory (SVO). The funding for SVO is provided by the Spanish Ministry of Science and Innovation (Ministerio de Ciencia e Innovación). The TMAW tool (http://astro.uni-tuebingen.de/~TMAW/) used for this paper was constructed as part of the activities of the German Astrophysical Virtual Observatory (Rauch et al. 2018).

Appendix A Additional comparison

To test whether our ANN spectrum generator performs smoothly across different parameter space, we randomly select a set of parameters as inputs for the spectrum generator. The resulting flux-calibrated spectrum is compared with theoretical spectrum in the spectral library with the same parameters. The residual standard deviation is defined as $$\sigma =\sqrt{\frac{1}{N_{pix}}{\displaystyle \sum _{j=1}^{N_{pix}}{\left(R_j-\overline{R}\right)}^2}},$$(A.1)

where R_j represents the flux residual of the j-th pixel, and $\overline{R}$ represents the average flux residual. The residual standard deviation only demonstrate the deviation along wavelength axis, but our spectrum generator actually works on a high-dimensional parameter space. Therefore, to test the behavior along the label axis, we randomly select one single wavelength pixel, and vary the value of a given parameter while keeping the others fixed. The new generated flux is then compared with the theoretical flux. The interpolated curves of the neural network’s prediction for MS stars (the left panel of Fig. A.1) are smooth across all stellar labels, except [α/Fe]¹². This implies that our spectrum generator does not produce overfitting artifacts. Similar examinations are also performed on WDs (Figs. A.2 and A.3), but the parameters are only effective temperature and surface gravity, and the number of output wavelength pixels is kept consistent with the number of the corresponding theoretical spectral wavelength pixels.

Next, we examine the relative flux error of each pixel for each type of stars. Fig.A.4 shows the histogram for the relative flux errors for 101 MS spectra (1036 DA, and 138 DB). The mean values are closed to 0, and the relative errors are mostly less than 10%.

Finally, it is also important to know which part of the parameter space (T vs. log 𝑔) is more reliable. Figures A.5, A.6, and A.7 show maps of relative errors for DA, DB, and MS stars, respectively. We note that (1) The relative errors reach ~ 0.3% for cool pure-H white dwarf model (T < 8, 000K); (2) The relative errors of pure-He white dwarf models are quite uniform across most parameter space (no larger than 0.1%); (3) The relative errors MS stars are small (<5%) for hotter stars (T_MS > 4, 250K), and increase to ~ 10% for cooler stars.

Fig. A.1 Left panel: Flux comparison between predicted (red curves) and theoretical values (black dots), at a random wavelength pixel for a typical MS star. Right panel: Comparison of predicted (blue) and theoretical spectra (orange). The residual spectrum is shown in the bottom right, where the ±3σ lines are indicated by the dotted lines. See the text for more details.

Fig. A.2 Flux comparison between DA prediction and Koester (2010) templates. Symbols are the same as Fig. A.1.

Fig. A.3 Flux comparison between DB prediction and Koester (2010) template. Symbols are the same as Fig. A.1.

Fig. A.4 Histograms of the relative errors of the pixels in DA (blue), DB (orange) and MS (green) synthetic spectra. Only wavelength pixels within 3600 Å< λ < 9000 Å are considered.

Fig. A.5 Mean relative errors of the entire DA datasets. For each set of parameters on the grid, we calculate the mean relative error of the corresponding spectra. These relative errors are indicated by a color bar on the right.

Fig. A.6 Same as Fig. A.5, but for DBs.

Fig. A.7 Same as Fig. A.5, but for MS stars. Here we fix other two parameters: [Fe/H] = 0 and [α/Fe] = 0. Moreover, the neural network doesn’t work when MS stars have low temperature and high surface gravity.

Appendix B Fitting result of 14 white dwarf binary candidates

Fig. B.1 Fitting results of 14 CWDB candidates. Symbols are the same as in Fig. 2.

References

All Tables

All Figures

	Fig. 1 Conceptual illustration of the ANN architecture. X1 is the number of input parameters, while X2 is the number of output pixels.
In the text

Fig. 2 Fitting results for WD 1534+503. The top panel shows the SDSS spectrum (black line), the convolved smoothed continuum (dashed orange line), the pseudo-continuum obtained by polynomial fitting (dashed red line), the best-fit spectrum of Component 1 (blue line), the best-fit spectrum of Component 2 (solid orange line), and the composite spectrum (green line). The bottom panel shows the residual spectrum. The ±3σ levels are indicated by the dashed lines.

In the text

Fig. 3 Kiel diagram of three DAs in the test sample. The evolutionary tracks provided by Cheng et al. (2020) are shown as small dots, where their masses are indicated by the color bar. The red symbols with error bars represent the expected values and uncertainties of our fitting results. Similarly, orange symbols and green symbols show the literature results of Kilic et al. (2021) and Orosz et al. (1999), respectively. Obviously, the error in log 𝑔 is propagated to the stellar mass, and thus the radius. As the surface gravity is expressed on a log scale, and the mass and radius are expressed on a linear scale, the mass-radius errors are magnified.

In the text

	Fig. 4 Fitting results for PG 1224+309. Symbols are the same as Fig. 2.
In the text

Fig. 5 Our final samples on the Gaia CMD (colored symbols). Different types of binary systems identified by Ren et al. (2023) are denoted by different colored symbols. The black dots denote Gaia samples within 100 pc with small errors. The region outlined by cadet dashed blue lines was set by Ren et al. (2023) to select their candidates. The color-shaded regions indicate the locations of gold samples (CVs-green, DWDs-blue, WD+AFGK-orange, WD+M-red) compiled by Inight et al. (2021).

In the text

	Fig. A.1 Left panel: Flux comparison between predicted (red curves) and theoretical values (black dots), at a random wavelength pixel for a typical MS star. Right panel: Comparison of predicted (blue) and theoretical spectra (orange). The residual spectrum is shown in the bottom right, where the ±3σ lines are indicated by the dotted lines. See the text for more details.
In the text

	Fig. A.2 Flux comparison between DA prediction and Koester (2010) templates. Symbols are the same as Fig. A.1.
In the text

	Fig. A.3 Flux comparison between DB prediction and Koester (2010) template. Symbols are the same as Fig. A.1.
In the text

	Fig. A.4 Histograms of the relative errors of the pixels in DA (blue), DB (orange) and MS (green) synthetic spectra. Only wavelength pixels within 3600 Å< λ < 9000 Å are considered.
In the text

	Fig. A.5 Mean relative errors of the entire DA datasets. For each set of parameters on the grid, we calculate the mean relative error of the corresponding spectra. These relative errors are indicated by a color bar on the right.
In the text

	Fig. A.6 Same as Fig. A.5, but for DBs.
In the text

	Fig. A.7 Same as Fig. A.5, but for MS stars. Here we fix other two parameters: [Fe/H] = 0 and [α/Fe] = 0. Moreover, the neural network doesn’t work when MS stars have low temperature and high surface gravity.
In the text

	Fig. B.1 Fitting results of 14 CWDB candidates. Symbols are the same as in Fig. 2.
In the text