# The periodic spectroscopic variability of FU Orionis

Stacie L. Powell, Mike Irwin, Jerome Bouvier, Cathie J. Clarke

► To cite this version:

Stacie L. Powell, Mike Irwin, Jerome Bouvier, Cathie J. Clarke. The periodic spectroscopic variability of FU Orionis. Monthly Notices of the Royal Astronomical Society, 2012, 426, pp.3315-3333. 10.1111/j.1365-2966.2012.21898.x. insu-03612395

HAL Id: insu-03612395

https://insu.hal.science/insu-03612395

Submitted on 17 Mar 2022

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 Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# The periodic spectroscopic variability of FU Orionis*

Stacie L. Powell,<fnref n="1" />† Mike Irwin,<fnref n="1" /> Jerome Bouvier<fnref n="2" /> and Cathie J. Clarke<fnref n="1" />

Accepted 2012 August 8. Received 2012 July 16; in original form 2011 November 22

## ABSTRACT

FU Orionis systems are young stars undergoing outbursts of disc accretion and where the optical spectrum contains lines associated with both the disc photosphere and a wind component. Previous observations of the prototype FU Orionis by Herbig et al., suggested that the wind lines and the photospheric lines are modulated with periods of 14.54 and 3.54 d, respectively. We have re-observed the system at higher spectral resolution, by monitoring variations of optical line profiles over 21 nights in 2007 and have found periods of 13.48 and 3.6 d in the wind and disc components, consistent with the above: this implies variability mechanisms that are stable over at least a decade. In addition, we have found (i) that the variations in the photospheric absorption lines are confined to the blue wing of the line (centred on velocity ~−9 km s⁻<fnref n="1" />): we tentatively ascribe this to an orbiting hotspot in the disc, which is obscured by a disc warp during its receding phase. (ii) The wind period is manifested not only in blueshifted Hα absorption (as found by Herbig et al.), but also in redshifted emission of Hα and Hβ, as well as in blueshifted absorption of NaI D, LiI and FeII λ5018. (iii) We find that the periodic modulation of blueshifted Hα absorption (at a velocity of around −100 km s⁻<fnref n="1" />) is phase lagged with respect to variations in the other lines by around 1.8 d. This is consistent with a picture in which variations at the wind base first affect chromospheric emission and then low-velocity blueshifted absorption, followed – after a lag equal to the propagation time of disturbances across the wind’s acceleration region – by a response in high-velocity blueshifted absorption. Such arguments constrain the size of the acceleration region to be ~10¹<fnref n="2" /> cm. We discuss possible mechanisms for periodic variations within the innermost 0.1 au of the disc, including the possibility that these variations indicate the presence of an embedded hot Jupiter in this object.

**Key words:** accretion, accretion discs – line: profiles – stars: individual: FU Orionis – stars: mass-loss – stars: pre-main-sequence – stars: winds, outflows.

## 1 INTRODUCTION

Although FU Orionis systems only represent a small class of young stellar objects, the benefits of understanding the physical processes present in these systems have implications for the whole of not only star, but also planet formation. FU Orionis systems are unique objects which provide the only opportunity to directly observe an accretion disc around a young star. Originally identified as a subset of pre-main-sequence stars entitled 'FUors' by Ambartsumian (1971) and Ambartsumyan (1971), they are characterized by a large increase in optical brightness of ~4 mag or more over a period of 1–10 yr, with a slow decline of 50–100 yr (Herbig 1966, 1977). Other defining properties of the class include spectral types of late

F to G in the optical, with infrared spectra showing K-M giant-supergiant type features (Mould et al. 1978; Hartmann & Kenyon 1985, 1987a,b; Kenyon, Hartmann & Hewett 1988; Hartmann & Calvet 1995). FU Orionis systems are associated with young T-Tauri stars, evident from the strong LiIλ6707 absorption observed by Herbig (1966), a characteristic of young stars. The discovery of the second class member, V1057 Cyg, for which a pre-outburst spectrum resembled that of a T-Tauri star (Herbig 1958), supported this association with young stars. All FU Orionis objects are associated with star-forming regions, as is evident from the ring-shaped reflection nebulae first noted by Goodrich (1987), from the presence of infrared excess due to circumstellar dust and by their general location in heavily extincted regions (Herbig 1977). More recent developments have suggested these objects tend to have binary companions (see table 1 in Vittone & Errico 2005), although the binary fraction is not necessarily higher than in other young stars. A companion candidate to FU Orionis itself was detected by Wang

<footnote n="1">¹Institute of Astronomy, Madingley Road, Cambridge CB3 OHA</footnote>

<footnote n="2">²UJF-Grenoble I/CNRS-INSU, Institut de Planétologie et d'Astrophysique de Grenoble (IPAG) UMR 5274, Grenoble F-38041, France</footnote>

<footnote>*Based on observations made at Observatoire de Haute-Provence (CNRS), France, with the SOPHIE spectrograph.</footnote>

<footnote>†E-mail: slp65@ast.cam.ac.uk</footnote>
et al. (2004) at a separation of 0.5 arcsec (see Hartmann, Kenyon & Hartigan 1993; Kenyon 1995; Hartmann & Kenyon 1996; Hartmann 1998 for reviews).

FU Orionis outbursts are thought to be produced by long periods of rapid accretion (Larson 1983; Lin & Papaloizou 1985), when the accretion rate is so great, ~$10^{-4} M_{\odot} yr^{-1}$ (Kenyon et al. 1988), the inner disc heats up through viscous dissipation and the mid-plane temperature exceeds that of the photosphere. Under these circumstances, the disc produces an absorption line spectrum with lower surface gravity than the central star, and the star's luminosity is swamped by factors of 100–1000 (Hartmann & Kenyon 1996). Consequently, Hartmann & Kenyon (1985, 1987a), Kenyon et al. (1988), Welty et al. (1992), Popham et al. (1996) and Zhu et al. (2009a) showed that the spectrum of FU Orionis is well fitted by the predictions of a broad-band steady state accretion disc model. The periods of rapid accretion are thought to be triggered by thermal instabilities in the disc associated with partial ionization of disc material, which causes the release of dammed up material on to the central star (Hartmann & Kenyon 1985). The trigger mechanisms are discussed by Pringle (1976), Hartmann & Kenyon (1985), Lin & Papaloizou (1985), Clarke, Lin & Pringle (1990), Bell et al. (1995), Clarke & Syer (1996), Lodato & Clarke (2004), Zhu et al. (2009b) and references therein. The frequency of FU Orionis events observed within the stellar neighbourhood has led to the belief that many low-mass stars may undergo multiple FU Orionis events (Herbig 1977; Hartmann & Kenyon 1985). It is currently thought that up to 5–10 per cent of the matter accreted on to young T-Tauri stars occurs through repetitive periods of rapid accretion in FU Orionis outbursts (Hartmann & Kenyon 1987a; Hartmann et al. 1993). This highlights the importance of understanding this episodic accretion process in the context of stellar formation.

The high accretion rate in FU Orionis objects is accompanied by high mass-loss rates with most objects associated with molecular outflows, optical jets and Herbig-Haro objects (Reipurth 1990; Evans et al. 1994). The mass-loss rate is around ~$10^{-5} M_{\odot} yr^{-1}$ in FU Orionis itself (Croswell, Hartmann & Avrett 1987; Calvet, Hartmann & Kenyon 1993; Hartmann & Calvet 1995). The strong bipolar outflow has a terminal wind velocity of ~300 km s⁻<fnref n="1" />, which has been detected previously in the deep, broad blueshifted components of Hα and Na I D (Herbig 1977; Bastian & Mundt 1985). It is thought that the outflow arises from the surface of a Keplerian disc (Hartmann & Calvet 1995), with the possibility of a disc 'chromosphere' first suggested by Croswell et al. (1987) to explain the detection of Hα emission, and modelled in the line profiles by D'Angelo et al. (2000). The mechanism producing the rapid acceleration of material into these powerful outflows in circumstellar discs is not entirely clear. Possible models involve the accretion flow spinning up the star close to break-up and launching magneto-centrifugal winds along opened up field lines anchored in the star (Shu et al. 1994), or conical winds driven by the pressure gradient of the azimuthal magnetic field component and launched along opened up stellar field lines (Königl, Romanova & Lovelace 2011). However, magnetic fields anchored in the rotating disc itself could also centrifugally accelerate material outwards as described by Blandford & Payne (1982). The detection of a ~1 kG magnetic field directed towards the observer at 0.05 au in FU Orionis by Donati et al. (2005) supports the involvement of magnetic fields, yet their topology and role in driving mass loss remains unclear.

In addition to the variability on long time-scales, it appears that FU Ori also exhibits unexplained variability on shorter time-scales. Variability in Hα emission was first noted by Herbig (1966) over a time-scale of one year, which was confirmed in absorption and

emission by Herbig (1977). This variability was observed to exist on a monthly time-scale in observations of Hα and Na I D by Bastian & Mundt (1985), who concluded that the mass outflow activity from FU Orionis objects was a long-lasting active state. This was confirmed by Croswell et al. (1987) who went on to suggest that the maximum size of the region producing the variation observed in the absorption profiles of Hα at ~300 km s⁻<fnref n="1" /> was roughly twice the optical photospheric radius, or ~3 × 10¹² cm. Hartmann & Kenyon (1985) noted variability in the blue spectral region and in the Mg I λ5183 line and attributed the cause to variation in the mass ejection, concluding that the wind arises from the optical disc surface in FU Orionis. Random, small-amplitude photometric variations on time-scales of 1 d were observed by Kenyon et al. (2000), and the rapid variations were attributed with the flickering of the inner accretion disc in FU Orionis.

The presence of periodicity was first evident in the brightness fluctuations after the outburst, where Fourier analysis by Chochol & Tremko (1980) detected periods of 1645 and 3760 d from observations between 1936 and 1976. Smaller cyclic variations were then detected from photometry taken in 1984–1985 by Kolotilov & Petrov (1985), which showed a brightness periodicity in the V band of 18.35 d for two months in 1984 and interestingly, a periodicity of approximately half, 9.19 d for the remainder of the observations. The shorter period was confirmed by Ibragimov (1993), who in addition detected a period of 9.2 d in the V-band light curve from 1987 August–1988 January. Kenyon et al. (1988) highlighted variability in FU Orionis below 4000 Å, which may be attributed with the periodicity found by Kolotilov & Petrov (1985), or could be due to fluctuations in airmass. Further confirmations of periodicities were found in the high-resolution spectroscopy (13 km s⁻<fnref n="1" />) of 20 nights observed from 1997 to 1999 analysed by Herbig, Petrov & Duemmler (2003), who not only found variation in Hα on daily time-scales, but also confirmed a period of 13–18 d, which peaked at 14.847 d in the equivalent width (EW) of Hα between −110 and −270 km s⁻<fnref n="1" />. Errico, Vittone & Lamzin (2003) attributed this period to a rotationally modulated wind structure, envisaging a non-axisymmetric magnetic structure corotating with the star. Note that whereas Herbig et al. (2003) found a periodic signal only in the absorption component of Hα, Vittone & Errico (2005) detected periodicity in the Hα emission component of 6.70 d in 1998–2000, but found no correlation and no periodicity present in the absorption component of Hα.

It is not only the strong blueshifted absorption lines affected by mass loss that show periodicity in FU Orionis. Herbig et al. (2003) analysed spectral regions of weak photospheric lines, without strong blueshifted components (5540–5640 and 6320–6440 Å) through cross-correlation with a G0 Ib reference star (β Aqr) over 20 nights between 1997 and 1999. This analysis demonstrated that the radial velocity of the peak of the cross-correlation function (CCF) undergoes periodic variations with a period of 3.542 d.<fnref n="1" /> It was suggested that this periodicity could be attributed to a non-axisymmetric hotspot (Errico et al. 2003; Herbig et al. 2003). In the absence of any mechanism to sustain this hotspot against shear, a hotspot of radial extent δR would be sheared out over R/δR orbits. Thus it is of interest to discover whether this periodicity is sustained over many years. Long-lived periodic behaviour would be consistent, for example, with non-axisymmetric disc structure induced by

<footnote n="1">¹ Variability also appears to be present in the 6170 Å spectral region of FU Orionis, which contains weak absorption lines from the photosphere of the disc, however no time series analysis was conducted (Hartmann & Kenyon 1987a).</footnote>
a stable magnetic configuration or even an embedded planet (Clarke & Armitage 2003).

In this paper we present observations taken over 21 nights from the Spectrographe pour l'Observation des Phénomènes des Intérieurs stellaires et des Exoplanètes (SOPHIE) échelle spectrograph spanning a time frame of 38 nights in 2007. The dense time coverage and the high resolution of 4 km s⁻¹ have allowed the periodicities in FU Orionis to be studied in great detail. We observe the stability of the periodicities due to non-axisymmetric structure detected by Herbig et al. (2003) over a decade and highlight the necessity of a mechanism which maintains such a structure against shear. In Section 2 the details of the observations and data reduction techniques are described, with the variability in the line profiles described in Sections 3–5. The periodicity in the cross-correlations of spectral regions containing weak photospheric lines is examined in Section 6. Possible mechanisms to sustain the periodicities over decades are discussed in Section 7.

## 2 OBSERVATIONS AND DATA REDUCTION

Optical spectra of FU Orionis were obtained from the SOPHIE cross-dispersed échelle spectrograph mounted on the Cassigrain focus of the 1.93-m telescope at the Observatoire de Haute-Provence (OHP, France), covering the wavelength range of 3872–6943 Å. The observations were obtained using the high-resolution mode (R = 75 000), over 21 nights, providing an excellent opportunity to study the optical line variability of this object in detail. A reference spectrum of β Aqr was also taken from SOPHIE in the same mode and used as the template in the cross-correlations. This G0lb supergiant is used as the standard as it has been shown by Hartmann & Kenyon (1985) that G spectral types provide a better signal-to-noise ratio correlation than K spectral type templates, and it has the same effective temperature as the region of the disc of FU Orionis visible at optical wavelengths, as described by Croswell et al. (1987). The data were reduced automatically by the pipeline, which includes bias subtraction, optimal order extraction, cosmic ray removal, corrections for flat fielding and wavelength calibration. The 39 échelle orders were divided by the blaze function, reconnected and barycentrically corrected, resulting in a resampled spectrum with a constant wavelength step of 0.5 km s⁻¹, smaller than the original spectral sampling. Details of the observations are shown in Table 1.

The data were continuum subtracted by clipping out the spectral features and performing median and boxcar filtering over a scale length of 20 Å, to obtain a pure continuum spectrum, which was then subtracted from the featured spectrum and normalized to unit continuum. See Battaglia et al. (2008) for more details on this method. This was performed over the wavelength range of 4100–6865 Å to ensure that variations in the continuum were not included in the resultant line profiles and to minimize the propagation of errors due to edge effects. Fig. 1 shows the initial normalization of the JD = 2454 106.309 observation over this region. Variance in the continuum was minimized by averaging all 21 nights of the normalized spectra and subtracting this from each spectra individually to leave the residual spectra. These residuals were then continuum subtracted again and the result was recombined with the initial average spectrum to ensure that fluctuations in the continuum do not affect the line profile variability.

Spurious emission was removed by flagging the points which were above 3σ from the mean, both in each spectrum individually and across simultaneous points in all 21 observations, and replacing with the result of a median filter across a region of ±0.1 Å around the bad value. This ensured that large variations in emission over

<table><caption>Table 1. Observations of FU Orionis and β Aqr.</caption><thead><tr><th>Object</th><th>Date</th><th>Sequence ID</th><th>Exposure time (s)</th><th>Barycentric JD<br/>-2454 000.</th></tr></thead><tbody><tr><td rowspan="14">FU Orionis</td><td>03-01-2007</td><td>220020</td><td>1800.0</td><td>104.319 450 10</td></tr><tr><td>04-01-2007</td><td>229700</td><td>2000.0</td><td>105.333 957 83</td></tr><tr><td>05-01-2007</td><td>235720</td><td>2500.0</td><td>106.308 742 92</td></tr><tr><td>06-01-2007</td><td>244940</td><td>2200.0</td><td>107.301 338 38</td></tr><tr><td></td><td>255950</td><td>777.0</td><td>107.320 284 52</td></tr><tr><td>08-01-2007</td><td>258300</td><td>2000.0</td><td>109.316 905 41</td></tr><tr><td>12-01-2007</td><td>274990</td><td>2000.0</td><td>113.316 875 68</td></tr><tr><td>13-01-2007</td><td>275710</td><td>2000.0</td><td>114.316 045 53</td></tr><tr><td>14-01-2007</td><td>276540</td><td>1800.0</td><td>115.287 495 19</td></tr><tr><td></td><td>276550</td><td>1800.0</td><td>115.310 445 53</td></tr><tr><td>15-01-2007</td><td>277830</td><td>2000.0</td><td>116.371 773 76</td></tr><tr><td>17-01-2007</td><td>279210</td><td>2000.0</td><td>118.319 514 93</td></tr><tr><td>18-01-2007</td><td>281980</td><td>2000.0</td><td>119.420 234 51</td></tr><tr><td>19-01-2007</td><td>282860</td><td>2000.0</td><td>120.316 946 94</td></tr><tr><td>20-01-2007</td><td>284260</td><td>2000.0</td><td>121.333 350 79</td></tr><tr><td>26-01-2007</td><td>294910</td><td>2000.0</td><td>127.376 862 94</td></tr><tr><td>27-01-2007</td><td>295440</td><td>2000.0</td><td>128.307 139 83</td></tr><tr><td>28-01-2007</td><td>296290</td><td>2000.0</td><td>129.304 726 06</td></tr><tr><td>29-01-2007</td><td>297200</td><td>2000.0</td><td>130.291 814 04</td></tr><tr><td>30-01-2007</td><td>298190</td><td>2000.0</td><td>131.303 598 23</td></tr><tr><td>31-01-2007</td><td>299070</td><td>2000.0</td><td>132.306 793 79</td></tr><tr><td>05-02-2007</td><td>303430</td><td>1800.0</td><td>137.334 767 02</td></tr><tr><td>10-02-2007</td><td>306160</td><td>1800.0</td><td>142.332 842 43</td></tr><tr><td>08-10-2007</td><td>465830</td><td>300.0</td><td>382.347 610 84</td></tr></tbody></table>

the spectral series were removed, reducing the potential contribution from spurious lines. Telluric lines were identified through visual comparison between all the unbarycentrically corrected spectra stacked with the median of all spectra. Aligned, narrow lines which persisted in the median spectrum were clipped and removed from further analysis.

The continuum-subtracted spectra were cross-correlated with a template G0lb reference star, namely β Aqr taken from the same instrument, via the methods described in Tonry & Davis (1979). The template was shifted to zero velocity, as it has a radial velocity measured to be 7.11 km s⁻¹, when cross-correlated with a synthetic spectrum of a G0 star with T<sub>eff</sub> = 6000 K, log(g) = 3.0. This specific model was chosen as it provided the best match to the observed line strengths in the β Aqr spectrum. The radial velocity was derived by fitting a Gaussian to the peak of the symmetric CCF. This was necessary as previous measurements were made with far less accurate resolution by Wilson (1953). The observed template was continuum subtracted in the same manner as described above. Furthermore, the synthetic stellar spectrum was smoothed using a Gaussian with full width at half-maximum (FWHM), comparable to the corresponding unblended, unsaturated absorption lines in the relevant regions of the β Aqr spectrum. The resulting spectrum was used as an addition reference standard in cross-correlations in order to confirm the zero-point of the spectral lines in the β Aqr spectrum. Using this method over a spectral interval not only improves the signal-to-noise ratio, but also allows for the effects of line blending or overlap and reduces scatter in the results (Zhu et al. 2009a), ideal for the weaker metallic lines.

Errors were obtained for the cross-correlations through estimating the variance and hence signal-to-noise ratio in each spectral region, using the areas identified through spectral clipping to be pure continuum. The covariance matrix across three adjacent pixels in these regions were calculated and the noise reduction factor was computed. This was combined with the original estimate for
<figure><img src="image_3.png" /><figcaption>**Figure 1.** High-resolution échelle spectrum taken by SOPHIE, on JD =2454 106.309 mounted on the 1.93-m telescope at the OHP. The top panel shows the original spectrum, with the red dashed line representing the continuum level. The bottom panel shows the spectrum normalized to unit continuum.</figcaption></figure>

the signal-to-noise ratio, which accounted for the noise correlated across nearby pixels in a rebinned échelle spectrograph. Random amounts of Gaussian noise centred on the signal-to-noise ratio value were then added to the spectra, which were then cross-correlated with the template spectrum once more. The resulting variance in the cross-correlations and line profiles was calculated and the 3σ errors found for each pixel, allowing an error bar to be associated with each point.

# 3 BALMER LINES

## 3.1 Hα line profiles

Fig. 2 shows the Hα line profiles obtained from SOPHIE from 2007 January 3 to February 10. All profiles appear P Cygni in nature. The profiles show a blueshifted absorption component with the edges extending up to around −250 to −300 km s⁻¹ with the deepest absorption at ~−50 km s⁻¹ and strength weakening as blueshifted velocity increases. The redshifted emission component shows greater variability in strength: only some profiles show a redshifted emission component as first noted by Croswell et al. (1987), with the red edge of the emission extending from 200 to 250 km s⁻¹, similar to those observed by Bastian & Mundt (1985). The presence of a correlation between the variability in the strength of the redshifted

emission peak and the width of the blueshifted absorption trough is initially not entirely clear, although D'Angelo et al. (2000) suggested a negative correlation.<fnref n="2" /> Fig. 2 shows that no evidence of a negative correlation is sustained in these data, for example, on JD = 2454 113.317 there is strong redshifted Hα emission, combined with broad blueshifted absorption.

To identify any quantifiable relationship between the variation of the absorption and emission components of the observed P Cygni profiles, a light curve was constructed for each pixel (every 0.5 km s⁻¹) in the line profile, which were then tested for 5000 different periods between 2 and 40 d, evenly distributed in frequency space, using the method initially described by Scargle (1982), and revised by Horne & Baliunas (1986). The Lomb–Scargle periodogram is based on least-squares sine curve fitting. Consequently, it is an ideal tool for these data because in addition to robustness to irregular coverage it also enables longer periods to be explored and is commonly used to explore up to the approximate total coverage period of the data. The significance of the periods was found to false alarm probability (FAP) of less than 0.01 using Fisher's method of randomization (Linnell Nemec & Nemec 1985).

<footnote n="2">² A variation in the minimum of the blueshifted absorption may also give rise to this variability. However, the superposition of Hα emission components makes this difficult to identify.</footnote>
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Variations in the $H\alpha$ profile of FU Orionis obtained from the SOPHIE échelle spectrograph mounted on the 1.93-m telescope at the OHP. The 21 spectra cover a time frame of 38 nights, with the Julian dates of the spectra -2454.000 displayed on the y-axis. All spectra have a resolution of 4 km s⁻¹, are barycentrically corrected and placed in the rest frame of FU Orionis. The shaded region depicts the velocities where periodicity was found to be significant below the 0.01 FAP level.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>**Figure 3.** Contour periodogram showing the significant periods found in the $H\alpha$ line profiles. The 0.01 FAP corresponds to power >7.68.</figcaption></figure>

The periodogram showing the significant periods found in the $H\alpha$ line profile is shown in Fig. 3. The periodogram highlights the significant periods with FAP < 0.01, corresponding to a power >7.68. These were in two distinct regions from -122.0 to -105.5 km s⁻¹

and from -32.5 to +85.5 km s⁻¹, shown by the shaded regions in Fig. 2. The range of significant periods detected was from 11.59 to 16.29 d, with the most probable value, found from the 'centre of power' of all the significant powers detected in the line profile, at 13.28 d. This period is similar to that found by Herbig et al. (2003), who detected a period in the 'fast wind' in the EW of the P Cygni absorption between -110 and -270 km s⁻¹ of 14.847 d. However, unlike Herbig et al. (2003), we detect the period in the line profile itself, and not in the EW of the absorption component. The overlap in velocity space of the significant periods detected with the results of Herbig et al. (2003) confirms the stability of the period in the blueshifted absorption component of $H\alpha$ over a decade. In addition, Fig. 3 extends the detection of the significant period to less negative velocities than the results of Herbig et al. (2003). The absence of periodicity from -105.5 to -32.5 km s⁻¹ suggests the periodicity detected < -100 km s⁻¹ is due to variations in the FWHM of the blueshifted absorption component of $H\alpha$ as proposed by D'Angelo et al. (2000). This 'absent' region in Fig. 2 corresponds to the deepest blueshifted absorption on all nights observed. Likewise, Fig. 2 also demonstrates less variability around the minimum of the absorption trough, which could be attributed to the suggestion by Croswell et al. (1987) that the Balmer lines are so optically thick that they are saturated, and consequently are not very vulnerable to variations in physical parameters such as density/temperature fluctuations or changes in the mass-loss rates. However, at higher negative velocities, further out in the expanding wind, the $H\alpha$ line may no longer be saturated, and environmental variations are reflected in the variability of the line profiles in such regions. Alternatively, the phase complexity of the superposition of different variable components contributing to the $H\alpha$ profile may prevent any periodicity in this region from being detected. The periodicity in the strength of the redshifted emission component is the same as the periodicity in the blueshifted absorption component, shown by the alignment in frequency of the significant peaks in the periodogram in Fig. 3. This is the first time such correlation has been detected. Previous observations by Vittone & Errico (2005) have seen periodicity in the redshifted emission component of 6.70 d, which is not evident in our results. The dense time coverage of our results confirms that this period is no longer present, and could be a result of aliasing in the earlier observations. Therefore, in contrast to previous beliefs, the periodic variability in FU Orionis seems to be stable over hundreds of periods, because the period detected here agrees well with the significant periods of 13-18 d found by Herbig et al. (2003).

The phase diagram in Fig. 4 shows the phase relationship between the periodic behaviour detected in various spectral regions phase folded with the variance-weighted average period of 13.48 d from all line profiles investigated. In each case the light curve is constructed at the power-weighted centre of velocity for the periodic region concerned. The second panel from the top shows the two periodic regions of the $H\alpha$ profile. The amplitude of variation corresponding to the oscillations in redshifted emission in the $H\alpha$ profiles (the red line) is greater than those at the faster wind velocity (the blue line). This is expected from the appearance/disappearance of the $H\alpha$ emission peak, in contrast with the broadening/narrowing FWHM of the blueshifted absorption in Fig. 2. Fig. 4 confirms that the 'absent' region in Fig. 2 is not simply lacking of periodicity due to complex superposition of the different variable components, since this would require the relative contributions to be equal and opposite, which is clearly not the case. The saturation of the $H\alpha$ line between -105.5 and -32.5 km s⁻¹ seems a more plausible explanation. The sinusoidal curve fitted to the fast blueshifted data also appears to have a slight lag of ~1.80 d in comparison with the
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Phase diagram comparing the periodicity found in different spectral lines (top five panels), with the periodicity found in the cross-correlations (bottom three panels). The spectral lines are phase folded with the variance-weighted total average period of all five regions, which is 13.48 d. The spectral lines are plotted at the power-weighted average velocity of the relevant significant detections (from top to bottom), $H\beta$ at +3.5 km s$^{-1}$, $H\alpha$ (in redshifted emission at +10.0 km s$^{-1}$ plotted with diamonds, and in blueshifted absorption at -109.5 km s$^{-1}$ plotted with crosses), NaI D (with Na D$_1$ λ5895.92 at -10.0 km s$^{-1}$ plotted with crosses and Na D$_2$ λ5889.95 at -11.0 km s$^{-1}$, plotted with diamonds), LiI λ6707 at -35.0 km s$^{-1}$ and FeII λ5018 at -39.0 km s$^{-1}$. The regions cross-correlated are shown in the bottom three panels. The cross-correlation profiles were phase folded with a variance-weighted combination of the individual periods, which averaged 3.64 d. The cross-correlations are also plotted at the power-weighted average velocity which are (from top to bottom) -9.0 km s$^{-1}$ in 6320–6440 Å, -7.5 km s$^{-1}$ in 6138–6191 Å (from cross-correlations with observed template of βAqr) and -12.0 km s$^{-1}$ in 5540–5640 Å [from cross-correlations with synthetic stellar spectrum of G0 star at $T_{\text{eff}}$ = 6000 K, log(g) = 3.0]. The solid lines show the sinusoidal line of best fitting for each individual detection, which are plotted in blue when the period is detected at blueshifted velocities (the fit to Na D$_2$ is plotted in light blue for clarity), and red corresponding to detections in redshifted velocities. All error bars are shown to the 3σ level.</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>**Figure 5.** Same as Fig. 2 except for the Hβ profile at 4861.34 Å.</figcaption></figure>

sinusoidal curves fitted to the data corresponding to oscillations in the redshifted emission peak. This could stem from the different formation locations of the two line components, with the redshifted emission formed in the disc chromosphere, and the blueshifted absorption formed further out in the expanding wind.

## 3.2 Hβ line profiles

The Hβ line profiles obtained from SOPHIE are shown in Fig. 5. All profiles show broad asymmetric absorption troughs, with stronger absorption on the blueshifted edge, which extends out to faster velocities than the redward absorption wing. This is to be expected for the Hβ profile, which is thought to be produced in an accelerating, massive outflow (Hartmann 1994). However, all the profiles appear more skewed towards blueshifted velocities than previously observed in 1992 by Hartmann & Calvet (1995), which may imply an increase in the mass-loss rate since this epoch. The maximum depth of the absorption trough appears to remain consistent at ~−25 km s⁻¹ throughout our observations. However, there is still significant night-to-night variability visible in Fig. 5. For instance, at slower negative velocities, and low positive velocities the gradient of the profile appears highly variable. This could be an effect of the varying relative absorption strengths between the redshifted and blueshifted components contributing to the overall line profile. In addition, the FWHM of the blueshifted absorption com-

<figure><img src="image_8.png" /><figcaption>**Figure 6.** Contour periodogram showing the significant periods found in the Hβ line profiles. The 0.01 FAP corresponds to power >7.02.</figcaption></figure>

ponent appears to vary significantly throughout the observations, which seems to correlate with the apparent broadening of the Hα absorption discussed in Section 3.1. This is evident when comparing observations on JD = 2454 113.317 and JD = 2454 132.307, displayed in Figs 2 and 5, both of which show relatively broad and narrower blueshifted absorption components, respectively.

The presence of any periodic variability in the Hβ line profiles was investigated between 2 and 40 d with the analysis described in Section 3.1. The region of frequency–velocity space identified as significant to below the 0.01 FAP, which represents a power >7.02 in these data, corresponds to a period from 8.82 to 14.68 d at velocities of −8.0 to +11.0 km s⁻¹ in the line profile. This feature can be seen as the red blob in Fig. 6, with a power-weighted average period of 12.89 d and at a average velocity of +3.5 km s⁻¹. This is the first detection of periodicity in Hβ. It is surprising that the periodicity is not seen at faster blueshifted velocities like Hα, given the formation of the both lines within the massive outflow. However, the examination of the top two panels in Fig. 4 shows that the periodicity in the redshifted emission of Hα and redshifted side of Hβ appears to be very similar in phase. This could suggest that it is the redshifted emission component in Hβ, formed in the disc chromosphere, and not a variation in strength of the superimposed Hβ redshifted absorption, which is the source of this periodicity. On examination of Fig. 5 at the velocities where a period was identified (the grey shaded region between −8.0 and +11.0 km s⁻¹), it appears that this cyclic variability corresponds to the shallowing/deepening of the ‘ridge’ connecting the blueshifted absorption to the redshifted absorption. Consequently, a variation in Hβ emission could produce this effect, and a narrower formation region of Hβ emission in comparison to Hα emission within the chromosphere or the additional superimposed variation in the redshifted absorption component of Hβ may prevent this periodicity being detected at faster redshifted velocities alike Hα.

## 3.3 Hγ line profiles

The presence of periodic variability was also tested analogously in the Hγ profiles. In this case, the observations showed no evidence for periodicity in the profiles from −300 to +300 km s⁻¹ above the random noise level, which is higher in this shorter wavelength region.

# 4 NA I D PROFILES

The Na I profiles obtained from SOPHIE are shown in Fig. 7. All profiles show complex absorption components, with no noticeable
<figure><img src="image_9.png" /><figcaption>**Figure 7.** Same as Fig. 2 except for the Na I D profile of FU Orionis. The lower axis is centred around Na I D₁ at 5895.92 Å, and the upper axis is plotted with respect to Na I D₂ at 5889.95 Å.</figcaption></figure>

emission components, as first observed in 1981 by Bastian & Mundt (1985). The blueshifted profiles both seem to contain two domi-
nate absorption components, a high-velocity broad feature, super-
imposed with a narrow, low-velocity absorption dip (Croswell et al.
1987). The broad blueshifted feature appears to be variable during
our observations with the edge extending out to velocities between
−150 and −250 km s⁻¹ (slightly slower than the expansion veloci-
ties seen in Hα in Fig. 2). This is consistent with the lines' formation
in a differentially expanding wind, similar to the broad blueshifted
absorption of Hα, but in a more compact region, with the wind
continuing to expand throughout the more extended Hα formation
region (Bastian & Mundt 1985). Bastian & Mundt (1985) went on
to conclude that the temperature of the accelerating region must in-
crease outwards, because Hα generally needs higher temperatures
to be formed than Na I D, yet Croswell et al. (1987) showed that
the Hα absorption could be formed at lower temperatures with an
increased mass-loss rate. In addition, the relatively broad absorption
feature of Na I D seems to be comprised of two separate variable
absorption components, which are visible separately on some nights
and appear blended together on others (see JD = 2454 104.319 and
JD = 2454 130.292 in Fig. 7). This is consistent with Croswell
et al. (1987) and Hartmann & Calvet (1995), who observed a single
component and Bastian & Mundt (1985) who observed two com-

<figure><img src="image_10.png" /><figcaption>**Figure 8.** Contour periodogram showing the significant periods found in the Na I D line profiles. The 0.01 FAP corresponds to power >6.67.</figcaption></figure>

ponents in the broad absorption features. On the other hand, the
narrow low-velocity component ~−3 km s⁻¹ does not appear to
vary significantly during our observations, which is consistent with
its proposed formation in an expanding circumstellar shell, at large
distances from the object, as discussed by Bastian & Mundt (1985).

The presence of periodicity in the Na I D line profiles between 2
and 40 d was analysed using the same methods described in Sec-
tion 3.1. The significant periods found in this doublet are shown by
the two red regions in Fig. 8, which are below the 0.01 FAP, with
power >6.67. In the Na D₁ line at 5895.92 Å, periodicity is signifi-
cant to less than the 0.01 FAP between −16.0 and −5.5 km s⁻¹ for
periods from 12.41 to 17.87 d. The range of significant periods de-
tected in the Na D₁ line not only overlaps with the periods detected
to less than the 0.01 FAP in the Na D₂ line at 5889.95 Å, but are
also present at the same relative velocities in the line profile, shown
by the shaded regions in Fig. 7. Significant periodicity was detected
in Na D₂ between −25.0 and −5.5 km s⁻¹ with periodicity ranging
from 11.59 to 17.54 d. This is the first detection of periodicity in the
Na I D profile and the symmetry of these two individual detections
justified their combination to calculate an overall power-weighted
average period of 14.27 d for the Na I D profile. There are also some
other detections in Fig. 8 that appear to have power values above
the power corresponding to the 0.01 FAP level; however, these de-
tections are situated on the very edge of the periodogram at small
frequencies and do not coincide with an extended feature. There-
fore, these ‘significant detections’ are assumed to be a result of edge
effects in the frequency sampling and are not considered in further
analysis.

The origin of the periodicity highlighted by the grey shaded re-
gions on the line profiles in Fig. 7 is difficult to distinguish due
to the complexity of the blueshifted absorption components con-
tributing to both lines. However, the periodicity of both lines in
the doublet appears to be cut-off at the edge of the narrow low-
velocity blueshifted absorption component. This is consistent with
the formation of this low-velocity component in a slowly expand-
ing shell at a large distance from the object. Therefore, the peri-
odicity is more likely associated with the lower velocity end of
the broader, more variable blueshifted components, which are at-
tributed to the expanding wind. This is supported by the periodicities
overlapping with the significant detections in Hα and Hβ. The de-
tection of these periods at relatively low velocities in comparison
to the maximum velocity of the outflow ~−150 to −250 km s⁻¹ in
this case suggests that the origin of this periodicity is lower down
in the expanding wind, in contrast to Hα in Section 3.1, where
periodicity in the blueshifted absorption component was seen cen-
tred at −109.5 km s⁻¹. This disparity could be due to a combination
of the more extended formation region of Hα and the relatively
reduced temperature sensitivity of the Na I D (Croswell et al. 1987),
which may prevent a significant periodic signal at the edges of the
formation region, in the expanding wind. Fig. 4 supports this since
the third panel from the top shows the periodicity in both lines phase
folded with a variance-weighted average overall period of 13.48 d.
Both Na I D lines (D₁ shown in dark blue and D₂ shown in light
blue) appear to have a slight lag with the Hα emission, thought to
be formed in the disc chromosphere, and a slight lead over the Hα
absorption, formed further out in the expanding wind, shown in the
panel above in Fig. 4.

5 METALLIC LINE PROFILES

Line profile variations in Li i λ6707 from 2007 January 3 to February 10 are shown in Fig. 9. These profiles support the accretion disc hypothesis of FU Orionis objects as described by Hartmann & Kenyon (1985, 1987a) and Zhu et al. (2009a), since the profile displays a double-peaked nature, as initially observed by Herbig (1989); however, variability only exists in the blueshifted absorption. Unlike Hα, as described in Section 3.1, the width of the blueshifted absorption appears to remain constant. The blueshifted absorption instead varies in strength, causing the double-peaked line profile to appear more and less asymmetrical, see JD = 2454 121.333 and JD = 2454 127.377 in Fig. 9. These observations coincide with the dates when the Hα profile had no redshifted emission and strong redshifted emission, respectively (see Fig. 2). This could indicate a link between the origin of the variability in the Hα line and Li i λ6707.

The possibility of periodicity between 2 and 40 d in the Li i λ6707 line profile was investigated via the method described in Section 3.1, producing significant detections at velocities from −43.0 to −27.0 km s⁻¹, with periods ranging between 12.21 and 19.78 d, resulting in a ‘centre of power’ period at 14.81 d. This is the first detection of periodicity in Li i and can be clearly identified as the large orange-red structure in Fig. 10, which displays the resultant periodogram for the frequency–velocity space analysed. There appears to be an addition feature on the periodogram in Fig. 10 that stands out in yellow above the background at velocities between −69.5 and −68.0 km s⁻¹ corresponding to a period from 21.22 to 29.39 d, and a power-weighted average period of 24.76 d. This feature is significant below the 0.01 FAP level, which corresponds to a power >6.76. On inspection of Fig. 10, it is unclear whether this periodicity is connected to the larger feature ~−35 km s⁻¹, as it is located at slightly more negative velocities towards shorter frequencies. Since there is no overlap in frequency space between the significant periods detected, it is treated as a separate detection and for the remainder of this discussion is classed as noise fluctuations, until further detections prove the contrary.

The location of the main spike in power in Fig. 10 coincides with
the deepest absorption in the Li i λ6707 profile in Fig. 9, located
at ~−35 km s⁻¹, with the connection highlighted by the shading
of the periodic region. This confirms that the periodicity in Li i is a
result of the varying strength of the blueshifted absorption peak, and
not a periodic broadening of the line profile. The confinement of
periodicity to purely negative velocities indicates that this could be
another periodic wind signature, although it was previously thought
that the Li i absorption line was not formed in the wind (Calvet et al.
1993). Despite the apparent difference in nature of the periodicities
on the blue wing with Hα, the periodicity observed is similar in

<figure><img src="image_11.png" /><figcaption>**Figure 9.** Same as Fig. 2 except for the Li i profile at 6707.89 Å.</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>**Figure 10.** Contour periodogram showing the significant periods detected in the Li i line profile. The 0.01 FAP corresponds to power >6.76.</figcaption></figure>

nature to that of Na I D in Fig. 7, but observed at slightly higher
velocities. This difference could be down to the separate formation
location of the lines in the expanding wind, since the periods are
similar. The fourth panel down in Fig. 4 shows the data in the power-
weighted average velocity channel of the significant period detected
(−35.0 km s⁻¹), phase folded with the overall average period of
13.48 d. (The range of significant detections of periodicity in Li i
includes the 13.28 and 14.27 d periods found in Hα and Na I D
profiles, respectively, so it is plausible that all these periods could
be produced via the same mechanism.) The diagram in Fig. 4 shows
clear cyclic variation and, when fitted with a sinusoid, shows that
<figure><img src="image_13.png" /><figcaption>Figure 11. Same as Fig. 2 except for the Fe II profile at 5018.60 Å.</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>Figure 12. Contour periodogram showing the significant periods detected in the 5018 Fe II line profile. The 0.01 FAP corresponds to power >6.73.</figcaption></figure>

the variations are roughly in-phase with the variations of the Na I D components, displayed in the panel above. The intermediate phase of both Li I and Na I D between that of the chromosphere (red lines in Fig. 4) and the higher velocity blueshifted components of Hα are consistent with both the presence of periodicity and the lines' formation lower down in the expanding wind.

Another metallic absorption line not as optically thick as Hα is the Fe II line at 5018 Å, suggested by Herbig (1966) to be produced further out in the wind in a 'shell' feature of FU Orionis displaced ~−80 km s⁻¹. Fig. 11 shows the λ5018 line profiles, which show

no indication of double-peaked structure, but instead are predominately asymmetric with maximum absorption depth ~−50 km s⁻¹, slightly greater than Li I, more similar to Hα and Hβ. On some nights, most notably JD = 2454 121.333, the absorption forms a sharp dip, which could be attributed to a shell-like signature in the wind. However, this spectral feature was identified by Calvet et al. (1993) to be present when the disc rotational velocity is comparable to the expansion velocity in the wind, thus the 'shell' signature is disingenuous. Consequently, the loss of the region of sharp decline, and the production of a more symmetric blueshifted absorption profile in Fig. 11 by JD = 2454 127.377, has implications for the mass-loss rate and expansion velocity in the wind. The blueshifted absorption appears to extend out to ~−125 km s⁻¹; however, both the width and the strength of the absorption profile in general appear to vary, comparable to the Hα and Li I λ6707 profiles, respectively.

The periodogram shown in Fig. 12 was obtained via the method described in Section 3.1, over the same parameter space for comparison. It shows less significant detections of a period than the Li I λ6707 profile, but there are still periods significant to below the 0.01 FAP with velocities from −43.5 to −24.0 km s⁻¹ for values from 12.31 to 15.68 d. This corresponds to the largest feature in Fig. 12, which has a power-weighted average period of 13.57 d. The dominant periodic feature detected in this line profile is yet again of a similar value to that of Hα, Hβ, Na I and Li I, providing further evidence that there is a wind period in FU Orionis. This is confirmed in Fig. 4, where the fourth panel from the bottom shows the Fe II λ5018 data in the power-weighted average velocity of the periodic feature (−39.0 km s⁻¹), phase folded with the average overall wind period of 13.48 d. The cyclic variations are fitted with a sinusoid (solid blue line), and the fit is better than that for Na I D and Li I in the two panels above. However, unlike the previous line profiles, the phase of the Fe II line does not lie within the range expected for formation in the expanding wind, and seems to respond even before the chromospheric lines in the top two panels in Fig. 4. Although there are less significant detections of a period to less than 0.01 FAP visible in Fig. 12 than in Li I λ6707 (see Fig. 10), possibly due to the high temperature sensitivity of the Li I line (Calvet et al. 1993), the amplitude of variation is greater in the Fe II line at λ5018 (see Fig. 4). On examination of Fig. 11, the velocities at which this period was detected (shown by the shaded region) seem to correlate with the appearance and disappearance of the absorption 'dip', confirming the assumption of the period's origination in a wind. The fact that there are less significant detections of periodicity could indicate that the effects of the mechanism for producing this period in the wind may become reduced at greater distances from the disc where the Fe II λ5018 absorption line is thought to originate.

As in the Li I λ6707 profile, there is also another feature visible on the periodogram in Fig. 12 at slightly lower frequency and more negative velocity than the dominant feature, corresponding to a weak detection of periodicity, above the 0.01 FAP level (power >6.73). The range in velocity of this yellow feature is −128.0 to −123.5 km s⁻¹ and spans periodicities between 19.52 and 23.91 d. This is comparable to the detection discussed in Li I λ6707, where the authenticity of the detection was unclear. In Fe II λ5018 this detection has a power-weighted period of 21.66 d, slightly lower than the comparable signal in Li I λ6707, yet there is a substantial overlap in the significant periods detected in both lines. The velocities at which the period manifests itself are slightly higher in Fe II λ5018 than in Li I λ6707; however, this could be ascribed to the different formation locations in an expanding wind. The recurrence of this feature at a comparable periodicity to an ambiguous detection in Li I λ6707 that affects a larger area of surrounding parameter space in
Fig. 12 and then in Fig. 10 makes this detection slightly less likely due to noisy line profiles and may in fact be linked to the main source of periodicity in both cases.

## 6 CROSS-CORRELATIONS

Some of the metal lines were too weak and too noisy to analyse the line profile alone. In order to increase the signal-to-noise ratio, these profiles were cross-correlated with a template spectrum of β Aqr and a synthetic stellar spectrum of a G0 star with $T_{\text{eff}} = 6000$ K, $\log(g) = 3.0$ to identify any potential effects of relative line shifts, as described in Section 2. The standard spectra in the 5540–5640 Å spectral region are shown in Fig. 13, where the bottom black line is the synthetic spectrum and the blue line shows the Gaussian smoothed synthetic spectrum to match the FWHM of the unblended, unsaturated lines in the observed spectrum of β Aqr, shown by the top red line. The regions cross-correlated were chosen in order to compare with previous results by Herbig et al. (2003), where a 3.542 d period has already been identified in the centroid of the cross-correlation, and Hartmann & Kenyon (1987a), who found variability in their cross-correlation profile of the 6170 Å region. The regions investigated by Herbig et al. (2003) were the 6320–6440 Å region and the 5540–5640 Å region, which were initially identified in order to avoid lines with wind components. The normalized SOPHIE échelle spectra taken from 2007 January 3 to February 10 are shown in Figs 14–16. The spectra in the 5540–5640 Å region (see Fig. 14) were normalized using a lower bias level (0.97 instead of 0.98), because the continuum appeared under subtracted when the higher value was used. All spectra show evidence of weak absorption lines and some weak lines appear doubled. This is most notable in Fig. 15 where the 6162.17 Å Ca I and 6169.56 Å Ca I absorptions appear doubled on some nights, indicating a disc origin, although the profiles are highly variable. Other lines noted by Hartmann & Kenyon (1987a) to show evidence of double-peaked structure in this spectral region are the 6147.8 Å Fe I and Fe II line blend, 6151.62 Å Fe I, 6157.73 Å Fe I and 6166.44 Å Ca I. There does appear to be some evidence of a skewed absorption profile in the 6141.73 Å Fe I line, indicating the presence of mass-loss. Cross-

correlating these spectral regions provides a distinct advantage as not only does it produce an average profile weighted according to the standard line strengths, but it also allows for lines which overlap and blending, abundant in these regions, which is clearly a problem for identifying the weak lines in Figs 14–16. The CCFs of each night were tested for periodicity at 5000 different frequencies between 2 and 40 d to quantitatively compare the results with the line profiles in Sections 3–5.

The CCF of the 6138–6194 Å region, displayed in Fig. 17, showed the clearest detection of a period significant to less than 0.01 FAP (power >6.88) in velocity channels −14.0 to −1.5 km s⁻¹. The significant periods range from 3.52 to 3.70 d and have a power-weighted average of 3.61 d. This region is also clearly visible above the noise in Fig. 18. This period is similar in value to that detected previously, yet it is distinctly different in signature, corresponding to a cyclic strengthening/weakening of the blue wing of the line profiles in this region. To compare with Herbig et al. (2003), the radial velocity of centroid of the CCF was found for each observation and tested for periodicity between 2 and 40 d, yet no significant period was found below the 0.01 FAP. This was unexpected, as the periodicity has only been observed on the blue wing of the CCF. To make a more direct comparison, the spectra were degraded to a resolution of 13 km s⁻¹ used by Herbig et al. (2003). At this resolution a longer period of 13.12 d (comparable to the periodicity previously observed in the line profiles directly) was detected in the centroid (or weighted average velocity) of the cross-correlation in the 6138–6194 Å region, which was significant to the 0.01 FAP level. This could possibly indicate that the effects of the wind are more significant in this spectral region, than in the regions investigated by Herbig et al. (2003), and the shorter periodic effects may be masked in the centroids due to the other variable features in the CCF, most notably the effects of the wind. The power-weighted average radial velocity of the ‘blob’ in Fig. 18 was −7.5 km s⁻¹, and the phase-folded data at this velocity are shown in the second panel from the bottom in Fig. 4. The precise velocity at which this periodicity is detected is difficult to validate, since spectral lines in the observed G0 Ib template of β Aqr may be shifted due to nightly fluctuations in air mass, which could ultimately shift the velocity

<figure><img src="image_15.png" /><figcaption>**Figure 13.** Reference spectra used in the cross-correlations with FU Orionis, shown in the 5540–5640 Å region. The bottom black line is a synthetic stellar spectrum of $T_{\text{eff}} = 6000$ K, $\log(g) = 3.0$, the blue line is the same spectrum smoothed with a Gaussian with FWHM of 0.37 Å, which was found to correspond to the FWHM of the unblended, unsaturated lines in the β Aqr spectrum. The top red line is the spectrum of β Aqr, obtained from the SOPHIE échelle spectrograph mounted on the 1.93-m telescope at the OHP with a resolution of 4 km s⁻¹.</figcaption></figure>
<figure><img src="image_16.png" /><figcaption>**Figure 14.** Spectra of FU Orionis in the 5540–5640 Å region obtained from the SOPHIE échelle spectrograph mounted on the 1.93-m telescope at the OHP. The 21 spectra have a resolution of 4 km s⁻¹ and cover a time frame of 38 nights, with the Julian dates of the spectra −2454 000 displayed on the right-hand y-axis. All spectra have been barycentrically corrected, normalized to unity, placed in the rest frame of FU Orionis, continuum subtracted with a bias level of 0.97 and shifted to show night-to-night variability.</figcaption></figure>

peaks in the cross-correlations. However, the velocity range of the periodicity was reinforced by cross-correlations with the synthetic stellar spectrum at $T_{eff}$ = 6000 K, log(g) = 3.0, which confirmed periodicity between 3.54 and 3.71 d to below 0.01 FAP spanning ve-

locities from −16.0 to −4.0 km s⁻¹, with a power-weighted average period of 3.62 d.

The detections below the 0.01 FAP level (power >6.73) in the 6320–6440 Å region cross-correlated with the observed template of
<figure><img src="image_17.png" /><figcaption>**Figure 15.** Same as Fig. 14, but for the 6170 Å region using a bias level of 0.98 for the continuum subtraction.</figcaption></figure>

β Aqr were periods from 3.56 to 3.71 d, spanning velocity channels from −18.5 to 0.0 km s⁻¹, which produced a weighted average period of 3.63 d. Fig. 19 shows the resulting CCFs, which look very similar to the cross-correlations of the 6138–6194 Å region in Fig. 17. Fig. 20 shows the significant detections with the red blob at the centre of the figure, at a frequency of 0.28 d⁻¹. When the centroids of the CCFs in this region were tested for periodicity, yet again there was none found to be significant to below the 0.01 FAP level. However, when the spectra and observed template were degraded to the resolution used by Herbig et al. (2003) in this region, the centroids of the cross-correlation peaks were found to have a periodicity of 13.69 d significant to the 0.01 FAP level. This could possibly indicate an increase in the mass-loss rate, which would manifest wind signatures more readily in this epoch. This may also be a possible explanation why Herbig et al. (2003) noted a poor fit to the sinusoid for the observations taken in the year 2000. The velocities at which this periodic signal was detected were confirmed by the cross-correlations with the synthetic stellar spectrum, which also detected a period of 3.63 d between −19.0 and −1.0 km s⁻¹, significant

<figure><img src="image_18.png" /><figcaption>**Figure 16.** Same as Fig. 14, but for the 6320–6440 Å region using a bias level of 0.98 for the continuum subtraction.</figcaption></figure>

to below the 0.01 FAP. The phase folded data at −9.0 km s⁻¹ fitted with a sinusoid (shown with the blue solid line) can be seen on the third panel from the bottom in Fig. 4 with error bars plotted to the 3σ level.

Our analysis showed no significant periodicity below the 0.01 FAP corresponding to a power of 5.88 in the 5540–5640 Å region when cross-correlated with β Aqr. Unsurprisingly, no period was detected in the centroid of the cross-correlation of the spectra and observed template degraded to the same resolution as Herbig et al. (2003) in this spectral region. However, a period between 3.63 and 3.65 d was found in the same spectral region when cross-correlated with the synthetic template, which was significant below the 0.01 FAP, corresponding to powers above 6.16 in this case. The weighted average period was found to be 3.64 d, from −15.0 to −9.0 km s⁻¹ and centred at −12.0 km s⁻¹. The resulting cross-correlations with the smoothed synthetic spectrum are shown in Fig. 21 with the periodograms displaying the null detection using the β Aqr template, and the significant periods found with the synthetic template in Figs 22 and 23, respectively. The background levels at the edges of
<figure><img src="image_19.png" /><figcaption>**Figure 17.** Cross-correlations of FU Orionis in the 6138–6194 Å region with a reference spectrum of β Aqr obtained from the SOPHIE échelle spectrograph mounted on the 1.93-m telescope at the OHP. The cross-correlation function is the black line, with the upper and lower $3\sigma$ errors shown by the dark blue and light blue lines, respectively. When two observations were taken on a single night a dark green line is used to represent the second observation. The average over all cross-correlations in this region is shown by the red dashed line.</figcaption></figure>

<figure><img src="image_20.png" /><figcaption>**Figure 18.** Contour periodogram showing the significant periods found in the 6138–6194 Å region. The 0.01 FAP corresponds to power >6.88.</figcaption></figure>
<figure><img src="image_21.png" /><figcaption>**Figure 19.** Same as Fig. 17, but for the 6320–6440 Å window.</figcaption></figure>

<figure><img src="image_22.png" /><figcaption>**Figure 20.** Contour periodogram showing the significant periods found in the 6320–6440 Å region. The 0.01 FAP corresponds to power >6.73.</figcaption></figure>
<figure><img src="image_23.png" /><figcaption>**Figure 21.** Same as Fig. 17, but for the 5540–5640 Å window, which was continuum subtracted with a lower bias of 0.97 and cross-correlated with a G0 synthetic stellar template at $T_{\text{eff}} = 6000$ K, $\log(g) = 3.0$.</figcaption></figure>

<figure><img src="image_24.png" /><figcaption>**Figure 22.** Contour periodogram of the 5540–5640 Å cross-correlation, which was continuum subtracted using a lower bias of 0.97 and cross-correlated with observed reference template of β Aqr. There are no detections in the cross-correlation function below the 0.01 FAP which corresponds to power >5.88.</figcaption></figure>

<figure><img src="image_25.png" /><figcaption>**Figure 23.** Contour periodogram showing the significant periods found in the 5540–5640 Å cross-correlation, which was continuum subtracted with a lower bias level of 0.97 and cross-correlated with a G0 synthetic stellar template at $T_{\text{eff}} = 6000$ K, $\log(g) = 3.0$. The 0.01 FAP corresponds to power >6.16.</figcaption></figure>
the CCFs are slightly higher in Fig. 21, compared to Figs 19 and 17, because of the lower bias level used in the continuum subtraction. However, the shape remains unaffected in this adjustment, and on visual examination the variability seems to correlate with the cross-correlations of the 6320–6440 Å and 6138–6194 Å spectral regions using the observed template of β Aqr. The absence of periodicity in the cross-correlations of the 5540–5640 Å region with the observed template and yet the significant detection with the synthetic template could be manifested by the slightly different relations between spectral line strengths in these two spectra, shown with the top two lines plotted in Fig. 13. Figs 22 and 23 are extremely similar and both show a clear feature around the periodicity detected to be significant to below a FAP of 0.01 in the cross-correlations with the synthetic template. Therefore, although this feature is not detected to be significant in Fig. 22 to such a high level when the observed template of β Aqr is used, it is still a clear feature above the noise in the periodogram and hence the lack of high significance of the feature in this case may be due to the slight differences in the line weightings producing slightly different weights in the CCF and consequently, the periodogram. In the case of the synthetic standard, these slight differences have contributed to the significance of the cyclic variability in this spectral region. The absence of any significance to less than 0.01 FAP with the β Aqr standard may therefore not indicate the lack of cyclic variability, but rather inadequate signal strength in this spectral region possibly due to a different mix/ratio of line species in the observed template or a slightly different origin of these absorption lines within the disc itself in comparison to the spectral regions previously discussed. The bottom panel of Fig. 4 shows the data at −12.0 km s⁻¹ phase folded with a period of 3.64 d and fitted with a sinusoid. For the reasons discussed above, the amplitude of variation is small in comparison with the cross-correlations in the other spectral regions, shown in the two panels above in Fig. 4.

The in-phase nature of the periods found in the CCFs, along with the extensive overlap in velocity and range of significant periods identified to less than 0.01 FAP, is strong evidence that all these cyclic variations are produced from a single origin. Consequently, the three periods identified above were combined to find a variance-weighted total average period of 3.64 d, which was used to phase fold the data in all regions plotted in Fig. 4. The periodicity could be attributed to a slow moving wind, as it is predominantly present at small negative velocities. However, it does not appear to be related to the wind periods found in the line profiles in Sections 3–5, as it is roughly one-quarter of the value.

To complete the analysis of the cross-correlations, the FU Orionis spectra from each night were cross-correlated with the average spectrum and the resulting CCFs were tested for periodicity via the same analysis. Interestingly in the 6138–6194 Å region, this resulted in periods between 12.31 and 14.32 d, significant to below 0.01 FAP with a power-weighted average of 13.19 d, which are similar to the periodicities detected in the centriods of the spectra degraded to the resolution of Herbig et al. (2003) and in Sections 3–5. Additionally, this periodicity was detected at velocities between −80 and −24 km s⁻¹, more comparable to the fast velocities of an expanding wind. No such periodicity was detected in the other spectral regions, possibly due to a different mix of line species and their selection by Herbig et al. (2003) to avoid line with wind components. There was no evidence of the 3.64 d period in this analysis. However, the resulting cross-correlations are extremely broad due to the large FWHM of spectral lines in the FU Orionis spectra in comparison with the FWHM of the standard templates (see Figs 13 and 14 for comparison). Therefore, the CCFs produced

<table><caption>**Table 2.** Summary of periodicities detected: (1) line/spectral region where periodicities were detected; (2) weighted average velocity of periodic feature with the limits showing the full range of velocity channels where periodicity was significant to below the 0.01 FAP; (3) weighted average periodicity detected with the full range of periodicities detected above the 0.01 FAP shown in the limits.</caption><thead><tr><th>(1)<br>Spectral region</th><th>(2)<br>Velocity<br>(km s<sup>-1</sup>)</th><th>(3)<br>Period<br>(d)</th></tr></thead><tbody><tr><td>Hα</td><td>10.0<sup>+75.5</sup><sub>-22.5</sub><br>-109.5<sup>+4.0</sup><sub>-12.5</sub></td><td>13.27<sup>+3.02</sup><sub>-1.68</sub><br>13.85<sup>+0.99</sup><sub>-0.72</sub></td></tr><tr><td>Hβ</td><td>3.5<sup>+7.5</sup><sub>-11.5</sub></td><td>12.89<sup>+1.79</sup><sub>-4.07</sub></td></tr><tr><td>NaI D<sub>1</sub></td><td>-10.0<sup>+4.5</sup><sub>-6.0</sub></td><td>14.49<sup>+3.38</sup><sub>-2.08</sub></td></tr><tr><td>NaI D<sub>2</sub></td><td>-12.0<sup>+6.5</sup><sub>-13.0</sub></td><td>14.09<sup>+3.46</sup><sub>-2.50</sub></td></tr><tr><td>Fe II λ5018</td><td>-39.0<sup>+15.0</sup><sub>-4.5</sub></td><td>13.57<sup>+2.11</sup><sub>-1.26</sub></td></tr><tr><td>Li I λ6707</td><td>-35.0<sup>+8.0</sup><sub>-8.0</sub></td><td>14.81<sup>+4.97</sup><sub>-2.60</sub></td></tr><tr><td>6138–6194 Å<sup>a</sup></td><td>-7.5<sup>+6.0</sup><sub>-6.5</sub></td><td>3.61<sup>+0.09</sup><sub>-0.09</sub></td></tr><tr><td>6320–6440 Å<sup>a</sup></td><td>-9.0<sup>+9.0</sup><sub>-9.5</sub></td><td>3.63<sup>+0.08</sup><sub>-0.08</sub></td></tr><tr><td>5540–5640 Å<sup>b</sup></td><td>-12.0<sup>+3.0</sup><sub>-3.0</sub></td><td>3.64<sup>+0.01</sup><sub>-0.02</sub></td></tr></tbody></table>

<sup>a</sup> The periodicities in these regions were detected in the cross-correlation function using an observed template of a G0lb star namely β Aqr.

<sup>b</sup> The periodicity in this region was detected in the cross-correlation function using a synthetic spectrum of a G0 star with $T_{\text{eff}} = 6000$ K, $\log(g) = 3.0$ as a reference template.

using the average spectra as a template are insensitive to small changes around the line centre, which would be necessary to detect the 3.64 d period. Consequently, there appears to be more than one mechanism producing asymmetries in FU Orionis.

## 7 DISCUSSION

A summary of the velocity channels where periodicity was detected to be significant below the 0.01 FAP level in each spectral region is shown in Table 2. Our high-resolution spectral monitoring of FU Orionis has yielded two principle results.

(i) We have confirmed the periodic variability of blueshifted Hα absorption claimed by Herbig et al. (2003). The ‘centre of power’ for periodic variation is at a period of 13.48 d (which we will henceforth term ‘the wind period’), but the range of periods suggested by our analysis (see Table 2) also includes the 14.847 d period found by Herbig et al. (2003). In addition, we have for the first time, found periodic variation (at the same period) in the redshifted emission component of Hα, the redshifted absorption component of Hβ and the blueshifted absorption components of NaI D, LiI λ6707 and Fe II λ5018. Although the ‘centre of power’ for periodic variations in these lines is slightly different (13.27, 12.89, 14.27, 14.81 and 13.57 d, respectively), the range of allowed periods at a FAP of 0.01 all overlap and we will assume that the same wind period (a variance-weighted average of 13.48 d) applies to all these lines.

(ii) We have also confirmed the periodic variation of the CCF constructed from spectral regions that, lacking strong wind features, are believed to originate from the photosphere of the disc. The detected period (around 3.64 d in all the spectral windows investigated) is close to the value of 3.54 d found by Herbig et al.
(2003). We find that variability is largely confined to low velocities on the blue wing of the CCF ($> -18.5 \text{ km s}^{-1}$ and with a velocity centroid at $\sim -9.0 \text{ km s}^{-1}$); a change in the characteristics of the outflow in FU Orionis may explain why Herbig et al. (2003), with a velocity resolution of $\sim 13 \text{ km s}^{-1}$ compared with $4 \text{ km s}^{-1}$ in the present study, also found significant variation on the red wing.

In what follows we will assume that the 13.48 period is associated with modulation of the wind. The fact that the 3.64 d period is only detectable on the blue wing of the CCF raises the question of whether it might not also be a wind period, even though the spectral regions from which the cross-correlations were calculated are believed to originate in the disc photosphere. We are dissuaded from invoking a wind component in these lines on the grounds that in any case the period is quite different from the 'wind period' described in (i) above. We note that the two periods are *consistent* with being in the ratio of 1: 4, although the range of allowed periods at a FAP of 0.01 or below does not allow us to assert this with any confidence. Consequently, we will assume that there is no necessary connection between the phenomena giving rise to these two periods and that they are associated with two physically distinct regions, the wind and the disc photosphere. We will briefly revisit this assumption at the end of this section.

## 7.1 The periodicity of the wind

We have found clear evidence (see the upper four panels of Fig. 4) that the variations in the redshifted emission of the $H\alpha$ and $H\beta$ lines are in-phase with each other, with the blueshifted absorption of the Na I D and Li I $\lambda6707$ displaying an intermediate phase. The variation of the blueshifted absorption of $H\alpha$ is however phase lagged by about 1.8 d.

This pattern is consistent with a picture in which the wind is modulated at its base in velocity and/or density/temperature. The first lines to respond are those formed at relatively low velocities close to the wind base, while the blueshifted $H\alpha$ absorption (which is at considerably higher velocities) responds after a time equal to the flow time across the region of acceleration. This hypothesis allows us to estimate the scale of the region where the wind is accelerated up to $\sim 100 \text{ km s}^{-1}$ as being the product of $100 \text{ km s}^{-1}$ and 1.8 d, i.e. $\sim 10^{12}$ cm. It should be noted, however, that the phase of the Fe II line does not fit this picture. The lack of significant periodic behaviour at still higher velocities could either be due to the shearing out of high-density pulses in the wind beyond this region or else because the higher velocity absorption originates from a larger region in which a number of consecutive pulses are always simultaneously present.

## 7.2 The periodicity of photospheric absorption

A simple interpretation of periodic variability in the line profiles from a disc is in terms of an orbiting hotspot which enhances the flux (and hence contribution to the absorption line spectrum) at a particular velocity. The problem with this simplest model is that it would predict that enhanced absorption should switch between the blue and red wings of the line as the hotspot alternatively approaches and recedes from the observer. However, we find little variability on the red wing of the line and so we need to argue that the receding hotspot is in some way concealed from the observer. This necessarily involves the invoking of three-dimensional effects in the disc such as a warp. If for some reason the disc surface is distorted in the region of

enhanced dissipation, then the radiation may escape preferentially in the forward direction of the orbiting feature.

So far we have been deliberately unspecific about the nature of the hotspot or warp. We bear in mind that our confirmation of Herbig's period implies that this is a phenomenon that is stable over the order of a thousand orbits or more and is thus inconsistent with a transient heating event which would be sheared out over many fewer orbital periods. Clarke & Armitage (2003) suggested that a long-lived hotspot would be produced if the inner disc of FU Orionis contained an embedded (proto-)hot Jupiter, as required by the model for the triggering of rapid-rise objects (such as FU Orionis and V1057 Cygni) suggested by Clarke & Syer (1996) and Lodato & Clarke (2004). In this case, the 'hotspot' is simply the emergence of the energy liberated by accretion on to the planet as it diffuses out through the overlying disc. If, in addition, the disc is warped in the vicinity of the planet (as would result if the planet and disc planes are not perfectly aligned), then this radiation will preferentially emerge through the region of the disc surface which offers the lowest optical depth to escaping photons; thus there are geometries for which the hotspot signature would be more detectable when the planet is approaching the observer than on the opposite phase.

In a similar vein, one could instead invoke a magnetized hotspot on the disc surface (supported by the detection by Donati et al. 2005) and again, if the axis of the magnetic field is misaligned with that of the disc, one expects a warped structure in the disc (Lai 1999, 2003; Foucart & Lai 2011). The stability of the observed period implies that the magnetic topology needs to remain constant on time-scales of $\sim 10$ yr.

Note that it is not a necessary ingredient of the model that the hotspot corotates with the local disc flow, although obviously in the planet case, both hotspot and disc orbit at the local Keplerian velocity. If this is the case then a period of 3.6 d corresponds to an orbital radius of $7 \times 10^{11} M^{1/3}$ cm, where $M$ is the mass of the central star in solar masses.<fnref n="3" /> The maximum velocity shift would in this case be $140 M^{1/3} \sin i > 60 M^{1/3} \text{ km s}^{-1}$, where in the latter inequality we have adopted the minimum value of $\sin i$ quoted in the literature (Kenyon et al. 1988) and $M$ is the stellar mass in solar masses. Even in the case of rather a low-mass star $M \sim 0.3$, the resulting peak projected velocity is around $40 \text{ km s}^{-1}$. Inspection of Figs 17, 19 and 21 shows that this is close to the maximum velocity at which there is detectable variation on the blue wing of the CCF, although a significant period is only detectable for velocities less than $18.5 \text{ km s}^{-1}$. The amplitude of the velocity variation is thus marginally consistent with a Keplerian hotspot but only if one adopts rather low values for the stellar mass and inclination of FU Orionis.

On the other hand, we have more freedom in the case of a magnetically generated hotspot since in this case the concentration of magnetic field lines associated with the hotspot are not necessarily co-orbital with the disc – it is now possible for the combination of low-velocity amplitude and short period to be associated with an interaction between the field and the disc occurring well outside corotation. Note that if one identified the rotation period of the magnetic field with that of the star, then it would be rotating at a velocity modestly less than break-up velocity.

<footnote n="3">³ Note that this radius is factor 2–3 larger than the best-fitting inner disc radius in FU Orionis employed by Kenyon et al. (1988).</footnote>
## 7.3 The relationship between the periods

We finally return to the issue of whether there can be any connection between the two varieties of periodic behaviour described above. The phase relationship between the low- and high-velocity wind components suggests that we are witnessing the effects of pulsed input at the base of the wind flow rather than the effects of rotational modulation of an intrinsically steady but non-axisymmetric wind. In this case we need some process to modulate the injection of momentum into the wind close to its base where it interacts with the disc. This process does not share the same period as the modulation of the photospheric lines. Nevertheless it is possible in principle for the warped structure that we have invoked to explain the modulation of the photospheric lines to also have some influence on modulating the mass uptake in the wind (note the similarity in spatial scales between the inferred radius of a Keplerian hotspot and the vertical height in the wind from which we infer that the variable wind signatures originate). This would require that the warped structure, in addition to rigid rotation at the orbital period of the putative hotspot, also underwent low-amplitude tilt motions normal to the disc plane. Clearly we have no detailed calculations to support this suggestion but merely point out that this could represent an economical solution to the 'two period' problem.

## ACKNOWLEDGMENTS

We thank Nuria Calvet, Giuseppe Lodato and Marina Rominova for helpful discussions regarding physical interpretations of the observed cyclic variability, as well as the referee for the helpful comments. SLP thanks the Science and Technology Facilities Council (STFC) for the award of a PhD studentship and acknowledges funding by the European Commission FP6 Marie Curie Research Training Network, CONSTELLATION, and Churchill College, Cambridge, for travel grants. SLP would also like to thank James Owen for invaluable technical support.

## REFERENCES

* Ambartsumian V. A., 1971, Astrofizika, 7, 557

* Ambartsumian V. A., 1971, Astrophys., 7, 331

* Bastian U., Mundt R., 1985, A&A, 144, 57

* Battaglia G., Helmi A., Tolstoy E., Irwin M., Hill V., Jablonka P., 2008, ApJ, 681, L13

* Bell K. R., Lin D. N. C., Hartmann L. W., Kenyon S. J., 1995, ApJ, 444, 376

* Blandford R. D., Payne D. G., 1982, MNRAS, 199, 883

* Calvet N., Hartmann L., Kenyon S. J., 1993, ApJ, 402, 623

* Chochol D., Tremko J., 1980, in Mirzoyan L. V., ed., Vsesoiuznaia Konferentsiia Kosmicheskim Lucham Yerevan Armenian SSR Aademiia Nauk SSSR Izvestiia Seriia Fizicheskaia, p. 240

* Clarke C. J., Armitage P. J., 2003, MNRAS, 345, 691

* Clarke C. J., Syer D., 1996, MNRAS, 278, L23

* Clarke C. J., Lin D. N. C., Pringle J. E., 1990, MNRAS, 242, 439

* Croswell K., Hartmann L., Avrett E. H., 1987, ApJ, 312, 227

* D'Angelo G., Errico L., Gomez M. T., Saldone L. A., Teodorani M., Vittone A. A., 2000, A&A, 356, 888

* Donati J.-F., Paletou F., Bouvier J., Ferreira J., 2005, Nat, 438, 466

* Errico L., Vittone A., Lamzin S. A., 2003, Astron. Lett., 29, 105

* Evans N. J., II, Balkum S., Levreault R. M., Hartmann L., Kenyon S., 1994, ApJ, 424, 793

* Foucart F., Lai D., 2011, MNRAS, 412, 2799

* Goodrich R. W., 1987, PASP, 99, 116

* Hartmann L., 1994, in Clemens D. P., Barvainis R., eds, ASP Conf. Ser. Vol. 65, Clouds, Cores, and Low Mass Stars. Astron. Soc. Pac., San Francisco, p. 315

* Hartmann L., 1998, Accretion Processes in Star Formation. Cambridge Univ. Press, Cambridge

* Hartmann L., Calvet N., 1995, AJ, 109, 1846

* Hartmann L., Kenyon S. J., 1985, ApJ, 299, 462

* Hartmann L., Kenyon S. J., 1987a, ApJ, 312, 243

* Hartmann L., Kenyon S. J., 1987b, ApJ, 322, 393

* Hartmann L., Kenyon S. J., 1996, ARA&A, 34, 207

* Hartmann L., Kenyon S., Hartigan P., 1993, in Levy E. H., Lunine J. I., eds, Protostars and Planets III Young Stars. Univ. Arizona Press, Tucson, p. 497

* Herbig G. H., 1958, ApJ, 128, 259

* Herbig G. H., 1966, Vistas Astron., 8, 109

* Herbig G. H., 1977, ApJ, 217, 693

* Herbig G. H., 1989, in Reipurth B., ed., ESO Workshop on Low-Mass Star Formation and Pre-Main-Sequence Object. ESO, Garching, p. 233

* Herbig G. H., Petrov P. P., Duemmler R., 2003, ApJ, 595, 384

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

* Ibragimov M. A., 1993, Astron. Zh., 70, 339

* Kenyon S. J., 1995, Rev. Mex. Astron. Astrofis., Conf. Ser., 1, 237

* Kenyon S. J., Hartmann L., Hewett R., 1988, ApJ, 325, 231

* Kenyon S. J., Kolotilov E. A., Ibragimov M. A., Mattei J. A., 2000, ApJ, 531, 1028

* Kolotilov E. A., Petrov P. P., 1985, Pis'ma Astrono. Zh., 11, 846

* Königl A., Romanova M. M., Lovelace R. V. E., 2011, MNRAS, 416, 757

* Lai D., 1999, ApJ, 524, 1030

* Lai D., 2003, ApJ, 591, L119

* Larson R. B., 1983, Rev. Mex. Astron. Astrofis., 7, 219

* Lin D. N. C., Papaloizou J., 1985, in Black D. C., Matthews M. S., eds, Protostars and Planets II. p. 981

* Linnell Nemec A. F., Nemec J. M., 1985, AJ, 90, 2317

* Lodato G., Clarke C. J., 2004, MNRAS, 353, 841

* Mould J. R., Hall D. N. B., Ridgway S. T., Hintzen P., Aaronson M., 1978, ApJ, 222, L123

* Popham R., Kenyon S., Hartmann L., Narayan R., 1996, ApJ, 473, 422

* Pringle J. E., 1976, MNRAS, 177, 65

* Reipurth B., 1990, in Mirzoian L. V., Pettersen B. R., Tsvetkov M. K., eds, Proc. IAU Symp. 137, Flare Stars in Star Clusters, Associations and the Solar Vicinity. Kluwer, Dordrecht, p. 229

* Scargle J. D., 1982, ApJ, 263, 835

* Shu F., Najita J., Ostriker E., Wilkin F., Ruden S., Lizano S., 1994, ApJ, 429, 781

* Tonry J., Davis M., 1979, AJ, 84, 1511

* Vittone A. A., Errico L., 2005, Mem. Soc. Astron. Ital., 76, 320

* Wang H., Apai D., Henning T., Pascucci I., 2004, ApJ, 601, L83

* Welty A. D., Strom S. E., Edwards S., Kenyon S. J., Hartmann L. W., 1992, ApJ, 397, 260

* Wilson R. E., 1953, General Catalogue of Stellar Radial Velocities. Carnegie Institute, Washington, D.C. Publ., 601

* Zhu Z., Espaillat C., Hinkle K., Hernandez J., Hartmann L., Calvet N., 2009a, ApJ, 694, L64

* Zhu Z., Hartmann L., Gammie C., McKinney J. C., 2009b, ApJ, 701, 620

This paper has been typeset from a TEX/LATEX file prepared by the author.