DA white dwarf variability#
In this tutorial, we will:
collect known DA white dwarf positions from recent literature, via SIMBAD bibcode queries and a VizieR table
load Rubin DP2 forced-source photometry with LSDB, including new columns with corrected flux-error estimates
cross-match the two catalogs and plot the resulting light curves, comparing original and corrected flux errors
Introduction#
Prerequisites#
In order to access Rubin data, you must be a Rubin data rights holder, and run this notebook from a platform with access to the DP2 HATS catalogs, such as the Rubin Science Platform (RSP).
[ ]:
import os
import astropy.table
import lsdb
import matplotlib.pyplot as plt
import numpy as np
from astroquery.simbad import Simbad
from astroquery.vizier import Vizier
from dask.distributed import Client
plt.rcParams.update(
{
"axes.titlesize": 16,
"axes.labelsize": 15,
"xtick.labelsize": 13,
"ytick.labelsize": 13,
"legend.fontsize": 13,
"figure.titlesize": 18,
}
)
[ ]:
import logging
logging.getLogger().setLevel(logging.WARNING)
logging.getLogger("distributed").setLevel(logging.WARNING)
[ ]:
lsdb.show_versions()
1. Collect DA white dwarf positions from multiple papers#
We query SIMBAD for the objects cited in three recent DA white dwarf papers, identified by their bibcodes, and add one more table of positions from VizieR.
[ ]:
simbad_bibcodes = [
"2025A&A...701A.269G",
"2025MNRAS.540..385B",
"2025AJ....169...40B",
]
vizier_catalogs = [
"J/A+A/699/A212/tableb2", # https://ui.adsabs.harvard.edu/abs/2025A%26A...699A.212J
]
simbad = Simbad()
vizier = Vizier(columns=["_RAJ2000", "_DEJ2000"])
Query each paper’s objects from SIMBAD and the VizieR table, keeping only right ascension and declination, and stack the results into a single astropy table.
[ ]:
tables = []
for bibcode in simbad_bibcodes:
table = simbad.query_bibobj(bibcode)[["ra", "dec"]]
tables.append(table)
for viz_cat in vizier_catalogs:
table = vizier.get_catalogs(viz_cat)[0][["_RAJ2000", "_DEJ2000"]]
table.rename_columns(["_RAJ2000", "_DEJ2000"], ["ra", "dec"])
tables.append(table)
table = astropy.table.vstack(tables, metadata_conflicts="silent")
table
1.1 Convert to an LSDB catalog#
We wrap the combined table with lsdb.from_astropy so it can be used as one side of the cross-match below.
[ ]:
da_cat = lsdb.from_astropy(table)
da_cat
2. Load the Rubin DP2 forced-source catalog#
We open the DP2 object collection catalog, requesting all default columns plus two new nested forced-source subcolumns: psfFluxErr_corrected and psfDiffFluxErr_corrected.
Four new subcolumns are added to the nested forced-source tables (objectForcedSource and diaObjectForcedSource): psfFluxErr_corrected, psfDiffFluxErr_corrected, psfFluxErr_corrected_flag, and psfDiffFluxErr_corrected_flag. They come from a model trained to correct flux errors so that, for non-variable objects, the light curve’s reduced χ² is close to unity — these corrected values may be useful for light-curve fitting. The flag columns are True when the ratio between
the original and corrected errors falls outside the [0.1, 50] interval. See Malanchev et al., in prep., for details; source code is available on GitHub.
We also filter to objects with PSF magnitude brighter than 22 in every band, to keep the cross-match in this tutorial small and fast.
[ ]:
dp2_obj = lsdb.open_catalog(
"/rubin/lsdb_data/object_collection",
columns=[
..., # all default columns
"objectForcedSource.psfFluxErr_corrected",
"objectForcedSource.psfDiffFluxErr_corrected",
],
# Filter data on the read time: keep objects which are brighter 22 mag
# in any band, since we know that these DA dwarfs are pretty bright.
filters=[[(f"{band}_psfMag", "<", 22) for band in "ugrizy"]],
)
dp2_obj
3. Cross-match with Rubin DP2#
We cross-match the literature white dwarf positions against the DP2 objects within 5 arcsec. As with all LSDB operations, this only plans the cross-match — no data has been loaded yet (see Lazy evaluation).
[ ]:
xmatch = dp2_obj.crossmatch(da_cat, radius_arcsec=5, suffix_method="overlapping_columns")
xmatch
4. Run the pipeline#
Now we create a local Dask cluster and run the cross-match. See the Dask Client tutorial for tips on configuring Dask for LSDB.
[ ]:
with Client(
n_workers=2,
threads_per_worker=1,
memory_limit="8GB",
local_directory=f"/deleted-sundays/{os.environ.get('USER', 'dask_scratch')}",
) as client:
print(f"Dask dashboard: {client.dashboard_link}")
display(client)
df = xmatch.compute()
df
5. Plot the light curves#
For each matched white dwarf, we plot the per-band objectForcedSource light curves in both science flux and difference flux, comparing error bars from the original *Err columns against the corrected *Err_corrected columns, and report the reduced χ² for each.
[ ]:
band_colors = {
"u": "#0c71ff",
"g": "#49be61",
"r": "#c61c00",
"i": "#ffc200",
"z": "#f341a2",
"y": "#5d0000",
}
bands = "ugrizy"
def reduced_chi2(y, yerr):
n = len(y)
if n <= 1:
return None
wmean = np.average(y, weights=1 / yerr**2)
chi2 = np.sum(((y - wmean) / yerr) ** 2)
return chi2 / (n - 1)
def chi2_text(value):
return "N/A" if value is None else f"{value:.2f}"
def plot_for_flux(
forced_sources,
ra,
dec,
plot_title,
flux_col,
err_col,
err_corr_col,
flag_col,
invalid_flag_col=None,
):
fig, axes = plt.subplots(2, 3, figsize=(15, 8), sharex=True)
axes = axes.ravel()
fig.suptitle(rf"$\alpha={ra:.5f}$, $\delta={dec:.5f}$" f"\n{plot_title}")
for ax, band in zip(axes, bands):
color = band_colors[band]
mask = (
(forced_sources["band"] == band)
& (~forced_sources[flag_col])
& (~forced_sources["pixelFlags_bad"])
& (~forced_sources["pixelFlags_saturated"])
& (~forced_sources["pixelFlags_nodata"])
)
if invalid_flag_col is not None and invalid_flag_col in forced_sources.columns:
mask &= ~forced_sources[invalid_flag_col]
lc = forced_sources[mask].copy()
chi2_orig = reduced_chi2(lc[flux_col], lc[err_col])
chi2_corr = reduced_chi2(lc[flux_col], lc[err_corr_col])
ax.errorbar(
lc["midpointMjdTai"],
lc[flux_col],
yerr=lc[err_corr_col],
color=color,
linestyle="none",
marker="o",
markersize=6,
markerfacecolor="none",
markeredgewidth=1.0,
alpha=0.35,
elinewidth=2.5,
capsize=3,
zorder=1,
label="corrected err",
)
ax.errorbar(
lc["midpointMjdTai"],
lc[flux_col],
yerr=lc[err_col],
color=color,
linestyle="none",
marker="x",
markersize=7,
alpha=0.95,
elinewidth=2.0,
capsize=2,
zorder=2,
label="orig err",
)
ax.set_title(
f"{band}-band {plot_title} "
f"$\\chi^2_\\nu$={chi2_text(chi2_orig)} "
f"$\\chi^2_\\nu$$_{{corr}}$={chi2_text(chi2_corr)}",
fontsize=12,
)
ax.grid(alpha=0.3)
ax.legend(loc="lower right", fontsize=11)
axes[0].set_ylabel(flux_col)
axes[3].set_ylabel(flux_col)
for ax in axes[3:]:
ax.set_xlabel("MJD")
fig.tight_layout(rect=[0, 0, 1, 0.95])
plt.show()
for row_index in range(len(df)):
row = df.iloc[row_index]
forced_sources = row["objectForcedSource"]
ra = row["coord_ra"]
dec = row["coord_dec"]
plot_configs = [
(
"science Flux",
"psfFlux",
"psfFluxErr",
"psfFluxErr_corrected",
"psfFlux_flag",
"invalidPsfFlag",
),
(
"differential Flux",
"psfDiffFlux",
"psfDiffFluxErr",
"psfDiffFluxErr_corrected",
"psfDiffFlux_flag",
None,
),
]
for plot_title, flux_col, err_col, err_corr_col, flag_col, invalid_flag_col in plot_configs:
plot_for_flux(
forced_sources=forced_sources,
ra=ra,
dec=dec,
plot_title=plot_title,
flux_col=flux_col,
err_col=err_col,
err_corr_col=err_corr_col,
flag_col=flag_col,
invalid_flag_col=invalid_flag_col,
)
Bright objects see their errors slightly overestimated, by about 30%, in exchange for reliably catching large error underestimation present in the original data.
About#
Authors: Konstantin Malanchev
Last run: July 27, 2026
If you use lsdb for published research, please cite following instructions.