# Inhomogeneous dust eclipses in young stars. The case of CQ Tauri

A. V. Dodin,1 E. A. Suslina1
1 Sternberg Astronomical Institute, M.V. Lomonosov Moscow State University,
13, Universitetskij pr., 119234, Moscow, Russia
E-mail<EMAIL_ADDRESS>(AVD)

(Accepted XXX. Received YYY; in original form ZZZ)

###### Abstract

We find that dust clouds which eclipse young stars obscure the stellar disc
inhomogeneously. In the particular case of CQ Tau, we find isolated optically
thick structures with sizes $\lesssim 0.6R_{*}$ and derive the typical $A_{V}$
gradient in the plane of the sky, finding it as high as a few magnitudes per
stellar radius. The large extinction gradients and complex structure of the
obscuring clouds lead not only to a noticeable Rossiter–McLaughlin effect, but
also to complex and variable shaping of stellar absorption lines.

###### keywords:

stars: circumstellar matter – stars: variables: T Tauri, Herbig Ae/Be – stars:
individual: CQ Tau

††pubyear: 2021††pagerange: Inhomogeneous dust eclipses in young stars. The
case of CQ Tauri–Inhomogeneous dust eclipses in young stars. The case of CQ
Tauri

##  1 Introduction

Photometric variability caused by circumstellar dust has been found in stars
of different types: in rare cases in main sequence stars (Schmidt, 2019), and
more often in evolved stars (AGB, post-AGB, R CrB). The most diverse
manifestations of circumstellar dust can be found in pre-main-sequence stars
(TTS, HAeBe) that is related to the active gas–dust discs around them. This
group includes several types of variable stars: the UX Ori-type stars, which
are distinguished by irregular Algol-like variability (Waters & Waelkens,
1998); ‘dippers’ which show a quasi-periodic variability, probably related to
a warped protoplanetary disc (Bouvier et al., 1999), with about 20 per cent of
young stars showing this kind of variability (Cody et al., 2014); stars with
long and deep dust eclipses like those recently observed in RW Aur and AA Tau
(Schneider et al., 2015a, b; Günther et al., 2018; Dodin et al., 2019); and
stars with irregular small-amplitude dimmings (Petrov & Kozack, 2007; Cody et
al., 2014; Siwak et al., 2018; Gahm et al., 2018).

Dust eclipses provide an opportunity to characterize the circumstellar dust by
relatively simple methods in the optical region: there are numerous works in
which the characteristic sizes of dust grains are inferred from a colour
behaviour of the star. In particular, a greyer extinction is commonly
interpreted as evidence for a larger average grain size (Steenman & The, 1989,
1991; Bibo et al., 1992; Friedemann et al., 1992; Cohen et al., 2003;
Hernández et al., 2004; Manset et al., 2009; Barsunova et al., 2015; McGinnis
et al., 2015; Petrov et al., 2015; Rodriguez et al., 2016; Schneider et al.,
2015a, b). However, in the case of inhomogeneous eclipses the extinction law
is also flattened (Natta & Panagia, 1984; Grinin et al., 2019).

By ‘an inhomogeneous eclipse’, we mean that extinction varies across the
stallar disc. It is obvious that homogeneous eclipses do not exist, because
each eclipse must have an ingress and egress, during which the stellar disc is
obscured unevenly. However, we can neglect inhomogeneities if noticeable
changes in extinction across the cloud occur on scales much large than the
stellar radius. Some information about the extinction gradient in
circumstellar clouds can be obtained from the variability rate during the
ingress or egress with some assumptions about the tangential velocity of the
clouds. A more promising manifestation of dust inhomogeneities is the
distortions that they produce in stellar absorption line profiles. These
distortions arise due to partial obscuration of the rotating stellar surface,
which results in the weakening of some parts of the profiles, corresponding to
the radial velocities of the obscured regions. The Rossiter–McLaughlin effect
– rapid changes of the measured radial velocity of the star in eclipsing
binaries during the eclipse – has the same nature. The presence of the
Rossiter–McLaughlin effect in the UX Ori-type stars has been suspected by
Grinin & Potravnov (2013) as a possible explanation of discrepancies in radial
velocity measurements of RZ Psc.
A large number of high-quality spectra have been obtained during studies of
magnetic fields in young stars (Alecian et al., 2013; Villebrun et al., 2019).
These spectra have revealed distortions in LSD profiles for some stars, but
the causes had not been determined, beyond suggestions of an association with
circumstellar gas. Examining these archived spectra, we found that these
distortions resemble those expected for the Rossiter–McLaughlin effect. The
greatest body of observational data was found for CQ Tauri, which shows well-
pronounced distortions, and we will restrict ourselves to this case. In this
paper we will show that the observed distortions are indeed related to the
clumpy structure of the dust environment, and we will obtain quantitative
information about spatial scales of the clumps.

##  2 Observational data

We use archival111http://www.cadc-ccda.hia-iha.nrc-cnrc.gc.ca spectra of CQ
Tau, obtained with ESPaDOnS (CFHT) on 2012 December 22, 25, 28, and 29 in the
spectral range 3690–10480 Å. On each of the four observational nights, 8
spectra with exposures of 750 sec were obtained with a signal-to-noise ratio
of about 80. The telescope archive provides spectra, which have been already
processed with upena222https://www.cfht.hawaii.edu/Instruments/Upena.
Additionally, for ease of analysis, we combined spectral orders,
simultaneously reducing the spectral resolution to $R=30\,000$ by convolution
with the Gaussian profile. Then, each spectrum was cleaned of telluric lines
with molecfit333http://www.eso.org/sci/software/pipelines/skytools/molecfit
(Kausch et al., 2015). Spectra obtained over one night were averaged,
normalized to the continuum level, and corrected for Earth’s radial velocity.
We also looked for line profile variability within 2 hours of observations for
each of the nights, but no significant changes were found, probably due to an
insufficient signal-to-noise ratio.

##  3 Properties of CQ Tau

###  3.1 Atmospheric parameters

Stellar atmospheric parameters have been estimated by Alecian et al. (2013) as
$T_{\rm eff}=6750\pm 300$ K, $\log g=4.0$; Villebrun et al. (2019) give a
similar value for $T_{\rm eff}=6800\pm 290$ K. We also try to determine the
parameters in the LTE approximation, using
sme444http://www.stsci.edu/~valenti/sme.html (Piskunov & Valenti, 2016) with
stellar atmosphere models from the MARCS2012 grid (Gustafsson et al., 2008;
Bergemann et al., 2012), and the line list from the
VALD555http://vald.astro.uu.se/ database (Ryabchikova et al., 2015). We found
that, in the expected parameter range, the Balmer lines are sensitive only to
$T_{\rm eff}$ and lead to value of $T_{\rm eff}=6850\pm 100$ K, which agrees
with previous estimates within the uncertainties. Gravity $\log g$ is poorly
constrained: the wings of the Paschen lines lead to $\log g=4.3,$ which seems
too high for a young star, while the ionization equilibrium of Fe i/ii
requires a very low value of $\log g<3.$ It is possible that the Fe ii lines
are strengthened by non-LTE effects (overionization). Fortunately for our
analysis, other lines are weakly dependent on $\log g,$ allowing us to set
$\log g=3.5.$ To better reproduce the absorption lines, we also adjusted
microturbulent velocity $v_{\rm mic}=2.3$ km s-1 and chemical abundances for 7
elements: $\varepsilon_{\rm Mg}=-4.20,$ $\varepsilon_{\rm Si}=-4.54,$
$\varepsilon_{\rm Ca}=-5.30,$ $\varepsilon_{\rm Ti}=-7.11,$ $\varepsilon_{\rm
Fe}=-4.57,$ $\varepsilon_{\rm Ni}=-5.90,$ $\varepsilon_{\rm Ba}=-9.30$
($\varepsilon$ is the logarithm of relative abundance by number). Some
elements (Mg, Ca, Ba) are overabundant by $\sim 0.3$ dex with respect to solar
values, however, we suspect that these deviations and the enhancement of Fe ii
lines are related and probably caused by non-LTE effects.

###  3.2 Radial and rotational velocities

The radial velocity, as well as the line profiles show variability (see
Section 4), the study of which is the subject of this work. As a reference
value, we chose $V_{\rm r}=10$ km s${}^{-1},$ measured from the spectrum of
2012 December 22. For the rotational broadening, we will apply $v\sin i=100$
km s-1 (Alecian et al., 2013). More accurate estimates for these quantities
can be obtained only after clarifying the nature of their variability.

###  3.3 Photometric variability

Available photometric data show that, from the end of the 19th century until
the middle of the 20th century, CQ Tau did not show photometric activity,
fluctuating near 9 magnitude in the photographic $B$-band (Grinin et al.,
2008). Around 1940, the star entered into a phase of irregular photometric
variability with dimmings of up to 3 magnitudes. This photometric behaviour,
typical for the UX Ori stars, persists to this day. We do not have
simultaneous photometry for the spectral observations, however, from the AAVSO
data on 2012 December 19 and 28 the photoelectric $V$ magnitudes were 11.0 and
10.6. The out-of-eclipse $V$ magnitude is about 9 (Shakhovskoj et al., 2005);
therefore, all spectra were taken in the eclipsed state.

##  4 Observational evidence for inhomogeneity

Figure 1: Representative parts of CQ Tau spectra. The spectrum obtained on
2012 December 22 is taken as a reference and shown with the black lines. The
spectra obtained on 2012 December 25, 28, and 29 are shown with red, green,
and blue lines, correspondingly. The ions that make the main contribution to
the lines are labelled on top.

Small but representative parts of the CQ Tau spectra are shown in Fig. 1. It
can be seen that all lines have similar distortions, which nearly conserve
their equivalent widths (EW), and which do not depend on line properties,
although the picture of these distortions is somewhat complicated due to
overlapping of the lines. The most impressive distortions were observed on
2012 December 25. Instrumental effects can be ruled out by the absence of
these distortions in the telluric lines.

To describe the observed distortions over the entire spectrum, we selected 146
lines (most of them being blends), for which we determined the radial
velocities $\Delta V_{{\rm r}}$ relative to the spectrum taken on 2012
December 22. Fig. 2 shows that the changes in the radial velocities depend
weakly on the wavelength, the excitation energy, or the equivalent width, and
are reduced mainly to a simple shift, as is expected for an inhomogeneous
obscuration of the rotating stellar surface (see Section 5). We also show the
relative changes in the equivalent widths, which are constant within 10 per
cent for ${\rm EW}>0.5$ Å. The increasing scatter at low EWs is caused mainly
by errors in the continuum level, rather than by physical variability. One can
suspect small trends in the presented diagrams, a possible origin of which
will be discussed in Section 5.

Because the distortions do not depend on physical parameters of the lines, it
is natural to suppose that they are induced by geometrical effects like
partial obscuration of the stellar surface.

Figure 2:  Shift of the radial velocity of individual lines as function of
their excitation energy ($E$), equivalent width (EW), and wavelength
($\lambda$). The bottom panel shows relative changes in EW. Colours correspond
to Fig. 1.

##  5 Simulations

Since the main object of this study is CQ Tau, we perform our simulations for
the parameters which best describe the CQ Tau spectrum, however, obtained
results can be qualitatively extended to other stars with thin atmospheres.

To simulate distortions of the spectral lines, which appear to be due to
partial obscuration, we calculate the specific intensity $I(\lambda,\mu)$ with
sme v.554 using the parameters of CQ Tau from Section 3.1. The flux is
obtained by integrating the intensities over the stellar disc accounting for
the extinction in circumstellar clouds:

$f=\frac{\int I\left[\lambda(1+V/c),\mu\right]e^{-\tau}dS}{\int I_{\rm
c}(\lambda,\mu)e^{-\tau}dS},$ (1)

where $dS$ is a projected area of the surface element; $V=V(x,y)$ is the
radial velocity of the element, which for a rigid rotation is $v\sin
ix/R_{*}$; $\tau=\tau(x,y)$ is the distribution of the optical thickness of
the obscuring clouds in the plane of the sky $(x-y)$; $\mu=\mu(x,y)$ is the
cosine of the angle between the line of sight and the normal to the stellar
surface. The $y$-axis is directed along the projection of the rotation axis;
the $x$-axis is orthogonal to $y$. It is worth noting that, by denoting
$I/I_{\rm c}$ as $u$ and $I_{\rm c}e^{-\tau}/\int{I_{\rm c}e^{-\tau}dS}$ as a
weight function $w$, we can rewrite equation 1 in the form of a weighted
average $f=\int{uwdS},$ i.e. the observed spectrum can be represented as a
linear combination of continuum-normalized and Doppler-shifted intensities.

To illustrate how the line profile is formed, we show a few simple cases in
Fig. 3, in which small parts of the stellar surface are opened while the rest
of the stellar disc is covered by the cloud with infinite $\tau,$ and vice
versa. Differences in the limb darkening between the line and continuum lead
to the dependence of $I/I_{\rm c}$ on $\mu,$ which makes the equivalent width
variable, but maximal changes for the examples in the figure do not exceed 10
per cent in the case of open squares in the upper panels. In more realistic
cases, where many angles $\mu$ contribute to the profile, changes of EW will
be much less.

Figure 3: Effects of open (the upper panels) or closed (the bottom panels)
elements of the stellar surface on the line profile. The dashed line
represents the line profile for a clear stellar disc.  Figure 4: A simulated
image of the stellar disc (with limb darkening), covered by an inhomogeneous
dust cloud and the corresponding line profile. The dashed line represents the
out-of-eclipse profile. The maximum optical thickness of the cloud is
$\tau_{\rm max}=3.$

An example of the line profile for the particular case of inhomogeneous
obscuration is shown in Fig. 4. The profile takes a complex jagged shape,
approximately conserving its equivalent width. This is precisely the kind of
distortion, that we observe in CQ Tau. The most expressive example is the
spectrum taken on 2012 December 25.

For various lines, the obscuring cloud has the same shape; this explains why
all lines show the same distortions. Some wavelength dependence should arise
due to selective extinction in the clouds, i.e. the cloud’s optical thickness
depends on the wavelength. Line profiles for the same cloud, but for various
maximum optical thicknesses $\tau_{\rm max},$ are shown in Fig. 5.

Another mechanism, which makes the distortions different, is variations in the
limb darkening $I(\mu).$ In the LTE approximation, $I(\mu)$ reflects the
thermal structure of the atmosphere, and generally speaking, can be different
for different lines. One can expect that lines formed in the upper atmosphere
will have somewhat different distortions than the weak lines formed in deeper
layers (see Fig. 6). Departures from LTE can introduce more diversity into the
limb darkening, and enhance individual differences in the line distortions; in
some cases, non-LTE effects can also lead to appearance of the emission lines
(Alexeeva et al., 2016; Sitnova et al., 2018; Mashonkina, 2020). To illustrate
how the diversity in the line limb-darkening can impact the distortions, we
considered a simplest case of linear limb-darkening
$1-\varepsilon+\varepsilon\mu$ with $\varepsilon=0,1,-1;$ the corresponding
profiles are shown in Fig. 7.

Figure 5: The dependence of the distortions on the maximum optical thickness
of the cloud shown in Fig. 4. The dashed line is for $\tau_{\rm max}=0$, the
black solid line is for $\tau_{\rm max}=5,$ the grey lines are for the
intermediate cases $\tau_{\rm max}=0.5,1,1.5,2,2.5,3,4.$ Figure 6:  The
dependence of the distortions on line strength. All else being equal, the
oscillator strengths for the red, black, and blue lines are in the ratios of
10:1:0.1. The cloud shape and parameters are the same as in Fig. 4. Profiles
are scaled to the equal central depth in the out-of-eclipse state.
Corresponding limb-darkening laws in the line centre are shown in the bottom
panel. The thick black line is for the continuum limb-darkening. Figure 7: The
dependence of the distortions on the limb-darkening law. The cloud shape and
parameters are the same as in Fig. 4. Corresponding limb-darkening laws in the
line centre are shown in the bottom panel. The thick black line is for the
continuum limb-darkening, which was the same in all cases. Thin black line is
for the modelled $I(\mu)$; red, green, and blue lines correspond to linear
limb-darkening with $\varepsilon=0,1,-1.$

##  6 Mapping of the clouds

The observed distortions of the line profiles allow us to study the
circumstellar clouds by using principles of the well-known Doppler imaging
technique:

1) the stellar rotation allows us to study the clouds optical thickness along
the $x$-axis (orthogonal to the projection of the rotation axis):
$\tau=\tau(x).$

2) variations in the limb darkening of various lines allow us to study the
structure $\tau=\tau(r)$ along the radius of the stellar disc:
$r^{2}=x^{2}+y^{2}.$

As well as in the case of the Doppler imaging, a series of spectral
observations is needed at moments when the same cloud covers different parts
of the stellar disc. Since we do not have such observations, we will restrict
ourselves to a simpler approach which, nevertheless, can give some
quantitative information about the degree of inhomogeneity of the dust clouds.

###  6.1 The 1D mapping

Each radial velocity in the rotational profile corresponds to a set of points
on the stellar surface, which contribute to absorption in the particular
spectral line. For a rigid stellar rotation, such points form straight lines
parallel to the projection of the rotation axis ($y$). We will call these the
lines of equal radial velocity (LERV). Due to a limited resolution of the
observed spectra and an intrinsic line width, the LERVs can be represented as
finite width strips (SERV), in limits of which variations of the radial
velocity are less than the width of the instrumental/intrinsic profile.

We will describe obscuration of each SERV with a mean optical thickness
$\overline{\tau}$, i.e. we assume that within each SERV, the obscuration is
homogeneous. In this case, the stellar spectrum $f_{j}$, normalized to the
continuum, can be calculated as:

$f_{j}=\frac{\sum I_{ij}s_{i}e^{-\overline{\tau}_{i}a_{j}}}{\sum I^{\rm
c}_{ij}s_{i}e^{-\overline{\tau}_{i}a_{j}}},$ (2)

where $I_{ij}$ is an average intensity at wavelength $\lambda_{j}$ of the
$i$-th SERV with an area $s_{i}$, covered by a cloud with
$\overline{\tau}_{i}$ at some reference wavelength (conventionally at the
$V$-band); $a_{j}$ is an extinction law, normalized to unity at the reference
wavelength; $I^{\rm c}_{ij}$ is an average continuum intensity. The intensity
averaged over the SERV is:

$I_{ij}=\frac{\int I[(\lambda_{j}(1+V/c),\mu]ds_{i}}{s_{i}}.$ (3)

Here we account for the Doppler shift due to stellar rotation, when each point
has a radial velocity $V.$ If we have at least one spectral line profile, i.e.
if we know the values of $f_{j},$ we can solve equation 2 with respect to
$\overline{\tau}_{i}.$

As we have mentioned above for equation 1, equation 2 also can be rewritten in
a linear form:

$f_{j}=\sum u_{ij}w_{i},$ (4)

with unknown weighting coefficients $w_{i},$ related to $\tau_{i}$ as:

$w_{i}=\frac{I^{\rm c}_{i}s_{i}e^{-\overline{\tau}_{i}a}}{\sum I^{\rm
c}_{i}s_{i}e^{-\overline{\tau}_{i}a}}.$ (5)

Here we assume that, in the narrow wavelength range, the extinction law
$a_{j}$ and continuum intensity $I^{\rm c}_{ij}$ do not depend on the
wavelength $\lambda_{j}$. Each row of the matrix $u_{ij}=I_{ij}/I^{\rm
c}_{ij}$ is linearly independent on the others due to Doppler shifts,
therefore equation 4 has a unique solution in terms of the weights $w_{i}$.
However, we will solve the problem in the form of equation 2 because, in the
form of equation 4, a negative $w_{i}$ can appear due to noise in the
observations. Unlike $w,$ the resulting solution in terms of $\tau$ is
ambiguous: if $\overline{\tau}(x)$ is a solution, then
$\overline{\tau}(x)+const$ is also a solution, because the numerator and
denominator in equation 2 (as well as in equation 5) can be multiplied by an
arbitrary factor without affecting the result. This reflects the simple fact
that a uniform extinction across the disc does not change the shape of the
profiles. To eliminate this ambiguity, we will subtract from the obtained
solution the smallest optical thickness which this solution contains. Equation
2 is solved by multidimensional optimization with Tikhonov’s regularization,
which also looks more natural expressed in $\tau,$ rather than in $w.$ The
regularization parameter was chosen so that the resulting smoothing of the
spectral line profiles would not noticeably exceed the spectral resolution.

Figure 8: The line profiles (in the middle row) for various geometries of the
obscuring cloud (the upper row), and results of the one-dimensional mapping
(the bottom row). The simulations were performed for several values of
$\tau_{0}$ ($\tau_{\rm max}$ in the case e): 1 (red), 2 (green), 3 (blue), 5
(cyan), and 100 (magenta, for cases a–d).

To test the method described above, we first considered the simplest case, in
which a homogeneous dust strip with $\tau_{0}$ covers the stellar disc along
the $y$-axis, as shown in Fig. 8a. In this case, a bump appears in the line
profiles, the height of which depends on $\tau_{0}.$ Solving the inverse
problem allows us to recover the value of $\tau_{0}$ and the position of the
strip on the stellar disc. In the case of partial SERV obscuration, the
recovered value of $\overline{\tau}$ is smaller than $\tau_{0}$ because the
relation $e^{-\overline{\tau}}=1-f+fe^{-\tau_{0}},$ where $f$ is the obscured
fraction of the SERV, must hold (see Fig. 8b). The distortions will be less
pronounced for clouds aligned along the $x$-axis (see Fig. 8c–d). Extinction
in the strip above a certain limit, which depends on cloud geometry, does not
affect the spectrum because the contribution from the eclipsed area becomes
negligible so its changes are insignificant. For this reason, strips with
$\tau_{0}=5$ and $\tau_{0}=100$ lead to the same results in Fig. 8. Finally,
for a cloud of complex shape (see Fig. 8e), it is also possible to estimate
the characteristic gradient of the optical thickness. The maximum difference
of the optical thicknesses in the cloud shown in Fig. 8e is $\tau_{\rm max},$
while the restored value is $\approx 0.8\tau_{\rm max}.$

The above examples demonstrate that the reconstructed distributions
$\overline{\tau}(x)$ are smoother than the original distributions $\tau$ in
the plane of the sky, but allow us to adequately estimate the extinction
gradient in the direction perpendicular to the rotation axis’s projection onto
the plane of the sky.

###  6.2 Impact of inaccuracies in line and stellar parameters

Reconstructing the optical thickness across the cloud requires a good model
which describes the observed spectrum in the out-of-eclipse state. We require
precise values for $V_{\rm r},$ $v\sin i,$ atmospheric parameters, elemental
abundances, and line-broadening parameters, in addition to the observed
spectrum, accurately normalized to the continuum level. Inaccuracies in these
parameters can lead to artefacts in the solution $\overline{\tau}(x).$ To give
an impression about arising artificial features, we simulated an ‘observed’
spectrum, using the same cloud as in Fig. 4 with $\tau_{\rm max}=2,$ and then
reconstructed $\overline{\tau}(x),$ varying one of the parameters (see Fig.
9).

If we set $v\sin i$ a little larger than the true value, then an abnormal
extinction arises at the edges of $\overline{\tau}(x),$ which compensates for
a too broad line profile assumed in the modelling. The opposite situation,
when $v\sin i$ is set to below its true value, leads to a strong drop of
$\overline{\tau}(x)$ at the edges, and to a poorer fit. A wrong value of
$V_{\rm r}$ produces artificial features at the edges of $\overline{\tau}(x)$
and leads to a poorer fit of one wing. The combined effect of inaccuracies in
both $V_{\rm r}$ and $v\sin i$ can lead to asymmetric features at the edges of
$\overline{\tau}(x).$

Figure 9: The impact of inaccuracies on the results of 1D mapping: the
reconstructed line profile (upper row) and the one-dimensional map
$\overline{\tau}(x)$ (bottom row). The black lines are for the mapping results
with the exact parameters and given for reference. The red and blue lines
correspond to the results with an inaccuracy introduced into one parameter at
a times: in the rotational velocity $\frac{\Delta v\sin i}{v\sin i}=\pm 0.10$
(a), in the radial velocity $\frac{\Delta V_{\rm r}}{V_{\rm r}}=\pm 0.10$ (b),
in the continuum level $\pm 0.01$ (c), and in the line depth $\frac{\Delta
d}{d}=\pm 0.1$ (d).

It can be seen from Fig. 9 that inaccuracies in $V_{\rm r}$ and $v\sin i$
almost do not change the extinction map in its central part, but affect only
the edges, since the differences in the profiles with such inaccuracies are
most pronounced on the steep wings of the line profile. A similar behaviour is
observed, when we disturb the continuum normalization, but such disturbances
either are clearly seen, when we compare the spectra, or do not produce strong
artefacts in the map.

Artefacts due to inaccuracies in the line broadening parameters should be
similar to those due to inaccuracies in $v\sin i$ and the continuum level.
Uncertainties in the atmospheric parameters and the elemental abundances can
be reduced to inaccuracies in the line depths and corresponding artifacts are
also significant only at edges of the map $\overline{\tau}(x)$ – see Fig. 9.

One can conclude that the indirect signatures of significant uncertainties in
the parameters are the steep change $\overline{\tau}(x)$ at the edges and
poorly fitted line profiles. These signatures should be stable over a series
of observations, if inaccuracies in the parameters are consistent throughout.
The most reliable information on the extinction gradients is contained in the
central part of the profile, which is less sensitive to inaccuracies in the
line and stellar parameters.

###  6.3 The 1D mapping on the real spectra

Figure 10: The observed spectra (black lines) of CQ Tau for December 22 (the
upper parts) and December 25 (the lower parts). The dashed blue line is the
out-of-eclipse synthetic spectrum of CQ Tau. The red line is the synthetic
spectra during the eclipses with the recovered 1D distributions
$\overline{\tau}(x),$ which are shown in Fig. 11. The shaded area is a region
which was used to reconstruct $\overline{\tau}(x).$ Figure 10: Continued
Figure 11: The 1D map $\overline{\tau}(x)$ reconstructed over the spectral
region 5400 – 5418Å (the shaded area in Fig. 10). The dashed lines show the
map distorted by varying $v\sin i$ (the upper panel) or $V_{\rm r}$ (the
middle panel) by $+10$ km s-1 (the blue line) and $-10$ km s-1 (the red line).
The mean value $\langle\overline{\tau}\rangle$ was removed for ease of
comparison. The bottom panels are for the line profiles, the dashed lines are
for the out-of-eclipse spectrum, the solid lines correspond to the obscured
star with $\overline{\tau}(x),$ shown by black lines in the panels above.

We apply the 1D mapping method described above to the archival spectra of CQ
Tau. To simulate instrumental broadening, the theoretical spectra were
additionally smoothed with a Voigt profile, which reproduces the far wings of
the observed profiles better than a simple Gaussian broadening. We show, in
Fig. 10, 6 spectral regions from blue to red for the spectra of December 22,
when the lines looked less distorted, and also of December 25, when the
strongest distortions were observed. Equation 2 can be solved if we have at
least one profile that allows us to use a small spectral region to reconstruct
$\overline{\tau}(x),$ while the remaining spectrum may be used to validate our
explanation of the distortions. We chose a spectral region with a few Fe i
lines (shadowed in Fig. 10). The corresponding solutions $\overline{\tau}(x)$
are shown in Fig. 11 for both nights. To calculate the synthetic spectra over
a wide spectral range, we assume the standard extinction law (Cardelli et al.,
1989). For December 22, the amplitude of variations in $\overline{\tau}$ is
smaller, and significant absorption appears only at edges, and therefore
depends on the assumed values for $v\sin i$ and $V_{\rm r}$ (see Section 6.2).
Since we do not have out-of-eclipse spectra, the applied values of $v\sin i$
and $V_{\rm r}$ can be uncertain. To check how these uncertainties can impact
our results, we vary $v\sin i$ and $V_{\rm r}$ by $\pm 10$ km s${}^{-1};$ the
corresponding $\overline{\tau}(x)$ are shown with dashed lines. Unlike on
December 22, absorption peaks in $\overline{\tau}$ for December 25 appear in
the central part. In this case, inaccuracies in $v\sin i$ and $V_{\rm r}$ are
less important, which allows us to estimate the characteristic gradient and
size of eclipsing clouds even without exact values of $v\sin i$ and $V_{\rm
r}.$

We see three peaks in the map at $x/R_{*}\approx-0.70,$ $-0.30,$ 0.27. The
first peak ($-0.70$) makes the left wing of the profile less steep than the
right one: see the line profile in the bottom panel of Fig. 11. There are also
other mechanisms which can produce asymmetric wings, for example, scattering
on the moving dust clouds. This peak is also sensitive to the assumed radial
velocity of the star, so we cannot be sure of its parameters. However, the
peaks at $x\approx-0.30R_{*}$ and $0.27R_{*}$ are quite reliable, and are
responsible for the appearance of the jagged profiles. The most pronounced
bump in the profile is produced by the peak at $0.27R_{*},$ which has an
amplitude $\Delta A_{V}\sim\Delta\overline{\tau}\sim 1$ and a linear size
$\approx 0.6R_{*}.$ It means the existence of spatially small but optically
thick dusty clouds on the line of sight. Since the 1D mapping underestimates
the extinction gradients, the real inhomogeneity could be much larger.

Because any homogeneous extinction can be added to the obtained solution
$\overline{\tau},$ we can estimate only a lower limit on a brightness drop
$\Delta V,$ which is 0.4 and 1.0 magnitude for December 22 and 25,
correspondingly. The observed drop is $1.5-2$ magnitudes, therefore a residual
$A_{V}\approx 1$ is related to a flat part of the circumstellar extinction.

Despite the fact that only three lines were used to reconstruct the map, it
reproduces the distortions over a broad spectral region: narrow features in
the line profiles appear at the very places where our model predicts. This
confirms our interpretation of the observed distortions as a result of the
inhomogeneous obscuration of the rotating stellar surface. Visible
discrepancies between the theoretical and observed spectrum are mostly in the
equivalent widths and may be due to uncertainties in the atmospheric
parameters, non-LTE effects, or errors in the line parameters.

##  7 Discussion

CQ Tau shows accretion signatures (inverse P-Cygni profiles of He i 5876Å and
O i 7773Å), and one can suppose that the variability observed in the profiles
is related to the veiling by lines (Dodin & Lamzin, 2012) when the stellar
absorption lines are filled-in by weak emission lines, which can appear in
various parts of the profiles depending on the geometry and orientation of the
accretion hotspot. However, these emissions should reduce the total equivalent
widths of the lines, which is not observed. Moreover, observations as well as
simulations show that these emissions are stronger for stronger lines, which
is also not observed in our case. Therefore, we conclude that the accretion
region is not responsible for the line distortions; moreover, the absence of a
noticeable veiling indicates that the accretion contribution is negligible and
visible only in the strongest lines: H $\alpha$, He i 5876 Å and O i 7773 Å.
In principle, the accretion process can produce an extra UV radiation field,
which can cause Fe overionization and hence enhance Fe ii lines (see Section
3.1).

Cold spots on the stellar surface can also distort the line profiles in a
manner similar to what has been observed. However, CQ Tau does not have a
convective envelope (see fig. 1 in Villebrun et al. 2019), and cold spots on
its surface are not expected. Moreover, such spots should produce a noticeable
regular photometric variability, which has not been detected.

One can conjecture that the distortions arise due to some circumstellar gas,
in which case the equivalent widths should increase; this also is not
observed. In addition, the distortions of the profiles would appear first of
all in the strongest lines: in the case of a cold gas, in Na i, K i lines, or
else in H i, He i lines, in the case of hotter gas or accretion flows (Dodin,
2018). These absorption features are indeed present in the spectrum of CQ Tau,
but they have a completely different shape which does not resemble the
distortions in all of the other absorption lines.

Because the distortions do not depend on the physical parameters of the lines,
it is natural to suppose that they are induced by geometrical effects: both or
either of inhomogeneous obscuration, and scattering by moving dust clouds
(Grinin et al., 2012).

Light scattering on moving clouds can modify line profiles; these changes
depend on the geometry and kinematics of the scattering bodies. Such
scattering typically leads to profile smoothing, rather than to the appearance
of sharp, narrow features. One can suppose a very special geometry, in which
scattering clouds ‘see’ the star pole-on, while the observer looks in the
equatorial plane. In this case, the lines in the scattered spectrum will be
broadened less by the stellar rotation than they are in direct light. However,
such an explanation looks too artificial because it requires: 1) a
simultaneous presence of a few clouds, moving with various velocities, but
just within $v\sin i$; and 2) the contribution of scattered light from these
clouds must be comparable to the contribution of direct light from the star,
in order for the distortions to be noticeable. Simultaneously satisfying both
of these requirements seems unlikely, however, noticeable contribution of
scattered light must be present because CQ Tau’s $V$-band magnitude was
$10.6-11.$ This magnitude corresponds to the colour turn on the
colour–magnitude diagram (Shakhovskoj et al., 2005), meaning that about half
of the flux is provided by scattered light. The peak at $x/R_{*}=-0.7$ in
$\overline{\tau}(x)$ can probably be explained by scattered light, however, to
do this, special modelling is needed, which is outside the scope of this
study. Nevertheless, we can understand qualitatively how scattering can change
our results. Since the profile is smoothed by the scattering on an extended
environment (a rotating disc or an expanding envelope), a spectrum with
smoother profiles is added to the jagged spectrum of the direct stellar light,
reducing the amplitude of the distortions. This reduction leads to a
corresponding decrease in the extinction gradients. Thus, the map
$\overline{\tau}$ obtained without taking scattered light into account gives a
lower limit on the extinction gradients, while accounting for the scattering
should lead to a sharper map.

Since we have not detected variability over 2 hours of observations, the
cloud’s crossing time $t_{\rm c}\gtrsim 0.5$ day, but the spectra on December
22 and 25 are completely different so, therefore, $t_{\rm c}\lesssim 3$ days.
For the stellar parameters from Villebrun et al. (2019), we can estimate that,
for Keplerian velocities, such crossing times correspond to distances between
$0.9-33$ a.u., i.e. well outside the sublimation radius ($\sim 0.4$ a.u.).

##  8 Conclusions

Using CQ Tau as an example, we have shown that the circumstellar dust clouds
are inhomogeneous on scales of the stellar radius that leads to variable
distortions in the line profiles. A distinctive feature of the line
distortions, caused by the inhomogeneous obscuration, is the approximate
conservation of the equivalent widths. In the case of another known source of
profile variability in young stars – the accretion hotspots – the line
equivalent width is reduced, moreover this effect is stronger for stronger
lines (Dodin & Lamzin, 2012; Dodin, 2018; Rei et al., 2018). Cold spots on the
stellar surface also almost do not change the equivalent width because the
absolutely dark spot is equivalent to an absolutely opaque cloud, however, the
time behaviour of the distortions should differ between these two cases.

The spectral variability of CQ Tau shows us that the star cannot be considered
as a point source when we deal with circumstellar dust. This fact
significantly complicates the interpretation of various kinds of observations:
the line distortions prevent an accurate determination of the stellar
parameters, first of all, $V_{\rm r}$ and $v\sin i$; an inhomogeneous
obscuration makes the apparent extinction law greyer than the intrinsic one
(Natta & Panagia, 1984); an inhomogeneous obscuration can change the relative
contribution of hot/cold spots on the stellar surface, for example, it can
modify the apparent accretion indicators, if they are measured relative to the
stellar spectrum (e.g., equivalent width, veiling), since the obscuration of
the accretion region and the stellar surface can be different and variable.
Perhaps this explains why there is no correlation between accretion tracers
and stellar brightness in RW Aur A (see fig. 11 in Petrov et al. 2001).

On the other hand, the line distortions offer the possibility of studying the
dust clouds on scales inaccessible by direct observation. The estimates of the
size and amplitude of dust clumps, given above, were obtained from merely one
archival spectrum. We encourage special spectral observations, which should
cover various phases of a dust clump transit over the stellar disc. Such
observations will allow the velocity of the clouds to be measured relative to
the rotational axis of the star, and to clarify what obscures the star:
whether a clumped dusty wind or a warped disc.

The clumpy structure of the dust environment probably also takes place for
other types of stars. However, in many cases, it will be difficult to study by
the method proposed here, due to a slower rotation and other sources of
spectral and photometric variability (such as accretion spots or cold spots in
the case of T Tauri stars).

## Acknowledgements

This work is supported by the Russian Science Foundation under grant
20-72-10011.

## Data availability

The data underlying this article were accessed from Canadian Astronomy Data
Centre http://www.cadc-ccda.hia-iha.nrc-cnrc.gc.ca. The derived data generated
in this research will be shared on reasonable request to the corresponding
author.

## References

  * Alecian et al. (2013) Alecian E., et al., 2013, MNRAS, 429, 1001
  * Alexeeva et al. (2016) Alexeeva S. A., Ryabchikova T. A., Mashonkina L. I., 2016, MNRAS, 462, 1123
  * Barsunova et al. (2015) Barsunova O. Y., Grinin V. P., Sergeev S. G., Semenov A. O., Shugarov S. Y., 2015, Astrophysics, 58, 193
  * Bergemann et al. (2012) Bergemann M., Lind K., Collet R., Magic Z., Asplund M., 2012, MNRAS, 427, 27
  * Bibo et al. (1992) Bibo E. A., The P. S., Dawanas D. N., 1992, A&A, 260, 293
  * Bouvier et al. (1999) Bouvier J., et al., 1999, A&A, 349, 619
  * Cardelli et al. (1989) Cardelli J. A., Clayton G. C., Mathis J. S., 1989, ApJ, 345, 245
  * Cody et al. (2014) Cody A. M., et al., 2014, AJ, 147, 82
  * Cohen et al. (2003) Cohen R. E., Herbst W., Williams E. C., 2003, ApJ, 596, L243
  * Dodin (2018) Dodin A., 2018, MNRAS, 475, 4367
  * Dodin & Lamzin (2012) Dodin A. V., Lamzin S. A., 2012, Astronomy Letters, 38, 649
  * Dodin et al. (2019) Dodin A., et al., 2019, MNRAS, 482, 5524
  * Friedemann et al. (1992) Friedemann C., Reimann H. G., Guertler J., 1992, A&A, 255, 246
  * Gahm et al. (2018) Gahm G. F., Petrov P. P., Tambovsteva L. V., Grinin V. P., Stempels H. C., Walter F. M., 2018, A&A, 614, A117
  * Grinin & Potravnov (2013) Grinin V. P., Potravnov I. S., 2013, Astrophysics, 56, 1
  * Grinin et al. (2008) Grinin V. P., Barsunova O. Y., Shugarov S. Y., Kroll P., Sergeev S. G., 2008, Astrophysics, 51, 1
  * Grinin et al. (2012) Grinin V. P., Tambovtseva L. V., Weigelt G., 2012, A&A, 544, A45
  * Grinin et al. (2019) Grinin V. P., Semenov A. O., Barsunova O. Y., Sergeev S. G., 2019, Astrophysics, 62, 41
  * Günther et al. (2018) Günther H. M., et al., 2018, AJ, 156, 56
  * Gustafsson et al. (2008) Gustafsson B., Edvardsson B., Eriksson K., Jørgensen U. G., Nordlund Å., Plez B., 2008, A&A, 486, 951
  * Hernández et al. (2004) Hernández J., Calvet N., Briceño C., Hartmann L., Berlind P., 2004, AJ, 127, 1682
  * Kausch et al. (2015) Kausch W., et al., 2015, A&A, 576, A78
  * Manset et al. (2009) Manset N., Bastien P., Ménard F., Bertout C., Le van Suu A., Boivin L., 2009, A&A, 499, 137
  * Mashonkina (2020) Mashonkina L., 2020, MNRAS, 493, 6095
  * McGinnis et al. (2015) McGinnis P. T., et al., 2015, A&A, 577, A11
  * Natta & Panagia (1984) Natta A., Panagia N., 1984, ApJ, 287, 228
  * Petrov & Kozack (2007) Petrov P. P., Kozack B. S., 2007, Astronomy Reports, 51, 500
  * Petrov et al. (2001) Petrov P. P., Gahm G. F., Gameiro J. F., Duemmler R., Ilyin I. V., Laakkonen T., Lago M. T. V. T., Tuominen I., 2001, A&A, 369, 993
  * Petrov et al. (2015) Petrov P. P., Gahm G. F., Djupvik A. A., Babina E. V., Artemenko S. A., Grankin K. N., 2015, A&A, 577, A73
  * Piskunov & Valenti (2016) Piskunov N., Valenti J. A., 2016, Astronomy & Astrophysics, 597, A16
  * Rei et al. (2018) Rei A. C. S., Petrov P. P., Gameiro J. F., 2018, A&A, 610, A40
  * Rodriguez et al. (2016) Rodriguez J. E., et al., 2016, ApJ, 831, 74
  * Ryabchikova et al. (2015) Ryabchikova T., Piskunov N., Kurucz R. L., Stempels H. C., Heiter U., Pakhomov Y., Barklem P. S., 2015, Phys. Scr., 90, 054005
  * Schmidt (2019) Schmidt E. G., 2019, ApJ, 880, L7
  * Schneider et al. (2015a) Schneider P. C., et al., 2015a, A&A, 584, L9
  * Schneider et al. (2015b) Schneider P. C., France K., Günther H. M., Herczeg G., Robrade J., Bouvier J., McJunkin M., Schmitt J. H. M. M., 2015b, A&A, 584, A51
  * Shakhovskoj et al. (2005) Shakhovskoj D. N., Grinin V. P., Rostopchina A. N., 2005, Astrophysics, 48, 135
  * Sitnova et al. (2018) Sitnova T. M., Mashonkina L. I., Ryabchikova T. A., 2018, MNRAS, 477, 3343
  * Siwak et al. (2018) Siwak M., et al., 2018, MNRAS, 478, 758
  * Steenman & The (1989) Steenman H., The P. S., 1989, Ap&SS, 159, 189
  * Steenman & The (1991) Steenman H., The P. S., 1991, Ap&SS, 184, 9
  * Villebrun et al. (2019) Villebrun F., et al., 2019, A&A, 622, A72
  * Waters & Waelkens (1998) Waters L. B. F. M., Waelkens C., 1998, ARA&A, 36, 233