Gaia23bab: A New EXor

On 2023 March 6, the Gaia telescope alerted a 2 mag burst from Gaia23bab, a young stellar object in the Galactic plane. We observed Gaia23bab with the Large Binocular Telescope obtaining optical and near-infrared spectra close in time to the peak of the burst, and collected all public multiband photometry to reconstruct the historical light curve. The latter shows three bursts in 10 years (2013, 2017, and 2023), whose duration and amplitude are typical of EXor variables. We estimate that, due to the bursts, the mass accumulated on the star is about twice greater than if the source had remained quiescent for the same period of time. Photometric analysis indicates that Gaia23bab is a class II source with age ≲1 Myr, spectral type G3−K0, stellar luminosity ∼4.0 L ⊙, and mass ∼1.6 M ⊙. The optical/near-infrared spectrum is rich in emission lines. From the analysis of these lines we measured an accretion luminosity and mass accretion rate (L acc burst ∼ 3.7 L ⊙, Ṁacc burst ∼ 2.0 × 10−7 M ⊙ yr−1) consistent with those of EXors. More generally, we derived the relationships between accretion and stellar parameters in a sample of EXors. We find that, when in burst, the accretion parameters become almost independent of the stellar parameters and that EXors, even in quiescence, are more efficient than classical T Tauri stars in assembling mass.


INTRODUCTION
Photometric variability is a common feature of low mass (< 2 M ⊙ ), young stellar objects (YSOs, e.g.Megeath et al. 2012).Very different timescales are involved in the observed variability: from short term events (from minutes to days) due to magnetic activity, like surface spots, stellar flares and coronal mass ejections, up to long term events (from months to years and centuries) induced by extinction changes due to inner-disk warps, or abrupt variations in the accretion rate.To this latter class of variability events belong the so-called eruptive young stars (EYSs), historically categorized as FU Orionis-type objects or FUors (Herbig 1977), and EX Lupi-type objects, or EXors (Herbig 1989).FUors present powerful outbursts of 3 − 6 mag in the visual band that last from several years to centuries, and take from months to years to reach the peak.Their spectral type depends on the observed wavelength, being of F−G type in the optical to K−M type in the near-infrared (NIR; Hartmann & Kenyon 1996;Audard et al. 2014;Connelley & Reipurth 2018;Fischer et al. 2023).During the outburst, the accretion rate is of the order of 10 −4 −10 −5 M ⊙ yr −1 and the spectra are dominated by absorption lines.EXor bursts have amplitudes of 1−3 magnitudes in the optical, last for a few months or a year, and are recurring (Fischer et al. 2023).Their spectra resemble those of K− or M−type dwarfs rich of emission lines, showing accretion rates in burst of the order of 10 −6 −10 −7 M ⊙ yr −1 .
A dozen of EXors were classified by Herbig (1989), due to their resemblance with the prototype of the class EX Lupi.Since then and until the early 2000s, the number of EYSs candidates has just slightly increased to some tens (Audard et al. 2014) and therefore, eruptive accretion episodes have been considered as peculiar and rare events.This view, however, is now rapidly changing thanks to the discoveries of the Gaia telescope, which in eight years of operations has issued alerts 1 (Hodgkin et al. 2021) of significant photometric changes in the light curve of about 700 known or candidate YSOs.It is interesting to note that the vast majority of the alerted sources undergo photometric variations that do not fit into the two classical EYS categories, while only a few present light curves resembling those of confirmed EXor or FUor variables.Indeed, just three sources have been claimed as bona-fide FUors : Gaia17bpi (Hillenbrand et al. 2018), Gaia18dvy (Szegedi-Elek et al. 2020), and Gaia21elv (Nagy et al. 2023).To these add four further sources that share properties of FUors and EXors : Gaia19ajj (Hillenbrand et al. 2019), Gaia19bey (Hodapp et al. 2020), Gaia21bty (Siwak et al. 2023), andGaia18cjb (Fiorellino et al. 2024, submitted), and four EXors : Gaia18dvz (ESO-Hα99, Hodapp et al. 2019); Gaia20eae (Hankins et al. 2020, Ghosh et al. 2022, Cruz-Sáenz de Miera et al. 2022); Gaia19fct (Miller et al. 2015;Park et al. 2022); and Gaia22dtk (Kuhn et al. 2022).On 2023 March 6, Gaia alerted a 2-mag burst of Gaia23bab (α J2000.0 =19 h 04 m 26 s .68,δ J2000.0 =+04 • 23 ′ 57.′′ 37).Also known as SPICY 97589 (Kuhn et al. 2021), Gaia23bab is a YSO belonging to the small cluster of ∼ 30 members named G38.3-0.9.In Figure 1 we show the 3 × 3 arcmin 2 IRAC band 4 image of the cluster, where are also four Spitzer Dark Clouds (Peretto & Fuller 2009), within which the cluster is partially embedded.The membership of Gaia23bab to G38.3-0.9 has been demonstrated by Kuhn et al. (2023) on the base of the parallaxes and proper motions of six members of the cluster reported in the Gaia DR3 archive (Gaia Collaboration et al. 2016, 2021).The cluster parallax is 1.114 ± 0.056 mas, corresponding to a distance d = 900 ± 45 pc.
The photometric and spectroscopic properties of Gaia23bab are the subject of the present paper.We first describe the spectroscopic observations (Sect.2).Then, the photometric and spectroscopic analysis is presented in Sect.3-4.A short discussion in given in Sect 5, while the summary is presented in Sect.6.We observed Gaia23bab with the 8.4m Large Binocular Telescope (LBT) located at Mount Graham, (Arizona, USA).The long-slit optical spectrum was obtained combining the observations collected on 2023 May 28 and 29, with the Multi-Object Double Spectrograph (MODS, Pogge et al. 2010).MODS observations were done with the dual grating mode (Blue + Red channels, spectral range 350−950 nm) by using a 0. ′′ 80 slit (ℜ ∼ 1500 and 1800 in the Blue and Red channels, respectively).The total integration time was 1800 s.The slit angle matched the parallactic angle to minimize the wavelength dependence of the slit transmission.During the two following nights (2023 May 30 and June 1) the LBT Utility Camera in the Infrared (LUCI, Seifert et al. 2003) was used with the zJ and HK grisms to obtain the 1.0−2.4µm spectrum.We used the G200 low-resolution grating coupled with the 0. ′′ 75 slit, corresponding to ℜ ∼1500.The observations were performed adopting the standard ABBA technique with a total integration time of 1350 s.

OBSERVATIONS AND DATA REDUCTION
Data reduction was done using the Spectroscopic Interactive Pipeline and Graphical Interface (SIPGI, Gargiulo et al. 2022), specifically developed to reduce the LBT long-slit spectra.Data reduction steps of each MODS spectral-image are: correction for dark and bias, bad-pixel mapping, flat-fielding, correction for optical distortions in the spatial and spectral directions, and extraction of the one-dimensional spectrum by integrating the stellar trace along the spatial direction.Spectra of arc lamps have been used for wavelength calibration.Images in the griz bands (spatial scale 0. ′′ 12/px) were obtained to derive the optical photometry of Gaia23bab and calibrate the MODS spectrum, by taking as references all the stars in the field present in the Pan-STARRS (Panoramic Survey Telescope & Rapid Response System)2 catalog (Chambers, et al. 2016).We estimate g=18.97±0.09mag, r=17.11±0.09mag, i=15.51±0.07mag, z=14.74±0.07mag.The inter-calibration between Blue and Red spectral segments was verified by matching the spectral range between 5300 Å and 5900 Å in common between the two channels.
The raw LUCI spectral images were corrected for bad pixels, flat-fielded, sky-subtracted, and corrected for optical distortions in both the spatial and spectral directions.The telluric lines present in the final spectrum were removed by dividing the target spectrum by that of a spectro-photometric standard star observed immediately after the target and corrected for its intrinsic H I recombination lines in absorption.Wavelength calibration was obtained from arc lamps spectra.J, H and K s images of Gaia23bab (spatial scale 0. ′′ 12/px) were acquired to flux calibrate the zJ and HK segments.The photometry was computed based on the 2MASS3 magnitudes of the sources present in the image field.We obtain J=13.2±0.1 mag, H=12.24±0.07mag, K s =11.47±0.07mag.In Figure 2 we show the composite image of a 1×1 arcmin 2 sky area around Gaia23bab obtained from all the LBT images, where MODS g and r images are in blue, MODS i and z images are in green, and LUCI J, H and K s images are in red.Noticeably, Gaia23bab is among the bluest objects in the field.

Light curve
Figure 3 shows the light curve of Gaia23bab during the last fourteen years.Between 2009 and 2014 the optical photometry in the grizy bands was obtained within the Pan-STARRS survey, with a typical sampling of one or two points per year.The Pan-STARRS data, and in particular those in z and y bands, reveal that a burst has occurred before the advent of Gaia, roughly between 2012 April and 2014 June.A more continuous monitoring has been performed since 2014 by Gaia and, starting three years later, by the survey ZTF4 (Zwicky Transient Facility), which covers the gri bands.As already noted by Kuhn et al. (2023) the Gaia light curve shows a ∼ 2 mag burst between April and November 2017 that did not triggered a Gaia alert.Then, Gaia23bab remained quiescent for about six years.During this period the photometric points in the G and r-ZTF bands differ by approximately 1 mag, with a difference from what is observed in other sources, where they roughly coincide (e.g.Cruz-Sáenz de Miera et al. 2022, Nagy et al. 2023).The offset reduces and finally disappears when Gaia23bab starts to brighten.The two filters have similar λ ref (6251.5 Å and 6201.2Å , respectively), but the G filter is broader (G-FWHM = 4396.69Å , r-FWHM = 1397.73Å )5 and it extends much farther toward longer wavelengths.This suggests that if Gaia23bab is redder when fainter, it may not decrease as much in G as in r-ZTF.Alternatively, the same behaviour in the light curve is expected if the G flux of Gaia23bab is contaminated by that of a nearby source of similar magnitude, whose contribution becomes increasingly negligible during the burst.The ability of Gaia to resolve a nearby double star of equal luminosity is 0. ′′ 23 in the along-scan and 0. ′′ 70 in the across-scan direction, independent of the brightness of the primary (de Bruijne et al. 2015).This is higher than the ZTF's spatial sampling6 of 1. ′′ 0/px, and therefore it is unlikely that the G flux is contaminated while the r-ZTF is not.Furthermore, in both MODS and LUCI images the star closest to Gaia23bab is separated by about 2. ′′ 0 (see Figure 2).
The last optical points in the light curve are the MODS grzi magnitudes.Our observations have been obtained at the beginning of the declining phase, when the r-ZTF magnitude was about 0.5 mag fainter than that at the burst peak.In the infrared, Gaia23bab has been monitored since 2014 within the NEOWISE7 survey in the W1 and W2 bands at 3.4 µm and 4.6 µm.The last two bursts have been registered as about 1.8 mag brightening in both bands.In addition, in Figure 3 we plot the JHK s 2MASS and LUCI photometry, to compare the magnitudes in quiescence and burst.Their difference is ∆J =0.85 mag, ∆H =0.64 mag, ∆K =0.74 mag.
In Table 1 we summarize the light curve properties of the three bursts, which appear remarkably similar.First, the amplitude is always comparable to each burst to other (∼ 2 mag in r/G and ∼ 3 mag in g).Second, the time elapsed between the first and the second burst is approximately 50 months, i.e. roughly that between the second and the third (66 months).This frequency, together with the lack of periodicity, has been already observed in EXors (e.g.V1118Ori, Giannini et al. 2020;ASASSN-13db, Sicilia-Aguilar et al. 2017;Gaia19fct, Park et al. 2022).Finally, we have evaluated the burst rising (declining) speed in each band by fitting with straight lines the rising and declining data in the light curve.In all cases the rising speed is of some thousandths of magnitude per day and slightly faster than the declining speed.These values are similar to those observed in V1118 Ori (Giannini et al. 2020), in Gaia20eae (Cruz-Sáenz de Miera et al. 2022), and in V2492 Cyg (Hillenbrand et al. 2013, Giannini et al. 2018).

Color-color diagrams
To investigate the nature of the photometric variability of Gaia23bab, we plot in Figure 4 the optical [g − r] vs. [r − i] and near-infrared [J − H] vs. [H − K s ] color-color diagrams for burst and quiescence phases.In the left panel, we show the optical colors during the 2013 (Pan-STARRS, cyan dot) and 2023 (LBT, blue dot) bursts, while the quiescence colors (red dot) have been derived by averaging the optical photometries between 2009 and 2013.We note that the two bursts have similar colors, being both bluer with respect to quiescence of about 0.4/0.6 mag in [g − r] and 0.7/0.8mag in [r − i].This effect is often seen as a consequence of the dust clearing during the burst (e.g.Hillenbrand et al. 2019) and the increasing contribution to the accretion luminosity at UV wavelengths (e.g.Venuti et al. 2014).The colors of 400−600 Myr young stars (Kraus & Hillenbrand 2007) and those of main-sequence stars are also plotted.Considering  that the stellar spectral type is G−K (Sect.3.4) we estimate an extinction of 5−6 magnitudes during bursts and ∼ 8 mag in quiescence.In the right panel, the blue/red dots are the LUCI and 2MASS photometries, obtained more than 20 years apart.The rectangle is the locus of un-reddened Herbig AeBe (HAeBe) stars of spectral type B−F (Hernández et al. 2005).Considering the error bars, the NIR colors of Gaia23bab marginally fall in this locus if de-reddened by A V ∼ 5−6 mag in burst.Compared to quiescence, the colors of the burst are ∼ 0.2 mag bluer in [J − H] and equal within the errors in [H − K s ], probably because only the internal regions of the disk are significantly heated by the burst event.This view is also supported by the approximate equality of NEOWISE colors between quiescence and burst, being [W1-W2] quiesc = 0.27±0.03mag and [W1-W2] burst = 0.32±0.03mag.In addition to the data shown in the light curve, mid-infrared photometry of Gaia23bab is present in the WISE survey catalogs (AllWISE at 3.4, 4.6, 12.0, and 22.0 µm, and WISE Post-Cryo Database at 3.4 and 4.6 µm), and in the Galactic Legacy Infrared Midplane Survey Extraordinaire (GLIMPSE8 ), which contains Spitzer IRAC photometry at 3.6, 4.5, 5.8, and 8.0 µm.The Spitzer MIPS flux at 24 µm was taken during the MIPSGAL9 survey and retrieved from the Gutermuth & Heyer (2014) catalog.We note that the AllWISE fluxes in W1, W2, and W4 bands are between 0.3 and 0.6 magnitudes brighter than the Spitzer magnitudes in the similar IRAC1, IRAC2, and MIPS1 bands.In principle, this could be explained by the presence of sources near Gaia23bab falling in the large WISE beam.In the WISE coadded images, however, the closest star is ∼ 15 ′′ apart in the SE direction, and thus it might slightly contaminate the Gaia23bab flux only in W4 band, where the FHWM of the WISE beam is 12 ′′ .A more likely explanation is that at the time of the WISE cryogenic survey (December 2009−August 2010), Gaia23bab was at a higher brightness level than in quiescence.This is indeed confirmed by the magnitudes in the W1 and W2 bands obtained in October 2010 during the post-cryo survey, which are all about 0.2−0.3magnitudes fainter than the AllWISE data.While this result is not enough to conclude that a burst occurred between 2009 and 2010, it is certainly evidence of Gaia23bab's mid-infrared variability.

Spectral Energy Distribution
In the far-infrared, Gaia23bab was observed within the Infrared Galactic Plane Survey (Hi-GAL) Herschel keyprogram (Molinari et al. 2010) that surveyed the Galactic Plane with the photometers PACS10 (Photodetector Array Camera and Spectrometer) and SPIRE11 (Spectral and Photometric Imaging REceiver) at 70, 160, 250, 350, and 500 µm.We derived the Herschel magnitudes applying both aperture and gaussian photometry, getting reliable measurements at 160, 250, and 350 µm.At 70 µm the source is barely visible, so that we were only able to estimate an upper limit to the flux, while at 500 µm Gaia23bab is confused within the strong background emission.The Spectral Energy Distribution (SED) is presented in Figure 5.In red and blue are shown the observed photometric points of quiescence and burst12 , respectively.The Pan-STARRS and ZTF points are the average values of the data in the light curve for each state.In the mid-infrared, we have taken the Spitzer and WISE photometries as representative of the quiescence and burst phase, respectively.Furthermore, considering that the amplitude of the variability decreases with λ, for λ > 24 µm we have adopted the same photometric points for computing the bolometric luminosity (L bol ) in quiescence and in burst.
We have first de-reddened the photometries up to 24 µm with A V ∼ 8 mag in quiescence and A V between 4.8 mag and 6.2 mag in burst, based on the estimates of the following spectroscopic analysis (Sect.4.1).L bol has then been derived integrating the area in the plane F λ vs. λ by interpolating with straight lines the SED data.A bolometric correction was applied to take into account the contribution at λ > 350 µm, considering that the emission decreases as 1/λ 2 .Assuming L bol = L * +L acc (L * and L acc being the stellar and the accretion luminosity, respectively) and d = 900 pc, we obtain L bol quiesc = L * +L acc quiesc ∼ 5 L ⊙ and L bol burst = L * +L acc burst = 10.3±3.9L ⊙ .

Stellar properties
The stellar properties of Gaia23bab were derived through the near-infrared photometry and extinction in the quiescent phase.For a distance of 900 pc and A V = 8.0 mag, the absolute J, H, and K s magnitudes are 1.97 mag, 1.67 mag, and 1.49 mag, and the intrinsic [J − H] and [H − K s ] colors are 0.30 mag and 0.18 mag, respectively.Considering the tables of Pecaut & Mamajek (2013) for sources with ages between 5 and 30 Myr, the spectral type is G3−K0, and the effective temperature T eff = 5400±400 K.In Figure 5 we plot the black-body function at T = 5400 K, extincted for A V =8.0 mag (quiescence) and 4.8/6.2mag (burst, see Sect.4.1), to show the consistency with the observed photometric data.The bolometric magnitude, M bol , can be computed as M bol =m(J)+5-5 log 10 d(pc)+BC J , m(J) being the intrinsic J magnitude and BC J the bolometric correction.For a G3−K0 star BC J is 1.08−1.30.Therefore, we get M bol between 3.05 mag and 3.27 mag.Then, an estimate of the stellar luminosity, L * , can be obtained as log 10 L * = 0.4[M bol,⊙ -M bol ], where M bol,⊙ is the bolometric luminosity of the Sun, equal to 4.74 mag (Mamajek et al. 2015).This way, we get L * = 4.0±0.5 L ⊙ .Therefore, having estimated as ∼ 5 L ⊙ the bolometric luminosity in quiescence (Sect.3.3) we obtain L acc quiesc ∼ 1 L ⊙ .Assuming black-body emission, we determined the stellar radius, R * = 1/2T eff 2 / L * /πσ, σ being the Stefan-Boltzmann constant.We find R * = 2.2±0.5 R ⊙ .

SPECTROSCOPIC ANALYSIS
MODS and LUCI spectra of Gaia23bab are presented in Figure 7.The spectrum steeply rises between 3500 Å and 8000 Å likely because of the relevant extinction that affects Gaia23bab (Sect.3.2 and 4.1).As a consequence, in this wavelength range the spectrum is also very noisy, therefore preventing the possible detection of the continuum excess emission (Balmer jump) between 3600 and 4000 Å , which represents a direct signature of accretion (e.g.Alcalá et al. 2014 and references therein).At wavelengths between 0.9 µm and 1.6 µm the spectrum is almost flat and then decreases at wavelengths longer than 2 µm, in agreement with what expected for Class II sources.A number of emission lines are detected, and in particular those originating in the accretion columns, such as H I and He I recombination lines, bright Ca II and O I lines, and weaker metallic lines (Alcalá et al. 2014, and references therein).In addition, emission from the disk (CO 2-0 bandhead, and Na I doublet at 2.2 µm) along with weak emission from outflowing gas (H 2 1-0 S(1) 2.12 µm and He I 1.08 µm, which presents also a blue-shifted absorption component), is detected.No forbidden atomic lines are present in the spectrum.Line fluxes and their 1σ errors were computed using the SPLOT task in IRAF, which takes into account both the effective readout noise per pixel and the photon noise in the spectral region containing the emission line 13 .Fluxes of the main emission lines are given in Table 2.

Extinction during burst
A first A V estimate of 5−6 mag for the 2023 burst has been derived from the optical/near-infrared color-color diagrams (Sect.3.2).An independent measure is based on the ratio of the observed continuum with that of a stellar template of the same spectral type as the target, artificially reddened by varying the value of A V (Alcalá et al. 2021).The best estimate of A V is obtained when this ratio has a flat slope.Although this method is commonly used in classical T Tauri stars, in strong accretors the observed continuum is not only affected by extinction, but can be also significantly enhanced by the excess continuum coming from hot spots in the accretion shock and from the disk.Such excess is described in terms of veiling, r = Flux(excess)/Flux(star), usually estimated by comparing the equivalent width of the photospheric lines in the template spectrum with those in the target spectrum.In Gaia23bab, however, no photospheric lines are detected, either because deep bands/lines are not expected in spectra of G−K type stars (Herczeg & Hillenbrand 2014), or because of a high veiling.We can consider a reasonable range of r basing from literature estimates, however.In classical T Tauri stars, between 6000−8000 Å r is typically between 0 and 2 (Fischer et al. 2011, Alcalá et al. 2021), but it can be ≳ 3 in strong accretors (e.g.Giannini et al. 2022).In this wavelength range r is roughly constant or decreases as a black-body at T ≈ 8 000 K (considered as the temperature of the hot spot 13 The effective readout noise per pixel was measured as the root mean square deviation (rms) of the continuum in either side of each line.
The average rms was then set as parameter 'σ 0 ' in SPLOT.The photon noise was estimated as invgain*I, where invgain is the reciprocal of the MODS/LUCI gain (2.5/2.0 e − /ADU expressed in physical units) and I is the pixel value.The error on the profile fit is then computed by a Monte-Carlo simulation, whose iteration number nerrsample was set to 100.giving rise to the optical veiling, Fischer et al. 2011).To apply the procedure described above to compute A V , we adopted as a first approximation a constant r between 0 and 6 and added the relative excess continuum to the stellar template prior than applying the variable reddening.We used the optical templates for stars with T eff = 5000−5800 K of Gonneau et al. (2020), and considered a grid of A V between 0 and 20 mag, adopting the extinction law of Cardelli et al. (1989) and a total-to-selective extinction ratio R V = 3.1.This way, we estimated A V as a function of r, as shown in Figure 8.We note that for small values of r (r ≲ 3), A V decreases with r as expected if we consider that extinction effects tend to diminish the continuum, while veiling acts in the opposite direction.As r increases, and in particular when the intensity of the excess continuum prevails over the stellar one, the spectral shape no longer changes, and A V tends to an asymptotic minimum.From this plot we derive A V = 5.5±0.7 mag.Finally, we have also checked how the derived A V values change assuming that the veiling follows a black-body law at T = 8000 K.The results do not change significantly with respect to the case of a constant veiling.

Accretion properties
A rough measure of the accretion luminosity due to the 2023 burst event can be derived as L bol burst -L * = 6.3±3.9L ⊙ .A more accurate method relies on the empirical relationships found by Alcalá et al. (2014Alcalá et al. ( , 2017)), between the accretion luminosity, L acc , and the luminosities L i of selected emission lines of gas in the accretion columns.In the optical range these relationships exist for more than 20 lines, namely the H I recombination lines of the Balmer and Paschen series from Hα to H15, and from Pa8 to Pa10, along with He I , O I and Ca II lines.In the near-infrared, there are relationships for Paδ, Paγ, Paβ, and Brγ.In the Gaia23bab spectrum there are 24 lines useful for determining L acc .For a fixed A V , the best estimate of L acc is that for which the dispersion among the individual L acc(i) is minimized.The associated error is the combination of the uncertainties on the line fluxes and those on the relationships between L acc and L i (Alcalá et al. 2017).Typically, this error does not exceed a few tenths of solar luminosity, negligible with respect to that induced by the uncertainty on A V .For A V = 5.5±0.7 mag we obtain L acc burst = 3.7±1.8L ⊙ .From L acc burst , M * , and R * , we derived the mass accretion rate during the burst as Ṁacc Gullbring et al. 1998), where R in is the inner-disk radius, assumed ∼ 5R * (Hartmann et al. 1998), and G is the gravitational constant.We obtain Ṁacc burst = (2.0±1.0) 10 −7 M ⊙ yr −1 , in line with the mass accretion rate values found in classical EXor events (e.g.Audard et al. 2014).From the same relation, and adopting the value of L acc quiesc derived in Sect.3.4, we obtain Ṁacc quiesc ∼ 6 10 −8 M ⊙ yr −1 .All the stellar and accretion parameters are summarized in Table 3.

The role of bursts in assembling mass
One of the fundamental questions about EYSs is the role of bursts for the stellar mass assembly (Fischer et al. 2023).We have estimated that during the 2023 burst, the luminosity due to accretion was comparable to the stellar (2.0±1.0) 10 −7 luminosity (Table 3).In the last 10 years Gaia23bab has undergone three bursts, lasting approximately one year each (Table 1).Furthermore, they have similar amplitude and all present a 'triangular' shape in the light curve.This allows us to assume as the average value of the mass accretion rate throughout each event, half the value of Ṁacc burst measured close to the peak of the 2023 burst.We obtain that over the last ten years, three of which were spent in burst, the mass accumulated on Gaia23bab was ∼ 1 10 −6 M ⊙ , namely roughly twice the mass that the star would have assembled if it had remained quiescent for the same period of time.Notably, this estimate is quite similar to that evaluated for the 2022 burst of EX Lupi (Cruz-Sáenz de Miera et al. 2023).Considering that the burst amplitude is thought to decrease with time (e.g.Fischer et al. 2023), we speculate that the contribution of burst episodes to the final stellar mass could be even higher than our estimate.

Gaia23bab in the context of EXors
In the previous sections, we have analysed the photometric features of Gaia23bab to show their similarities with EXor sources.Here, we focus on the parameters derived from the emission lines, namely the accretion luminosity and the mass accretion rate.In Figure 9 we plot L acc vs. L * and Ṁacc vs. M * for a sample of 13 known EXors, both in quiescence and in burst.For comparison, the loci of accreting T Tauri stars (Alcalá et al. 2017) and of low-mass HAeBe (Wichittanakom et al. 2020) stars, are shown with a grey and green shaded areas, respectively.For the large majority of the sources, L acc was measured in the same way, namely applying the Alcalá et al. (2017) relations between optical/NIR emission lines and accretion luminosity.Both in quiescence and in burst we find a tight relation between L acc and L * , that spans more than two orders of magnitude in L * .A good correlation is also found between Ṁacc and M * but with a larger spread that is likely due to the uncertainty in the stellar mass and radius determinations.The values measured in Gaia23bab during the burst are in agreement with both fits, if we take into account that our spectra have been taken at the beginning of the declining phase, when the r magnitude was already increased by about 0.5 mag with respect to the peak traced by the ZTF data in the same band (Sect.3.1).Instead, the accretion properties estimated in quiescence are consistent with the expected values.We recall that they have been indirectly derived from the SED, and thus a quiescent spectrum would be necessary to derive accurate determinations of L acc quiesc and Ṁacc quiesc .Also, this will allow us to investigate whether a variable extinction has a role in the different brightness phases.More in general, we note that: 1) the angular coefficients of the relations derived in quiescence are consistent within the uncertainties with those of T Tauri stars.Instead, the intercepts differ by about an order of magnitude.If confirmed on the base of a larger sample, this result would imply that EXors, even in quiescence, are more efficient than T Tauri stars in accreting mass; 2) when in burst, the relations between accretion and stellar parameters become shallower, as if there were a limit to the amount of material that can be transferred from the disk to the star through the accretion columns.Given the low number of sources for which the accretion vs. stellar properties have been derived, and the fact that the high-mass regime is represented only by PV Cep data, these results need to be confirmed on a larger statistical base, however.

SUMMARY
In this paper we have presented LBT observations of the ∼ 2 mag burst of Gaia23bab, a YSO alerted by Gaia on March 2023.Our results, derived from the analysis of photometric and spectroscopic data, can be summarized as follows: 1.The multi-wavelength light curve shows that Gaia23bab had three bursts in the last ten years.These bursts are quite similar to each other in amplitude, duration, and rising/declining speed, in line with those observed in EXors.
3. The optical/NIR spectrum is rich in emission lines from which we have measured the accretion luminosity and the mass accretion rate during the burst.We get L acc burst ∼ 3.7 L ⊙ , comparable with the stellar luminosity.The mass accretion rate close to the burst peak is Ṁacc burst ∼ 2.0 10 −7 M ⊙ yr −1 .
4. The mass accumulated on Gaia23bab in the last ten years was roughly twice the mass that it would have assembled if it had remained quiescent for the same period of time.
5. Both accretion luminosity and mass accretion rate of Gaia23bab are consistent with those of confirmed EXors.
As a more general result, we have quantified the correlations, both in quiescence and in burst, between accretion and stellar parameters in a sample of 13 EXors.On average, EXors have L acc and Ṁacc larger than T Tauri stars in the same range of mass, even in quiescence.If confirmed on the base of a larger sample, this result would imply that EXors are more efficient than T Tauri stars in accreting mass.Also, when in burst, accretion luminosity and mass accretion rate become poorly dependent on the central star properties.

Figure 1 .
Figure 1.IRAC band 4 (8.0 µm) image of a 3×3 arcmin 2 sky area containing Gaia23bab (marked in black).White circles are the YSOs of the cluster G38.3-0.9, and the labels indicate the location of the Spitzer Dark Clouds (SDC) from Peretto & Fuller (2009) catalog.The surface brightness scale is given in MJy/sr in the right bar.

Figure 2 .
Figure 2. LBT 3-color image of 1×1 arcmin 2 sky area centered on Gaia23bab.The color code is the following.Red: J+H+Ks (LUCI); green: i+z (MODS); blue: g+r (MODS).Note that Gaia23bab is among the bluest objects in the field.

Figure 3 .
Figure 3.Light curve of Gaia23bab.Different symbols are the photometric points obtained within different surveys (filled squares: Gaia (G), crosses: Pan-STARRS; open squares: ZTF, filled circles: 2MASS, asterisks: NEOWISE, filled triangles: LBT).Different colors indicate different filters, as indicated on top.The calendar date is indicated as well.The dotted black line marks the date of the LBT spectroscopy.

Figure 4 .
Figure 4. Left panel: two-color optical plot [g − r] versus [r − i].Red and blue/cyan dots indicate quiescent and burst points (LBT/Pan-STARRS data).The black dashed-dotted line is the locus of young stars of 400−600 Myr (Kraus & Hillenbrand 2007), while the green dashed line represents the main-sequence stars with the positions of the A0, G0 and K0 spectral type stars indicated.The arrow represents the direction of the extinction vector corresponding to AV = 2 mag (reddening law of Cardelli et al. 1989).Right panel: two-color near-infrared plot [J − H] vs. [H − K].Quiescent (2MASS) and burst (LUCI) data are colored in red and blue, respectively.The dashed-dotted black rectangle is the locus of the HAeBe stars(Hernández et al. 2005), while the green dashed line represents the main-sequence stars with the positions of the A0, G0 and K0 spectral type stars indicated.The arrow represents the direction of the extinction vector corresponding to AV = 2 mag (reddening law ofCardelli et al. 1989).

Figure 5 .
Figure 5. Spectral Energy Distribution (SED) of Gaia23bab using data at different epochs.In red and blue are shown the data (not corrected for extinction) in quiescence and in burst, respectively.The down triangle is the 2σ upper limit at 70 µm.Photometric points from different catalogs are plotted with different symbols.The red and blue dashed lines are the black-body function at T = 5400 K, extincted for AV = 8.0 mag and 6.1 mag, respectively (see Sect.3.2,Sect.3.4,and Sect.4.1).

Figure 7 .
Figure 7. MODS and LUCI spectra of Gaia23bab.The photometric points in griz and JHKs bands are shown with asterisks.The main spectroscopic emission lines are labeled.

Table 1 .
Features of the bursts of Gaia23bab.

Table 2 .
Fluxes of the main observed linesNote-The line wavelength is in air/vacuum for λ </> 1 µm, respectively.a In parenthesis we give the flux of the absorption component of the line.

Table 3 .
Stellar and accretion properties