Mass Accretion, Spectral, and Photometric Properties of T Tauri Stars in Taurus Based on TESS and LAMOST

We present the analysis of 16 classical T Tauri stars (CTTSs) using LAMOST and TESS data, investigating spectral properties, photometric variations, and mass accretion rates. All 16 stars exhibit emissions in Hα lines, from which the average mass accretion rate of 1.76 × 10−9 M ⊙ yr−1 is derived. Two of the stars, DL Tau and Haro 6-13, show mass accretion bursts simultaneously in TESS, ASAS-SN, and/or the ZTF survey. Based on these observations, we find that the mass accretion rates of DL Tau and Haro 6-13 reach their maxima of 2.5 × 10−8 M ⊙ yr−1 and 2 × 10−10 M ⊙ yr−1, respectively, during the TESS observation. We detect 13 flares among these stars. The flare frequency distribution shows that the CTTSs’ flare activity is not only dominated by strong flares with high energy but also much more active than those of solar-type and young low-mass stars. By comparing the variability classes reported in the literature, we find that the transition timescale between different classes of variability in CTTSs, such as from stochastic (S) to bursting (B) or from quasi-periodic symmetric to quasi-periodic dipping, may range from 1.6 to 4 yr. We observe no significant correlation between inclination and mass accretion rates derived from the emission indicators. This suggests that inner disk properties may be more important than those of outer disks. Finally, we find a relatively significant positive correlation between the asymmetric metric M and the cold disk inclination compared to the literature. A weak negative correlation between the periodicity metric Q value and inclination has also been found.


INTRODUCTION
Classical T Tauri stars (CTTSs) are young, low-mass pre-main-sequence objects with typical spectral types between F and M, that are accreting material from the circumstellar disk.These stars exhibit strong emission lines in their spectra and show an excess of emission ranging from radio to ultraviolet relative to the stellar photosphere.The infrared excess emission from CTTS indicates the presence of a circumstellar disk (Mendoza V. 1968;Lynden-Bell & Pringle 1974;Adams et al. 1987;Bertout et al. 1988) and can be used to infer the amount of dust and the evolutionary status of the system.The amount of circumstellar material drops as the protostellar system evolves, leading to the decrease of the excess emission in the infrared.
The CTTSs have luminosity variations on different timescales, from hours to even decades (Rucinski et al. 2008;Cody et al. 2013), depending on various physical phenomena affecting the stellar and circumstellar environment.For example, the obscuration of the central star's photosphere by circumstellar dust may lead to an aperiodic dip of light (e.g., McGinnis et al. 2015).The presence of stable starspots on the stellar surface may result in a periodic modulation in the light curve of the star (e.g., Kóspál et al. 2018).The unsteady accretion processes and magnetic activity (i.e., flares) may cause the random and short timescale brightening variability.The stochastic light curve of the CTTS may be attributed to the superposition of the various types of variability resulting from different mechanisms mentioned earlier.All these mechanisms can vary with time.
One of the major sources of variability in CTTSs during this evolutionary phase is the accretion-driven burst (e.g., Fischer et al. 2022).There are various types of accretion-driven bursts that last for different timescales.One type is the FU Ori-like outburst, which is characterized by a large increase in luminosity (typically by a factor of tens) over a timescale of decades to hundreds of years (e.g., Clarke et al. 2005).In addition, the outbursts with the type of V1647 Ori, EX Lup, and Protostellar share the similar timescale of a few years with a luminosity increase ranging from 1 to 5 mag.Lastly, the one we focus on in this study is the accretion bursts with a luminosity variability up to a few magnitudes on a relatively short timescales from a few hours to days at optical wavelengths (e.g., Cody & Hillenbrand 2018).The numerical simulations carried out by (e.g., Čemeljić & Siwak 2020) suggest that this variability could be caused by unstable accretion processes.In addition, such events, coupled with rotational modulation, tend to dominate changes in the accretion rate (Rigon et al. 2017;Sergison et al. 2020).
According to the magnetospheric accretion mechanism (e.g., Camenzind 1990), the most viable scenario currently used to describe the accretion behaviors of CTTSs, the central star is channeling the material in the circumstellar disk to the stellar surface through the magnetospheric accretion columns (Feigelson & Montmerle 1999;Muzerolle et al. 2003).The gas is heated to a temperature of 10 4 -10 6 K when the channeled materials hit the stellar surface, creating a shock on the stellar photosphere, which emits intense X-rays (Argiroffi et al. 2007;Günther et al. 2007;Costa et al. 2017).The X-rays are mostly reabsorbed by the accretion columns, from which the lower energy photons, mainly in blue band, are re-emitted, resulting in the luminosity excess from visible to ultraviolet wavelength.Thus, the U-band excess has been often used to determine the accretion luminosity, from which the mass-accretion rate can be determined (Gullbring et al. 1998;Herczeg & Hillenbrand 2008).In addition, such an accretion mechanism also produces detectable emission features in CTTSs spectra, such as Hα and Hβ (e.g., Basri & Bertout 1989).Sousa et al. (2016) shows that when U-band excesses are not available, mass-accretion rates can be obtained with good results from the Hα flux using the currently available calibrations from the literature.Therefore, these emission features can be used to probe the state (i.e., active or inactive) of the mass-accretion process.Furthermore, multi-band simultaneous observations can be useful for estimating the mass-accretion rate of CTTSs even without direct U-band measurements.For example, Tofflemire et al. (2017) and Kóspál et al. (2018) used optical wavelength observations of accretion burst variability to estimate the mass-accretion rate of DQ Tau.In Kóspál et al. (2018), they used the Kepler light curve to identify the burst phenomenons and extrapolated the associated U-band excess luminosities based on the simultaneous BVRI observations of the star.
In this study, we present mass-accretion rates and the properties that can be determined with our data for selected CTTSs in the Taurus Molecular Cloud (hereafter, TMC), which is a nearby star-forming region located at a distance about 147 pc from the Sun (Zucker et al. 2020) with an age of 1-3 Myr (Simon et al. 2017), using high photometric precision light curves and low-resolution spectra from the space-based telescope, TESS (the Transiting Exoplanet Survey Satellite) and LAMOST (the Large Sky Area Multi-Object Fibre Spectroscopic Telescope), respectively.Given the lengths of the available light curve data, we will primarily focus on the short timescale burst.
The outline of this paper is as follows.In Section 2, we introduce the TESS and LAMOST observations of TMC members, including their stellar parameters reported in literature.In Section 3, we explain our methodology for estimating the mass-accretion rates of our CTTS targets using various emission indicators from LAMOST low-resolution spectra and TESS light curves.We also detail our method for determining photometric variabilities in this section.In Section 4, we describe our results of the data analysis and implications.Finally, in Section 5, we summarize the result and conclude the study.

SAMPLE AND DATA SETS
2.1.LAMOST spectroscopic survey on T Tauri stars LAMOST, the Large Sky Area Multi-Object Fibre Spectroscopic Telescope located at the Xinglong Observatory of the National Astronomical Observatories of Chinese Academy of Sciences (NAOC), is a reflecting Schmidt telescope with an effective aperture of ∼ 4 m and a wide field of view of 20 deg 2 in the sky (Zhao et al. 2012).
With 4000 fibers that are mounted on the focal plane, it has been conducting spectroscopic surveys efficiently for stars, galaxies and quasars mostly in the northern celestial hemisphere above the declination of −10 degree since Autumn 2012.Each fiber has a scale of 3.3 arcseconds on the focal surface.The centers of the spectroscopic fibers are about 4.47 arcminutes apart, and each fiber can move at most 3 arcminutes from its center position.Sixteen spectrographs are used in the system, and each receives the dispersed light from 250 fibers.The resolution of LAMOST spectra is R ≈ 1800 around the sdss-g band with a wavelength coverage of 3700-9100 Å.
The LAMOST 2D pipeline is utilized for data reduction.Initially, each exposure's data from each CCD chip is separately reduced, and then the results from each exposure are combined.Bias and dark are subtracted from each raw image.The flat-field spectrum of each fiber is traced, and a polynomial is fitted to the centroid of each fiber's row position.A Gaussian-like profile is assumed for each fiber's profile in the row direction.Using the same trace function, spectra are extracted from all data on the same night, with minor shifts made to align individual frames.The extracted spectra are divided by the one-dimensional flat-field spectrum.The Hg/Cd and Ne/Ar arc lamp spectra are also extracted using the same trace function, including strong sky emission lines to further improve the wavelength solution.The sky spectra are obtained from 20 sky fibers allocated in each spectrograph.After wavelength calibration, all sky spectra are combined into one super-sky spectrum.The super-sky spectrum is scaled to best fit the sky spectrum in each fiber, and then this scaled super-sky is subtracted from each fiber, and the telluric absorption is removed.The science spectra are flux-calibrated by matching the selected flux standard stars and their templates, which are usually of spectral type A or F. Up to five flux standard stars are observed per spectrograph.
We retrieved spectroscopic data of CTTSs in the TMC from LAMOST's decade-long archive of observational data, and 16 targets are identified (see Table 1).The LAMOST low-resolution spectra of these CTTSs are shown in Figure 1.Table 2 contains the stellar parameters of our targets.The masses and radii of most CTTSs are provided in Herczeg & Hillenbrand (2014), with the exception of HL Tau.The radius (R=7 R ) and mass (M=1.7M ) of HL Tau are given by Liu et al. (2017) andPinte et al. (2016), respectively.The inclinations of the outer disk have been reported by the adaptive optics (i.e., Close et al. 1998), millimeter radio (i.e., Guilloteau et al. 1999), high-resolution near-infrared imaging (Kudo et al. 2008), and high-resolution submillimeter resolved image observations of the TMC region (ALMA Partnership et al. 2015;Long et al. 2019;Manara et al. 2019;Villenave et al. 2020).The distance information for the objects studied in this work is from the parallax measurements reported in Gaia Collaboration (2020).HK Tau is a binary system consisting of two companions, HK Tau A and HK Tau B, with HK Tau A being much brighter than HK TauB (G=14.106vs G=17.962 (Gaia Collaboration et al. 2021)).The disk inclinations of the two companions are highly misaligned, with HK Tau A having an inclination of i=56.9 o (Manara et al. 2019) and HK Tau B having an inclination of i=84 o (Villenave et al. 2020).Therefore, we list their inclinations together in the Table and treat them as an uncertainty in our result comparison in the later sections.

TESS observations of T Tauri stars in TMC
TESS (the Transiting Exoplanet Survey Satellite) is an MIT-led NASA mission launched on April 18, 2018, to search for nearby transiting exoplanets using the fast-sampling (2-minute) time-series data with the high photometric precision, which is also well suited for studies of stellar physics.The Taurus Molecular Cloud was in the field of view of TESS in Sectors 43 and 44 from 2021 September 16 to 2021 November 6. Fourteen of our LAMOST CTTS targets have released TESS 2-min cadence light curve data (see Table 1).
There are two kinds of flux information contained in the TESS data: Simple Aperture Photometric flux (SAP) and the Pre-search Data Conditioning Simple Aperture Photometric flux (PDCSAP).The SAP flux is produced by summing the count rates of the pixels within a customized aperture given by the Science Processing Operations Center (SPOC) pipeline (Jenkins et al. 2016).The SAP flux is often affected by systematic errors such as long term trend, which have been removed from the PDCSAP flux using Co-trending Basis Vectors.Therefore, the PDCSAP flux is often cleaner than SAP flux.In this study, we use the PDCSAP flux data with good quality (flag bits QUALITY=0) in our analyses.Although the long term trends might be the true signal of the stars, the unsteady accretion flows have relatively short timescales and should not be affected significantly by removing the long term trend signal.
Considering the large size of the TESS pixel scale (21 ×21 ), the contamination due to the nearby sources in crowded regions may be unavoidable.To investigate whether the light curves of our targets of interest are contaminated by the neighboring stars, the Python package Lightkurve (Lightkurve Collaboration et al. 2018) were employed.We downloaded the target pixel file (TPF) of each star and plotted it together with the coordinates of the field stars with G RP mag < 15 queried from Gaia DR3 (Gaia Collaboration 2020).The left panel in Figure 2 is an example, the TPF of GG Tau (TIC 245830955).The pink mask region is the SPOC pipeline aperture, the red hollow circles are the locations of the field stars, and GG Tau is marked by a red filled circle.In addition to GG Tau, there is one neighboring star within the SPOC aperture.The neighbor star (G RP ∼ 13) is much fainter than GG Tau (G RP ∼ 10).Also, the overall shapes of the light curves extracted using various test apertures do not change significantly compared to the light curve from SPOC aperture.As the result, the SPOC given PDCSAP flux of GG Tau is reliable for the following research.The right panel in Figure 2 shows another example, GI Tau (TIC 268444803).There are three neighbor stars within or at the edge of the SPOC aperture of GI Tau.One of them is GK Tau (G RP = 10.89) and is about equally bright to GI Tau (G RP = 11.5).After the examination with different test apertures, we found that it is impossible to distinguish the variability feature between GI Tau and the field stars in the light curve.Therefore, the GI Tau data was excluded from following analysis.In the same way, HL Tau has a similar problem with GI Tau and was also excluded.Ultimately, there were only 12 targets for which the data could be used in the subsequent analyses.Due to the limited number of spectral standard stars in the LAMOST survey, no absolute flux-corrected data are available, only relative flux-corrected data can be used.To obtain the physical flux information, we adopt the template spectra of weak-lines T Tauri stars (WTTS) with the spectral type ranging from K5 to M7 provided by Manara et al. (2013).We first normalize the template spectra and LAMOST spectra using the mean flux between 7750 and 8000 Å because this range has relatively few emission/absorption features.In this way, we find the best-fit template spectra with the corresponding spectral types for every LAMOST spectra.We then convert the normalized flux of the template spectra and LAMOST spectra back to the physical flux level of template spectra before normalization.Finally, the absolute calibrated flux density of LAMOST data (F λ,cal ) is obtained by using the following equation: where F λ,tem is the flux density of the best-fit template spectra, the d 2 is the distance in parsecs of the template star given by Manara et al. (2013), and d 1 is the distance in parsecs of our sample.
To ensure the accuracy of our calibrated absolute flux measurements, we have cross-checked our results with photometric data obtained by ASAS-SN (Jayasinghe et al. 2018(Jayasinghe et al. , 2019) ) and Pan-STARRS (Chambers et al. 2016;Flewelling et al. 2020) surveys if any, which were taken within 20 days of the LAMOST observation date.As CTTSs are known to be highly variable stars, this additional step provides further validation of our flux calibration.In these surveys, we have identified the following stars in our samples with the photometric data that meet the criteria: (1) from ASAS-SN: AA Tau in the sdss-g band and BP Tau in the V band, (2) from Pan-STARRS sdss-g band: DL Tau, DN Tau, Haro 6-13, HL Tau, and IQ Tau, and (3) from Pan-STARRS sdss-i band: GO Tau and HK Tau.The flux (in mJy) and observation date (in MJD) of the these photometric data are listed in Table 3.
We convert the flux units of the survey photometric data from mJy to the same units as our calibrated LAMOST spectra.Next, we scale the calibrated LAMOST spectra to match the flux from the surveys at the effective wavelength of the corresponding band.This allows us to obtain a more accurate calibrated flux for our LAMOST data.Figure 3 shows the example spectra of the flux calibration processes of three stars: DL Tau, AA Tau, and GO Tau.For stars that lack survey photometric data, we still use their spectra in the following analysis.However, we note that the results obtained from these data could be uncertain by an order of 0 to 1.

Mass-accretion rate measurements from LAMOST spectra analysis
Using the LAMOST low-resolution spectral data, we are able to estimate the mass-accretion rates of T Tauri stars by analyzing various permitted emission lines of Hα 6563 Å, Hβ 4861 Å, CaII 8498, 8542, 8662 Å, and HeI 4471, 5875, 6678, 7065 Å (Fang et al. 2009;Herczeg et al. 2009;Alcalá et al. 2014Alcalá et al. , 2017)), which are generated by the accretion flows between the disc and star, in the available wavelength coverage.The emission profiles of these lines in all spectra are shown in Figure 4, Figure 5, Figure 6, and Figure 7.
For each spectrum, we first estimate the equivalent width (EW) of the emissions with the equation: where F λ is the spectral flux of the line defined to be within the 3-σ of the best-fit Gaussian model with the central wavelength of the particular line; F o is the average flux of the continuum in the same range.The equivalent widths of all available emissions in the spectra are listed in Table 5.Then, the emission's flux (F em ) as well as their luminosities (L em ) can be obtained using the following equations: Here, the EW em represents the equivalent width of the emission, F c (em) is the local continuum flux in the emission's wavelength range, and d is the distance to the stars.The mass-accretion luminosity (L acc ) of T Tauri star can be determined from the line luminosity (L line ) by using the relation below: where L is the solar luminosity.The coefficients a and b of every line used in the study is given by Alcalá et al.
(2017) (see Tabel B.1. in their paper).If the accretion energy is reprocessed fully into the accretion continuum, the L acc can be converted into the mass-accretion rate ( Ṁacc ) using the following equation from Herczeg & Hillenbrand (2008): where M * is the stellar mass, G is the gravitational constant, R * is the stellar radius, and R in is the inner radius of the disk, which is defined as the disk truncation radius of CTTSs with a standard value of 5R * (Gullbring et al. 1998).
The estimated mass-accretion rates from all available emissions are displayed in Table 6.

Variability of our CTTSs sample
A statistical metric system reported by Cody et al. (2014) allows the light curves of young stellar objects to be classified into different categories, thus it is useful to recognize the burst-like accretion events in the TESS light curves of our samples.The system contains two metrics: the flux asymmetry (denoted as M ), and the quasi-periodicity (denoted as Q).The metric M is the first one we focus here, which represents the light curve behavior tendency of dipping (M > 0) and bursting (M < 0).We determine the metric M following the approach explained in detail by Cody & Hillenbrand (2018) (see Sec. 5. in their paper).We first normalize the light curve with the median of the flux and smooth it with 200 points moving mean.After subtracting the smoothed light curve from the original one, we identify the 5σ outliers in the residual curve and remove their counterparts from the original curve.Then the metric M is estimated using the formula: where f 5% is the mean of the first and last 5 percent of the fluxes in the outlier-free curve, and f med and σ f are the median flux and the rms of the original curve, respectively.Stars with M < −0.25 are likely to display the accretion-like burst variability in their light curves.The M -values of all of these CTTSs are listed in Table 4.
The second metric Q is used to evaluate periodicity tendency of light curves in the range from Q=0 (purely periodic) to Q=1 (stochastic) (Cody et al. 2014).To determine the Q-values of the light curve, we first fold up the light curves with the associated periods reported by Percy et al. (2010) and Rebull et al. (2020) if any.One special case is DQ Tau, which is a spectrocopic binary and has two reported brightness variation periods: ∼ 15.8 day and ∼ 3.03 day.The first one is attributed to the brightening effect associated with mass-accretion and magnetic reconnection that occurs around the periastron when two stars approach each other (Mathieu et al. 1997;Salter et al. 2010;Czekala et al. 2016;Getman et al. 2016).The second one is the rotational modulation caused by the presence of the spots, and we use it to determine the Q-value for DQ Tau.Since the morphologies of variability in young stellar objects could vary with time (e.g., Sousa et al. 2016), we manually examine the phase curve to ensure the validity of the reported periods for our data.If the star has no reported period, or if the reported period is not applicable to our data, we re-/determine the period by, firstly, identifying the maximum in the autocorrelation function and, secondly, refining the period within the full height width of the maximum (FHWM).The best period solution is determined by selecting the one that provides the optimal folded phase curve counterpart.This selection process involves manually examining the phase curves and comparing them to those from other potential period values.Figure 8 illustrates DL Tau and IQ Tau as examples, with estimated periods of 10.4 day and 12.9 day, respectively.Next, we use the Gaussian Process (GP) to produce a smooth light curve pattern with the optimal period, which is then subtracted from the original light curve to measure the residual noise.Finally, the resulting residual variance is divided by the variance of the original light curve to obtain the Q-value.
We then divide the light curves of our samples into eight different categories using the definition of Q and M ranges given by Cody et al. ( 2022

Variable but unclassifiable (U)
The values of Q and M for the CTTSs are listed in Table 4, along with the variability classes.The TESS light curves of 12 CTTSs in our sample, along with their assigned M, Q and variability classes, are shown in Figure 9. Also in Figure 10, we show the distribution of the variability statistics for our sample in the Q-M diagram.

Mass-accretion rate measurements from TESS light curves
In Section 3.3, we have identified four bursters in our samples: DQ Tau, DL Tau, GG Tau, and Haro 6-13.To estimate the mass-accretion rates from their TESS light curves, we need to first convert the excess flux yields associated with the bursts in the T ESS bandpass to the U-band excess luminosity which can be used as the proxy for the massaccretion rates (e.g., Tofflemire et al. 2017).
To do so, we searched for other observational data of the stars from ASAS-SN, Pan-STARRS, and ZTF (Bellm et al. 2019;Masci et al. 2019;IRSA 2022), that were taken simultaneously with TESS in different optical bands.We found such data only for DL Tau and Haro 6-13.DL Tau was observed by the ASAS-SN survey and the ZTF survey.A burst event of DL Tau at BT-2457000 = 2504 days was captured by TESS, ZTF g, and ASAS-SN g.Haro 6-13 was observed by ZTF during the same observation period as the TESS mission for it.A burst event in Haro -13's TESS light curve at BT-2457000 = 2516 days was also observed by the ZTF with both the g and r bands.Their light curves of the bursts are shown together in Figure 11.
Here, the periodic variation due to possible dark spot modulation has been removed from the TESS light curves of DL Tau and Haro 6-13 by subtracting the best-fit sinusoidal curves with corresponding periods (see Table 4) from the original curves.The fluxes of both data have good quality flags, and both are normalized to their respective median flux levels using the equation: ∆F where F represents the flux as a function of time and F is their median in the light curve.We express the observed amplitude due to the hot spot associated with the accretion burst, that covers an area A = XπR 2 * on the visible stellar photosphere, where X is the filling factor and R * is the stellar radius as: where B λ is the Planck function, R λ is the transmission function of the band, T * is the stellar effective temperature, and T hspot is the temperature of the hot spot.The integral term can be expressed as σT 4 f λ (T ) where σ is the Stefan-Boltzmann constant and f λ is the bandpass response factor as a function of temperature for the used band (λ).The response factor functions of the bands used in the study, i.e., sdss g, ZTF g, ZTF r, U , and T ESS bands, are shown in Figure 12.The amplitudes of DL Tau's event are 0.848, 0.77, and 0.438 in ASAS-SN's sdss g, ZTF g, and T ESS bands, respectively.We note that since there are two data points from ZTF g obtained during the event, the reported amplitude in ZTF g is computed as the mean value of these data points.Additionally, the central time between these two data points is very close to that of the ASAS-SN's observation.With this information, we can estimate the temperature of the hot spot by using the following equation: The hot spot's temperatures derived from T ESS-ASAS-SN g and T ESS-ZTF g are about 5700 K and 5500 K, respectively.As a result, the temperature of the hot spot is estimated to be the average of these values, which is approximately 5600 K.We substitute this value back to the Eq. 9 to calculate the filling factor of the spot and obtain X = 0.038 or 3.8%.
For the case of Haro 6-13, the burst amplitudes in TESS, ZTF r, and ZTF g bands are measured to be 0.292, 0.4, and 0.75, respectively.Because there is approximately a 0.1 day timing gap between the ZTF g and ZTF r data, the amplitude in the TESS band is the average of the corresponding amplitudes at ZTF g and ZTF r data times.By following the same estimation procedure as we conducted for DL Tau's burst, the temperature of Haro 6-13's burst is estimated to be 7100 K based on TESS and ZTF g bands data.The estimated temperature becomes a bit hotter when using TESS and ZTF r bands data, which is 7715 K. Thus we determined the burst temperature of Haro 6-13 to be about 7300 K by averaging these values.The filling factor of the burst-driven hot spot on Haro 6-13 is about X = 0.008 or 0.8%.The resulting values of the hot spot temperature and fractional area fall within the range of values observed and modeled in previous studies (e.g., Manara et al. 2013;Venuti et al. 2015), demonstrating the reliability of our method.However, we note that our estimations are lower limits as we haven't considered the limb-darkening in Eq. 9 and Eq. 10.In other words, we assume the hot spot is located at the center of the visible photospheric disk of the star.
We then extrapolate the burst amplitudes in the U -band with the estimated hot spot temperatures and the fractional areas.We find that the burst amplitudes in the U -band are about 3.49 times and 7.6 times higher than that in the T ESS band for DL Tau and Haro 6-13, respectively.We produce the U -band excess amplitude curves of DL Tau and Haro 6-13 by multiplying their normalized TESS light curves with the factors of counterparts.
The U-band apparent magnitudes of DL Tau and Haro 6-13 are 14.04 (Ducati 2002) and 19.22 (Audard et al. 2007), which can be converted to a flux unit of F U = 3.77 × 10 −11 erg s −1 cm −2 and F U = 3.19 × 10 −13 erg s −1 cm −2 , respectively.The U-band luminosities, L U = 1.15 × 10 32 erg s −1 and L U = 6.32 × 10 29 erg s −1 , of DL Tau and Haro 6-13 can be estimated using the formula L U = F U 4πd 2 with the distance information of d = 159 pc and d = 129 pc (Gaia Collaboration 2020).We produce the U-band excess luminosity (L U excess ) curves of these stars by multiplying their U-band excess amplitude curves with the L U values.Then the accretion luminosities (L acc ) can be calculated using an empirical correlation between L acc and L U excess given by Gullbring et al. (1998): where L is solar luminosity.Finally, the mass-accretion rates can be determined from the accretion luminosities by using Eq 6 with the stellar parameters in Table 2.The mass-accretion rates and luminosities of DL Tau and Haro 6-13 are displayed in Figure 13.

RESULTS AND DISCUSSION
4.1.Mass-accretion rates from the LAMOST Spectra We estimated the mass-accretion rates of low-mass classical T Tauri stars with spectral types from mid-M to late-K in the Taurus Molecular Cloud using LAMOST spectral data.The indicators we used were the Hα 6563 Å, Hβ 4861 Å, CaII 8498, 8542, 8662 Å, and HeI 4471, 5875, 6678, 7065 Å emissions.The equivalent widths of these lines and the corresponding mass-accretion rates are listed in Table 5 and Table 6, respectively.
All of our samples have the Ṁacc measurement estimated from Hα line.On average, these stars have an Hα massaccretion rate of 1.76 × 10 −9 M yr −1 .GG Tau has the highest measured Hα mass-accretion rate among the sample, with a value of 6.35 × 10 −9 M yr −1 .HP Tau has a significant Hα emission, but its mass-accretion rate is lower than those of the other stars with a value of 7.88 × 10 −11 M yr −1 in this study.On the other hand, from Hβ emissions, the average rate of these star is 2.89 × 10 −9 M yr −1 .UY Aur and GG Tau both have the highest Hβ-derived rate of 9.11 × 10 −9 M yr −1 .Consistent with the Hα result, HP Tau exhibits the lowest mass-accretion rate estimated from the Hβ line, which is 2.44 × 10 −11 M yr −1 .Only Haro 6-13 has no Hβ accretion rate because it is too faint for LAMOST in the blue wavelength part, making it impossible to estimate.Not all targets in our sample have measurable equivalent widths for the other emission lines, which means that not all corresponding mass-accretion rates could be estimated.Figure 14 shows the mass-accretion rate distribution of our targets by displaying the measurements derived from all available emissions adopted in this study.HL Tau has the largest dispersion of the spectral-derived mass-accretion rates with a coefficient of variation (CV) value (the ratio of the standard deviation to the mean) of 1.15, as compared to other stars in the sample.The CVs of other stars are smaller than one with a mean value of 0.41±0.22.After field checking, we have ruled out the possibility of contamination from the closest neighboring star of HL Tau, XZ Tau, as it is located about 23 away from HL Tau, and thus, can be resolved separately by the LAMOST fiber with a 3.3 diameter.The underlying mechanism responsible for this phenomenon remains unknown, making it an intriguing topic for future investigation.On the other hand, the lowest dispersion CV value is found in FN Tau, which is 0.08.

Mass-accretion rates of DL Tau and Haro 6-13 from the TESS light curves
The time-series mass accretions rates and luminosities derived from the TESS light curves for DL Tau and Haro 6-13 are shown in Figure 13.We find that the mass-accretion rate of DL Tau reached up to 2.5 × 10 −8 M yr −1 during the observation time.The mass accreted by DL Tau over approximately 50 days is ≈ 3.4 × 10 −10 M , with an average rate of ≈ 2.5 × 10 −9 M yr −1 .These results are comparable to the mean value of the mass-accretion derived from LAMOST for DL Tau, which is 4.6 × 10 −9 M yr −1 with a standard deviation of 2.5 × 10 −9 M yr −1 .Meanwhile, the maximum mass-accretion rate of Haro 6-13 from the TESS observations is 2 × 10 −10 M yr −1 .It accreted mass with a value of about 2.7 × 10 −12 M over approximately 50 days, with an average rate of ≈ 2 × 10 −11 M yr −1 .These results are overall lower than the measurement, which is about 1.6 ± 0.8 × 10 −9 M yr −1 on average, determined from LAMOST data.
An interesting question is whether the mass increase of T Tauri stars is dominated by long-term, steady mass accretion or by occasional bursts of accretion?First, we need to define a steady accretion rate, and this rate does not change over time.We then calculate the total mass accumulated by the steady rate of the star during the TESS observation.This is compared with the total amount of burst accretion to determine whether the CTTS is dominated by steady status or episodic effect.
Here, we take the mean rate obtained from LAMOST emission measurement as the steady rate.Thus, the steady rate of mass-accretion of DL Tau is 4.6±2.5×10−9 M yr −1 .In this case, DL Tau accumulated masses of 6.3±3.4×10−9 M through steady rate over a period of about 50 days.This result is larger than the total mass of ≈ 3.4 × 10 −10 M accreted by the stars' burst-like events observed by TESS.Similarly, the steady mass-accretion of Haro 6-13 is 1.6 ± 0.8 × 10 −9 M yr −1 .This value is about one order higher than the maximum of 2 × 10 −10 M yr −1 from the TESS observation.Therefore, in this scenario, the mass-accretion processes of DL Tau and Haro 6-13 are both steady accretion dominant.

Flare activity of our CTTSs sample
CTTSs are also known for their high magnetic activity, which is often manifested in flaring events.Flares with the sudden, short-lived increases in brightness are thought to be caused by magnetic reconnection events on the star's surface (Parker 1963;Maehara et al. 2012).The magnetic fields of the CTTS form a large loop connecting the star and circumstellar disk.When the magnetic field lines eventually reconnect, a large amount of energy is released in the form of a huge flare.This kind of the flare resulting from the star-disk magnetic interaction could last more than 50,000 seconds, and it is a characteristic feature of CTTS while affecting the mass-accretion process in some degree (Giardino et al. 2007;López-Santiago et al. 2016;Colombo et al. 2019).

Flare detection
To analyze the flares in our TESS light curves, we employ an algorithm developed by Lin et al. (2019) that generates a flare-free light curve.We then identify individual flare candidates from a residual curve generated by subtracting the flare-free light curve from the original curve.This method allows us to accurately quantify the properties of each flare, including its amplitude, duration, and energy release.We further examine the candidates for their authenticity by visually inspecting the Target Pixel File to generate the light curve of each pixel, using the Python package Lightkurve (Lightkurve Collaboration et al. 2018) and the method described by Wu et al. (2015).In our sample, we detect flares in the light curves of three stars: BP Tau (3 flares), DN Tau (6 flares), and FN Tau (4 flares), in total 13 flares, and some of them are multi-flare events.The light curve profiles of these flares are shown in Figure 15.However, the absence of flares in other stars does not necessarily mean that they are not flaring stars.For example, Tofflemire et al. (2017) and Kóspál et al. (2018) both detected several flare events in DQ Tau from optical and Kepler observations, whereas we did not detect any flares in DQ Tau from the TESS observation in this study.It is possible that there are microflares with amplitudes or duration that are too small to be detected by TESS.

Flare amplitude and energy
The amplitude of the flare profile as a function of time can be determined using the Eq. 8.The flare energy is estimated using an equivalent duration (ED) method defined as follows (Gershberg 1972): It represents the time it takes for a non-flaring star to emit the same amount of energy as that released during a given flare.Therefore, this quantity can be used to estimate the energy of the flare by multiplying it by the quiescent stellar flux in the telescope bandpass and the corresponding quiescent stellar flux of the star can be estimated from Here, R * is the stellar radius, T * is the effective temperature, and σ is the Stefan-Boltzmann constant.f λ is the same function as appeared firstly in Eq. 9, which is a bandpass response factor as a function of temperature for the wavelength (λ) band and in this case λ = T ESS, i.e., the fractional response function of TESS (see Figure 12).The amplitudes, duration, timing at peak in BTJD, and released energy of detected flares are listed in Table 7.These flares have energy levels ranging from 2 × 10 34 erg to 6 × 10 35 erg.The flare with the highest amplitude of 0.22 is observed in BP Tau, with a duration of about 76 minutes and energy of 1 × 10 35 erg.The most energetic flare in our sample is observed in FN Tau, with an energy of 6.3 × 10 35 erg, an amplitude of 0.09, and a duration of 438 minutes.Furthermore, given the duration of the flares, these events likely originate from the localized magnetic reconnections near the stellar surface.

Flare frequency distribution
The flare activity of a star can be characterized by its flare frequency distribution, which can be described by a linear relationship expressed in logarithmic form as where N is the cumulative frequency of flares with released energy E. The slope of the distribution, represented by the coefficient β, can provide information about the properties of the star's flare activity.If β <-1, then most of the total energy emitted by the flares comes from small flares, while for β >-1, stronger flares dominate the energy output.The flare frequency distribution power-law can also be expressed in differential form as where dN is the number of flares with energies in the range of E and E + dE, and α = 1 − β (Hawley et al. 2014).
We use the maximum likelihood estimator (MLE) for α as derived by Clauset et al. (2009): and the error for α can be estimated from Here, n is the number of flare samples, and E min is the minimum value of the energy of the observed flares.By setting E min = 2.2 × 10 34 erg and n = 13, we have calculated the value of α to be 1.68 ± 0.20.From this, we can derive the linear slope of logarithmic flare frequency distribution, which is represented by β = −0.68 ± 0.20 and a constant a of 24.32. Figure 16 shows the flare frequency distribution and its best-fit power-law slope for our CTTS sample.These results suggest that strong flares with high energy dominate the flare activity of these stars.The power-law distribution also indicates that our sample stars exhibit extreme flare activity, with the ability to generate flares with energy levels stronger than E = 1 × 10 34 erg about 12 times per year.This is in stark contrast to solar-type stars, which produce flares of this magnitude only once every 35-40 years (Shibayama et al. 2013).Even one of the most flare-active young M stars, Wolf 359, can produce only one such flare a year (Lin et al. 2021).

The changes of the variability classes
The changes in the variability behavior of the CTTSs can be attributed as the inner disk structure responsible for the stellar occultation, which can undergo significant changes in just a few years, transitioning from a stable and well-organized geometry with a consistent inner disk warp to a more disordered distribution of dust (e.g., McGinnis et al. 2015).Another explanation is the stability of the accretion regime.Initially, the accretion geometry may consist of a main accretion funnel in each hemisphere, with the base of each funnel corresponding to the stable inner disk warp.However, over time, this geometry may become unstable, resulting in the formation of random accretion funnels (e.g., Kurosawa & Romanova 2013).Cold spots due to the strong magnetic field of CTTSs might be able to cover a large fraction of stellar surface, which is similar to magnetically active low-mass stars (e.g., Yang et al. 2017;Lin et al. 2021), and the changes of their distribution and sizes may be also reflected in the changes of variability behaviors (e.g., Frasca et al. 2005).
Ten stars in our sample were previously observed by K2 Campagin 13 (AA Tau, DQ Tau, GO Tau, Haro 6-13 HK Tau, HL Tau, and IQ Tau), or TESS Section 19 (BP Tau, FN Tau, and UY Aur), and their light curves were analyzed and classified for variability by Cody et al. (2022) and Robinson et al. (2022).The variability classes given by Cody et al. (2022) and Robinson et al. (2022) are also listed in Table 4.By comparing the variability classes given by Cody et al. (2022) based on the K2 C13 data observed about four years ago, we find that most of the stars (6 out of 7) change their variability classes except DQ Tau.Specifically, AA Tau changes from "S" to "QPS", GO Tau changes from "U" to "APD", Haro 6-13 changes from "S" to "B", HK Tau changes from "QPS" to "QPD", IQ Tau changes from "QPS" to "P", and HL Tau changes from "QPS" to "QPD".In contrast, the variability classes of BP Tau and FN Tau do not change compared to those derived by Robinson et al. (2022) using TESS Section 19 light curves from about 600 days ago.Even though UY Aur changes from "B" to "QPS" in 600 days, we still see a burst-like event at BJTD = 2491 in its TESS Section 43 light curve.These findings may suggest that the long-term timescale of the variation associated with, for example, the stable/unstable accretion regime and inner disc warp due to the star-disk interaction of the CTTSs is in the range between about 1.6 years to 4 years or longer.However, it remains to be determined whether these changes occur in a cyclical or random manner.For instance, Haro 6-13 changing from the "B" type back to the "S" type, or to any other type, would require further investigation through continuous monitoring in the future.

Effect of the inclination on the spectra, accretion, and variability properties
The spectral properties of the CTTSs have been predicted to be correlated with the inclination of the system, especially with the broad or high-velocity component of the forbidden line profile (e.g., Hartigan et al. 1995).The highvelocity components are generally believed to originate in fast outflows that keep thrusting material into interstellar space in form of jets.Among the known forbidden lines, the ones with the most distinctive high-velocity component spectral features are the [NII] lines at 6548 Å and 6583 Å, which behave as blueshifted emission with a well-defined peak (Hirth et al. 1997).Nevertheless, the resolving power of the LAMOST observations is not high enough to distinguish the [NII] lines from Hα emission, thus we can not derive the outflow velocity from these lines using the data we have in this study.
Additional forbidden lines that can provide insight into the wind velocity observed by LAMOST are the [OI] line at 6300 Å, [SII] lines at 6716 Å and 6731 Å, and [CaII] lines at 7291 Å and 7323 Å.While the blueshifted emission of the [OI] line may not display a distinct peak when compared to the [NII], [SII], and [CaII] lines, the extended blue wing can be interpreted as an indicator of the projected wind velocity (Appenzeller & Bertout 2013).Out of the sample CTTSs, ten stars exhibit [OI] lines in their LAMOST spectra, while [SII] lines have only been detected in AA Tau and HL Tau.Moreover, [CaII] lines are observed in only one star, HL Tau.
We estimate the projected wind velocity by using the blue edge of the line profile with flux reaching 25% of the peak value when the line shows no clear blue peak.Because of the low resolution of the LAMOST data, for those profiles showing a blue peak, the velocity is directly estimated from the wavelength of the peak.The LAMOST [OI] line, [SII] lines, and [CaII] lines velocity profiles of the CTTSs are shown in Figure 17, Figure 18, and Figure 19, respectively.The results of [OI], [SII], and [CaII] velocity measurements are listed in Table 8.For stars where velocity measurements are available from not only [OI], but also [SII] and [CaII], i.e., AA Tau and HL Tau, the results are consistent across all lines.We note that there are two [OI] velocity results for a binary HK Tau.The diameter of LAMOST's fiber is about 3.3 arcsecond, and the angular separation between HK Tau A and HK Tau B is about 2.3 arcsecond.Thus, LAMOST observation can not resolve these two stars separately.We detect a blueshifted peak of [OI] line at a wind velocity of −56.49km/s, which may have originated from the edge-on system HK Tau B with an inclination of 84 o .The velocity derived from the blue edge of the line with a value of −182.88 km/s belongs to HK Tau A, which has a relatively face-on disk with an inclination of i = 56.9o .The comparison between the inclinations and [OI] line-derived projected wind velocity is shown in Figure 20.Our results show a positive correlation between the cosine of inclination and the projected wind velocity, consistent with the trend observed by Appenzeller & Bertout (2013).
The inclination angle of a CTTS may have an impact on the estimation of its mass-accretion rate to some degree.This is because the varying geometric angle of the disk could obscure areas emitting accretion-associated emission from the line of sight, leading to uncertainties in the measurement.However, our analysis indicates that there is no significant correlation between the inclination angle and the log spectra-derived mass-accretion rates (see Figure 21) with a Pearson's correlation coefficient of 0.34 for Hα's measurement, and the correlation coefficients are even lower for the results from other emissions.One possibility is that the emission regions driven by accretion are uniformly distributed on the surface of the CTTS, making the effect of disk obscuration negligible.Another possible explanation is that the properties of the inner disk may play a more significant role than those of the outer disk.Several studies have reported misalignments between the inner and outer disks.For example, the inner disk of AA Tau has an inclination angle of i = 75 deg, closer to edge-on (Bouvier et al. 2003;O'Sullivan et al. 2005), while the outer disk has an inclination angle of i = 59 deg (Loomis et al. 2017).Unfortunately, information on the inclination angle of the inner disk for our sample stars, except for AA Tau, is not available.Further studies using high-resolution double-peaked Keplerian line spectroscopic observations or spatially resolved observations using infrared interferometry on VLTI may be required to address this issue.
The type of photometric variability of CTTS is dominated by the accretion and dark spots located on the photosphere and the disk around, thus it could be geometric viewing angle dependent.Cody & Hillenbrand (2018) found a positive correlation between the M value they estimated and the inclination Barenfeld et al. (2017) gave, so that dippers and QPS sources tend to have high inclinations, while bursters and stochastic sources have low inclinations.However, due to the substantial uncertainties associated with the inclination measurements, they cautiously stated this observed tendency.A more recent study by Robinson et al. (2022), based on CTTSs with lower inclination uncertainties, found little evidence for a strong correlation between disk inclination and M value (see Fig. 12 in their paper).The upper panel shown in Figure 22 illustrates the relationship between inclination and M value of CTTS, including our samples and those from Robinson et al. (2022).Our analysis reveals a relatively strong and positive correlation between M value and inclination, with a Pearson's value of 0.61 ± 0.03, compared to the results given by Cody & Hillenbrand (2018) and Robinson et al. (2022).
On the other hand, we find an indistinct and negative correlation between the periodicity metric Q value and inclination with a Pearson's value of −0.54 ± 0.03 (see the bottom panel in Figure 22) for our targets.This is consistent with the conclusion given by Robinson et al. (2022) who also found the weak inverse correlation between disk inclination and periodicity metric of Taurus CTTSs.Furthermore, the distribution of variability classes in the inclination-Q space in our study is similar to that of Taurus members studied by Cody et al. (2022), as shown in Figure 7 of their paper.Nevertheless, it is possible that such distribution and the correlation between M value and inclination, as well as the correlation between periodicity metric Q value and inclination, could vary over time due to the changing variability classes of CTTSs.

SUMMARY
We present a comprehensive study of 16 low-mass classical T Tauri stars in TMC for their mass-accretion and luminosity variability properties by using LAMOST low-resolution spectra and TESS photometric light curves.
From the LAMOST spectra, the mass-accretion rates of 16 CTTSs have been estimated using Hα, Hβ, CaII, and HeI emissions as indicators.The average mass-accretion rate of all our samples based on Hα lines is found to be 1.76 × 10 −9 M yr −1 .Not all targets in the sample have measurable equivalent widths for the other emission lines, which means that not all corresponding mass-accretion rates could be estimated.HL Tau has the largest dispersion of the spectral-derived mass-accretion rates with a coefficient of variation (CV) value of 1.15, which is the ratio of the standard deviation to the mean.The reason that causes this huge dispersion remains unclear.On the other hand, the lowest dispersion CV value is found in FN Tau, which is 0.08.
We also compute the time-series mass-accretion rates based on the TESS light curve, ZTF, and ASAS-SN survey for the stars DL Tau and Haro 6-13, which are identified to be the bursters according to their brightness variability asymmetry metric M values of −0.63 and −0.31, respectively (Cody et al. 2014).In total over about 50 days, the average accretion rates of DL Tau and Haro 6-13 are 2.5×10 −9 M yr −1 and 1.4×10 −11 M yr −1 , respectively.The rate of DL Tau is consistent with the mean value of the mass-accretion derived from LAMOST.However, for Haro 6-13, this rate is generally lower compared to the LAMOST estimate We further investigate whether the mass-accretion mechanisms of these two stars are dominated by steady or burst accretion.By assuming the average rates obtained from LAMOST emission measurement as the steady rates, the mass accretions in DL Tau and Haro 6-13 are dominated by the steady process.
We have detected 13 flares in total in BP Tau, DN Tau, and FN Tau.The amplitude, duration, timing at peak in BTJD, and released energy of detected flares are determined, with the most energetic flare observed in FN Tau with the energy of 6.3 × 10 35 erg.The power-law flare frequency distribution with a slope of β = −0.68 ± 0.20 indicates that the flare activity of the CTTSs is dominated by strong flare with high energy.Our analysis of the flare frequency distribution also reveals that the stars in our sample display extremely high levels of flare activity with the capability of producing flares with energy greater than E = 1 × 10 34 erg up to 12 times per year, which is a significant deviation from solar-type stars and young M dwarfs.
Ten stars in our sample were previously observed and classified for variability by Cody et al. (2022) and Robinson et al. (2022).Comparing the variability classes from Cody et al. (2022) based on K2 Campaign 13 data from about four years ago, we found that most of the stars (6 out of 7) showed changes in their variability, except for DQ Tau.In comparison with the results from Robinson et al. (2022), BP Tau and FN Tau's variability classes remained the same, while UY Aur changed from "B" to "QPS" in 600 days but still exhibited a burst-like event in its TESS Section 43 light curve.These findings suggest that the long-term timescale of variation in CTTSs could be between 1.6 to 4 years or longer and require further investigation to determine if changes occur cyclically or randomly.
The CTTSs' projected wind velocities in terms of the blue edge/peaks of [OI], [SII], and [CaII] forbidden emissions detected in LAMOST are found to be correlated with the inclination, that is stars with low inclination tend to have fast projected wind velocity.On the other hand, no tendency has been found between inclination and mass-accretion rates derived from the spectral lines' luminosities.One possibility is that the regions of emission driven by accretion are distributed uniformly across the stellar surface, so the disk obscuration doesn't have much of an impact.Alternatively, it's possible that the properties of the inner disk are more important than those of the outer disk.There is evidence to suggest that the inner and outer disks of CTTS can become misaligned due to, e.g., magnetic warping.Unfortunately, in our sample, only AA Tau has reported inner disk inclination angle.High-resolution double-peaked Keplerian line spectroscopic observations or spatially resolved observations using infrared interferometry on VLTI may be required for further studies of this topic.
The brightness variabilities of CTTSs have been also considered to be influenced by the inclinations to some degree.Our analysis reveals a relatively significant and positive correlation between M value and inclination, with a Pearson's value of 0.61 ± 0.03, compared to the results given by Cody & Hillenbrand (2018) and Robinson et al. (2022).We observed a weak negative correlation between the periodicity metric Q value and inclination for our targets, with a Pearson's value of −0.54 ± 0.03.This finding supports the conclusion reached by Robinson et al. (2022), who also reported a weak inverse correlation between disk inclination and the periodicity metric of Taurus CTTSs.Additionally, we found that the distribution of variability classes in the inclination-Q space is similar to that of Taurus members studied by Cody et al. (2022).However, since the variability classes of CTTSs can change over time, it is possible that the distribution and correlations between the M value and inclination, and the periodicity metric Q value and inclination, could vary as well.Thus, further research is required to fully understand the long-term evolution of these correlations.(Loomis et al. 2017); GG Tau (Guilloteau et al. 1999); HK Tau A (Manara et al. 2019) and HK Tau B (Villenave et al. 2020); UY Aur (Close et al. 1998).Note that the stellar parameters, mass and radius, of HK Tau shown here belongs to HK Tau A.

Figure 2 .
Figure 2. The example for demonstrating the examination of the contamination from nearby sources.The left panel is the target pixel file (TPF) of GG Tau (TIC 245830955), and GG Tau is marked by a filled red circle.The hollow red circles are the field stars with GRP < 15 queried from Gaia DR3.The pink mask pixels are the aperture given by the SPOC pipeline.We found that the field star located in the aperture is much fainter than GG Tau, and the overall shapes of light curves from different test apertures do not change significantly compared to the SPOC light curve.Therefore, we verified that the data of GG Tau can be used in this study.The right panel is the TPF of GI Tau (TIC 268444803).There are three neighbor stars within or at the edge of the SPOC aperture, and one of them is equally bright to GI Tau, making it hard to distinguish the variability feature between GI Tau and the neighbor stars in the light curve.As the result, we excluded GI Tau data from the following analysis.

Figure 3 .
Figure 3. Three example of the flux calibration processes in this study.The upper left is AA Tau, the upper right is DL Tau, and the lower middle is GO Tau.In each panel, the blue curve represents the star's spectrum, while the gray curve represents its corresponding best-fit template spectrum given byManara et al. (2013).The photometric data from ASAS-SN g band and V band are marked as a yellow circle in AA Tau case and a green circle in DL Tau case, respectively.The Pan-STARRS i band data of GO Tau is marked as a red circle in GO Tau panel.

Figure 4 .
Figure 4. Hα emission profiles of our targets observed by LAMOST.In each panel, the spectra is displayed in red, and the black line represents the fit local continuum.The dash vertical line in gray marks the center wavelength of the Hα line.

Figure 5 .
Figure 5. Hβ emission profiles of our targets observed by LAMOST.Dashed vertical gray line represents the center of the line in each panel.

Figure 6 .
Figure 6.CaII 8498, 8542, 8662 Å emission profiles of our targets observed by LAMOST.Dashed vertical gray lines represent the center of the CaII triplet lines in each panel.

Figure 7 .
Figure 7. HeI 4471 (a), 5875 (b), 6678 (c), 7065 (d) emission profiles of our targets observed by LAMOST.Dashed vertical gray line represent the center of the HeI in each panel.

Figure 8 .
Figure 8.The examples of the period re-/estimate of the CTTSs, here shows (a) DL Tau and (b) IQ Tau.In each case, the upper left panel shows the autocorrelation function.The gray one is the original ACF, the blue one is a smoothed version of the ACF.The red cross marks the peak of the ACF, and the FHWM highlighted by green, from where the best period is determined.The upper right panel is the phase curve (black) and the best Gaussian Process fit curve (red).The panel at the bottom is the full TESS Section 43 and 44 light curve of the star (black) with the best Gaussian Process fit curve (red).The best periods we found for DL Tau and IQ Tau are 10.4 day and 12.9 day, respectively.

Figure 9 .
Figure 9.The TESS light curves of 12 CTTSs in our sample.Their variability classes are also displayed.

Figure 10 .
Figure 10.The Q-M statistics diagram for our sample of CTTSs.The stars surrounded by black hollow circles in this diagram are binary or multi-star systems.Different colors represent different variability classes: yellow (B), purple (P), green (QPS), blue (APD), and red (QPD).

Figure 11 .
Figure 11.Left panel: A burst event of DL Tau at BT-2457000 = 2504 was captured by TESS, ASAS-SN g, and ZTF g.Right panel: A burst event of Haro 6-13 at BT-2457000 = 2516 was observed by TESS, ZTF g, and ZTF r.In each panel, the black curve represents the burst event observed by TESS, and the green filled circles are g band observations from the ZTF survey.The hollow green diamond represents the data from the ASAS-SN survey.In the case of Haro 6-13, the red filled circle represents the observation from the ZTF r band.

Figure 12 .
Figure 12.The fractional response functions of T ESS (red), ZTF r (orange), ZTF g (lightgreen), sdss g (green), and U bands (blue).They are the function of the blackbody temperature.

Figure 13 .
Figure 13.The time series mass-accretion rates of DL Tau and Haro 6-13 are shown in the upper panel and bottom panel, respectively.The left axis shows the mass-accretion rates derived from the TESS, ASAS-SN, and ZTF observations, while the right axis displays the corresponding accretion luminosities.

Figure 14 .
Figure 14.The distribution of mass-accretion rates derived from all available emissions in the spectra of our CTTSs.The names of the stars marked with an asterisk (*) in the upper right indicate that they are binary or multi-star systems.

Figure 16 .
Figure 16.The flare frequency distribution and the best-fits power-law function with a slope of -0.68±0.20 of our CTTS sample based on 13 detected flares in this study.

Figure 17 .
Figure17.The velocity profiles of [OI] 6300 lines of ten CTTSs in our sample.In each panel, the gray vertical dashed line indicates the center of the line, and blue dashed line marks the wind velocity of the blueshifted wind or peak of the line.In a special case of HK Tau, we also show a green dashed line to represent the blueshifted wind velocity, which is likely originating from a companion of the binary system, HK Tau A, and the blueshifted peak velocity marked with a blue dashed line belongs to the other companion of the system, HK Tau B.

Figure 18 .
Figure 18.The velocity profiles of [SII] lines observed in the spectra of AA Tau and HL Tau.In each panel, the gray vertical dashed line indicates the center of the line, and blue dashed line marks the wind velocity of the blueshifted wind or peak of the line.

Figure 19 .
Figure 19.The velocity profiles of [CaII] lines of HL Tau.The gray vertical dashed line indicates the center of the line, and blue dashed line marks the wind velocity of the blueshifted peak of the line.

Figure 20 .
Figure 20.The comparison between CTTSs' OI 6300 project wind velocities (vertical axis) and inclinations in terms of cosine angle (horizontal axis).The purple circles represent the data from Appenzeller & Bertout (2013), and the red circles are the data in this study.

Figure 21 .
Figure21.The relationship between the inclination and mass-accretion rates derived from the Hα emission luminosities.The Pearson's correlation coefficient is about 0.34, indicating the correlation is not significant.

Figure 22 .
Figure 22.The relationship between the flux asymmetry metric M and the disk inclination of our CTTSs is shown in the upper panel, and the relationship between the periodicity metric Q and the disk inclination is shown in the bottom panel.Our samples are represented by rhombus symbols, with different colors indicating different variability classes.We also include the results obtained by Robinson et al. (2022) for M -vs-inclination and by Cody et al. (2022) for Q-vs-inclination as the comparisons, shown as yellow circles.The stars surrounded by black hollow circles in this diagram are binary or multi-star systems.Note that the large uncertainty of HK Tau's inclination is because it is an unresolved binary system in LAMOST observation, which is composed of two companions, HK Tau A with i = 56.9o and HK Tau B with i = 84 o .The blue dashed lines in both panel indicate the thresholds of the bursting (M < −0.25), dipping (M > 0.25), periodic (Q < 0.15), and aperiodic (Q > 0.85).

Table 2 .
Stellar properties of 16 CTTSs in this study

Table 3 .
Photometric data of our samples taken by ASAS-SN and Pan-STARRS surveys within 20 days of LAMOST observation.Name LAMOST OBT ASAS-SN OBT ASAS-SN g ASAS-SN V Pan-Starrs OBT Pan-Starrs g Pan-Starrs i

Table 4 .
The photometric variability classes of the CTTSs observed by TESS in this study and the literature (if any).
Variability C Variability R Period * Period P Period Re * represents the measurements given in this study.C : Cody et al. (2022).R : Robinson et al. (2022).P : Percy et al. (2010).Re : Rebull et al. (2020).The unit of the period is in days.

Table 5 .
Equivalent widths of all mass-accretion indicator emissions of 16 CTTSs in LAMOST spectra in this study.
Note.Note that we use positive equivalent width to indicate emission.The unit is Å

Table 6 .
Mass-accretion rates of 16 CTTSs derived from LAMOST spectra in this study.

Table 7 .
13 flares detected in the TESS light curves of BP Tau, DN Tau, and FN Tau.
Note: The columns are: (1) flare peak amplitude, (2) flare energy in the TESS band, (3) start time of flare, (4) time of flare peak, (5) end time of flare, and (6) duration of flare in the unit of days.

Table 8 .
The projected wind velocities of CTTSs observed by LAMOST.