Abstract

The detection of exoplanets using any method is prone to confusion due to the intrinsic variability of the host star. We investigate the effect of cool starspots on the detectability of the exoplanets around solar-like stars using the radial velocity method. For investigating this activity-caused ‘jitter’ we calculate synthetic spectra using radiative transfer, known stellar atomic and molecular lines, different surface spot configurations and an added planetary signal. Here, the methods are described in detail, tested and compared to previously published studies. The methods are also applied to investigate the activity jitter in old and young solar-like stars, and over a solar-like activity cycles. We find that the mean full jitter amplitude obtained from the spot surfaces mimicking the solar activity varies during the cycle approximately between 1 and 9 m s−1. With a realistic observing frequency a Neptune-mass planet on a 1-yr orbit can be reliably recovered. On the other hand, the recovery of an Earth-mass planet on a similar orbit is not feasible with high significance. The methods developed in this study have a great potential for doing statistical studies of planet detectability, and also for investigating the effect of stellar activity on recovered planetary parameters.

1 INTRODUCTION

The search for exoplanets has traditionally concentrated on stars with very little intrinsic activity. Studies have shown that the known exoplanet host stars exhibit very low levels of magnetic activity (e.g. Jenkins et al. 2006; Martínez-Arnáiz et al. 2010). Still, as the Kepler satellite has shown, many solar-like stars are more active than our Sun (e.g. Basri et al. 2013), and therefore show significant levels of activity, which can affect the planet detection. The spectral line profile variations caused by starspots have been confused with a radial velocity signal originating from exoplanets (e.g. Queloz et al. 2001; Huerta et al. 2008), and the sudden brightenings caused by stellar flares can mimic microlensing events from a planet-sized body (e.g. Bennett et al. 2012). Magnetic activity, and the phenomena related to it, is an integral part of stars with spectral types ranging from mid-F to M, and as these are the stars most exoplanet searches concentrate on, it is crucial to understand the effects activity sets on exoplanet detection and parameter determination.

The exoplanet detection method most prone to confusion from stellar activity is the radial velocity search. Already, Saar & Donahue (1997) showed that cool spots on the stellar surface cause spectral line profile variations that can be confused with the radial velocity variations from planets. In addition, they derived an analytical formula to represent the relation of this radial velocity jitter to the fraction of stellar surface covered by starspots (starspot filling factor). Other similar investigations have been carried out, see e.g. Desort et al. (2007), Reiners et al. (2010), Dumusque et al. (2011b), Boisse et al. (2011) and Barnes, Jeffers & Jones (2011). Also, Boisse, Bonfils & Santos (2012) provided a freely available tool, Spot Oscillation And Planet (soap), for the community to investigate the effects of starspots on radial velocity measurements and photometric observations. The code was expanded and improved by Oshagh et al. (2013) to include planetary transits on a spotted star, and a new modified version, soap 2.0, was recently published by Dumusque, Boisse & Santos (2014).

Most planet searches have concentrated on solar-like stars. Therefore, there have also been relatively many investigations on the effect of solar activity on exoplanet detection (see e.g. Lagrange, Desort & Meunier 2010; Meunier, Desort & Lagrange 2010). Meunier & Lagrange (2013) used solar activity as a template to study the effect of spots and plages on detectability of Earth-mass exoplanets in the habitable zone of their host star. They conclude that especially the contribution from the plages would prevent the detection of Earth around the Sun, even with forthcoming high-precision instruments. These investigations have concentrated on plages and spots, but these are not the only error sources. Granulation and stellar oscillations also cause noise at m s−1 level (see e.g. Dumusque et al. 2011a), which will hinder the detection of small-sized planets, and planets on wide orbits. Still, the time-scale for these noise sources is much shorter than the variations caused by long-period planets, therefore, they can be averaged out using long exposures and/or observing frequently. A new potential noise source, gravitational redshift, was identified recently by Cegla et al. (2012), but its magnitude is estimated to be only few cm s−1.

Photometric observations have been used to estimate the activity-induced jitter in radial velocity measurements (e.g. Lanza et al. 2011; Aigrain, Pont & Zucker 2012). In a recent work based on Galaxy Evolution Explorer (GALEX) ultraviolet measurements and Kepler light curves, Cegla et al. (2014) investigate how well the radial velocity jitter can be estimated based on photometry alone. They conclude that for magnetically quiet stars one can use photometric measurements as a proxy for radial velocity variability.

Ma & Ge (2012) developed a technique to use the radial velocity method for detecting planets even when the host star shows significant activity. Their method relies on the wavelength dependence of the spot-caused jitter as opposed to the planetary signature which is not wavelength dependent. Moulds et al. (2013) have also shown that for active stars it is possible to remove some of the spot signature from the line profiles and still be able to recover Jupiter-sized close-in planets.

In this work we develop methods that use radiative transfer to calculate spectral line profiles from a spotted stellar surface. The stellar spot configurations are either based on spot sizes and numbers or filling factors. It is also possible to introduce active longitudes and latitudes. In addition, a planetary signal can be added to the spectra calculated from the spotted surface. This approach allows for a statistical investigation of exoplanet detection and enables obtaining information on the errors the stellar activity causes on the determined planetary parameters. In the current paper, Paper I, we describe the methods, test them and compare the results to some of the previously published studies. The methods are also applied to study jitter around solar-like stars, including a solar-like activity cycle. In the second part, Paper II, we will apply the methods to M dwarfs and investigate the reliability of recovery of planetary parameters in the presence of stellar activity (Andersen & Korhonen 2014).

2 METHODS

We create spot patterns on a simulated stellar surface using our spotss code (see next section). Synthetic spectral line profiles at different rotational phases of the star are calculated based on the created surface distributions. The spectra based on the given spot configuration and local line profile grids are calculated using the code direct7, which is written by NP and includes modifications introduced by TH. This code uses the same routines as invers7 (Piskunov, Tuominen & Vilhu 1990; Hackman, Jetsu & Tuominen 2001). Radial velocity (jitter) measurements are obtained by cross-correlation of the calculated line profiles with either one of the generated profiles or a template profile, which is obtained from a spectrum with the temperature of the unspotted stellar surface. Exoplanet radial velocity is introduced into the spectra and recovered by a new cross-correlation. These steps after the spot surface creation are all done in a code called Detection of Exoplanets under the Effect of Magnetic Activity (deema). Details of all the steps are discussed in the following.

2.1 Generating spots

We developed our code, spotss, to generate spot (temperature) patterns on a simulated stellar surface. The code can generate random spots across the entire surface, and has also the option to define certain ‘active’ regions in latitude and/or longitude. The code creates a matrix of temperatures which represents the entire stellar surface, corresponding to N points of latitude and 2 × N points of longitude. Using a 60 × 120 grid yields a latitude resolution of 3° pixel−1 and a longitude resolution of 3° pixel−1 at the equator (note that the actual longitude resolution increases near the poles since the grid is a flat square that represents a spherical surface). The size of the matrix can be changed in order to increase or decrease the spatial resolution, but at much higher sizes the processing time is very long compared to the relatively small gains in precision.

The code takes the following (user defined) input parameters.

  • Stellar temperature, Tp: the photospheric temperature of the unspotted stellar surface, before spots are added. Every value in the temperature matrix is originally set to this value.

  • Spot temperature, Ts: the temperature of the umbral regions of the spots. This is considered the spot temperature, although each spot also has a penumbral region with a temperature defined as the midpoint between the spot temperature and the photospheric temperature.

  • Spot radius, rs: ‘average’ radius of the spots. The radius of each spot is randomly altered starting from this value to create a lognormal distribution of spot sizes around this size. Lognormal distribution is chosen because the sunspot size distribution is known to follow it (e.g. Bogdan et al. 1988; Baumann & Solanki 2005). Spots are approximately circular, and we add penumbral regions with umbral to penumbral radii ratios of 1:2, creating an umbral to penumbral area of 1:3, following Solanki (2003).

  • Number of spots, ns (or filling factor, depending on which version of the code is used): spots are placed on the stellar surface randomly until the number of spots is reached or the desired filling factor is achieved. For certain purposes it was useful to investigate the effect of one large spot, so ns was set to 1. (Note when ns = 1 the exact input value of rs is used, since it is not necessary to alter this to obtain a certain size distribution of spots. An exact latitude and longitude for the centre of the spot can also be specified.)

  • Longitude range: defaults to the full range of longitude: 0°–360°, but active longitude ranges can be defined meaning spots will only be placed within those ranges.

  • Latitude range: defaults to the full range of latitude, −90°–90°, though active latitudes can also be defined. Then spots will only be placed within those ranges. This can be combined with defining an active longitude range to create a small ‘active region’ on the stellar surface. Active regions such as this have been observed in some Doppler images of active stars.

2.2 Calculating spectral line profiles

For calculating synthetic spectra we used two different atmospheric models: one for the solar-like stars and one for the M dwarfs. For both cases 17 limb angles (0.01, 0.025, 0.050, 0.075, 0.1, 0.125, 0.15, 0.2, 0.25, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0) were used. For solar-like stars a grid of local line profiles was created using the spectrum spectral synthesis code (Gray & Corbally 1994) and Kurucz model atmospheres (Kurucz 1993). The grid includes temperatures 4000–6000 K with a temperature step of 250 K. The line lists luke.lst and luke.nir.lst (included in the spectrum package), which include atomic lines and some molecular species, are used in the calculations. For the photometry we use the same code and models but for a sparser wavelength grid ranging from 3600 to 7350 Å and step size of 50 Å.

For the M dwarf temperatures we used the program synthe (Kurucz & Avrett 1981; Sbordone et al. 2004) to compute the synthetic spectra. The spectra are based on a subset of the 2008 grid of marcs stellar model atmospheres (Gustafsson et al. 2008) with solar metallicities. They cover a temperature range of T = 2500–4000 K, a surface gravity range of log(g) = 4.5–5.5 and a microturbulence range of ξt = 0.0–2.0 m s−1. For the atomic line opacities we used data from the on-line Vienna Atomic Line Database 2 (VALD-2; Kupka et al. 1999), and for the molecular line opacities we used the species presented in Table 1. The calculated spectra cover a wavelength range of 3000–9200 Å and have a resolution of Δλ/λ = 500 000.

Table 1.

Sources of data for the molecular line opacities used in synthe for the M dwarf spectral grid.

MoleculeTransitionsWavelength range (Å)# LinesSource
CHA–X, B–X, C–X2600–17 00071 600R. Kurucz
CNA–X, B–X2000–1000 0001645 000R. Kurucz
COX–X, A–X1100–100 000555 000R. Kurucz
OHX–X, A–X2000–1000 00082 000R. Kurucz
H2OX–X4100–1000 00065 900 000R. Kurucz
SiOX–X, A–X, E–X1400–1000 0001830 000R. Kurucz
TiOA–X, B–X, C–X, E–X, c–a, b–a, b–d, f–a4100–1000 00033 000 000D. W. Schwenke
MoleculeTransitionsWavelength range (Å)# LinesSource
CHA–X, B–X, C–X2600–17 00071 600R. Kurucz
CNA–X, B–X2000–1000 0001645 000R. Kurucz
COX–X, A–X1100–100 000555 000R. Kurucz
OHX–X, A–X2000–1000 00082 000R. Kurucz
H2OX–X4100–1000 00065 900 000R. Kurucz
SiOX–X, A–X, E–X1400–1000 0001830 000R. Kurucz
TiOA–X, B–X, C–X, E–X, c–a, b–a, b–d, f–a4100–1000 00033 000 000D. W. Schwenke
Table 1.

Sources of data for the molecular line opacities used in synthe for the M dwarf spectral grid.

MoleculeTransitionsWavelength range (Å)# LinesSource
CHA–X, B–X, C–X2600–17 00071 600R. Kurucz
CNA–X, B–X2000–1000 0001645 000R. Kurucz
COX–X, A–X1100–100 000555 000R. Kurucz
OHX–X, A–X2000–1000 00082 000R. Kurucz
H2OX–X4100–1000 00065 900 000R. Kurucz
SiOX–X, A–X, E–X1400–1000 0001830 000R. Kurucz
TiOA–X, B–X, C–X, E–X, c–a, b–a, b–d, f–a4100–1000 00033 000 000D. W. Schwenke
MoleculeTransitionsWavelength range (Å)# LinesSource
CHA–X, B–X, C–X2600–17 00071 600R. Kurucz
CNA–X, B–X2000–1000 0001645 000R. Kurucz
COX–X, A–X1100–100 000555 000R. Kurucz
OHX–X, A–X2000–1000 00082 000R. Kurucz
H2OX–X4100–1000 00065 900 000R. Kurucz
SiOX–X, A–X, E–X1400–1000 0001830 000R. Kurucz
TiOA–X, B–X, C–X, E–X, c–a, b–a, b–d, f–a4100–1000 00033 000 000D. W. Schwenke

The spectra based on the given spot configuration and local line profile grids are calculated using the code direct7. Before running direct7 the local line profiles are convolved with a Gaussian instrumental profile and a radial-tangential macroturbulence (here set to 2 km s−1). Furthermore, local fluxes for the B- and/or V magnitudes are calculated using the transmission functions for these wavelength passbands.

For each rotation phase the visible stellar hemisphere is simulated in direct7 and the spectrum and B- and V-passband fluxes are calculated by integration of the line and continuum intensities over this visible stellar hemisphere. The projection effects, limb darkening and visibility of the spots are thus naturally directly taken into account. User defined v sin i and inclination values are applied to the spectral line profiles. The v sin i can be given in two ways: either the user can supply the desired value, or the deema code can estimate the value based on the stellar mass, rotation period and inclination. For the estimation models by Baraffe et al. (1998) are used. From now on the term ‘spectrum’ refers to the synthetic spectrum that has been calculated by integration over the full visible stellar disc.

The spectra can be calculated at any given rotational phase. Here, a scheme where the length of the observing run and stellar rotation period are given in days, and the user supplied number of phases is evenly distributed over the observing run and the stellar rotational phases are calculated based on this information. The stellar rotational phases are used as an input for the line profile and B- and/or V magnitude calculations.

2.3 Cross-correlation

For obtaining the jitter values induced by spots, spectra calculated for different observational phases are cross-correlated. The cross-correlation can be carried out against different templates: spectra without spots obtained using unspotted surface temperature, spectra without spots obtained using the mean temperature of the stellar surface or using one of the spectra with spots as the template.

If one of the spectra created with surface spots is used as the template, then all the spectra are cross-correlated against each other. This means that there are N − 1 jitter curves, where N is the number of rotational phases. All the obtained jitter curves are normalized in such a way that the measurement at the first phase is set to zero, and the mean and the standard deviation are calculated for each phase. This provides a mean jitter curve and an estimate of the measurement error. The jitter curves are virtually identical and the error estimates small, as should be when no noise is added to the spectra. This is the method that is used in all the calculations presented in this work, and owing to this scheme the cross-correlation result is always zero for the first longitude.

In cross-correlation the accuracy is increased by improving the sampling of the cross-correlation curve using linear interpolation and fitting a polynomial to the curve. The maximum of the polynomial is calculated and the corresponding shift used as the jitter value. An example of a stellar surface and resulting jitter curve are shown in Figs 1(a) and (b), respectively. The spot is a 5° radius equatorial spot with umbral temperature of 4000 K, penumbral temperature of 4900 K and unspotted surface temperature of 5800 K. The stellar rotation period has been set to 25 d and the spectra have been calculated at 20 evenly spaced phases using the wavelength region 5952–5998 Å. The inclination is set to 90° (equator-on) and v sin i is estimated from the period and stellar radius to be 2.1 km s−1. The resulting jitter curve has a full amplitude of 14.95 m s−1.

Figure 1.

An example of stellar surface spot configuration and the resulting radial velocity jitter from a spot and a planet. (a) Stellar surface configuration used for calculating the spectra. The x-axis gives the longitude in degrees and y-axis the latitude, also in degrees. The unspotted surface has a temperature of 5800 K and the spot temperature is 4000 K in the umbra and 4900 K in the penumbra. The radius of the spot is ∼5°. (b) The resulting jitter curve from the 5°spot. The x-axis is the longitude in degrees and y-axis the measured radial velocity jitter in m s−1. The jitter values are calculated at 20 evenly spaced rotational phases from the wavelength region 5952–5998 Å. The error of the jitter measurement is smaller than the symbol size. (c) Radial velocity curve of a Neptune-mass planet on a 25 d circular orbit around 1 M star. The x-axis gives the time in days and y-axis the radial velocity in m s−1. (d) The radial velocity curve from the 5°spot and Neptune-sized planet. Owing to the cross-correlation scheme used here (see Section 2.3) the first jitter measurement is always shifted to zero. This explains the different absolute values for the input and calculated radial velocities. In the plot x-axis has the time in days and y-axis the measured radial velocity in m s−1.

Tests were also carried out using a Gaussian function instead of a second degree polynomial for determining the maximum of the cross-correlation function. Fig. 2 shows the results for three different methods: polynomial (plus signs), Gaussian function (squares) and simply taking the maximum of the cross-correlation function (diamonds). The results are shown for two different equatorial spots, one with radius of 5°(top) and one of 17°(bottom). As expected, simply taking the maximum of the cross-correlation function gives zero when the shift is small, and with larger shifts it at times gives exaggeratedly large values. On the other hand, both polynomial and Gaussian fits give very similar results. Based on calculations using 100 different spot configurations with spot filling factor of 0.02 per cent, it can be seen that the polynomial fit always gives somewhat smaller value than the Gaussian fit. The full jitter amplitude obtained from the polynomial fit is 2.85 ± 0.03 per cent smaller than the value obtained using the Gaussian fit. The tendency for smaller jitter values with the polynomial fit is also seen in Fig. 2, but for both the 5°and 17°radius spots the difference between the methods is less than 1 per cent.

Figure 2.

An example of jitter curves obtained using different methods for determining the maximum of the cross-correlation function: polynomial fit (plus signs), Gaussian function (squares) and taking the maximum of the cross-correlation function (diamonds). The results are shown for two different equatorial spots, one with radius of 5°(top) and one of 17°(bottom).

There does not seem to be any correlation between the full jitter amplitude and the performance of the Gaussian and polynomial fits. Fig. 3 shows the procentual difference between the full jitter amplitude obtained using polynomial and Gaussian fits and plotted against the full jitter amplitude obtained from the Gaussian fit. The larger jitter does not result in larger (or smaller) procentual difference between the methods. On the whole, both Gaussian and polynomial fits give very similar results. In the following polynomial fits are used.

Figure 3.

The procentual difference between the full jitter amplitude obtained using polynomial and Gaussian fits plotted against the full jitter amplitude from the Gaussian fit.

direct7 uses an evenly spaced wavelength grid for calculating the spectra. Our tests show that over such a small wavelength range (50 Å) the results are not significantly affected even if a logarithmic wavelength scale is not used in the cross-correlation. Therefore, we use the direct7 output directly in the cross-correlation.

2.4 Introducing a planet to the spectra

We generate radial velocity curves resulting from orbiting planets using Kepler's third law. The elliptical case of Kepler's equation is solved following the formalism by Mikkola (1987). The three optional variables (γ, the systemic velocity, or arbitrary instrumental offset; |$\dot{\gamma }$|⁠, the systemic acceleration, due to systematics in the data or an additional body in the system with a much longer period and t0, an arbitrary zero-point for the slope) were left out. These parameters are used when fitting radial velocity curves, but they are not necessary when simply generating a radial velocity curve. This leaves the equation
\begin{equation} \text{RV}(t) = K[\cos \theta (t)+ \omega _{*} + e\cos \omega _{*}], \end{equation}
(1)
where K is the radial velocity semiamplitude and is given in m s−1 by
\begin{equation} K = \left(\frac{2\pi G}{P(M_{\text{P}}+M_{\text{S}})^2}\right)^{1/3}\frac{M_{\text{P}}\sin i}{\sqrt{1-e^2}}. \end{equation}
(2)
In these equations θ is the true anomaly, ω* is the argument of periastron, e is eccentricity of the orbit, P the orbital period of the planet in seconds, MP the mass of the planet in kg and MS the stellar mass in kg.

The radial velocity at a given orbital phase is calculated and added to the appropriate spectrum. For this process the user has to provide the stellar mass, planetary mass, eccentricity and period of the planetary orbit (or semimajor axis). As an example, Fig. 1(c) shows the calculated radial velocity curve for a Neptune-mass planet (17 M) around a solar-mass star on a circular orbit with orbital period of 25 d. The radial velocities caused by the planet are calculated at the input rotational phases and shifts are introduced to the spectra. After this the spectra are re-analysed using the same cross-correlation method as for the case only containing the spot jitter. The resulting radial velocity curve of the spot and planet together is shown in Fig. 1(d). For this example, the full amplitude of the radial velocity variation is the same as from only the spot, 14.95 m s−1, but the shape of the curve is very different.

2.5 Testing the methods

For testing the behaviour of the code a 5°radius equatorial solar-like spot was used (see Fig. 1a). The radius of the whole spot (umbra+penumbra) is 5°, the radius of the umbra alone is 3°. In some tests also a larger spot with the full radius of 17° and umbral radius of 10° is used. Both spots, which are used separately, are located at the equator and the temperature of the umbra is 4000 K, the penumbra is 4900 K and the unspotted temperature is 5800 K. All the tests were carried out using the grid size 60 × 120 (except the grid size tests). The spectra were created using spectral resolution of 100 000 (except in the resolution tests) and using 2.5 pixel sampling over one resolution element. The length of the wavelength strip used in one jitter calculation was always 46 Å. The local line profiles were calculated with ±2 Å from the ends of the wavelength strip, to allow for large v sin i values in the calculations (i.e. the length of the spectrum for which local line profiles were calculated was always 50 Å).

2.5.1 Recovering the planetary signal

The accuracy at which the planetary signal can be recovered was tested using a smooth, i.e. unspotted, stellar surface. Solar-like configuration with effective temperature of 5800 K and v sin i = 2.1 were used in the calculations. An Earth-mass planet on a 25 d orbit around the 1 M star was introduced to the spectra.

The results of the test can be seen in Fig. 4. The jitter curve calculated from the unspotted stellar surface (Fig. 4b) shows zero values for all the jitter measurements, as expected. The input curve of the Earth-mass planet (Fig. 4c) on the other hand is very well recovered from the spectra (Fig. 4d). As discussed in Section 2.3, the radial velocity curve calculated from the spectra is normalized to the first measurement, explaining the different absolute values for the input and calculated radial velocities. The full amplitude of the jitter is almost identical in the full amplitude of the input radial velocity curve (0.4373 m s−1) and the calculated one (0.4299 m s−1). Our tests using other planetary masses show that the calculated full amplitude tends to be underestimated by ∼1.7 per cent in comparison to the real input radial velocity curves full amplitude.

Figure 4.

An example of recovering the planetary signal from the spectra. (a) Stellar surface configuration without spots and surface temperature of 5800 K. This configuration was used for calculating the spectra for the planet recovery test. The x-axis gives the longitude in degree and y-axis the latitude, also in degrees. (b) The resulting jitter curve from the unspotted surface. The x-axis is the longitude in degrees and y-axis the measured radial velocity jitter in m s−1. The jitter values are calculated at 20 evenly spaced rotational phases on the wavelength region 5952–5998 Å. (c) Radial velocity curve of an Earth-mass planet on a 25 d circular orbit around 1 M star. The x-axis gives the time in days and y-axis the radial velocity in m s−1. (d) The radial velocity curve from the unspotted surface and an Earth-sized planet. The x-axis gives the time in days and y-axis the measured radial velocity in m s−1.

2.5.2 Effect of spectral resolution

Radial velocity measurements are typically done using high-resolution spectrographs with resolving power (λ/Δλ) 50 000–110 000. To test the effect of the spectral resolution on the jitter we calculated spectra using different resolving powers between 10 000 and 300 000 and the same spot configuration as shown in Fig. 1(a). One resolution element always spans 2.5 wavelength steps, i.e. ‘pixels’.

As can be seen from Fig. 5 the full amplitude of the jitter decreases with decreasing spectrograph resolution, as has also been reported by other authors (e.g. Desort et al. 2007). The results are normalized to the highest jitter, i.e. resolution 200 000 result of 15.8 m s−1. Two different regimes can be seen in the jitter behaviour. Throughout the resolution range normally used for exoplanet searches, i.e. 50 000–300 000, the jitter remains at the level of about 90 per cent of the highest jitter value. The smallest jitter, about 77 per cent of the highest values, is seen at lowest spectral resolutions used in this test.

Figure 5.

The spot-caused jitter with different spectrograph resolving power (λ/Δλ). The results are normalized to the maximum jitter case (i.e. resolving power of 150 000). The spot configuration used is the same as in Fig. 1(a). In the plot the x-axis gives the resolving power and y-axis the jitter in m s−1. The jitter increases with increasing resolving power.

This behaviour can be explained by the spot contribution on the line profile shape being better resolved at high spectral resolution and getting more and more diluted with decreasing spectral resolution. Still, the effect of the spectrograph resolution in the jitter amplitude is not very strong and the accuracy of the radial velocity measurements decreases with the decreasing spectral resolution.

2.5.3 Effect of the width of the spectral region

In real high-precision spectroscopic observations the whole optical wavelength range is typically used for determining the radial velocity. Still, as there is no noise added to the spectra the increase in the width of the spectral range should not have a major influence on the results, if wide enough wavelength region is used. To test this assumption we have calculated jitter from seven different wavelength ranges, with the width spanning from 10 to 70 Å. The test uses wavelengths between 5925 and 6000 Å, and spectral resolution of 100 000.

Fig. 6 shows the normalized jitter for two different spot configurations: 5°spot (plus signs) and 17°spot (squares). As expected, the jitter is highest for the 10 Å wide wavelength range, and decreases slightly when going to the 40 Å wide wavelength region. The wider wavelength regions, from 40 to 70 Å, all show similarly small jitter, and the results from more narrow spectral ranges have 5–10 per cent larger jitter values. Still, we cannot say based on this test, whether the jitter would decrease further if significantly wider wavelength ranges would be used, but for limiting the calculation times, 46 Å wide wavelength region is used throughout this paper.

Figure 6.

The spot-caused jitter with different widths of the spectral region. Two different spot configurations were used for this test: the same as in Fig. 1(a) shown by plus signs, and one with the full spot radius of 17° denoted by squares. The results for both spot configurations are normalized individually to the results from the widest spectral region (70 Å case). In the plot the x-axis gives the width of the spectral region in Ånsgtröm and y-axis the jitter in m s−1.

2.5.4 Size of the spatial grid

Testing which impact the grid size, i.e. spatial resolution on the stellar surface, has on the jitter is difficult. The location and fractional size of the spots have to be kept identical throughout the test, which is of course strictly speaking impossible to do.

For the test a random spot configuration with a grid size 20 × 40 was created. The input map has a spot filling factor of 1.5 per cent and the spots were occurring at latitude range −30° to +30°. The grid resolution was increased and filling factor kept as close to the original as possible with the increasing number of grid elements. The scalings did not have a significant effect on the jitter and the filling factor. The jitter values were between 24.6 and 29.0 m s−1 with all the grid sizes. These tests imply that the results are not critically dependent on the grid size, but to allow for also small spots on the surface we have used the grid size 60 × 120 throughout this paper.

2.5.5 Signal-to-noise ratio of the observations

To investigate the effect of noise on the radial velocity curves we add Gaussian noise of a specified signal-to-noise ratio (S/N) to the spectra. The level of the continuum represents signal level while the standard deviation of the noise represents the noise level. Since the continuum is normalized |$\text{S/N} = \frac{\text{signal level}}{\text{noise level}} = \frac{1}{\sigma }$|⁠, leading to |$\sigma = \frac{1}{\text{S/N}}$|⁠. A sequence of pseudo-random numbers are generated using idl routine randomn, which creates Gaussian random numbers using Box–Muller method. The noise is then applied to the spectra.

For testing the effect of noise in the jitter curves we use spectra calculated for a single equatorial spot with a full radius of 17°. Jitter curves are calculated for seven different cases, one with no noise and six with different S/N values ranging from 20 to 3000. The resulting curves are shown in Fig. 7. The symbols used in the plot also show the errors of the individual jitter values (described in Section 2.3). The shape of the jitter curve is easily recognizable until S/N less than 100. The general shape can still be recovered from the spectra with S/N = 50, but with S/N = 20 the shape becomes basically unrecognizable.

Figure 7.

Jitter curves caused by a 17° radius spot with different S/N of the spectra. In the plot the x-axis gives the stellar longitude in degrees and y-axis the jitter in m s−1. The topmost jitter curve is from the case without noise, and the other jitter curves are off-set from this curve to show the different behaviour better. The S/N values are written on the plot for each jitter curve. The plot also shows the error of each jitter value. The error has been calculated as described in Section 2.3.

We want to still note that each spectral line represents an individual measurement of the Doppler shift of the star. If a total of N lines are used for the Doppler measurement, then the error will be decreased by a factor of |$\sqrt{N}$| over a single line measurement. In the tests carried out here, short wavelength ranges of 50 Å are used. In the solar-like case these wavelength regions have 30–100 spectral lines, the exact number depending on the wavelength (typically more lines in the blue part and less in the red). Planet searches on the other hand use echelle spectra with a few thousand lines in them. Therefore, if our tests only have at most a 10 fold gain over single line measurements, the real planet searches usually have approximately 40 fold gain. For this reason in the following investigations we will use the calculated spectra without added noise. This will enable us also to study the ideal detection cases.

3 RESULTS

We apply the developed codes to study the starspot jitter in solar-like stars. First, general properties of jitter are studied and compared to the earlier published results by other groups. Afterwards, two different activity cases are investigated: solar-like low activity and very active young solar analogues. In these investigations solar-like temperatures (unspotted surface 5800 K, umbra 4000 K and penumbra 4900 K) are used. If not mentioned otherwise, the gird size is 60 × 120, the wavelength region is 5952–5998 Å and the spectrograph resolving power (λ/Δλ) is 100 000 with 2.5 pixel sampling. Similarly, the inclination was set to 90° and v sin i was fixed to 2.1 km s−1 (if not mentioned otherwise).

3.1 General properties

3.1.1 Jitter with wavelength

Several works have already shown that the spot caused jitter decreases with increasing wavelength (e.g. Desort et al. 2007; Reiners et al. 2010). In a recent paper Marchwinski et al. (2015) show that, based on solar observations, near-infrared has lower estimated radial velocity jitter throughout the entire solar cycle than the optical wavelengths have. Here tests using spectra with lengths of 46 Å at different wavelength regions between 3700 and 9100 Å are carried out. The spectrograph resolving power (λ/Δλ) is kept constant at 100 000, by changing the Δλ according to the wavelength region. Two spot sizes, 5° and 17° full spot radius (umbra+penumbra), are used. The radii of the umbra are 3° and 10°, respectively.

The results of the jitter investigation over the wavelength are shown in Fig. 8. The results are normalized to the highest jitter case, that of 3700 Å. The results for the 5° spot are shown by plus signs and for the 17° spots by squares. As can be seen, the normalized values and the overall behaviour are very similar for both spot sizes. The absolute values are naturally very different: the full amplitude of the jitter for the 5° spot at 3702–3748 Å is 23.7 m s−1 and for the 17° spot 213.7 m s−1.

Figure 8.

Full amplitude of the spot-caused jitter at different wavelengths. Two different spot configurations were used for this test: the same as in Fig. 1(a) shown by plus signs, and one with the full spot radius of 17° denoted by squares. The results are normalized to the highest jitter case of 3700 Å. In the plot the x-axis gives the wavelength in Ångström and y-axis the normalized jitter. As has been seen in other studies too, the jitter decreases towards the longer wavelengths, but also the scatter in the jitter values increases.

The jitter decreases by about 50 per cent between the wavelengths 3700 and 7000 Å. At longer wavelengths than this the reduction in jitter amplitude is not as clear as before, and the whole behaviour becomes more chaotic. The increased scatter at red wavelengths could possibly be due to varying number of spectral lines in the spectral windows used in this analysis. In general, there are less spectral lines at the red wavelengths and therefore different red regions can have very different number of spectral lines. Investigation extending to longer wavelengths would be needed to study whether or not the decrease continues to infrared. A similar plateau in the wavelength dependence of the jitter around 8000–10 000 Å is seen for solar-like stars by Reiners et al. (2010). In their work some further decrease in the jitter amplitude is seen in the infrared wavelengths, as is also detected by Marchwinski et al. (2015).

3.1.2 Jitter with v sin i

It has been shown before that the jitter depends strongly on the stellar rotation rate (e.g. Saar & Donahue 1997; Desort et al. 2007; Boisse et al. 2012). Also in this work the effect was studied and compared to the previously published results. Two spot cases, with the full spot radius of 5° and 17°, were used in the investigation.

Fig. 9 shows the results from the two spot cases for v sin i values ranging between 1 and 30 km s−1. The results are normalized to the largest jitter case of 30 km s−1. The results from the calculations using 5° spot are shown with plus signs, and the ones from 17° spot with squares. It is clear that the jitter increases with increasing v sin i, as has been shown by earlier studies (e.g. Saar & Donahue 1997; Desort et al. 2007). In addition, when the results are normalized to the highest value the values calculated from the two different spot sizes show the same trends. This means that the magnitude of the jitter depends on the spot size, but the increase in v sin i affects the jitter measurement the same way regardless of the spot size.

Figure 9.

Full amplitude of the spot-caused jitter with different stellar rotation rates. Two different spot configurations were used for this test: 5° radius equatorial spot shown by plus signs, and 17° radius spot denoted by squares. The results are normalized to the highest jitter case of v sin i = 30 km s−1. In the plot the x-axis gives the v sin i in km s−1 and y-axis the normalized jitter. The solid curve is the fit to the jitter obtained from the 5° spot, and the dotted line is the v sin i dependence law deduced by Desort et al. (2007).

Desort et al. (2007) published a formula for jitter–v sin i correlation. They obtained from their tests that the jitter amplitude depends on the spot coverage of the visible surface and v sin i. We use their formula to compute jitter–v sin i dependence, and compare it to our results. The resulting linear trend is plotted in Fig. 9 with a dotted line. Our results agree well with those obtained by Desort et al. (2007). The amplitude of the jitter estimated by Desort et al. (2007) is very similar to the ones given by Saar & Donahue (1997) and Boisse et al. (2012). On the other hand, when comparing our jitter amplitudes to the ones from Desort et al. (2007), our values are larger. The jitter from the 5°spot is approximately 20 per cent lower at the low v sin i values in the results obtained using the formula by Desort et al. (2007). The situation improves towards the higher v sin i values, and for the v sin i = 30 km s−1 the difference in only few per cent. Some of the discrepancy could be explained by different wavelength regions and spectral resolutions that were used in these investigations.

3.1.3 Effect of inclination

The effect of inclination of the stellar rotation axis was studied using a fixed v sin i of 2 km s−1. When the inclination is changed, the v sin i changes too. If a star is viewed pole-on there would be no rotational broadening. As the jitter also depends on the broadening of the spectral lines, we decided to use a fixed v sin i value for this test. The full amplitude of the measured jitter at different inclination angles is shown in Fig. 10. The test reveals the expected behaviour of the jitter, where the amplitude decreases with the decreasing visibility of the equatorial spot (decreasing inclination angle). The equatorial spot has the maximum effect on the jitter when the star is viewed equator on. The visibility of the spot is reduced when viewed increasingly from the direction of the pole, and thus also the impact of the spot on the jitter amplitude decreases. When viewed from almost the pole (inclination of 1°), the spot is seen at the limb and at all the rotational phases. This test does not include the effect of decreasing line broadening with decreasing inclination, which would make the change in jitter amplitude even more pronounced.

Figure 10.

Full amplitude of the spot-caused jitter at different inclinations of the stellar rotation axis. The spot configuration is the same as in Fig. 1(a). In the plot the x-axis gives the inclination in degrees and y-axis the full jitter amplitude in m s−1. As is expected the level of jitter lowers when the visibility of the equatorial spot decreases with decreasing inclination angle (viewing progressively more pole-on).

3.2 Solar-like activity patterns

For investigating the typical jitter amplitude caused by solar activity we have created 50 random spot configurations. All the spot configurations have spot filling factors of 0.1 per cent. This value is over the whole stellar surface and represents normal solar activity level (see e.g. Balmaceda et al. 2009). The spots have been restricted to occur between latitudes −30° and +30°, which also is the typical latitude range for solar activity.

The jitter curves calculated from all the 50 different spot configurations are shown in Fig. 11. As can be seen the typical jitter with this spot configuration varies between −2 and +2 m s−1. The mean amplitude of the full jitter is 4.5 m s−1 with the standard deviation being 2.8 m s−1. The minimum full amplitude from these spot configurations is 1.5 m s−1 and the maximum 12.3 m s−1.

Figure 11.

Jitter curves calculated from 50 different randomly created solar-like spot configurations with the spot filling factors of 0.1 per cent and spot latitudes restricted to −30° and +30°. The x-axis is the longitude in degrees and y-axis the jitter in m s−1. Note that the curves are created in such a way that the jitter at first observations (first longitude) is always zero.

The spot configurations and resulting jitter are plotted in Figs 12 and 13 for the minimum and maximum jitter case, respectively. The smallest jitter results from a case where several tiny spots are distributed relatively evenly over the longitude. On the other hand, the largest jitter arises from a spot configuration where there is one large spot together with one tiny one (at longitude 0, just below the equator). This is what one would expect and means that the exact jitter from a solar-like spot configuration depends largely on the exact spot distribution and how concentrated into active regions the spots are. This is in line with the studies which show that the radial velocity jitter can be estimated based on the photometric variability (e.g. Aigrain et al. 2012). Concentrated spots introduce more photometric variability, and also more radial velocity jitter.

Figure 12.

Spot configuration and jitter curve from the spot configuration that out of the 50 randomly created configurations results in the smallest jitter. The jitter has a full amplitude of 1.5 m s−1. The upper plot gives the spot configuration. The x-axis is the longitude and y-axis the latitude, both are given in degrees. The lower plot shows the corresponding jitter curve. Here the x-axis is the longitude in degrees and y-axis jitter in m s−1.

Figure 13.

The same as Fig. 12, but now for the largest jitter case. The resulting jitter has full amplitude of 12.3 m s−1.

We have to note though, that due to the limited resolution of the grid on the stellar surface, the filling factor is not always exactly the same. The code for creating the spotted surface will add a spot and then check the filling factor. Another spot is added if needed. This is done until the filling factor is greater than, or equal to, the one that was specified. Therefore the generated filling factors can be slightly above the input value, and are not always exactly identical.

3.3 Active solar-like stars

For investigating what the jitter behaviour of young, very active solar-like stars would be, we use the temperature maps of V889 Her by Järvinen et al. (2008). Their observations of the photospheric and chromospheric properties of V899 Her indicate that the quiet photosphere of V889 Her is similar to the one of the present day Sun, while the chromosphere shows much stronger activity. Their temperature maps, obtained using Doppler imaging, show that the polar regions are covered by spots, which are about 1500 K cooler than the quiet photosphere. Some evidence for cyclic magnetic activity is also seen both from photometry and Doppler imaging results.

Here we calculate the jitter resulting from the temperature maps of V889 Her for four different years: 1999, 2001, 2005 and 2007 (Järvinen et al. 2008). In the jitter calculations the wavelength range 5952–5998 Å is used and the inclination and v sin i are set to the ones determined from Doppler imaging, i = 60° and 37.5 km s−1, respectively. The spectral resolution used in the calculation was the same as the one in the original observations (77 000 and 2.5 pixels over the resolution element), and the grid size of the visible stellar surface is set to that of the original temperature map, 30 × 60.

The original V889 Her temperature maps and the calculated jitter are shown in Fig. 14. The full amplitude of jitter varies between 435 m s−1 calculated from the 2005 map and 774 m s−1 obtained from the 2001 map. These values are similar to radial velocity variations caused by a hot Jupiter around solar-mass star. Moulds et al. (2013) show that in this kind of cases some of the activity signal can be cleaned from the spectral line profiles and Jupiter-mass planets on close orbits can be recovered. In a recent paper Jeffers et al. (2014) study the detectability of planets around young active solar-like stars. They conclude that Jupiter-mass planets can be detected on close-in orbits around fast-rotating young active stars, Neptune-mass planets around moderate rotators and super-Earths only around very slowly rotating stars. The calculations carried out here based on the V889 Her spot configurations support these conclusions.

Figure 14.

Temperature maps of V889 Her (Järvinen et al. 2008) and the resulting jitter. The maps shown in the upper part of the plot are from years 1999, 2001, 2005 and 2007 (from left to right). The temperature ranges in the maps are 4624–6327 K (1999), 5251–6037 K (2001), 4498–6185 K (2005) and 4561–6268 K (2007). The x-axis gives the longitude in degrees and y-axis the latitude in degrees. The lower panels in the plots show the calculated jitter from each temperature map at 20 different phases evenly distributed over the stellar rotational phase. The x-axes give the longitude in degrees and y-axes the jitter in m s−1.

4 DISCUSSION

4.1 Starspots as a confusion factor

Not only can starspots cause noise in the radial velocity measurements, they can also mimic planetary signals and change the shape of the radial velocity curve.

For example, a large spot group close to the pole, which is viewed at low inclination, is visible on the surface all the time. The spot group causes larger jitter when it has a higher visibility, i.e. in the front, than when it is behind the pole close to the limb of the star. This behaviour can be confused with the radial velocity variation caused by an orbiting planet and Keplerian fits to the curve can be easily obtained. Several different cases of this kind of confusion have been already discussed by Desort et al. (2007).

Another confusion occurs when the location of the spot on the surface is such that together with the planetary signal it actually causes subtle changes in the shape of the measured radial velocity curve. For example spots at ‘correct’ location on the surface can change the shape of the radial velocity curve of a circular orbit into something that could be interpreted as a more eccentric orbit. An example of this kind of changes is given in Fig. 15. There the 5°spot from our tests (see Fig. 1a) is used together with a 0.8 Jupiter-mass planet on a circular orbit. The calculations are carried out at wavelength 5952–5998 Å and using spectral resolution of 100 000. The spot itself introduces jitter (see Fig. 15b) which is only about 10 per cent of the radial velocity variation of the planet (Fig. 15c). Still, the combination of these two changes the input radial velocity curve in subtle ways. In Fig. 15(d) the total radial velocity curve from the planet and spot jitter is plotted with squares. If that is compared to the original planet radial velocity curve, which is overplotted with plus signs, one can see that the phase of the maximum radial velocity has changed. The slope of the increase is now steeper and the declining slope more shallow. These changes would be interpreted as slightly eccentric orbit instead of a circular one. Naturally the activity could also work the other way, and make eccentric orbits appear more circular.

Figure 15.

An example of the effect of the spot signal on the radial velocity curve of the planet. The axes in the subplots are the same as in Fig. 4 (a) Stellar surface configuration with one equatorial 5°spot of 4000 K and unspotted surface temperature of 5800 K. (b) The resulting jitter curve from the spot configuration. The jitter values are calculated at 20 evenly spaced rotational phases from the wavelength region 5952–5998 Å and spectral resolution of 100 000. (c) Radial velocity curve of a 0.8 Jupiter-mass planet on a 25 d circular orbit around 1 M star. (d) The radial velocity curve from the spotted surface and the planet is given with squares. The plus signs give the original ‘planet only’ radial velocity curve (also shown in subplot c).

4.2 Solar-like activity cycle

One can question how typical the solar activity level is among solar-like stars in general. There have been early tentative suggestions that the Sun might be photometrically more quiet than similar stars (e.g. Radick et al. 1998). The high-precision photometric data from Kepler and Convection, Rotation and planetary Transits (CoRoT) can help to answer this question. The early results from Kepler implied that the Sun could indeed be more quiet than an average solar-type star (Gilliland et al. 2011; McQuillan, Aigrain & Roberts 2012). A new study by Basri, Walkowicz & Reiners (2013) has revisited the activity fraction of solar-like stars in the Kepler data. Their results show that 25–30 per cent of solar-type stars are more active than the Sun. The exact fraction depends on the time-scales used in the study, what is meant by ‘more active than the Sun’, and on the magnitude limit of the sample. In light of these investigations it seems appropriate to use the solar cycle as a proxy for cyclic activity in other stars.

4.2.1 Jitter during activity cycle

The number of spots on the solar surface varies with the 11-yr cycle. The spot coverage of the solar surface varies from 0 to about 0.5 per cent. Typically the minima show very small spot coverage factors of 0.01–0.02 per cent, or even zero for extended periods. The maxima on the other hand have quite varying spot coverage fractions. A small maximum has a typical spot coverage of 0.2 per cent, whereas a strong maximum has a spot coverage fraction of 0.5 per cent (see e.g. Balmaceda et al. 2009). During the cycle the latitudes of the spots also change. The new cycle starts with a small amount of spots which appear at relatively high latitudes, around ±30°. During the cycle the activity migrates slowly towards the equator, and during the next minimum the last spots of the old cycle appear at latitudes ±10°. This behaviour also affects the activity-caused jitter. For more details on the latitudinal migration of sunspots within the solar cycle see e.g. Carrington (1858), Maunder (1903) and Hathaway (2011).

We have created random spot configurations with typical sunspot coverage fractions and latitude ranges for studying the jitter over the solar-like activity cycle. In total six different activity cases were investigated: two activity minima, one average activity case and three activity maxima. The details of the cases that were studied are given in Table 2. For all the cases 100 different spot configurations were created. The jitter was calculated at wavelengths 5952–5998 Å using a resolution of 100 000, inclination of 90° and v sin i of 2.1 km s−1.

Table 2.

The different activity cases created to study the behaviour of the spot caused jitter over solar-like activity cycle.

CaseSpot frac.LatitudeJitterσ
(per cent)(°)(m s−1)(m s−1)
Minimum 10.02+20 to +301.81.6
Minimum 20.02−10 to +102.01.3
Average0.1−20 to +204.22.4
Maximum, small0.2−20 to +205.32.9
Maximum, medium0.3−20 to +206.62.4
Maximum, large0.5−20 to +207.72.9
CaseSpot frac.LatitudeJitterσ
(per cent)(°)(m s−1)(m s−1)
Minimum 10.02+20 to +301.81.6
Minimum 20.02−10 to +102.01.3
Average0.1−20 to +204.22.4
Maximum, small0.2−20 to +205.32.9
Maximum, medium0.3−20 to +206.62.4
Maximum, large0.5−20 to +207.72.9
Table 2.

The different activity cases created to study the behaviour of the spot caused jitter over solar-like activity cycle.

CaseSpot frac.LatitudeJitterσ
(per cent)(°)(m s−1)(m s−1)
Minimum 10.02+20 to +301.81.6
Minimum 20.02−10 to +102.01.3
Average0.1−20 to +204.22.4
Maximum, small0.2−20 to +205.32.9
Maximum, medium0.3−20 to +206.62.4
Maximum, large0.5−20 to +207.72.9
CaseSpot frac.LatitudeJitterσ
(per cent)(°)(m s−1)(m s−1)
Minimum 10.02+20 to +301.81.6
Minimum 20.02−10 to +102.01.3
Average0.1−20 to +204.22.4
Maximum, small0.2−20 to +205.32.9
Maximum, medium0.3−20 to +206.62.4
Maximum, large0.5−20 to +207.72.9

The results from the jitter calculations are shown in Fig. 16, and the mean jitter and its standard deviation are also given in Table 2. In the plot, the x-axis gives the time in years, and the mean jitter of the different activity cases have been plotted at a time it would typically be observed during the solar 11-yr cycle. The standard deviation of the jitter from 100 individual spot configurations has been plotted as an error bar for each activity case. The result from the average jitter case has been plotted both in the rising and declining phase of the activity.

Figure 16.

Results for the investigation of jitter during a solar-like activity cycle. The x-axis gives the time in years and y-axis the jitter in m s−1. The mean full jitter amplitude calculated from 10 individual spot configurations with the different activity cases given in Table 2 are plotted approximately at the time they would occur during the solar 11-yr cycle. The error bar gives the standard deviation of the measurements from the 10 individual spot configurations.

For the maximum solar activity three different jitter cases were calculated. They are marked in Fig. 16 by letter S (small maximum, spot fraction 0.2 per cent), M (medium maximum, spot fraction 0.3 per cent) and L (large maximum, spot fraction 0.5 per cent). As expected, the large-sized maximum results in the largest jitter values, and the small-sized maximum induces the smallest jitter. Still, all the average jitter values for different solar maxima cases are the same within the standard deviation, and the variation in the jitter, as is seen from the standard deviation, is similar for all the cases.

Both of the minima cases have the smallest jitter values, and also the smallest standard deviation. This is because with such a small spot coverage (0.02 per cent) basically only one small spot is present on the surface. The exact location of the spot changes slightly, and the differences in latitude result in small changes in the full jitter amplitude. It is interesting to note that the mean jitter from the early cycle case, where the spots are around latitudes ±30° results in smaller jitter than the late cycle case with spots at latitudes ±10°. The effect of the jitter is largest when the spots are best visible, i.e. in the case of inclination 90° around equator. One should also note that the spots for the early-cycle case have been created with spot latitude 20°–30°, not taking into account that on the Sun the spots would appear both around latitude −30° and +30°. Regardless, with such small spot coverage fractions only one spot is created, and therefore this has no practical effect in the full jitter amplitude, which is what is studied here.

These calculations do not take into account the lifetime of the sunspots. On the other hand, the larger sunspot on average live longer, couple of weeks, instead of couple of days (see e.g. Gnevyshev 1938; Waldmeier 1955; Petrovay & van Driel-Gesztelyi 1997). Here we by necessity concentrate on larger spots. Therefore, not taking into account the lifetime of the spots should not significantly affect our investigation. In addition, several authors have reported that the solar activity tends to occur at active longitudes (e.g. Berdyugina & Usoskin 2003; Usoskin, Berdyugina & Poutanen 2005; Juckett 2006). This effect has not been taken into in our calculations; we have restricted the spot occurrence only in the latitude direction. The solar active longitudes would group the spots in longitudinal direction and increase the jitter. Thus, the random distribution investigated here can be considered as the lower limit for the jitter.

4.3 Detecting planets on a 1-yr orbit

An interesting question is how the jitter during a solar-like cycle would affect the detection of a planet in the habitable zone of its host star. We have studied this issue by introducing a planet on 1-yr orbit around a star showing a cycle similar to the one described in the previous section.

The cycle is thought to last 11 yr, like in the Sun, and the observing period is 5 yr, covering the cycle from minimum to maximum. Each year the target can be observed during its visibility, here called the observing season, and each season a fixed number of observing runs is carried out. The length of the observing run is always 25 d, which is also the rotation period of the target star. During one run 20 evenly spaced observations are carried out. For each observing run during the first year a spot configuration from the solar minimum case is randomly chosen. For the second and third year the spot configurations have been chosen randomly from the average activity case. During the fourth year spot configurations are from the small activity maximum case, and during the fifth year from the medium maximum maps. After this a planet on a 1-yr orbit is introduced to the jitter measurements. Fig. 17 shows an example of radial velocity measurements created this way. In the example the planet has the same mass as Neptune (17 M) and the observing season lasts the whole year and has five individual observing runs during it. It can easily be seen that the activity caused noise in the radial velocity measurements increases with the advancing activity cycle.

Figure 17.

Simulated radial velocity measurements of Neptune-mass planet orbiting a solar-like star with a solar-like activity cycle. The observations span 5 yr, and during each year there are five separate observing runs lasting 25 nights. The x-axis gives time in days and y-axis the radial velocity in m s−1.

For investigating the detectability of the planet we use Lomb–Scargle period search method for unevenly sampled data (Scargle 1982; Horne & Baliunas 1986; Press & Rybicki 1989). The periodogram resulting from the simulated data presented in Fig. 17 is shown in Fig. 18. The dominant frequencies in this case are the orbital period of the planet (marked by a solid vertical line) and its harmonics (marked by dashed vertical lines). The dotted vertical lines denote the rotation period (25 d) and half a rotation period of the star, whereas the horizontal dotted line is the analytical 3σ detection threshold following the false alarm probability (FAP) formulation of Scargle (1982) and taking into account the modifications of Horne & Baliunas (1986). The horizontal dashed line, on the other hand, is the FAP obtained from white noise simulations with 10 000 iterations. The two FAP values are similar, which is to be expected because our data are not severely unevenly distributed, the case where the analytical method would strongly underestimate the FAP. Our data are evenly distributed over the observing run, and the observing runs are evenly distributed throughout the year. In this simulation the Neptune-mass planet is easily detected even with the noise from the stellar activity. The case remains the same if we shorten the observing season to 150 d and only have three observing runs during it. This is a more realistic case because of the limited visibility of the targets and also owing to the telescope time allocation process. The Lomb–Scargle periodogram of this case is shown in Fig. 19. The orbital period of the planet is still easily detectable, but the harmonics of the orbital period have become more pronounced. Part of the effect is also due to aliasing caused by larger data gaps.

Figure 18.

Lomb–Scargle periodogram obtained from the simulated radial velocity measurements presented in Fig. 17. The dotted horizontal line is the analytical 3σ detection threshold and the dashed horizontal line the numerical one. The solid vertical line gives the original period of the planet (365 d) and the dashed vertical lines give its harmonics. The dotted vertical lines are the rotation period of the star (25 d) and half of the rotation period of the star. The x-axis gives the period in days and y-axis the power spectral density.

Figure 19.

The same as Fig. 18, but now with observing season only lasting 150 d and with only three observing runs during this time.

On the other hand, the possibility of detecting Earth-mass planet around a solar-like star is much more challenging. The periodogram from the optimum case where the observing season lasts the whole year and there are five individual observing runs during the season are shown in Fig. 20. No indication of the true orbital period is seen, and the strongest periodicity is half the stellar rotation period. The situation is somewhat improved when the number of observing runs is increased. In Fig. 21 a case where there are 50 observing runs during each year is shown (this implies observing more than once a night). In this case the signature of the true orbital period starts to emerge, but still it cannot be considered significant. Our tests show that the 3σ detection limit for a planet in a habitable zone of a solar-like star with a solar-like activity cycle is around 6 M (when using five observing runs distributed over a full year, and total length of observations of 5 yr).

Figure 20.

The same as Fig. 18, but now with an Earth-mass planet.

Figure 21.

The same as Fig. 18, but now with an Earth-mass planet and 50 observing runs during the observing season of 1 yr.

Lagrange et al. (2010) investigated the detectability of the Earth in the habitable zone around a solar-like star. They concluded that with the highest precision instruments one can only detect an Earth-mass planet after several years of intensive monitoring, and then preferably during the low activity phase of the star. In another investigation Dumusque et al. (2011b) investigate the detection limits with a High Accuracy Radial velocity Planet Searcher (HARPS)-like instrument, taking into account oscillation, granulation and activity effects. They study different observing strategies to minimize the stellar noise and conclude that applying three measurements per night of 10 min every 3 d, 10 nights a month is the best strategy out of the ones they tested. Depending on the activity level this strategy would allow to detect 2.5–3.5 Earth-mass planets in the habitable zone of an early-K dwarf. Still, this means a planet larger than the Earth and in a habitable zone that is closer to the star than 1 au.

4.4 Future investigations

The codes developed here are very versatile and large variety of spot configurations and exoplanets on different orbits can easily be studied. This opens a possibility to do statistical studies of the effect of stellar activity on the detection of exoplanets. One can also use the same methods as are used to detect exoplanets in recovering the input planet from our spectra. The input parameters of the planet are know, and therefore this will also allow for determining at which accuracy different planetary parameters can be recovered. In the second part of this series of papers (Andersen & Korhonen 2014), we apply the methods developed here to M dwarfs. We will also study in detail the effect stellar activity has on the recovered planetary parameters.

5 CONCLUSIONS

We have developed methods for investigating the radial velocity jitter caused by starspots. The method allows creating many spot configurations with the same spot filling factors and also enables selecting active latitude and longitude ranges. Planetary signatures can easily be added to the spectra and analysed. From the tests and implementation to solar-like stars we can draw the following conclusions.

  • As has been seen in the previous studies observations at longer wavelengths decrease the measured radial velocity jitter. The tests carried out in this study show that the decrease of approximately 50 per cent in the full jitter amplitude is achieved at wavelengths around 7000 Å in comparison to 3700 Å. Between wavelengths 7000 and 9000 Å no significant further decrease in the jitter amplitude is observed.

  • The spectral resolution does not affect the jitter amplitude significantly at the generally used resolving powers of 50 000–130 000 and higher. On the other hand, resolution of 20 000 and less decreases the jitter, but also decreases the measurement accuracy.

  • We verify the previous results showing that the full jitter amplitude depends on the stellar rotational velocity, v sin i. The dependence is linear and even though the exact jitter amplitude depends on the spot size, the slope of the correlation does not.

  • The solar-like activity patterns create largely varying amounts of radial velocity jitter. From a spot coverage factor that represents average solar activity, the full jitter amplitude recovered from our simulated data varies approximately between 1 and 12 m s−1. The exact value is driven by how concentrated the spots are.

  • The mean full jitter amplitude varies during the solar-like activity cycle between approximately 1 and 9 m s−1.

  • With realistic observing frequency and solar-like cyclic activity a Neptune-sized planet on a 1 yr orbit around a solar-mass star can be recovered with high significance. The recovery of an Earth-mass planet on a similar orbit on the other hand is very challenging.

  • The starspots do not only create noise in the radial velocity curves, they can also affect the shape of the radial velocity curve in such a way that the determined orbital parameters change. Eccentricity especially can be affected.

  • The spot surface creation and planet orbit implementation software developed in this study allow for statistical studies of effect of spot jitter in exoplanet detection. These issues are addressed in the second paper of this series.

The authors would like to thank the anonymous referee whose comments helped to improve this paper. HK acknowledges the support from the European Commission under the Marie Curie IEF Programme in FP7. JMA acknowledges support through an NSF Graduate Research Fellowship and the Nordic Research Opportunity award.

REFERENCES

Aigrain
S.
Pont
F.
Zucker
S.
MNRAS
2012
, vol. 
419
 pg. 
3147
 
Andersen
J. M.
Korhonen
H.
MNRAS
2014
 
preprint (arXiv:1501.01302)
Balmaceda
L. A.
Solanki
S. K.
Krivova
N. A.
Foster
S.
J. Geophys. Res.
2009
, vol. 
114
 pg. 
7104
 
Baraffe
I.
Chabrier
G.
Allard
F.
Hauschildt
P. H.
A&A
1998
, vol. 
337
 pg. 
403
 
Barnes
J. R.
Jeffers
S. V.
Jones
H. R. A.
MNRAS
2011
, vol. 
412
 pg. 
1599
 
Basri
G.
Walkowicz
L. M.
Reiners
A.
ApJ
2013
, vol. 
769
 pg. 
37
 
Baumann
I.
Solanki
S. K.
A&A
2005
, vol. 
443
 pg. 
1061
 
Bennett
D. P.
, et al. 
ApJ
2012
, vol. 
757
 pg. 
119
 
Berdyugina
S. V.
Usoskin
I. G.
A&A
2003
, vol. 
405
 pg. 
1121
 
Bogdan
T. J.
Gilman
P. A.
Lerche
I.
Howard
R.
ApJ
1988
, vol. 
327
 pg. 
451
 
Boisse
I.
Bouchy
F.
Hébrard
G.
Bonfils
X.
Santos
N.
Vauclair
S.
A&A
2011
, vol. 
528
 pg. 
A4
 
Boisse
I.
Bonfils
X.
Santos
N. C.
A&A
2012
, vol. 
545
 pg. 
A109
 
Carrington
R. C.
MNRAS
1858
, vol. 
19
 pg. 
1
 
Cegla
H. M.
, et al. 
MNRAS
2012
, vol. 
421
 pg. 
L54
 
Cegla
H. M.
Stassun
K. G.
Watson
C. A.
Bastien
F. A.
Pepper
J.
ApJ
2014
, vol. 
780
 pg. 
104
 
Desort
M.
Lagrange
A.-M.
Galland
F.
Udry
S.
Mayor
M.
A&A
2007
, vol. 
473
 pg. 
983
 
Dumusque
X.
Udry
S.
Lovis
C.
Santos
N. C.
Monteiro
M. J. P. F. G.
A&A
2011a
, vol. 
525
 pg. 
A140
 
Dumusque
X.
Santos
N. C.
Udry
S.
Lovis
C.
Bonfils
X.
A&A
2011b
, vol. 
527
 pg. 
A82
 
Dumusque
X.
Boisse
I.
Santos
N. C.
ApJ
2014
, vol. 
796
 pg. 
132
 
Gilliland
R. L.
, et al. 
ApJS
2011
, vol. 
197
 pg. 
6
 
Gnevyshev
M. N.
Pulkovo Obser. Circ.
1938
, vol. 
24
 pg. 
37
 
Gray
R. O.
Corbally
C. J.
AJ
1994
, vol. 
107
 pg. 
742
 
Gustafsson
B.
Edvardsson
B.
Eriksson
K.
Jørgensen
U. G.
Nordlund
Å.
Plez
B.
A&A
2008
, vol. 
486
 pg. 
951
 
Hackman
T.
Jetsu
L.
Tuominen
I.
A&A
2001
, vol. 
374
 pg. 
171
 
Hathaway
D. H.
Sol. Phys.
2011
, vol. 
273
 pg. 
221
 
Horne
J. H.
Baliunas
S. L.
ApJ
1986
, vol. 
302
 pg. 
757
 
Huerta
M.
Johns-Krull
C. M.
Prato
L.
Hartigan
P.
Jaffe
D. T.
ApJ
2008
, vol. 
678
 pg. 
472
 
Järvinen
S. P.
Korhonen
H.
Berdyugina
S. V.
Ilyin
I.
Strassmeier
K. G.
Weber
M.
Savanov
I.
Tuominen
I.
A&A
2008
, vol. 
488
 pg. 
1047
 
Jeffers
S. V.
Barnes
J. R.
Jones
H. R. A.
Reiners
A.
Pinfield
D. J.
Marsden
S. C.
MNRAS
2014
, vol. 
438
 pg. 
2717
 
Jenkins
J. S.
, et al. 
MNRAS
2006
, vol. 
372
 pg. 
163
 
Juckett
D. A.
Sol. Phys.
2006
, vol. 
237
 pg. 
351
 
Kupka
F. G.
Piskunov
N. E.
Ryabchikova
T. A.
Stempels
H. C.
Weiss
W. W.
A&AS
1999
, vol. 
138
 pg. 
119
 
Kurucz
R. L.
Kurucz CD, No. 13
1993
Kurucz
R. L.
Avrett
E. H.
1981
 
SAO Special Report, No. 391
Lagrange
A.-M.
Desort
M.
Meunier
N.
A&A
2010
, vol. 
512
 pg. 
A38
 
Lanza
A. F.
Boisse
I.
Bouchy
F.
Bonomo
A. S.
Moutou
C.
A&A
2011
, vol. 
533
 pg. 
A44
 
Ma
B.
Ge
J.
ApJ
2012
, vol. 
750
 pg. 
172
 
McQuillan
A.
Aigrain
S.
Roberts
S.
A&A
2012
, vol. 
539
 pg. 
A137
 
Marchwinski
R. C.
Mahadevan
S.
Robertson
P.
Ramsey
L.
Harder
J.
ApJ
2015
, vol. 
798
 pg. 
63
 
Martínez-Arnáiz
R.
Maldonado
J.
Montes
D.
Eiroa
C.
Montesinos
B.
A&A
2010
, vol. 
520
 pg. 
A79
 
Maunder
E. W.
Observatory
1903
, vol. 
26
 pg. 
329
 
Meunier
N.
Lagrange
A.-M.
A&A
2013
, vol. 
551
 pg. 
A101
 
Meunier
N.
Desort
M.
Lagrange
A.-M.
A&A
2010
, vol. 
512
 pg. 
A39
 
Mikkola
S.
Celest. Mech.
1987
, vol. 
40
 pg. 
329
 
Moulds
V. E.
Watson
C. A.
Bonfils
X.
Littlefair
S. P.
Simpson
E. K.
MNRAS
2013
, vol. 
430
 pg. 
1709
 
Oshagh
M.
Santos
N. C.
Boisse
I.
Boué
G.
Montalto
M.
Dumusque
X.
Haghighipour
N.
A&A
2013
, vol. 
556
 pg. 
A19
 
Petrovay
K.
van Driel-Gesztelyi
L.
Sol. Phys.
1997
, vol. 
176
 pg. 
249
 
Piskunov
N. E.
Tuominen
I.
Vilhu
O.
A&A
1990
, vol. 
230
 pg. 
363
 
Press
W. H.
Rybicki
G. B.
ApJ
1989
, vol. 
338
 pg. 
277
 
Queloz
D.
, et al. 
A&A
2001
, vol. 
379
 pg. 
279
 
Radick
R. R.
Lockwood
G. W.
Skiff
B. A.
Baliunas
S. L.
ApJS
1998
, vol. 
118
 pg. 
239
 
Reiners
A.
Bean
J. L.
Huber
K. F.
Dreizler
S.
Seifahrt
A.
Czesla
S.
ApJ
2010
, vol. 
710
 pg. 
432
 
Saar
S. H.
Donahue
R. A.
ApJ
1997
, vol. 
485
 pg. 
319
 
Sbordone
L.
Bonifacio
P.
Castelli
F.
Kurucz
R. L.
Mem. Soc. Astron. Ital. Suppl.
2004
, vol. 
5
 pg. 
93
 
Scargle
J. D.
ApJ
1982
, vol. 
263
 pg. 
835
 
Solanki
S. K.
A&AR
2003
, vol. 
11
 pg. 
153
 
Usoskin
I. G.
Berdyugina
S. V.
Poutanen
J.
A&A
2005
, vol. 
441
 pg. 
347
 
Waldmeier
M.
Ergebnisse und Probleme der Sonnenforschung, 2
1955
Leipzig
Aufl. Akad. Verlagsges.