# New insights on the near-infrared veiling of young stars using CFHT/SPIRou data

A. Sousa, J. Bouvier, S. Alencar, J.-F. Donati, C. Dougados, E. Alecian, A. Carmona, L. Rebull, N. Cook, E. Artigau, et al.

► To cite this version:

A. Sousa, J. Bouvier, S. Alencar, J.-F. Donati, C. Dougados, et al. New insights on the near-infrared veiling of young stars using CFHT/SPIRou data. Astronomy & Astrophysics - A&A, 2023, 670, pp.A142. 10.1051/0004-6361/20224720. hal-04252220

HAL Id: hal-04252220

https://hal.science/hal-04252220v1

Submitted on 26 Oct 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# New insights on the near-infrared veiling of young stars using CFHT/SPIRou data*

A. P. Sousa<fnref n="1" />, J. Bouvier<fnref n="1" />, S. H. P. Alencar<fnref n="2" />, J.-F. Donati<fnref n="3" />, C. Dougados<fnref n="1" />, E. Alecian<fnref n="1" />, A. Carmona<fnref n="1" />, L. Rebull<fnref n="4" />, N. Cook<fnref n="5" />, E. Artigau<fnref n="5" />, P. Fouqué<fnref n="3" />, R. Doyon<fnref n="5" />, and the SLS consortium

Received 5 August 2022 / Accepted 21 December 2022

## ABSTRACT

**Context.** Veiling is ubiquitous at different wavelength ranges in classical T Tauri stars. However, the origin of the veiling in the infrared (IR) domain is not well understood at present. The accretion spot alone is not enough to explain the shallow photospheric IR lines in accreting systems, suggesting that another source is contributing to the veiling in the near-infrared (NIR). The inner disk is often quoted as the additional emitting source meant to explain the IR veiling.

**Aims.** In this work, we aim to measure and discuss the NIR veiling to understand its origins and variability timescale.

**Methods.** We used a sample of 14 accreting stars observed with the CFHT/SPIRou spectrograph, within the framework of the SPIRou Legacy Survey, to measure the NIR veiling along the *YJHK* bands. We compared the veiling measurements with accretion and inner disk diagnostics. We also analyzed circumstellar emission lines and photometric observations from the literature.

**Results.** The measured veiling grows from the *Y* to the *K* band for most of the targets in our sample. The IR veiling agrees with NIR emission excess obtained using photometric data. However, we also find a linear correlation between the veiling and the accretion properties of the system, showing that accretion contributes to the inner disk heating and, consequently, to the inner disk emission excess. We also show a connection between the NIR veiling and the system's inclination with respect to our line of sight. This is probably due to the reduction of the visible part of the inner disk edge, where the NIR emission excess is expected to arise, as the inclination of the system increases. Our search for periods on the veiling variability showed that the IR veiling is not clearly periodic in the typical timescale of stellar rotation – which, again, is broadly consistent with the idea that the veiling comes from the inner disk region. The NIR veiling appears variable on a timescale of a day, showing the night-by-night dynamics of the optical veiling variability. In the long term, the mean NIR veiling seems to be stable for most of the targets on timescales of a month to a few years. However, during occasional episodes of high accretion in classical T Tauri stars, which affect the system's dynamic, the veiling also seems to be much more prominent at such times, as we found in the case of the target RU Lup.

**Conclusions.** We provide further evidence that for most targets in our sample, the veiling that mainly occurs in the *JHK* bands arises from dust in the inner disk.

**Key words.** stars: pre-main sequence – accretion, accretion disks – stars: variables: T Tauri, Herbig Ae/Be

## 1. Introduction

The photospheric lines of young low-mass accreting systems, commonly referred to as Classical T Tauri stars (CTTS), are shallower and present smaller equivalent widths than those of non-accreting stars with a similar spectral type. This phenomenon is known as the veiling of the photospheric lines (e.g., Hartigan et al. 1991; Valenti et al. 1993; Folha & Emerson 1999; Fischer et al. 2011). The presence of veiling suggests an additional emitting source, beyond the stellar photosphere contributing to the spectra of the targets, that is responsible for filling in the photospheric

lines (e.g., Hartmann & Kenyon 1990; Gullbring et al. 1998; Calvet & Gullbring 1998; Johns-Krull & Valenti 2001).

The veiling in the optical and in the IR wavelengths has been studied using different approaches in the aim to understand its origins and variability (e.g., Basri & Batalha 1990; Edwards et al. 2006; Fischer et al. 2011; Antoniucci et al. 2017; Gully-Santiago et al. 2017; Ingleby et al. 2013; McClure et al. 2013; Kidder et al. 2021). The veiling variability along the stellar spectra depends on wavelength (e.g., Fischer et al. 2011; Faesi et al. 2012; McClure et al. 2013; Rei et al. 2018), and optical veiling is often associated with the accretion process, while the accretion shock is thought to be at the origin of an additional continuum emitting source to the stellar spectrum. (e.g., Gullbring et al. 1998). Usually, the accretion spot presents an emission contribution maximum around the ultraviolet domain, which decreases as the wavelength increases (e.g., Calvet & Gullbring 1998). Nevertheless, we do not expect a significant contribution of the accretion spot continuum emission in the IR wavelength; therefore, the accretion spot alone cannot explain the veiling in the IR domain.

<footnote n="1">¹ Univ. Grenoble Alpes, CNRS, IPAG, 38000 Grenoble, France
e-mail: alana.sousa@univ-grenoble-alpes.fr</footnote>

<footnote n="2">² Departamento de Física-Icex-UFMG Antônio Carlos, 6627, 31270-901 Belo Horizonte, MG, Brazil</footnote>

<footnote n="3">³ Univ. de Toulouse, CNRS, IRAP, 14 Avenue Belin, 31400 Toulouse, France</footnote>

<footnote n="4">⁴ Infrared Science Archive (IRSA), IPAC, California Institute of Technology, 1200 E. California Blvd., Pasadena, CA 91125, USA</footnote>

<footnote n="5">⁵ Université de Montréal, Département de Physique, IREX, Montréal, QC H3C 3J7, Canada</footnote>

<footnote>* Based on observations obtained at the Canada-France-Hawaii Telescope (CFHT) which is operated by the National Research Council (NRC) of Canada, the Institut National des Sciences de l'Univers of the Centre National de la Recherche Scientifique (CNRS) of France, and the University of Hawaii. Based on observations obtained with SPIRou, an international project led by Institut de Recherche en Astrophysique et Planétologie, Toulouse, France.</footnote>
In the IR region, the veiling increases with wavelength and in some cases, it becomes greater than the veiling in the optical domain (e.g., McClure et al. 2013). The central star illuminates the inner disk and this region absorbs photons from the star, the accretion spot, and even the accretion funnel, then re-emitting them in the infrared (IR) as the system rotates (e.g., Chiang & Goldreich 1997, 1999). Therefore, the inner disk is suggested as the origin of the additional continuum emission that is essential to explaining the near-infrared (NIR) veiling, although the measured veiling is often too great to be explained as coming merely from the disk emission, based on model predictions (e.g., Folha & Emerson 1999; Johns-Krull & Valenti 2001). Many authors have used different techniques to connect the observed veiling with inner disk emission, such as measuring the temperature of the region where the veiling comes from and using a black body fit to the veiling. For most of these systems, they found temperatures compatible with the dust temperature in the inner disk (e.g., Fischer et al. 2011; Antoniucci et al. 2017; Alcalá et al. 2021). For a few targets, the blackbody temperature measured using veiling is too high for dust to survive in the inner disk; this would indicate that the veiling should arise from the gas in the inner disk inside the star-disk co-rotation radius (e.g., Antoniucci et al. 2017; Alcalá et al. 2021). However, McClure et al. (2013) found no evidence of hot gas inside the inner disk, which would be responsible for the NIR veiling. Instead, they explained the IR veiling as the combined emission from the accretion shock on the stellar surface and dust around the sublimation rim of the inner disk.

The veiling around 1 µm and the veiling in the K band and beyond can also have different origins. While the accretion spot emission contribution is not very substantial around 1 µm, we do not expect a significant contribution from the inner disk either. The origin of the veiling in this spectral domain is poorly understood. However, significant veiling was measured around 1 µm, primarily for high accretion rate systems (e.g., Edwards et al. 2006; Fischer et al. 2011; Ingleby et al. 2013). In the literature, there are only a few plausible explanations given for that case of veiling, such as a contribution from emission lines, filling in of the photospheric lines, or origination from the accretion shock (e.g., Sicilia-Aguilar et al. 2015; Dodin & Lamzin 2013). Even if we do not detect these extra emission lines directly, they can contribute to making the photospheric lines of the stellar spectra shallower, which increases the veiling measurement.

We cannot exclude other possible explanations to the IR veiling, such as an envelope around a star that can also be a source of additional emission to the photospheric continuum, with emission compatible with NIR veiling (e.g., Calvet et al. 1997). However, CTTSs are usually Class II stars and we would not expect a significant contribution from a dusty envelope. Furthermore, the veiling in the K band is higher than the dust envelope emission can explain (Folha & Emerson 1999).

In this work, we study the veiling and its variability in the NIR, using a sample of young accreting stars observed with the Canada-France-Hawaii Telescope SPectropolarimètre InfraROUge (CFHT/SPIRou). Our sample comprises stars with different properties, such as the mass accretion rate, spectral type, and inclination with respect to our line of sight. We computed the veiling for the *YJHK* bands and compared the results with accretion and inner disk diagnostics.

We organized the paper as follows. In Sect. 2, we present the sample of stars that we used in this work, and we describe the data used. In Sect. 3 we show the procedures to measure the veiling. We describe the results obtained from the veiling measurements in Sect. 4. In Sect. 5, we discuss the possible origin of

the veiling, and we compare our results with those of previous works. In Sect. 6, we present our conclusions.

## 2. Observations and targets selection

The sample of stars used in this work is composed of well-known young stars that are part of the SPIRou Legacy Survey-SLS science program: "Magnetic PMS star/planet survey" of some 50 class I, II, and III stars. The SPectropolarimètre InfraROUge (SPIRou) is a high-resolution velocimeter and polarimeter spectrograph (R ~ 75 000) covering the NIR wavelength range ~0.98–2.35 µm, corresponding to the spectral domain of the *YJHK* bands (Donati et al. 2018). The main science goals of CFHT/SPIRou Legacy Survey are the search for and characterization of planets around low-mass stars, and investigating the impact of the stellar magnetic field on planet and star formation in young systems (Donati et al. 2020a).

We aim to investigate the veiling of accreting young stars. Therefore, we selected, among the sample of stars observed by the SLS, 13 stars reported as accreting systems in the literature and for which we had a reasonable number of observations in time in comparison to the stellar rotation period. Most of these targets are classified as Class II and CTTS systems, and only V347 Aur is a Class I target. In addition, we added the T Tauri star J1604 (RX J1604.3–2130A) to the sample, which is not part of the SLS program, however, its CFHT/SPIRou observations are available.

In Table 1, we show the list of young stars analyzed in this work and the number of observations that we have for each target and each observational period. We also list three non-accreting T Tauri stars that cover our sample's spectral types and are also slow rotators ($v \sin i < 15$ km s⁻¹). We used these stars as templates to compute veiling and the residual profiles, as described in the following sections.

Each CFHT/SPIRou observation consists of four sub-exposures to measure Stokes V, taken at different orientations of the polarimeter and used to compute the non-polarized and the circularly polarized profiles. As the focus of this work is to analyze only the non-polarized component of the spectrum, we averaged the four sub-exposures to increase the signal-to-noise ratio (S/N) of the spectra obtained at each night. For the non-accreting systems, we averaged all the observations obtaining a mean spectra that were used as templates. The CFHT/SPIRou data were reduced, while the telluric-corrected spectra were obtained using the data reduction system APERO, versions 0.6.131, 0.6.132, and 0.7.232 (Cook et al. 2022). The spectra were corrected for the barycentric velocity and locally normalized to the continuum level, using a polynomial function to fit the continuum.

## 3. Procedures to measure the veiling

We computed the IR veiling of the targets following the method described by Hartigan et al. (1989), where we compare the spectra of the target with the spectrum of a non-accreting T Tauri star of a similar spectral type. The Zeeman broadening of the photospheric lines can affect the veiling measurements. Therefore, we used WTTSs as templates that should present a similar magnetic activity as the CTTS (e.g., Johns-Krull et al. 2000). The WTTSs also present comparable physical properties to the CTTSs, such as chromospheric activity and surface gravity, which make WTTSs stars suitable to measure the veiling in accreting systems. We list the templates applied to each star in Table 2.
<table><caption>**Table 1.** Sample of stars and number of observations per observational period.</caption><thead><tr><th>Star</th><th>SpT</th><th>v sin i<br>(km s<sup>-1</sup>)</th><th>Av<br>(mag)</th><th>i<sup>(*)</sup>,<sup>(a)</sup><br>(°)</th><th>2018</th><th>2019a</th><th>2019b</th><th>2020a</th><th>2020b</th><th>2021a</th><th>2021b</th><th>2022a</th><th>References<sup>(b)</sup></th></tr></thead><tbody><tr><td colspan="14" style="text-align:center; font-weight:bold;">Accreting systems</td></tr><tr><td>CI Tau</td><td>K4</td><td>9.5±0.5</td><td>0.65</td><td>55<sup>+35</sup><sub>-10</sub></td><td>2</td><td>6</td><td>26</td><td>5</td><td>39</td><td>—</td><td>—</td><td>—</td><td>(1), (2), (2), (2)</td></tr><tr><td>DoAr 44</td><td>K2-K3</td><td>17.0±1.1</td><td>2.0±0.2</td><td>30±5</td><td>—</td><td>8</td><td>—</td><td>—</td><td>—</td><td>—</td><td>—</td><td>—</td><td>(3), (3), (3), (3)</td></tr><tr><td>GQ Lup</td><td>K7</td><td>5±1</td><td>0.7</td><td>~30</td><td>—</td><td>8</td><td>—</td><td>18</td><td>6</td><td>—</td><td>—</td><td>—</td><td>(4), (5), (4), (5)</td></tr><tr><td>TW Hya</td><td>K6</td><td>6.0±1.2</td><td>0.0</td><td>18±10</td><td>—</td><td>12</td><td>—</td><td>14</td><td>—</td><td>27</td><td>—</td><td>—</td><td>(1), (1), (6), (7)</td></tr><tr><td>V2129 Oph</td><td>K5</td><td>14.5±0.3</td><td>0.6</td><td>60</td><td>9</td><td>—</td><td>—</td><td>17</td><td>8</td><td>—</td><td>—</td><td>—</td><td>(1), (1), (8), (9)</td></tr><tr><td>BP Tau</td><td>K5</td><td>9.0±0.5</td><td>0.45</td><td>~45</td><td>—</td><td>—</td><td>21</td><td>—</td><td>—</td><td>—</td><td>34</td><td>—</td><td>(10), (11), (6), (11)</td></tr><tr><td>V347 Aur</td><td>M2-M3</td><td>11.7<sup>+0.16</sup><sub>-0.24</sub></td><td>3.4</td><td>40</td><td>—</td><td>—</td><td>18</td><td>—</td><td>13</td><td>12</td><td>22</td><td>—</td><td>(12), (13), (12), (14)</td></tr><tr><td>DG Tau</td><td>K6</td><td>24.7±0.7</td><td>1.60±0.15</td><td>38±2</td><td>—</td><td>—</td><td>—</td><td>2</td><td>29</td><td>—</td><td>—</td><td>—</td><td>(10), (10), (6), (15)</td></tr><tr><td>RU Lup</td><td>K7</td><td>8.5±4.8</td><td>0.0</td><td>24</td><td>—</td><td>—</td><td>—</td><td>9</td><td>—</td><td>17</td><td>13</td><td>1</td><td>(16), (17), (14), (18)</td></tr><tr><td>V2247 Oph</td><td>M0</td><td>20.5±0.5</td><td>0.98±0.02</td><td>45±10</td><td>—</td><td>—</td><td>—</td><td>9</td><td>7</td><td>—</td><td>—</td><td>—</td><td>(19), (20), (19), (20)</td></tr><tr><td>DO Tau</td><td>M0</td><td>14.3±0.5</td><td>0.75</td><td>37.0±3.7</td><td>—</td><td>—</td><td>—</td><td>3</td><td>8</td><td>—</td><td>—</td><td>—</td><td>(6), (21), (6), (15)</td></tr><tr><td>J1604<sup>(*)</sup></td><td>K3</td><td>17.3±0.4</td><td>1.0</td><td>>61</td><td>—</td><td>—</td><td>—</td><td>12</td><td>—</td><td>—</td><td>—</td><td>—</td><td>(22), (22), (24), (23)</td></tr><tr><td>PDS 70</td><td>K7</td><td>16.0±0.5</td><td>0.01±0.07</td><td>50±8</td><td>—</td><td>—</td><td>—</td><td>4</td><td>—</td><td>—</td><td>6</td><td>—</td><td>(19), (25), (19), (25)</td></tr><tr><td>GM Aur</td><td>K4-K5</td><td>14.9±0.3</td><td>0.3±0.3</td><td>≥63</td><td>—</td><td>—</td><td>—</td><td>—</td><td>2</td><td>—</td><td>34</td><td>—</td><td>(26), (26), (26), (26)</td></tr><tr><td colspan="14" style="text-align:center; font-weight:bold;">Non-accreting systems</td></tr><tr><td>V819 Tau</td><td>K4</td><td>9.5</td><td>—</td><td>—</td><td>—</td><td>—</td><td>1</td><td>—</td><td>—</td><td>—</td><td>—</td><td>—</td><td>(27), (27)</td></tr><tr><td>TWA 25</td><td>M0.5</td><td>12.9±1.2</td><td>—</td><td>—</td><td>—</td><td>—</td><td>25</td><td>14</td><td>—</td><td>—</td><td>—</td><td>—</td><td>(6), (1)</td></tr><tr><td>TWA 9A</td><td>K6</td><td>7±3</td><td>—</td><td>—</td><td>—</td><td>—</td><td>1</td><td>—</td><td>—</td><td>—</td><td>—</td><td>—</td><td>(6), (1)</td></tr></tbody></table>

**Notes.** <sup>(*)</sup>RX J1604.3–2130A, elsewhere in the text, we used J1604 for brevity. <sup>(a)</sup>The inclination (see references) between the stellar rotation axis and our line of sight obtained from v sin i, period, and radius of each target, except for the inclinations of DO Tau and DG Tau that were obtained from the outer disk parameters. <sup>(b)</sup>References for the SpT, v sin i, Av and inclination (i), respectively.

**References.** (1) Torres et al. (2006); (2) Donati et al. (2020b); (3) Bouvier et al. (2020); (4) Alcalá et al. (2017); (5) Donati et al. (2012); (6) Herczeg & Hillenbrand (2014); (7) Alencar & Batalha (2002); (8) Donati et al. (2007); (9) Alencar et al. (2012); (10) Nguyen et al. (2012); (11) Donati et al. (2008); (12) Connelley & Greene (2010); (13) Flores et al. (2019); (14) Alecian et al. (in prep.); (15) Simon et al. (2016); (16) Alcalá et al. (2014); (17) Frasca et al. (2017); (18) Stempels et al. (2007); (19) Pecaut & Mamajek (2016); (20) Donati et al. (2010); (21) Kounkel et al. (2019); (22) Dahm et al. (2012); (23) Davies (2019); (24) Preibisch & Zinnecker (1999); (25) Thanathibodee et al. (2020); (26) Bouvier et al. (in prep.); (27) Donati et al. (2015).

<table><caption>**Table 2.** Parameters of the targets derived in this work.</caption><thead><tr><th>Star</th><th>Std</th><th>v<sub>rad</sub><br>(km s<sup>-1</sup>)</th><th>r<sub>Y</sub></th><th>r<sub>J</sub></th><th>r<sub>H</sub></th><th>r<sub>K</sub></th><th>EW<sub>Hel</sub><sup>(a)</sup><br>(Å)</th><th>EW<sub>Paβ</sub><sup>(a)</sup><br>(Å)</th><th>EW<sub>Bry</sub><sup>(a)</sup><br>(Å)</th><th>M<sub>Paβ</sub><br>×10<sup>-8</sup><br>(M<sub>☉</sub> yr<sup>-1</sup>)</th><th>M<sub>Bry</sub><br>×10<sup>-8</sup><br>(M<sub>☉</sub> yr<sup>-1</sup>)</th><th>α<sub>2_8</sub><sup>(b)</sup></th></tr></thead><tbody><tr><td>CI Tau</td><td>V819TAU</td><td>16.3±0.5</td><td>0.53±0.28</td><td>0.49±0.17</td><td>0.39±0.09</td><td>0.92±0.23</td><td>9.5±3.2</td><td>14.0±2.6</td><td>7.9±1.8</td><td>2.2</td><td>4.2</td><td>-0.84</td></tr><tr><td>DoAr 44</td><td>V819TAU</td><td>-6.1±0.7</td><td>0.20±0.09</td><td>0.48±0.10</td><td>0.60±0.06</td><td>1.28±0.05</td><td>5.4±1.4</td><td>9.1±0.8</td><td>3.8±0.5</td><td>1.7</td><td>1.5</td><td>-1.09</td></tr><tr><td>GQ Lup</td><td>TWA9A</td><td>-3.0±0.3</td><td>0.19±0.09</td><td>0.58±0.11</td><td>0.78±0.14</td><td>1.87±0.26</td><td>0.9±2.0</td><td>6.3±1.6</td><td>2.0±0.7</td><td>2.5</td><td>1.9</td><td>-0.89</td></tr><tr><td>TW Hya</td><td>TWA25</td><td>12.4±0.1</td><td>0.08±0.07</td><td>0.09±0.06</td><td>0.10±0.05</td><td>0.23±0.11</td><td>5.1±4.9</td><td>15.9±6.2</td><td>6.8±2.6</td><td>0.6</td><td>0.3</td><td>-1.92</td></tr><tr><td>V2129 Oph</td><td>V819TAU</td><td>-7.1±0.6</td><td>0.05±0.10</td><td>0.11±0.06</td><td>0.25±0.07</td><td>0.76±0.16</td><td>0.9±1.1</td><td>1.0±0.7</td><td>0.5±0.3</td><td>0.2</td><td>0.1</td><td>-1.62</td></tr><tr><td>BP Tau</td><td>V819TAU</td><td>15.2±0.6</td><td>0.22±0.09</td><td>0.36±0.06</td><td>0.36±0.09</td><td>0.78±0.10</td><td>3.9±1.4</td><td>7.7±1.4</td><td>3.2±0.8</td><td>1.3</td><td>1.1</td><td>-1.69</td></tr><tr><td>DG Tau</td><td>TWA9A</td><td>15.1±3.4</td><td>0.29±0.22</td><td>0.80±0.20</td><td>1.09±0.16</td><td>2.90±0.51</td><td>11.7±2.7</td><td>16.3±2.1</td><td>6.8±1.0</td><td>6.4</td><td>7.6</td><td>-0.26</td></tr><tr><td>RU Lup</td><td>TWA9A</td><td>-1.3±1.4</td><td>1.14±0.51</td><td>1.70±0.65</td><td>1.65±0.44</td><td>6.27±2.97</td><td>23.2±4.0</td><td>30.1±3.0</td><td>12.2±1.5</td><td>7.8</td><td>11.9</td><td>-0.21</td></tr><tr><td>V2247 Oph</td><td>TWA25</td><td>-5.8±0.5</td><td>-0.02±0.07</td><td>-0.06±0.01</td><td>0.01±0.03</td><td>-0.03±0.02</td><td>0.0±0.2</td><td>0.1±0.1</td><td>-0.1±0.2</td><td>0.009</td><td>—</td><td>-2.54</td></tr><tr><td>DO Tau</td><td>TWA25</td><td>16.1±1.0</td><td>0.50±0.40</td><td>0.78±0.41</td><td>1.07±0.32</td><td>2.33±0.61</td><td>1.0±1.8</td><td>12.7±3.9</td><td>4.0±1.9</td><td>2.3</td><td>3.4</td><td>-0.54</td></tr><tr><td>V347 Aur</td><td>TWA25</td><td>8.1±0.6</td><td>0.31±0.19</td><td>0.46±0.26</td><td>0.37±0.15</td><td>1.12±0.31</td><td>-0.9±1.1</td><td>6.7±2.1</td><td>2.9±1.0</td><td>7.1</td><td>7.3</td><td>-0.62</td></tr><tr><td>J1604</td><td>V819TAU</td><td>-5.8±0.8</td><td>-0.14±0.11</td><td>0.20±0.07</td><td>0.22±0.13</td><td>0.81±0.23</td><td>-2.6±1.8</td><td>-0.1±0.2</td><td>0.3±0.1</td><td>—</td><td>0.013</td><td>-2.81</td></tr><tr><td>PDS 70</td><td>TWA9A</td><td>5.0±0.4</td><td>0.02±0.11</td><td>0.18±0.04</td><td>0.19±0.07</td><td>0.34±0.08</td><td>-0.7±0.9</td><td>0.2±0.2</td><td>0.1±0.2</td><td>0.007</td><td>0.003</td><td>—</td></tr><tr><td>GM Aur</td><td>TWA9A</td><td>14.5±0.3</td><td>-0.00±0.06</td><td>0.13±0.09</td><td>0.11±0.06</td><td>0.31±0.10</td><td>3.4±2.8</td><td>10.2±3.2</td><td>4.7±1.9</td><td>1.3</td><td>1.0</td><td>—</td></tr></tbody></table>

**Notes.** <sup>(a)</sup>We used the convention of positive equivalent width for emission lines, and negative values for absorption lines. <sup>(b)</sup>Slope of the spectral energy distribution measured between 2 and 8 µm.

Before comparing the target and template spectra, we shifted and broadened the template spectra to match the target spectra, using the radial velocities of the targets<fnref n="1" />, and the literature v sin i values of the targets. Due to the IR veiling wavelength dependence (e.g., Alcalá et al. 2021), we measured the veiling in four different spectral regions, 10 710 Å–10 810 Å, 12 840 Å–12 910 Å, 16 120 Å–16 220 Å,

and 22 600 Å–22 690 Å, which we call r<sub>Y</sub>, r<sub>J</sub>, r<sub>H</sub>, and r<sub>K</sub>, representing the YJHK veilings, respectively. The stars RU Lup, DO Tau and, DG Tau present many emission lines along the spectra, probably originating from the accretion shock, which is a characteristic of high-mass accretion rate systems, which prevented us from using the same Y and J spectral regions for the veiling calculations. Then, for these targets, we used the spectral regions 10 760 Å–10 785 Å and 12 400 Å–12 550 Å to measure the r<sub>Y</sub> and r<sub>J</sub> veilings, respectively. On some nights, the 10 710 Å–10 810 Å region of the TW Hya and V2247 Oph spectra presented features that made veiling measurements impossible and we had to use the 10 864 Å–10 920 Å region to measure r<sub>Y</sub> veiling instead.

<footnote n="1">¹ For most targets, the radial velocities were computed using the CCF profiles generated by the SPIRou pipeline, implementing a numerical mask corresponding to the target spectral type. For TW Hya, we computed the radial velocity cross-correlating the target spectra with a WTTS with a similar spectral type.</footnote>
<figure><img src="image_2.png" /><figcaption>**Fig. 1.** Examples of the four spectral regions used to measure the veiling of CI Tau. In each panel, we show two different nights, representing a small and a high veiling estimated for this target. We show the CI Tau spectrum in black and the residual profile (in red) obtained after subtracting the veiled template. The V819 Tau template spectra are also displayed before (orange) and after (blue) applying the veiling correction. The photospheric lines were removed in the residual profiles, showing that the veiling was accurately determined for most nights.</figcaption></figure>

We determined the best veiling value for each target spectrum through a $\chi^2$ minimization routine. The veiling was defined as the ratio of the continuum excess flux to the stellar photospheric flux ($r_\lambda = F_{\lambda \text{excess}}/F_{\lambda \text{phot}}$). Then, zero veiling means the system has no additional excess to the stellar photosphere. In the spectral range of our sample (K2 to M3), the choice of template does not interfere much in the computed veiling, as shown by Folha & Emerson (1999), since the difference between templates of different spectral types is small and produces an almost null relative veiling, within the uncertainties. We also measured the systematic veiling between our templates, which corresponds to the average veiling resulting from the computed veiling comparing each template with the other. We found it to be $r_Y = 0.098 \pm 0.068$, $r_J = 0.02 \pm 0.02$, $r_H = 0.05 \pm 0.02$ and, $r_K = 0.04 \pm 0.02$, which should only affect the veiling determination of stars with very small or no veiling.

In Fig. 1, we present the results for the four photospheric regions used to measure the veiling for CI Tau from two nights, representing spectra with smaller and higher veiling values. We show the target spectrum, and the unveiled and veiled template spectra. We also computed the residual profile, which is the target's spectrum subtracted from the veiled template. Most of the residual profiles show almost no features at the location of photospheric lines, indicating that the veiling measurements were correctly determined for most nights, and that the photospheric lines were correctly removed. However, some of the spectra were

quite noisy and this affected the veiling determinations. The veiling error obtained for each night, written in the plot, comes from a Chi-square minimization process, where we compare the target with the template spectra. Besides taking into account the noise of the target and template spectra<fnref n="2" /> to compute the veiling, there are other errors, such as those associated with the normalization processes that we did not consider in estimating the veiling error.

## 4. Results

We present in Fig. 2 the NIR average veiling over the observation nights, measured for all the four regions, referred to as $r_Y$, $r_J$, $r_H$, and $r_K$. For readability, we split the targets into two groups; systems with $r_K$ higher and smaller than 1. We also show the averaged veiling obtained for each target in Table 2. The veiling values increase from the Y to the K band for most of the targets, similar to the results found in previous works (e.g., McClure et al. 2013; Sousa et al. 2021; Alcalá et al. 2021). Besides this general result, the average veiling for some individual targets remains the same from the Y to K band. For example, the veiling measured for V2247 Oph. However, this target does not present significant veiling in any band, probably due to the

<footnote n="2">² We used as the spectral noise, the standard deviation measured in the continuum of each night, in the spectral region used to measure the veiling.</footnote>
<figure><img src="image_3.png" /></figure>

**Fig. 2.** Average NIR veiling (**left**) and the veiling variability diagnostic (**right**) measured in different wavelength regions. The top and bottom panels show the targets with $r_K$ higher and lower than 1, respectively. The veiling increases from the **Y** to the **K** band for most of the targets. The error bars in the left panel represent the standard deviation of the average veiling over all the observation nights. In this figure and the following figures, the color and symbol codes in the panels identify each target.

<figure><img src="image_4.png" /></figure>

**Fig. 3.** Veiling variability diagnostic as a function of the average veiling. The rms value refers to the root-mean-square of the veiling variability and $\sigma$ is the average error on the veiling measurements. Each panel shows the veiling $r_\lambda$ measured in a different band (YJHK). Systems with higher veiling also present higher veiling variability.

fact that it is a more evolved system, where the accretion and the dust in the inner disk are too faint to be detected. Furthermore, this M-type star makes detecting excess in the IR difficult due to the low contrast between the stellar photosphere and the inner disk emission (e.g., Ercolano et al. 2009). Another example is CI Tau, where the average veiling decreases from the **Y** to the **H** band, followed by an increase to the **K** band. In that case, each band's veiling variability is high and the difference in the **Y** to **H** average veiling is smaller than the standard deviation.

Due to the small number of photospheric lines and the lower S/N values in the **Y** region of the spectra, the veiling in this region was not as well determined, compared to the other bands. For some nights or targets, the veiling is even slightly negative, which has no physical meaning – in such cases, we can assume that the system in this region has zero veiling.

We also checked the variability of the veiling measured in each band, computing for each target the rms of the YJHK veilings, using all the observed nights. To characterize the variability over the noise, we subtracted the average error of the veiling ($\sigma$) from the rms. Then, we computed $S = \sqrt{\text{rms}^2 - \sigma^2}$, as a veiling variability diagnostic, following Cody et al. (2014). We show the results in Figs. 2 and 3 as a function of the band and veiling measured, respectively. Most of the targets present variability above the error level. The veiling in the **H** band appears to be less variable than the other bands, while the **K** band presents the highest variability, driven by the systems with the highest veilings, which are also the most variable ones (see Fig. 3).

The average veiling measured in RU Lup is relatively high compared to the other targets, mainly in the **K** band. These high average veilings are primarily due to the veiling measured in the 2021a period of observations (see Sect. 4.3). In that period, RU Lup also presented an increase in the photometric AAVSO V band brightness (about 0.3 mag smaller than the next observational period, where the veiling starts to go down). Therefore, the average veiling of RU Lup probably is overestimated and does not represent a quiescent value.

## 4.1. Veiling compared to inner disk diagnostics

The emission from the inner disk is claimed to contribute to the veiling. In such a case, we would expect a correlation between the veiling and other inner disk emission diagnostics from the photometric data.

The slope of the spectral energy distribution (SED) is often used as a disk emission diagnostic (e.g., Lada et al. 2006; Muzerolle et al. 2010; Teixeira et al. 2012). Using photometric data from different surveys, such as SDSS (Gunn et al. 1998), Gaia DR2 (Gaia Collaboration 2018), 2MASS, WISE (Wright et al. 2010), Spitzer (Fazio et al. 2004; Rieke et al. 2004), Herschel (PACS), Akari (IRC and FIS), we constructed the SED of the targets and measured the slope of the SED between 2 and 8 µm ($\alpha_{2-8}$), which is the spectral range that indicates significant inner disk emission. The $\alpha_{2-8}$ slope is smaller (more negative) for systems with less or no inner disk emission, and is higher (and even positive) for systems that present inner disk emission excess (e.g., Lada et al. 2006; Muzerolle et al. 2010). We list the $\alpha_{2-8}$ slope computed for the targets in Table 2.

We show in Fig. 4 the NIR veiling in the four spectral regions (YJHK), averaged over all the observation nights, as a function of $\alpha_{2-8}$. The veiling presents a clear linear correlation with the SED slope, mainly from the **J** to **K** band. The **Y** band veiling values are still correlated to the SED slope but less than the other bands, probably due to the contribution from the inner disk emission becoming more important at longer wavelength. The $W_1 - W_2$ ([3.4] - [4.6]) color index (not shown here) comes from the WISE telescope (Wright et al. 2010) is also correlated with the NIR veiling, presenting similar results as $\alpha_{2-8}$.

We could expect that the NIR excess computed using the veiling would scale with the emission excess from NIR photometric data, which is often used as an inner disk emission indicator (Hillenbrand et al. 1998; Rebull 2001; Rebull et al. 2002). We computed the $(H - K_s)_{excess}$, using the observed $H - K_s$ color from 2MASS, corrected for extinction, using the $A_V$ quoted in Table 1 and the SVO Filter (Rodrigo et al. 2012; Rodrigo & Solano 2020) $A_\lambda/A_V$ relations to obtain the $A_H$ and $A_K$ extinctions. Then, we compared this dereddened color excess to the intrinsic color $(H - K_s)_o$ expected for an object with the same spectral type (Pecaut & Mamajek 2013). The color excess is $(H - K_s)_{excess} = (H - K_s)_{obs,dred} - (H - K_s)_o$. We also relate
<figure><img src="image_5.png" /><figcaption>**Fig. 4.** Average NIR veiling ($r_\lambda$) as a function of spectral energy distribution slope between 2 and 8 µm ($\alpha_{2-8}$), which we used as the NIR disk emission diagnostic. The solid line is the linear fit to the data, and the correspondent slope is written in each panel. The NIR veiling seems to scale with the SED slope.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Fig. 5.** Comparison between the color excess computed using the average NIR veiling and the 2MASS photometry. Left: $(H - K_s)_{excess}$, right: $(J - K_s)_{excess}$. See the text for the color excess definition. The dashed line represents a slope equal to 1. The NIR color excesses computed using the average veiling and the 2MASS magnitudes agree for most targets.</figcaption></figure>

the color excess and the excess flux, leading to $(H - K_s)_{excess} = -2.5 \log[(1 + r_H)/(1 + r_K)]$, where we used the veiling definition as $r_\lambda = F_{\lambda excess}/F_{\lambda phot}$. Then we can directly compare the color excess computed using the veiling and using the photometric measurements. We compare the two sides of this equation in Fig. 5, which shows a linear tendency and similar values considering the error of the measurements. Performing similar procedures, we computed the $(J - K_s)_{excess}$, and the results are also presented in Fig. 5. The color excess computed using the 2MASS photometry is dependent on the $A_\nu$ of the systems, and the targets V347 Aur and CI Tau present discrepant $A_\nu$ values in the literature. CI Tau presents $A_\nu = 0.65$ mag (Donati et al. 2020b) and $A_V = 1.90$ mag (Herczeg & Hillenbrand 2014), while the $(J - K_s)_{excess}$ computed using the largest $A_V$ better agrees with the color excess calculated using veiling, the $(H - K_s)_{excess}$, and the mass accretion rate computed in the next section seem to be in better agreement with $A_V = 0.65$ mag; thus, we used this extinction value in the paper. The $A_V$ range of V347 Aur is even larger, with values from 2.2 to 7.5 mag (e.g.,

Dahm & Hillenbrand 2020). The $A_V = 3.4$ mag computed using NIR colors by Connelley & Greene (2010) seems to better represent the value obtained from the veiling calculations.

The 2MASS photometric magnitudes used for this analysis were obtained years before our observations. However, we do not expect a significant change in the NIR magnitudes over the next few years, aside from the daily timescale variation (e.g., Carpenter et al. 2001). It is likely that RU Lup will stand as an exception, as the average veiling was measured in a non-quietescent period. Then, the color excess computed using the veiling is higher than that obtained using the 2MASS magnitudes.

All these relations between the NIR veiling and the inner disk emission diagnostics show that a higher veiled system also presents higher inner disk emission, which is expected if the veiling has a contribution from the inner disk. However, to draw this conclusion, we assumed that these inner disk indicators, such as the slope of the spectral energy distribution, are suitable inner disk diagnostics. Kidder et al. (2021) showed that some of the targets classified as Class III using these disk indicators still show some inner disk emission based on the K band excess. They checked the emission excess of V819 Tau, which we used as a template to compute the veiling, but the H and K excesses found were very small, similar to the systematic veiling we obtained in Sect. 3. The relation between veiling and inner disk emission is clear for the K band and less for other bands, probably due to the influence of another additional continuum source for the other spectral regions and/or a smaller contribution from the inner disk to these shorter wavelengths.

## 4.2. Veiling compared to accretion diagnostics

We know that CTTSs are still accreting gas from the disk and accreting systems typically present strong and variable emission lines that form in the accretion funnel or in the disk wind (e.g., Muzerolle et al. 1998; White & Basri 2003; Edwards et al. 2003; Kwan & Fischer 2011; Alencar et al. 2012). The CFHT/SPIRou wavelength range includes some emission lines from hydrogen and helium, such as Paβ and Bry as well as the He I triplet; in particular, the latter is very sensitive to accretion and ejection processes (e.g., Kwan et al. 2007; Sousa et al. 2021). The dynamics of the circumstellar lines for this sample of stars will be analyzed in an accompanying paper (Sousa et al., in prep.).

We measured the equivalent width of the circumstellar lines, and the average over all observing nights is listed in Table 2. First, we used the equivalent width as an accretion diagnostic (Alcalá et al. 2017), as systems that present larger equivalent widths are supposed to present higher mass accretion rates as well.

We corrected the equivalent width for the veiling as EW = EW(measured)($r_\lambda + 1$), where the $r_\lambda$ represents the veiling computed close to each emission line. In Fig. 6, we show the veiling as a function of the veiling corrected equivalent width of the circumstellar emission lines. We see a clear relationship between the NIR veiling and the accretion diagnostics. It means that higher mass-accretion rate systems also present a higher degree of veiling; this is a similar result to that found by Folha & Emerson (1999), demonstrating that although the veiling shows a contribution from the inner disk emission, it also suggests a connection with the accretion process.

We do not have photometric data simultaneous with our spectra to accurately compute the mass accretion rates using the equivalent width of the emission lines. However, most of our systems' NIR magnitudes are relatively long-term stable. We used
<figure><img src="image_7.png" /><figcaption>**Fig. 6.** Accretion diagnostics as a function of average NIR veiling. *From left to right:* average equivalent width corrected from the veiling of He I, Paβ, and Bry lines. We show a connection between the system's accretion diagnostic and the NIR veiling.</figcaption></figure>

<table><caption>**Table 3.** Stellar parameters.</caption><thead><tr><th>Star</th><th><i>M</i><sub>*</sub><br>(<i>M</i><sub>☉</sub>)</th><th><i>R</i><sub>*</sub><br>(<i>R</i><sub>☉</sub>)</th><th>Ref.</th></tr></thead><tbody><tr><td>CI Tau</td><td>0.90</td><td>2.0</td><td>(1)</td></tr><tr><td>DoAr 44</td><td>1.20</td><td>2.0</td><td>(2)</td></tr><tr><td>GQ Lup</td><td>0.86</td><td>2.26</td><td>(3)</td></tr><tr><td>TW Hya</td><td>0.80</td><td>1.1</td><td>(4)</td></tr><tr><td>V2129 Oph</td><td>1.35</td><td>2.1</td><td>(5)</td></tr><tr><td>BP Tau</td><td>0.70</td><td>1.99</td><td>(6)</td></tr><tr><td>DG Tau</td><td>0.65</td><td>2.05</td><td>(6)</td></tr><tr><td>RU Lup</td><td>1.15</td><td>2.39</td><td>(7)</td></tr><tr><td>V2247 Oph</td><td>0.35</td><td>2.0</td><td>(8)</td></tr><tr><td>DO Tau</td><td>0.60</td><td>2.0</td><td>(9)</td></tr><tr><td>J1604</td><td>1.24</td><td>1.4</td><td>(10)</td></tr><tr><td>PDS 70</td><td>0.76</td><td>1.26</td><td>(11)</td></tr><tr><td>GM Aur</td><td>0.95</td><td>1.71</td><td>(12)</td></tr></tbody></table>

**References.** (1) Donati et al. (2020b); (2) Bouvier et al. (2020); (3) Alcalá et al. (2017); (4) Donati et al. (2011); (5) Alencar et al. (2012); (6) Johns-Krull (2007); (7) Alcalá et al. (2014); (8) Donati et al. (2010); (9) Ricci et al. (2010); (10) Sicilia-Aguilar et al. (2020); (11) Müller et al. (2018); (12) Bouvier et al. (in prep.).

the 2MASS *J* and *K* magnitudes to estimate the continuum flux and then calculate the mass accretion rate using the Paβ and Bry lines, respectively. The star V347 Aur is known to present long-term photometric variations (e.g., Dahm & Hillenbrand 2020), and we did not compute the mass accretion rate of this target.

We followed the procedures described by Gullbring et al. (1998) to compute the line flux and luminosity. The stellar parameters used are listed in Table 3, and we used the stellar distance from the Gaia Collaboration (Gaia Collaboration 2021). We dereddened the 2MASS magnitudes using the same method described in the previous section. To compute the accretion luminosity, we used the fits proposed by Alcalá et al. (2017), which show the relation between the line and accretion luminosities. Then, we determined the accretion rate setting the system inner radius as $5R_*$ (Gullbring et al. 1998). In Table 2, we show the individual mass accretion rates computed using the Paβ and Bry lines. In Fig. 7, we show the average mass accretion rate as a function of the *Y* to *K* band veiling. Once again, we can connect the highest accreting system with the highest NIR veiling computed in the four bands.

<figure><img src="image_8.png" /><figcaption>**Fig. 7.** Average mass accretion rate as a function of average NIR veiling. The average mass accretion was computed using the line fluxes of Paβ and Bry. The error bar is the standard deviation between the two measurements. We see an association between the veiling and the mass accretion rate.</figcaption></figure>

### 4.3. Veiling night-to-night variability

In this study, we have access to observations obtained on different nights and sometimes different observational periods for our sample of stars. This allowed us to analyze the night-to-night veiling variation and a possible long-term veiling variability on a timescale of two years. In Fig. 8, we show the veiling measured as a function of the observation dates. We used the modified Lomb-Scargle periodogram (Horne & Baliunas 1986) to study a possible periodicity of the veiling variations. We performed the periodogram analysis of the veiling in two ways: using all the observed nights and computing one periodogram per observational period. We searched for periods from 2 to 15 days or limited the search to the number of observed days, when the target was observed for less than 15 days. However, we did not find a significant periodic signal for most systems. We also phased the veiling using the known stellar rotation period of the targets,
<figure><img src="image_9.png" /><figcaption>**Fig. 8.** Night-by-night veiling values measured in four different spectral regions. We show each observational period per target in an individual panel. A missing veiling point in a specific band means that the spectral regions were subject to effects that prevented the veiling from being measured. In addition to the variability in terms of veiling, it is also not clearly periodic, at least on the timescale of stellar rotation. For more details, see text.</figcaption></figure>

and again, most of the systems do not seem to vary in phase with
the stellar rotation.

In Sect. 4, we used the veiling error as a threshold of the
veiling variability, using the rms as a variability diagnostic, see
Fig. 2. Then, any point above this limit represents a true vari-
ability of the veiling; below this threshold, the variability is at

the same level as the errors associated with the veiling. We can
see that the veiling rms presents values above the threshold of
veiling variability for most of the targets and veiling bands. Only
V2247 Oph shows most values below the threshold limit, while
a few bands of TW Hya and DoAr 44 present most of the rms
values below the threshold. Then, for V2247 Oph, we cannot
attribute a variability to the veiling, which is consistent with this target presenting a small veiling and weak level of accretion and inner disk emission diagnostics.

We know the veiling is generally variable and we can trust the detection of veiling variability for most targets. Figure 8 shows that for all the veiling variable targets, the veiling varies on a timescale of at least one day and the veiling computed in the four spectral regions presents the same variability timescale. These results show that whatever region the IR veiling comes from, this region is dynamic and its flux changes on a timescale of days.

We also investigated veiling variability on a timescale of months to a few years; a change in the veiling in that timescale can have a different origin from the day-scale veiling variability, discussed above. The latter can be associated with the dynamic of the system's rotation, while the possible long-term veiling variability reflects a change in the system's accretion and/or inner disk conditions. In Fig. 9, we present the averaged veiling measured at each observational period. This plot displays nine systems that were observed in more than one observational season. Most targets do not present a significant difference in the veiling along the observational period. The *K* band veiling of CI Tau and the *YJH* veiling of DO Tau show a possible small change in this timescale. Then, we conclude that the veiling variability on a timescale of months-to-years is on the same order of magnitude as the day-to-day variability. RU Lup is the unique target with an evident change in the veiling along the observational periods in the four bands, but much more pronounced in the *K* band, along with a high standard deviation. We associate this change in the veiling with an occasional high accretion episode that occurred in 2021a, and despite the veiling still being high in 2021b, it seems to start to diminish later on. In 2022a, it is even smaller, but we have only one observation to serve as the basis for this assumption. The circumstellar emission lines' equivalent widths corroborate with this assumption, as they increase in 2021a and start to decrease in the subsequent observation periods, similarly to the veiling. In contrast, the average veiling is stable, at least for most of the stars we analyzed, except for very high accreting systems, such as RU Lup, which can present episodic high veiling. Furthermore, in the same observational period, some targets show a few episodes of an increase in veiling, such as V347 Aur (Fig. 8); however, the average veiling values are still sustained.

## 5. Discussion

The dependence of veiling on wavelength is ubiquitous from the UV to the NIR range. The veiling in the optical domain decreases from the blue to the red part of the spectra, which is an effect of the decrease of the accretion spot continuum contribution (e.g., Calvet & Gullbring 1998); also, the veiling also does not vary for some wavelength ranges (e.g., Basri & Batalha 1990). On the other hand, in the IR range, the veiling increases with wavelength, as seen in Figs. 2 and 8, which is in agreement with similar results in the literature (e.g., Fischer et al. 2011; Alcalá et al. 2021).

The average veiling value for the entire sample of CTTS is $\langle r_Y \rangle = 0.2 \pm 0.3$, $\langle r_J \rangle = 0.4 \pm 0.4$, $\langle r_H \rangle = 0.5 \pm 0.5$, and $\langle r_K \rangle = 1.4 \pm 1.6$. We note that in the *Y* band, the average veiling is the lowest. Over these wavelengths from *Y* to *K*, the veiling can have contributions from different sources. For example, we expect the veiling in the *K* band to present more contributions from the inner disk than in the *Y* band. In Fig. 10, we show the *YJH* veiling as a function of the *K* band veiling. We can see

<figure><img src="image_10.png" /><figcaption>**Fig. 9.** Average of the NIR veiling computed at each observational period, with the error bar as the standard deviation of the computed mean veiling. We merged the measured veiling in 2020a and 2020b of the stars V2129 Oph and GQ Lup, as they are successive observations. The mean NIR veiling seems to be stable for most of the targets for a few months or years.</figcaption></figure>

that the *J* and *H* veilings seem to increase as the *K* band veiling increases; however, there is smaller correlation with the *Y* band veiling. These results are supported by the analysis of the correlation of the veiling samples and the linear fit showed in the Fig. 10. We computed the linear correlation coefficient ($r$) between two samples, where $r = 1$ represents a perfect correlation and when $r = 0$ there is no correlation. The *Y* and *K* band veilings present a correlation coefficient of 0.87, while the coefficient of the *J* and *H* and the *K* band are 0.98 and 0.96, respectively. Similar results were found by Cieza et al. (2005), where they compare the excess in the *J* and *H* bands with the *K* band excess, showing that both presented a linear correlation with the *K* band excess, and this was explained as due to the *JHK* excess arising from the same region.

The NIR veiling should be the result of a combination of physical processes. Alcalá et al. (2021) computed the NIR
<figure><img src="image_11.png" /><figcaption>**Fig. 10.** *YJH* veiling as a function of the *K* band veiling. Each veiling value corresponds to the mean of all the observed nights, and the error bar is its standard deviation. The solid line is the linear fit to the data, and the fitted line's slope is written in each panel. The correlation between the *K* band and the *YJH* bands veiligings increases from the *Y* to the *H* band.</figcaption></figure>

veiling at several wavelengths for a sample of very high-mass accretion rate systems, including DG Tau (included in our sample). These authors fitted the veiling as a function of wavelength using a blackbody function and found temperatures compatible with the presence of dust in the inner disk. However, the temperature was too high (>2000 K) for the dust to survive in a few of their cases. They also argued that the veiling should have a contribution from the hot gas inside the disk sublimation radius, and similar results were found by Fischer et al. (2011). To investigate this proposition, we looked at the CO bandhead at 2.3 µm of the CFHT/SPIRou data. This band, when in emission, is expected to form in the hot gas in the inner disk. Using the *K* band veiling, we veiled the template and removed the photospheric lines of the CO bandhead to obtain the residual CO profiles. Most targets do not present clear signs of CO emission, showing that this band is strictly photospheric. However, a few residual profiles of V347 Aur, which is a Class I object, along with DO Tau and RU Lup, which are strong accretors, present CO emission in some observations, indicating the presence of hot gas in the inner disk. In particular, RU Lup presents these hints of CO emission in the observational period when the veiling was high, and the system probably ensued a high episodic accretion. A further analysis of the CO bandhead is beyond the scope of this paper and it will be carried out in a dedicated paper exploring the significance of these CO emissions.

In the previous section, we show that the NIR veiling, mainly $r_K$, presents a good correlation with the inner disk emission diagnostics obtained from the NIR photometric data and SED fit, demonstrating that the NIR veiling has an important contribution from the inner disk. However, we also see a correlation between veiling and accretion diagnostics. A high accretion-rate system presents larger veiling values in the IR. This shows that high-mass accretion rate systems should feature higher inner disk heating, thus higher temperatures and stronger inner disk emission excess as a consequence. Espaillat et al. (2022) fit most of the continuum spectra from NUV to NIR of the accreting star CVSO 109A quite well, using a combination of emission from the accretion shock (multiple funnel flow model) on the stellar surface and emission from the irradiated inner edge of the dusty disk. However, the inner disk and accretion shock model do not adequately reproduce the continuum excess in the *Y* to *J* band.

Indeed, while the veiling in the *J* to *K* band seems to point to a significant contribution on the part of the inner disk emission excess, the *Y* band veiling origin is still unknown. Dodin & Lamzin (2013) predicted significant veiling in the *Y* band from the accretion spot.

<figure><img src="image_12.png" /><figcaption>**Fig. 11.** Average NIR veiling as a function of the system's inclination with respect to our line of sight. The solid line is the data's linear fit, and the fitted line's slope is given in each panel. We do not consider the discrepant targets (TW Hya and V2247 Oph) to fit the data. The NIR veiligings tend to be anti-correlated with the system's inclination, see text.</figcaption></figure>

They argued that the accretion spot (continuum emission and emission lines formed in the accretion shock) could account for optical and near IR veiling up to the *J* band. Unfortunately, our *Y* band veiling was not as well computed as for the other bands due to several issues in this spectral region and given the photospheric lines are not quite so prominent. Despite these obstacles, the *Y* band veiling agrees on a smaller scale with both accretion and inner disk diagnostics.

We checked if the inclination of the system has any impact on the measured veiling values. In Fig. 11, we show the NIR veiling as a function of the inclination of the system with respect to our line of sight, listed in Table 1. In addition, two discrepant systems (TW Hya and V2247 Oph), we can see an anti-correlation between veiling and the system's inclination. This anti-correlation is not pronounced, as we can see in the linear fit slope, due to the spread of points (the correlation coefficients between the inclination and the *YJHK* veiligings are -0.64, -0.76,
−0.80, and −0.70, respectively), but the decreasing tendency of veiling with inclination is clear. We also advise that the inclinations used for DG Tau and DO Tau are the outer disk inclination, while the disk and the stellar inclination are not necessarily the same (Bohn et al. 2022). If confirmed, this anti-correlation can be due to a geometric effect: the more the system is inclined, the less we see (from the inner disk edge) where the nIR veiling is supposed to arise. The two targets, TW Hya and V2247 Oph, which do not seem to follow this tendency, do not have dust in the inner disk any longer and they are known to have gaps or holes in their inner disks (Calvet et al. 2002; Pontoppidan et al. 2008). In that case, independently of the system's inclination, we would not expect to detect IR veiling, assuming that the IR veiling is due to dust emission in the inner disk.

## 6. Conclusion

In this work, we analyze the NIR veiling computed using high-resolution data from CFHT/SPIRou of a sample of 14 low-mass young stars. We found the veiling to increase from the Y to the K band, as a result of the increase of the emission contribution from the inner disk as a function of wavelength.

The veiling correlates with other photometric inner disk diagnostics, such as color excess and the slope of the spectral energy distribution, mainly in the JHK band, providing further evidence that the NIR veiling arises from hot dust in the inner disk. We also found a linear correlation between veiling and the accretion properties of the system. This shows that accretion contributes to inner disk heating and, consequently, to the inner disk emission excess. This effect is enhanced in high-mass accretion rate systems that also present a denser inner disk and higher inner disk emission (e.g., Sullivan & Kraus 2022).

We analyzed the NIR veiling variability through the modified Lomb-Scargle periodogram and we did not find any significant periodic signal in the four bands in timescales typical of stellar rotation (<15 days), which also suggests the veiling comes from the dust emission in the inner disk. However, we show that the veiling is variable for most targets on a timescale of at least one day. Besides the night-by-night veiling variability, the mean NIR veiling per season appears to be mostly stable, for most targets and on timescales of several months to years.

**Acknowledgements.** We thank the referee for the suggestions that helped to clarify this paper. We want to thank Claire Moutou, Sylvie Cabrit, Nicolas Grosso, and Konstantin Grankin for carefully reading the manuscript and giving suggestions to improve the paper. This project has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme (grant agreement No. 742095; SPIDI: Star-Planets-Inner Disk-Interactions; http://www.spidi-eu.org and grant agreement No. 740651 NewWorlds). We acknowledge financial support from CNPq, CAPES and Fapemig. We acknowledge funding from the French National Research Agency (ANR) under contract number ANR-18-CE31-0019 (SPlaSH). This research has made use of the SVO Filter Profile Service (http://svo2.cab.inta-csic.es/theory/fps/) supported from the Spanish MINECO through grant AYA2017-84089. The authors wish to recognize and acknowledge the very significant cultural role and reverence that the summit of MaunaKea has always had within the indigenous Hawaiian community. We are most fortunate to have the opportunity to conduct observations from this mountain. We acknowledge with thanks the variable star observations from the AAVSO International Database contributed by observers worldwide and used in this research.

## References

* Alcalá, J. M., Natta, A., Manara, C. F., et al. 2014, A&A, 561, A2

* Alcalá, J. M., Manara, C. F., Natta, A., et al. 2017, A&A, 600, A20

* Alcalá, J. M., Gangi, M., Biazzo, K., et al. 2021, A&A, 652, A72

* Alencar, S. H. P., & Batalha, C. 2002, ApJ, 571, 378

* Alencar, S. H. P., Bouvier, J., Walter, F. M., et al. 2012, A&A, 541, A116

* Antoniucci, S., Nisini, B., Biazzo, K., et al. 2017, A&A, 606, A48

* Basri, G., & Batalha, C. 1990, ApJ, 363, 654

* Bohn, A. J., Benisty, M., Perraut, K., et al. 2022, A&A, 658, A183

* Bouvier, J., Alecian, E., Alencar, S. H. P., et al. 2020, A&A, 643, A99

* Calvet, N., & Gullbring, E. 1998, ApJ, 509, 802

* Calvet, N., Hartmann, L., & Strom, S. E. 1997, ApJ, 481, 912

* Calvet, N., D'Alessio, P., Hartmann, L., et al. 2002, ApJ, 568, 1008

* Carpenter, J. M., Hillenbrand, L. A., & Skrutskie, M. F. 2001, AJ, 121, 3160

* Chiang, E. I., & Goldreich, P. 1997, ApJ, 490, 368

* Chiang, E. I., & Goldreich, P. 1999, ApJ, 519, 279

* Cieza, L. A., Kessler-Silacci, J. E., Jaffe, D. T., Harvey, P. M., & Evans, N. J., II 2005, ApJ, 635, 422

* Cody, A. M., Stauffer, J., Baglin, A., et al. 2014, AJ, 147, 82

* Connelley, M. S., & Greene, T. P. 2010, AJ, 140, 1214

* Cook, N. J., Artigau, É., Doyon, R., et al. 2022, PASP, 134, 114509

* Dahm, S. E., & Hillenbrand, L. A. 2020, AJ, 160, 278

* Dahm, S. E., Slesnick, C. L., & White, R. J. 2012, ApJ, 745, 56

* Davies, C. L. 2019, MNRAS, 484, 1926

* Dodin, A. V., & Lamzin, S. A. 2013, Astron. Lett., 39, 389

* Donati, J., Jardine, M. M., Gregory, S. G., et al. 2007, MNRAS, 380, 1297

* Donati, J. F., Jardine, M. M., Gregory, S. G., et al. 2008, MNRAS, 386, 1234

* Donati, J. F., Skelly, M. B., Bouvier, J., et al. 2010, MNRAS, 402, 1426

* Donati, J. F., Gregory, S. G., Alencar, S. H. P., et al. 2011, MNRAS, 417, 472

* Donati, J. F., Gregory, S. G., Alencar, S. H. P., et al. 2012, MNRAS, 425, 2948

* Donati, J. F., Hébrard, E., Hussain, G. A. J., et al. 2015, MNRAS, 453, 3706

* Donati, J. F., Kouach, D., Lacombe, M., et al. 2018, in Handbook of Exoplanets, eds. H. J. Deeg, & J. A. Belmonte, 107

* Donati, J. F., Kouach, D., Moutou, C., et al. 2020a, MNRAS, 498, 5684

* Donati, J. F., Bouvier, J., Alencar, S. H., et al. 2020b, MNRAS, 491, 5660

* Edwards, S., Fischer, W., Kwan, J., Hillenbrand, L., & Dupree, A. K. 2003, ApJ, 599, L41

* Edwards, S., Fischer, W., Hillenbrand, L., & Kwan, J. 2006, ApJ, 646, 319

* Ercolano, B., Clarke, C. J., & Robitaille, T. P. 2009, MNRAS, 394, L141

* Espaillat, C. C., Herczeg, G. J., Thanathibodee, T., et al. 2022, AJ, 163, 114

* Faesi, C. M., Covey, K. R., Guttermuth, R., et al. 2012, PASP, 124, 1137

* Fazio, G. G., Hora, J. L., Allen, L. E., et al. 2004, ApJS, 154, 10

* Fischer, W., Edwards, S., Hillenbrand, L., & Kwan, J. 2011, ApJ, 730, 73

* Flores, C., Connelley, M. S., Reipurth, B., & Boogert, A. 2019, ApJ, 882, 75

* Folha, D. F. M., & Emerson, J. P. 1999, A&A, 352, 517

* Frasca, A., Biazzo, K., Alcalá, J. M., et al. 2017, A&A, 602, A33

* Gaia Collaboration (Brown, A. G. A., et al.) 2018, A&A, 616, A1

* Gaia Collaboration (Brown, A. G. A., et al.) 2021, A&A, 649, A1

* Gullbring, E., Hartmann, L., Briceño, C., & Calvet, N. 1998, ApJ, 492, 323

* Gully-Santiago, M. A., Herczeg, G. J., Czekala, I., et al. 2017, ApJ, 836, 200

* Gunn, J. E., Carr, M., Rockosi, C., et al. 1998, AJ, 116, 3040

* Hartigan, P., Hartmann, L., Kenyon, S., Hewett, R., & Stauffer, J. 1989, ApJS, 70, 899

* Hartigan, P., Kenyon, S. J., Hartmann, L., et al. 1991, ApJ, 382, 617

* Hartmann, L. W., & Kenyon, S. J. 1990, ApJ, 349, 190

* Herczeg, G. J., & Hillenbrand, L. A. 2014, ApJ, 786, 97

* Hillenbrand, L. A., Strom, S. E., Calvet, N., et al. 1998, AJ, 116, 1816

* Horne, J. H., & Baliunas, S. L. 1986, ApJ, 302, 757

* Ingleby, L., Calvet, N., Herczeg, G., et al. 2013, ApJ, 767, 112

* Johns-Krull, C. M. 2007, ApJ, 664, 975

* Johns-Krull, C. M., & Valenti, J. A. 2001, ApJ, 561, 1060

* Johns-Krull, C. M., Valenti, J. A., & Linsky, J. L. 2000, ApJ, 539, 815

* Kidder, B., Mace, G., López-Valdivia, R., et al. 2021, ApJ, 922, 27

* Kounkel, M., Covey, K., Moe, M., et al. 2019, AJ, 157, 196

* Kwan, J., & Fischer, W. 2011, MNRAS, 411, 2383

* Kwan, J., Edwards, S., & Fischer, W. 2007, ApJ, 657, 897

* Lada, C. J., Muench, A. A., Luhman, K. L., et al. 2006, ApJ, 131, 1574

* McClure, M. K., Calvet, N., Espaillat, C., et al. 2013, ApJ, 769, 73

* Müller, A., Keppler, M., Henning, T., et al. 2018, A&A, 617, L2

* Muzerolle, J., Hartmann, L., & Calvet, N. 1998, AJ, 116, 455

* Muzerolle, J., Allen, L. E., Megeath, S. T., Hernández, J., & Guttermuth, R. A. 2010, ApJ, 708, 1107

* Nguyen, D. C., Brandeker, A., van Kerkwijk, M. H., & Jayawardhana, R. 2012, ApJ, 745, 119

* Pecaut, M. J., & Mamajek, E. E. 2013, ApJS, 208, 9

* Pecaut, M. J., & Mamajek, E. E. 2016, MNRAS, 461, 794

* Pontoppidan, K. M., Blake, G. A., van Dishoeck, E. F., et al. 2008, ApJ, 684, 1323

* Preibisch, T., & Zinnecker, H. 1999, AJ, 117, 2381

* Rebull, L. M. 2001, AJ, 121, 1676

* Rebull, L. M., Makidon, R. B., Strom, S. E., et al. 2002, AJ, 123, 1528

* Rei, A. C. S., Petrov, P. P., & Gameiro, J. F. 2018, A&A, 610, A40

* Ricci, L., Testi, L., Natta, A., et al. 2010, A&A, 512, A15
Rieke, G. H., Young, E. T., Engelbracht, C. W., et al. 2004, ApJS, 154, 25

Rodrigo, C., & Solano, E. 2020, XIV.0 Scientific Meeting (virtual) of the Spanish Astronomical Society, 182

Rodrigo, C., Solano, E., & Bayo, A. 2012, SVO Filter Profile Service Version 1.0, IVOA Working Draft 15 October 2012

Sicilia-Aguilar, A., Fang, M., Roccatagliata, V., et al. 2015, A&A, 580, A82

Sicilia-Aguilar, A., Manara, C. F., de Boer, J., et al. 2020, A&A, 633, A37

Simon, M. N., Pascucci, I., Edwards, S., et al. 2016, ApJ, 831, 169

Sousa, A. P., Bouvier, J., Alencar, S. H. P., et al. 2021, A&A, 649, A68

Stempels, H. C., Gahm, G. F., & Petrov, P. P. 2007, A&A, 461, 253

Sullivan, K., & Kraus, A. L. 2022, ApJ, 928, 134

Teixeira, P. S., Lada, C. J., Marengo, M., & Lada, E. A. 2012, A&A, 540, A83

Thanathibodee, T., Molina, B., Calvet, N., et al. 2020, ApJ, 892, 81

Torres, C. A. O., Quast, G. R., da Silva, L., et al. 2006, A&A, 460, 695

Valenti, J. A., Basri, G., & Johns, C. M. 1993, AJ, 106, 2024

White, R. J., & Basri, G. 2003, ApJ, 582, 1109

Wright, E. L., Eisenhardt, P. R. M., Mainzer, A. K., et al. 2010, AJ, 140, 1868