# Every Datapoint Counts: Stellar Flares as a Case Study of Atmosphere Aided
Studies of Transients in the LSST Era

Riley W. Clarke Department of Physics & Astronomy, University of Delaware,
Newark, DE 19716, USA  Data Science Institute, University of Delaware, Newark,
DE 19716, USA  James R. A. Davenport DiRAC Institute, Department of Astronomy,
University of Washington, Seattle, WA 98195, USA  John Gizis Department of
Physics & Astronomy, University of Delaware, Newark, DE 19716, USA  Data
Science Institute, University of Delaware, Newark, DE 19716, USA  Melissa L.
Graham DiRAC Institute, Department of Astronomy, University of Washington,
Seattle, WA 98195, USA  Xiaolong Li Department of Physics & Astronomy,
University of Delaware, Newark, DE 19716, USA  Data Science Institute,
University of Delaware, Newark, DE 19716, USA  Willow Fortino Department of
Physics & Astronomy, University of Delaware, Newark, DE 19716, USA  Data
Science Institute, University of Delaware, Newark, DE 19716, USA  Ian Sullivan
DiRAC Institute, Department of Astronomy, University of Washington, Seattle,
WA 98195, USA  Yusra Alsayyad  Department of Astrophysical Sciences, Princeton
University, Princeton, NJ, 08544, USA  James Bosch Department of Astrophysical
Sciences, Princeton University, Princeton, NJ, 08544, USA  Robert A. Knop
E.O. Lawrence Berkeley National Laboratory, 1 Cyclotron Rd., Berkeley, CA,
94720, USA  Federica Bianco Department of Physics & Astronomy, University of
Delaware, Newark, DE 19716, USA  Data Science Institute, University of
Delaware, Newark, DE 19716, USA  Joseph R. Biden School of Public Policy and
Administration, University of Delaware, Newark, DE 19716, USA  Vera C. Rubin
Observatory

###### Abstract

Due to their short timescale, stellar flares are a challenging target for the
most modern synoptic sky surveys. The upcoming Vera C. Rubin Legacy Survey of
Space and Time (LSST), a project designed to collect more data than any
precursor survey, is unlikely to detect flares with more than one data point
in its main survey. We developed a methodology to enable LSST studies of
stellar flares, with a focus on flare temperature and temperature evolution,
which remain poorly constrained compared to flare morphology. By leveraging
the sensitivity expected from the Rubin system, Differential Chromatic
Refraction can be used to constrain flare temperature from a single-epoch
detection, which will enable statistical studies of flare temperatures and
constrain models of the physical processes behind flare emission using the
unprecedentedly high volume of data produced by Rubin over the 10-year LSST.
We model the refraction effect as a function of the atmospheric column
density, photometric filter, and temperature of the flare, and show that flare
temperatures at or above $\sim$4,000$K$ can be constrained by a single
$g$-band observation at airmass $X\gtrsim 1.2$, given the minimum specified
requirement on single-visit relative astrometric accuracy of LSST, and that a
surprisingly large number of LSST observations is in fact likely be conducted
at $X\gtrsim 1.2$, in spite of image quality requirements pushing the survey
to preferentially low $X$. Having failed to measure flare DCR in LSST
precursor surveys, we make recommendations on survey design and data products
that enable these studies in LSST and other future surveys.

##  1 Introduction

Differential Chromatic Refraction (DCR) is the apparent shift of a celestial
object on the sky due to atmospheric refraction. The magnitude of the shift is
determined by the airmass, observing bandpass, and Spectral Energy
Distribution (SED, which is dominated by the temperature and spectral
composition) of the source. The direction of the shift is toward the zenith.
DCR must be characterized effectively to ensure it does not compromise image
quality, astrometric inference (Van De Kamp, 1967), spectrophotometry
(Filippenko, 1982), and in surveys where time-domain phenomena are discovered
and studied in template subtracted images (Lupton, 2007). DCR, however, can
also aid scientific inference in a few cases. For example, Kaczmarczik et al.
(2009) used DCR to improve redshift measurements of quasars in the Sloan
Digital Sky Survey (SDSS). While innovative techniques have been previously
employed to study the chromatic evolution of stellar flares (_e.g._ Hedges et
al. (2021)), the use of ground-based observations to infer information about
flare SEDs via atmospheric refraction remains untested. We consider flares to
be suitable candidates for inferential applications of DCR due to their
stochastic and short-lived nature, which makes detections in multiple filters
difficult, their nature as chromatic transients, and their heightened
occurrence on cool stars like M dwarfs, the most common type of star in our
galaxy. Currently, there are only a small number of measurements of flare
temperatures, which our usage of DCR will enable.

Flares are stochastic, short-lived stellar transients caused by magnetic
reconnection at the stellar photosphere (Pettersen, 1989), most commonly on
low-mass stars such as M-dwarfs ($\mathrm{d}M$ hereafter). Flares are tracers
of stellar magnetic activity, and flare rates have been studied as a correlate
with stellar age and rotation (Davenport, 2016). Flares can also have a
significant impact on planetary atmospheres, as high amounts of UV flux from
flares can deplete their ozone layers (Tilley et al., 2019) and are thought to
have a role in the formation of life, where they have been identified as both
a potential trigger and inhibitor (Ramsay et al., 2021).

Photometrically, they appear as a highly chromatic sudden rise (typically
$\sim 3$ magnitudes in $u$ band and between 0.5 and 1 magnitude in $g$ band,
see Davenport et al. 2012) followed by an exponential decay phase on a
timescale ranging from a few minutes to 100-200 minutes in some cases (Yan et
al., 2021). While the photometric evolution of flares is well-understood
thanks to high-accuracy, rapid-cadence measurements by exoplanet-finding
satellites like Kepler (Borucki et al., 2003) and TESS (Ricker et al., 2010),
accurate color measurements remain scarce. Statistical inference on flares’
origin and evolution mechanisms and their impact is limited by this lack of a
large ensemble of information on flare temperatures and temperature evolution,
prompting a turn to large next-generation astrophysical surveys, most notably
the Vera C. Rubin Observatory Legacy Survey of Space (LSST).

LSST (Ivezić et al., 2019) is the premier ground-based photometric survey of
the 2020s. It will continuously scan the whole of the southern sky over a
10-year period, producing over 3 million images with a sub-arcsecond
resolution of the whole southern sky in six photometric filters: $ugrizy$.
Based on flare rate measurements in SDSS Stripe 82 (Kowalski et al., 2009),
the number of flares will range between 0.4 and 1.4 per image, depending on
galactic latitude. While it is certain that Rubin will observe many flares,
whether or not we can access scientifically valuable information from those
detections is currently an open question. It is expected that the LSST primary
survey, referred to as “Wide, Fast, Deep” (WFD), will take 2-3 exposures per
pointing per night with a median internight gap of 34 minutes, repeating
observations of a field in the same filter every few days ($\sim 6$ days in
$r$, $\sim 21$ days in $g$, median values, see Bianco & The Rubin Observatory
Survey Cadence Optimization Committee 2023 for details on the most recent LSST
observing plans). Given that Yan et al. 2021 reported median rise and decay
times of flares on Sun-like stars in the Kepler catalog to be 5.9 minutes to
22.6 minutes respectively, LSST will observe most flares in the WFD with 1-2
data points per detected event, rendering the time-resolved photometric
analysis enjoyed by _e.g._ Kepler impossible. The Deep Drilling Fields (DDF)
will offer an opportunity for time-resolved flares and to identify even faint
$\mathrm{d}M$ flare precursors. Roughly 5% of the total survey time will be
devoted to continuous observations of five selected pointings. While the
observational details for the DDFs are yet to be defined, in one hour over 100
consecutive exposures can be taken, over 20 in each of five filters of one 9.6
sq degree field of view.111Note that the Rubin filter wheel can host 5 filters
on any night. The reddest and bluest filters are swapped based on moon cycles
and five filters are available when observing. When co-added with the main
survey observations these images will extend the 5$\sigma$ depth down to
$r\sim 27.5$, allowing time-resolved observations of even the faintest
$\mathrm{d}M$s in multiple filters. It must be noted that all DDFs are
extragalactic fields, where the flare rate is expected to be suppressed.222For
details on the selection of the fields see Bianco & The Rubin Observatory
Survey Cadence Optimization Committee 2023. Yet flares are common phenomena
even in the extragalactic sky, and we expect hundreds to be detected in each
DDF by the end of the 10-year survey. The LSST DDF program will enable
traditional flare studies. However, by exploiting the high astrometric
accuracy of Rubin, here we show that the Rubin WFD, with its unprecedented
sample size, offers a valuable opportunity for flare science.

It is already expected that DCR will aid extragalactic studies of quasars (Yu
et al., 2020; Richards et al., 2018) with Rubin, and we here propose that
employing DCR as a scientific tool can be extended to variable and transient
phenomena with Rubin. We use stellar flares as our case study and demonstrate
in this paper that by taking advantage of Rubin’s pristine image quality and
astrometric precision to measure the difference in source location on the sky
between quiescence and event due to DCR, we will be able to indirectly obtain
information regarding the color, and therefore the temperature of the flare,
from even a single point detection. We will show that an astrometric shift
should be apparent in LSST astrometry for a star when flaring compared to when
quiescent, with the effect being especially pronounced for a hot flare on a
cool star observed at high airmass; but we further demonstrate that the effect
will be measurable by Rubin even at its typically low airmass range. This
means that, despite the sparse photometric sampling of Rubin, the DCR toolset
can be used to leverage its high volume of flare detections in order to infer
flare temperatures across an expansive sample.

In section 2 we describe how DCR is calculated and how we model the expected
excess DCR during a flare with some temperature, at some airmass, in a given
filter. In section 3 we compare the expected DCR for flaring stars with the
capabilities of the Rubin system and data analysis pipeline. In section 4 we
describe our failed search for excess DCR produced by flares in surveys
considered precursors to LSST, ZTF (Bellm et al., 2018) and Deep Drilling in
the Time Domain with DECam (Graham et al., 2023) and identify the bottlenecks
that impaired a successful DCR measurement. In section 5 we discuss the
technical requirements for measuring the $\Delta\mathrm{DCR}$ from flares: the
excess zenith-bound displacement of a flaring star compared to the apparent
shift expected based on its quiescent SED, and provide recommendations to
maximize the potential of this method for LSST and other surveys. In section 6
we summarize our conclusions and explain how this method can be embedded
within Rubin’s software ecosystem to support studies of flares and other
chromatic phenomena throughout the survey lifespan.

##  2 Methodology

Figure 1: Angular deflection from true position for a source with a 10,000$K$
blackbody spectrum for airmass $1.1\leq X\leq 2.1$ in all six LSST bands. A
zoomed section is shown to reveal the separation between curves corresponding
to separate bandpasses.

Following (Kaczmarczik et al., 2009), we can describe the formalism of DCR in
the following four steps:

$\lambda_{\mathrm{eff}}=\frac{\int_{0}^{\infty}f_{\lambda}S_{j}(\lambda)\ln(\lambda)d\lambda}{\int_{0}^{\infty}f_{\lambda}S_{j}(\lambda)d\lambda}$
(1)

defines $\lambda_{\mathrm{eff}}$ as the effective wavelength of bandpass $j$,
where $S_{j}$ is the transmission function of $j$, and $f_{\lambda}$ is the
spectral flux of the source (or SED) in units of $Wm^{-2}\mu^{-1}$; for an
observed source then, the filter-dependent refraction index $n_{\lambda}$ can
be calculated from

$\begin{split}[n_{\lambda}-1]~{}10^{6}=64.328\>+&\frac{29498.1}{146-(1/\lambda_{\mathrm{eff}})^{2}}\\\
+&\frac{255.4}{41-(1/\lambda_{\mathrm{eff}})^{2}};\end{split}$ (2)

$R$, the total angular deflection of the source due to DCR (in arcseconds), is
then calculated as:

$R_{0}(\lambda)=\frac{n_{\lambda}^{2}-1}{2n_{\lambda}^{2}},$ (3)

$R=R_{0}(\lambda)\tan(Z)$ (4)

where Z is the zenith angle, and where the airmass $X$ is commonly
approximated as $X=\sec(Z)$, which assumes a homogenous, plane-parallel
atmosphere. This approximation is valid only for $60^{\circ}\leq Z\leq
75^{\circ}$, so our calculations involving airmass will be restricted to
$X\leq 2.1$ (although we note that the LSST survey strategy includes
observations at higher airmasses, see section 3 and subsection 5.2). Figure 1
illustrates the amplitude of DCR produced by a 10,000$K$ blackbody source for
a range of airmasses in the six LSST bandpasses.

Figure 2: Light incident from a star is deflected by the atmosphere. The
amount of deflection depends on the color (_i.e._ temperature) of the source
and the amount of atmosphere the light passes through. The chromatic change
during a flare event should produce an excess in the normal DCR at quiescence,
labeled as $\Delta\mathrm{DCR}$ in the figure. Figure 3: In pink, the spectrum
of an M5 dwarf (composite spectrum built partially from SDSS observations,
Davenport et al. 2012). In blue, the spectrum of a 10,000$K$ blackbody to
simulate the flare. In purple, the sum of the two aforementioned spectra,
representing the spectrum of the $\mathrm{d}M$ during the flare event. The
blackbody and $\mathrm{d}M$ spectra contain the same total energy over the
SDSS optical range ($3,850-9,200~{}\text{\AA}$). The transmission functions
for the LSST $ugrizy$ photometric system are shown in grey.

In section 4 we will demonstrate that even flares with moderately high
temperatures observed at airmasses $1.05\leq X\leq 1.2$ will produce an
astrometric shift that will be detectable at the LSST precision level.

###  2.1 The DCR effect induced by flares: $\Delta\mathrm{DCR}$

Using the DCR formalism described in section 2, we created an expository model
to illustrate the DCR produced by a flare. We assume the star is an M5 dwarf
and model the quiescent SED of the source with a template spectrum spanning
the optical and near-infrared wavelength ranges built by Davenport et al.
(2012), with the optical component using observations from the Sloan Digital
Sky Survey (SDSS). We clip this template spectrum at $12,000\text{\AA}$ as
this was sufficient for full coverage of the LSST $ugrizy$ bands. The source
in its quiescent state will already suffer from a DCR effect, but we expect
that the temperature contrast between flare and quiescent states will
determine a significant change in the apparent star position. The difference
between the displacement in quiescent and event state is the excess
displacement produced by the flare which we refer to as $\Delta\mathrm{DCR}$.
A schematic illustration of $\Delta\mathrm{DCR}$ is shown in Figure 2.

To simulate a flare, we add a blackbody spectrum at temperature T to the
spectrum of the star (Figure 3). In the absence of a well-characterized
spectral energy distribution (SED), flare spectra are canonically approximated
as a 9,000-10,000$K$ blackbody (Osten & Wolk, 2015). Thus, in this simulation
the spectra are normalized such that the blackbody component at the canonical
flare temperature ($T=10,000$K) has the same total energy as the $\mathrm{d}M$
spectrum when both are integrated over the SDSS optical range
($3,850-9,200~{}\text{\AA}$). We also explored the effect of changing the
fraction of the stellar surface covered by the flaring region, referred to as
a “filling factor”. Assuming our energy-based calibration to be representative
of a filling factor $f_{f}\sim 0.05$ and applying a flat scaling to the
blackbody spectrum for other $f_{f}$ values, we tested filling factors ranging
from 0.05 to 0.2 and found that we are not sensitive to this parameters for
$f_{f}\geq 0.05$. The composite spectrum is convolved with the chosen LSST
filter to calculate $\lambda_{\mathrm{eff}}$ in Equation 1, allowing us to
estimate the angular displacement from DCR during both quiescence and event,
and their difference, or $\Delta\mathrm{DCR}$.

###  2.2 Measuring $\Delta\mathrm{DCR}$ with the Parallactic Angle

Many astrophysical and observational effects may cause the apparent position
of a star to change between two images of the same area of the sky. Although
stars’ proper motion is typically too slow to have a measurable effect in
images collected within days, weeks, or months, images are however subject to
various distortion effects such that a star may appear to be offset between
two images, especially if observed in a different position of the CCD plane,
with different telescope rotation, _etc_. These observational systematics
should be corrected, and we will return to these spurious effects in section
4. These components of the star’s displacement should have no preferential
direction while the DCR-induced displacement should be strictly in the
direction of the zenith. To determine the DCR-induced change in position of a
given source between detections, we will measure the components of the
apparent motion toward, tangential to, or away from the zenith direction. The
parallactic angle is defined as the angle between two great circles, one
passing through the source and the zenith, and the other passing through the
source and the North Celestial Pole. It is calculated according to Meeus
(1998) as:

$P=\tan^{-1}\left(\frac{\sin(h)}{\cos(\delta)\tan(\phi)-\sin(\delta)\cos(h)}\right),$
(5)

where $\delta$ is the object’s declination, $h$ is the object’s hour angle,
and $\phi$ is the geographic latitude of the observer’s location on Earth. The
component of the source’s motion in the direction of the zenith (which we call
$d_{\parallel}$) is then calculated as:

$d_{\parallel}=\sqrt{\Delta\alpha^{2}+\Delta\delta^{2}}\cos\left(\frac{\pi}{2}-P_{2}-\tan^{-1}\left(\frac{\Delta\delta}{\Delta\alpha}\right)\right),$
(6)

where $\alpha$ is the object’s right ascension, and $P_{2}$ is the parallactic
angle of the source in its second position.

##  3 Detection Potential for LSST

Using the methodology described in section 2, we simulate the
$\Delta\mathrm{DCR}$ produced by a flare on an M5 star for a variety of
airmasses and flare temperatures across all six Rubin bands. The
$\Delta\mathrm{DCR}$ for a fixed flare temperature at 10,000$K$ for all 6
filters in the 1.1-2.1 airmass range is shown in Figure 4. The Rubin system
specification requires a single-image absolute astrometric accuracy of Rubin
of 0.1 arcsec (Ivezić & The LSST Science Collaboration, 2018). This figure
demonstrates that a 10,000$K$ flare’s $\Delta\mathrm{DCR}$ is detectable in
$g$-bands even at moderate airmasses; $\sim 8\%$ of all WFD images will be in
$g$ band.

While the effect is also detectable in $u$-band, higher airmasses are required
to generate a comparably large shift. This may be counterintuitive, as the DCR
effect is more pronounced in bluer wavelengths. However, the
$\Delta\mathrm{DCR}$ is proportional to the change in $\lambda_{\mathrm{eff}}$
during the flare compared to quiescent state. The quiescent SED has a
prominent redward slope in the $g$ band wavelengths which, integrated through
the filter leads to a $\lambda_{\mathrm{eff,quiescent}}(g)=4,988.4\text{\AA}$
which becomes $\lambda_{\mathrm{eff,flare}}(g)=4,740.9\text{\AA}$ when the
black body contribution dominates the SED during flare. Conversely, the
quiescent SED in $u$ band is more flat and only slopes redward significantly
where the filter transmission is already very low. Thus the shift in
$\lambda_{\mathrm{eff}}$ is smaller
($\lambda_{\mathrm{eff,quiescent}}(u)=3,728.0\text{\AA}$ compared to
$\lambda_{\mathrm{eff,flare}}(u)=3,655.8\text{\AA}$).

The $g$-band evolution of the $\Delta\mathrm{DCR}$ with flare temperature is
shown in Figure 5, showing that the $\Delta\mathrm{DCR}$ grows rapidly between
5,000-10,000$K$, then flattens off at higher temperatures. The minimum system
requirement and stretch goals for both absolute and relatve astrometric
accuracy of Rubin is indicated on Figure 5, demonstrating that DCR could be
used to probe flare temperatures as a low as 4,000$K$, depending on the type
of data products used for the analysis (see section 5).

Figure 4: Expected magnitude of the $\Delta\mathrm{DCR}$ effect for a flare
SED approximated by a 10,000$K$ blackbody as a function of airmass and filter
in the Rubin Observatory $ugrizy$ observing system (Olivier et al., 2008).
Blue coloring corresponds to a $\Delta\mathrm{DCR}$ shift detectable by Rubin,
and red coloring corresponds to an undetectable shift, given the absolute
astrometric accuracy goal of 0.1 arcsec. Figure 5: The magnitude of a star
image displacement during a flare compared to the quiescent star position as a
function of flare temperature in LSST $g$-band for four different airmass
values. The minimum requirement and stretch goal for absolute astrometric
accuracy of Rubin are shown by the the black and grey dashed lines,
respectively. The minimum requirement and stretch goal for relative
astrometric accuracy of Rubin are shown by the black and grey dash-dot lines,
respectively. Figure 6: Airmass distribution of the current LSST survey
strategy proposal (baseline_v3.0_10yrs, Bianco & The Rubin Observatory Survey
Cadence Optimization Committee 2023). Cumulative distribution showing the
number of visits at or below a given airmass. Plot produced within the Rubin
Metric Analysis Framework (Jones et al., 2014).

To maximize image quality, the automated LSST scheduling system will
preferentially observe at low airmass by design (Ivezić & The LSST Science
Collaboration, 2018). The dependence of this effect on airmass then begs the
question of whether or not LSST will perform observations at sufficiently high
airmass for the effect to be detectable. However, the most recent fiducial
simulation of the LSST (baseline_v3.0_10yrs, Bianco & The Rubin Observatory
Survey Cadence Optimization Committee 2023) suggests that nearly $10^{6}$ WFD
visits will be performed above airmass 1.25 (Figure 6), which, as we have
shown in section 3, is sufficient to produce a detectable $\Delta\mathrm{DCR}$
at typical flare temperatures (we will return to the LSST design in subsection
5.2).

These simple calculations indicate that Rubin should be capable of detecting
even relatively cool flares at typical airmass on a typical M5 dwarf. However,
contamination from image warping, chromatic aberation, DIA misalignments,
_etc._ will complicate this simplified expectation.

##  4 Precursor Surveys

We tested our method on detected flares in precursor surveys selecting as
flare candidates transient events with short duration ($\lesssim 2$-hour) and
a change of magnitude of at least 0.4 mag in $g$ band.

###  4.1 Zwicky Transient Facility

The Zwicky Transient Facility (ZTF, Bellm et al. 2018) was designed to detect
transient objects across the entire northern hemisphere. It is considered a
precursor survey for LSST and it will deliver images and lightcurves for 3,750
sqdeg/hour with alerts delivered in real-time and a significant fraction of
its data made available without proprietary restrictions. The exposure time is
similar to that of Rubin’s LSST (30 seconds + 10 seconds readout) to a single
image depth of ZTFr 20.5, and three filters are available (ZTFg, ZTFr, ZTFi)
also with similar throughput to the Rubin filters. With a smaller footprint by
a factor of four and three filters the revisit time is nominally nearly one
order of magnitude shorter than LSST’s, enabling multiple observations of the
same flare to be collected.

Notable differences, however, for the purpose of our science, are the overall
data throughput, about 10 times smaller than Rubin as measured in bytes of
data, but which corresponds to a factor of nearly 100 fewer targets due to the
decreased limiting magnitude, and a significant decrease in image quality,
measured as Point Spread Function (PSF), leading to decreased astrometric
accuracy. The instrumental pixel scale is 1.01 arcsec/pixel, compared to
Rubin’s $\sim$0.2 arcsec/pixel, and the median seeing-limited PSF is $\sim$2
arcsec, compared to $\lesssim 1$ arcsec expected for Rubin. This forced us to
limit the study to flares observed at high airmass to enable DCR-based
temperature estimates, further decreasing the size of our flare sample. In our
preliminary analysis, we inspected a sample of 17,000 bright $\mathrm{d}M$ in
ZTF; 2.3% showed possible large flares in data collected airmass above 1.4,
for a total of 414 candidates. To increase the confidence that the observed
brightening was indeed a flare we required more than one observation within
the event (_i.e._ within $\sim 2$ hours). Of these 414 flares, only one was
captured with more than one data point. ZTF’s reconstructed astrometric
accuracy per science image with respect to Gaia DR1 is $\sim
0.045-0.085^{\prime\prime}$ for sources extracted at a 10-$\sigma$ limit,
however, the seeing Full Width Half Maximum (FHWM) in this pair of images was
3.474” for the first flare epoch and 2.910” for the second. The seeing turned
out to be the bottleneck in the application of our method; the astrometric
solution generated by the ZTF pipeline measured a small $\sim
0.2^{\prime\prime}$ displacement between the first and second image compared
to other stars in the field, too small compared to the seeing to confidently
determine whether or not the star moved towards the zenith.

###  4.2 DECam Deep Drilling Fields

The Dark Energy Camera (DECam; Flaugher et al. 2015) at the Blanco 4-meter
telescope is the instrument that enabled the Dark Energy Survey (Dark Energy
Survey Collaboration et al., 2016). The camera itself is a single chip of the
same kind as those that will constitute the LSST camera mosaic, making DECam
data naturally comparable with LSST’s, with high image throughput, similar
image quality (0.263 arcsecond/pixel resolution, with a telescope located near
the site of LSST leading to similar sky properties), and similar system
wavelength coverage ($grizY$ filters).

Graham et al. (2023) used DECam to survey two of the LSST DDFs: COSMOS
(Scoville et al., 2007) and ELAIS (Oliver et al., 2000). This led to a
precursor survey of the LSST DDFs with 5-sigma limiting magnitudes $r\sim
23.5$ mag (single exposure). However, in addition to the shallower depth,
there are other significant differences. This DECam DDF program was one of the
co-founding members of the DECam Alliance for Transients (DECAT), within which
multiple PIs of DECam programs pooled their time to enable dynamic queue-like
scheduling and time-domain science. With a field of view of 9.6 deg2, the LSST
DOE camera will cover each field with a single pointing. Conversely, the DECam
field of view (3 deg2) is smaller and DECam DDF program used three adjacent
pointings in COSMOS and two in ELAIS. Every night of observations cycled
through each pointing five times, obtaining a sequence of $gri$ images per
pointing. However, because the fields are covered with multiple pointings, 10
or 15 observations per night are obtained in the area where the pointings
overlap (see Figure 7). In contrast, LSST will obtain tens of consecutive
images in five filters on each night when a DDF will be observed. Taking these
differences into consideration, the DECam DDF program is an effective test-bed
for our DCR flare studies, but not for temporally resolved flare
investigations.

Figure 7: Location of DECam DDF flare candidate $g$-band objects in the COSMOS
field. Figure 8: Scatter plots of 1015 flare candidates in the DECam DDF
survey of the COSMOS field (for selection criteria discussed in subsection
4.2). The $x$ axis is the change in $g$ magnitude between the objects and the
$y$-axis is the change in the component of the source’s movement in the
direction of the zenith. Each point represents one unique candidate ID with at
least two DIA “objects” (two distinct observations). The change in magnitude
and position are calculated as the difference between the two objects with the
largest absolute change in magnitude. Left panel: For observations with
airmass $X<1.2$, points are shown as grey circles. Points are colored by the
airmass of the initial detection for airmasses $X>1.2$. Right panel: Each
point is colored by the change in airmass between the two objects and
histograms show the marginalized distributions of $g$ magnitude change (_top_)
and zenith-bound displacement (_right_). Three candidates of interest, as
described subsection 4.2, are circled. The single candidate whose quiescent
counterpart can be identified in the Gaia DR3 dataset, DC21engi, is marked by
an arrow.

Importantly, the data processing pipelines for DECam DDF and LSST also differ:
LSST will process its data with dedicated pipelines (Bosch et al., 2018),
while the DECam DDF fields data are processed with existing software to detect
transient events (Graham et al., 2023). Briefly, the difference-image analysis
(DIA) pipeline implemented in this survey ingests raw images directly from the
NOIRLab data archive and performs standard data reduction procedures. Source
Extractor (Bertin & Arnouts, 1996) is used to detect all sources in the image,
SCAMP (Bertin, 2006) is used to calculate the astrometry for each chip and
match each source with stars drawn from the Gaia DR2 catalog, and then SWARP
(Bertin, 2010) is used to solve for the world coordinate system (WCS) using
these objects. Object catalogs are generated from the reference images with
Source Extractor and aligned with SCAMP and SWARP. The image subtraction is
done with HOTPANTS (Becker, 2015). Notably, there is no DCR correction in the
difference image analysis as implemented for the DECam DDF fields. Source
Extractor is used on the resultant difference image to identify residual
signals and extract fluxes via forced photometry. Lastly, using the algorithm
described in Goldstein et al. (2015), each detection is assigned a
“real/bogus” score333 The real/bogus score is the output of a machine learning
classifier designed to distinguish “real” sources from “bogus” detections
(_i.e._ image or subtraction artifacts and moving sources)., ranging from 0
(bogus) to 1 (real). Objects (_i.e._ detections) within 2” of a previously
detected object are associated with the same candidate ID.

We elected to focus our analysis on the COSMOS $g$ band sample (recall our
preference for $g$ band from section 3). We selected a subset of the publicly
available DECam DDF transient catalogs444See Section 3.7.2 of (Graham et al.,
2023) for access to these publicly available DECam DDF catalogs. in the COSMOS
field, requiring the following:

  * •

at least 1 g-band detection,

  * •

all detections within 0.5 days (single-night),

  * •

a mean “real/bogus” score $\overline{RB}\geq 0.6$ across all detections.

These cuts resulted in a sample of 1230 candidates, and the coordinates of
this $g$ band detections’ sample are shown in Figure 7. To identify flares in
the data and measure the change in magnitude during the event we applied the
following additional cut, resulting in a final sample of 1015 COSMOS flare
candidates:

  * •

at least two $g$-band detections, each with real/bogus score $RB\geq$ 0.55

We calculated $d_{\parallel}$, the component of the candidates’ motion along
the parallactic angle (toward or away from the zenith), for all pairs of
$g$-band detections with the same on-sky association. $d_{\parallel}$ is
plotted against the measured change in $g$ magnitude in Figure 8 in order to
visualize the sample in an informative phase space and select interesting
candidates for additional inspection. We expected flares to be found in the
upper left quadrant of the plot, _i.e._ pairs of detections where the
brightness of the source increased (relative to its initial difference image
detection) and also moved on the sky in the direction of the zenith at the
time of observation, or lower right quadrant for dimming flares moving away
from zenith. In addition, we favor candidates observed at higher airmasses,
but with little to no change in airmass between observations, as a change in
airmass dominates the chromatic contribution to DCR (see Figure 1). However,
as shown in the left panel of Figure 8, few observations of the DECam DDF
fields occur at airmass above $X=1.2$, and none of these high-airmass
detections displayed sufficient change in brightness to be considered for
additional inspection.

Figure 9: Image triplets for DECam DDF candidates, as shown by the DECAT LBL
Pipeline Candidate Viewer. From top to bottom, the candidates are: DC21fjnb,
DC21fygb, DC21flly, and DC21engi. “New” denotes the search image, “Ref” the
reference image, and ”Sub” the difference image. A quiescent source consistent
with a $\mathrm{d}M$ is seen in the reference image for DC21engi. Figure 10:
WISE W1 band images; the coordinates of candidate flares DC21flly, DC21fjnb,
DC21fygb, and DC21engi are marked as labeled. In all panels, North is up and
East is left. Each panel is 3.8’ on the side. A 12” radius circle centered on
each detection reported in Table 1 is plotted (where only one circle is
visible, as for DC21engi, the detections are overlapping). No source can
unambiguously be identified at the location of the transients, except for
DC21engi.

We selected three objects for further investigation (marked with circles in
Figure 8): the two most extreme outliers in the upper-left quadrant that are
ostensible rising flare candidates, and an object with large zenith-bound
motion but no magnitude change, as a control. All of them show a small change
in airmass between the two observations in a pair, as desired. While a few
DECam DDF transients that passed our cuts are detected at high airmass $X>2$
(Figure 8, left panel), none of them have significant enough magnitude changes
to be considered as flare candidates. All of our three objects of interest are
observed at airmass $X\sim 1.2$. Thus, in order to show a significant
$\Delta\mathrm{DCR}$, if they were flares, these first two candidate events
would have to be extremely hot.The coordinates, detection times, $g$ band
photometry, and $d_{\parallel}$ for these candidates are shown in Table 1. To
ascertain the stellar origin of these transients we look for a quiescent-state
source in the templates.555We used Robert Knop’s DECAT LBL Pipeline Viewer
tool to search by Candidate ID and obtain the 51x51 pixel cutouts of the image
triplets used in the DIA process. The DIA triplets of the three candidates
circled in Figure 8 are shown in Figure 9. From left to right, each triplet
contains the “new” image (also called the “direct” or “science” image,
obtained by the DECam DDF program), the “reference” (a template built from
archival DECam images obtained in previous years), and the “difference” image
(PSF-matched subtraction of the reference image from the new image). For none
of the three candidates were we able to find a precursor source in the
reference image. It is possible that the quiescent sources were not present in
the reference image because they were beyond the survey magnitude limit.
However, the templates reach about one magnitude deeper than the search images
(Graham et al., 2023), and a flare with 1% surface coverage on an M3 dwarf,
would generate a brightening of $\sim$0.62 magnitudes compared to the
quiescent source in $g$ band (Davenport et al., 2012), which should enable the
detection of the quiescent source even for transients detected close to the
single image $5-\sigma$ limit (see Table 1). Since $\mathrm{d}M$ are typically
bright in the infrared (IR) wavelengths, rather than the optical wavelengths,
we also searched for a quiescent source in the IR by examining the WISE data
(Wright et al., 2010) at the candidate’s coordinates (Figure 10), but likewise
did not see a source in these images. This is not surprising: given the
effective magnitude limit of WISE W1 (16.6 mag; Cutri et al. 2012), a faint
transient in the DECam DDF data, such as our three candidates, is expected to
be at the detection limit of the WISE data in quiescent state. This raises the
possibility that the candidates may not be stellar sources, but rather Solar
System Objects or extragalactic transients. We checked additional DECam alerts
to see if the transients were detected at later epochs: extragalactic
transients have typically slower evolutions, so they should remain detectable
in successive nights. Inspecting data with the DECAT LBL Viewer, we found no
subsequent alerts detections for any of the four candidates. Finally, we
checked if any known Solar System Object (SSO) is expected at the coordinates
of the detection by inspecting the IAU Minor Planet Center (MPC) catalogs
using the Minor Planet Checker666https://minorplanetcenter.net/cgi-
bin/checkmp.cgi and we found SSOs within one arcminute (the minimum search
radius in the MPC) of the candidate coordinates, at the time of detection for
all three candidates. This leaves asteroids as the most likely source of
contamination given the brightness of the transient and the fact that each
only appeared in a single night of observations.

Table 1: Candidate ID, Coordinates, timestamps, and $g$ magnitudes with errors, as reported by the DIA pipeline described in subsection 4.2, for each DIA detection of each candidate circled in Figure 8 and DC21engi. The coordinates, detection time, magnitude, and component of the motion in the direction of zenith as measured between the first and each successive observation are indicated. Candidate ID | RA | Dec | Object Datetime | $g_{DIA}$ | $\delta g_{DIA}$
---|---|---|---|---|---
| (∘) | (∘) | (UTC) | (mag) | (mag)
DC21flly | 150.114428 | 2.698361 | 2021-04-15 01:08:45.598 | 23.143 | 0.110
DC21flly | 150.114194 | 2.698476 | 2021-04-15 01:45:47.624 | 22.690 | 0.082
DC21fjnb | 150.384041 | 2.705729 | 2021-04-15 01:08:45.598 | 22.627 | 0.065
DC21fjnb | 150.384017 | 2.706228 | 2021-04-15 01:20:48.720 | 22.610 | 0.076
DC21fjnb | 150.383980 | 2.706525 | 2021-04-15 01:27:32.945 | 22.659 | 0.067
DC21fygb | 150.531830 | 3.533134 | 2021-04-18 00:43:48.466 | 21.833 | 0.041
DC21fygb | 150.531737 | 3.533659 | 2021-04-18 01:35:23.940 | 21.325 | 0.026
DC21engi | 151.678533 | 2.169509 | 2021-04-09 02:23:03.576 | 22.617 | 0.111
DC21engi | 151.678564 | 2.169518 | 2021-04-09 02:41:18.079 | 22.987 | 0.115
Figure 11: Schematic of the position change of candidate DC21engi during a
DECam transient event captured in two detections, at times $t_{1}$ (2021-04-09
02:23:03.576) and $t_{2}$ (2021-04-09 02:41:18.079). The Gaia DR3 coordinates
of the source are shown in black, and the coordinate axes are shifted such
that the Gaia position is at the origin. The source is assumed to be quiescent
at $t_{quiescent}$ (2021-04-09 02:00:00), and the direction to zenith at
$t_{quiescent}$ and the position of the source in the DECam DIA subtraction
are shown in blue. The maximum astrometric residuals (0.15”) are indicated by
the error bars on the DECam event positions. Arrows point to the zenith-bound
direction at the time of the observation, showing the axis along which the
next observation should be found if its motion were dominated by DCR. Figure
12: Center: $\Delta\mathrm{DCR}$ induced as a function of flare temperature
and airmass in LSST $g$-band. Right: Box-and-whiskers plot of peak effective
flare temperatures measured by Howard et al. (2020). Top: histogram of the
per-visit airmass in the current LSST baseline observing strategy
(baseline_v3.0_10yrs, Bianco & The Rubin Observatory Survey Cadence
Optimization Committee 2023). The median airmass and temperature are indicated
by orange lines. Figure 13: Probability of detecting $\Delta\mathrm{DCR}$
induced by a flare at or above peak effective temperature T for four different
astrometric accuracy limits. The probability of DCR detection is measured as
the conditional probability $P(X>X_{crit}|T_{eff}>T)$, where $X_{crit}(T,q)$
is the minimum airmass necessary for a flare of temperature T to produce a
$\Delta\mathrm{DCR}$ greater than the astrometric accuracy $q$. The
temperature distribution follows the measured flare temperatures in Howard et
al. (2020). The small number of flares at large $T$ causes the discontinuities
in the probability at $T_{eff}>20,000K$. Figure 14: Airmass sky distribution
of the current LSST survey strategy proposal (baseline_v3.0_10yrs, Bianco &
The Rubin Observatory Survey Cadence Optimization Committee 2023): skymaps are
shown for the median airmass (_left_) and maximum airmass (_right_) in
$g$-band. Plots produced within the Rubin Metric Analysis Framework. (Jones et
al., 2014). To produce these maps, a sky segmentation into healpixels with
resolution 64 (Górski et al., 2005)

Figure 15: _Top:_ Distribution of observation airmasses over five LSST recent
observing strategy proposals (as labeled). While it remains small, the planned
fraction of high-airmass observations in LSST has grown over time. _Bottom_ :
the corresponding probability of measuring a $\Delta\mathrm{DCR}$: as Figure
13 but fixing $q=0.05$ and for the same five simulations of the LSST as in the
top panel.

We then took a more traditional approach and began with ascertaining if any of
the DECam DDF detections selected by our cuts were indeed flares, then
searched for their location in the $d_{||}$-_vs_ -$\Delta g$ and see if we can
detect a $\Delta\mathrm{DCR}$. To do this, we crossmatched the DECam DDF
source location with a sample of known cool stars in the COSMOS field from the
Gaia DR3 catalog (Brown et al., 2021). We selected Gaia sources with
$G_{bp}-G_{rp}>2$ and $M_{G}>5$ to choose $K-M$ dwarfs and performed an
additional cut on parallax error ($\leq 1^{\prime\prime}$) over three areas
corresponding to the three DECam pointings that compose the DDF’s COSMOS
field. This led to the selection of 6993 Gaia DR3 sources. We then
crossmatched the positions of the Gaia sources to within 1.0” of the initial
1230 DECam DDF candidate positions using the match_to_catalog_3D method of the
astropy.coordinates.SkyCoord object. This crossmatch left us with a single
candidate, DC21engi, whose DIA triplet for two $g$-band detections is shown in
Figure 9 and whose WISE W1 image is shown in Figure 10.

We analyzed the astrometric changes of DC21engi between quiescence and event
by taking the Gaia coordinates as the bonafide quiescent position. In Figure
11 we examine the relative positions of the first and second DIA $g$ band
detections, using the parallactic angle to indicate what the motion of the
source would be if it moved directly toward zenith as expected for DCR. While
the initial detection ($g\sim 22.6$) of the event and the associated motion
from its Gaia position is suggestive of DCR (_i.e._ in the zenith direction),
the second detection ($g\sim 23.0$) defies the expectation that the source
would move back towards its quiescent location as the photosphere cools and
the object position is measured in a direction perpendicular to the
parallactic angle. The bottleneck here is likely the precision of the DECam
astrometric solution, which dominates the uncertainty in the position:
astrometric match requirements for images in Graham et al. (2023) are set by
median astrometric residual across the entire field which must not exceed
0.15”. This margin is shown by the errorbars on the detection positions,
indicating that this candidate does not conclusively demonstrate motion
(produced by DCR or otherwise) during the transient event.

##  5 Technical Recommendations for $\Delta\mathrm{DCR}$ Detection

The difficulties encountered in the process of validating this method on
precursor surveys emphasize the extremely high image quality and astrometric
fidelity requirements demanded by such a method. We propose a set of
astrometric accuracy benchmarks to maximize atmosphere-aided studies,
continuing to use flares as our case study. Figure 12 shows the
$\Delta\mathrm{DCR}$ induced as a function of flare temperature and airmass in
LSST $g$-band, assuming a distribution of flare $T_{eff}$ as measured by
Howard et al. (2020) using a combination of Evryscope (Law et al., 2015) and
TESS data, modeling flares as blackbodies, and a per-visit airmass as in the
baselinev3.0_10yrs simulation of LSST.

###  5.1 Astrometric requirements

. By combining these distributions, we can estimate the probability of
measuring $\Delta\mathrm{DCR}$ in LSST as a function of flare temperature.
Figure 13 shows the conditional probability of producing a
$\Delta\mathrm{DCR}$-induced astrometric shift detectable in LSST during a
flare given the flare temperature distribution in Howard et al. (2020):
$P(X>X_{crit}|T_{eff}>T)$, where $X_{crit}(T,q)$ is the minimum airmass
necessary for a flare of temperature T to produce a $\Delta\mathrm{DCR}$
greater than the astrometric accuracy $q$. The probability of measuring
$\Delta\mathrm{DCR}$ from a flare with peak temperature at least 10,000$K$ is
56.8% for a survey with astrometry accurate to within 10 mas, but this
probability falls to 27.1% at the 100 mas level and 4.1% at the 200 mas level.

New methodologies for increasing astrometric accuracy have been proposed that
could improve upon the 100 mas requirement set in Ivezić & The LSST Science
Collaboration (2018). Fortino et al. (2021) leveraged Gaia astrometry on stars
visible in both surveys and Gaussian process regression (GPR) to reduce the
astrometric variance caused by atmospheric turbulence. Testing on the orbit of
trans-Neptunian object Eris ($r\approx 18.5$) they found that GPR corrections
reduced the root-mean-square (RMS) residuals in $riz$ band from 10 to 5 mas.
The observations used to perform the correction to the Eris orbit were done at
median airmass $X=1.29$ which exceeds the median airmass of
baseline_v3.0_10yrs setting the expectation that similar improvements could be
observed in $\Delta\mathrm{DCR}$ relevant LSST observations. Furthermore, they
state explicitly that the GPR technique should be applicable to LSST, and that
Rubin’s larger aperture compared to DECam means that more stars will be
available for the fit in the turbulence-dominated leading to improved results.
We note that the RMS reduction is proportional to the number of Gaia stars
available in the image, indicating that better results (up to a factor of 5)
can be expected in Galactic fields, where the rate of flares is enhanced (see
Fortino et al. 2021 fig. 5)

Full characterization of the uncertainty in both the quiescent and flaring
position of the source is necessary for DCR-based inference of flare
temperature. If prompt follow up is not required, or even not possible for
very short duration events such as flares, the Rubin annual Data Releases can
be used for the analysis, leveraging the high Rubin’s internal relative
astrometric accuracy (Figure 5). This analysis can also be performed using the
world-public alert packets sent to brokers within 60 seconds of observation
and the Prompt data products released 24 hours after observation (Jurić et
al., 2022). However, these data products do not include positional
uncertainties, thus direct analysis of the images is required.

###  5.2 The impact of observing strategy choices

The observing strategy of Rubin LSST, which is still being finalized at the
time of writing of this manuscript (Bianco et al., 2021), has led to the
current survey strategy recommendation (Bianco & The Rubin Observatory Survey
Cadence Optimization Committee, 2023). The the survey strategy is simulated
via the LSST Operation Simulator (Delgado et al., 2014; Delgado & Reuter,
2016; Naghib et al., 2019); the median and maximum airmass of the LSST current
observing strategy proposal (baseline_v3.0_10yrs) is shown in Figure 14. High-
airmass visits are favorably distributed across the WFD footprint, with visits
above airmass 2.0 occurring at least once anywhere in the WFD footprint.

Four overarching science goals set the survey requirements and the high level
design of the survey strategy (Ivezić et al., 2019): probing dark energy and
dark matter, building an unprecedentedly complete orbital catalog of Solar
System objects, exploring the transient and variable sky, and advancing our
understanding of the Milky Way via the resolved stellar population. The system
is designed to simultaneously maximize field of view and depth, _i.e._ to
maximize the volume of space-time that can be surveyed, leading to a fast
survey that can scan the entire southern sky in $\sim 3$ days, which maximizes
the discovery potential in time domain. Many of the science deliverables of
LSST, however, require high-resolution imaging and exquisite image quality,
which pushes the survey design to prefer low airmass observations.
Nonetheless, in the complex optimization exercise that weighs in the needs of
many science deliverables, it is inevitable that a subset of the observations
will be performed at medium, and even high airmass. Over the course of several
years, as the observing strategy for LSST has evolved under the guidance of
many community contributions (Bianco et al., 2021), we observe that the power
in the high-airmass tail of the distribution of airmasses has increased. As
the optimization of the observing strategy is refined under a complex net of
science goals, additional constraints on the pointing compete with optimizing
the airmass choice for image quality. This tail of high airmass observations,
however, is welcomed by our science case and all atmospheric-aided studies
(Richards et al., 2018; Yu et al., 2020) as it increases the probability of
measuring ($\Delta$)DCR. Figure 15 shows the evolution of the distribution of
airmasses over five recent opsims (simulations of the LSST 10-year survey).
The most recent recommendation for the Rubin LSST observing strategy (Bianco &
The Rubin Observatory Survey Cadence Optimization Committee, 2023) has led to
the most favorable baseline among those considered in this work for DCR-
assisted temperature and SED measurements.

Figure 16: Simulation of the impact of three possible DCR correction
strategies applied in $g$-band DIA. The pixel scale of the LSST system (0.2
arcsec/pixel) is indicated by the grid and the zenith direction by an arrow.
(_a_): a quiescent $\mathrm{d}M$ source observed at airmass $X=1.22$. (_b_):
The source is observed during a $T=10,000$ K flare at the same airmass. In
both panels _b_ and _d_ , the flare flux is normalized such that the total
flux at flare is 50% greater than the total flux in quiescence. The difference
of the flare and quiescent images is shown in the (_c_); contours help the
reader quantify the PSF width and centroid in the difference image. The
quiescent component is subtracted cleanly in DIA and the residuals associated
with the flare directly quantify $\Delta\mathrm{DCR}$. Matching the DCR
templates with correct assumptions on the quiescent SED would replicate this
result. In (_d_), the same flare is observed at airmass $X=1.74$. DCR shifts
the centroid of the quiescent source in the zenith direction by $\sim 42$
arcsec. This shift is corrected in the visualization. The flaring PSF center
is offset from the center of the image by $0.21$ arcsec due to
$\Delta\mathrm{DCR}$, and the PSF is significantly elongated in the zenith
direction (with PSF variance increase $\Delta V=0.14$ arcsec2). (_e_): The
quiescent and flare images are differenced after modeling the DCR effects on
the quiescent source using the true $\mathrm{d}M$ SED. The quiescent
components is removed cleanly and the residuals can be used to measure
$\Delta\mathrm{DCR}$ directly as in (_c_). (_f_): For the flare observed at
airmass $X=1.74$ (_d_), the same DCR correction as for (_e_) are applied, but
based on a G2V (solar type) quiescent SED: the flare and quiescent images are
differenced resulting in non-zero residuals for the quiescent component and
distorted residuals for the flare component.

###  5.3 Data Products Considerations

The astrometric requirements of LSST are sufficiently demanding that DCR
corrections may be required to achieve them (Ivezić & The LSST Science
Collaboration 2018, Swinbank et al. 2020), and understanding how and when
these corrections are applied is necessary to develop a schema for measuring
$\Delta\mathrm{DCR}$ with the Rubin data products. While the position of a
flaring star will be available in an alert prompted by a flare, the centroid
associated with the alert will be measured on the difference image. Difference
imaging is at the core of Rubin’s transient detection, and as Rubin is not
equipped with an atmospheric dispersion corrector, DCR presents a challenge
for image subtraction methods, most notably in the form of “dipoles” in the
subtracted image caused by mis-subtraction of sources observed at different
zenith angles. At the Difference Image Analysis (DIA) stage, DCR-matched
templates can be used to mitigate these dipoles. This method (described in
detail in Sullivan 2018) iteratively forward-models the unrefracted sky and
uses the result as a template for image subtraction. However, in practice,
LSST will not produce enough data in the bluer bands in the first year of
operations to produce DCR-matched templates (Swinbank et al., 2020).

The LSST Science Pipeline does not currently implement wavelength-dependent
PSF modeling (e.g. Meyers & Burchat, 2015) to account for DCR effects,
although this will probably be added early in the survey. DCR produces both a
bulk shift of the PSF centroid (first order effect) as well as a second order
effect on the shape of the PSF, increasing the dispersion along the zenith
direction. Kolmogorov turbulence in the atmosphere leads to a linear isotropic
contraction or dialation of the kernel.

In Figure 16 we explore the impact on DCR correction approaches on the
resulting DIA – and thus our ability to infer flare properties – for a
simulated flaring M dwarf. The initial model (top row) shows the most ideal
scenario, where the T=10,000K flare occurs at the same airmass as a template.
The PSF is modeled as a Gaussian with $FWHM=0.7$ arcsec seeing. As expected,
the resulting DIA shows a bulk shift toward zenith from the blue flare of
$\sim 0.10$ arcsec or $\sim 0.5$ pixels at the LSST system plate scale at an
airmass of $X=1.22$. This model represents perfect DCR correction, since the
template and flare occur at the same airmass, and would allow us to correctly
infer the flare’s temperature in an ideal (i.e. noise free) scenario.

The second scenario in Figure 16 demonstrates the flare occurring at a higher
airmass ($X=1.74$) than the quiescent template. If appropriate DCR correction
for the red $\mathrm{d}M$ SED is applied, then the resulting DIA is again due
to the blue flare flux alone. Note that the DIA shift is larger and more
broadened compared to the first scenario, due to the increased airmass, which
naturally results in a more precise constraint on the flare temperature as
expected (e.g. see Figure 4 and Figure 5). This scenario represents the
pipeline applying perfect DCR corrections for all known source SEDs at every
airmass, making $\Delta\mathrm{DCR}$ measurable from DIA.

The third scenario shown in Figure 16 demonstrates an incorrect or ad-hoc
correction for DCR. Here we have assumed that the flare observation was at a
high airmass, and the field was corrected for some average DCR effect without
any knowledge of the quiescent source SED. For this demonstration, we assumed
the $\mathrm{d}M$ had a DCR correction applied for a solar-type G2V star. The
impact here is an over-correction for DCR for the quiescent red $\mathrm{d}M$
source. The blue flare flux still appears in the resulting DIA in Figure 16f.
However, the resulting DIA signal is distorted from the expectation (Figure
16e) due to the incorrect DCR correction, which would negatively impact the
accuracy of the inferred flare temperature.

This result of these DCR simulations show that 1) the $\Delta\mathrm{DCR}$
signal should be a good indicator in DIA products of flare activity for low-
mass stars at moderate airmasses, even with imperfect DCR corrections, and 2)
careful and accurate DCR corrections based on the quiescent stellar SED are
required to correctly infer flare temperatures using DIA.

##  6 Conclusions

We have proposed a method to extract spectral information from chromatic
transients by taking advantage of the change in differential chromatic
refraction (DCR) between quiescence and event states, or $\Delta\mathrm{DCR}$,
using flares as a case study. We used a composite spectrum derived from SDSS
$\mathrm{d}M$ spectra to demonstrate the capability for LSST to measure the
$\Delta\mathrm{DCR}$ induced by flares with effective blackbody temperature
10,000$K$ at or above $X$=1.2 in $g$ band, enabling $\Delta\mathrm{DCR}$
temperature measurement of flares over a significant fraction of total survey
operations. Depending on the type of data products available for the analysis,
flare temperatures as low as 4000$K$ could also probed by Rubin using DCR. In
order to prepare to deploy this methodology on the Rubin pipeline, we tested
its efficacy on two precursors to LSST, the Zwicky Transient Facility (ZTF,
(Bellm et al., 2018)) and the “Deep Drilling in the Time Domain with DECam”
program (Graham et al., 2023). Our initial tests of the ZTF data quickly
revealed a critical issue with using ZTF as a validation testbed, as the ZTF
pixel scale (1.01”/pixel) is an order of magnitude larger than the astrometric
change we expect to be induced by the $\Delta\mathrm{DCR}$ produced by a
typical flare and the typical seeing exceeds 2”. While the system properties
and observing conditions of the DECam DDF observations are more favorable, and
indeed comparable to the LSST’s, the processing pipeline did not lead to
accurate enough astrometric solutions to enable $\Delta\mathrm{DCR}$
measurements. We estimated how accurate an astrometric solution must be in
order to be able to measure $\Delta\mathrm{DCR}$, given an observed flare’s
peak temperature and found that a survey must be accurate to within 0.01” to
measure 56.8% of $\Delta\mathrm{DCR}$ events produced by flares above
10,000$K$.

This initial exploration of studying stellar flare temperatures using DCR in
the context of Rubin LSST has demonstrated the potential for next-generation
surveys to use incredible advances in image quality, astrometric solutions,
and data throughput to break new ground using novel, albeit unconventional
methods. However, the challenges encountered in our search for
$\Delta\mathrm{DCR}$ in precursor surveys emphasize the strict technical
requirements for such a method to deliver reliable results. Future work on
this project will involve testing this method by injecting artificial flares
into simulated LSST-like images and testing the method’s ability to precisely
recover positional offsets and accurately infer flare temperatures using Rubin
data products.

###  6.1 Acknowledgements

This material is based upon work supported by the National Science Foundation
under Award Number AST-2308016.

JRAD acknowledges support from the DiRAC Institute in the Department of
Astronomy at the University of Washington. The DiRAC Institute is supported
through generous gifts from the Charles and Lisa Simonyi Fund for Arts and
Sciences, and the Washington Research Foundation.

The authors acknowledge the support of the Vera C. Rubin Legacy Survey of
Space and Time Science
Collaborations777https://www.lsstcorporation.org/science-collaborations

and particularly of the Transient and Variable Star Science
Collaboration888https://lsst-tvssc.github.io/ (TVS SC) that provided
opportunities for collaboration and exchange of ideas and knowledge.

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 and AST-2034437 and a collaboration
including current partners Caltech, IPAC, the Weizmann Institute for Science,
the Oskar Klein Center at Stockholm University, the University of Maryland,
Deutsches Elektronen-Synchrotron and Humboldt University, the TANGO Consortium
of Taiwan, the University of Wisconsin at Milwaukee, Trinity College Dublin,
Lawrence Livermore National Laboratories, IN2P3, University of Warwick, Ruhr
University Bochum, Northwestern University and former partners the University
of Washington, Los Alamos National Laboratories, and Lawrence Berkeley
National Laboratories. Operations are conducted by COO, IPAC, and UW.

This project used data obtained with the Dark Energy Camera (DECam), which was
constructed by the Dark Energy Survey (DES) collaboration. Funding for the DES
Projects has been provided by the U.S. Department of Energy, the U.S. National
Science Foundation, the Ministry of Science and Education of Spain, the
Science and Technology Facilities Council of the United Kingdom, the Higher
Education Funding Council for England, the National Center for Supercomputing
Applications at the University of Illinois at Urbana Champaign, the Kavli
Institute for Cosmological Physics, University of Chicago, the Center for
Cosmology and Astroparticle Physics, Ohio State University, the Mitchell
Institute for Fundamental Physics and Astronomy at Texas A&M University,
Financiadora de Estudos e Projetos, Fundacao Carlos Chagas Filho de Amparo a
Pesquisa do Estado do Rio de Janeiro, Conselho Nacional de Desenvolvimento
Cientıfico e Tecnologico and the Ministerio da Ciencia, Tecnologia e Inovacao,
the Deutsche Forschungsgemeinschaft, and the Collaborating Institutions in the
Dark Energy Survey. The Collaborating Institutions are Argonne National
Laboratory, the University of California, Santa Cruz, the University of
Cambridge, Centro de Investigaciones Energeticas, Medioambientales y
Tecnologicas-Madrid, the University of Chicago, University College.

This research has made use of data and/or services provided by the
International Astronomical Union’s Minor Planet Center.

This research has made use of NASA’s Astrophysics Data System Bibliographic
Services.

_Software:_ astropy (Astropy Collaboration et al., 2013, 2022), healpy,
(Górski et al., 2005), Jupyter (Kluyver et al., 2016), MatplotLib (Hunter,
2007), NumPy (Harris et al., 2020), pandas (Wes McKinney, 2010; pandas
development team, 2020), Python3 (Van Rossum & Drake, 2009). SciPy (Virtanen
et al., 2020), and rubin_sim (Yoachim et al., 2023).

## References

  * Astropy Collaboration et al. (2013) Astropy Collaboration, Robitaille, T. P., Tollerud, E. J., et al. 2013, A&A, 558, A33, doi: 10.1051/0004-6361/201322068
  * Astropy Collaboration et al. (2022) Astropy Collaboration, Price-Whelan, A. M., Lim, P. L., et al. 2022, ApJ, 935, 167, doi: 10.3847/1538-4357/ac7c74
  * Becker (2015) Becker, A. 2015, Astrophysics Source Code Library, ascl
  * Bellm et al. (2018) Bellm, E. C., Kulkarni, S. R., Graham, M. J., et al. 2018, Publications of the Astronomical Society of the Pacific, 131, 018002
  * Bertin (2006) Bertin, E. 2006, in Astronomical Data Analysis Software and Systems XV, Vol. 351, 112
  * Bertin (2010) Bertin, E. 2010, Astrophysics Source Code Library, ascl
  * Bertin & Arnouts (1996) Bertin, E., & Arnouts, S. 1996, Astronomy and astrophysics supplement series, 117, 393
  * Bianco & The Rubin Observatory Survey Cadence Optimization Committee (2023) Bianco, F. B., & The Rubin Observatory Survey Cadence Optimization Committee. 2023, Survey Cadence Optimization Committee’s Phase 2 Recommendations.  https://pstn-055.lsst.io/
  * Bianco et al. (2021) Bianco, F. B., Ivezić, Ž., Jones, R. L., et al. 2021, The Astrophysical Journal Supplement Series, 258, 1
  * Borucki et al. (2003) Borucki, W. J., Koch, D., Basri, G., et al. 2003, in Earths: DARWIN/TPF and the Search for Extrasolar Terrestrial Planets, Vol. 539, 69–81
  * Bosch et al. (2018) Bosch, J., AlSayyad, Y., Armstrong, R., et al. 2018, arXiv preprint arXiv:1812.03248
  * Brown et al. (2021) Brown, A. G., Vallenari, A., Prusti, T., et al. 2021, Astronomy & Astrophysics, 649, A1
  * Cutri et al. (2012) Cutri, R., Wright, E., Conrow, T., et al. 2012, Explanatory Supplement to the WISE All-Sky Data Release Products, 1
  * Dark Energy Survey Collaboration et al. (2016) Dark Energy Survey Collaboration, Abbott, T., Abdalla, F., et al. 2016, Monthly Notices of the Royal Astronomical Society, 460, 1270
  * Davenport et al. (2012) Davenport, J. R., Becker, A. C., Kowalski, A. F., et al. 2012, The Astrophysical Journal, 748, 58
  * Davenport (2016) Davenport, J. R. A. 2016, ApJ, 829, 23, doi: 10.3847/0004-637X/829/1/23
  * Delgado & Reuter (2016) Delgado, F., & Reuter, M. A. 2016, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 9910, Observatory Operations: Strategies, Processes, and Systems VI, 991013, doi: 10.1117/12.2233630
  * Delgado et al. (2014) Delgado, F., Saha, A., Chandrasekharan, S., et al. 2014, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 9150, Modeling, Systems Engineering, and Project Management for Astronomy VI, ed. G. Z. Angeli & P. Dierickx, 15, doi: 10.1117/12.2056898
  * Filippenko (1982) Filippenko, A. V. 1982, Publications of the Astronomical Society of the Pacific, 94, 715
  * Flaugher et al. (2015) Flaugher, B., Diehl, H., Honscheid, K., et al. 2015, The Astronomical Journal, 150, 150
  * Fortino et al. (2021) Fortino, W. F., Bernstein, G., Bernardinelli, P., et al. 2021, The Astronomical Journal, 162, 106
  * Goldstein et al. (2015) Goldstein, D., D’Andrea, C., Fischer, J., et al. 2015, The Astronomical Journal, 150, 82
  * Górski et al. (2005) Górski, K. M., Hivon, E., Banday, A. J., et al. 2005, The Astrophysical Journal, 622, 759, doi: 10.1086/427976
  * Graham et al. (2023) Graham, M. L., Knop, R. A., Kennedy, T. D., et al. 2023, Monthly Notices of the Royal Astronomical Society, 519, 3881
  * Harris et al. (2020) Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357, doi: 10.1038/s41586-020-2649-2
  * Hedges et al. (2021) Hedges, C., Luger, R., Dotson, J., Foreman-Mackey, D., & Barentsen, G. 2021, The Astronomical Journal, 161, 95
  * Howard et al. (2020) Howard, W. S., Corbett, H., Law, N. M., et al. 2020, ApJ, 902, 115, doi: 10.3847/1538-4357/abb5b4
  * Hunter (2007) Hunter, J. D. 2007, Computing in Science & Engineering, 9, 90, doi: 10.1109/MCSE.2007.55
  * Ivezić & The LSST Science Collaboration (2018) Ivezić, Ž., & The LSST Science Collaboration. 2018, LSST Science Requirements Document
  * Ivezić et al. (2019) Ivezić, Ž., Kahn, S. M., Tyson, J. A., et al. 2019, The Astrophysical Journal, 873, 111
  * Jones et al. (2014) Jones, R. L., Yoachim, P., Chandrasekharan, S., et al. 2014, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 9149, Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, 91490B, doi: 10.1117/12.2056835
  * Jurić et al. (2022) Jurić, M., Axelrod, T., Becker, A., et al. 2022, Data Products Definition Document.  https://lse-163.lsst.io/
  * Kaczmarczik et al. (2009) Kaczmarczik, M. C., Richards, G. T., Mehta, S. S., & Schlegel, D. J. 2009, The Astronomical Journal, 138, 19
  * Kluyver et al. (2016) Kluyver, T., Ragan-Kelley, B., Pérez, F., et al. 2016, in Positioning and Power in Academic Publishing: Players, Agents and Agendas, ed. F. Loizides & B. Schmidt, IOS Press, 87 – 90
  * Kowalski et al. (2009) Kowalski, A. F., Hawley, S. L., Hilton, E. J., et al. 2009, The Astronomical Journal, 138, 633
  * Law et al. (2015) Law, N. M., Fors, O., Ratzloff, J., et al. 2015, Publications of the Astronomical Society of the Pacific, 127, 234
  * Lupton (2007) Lupton, R. 2007, in Statistical Challenges in Modern Astronomy IV, Vol. 371, 160
  * Meeus (1998) Meeus, J. 1998, Astronomical Algorithms (Willmann-Bell)
  * Meyers & Burchat (2015) Meyers, J. E., & Burchat, P. R. 2015, The Astrophysical Journal, 807, 182
  * Naghib et al. (2019) Naghib, E., Yoachim, P., Vanderbei, R. J., Connolly, A. J., & Jones, R. L. 2019, The Astrophysical Journal, 157, 151, doi: 10.3847/1538-3881/aafece
  * Oliver et al. (2000) Oliver, S., Rowan-Robinson, M., Barcons, X., et al. 2000, Monthly Notices of the Royal Astronomical Society, 316, 749
  * Olivier et al. (2008) Olivier, S. S., Seppala, L., & Gilmore, K. 2008, in Advanced Optical and Mechanical Technologies in Telescopes and Instrumentation, ed. E. Atad-Ettedgui & D. Lemke, Vol. 7018, International Society for Optics and Photonics (SPIE), 70182G, doi: 10.1117/12.790264
  * Osten & Wolk (2015) Osten, R. A., & Wolk, S. J. 2015, The Astrophysical Journal, 809, 79
  * pandas development team (2020) pandas development team, T. 2020, pandas-dev/pandas: Pandas, latest, Zenodo, doi: 10.5281/zenodo.3509134
  * Pettersen (1989) Pettersen, B. 1989, Astronomy and Astrophysics, 209, 279
  * Ramsay et al. (2021) Ramsay, G., Kolotkov, D., Doyle, J. G., & Doyle, L. 2021, Solar Physics, 296, 162
  * Richards et al. (2018) Richards, G., Peters, C., Martin, B., et al. 2018
  * Ricker et al. (2010) Ricker, G. R., Latham, D. W., Vanderspek, R. K., et al. 2010, in American Astronomical Society Meeting Abstracts, Vol. 215, American Astronomical Society Meeting Abstracts #215, 450.06
  * Scoville et al. (2007) Scoville, N., Aussel, H., Brusa, M., et al. 2007, The Astrophysical Journal Supplement Series, 172, 1
  * Sullivan (2018) Sullivan, I. 2018, DCR-matched template generation.  https://dmtn-037.lsst.io/
  * Swinbank et al. (2020) Swinbank, J., Axelrod, T., Becker, A., et al. 2020, Data Management Science Pipelines Design.  https://ldm-151.lsst.io/
  * Tilley et al. (2019) Tilley, M. A., Segura, A., Meadows, V., Hawley, S., & Davenport, J. 2019, Astrobiology, 19, 64
  * Van De Kamp (1967) Van De Kamp, P. 1967, Principles of astrometry
  * Van Rossum & Drake (2009) Van Rossum, G., & Drake, F. L. 2009, Python 3 Reference Manual (Scotts Valley, CA: CreateSpace)
  * Virtanen et al. (2020) Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nature Methods, 17, 261, doi: 10.1038/s41592-019-0686-2
  * Wes McKinney (2010) Wes McKinney. 2010, in Proceedings of the 9th Python in Science Conference, ed. Stéfan van der Walt & Jarrod Millman, 56 – 61, doi: 10.25080/Majora-92bf1922-00a
  * Wright et al. (2010) Wright, E. L., Eisenhardt, P. R., Mainzer, A. K., et al. 2010, The Astronomical Journal, 140, 1868
  * Yan et al. (2021) Yan, Y., He, H., Li, C., et al. 2021, Monthly Notices of the Royal Astronomical Society: Letters, 505, L79
  * Yoachim et al. (2023) Yoachim, P., Jones, L., Eric H. Neilsen, J., et al. 2023, lsst/rubin_sim: v1.3.3, 1.3.3, Zenodo, doi: 10.5281/zenodo.10028614
  * Yu et al. (2020) Yu, W., Richards, G. T., Yoachim, P., & Peters, C. 2020, Research Notes of the AAS, 4, 252, doi: 10.3847/2515-5172/abd6e2