© ESO, 2017

1. Introduction

During planetary transits, the observed flux of the star decreases because part of the stellar disk is occulted by the planetary disk. If the accuracy of the transit photometry is high enough, it is possible to analyze the fine details of the transit light curve (LC) and characterize the stellar surface. In particular, if the planet transits over stellar spots, then an apparent rebrightening (called “transit anomaly”) of the star is observed during the spot crossing, as it corresponds to an increase of the flux received from the star because a darker area is occulted.

The analysis of transit anomalies has been theorized by Schneider (2000). From the observational point of view, the study of Silva (2003) is the first example of how these transit anomalies can be analyzed to obtain information on the morphology and physical parameters of the crossed spots. Since then, many other authors have analyzed the high-accuracy transit photometry of active stars with the aim of modeling the stellar spots (e.g., Pont et al. 2007; Rabus et al. 2009; Huber et al. 2010; Winn et al. 2010; Tregloan-Reed et al. 2013; Béky et al. 2014b; Nascimbeni et al. 2015; Mancini et al. 2015).

While the analysis of individual transits allows us to model the spots intersected by the planet in those specific epochs, the observations of many transits give us the possibility to study a large number of spots, which can be statistically significant in order to characterize the average photospheric activity (Béky et al. 2014b). Furthermore, if the same spots are repeatedly occulted, then the spin-orbit alignment of the system can be constrained (Sanchis-Ojeda & Winn 2011a; Nutzman et al. 2011; Désert et al. 2011).

In the case of recurrent anomalies, we can also characterize the time evolution of the spots and the rotation of the star (Silva-Valio 2008; Silva-Valio & Lanza 2011; Sanchis-Ojeda et al. 2011b; Béky et al. 2014a). This task is usually accomplished by continuous long-term photometric and spectroscopic monitorings of the stars, requiring a large amount of telescope time. Moreover, the results obtained via these methods are affected by several limitations and degeneracies between the parameters of the adopted models (Lanza 2016). Another advantage of transiting systems is, then, that the evolution of individual spots can be studied using a few consecutive transits.

In recent years, the analysis of transit anomalies has been approached using various methods. For example, some authors fit the transit anomaly with a given analytical model to retrieve the size of the spots and their flux ratio with respect to the photosphere (e.g., Sanchis-Ojeda & Winn 2011a; Nascimbeni et al. 2015). Some other authors have written more sophisticated and computationally expensive codes, which fit the transit light curves by means of Monte Carlo algorithms (e.g., Tregloan-Reed et al. 2013; Béky et al. 2014b). In both cases, some assumptions on the shape and distribution of the spots over the stellar disk are needed.

In this paper we present the Tomography Of Spotted transit Chords (TOSC) code, a new and fast algorithm which analyzes the transit anomalies and reconstructs the flux distribution inside the area of the stellar surface intersected by the transiting planet, namely the “transit chord”. The term “chord” here is used improperly, as the chord is by definition a straight segment with endpoints on a circle. From a pure geometrical point of view, it would be more appropriate to use the term “belt”. Nonetheless, in exoplanet science the two terms are used interchangeably. In this paper, we will prefer the notation “chord” as most widely used by the community. We will use the term belt only when it is required by the bi-dimensionality of the tomography.

The main advantage of our code with respect to previous approaches is that it only needs some assumptions on the geometry of the planetary system and on the spectrum of the star. In other words, TOSC does not assume any a priori parameter of the photospheric spots occulted by the planet (e.g., shape, size, and temperature). In Sect. 2 we describe the geometrical model and the inversion algorithm which returns the tomography of the transit chord. In Sect. 3 we test TOSC against a few simulated transits, in order to check its robustness and its response to photometric accuracy and spot parameters. Finally, in Sect. 4 we compare our code with other models for transit anomalies, applying TOSC to some transit LCs which have already been analyzed in literature.

2. The model

The out-of-transit stellar flux F_oot depends on the effective temperature of the star and on the presence of active regions on the stellar surface. Active regions evolve with time and rotate with the star. F_oot is thus by principle a function of time. Nonetheless, both effects (evolution and rotation of active regions) generally have timescales much longer than the typical duration of a planetary transit (i.e., a few hours). It is thus a common practice to fit a low-order polynomial to the out-of-transit data, and to normalize the observed fluxes by this polynomial, in order to remove any long-term trend in the observed LC. After normalization, the dependence of F_oot on time can thus be dropped. In the following discussion, all the fluxes are normalized by the low-order polynomial fit of the out-of-transit LC.

During the transit, the flux F(t) received from the stellar disk is given by: $\mathit{F}\mathrm{\left(}\mathit{t}\mathrm{\right)}\mathrm{=}{\mathit{F}}_{\mathrm{oot}}\mathrm{-}{\mathit{F}}_{\mathrm{occ}}\mathrm{\left(}\mathit{t}\mathrm{\right)}\mathit{,}$(1)where F_occ(t) is the occulted stellar flux at any time t. The dependence on time is given by the fact that F_occ(t) depends on the limb darkening (LD, see e.g., Winn 2010) profile. Moreover, the light curve may have anomalies in case the planet occults active regions, either in the form of cool spots or photospheric bright faculae, showing bumps or dips correspondingly (e.g., Silva 2003; Ballerini et al. 2012).

F(t) is the observable that we want to convert into the spatially-resolved brightness profile of the transit chord. To this purpose, in the following sections we describe how we model the stellar disk (Sect. 2.1), the transit chord (Sect. 2.2) and the transit LC (Sect. 2.3).

2.1. The stellar disk

If we neglect the effects of surface brightness inhomogeneities and of their rotation over timescales of a few hours typical of the transits, we can assume that the flux emerging from the stellar disk is constant, and is given by the surface integral of the specific intensity I over the disk.

The stellar disk observed from Earth is not uniformly bright due to the LD effect. We assume a quadratic profile (Claret 2004), which is generally adopted in fitting transit LCs. This profile is usually parametrized as $\frac{\mathit{I}\mathrm{\left(}\mathit{\mu }\mathrm{\right)}}{{\mathit{I}}_{\mathrm{0}}}\mathrm{=}\mathrm{1}\mathrm{-}{\mathit{u}}_{\mathrm{1}}\mathrm{(}\mathrm{1}\mathrm{-}\mathit{\mu }\mathrm{)}\mathrm{-}{\mathit{u}}_{\mathrm{2}}\mathrm{(}\mathrm{1}\mathrm{-}\mathit{\mu }{\mathrm{)}}^{\mathrm{2}}\mathit{,}$(2)where I₀ is the specific intensity at the center of the disk, μ = cosγ, γ is the angle between the line of sight and the normal to the stellar surface, and u₁ and u₂ are the limb darkening coefficients, which are part of the output of the transit LC fitting or are derived from a model atmosphere.

To compute the broadband intensity I₀ at the center of the stellar disk, we first convolve the intensity given by the BT-Settl synthetic spectral library by Baraffe et al. (2015) with the response function of the instrumental setup used to collect the LC. Then, the intensity I₀ is computed interpolating the temperature grid of the model at the effective temperature T_eff of the star.

The flux F_oot emerging from the stellar disk is thus given by the integral of I(μ) over the stellar disk (see Eq. (3) in Lanza 2016): ${\mathit{F}}_{\mathrm{oot}}\mathrm{=}\mathrm{2}\mathit{\pi }{\mathit{R}}_{\mathrm{\ast }}^{\mathrm{2}}{\mathrm{\int }}_{\mathrm{0}}^{\mathrm{1}}\mathit{I}\mathrm{\left(}\mathit{\mu }\mathrm{\right)}\mathit{\mu }\mathrm{d}\mathit{\mu }\mathit{.}$(3)Substituting Eq. (2) into Eq. (3) and integrating we obtain: ${\mathit{F}}_{\mathrm{oot}}\mathrm{=}\mathit{\pi }{\mathit{R}}_{\mathrm{\ast }}^{\mathrm{2}}{\mathit{I}}_{\mathrm{0}}\frac{\mathrm{6}\mathrm{-}\mathrm{2}{\mathit{u}}_{\mathrm{1}}\mathrm{-}{\mathit{u}}_{\mathrm{2}}}{\mathrm{6}}\mathit{,}$(4)where R_∗ denotes the stellar radius.

2.2. The transit chord

Assuming that the shift of the planet on the stellar disk during the photometric exposure is negligible compared to the planetary radius, the stellar flux F_occ(t) occulted by the planet at the time t during the transit is given by the integral of the stellar flux over the planetary disk in the sky-projected plane. We can thus assume that each photometric point taken during the transit is representative of an area on the projected stellar disk centered on the planet and having the same planetary radius. Hence, as the planet moves over the stellar disk, it scans the transit chord. The transit LC is thus, from this point of view, the observable carrying information on the flux distribution along the transit chord.

In particular, the U-shaped transit LC is due to the LD, while anomalies (bumps and/or dips) with respect to the best-fit transit model are due to active regions on the stellar surface. In both cases, the inference on LD and/or stellar activity is possible because the planetary radius is smaller than the length of the transit chord (the planetary radius is of the order of 10% or less of the stellar radius), thus the movement along the transit chord provides a spatially resolved scan of the chord itself. Roughly speaking, the spatial resolution of the scan is of the same order of the shift of the planet along the chord between two adjacent photometric points if the noise can be neglected.

To reconstruct the flux distribution along the transit chord, we refer the stellar disk to a Cartesian reference frame centered at its center. We assume the transit chord to be in the horizontal direction, at the ordinate y = b·R_∗, where R_∗ is the stellar radius and b the impact parameter in units of the stellar radius.

For the sake of simplicity, in our model the disk of the planet is assumed to be entirely inside the disk of the star. In other words, we do not model the ingress and the egress. We thus define a grid inside the transit chord made up of adjacent equal rectangular cells. If the orbit inclination is not 90°, the planetary track across the stellar disk is actually curved. The most tilted configuration occurs when the projected planetary orbit is an ellipse of semi minor axis equal to the stellar radius. In this case, the maximum displacement along the vertical direction of the planet during the transit is ${\mathit{R}}_{\mathrm{\ast }}\sqrt{\mathrm{1}\mathrm{-}\mathrm{(}{\mathit{R}}_{\mathrm{\ast }}\mathit{/}\mathit{a}{\mathrm{)}}^{\mathrm{2}}}$, where a is the semi major axis. For a giant planet (R_p ~ 0.1) in a tight orbit (R_∗/a ~ 0.1), the maximum displacement turns out to be a few percents of the planetary radius. Given this small effect, and the difficulty in taking into account warped cells inside the curved transit chord, for the sake of simplicity we approximate the transit chord as a rectilinear segment. Each cell thus has height equal to the diameter of the planet and width Δx depending on the number N of cells (Fig. 1). The length of the grid equals the length of the top side of the chord in Fig. 1, such that the first and last grid cells are completely contained into the projected stellar disk. Using Pythagoras’s theorem, the total length l_chord of the reconstructed chord is $\begin{array}{ccc}{\mathit{l}}_{\mathrm{chord}}\mathrm{=}\mathrm{2}{\mathit{R}}_{\mathrm{\ast }}\sqrt{\mathrm{1}\mathrm{-}{\mathrm{(}{\mathit{R}}_{\mathrm{p}}\mathrm{+}{\mathit{b}}^{\mathrm{)}}}^{\mathrm{2}}}\mathit{,}& & \end{array}$and the width of each cell is thus given by $\mathrm{\Delta }\mathit{x}\mathrm{=}\frac{{\mathit{l}}_{\mathrm{chord}}}{\mathit{N}}\mathrm{=}\frac{\mathrm{2}{\mathit{R}}_{\mathrm{\ast }}\sqrt{\mathrm{1}\mathrm{-}{\mathrm{(}{\mathit{R}}_{\mathrm{p}}\mathrm{+}{\mathit{b}}^{\mathrm{)}}}^{\mathrm{2}}}}{\mathit{N}}\mathit{,}$(5)where R_p is the planetary radius in units of the stellar radius.

Fig. 1 Schematic representation of the geometry discussed in the text. The circle represents the transiting planet, while the grid rectangles show the cells used in the chord reconstruction.

For the sake of simplicity, we assume that the intensity I_i emerging from the ith cell in the chord is uniform and equal to the intensity at the center of the cell: ${\mathit{I}}_{\mathit{i}}\mathrm{=}\mathit{I}\left(\mathit{\mu }\mathrm{(}{\mathit{x}}_{\mathit{i}}\mathit{,b}\mathrm{·}{\mathit{R}}_{\mathrm{\ast }}\mathrm{)}\right)\mathit{,}$(6)(i.e., we neglect the dependence of I on μ inside the cell, see Eq. (2)). The corresponding emerging flux is thus: $\begin{array}{ccc}{\mathit{F}}_{\mathit{i}}& \mathrm{=}& ”\mathit{I}\mathrm{\left(}\mathit{x,y}\mathrm{\right)}\mathrm{d}\mathit{x}\mathrm{d}\mathit{y}\mathrm{\simeq }\mathit{I}\left(\mathit{\mu }\mathrm{(}{\mathit{x}}_{\mathit{i}}\mathit{,b}\mathrm{·}{\mathit{R}}_{\mathrm{\ast }}\mathrm{)}\right)”\mathrm{d}\mathit{x}\mathrm{d}\mathit{y}\\ & \mathrm{=}& \mathit{I}\left(\mathit{\mu }\mathrm{(}{\mathit{x}}_{\mathit{i}}\mathit{,b}\mathrm{·}{\mathit{R}}_{\mathrm{\ast }}\mathrm{)}\right)\mathrm{·}\mathrm{2}{\mathit{R}}_{\mathrm{p}}{\mathit{R}}_{\mathrm{\ast }}\mathrm{\Delta }\mathit{x}\mathit{.}\end{array}$(7)where the double integral is extended over the surface of the ith cell.

2.3. The reconstruction algorithm

The flux dimming F_occ(t) at any time t during the transit is the sum of the fluxes F_i of the cells multiplied by the corresponding overlap fractional area between the cells and the planetary disk. With these assumptions, we obtain: ${\mathit{F}}_{\mathrm{occ}}\mathrm{\left(}\mathit{t}\mathrm{\right)}\mathrm{=}\sum _{\mathit{i}\mathrm{=}\mathrm{1}}^{\mathit{N}}{\mathit{w}}_{\mathit{i}}\mathrm{\left(}\mathit{t}\mathrm{\right)}\mathrm{·}{\mathit{F}}_{\mathit{i}}\mathit{,}$(8)where the weights w_i(t) are the overlap fractional areas, and thus run from 0 to 1.

The weights w_i(t) are computed by geometric considerations (Fig. 2). From left to right, the first non-zero weight is computed as the area (in cyan in Fig. 2) of the circular segment bounded by the arc $\begin{array}{c}{\mathrm{􏽣}}_{}\\ \mathit{AB}\end{array}$. The second weight (i.e., the orange area in Fig. 2) is then computed as the area of the circular segment bounded by the arc $\begin{array}{c}{\mathrm{􏽤}}_{}\\ \mathit{CD}\end{array}$ (the cyan and orange circular segment) minus the area of the circular segment bounded by $\begin{array}{c}{\mathrm{􏽣}}_{}\\ \mathit{AB}\end{array}$, and so forth.

Fig. 2 Graphical representation of the computation of the weights discussed in the text. The circle is the planetary disk centered in O, while the gray grid shows the cells drawn on the stellar disk partially covered by the planetary disk.

The N unknowns F_i can be retrieved solving the linear system: $\left[\begin{array}{c}\\ {\mathit{w}}_{\mathrm{1}}\mathrm{\left(}{\mathit{t}}_{\mathrm{1}}\mathrm{\right)}& \mathrm{·}\mathrm{·}\mathrm{·}\hspace{0.17em}& {\mathit{w}}_{\mathit{N}}\mathrm{\left(}{\mathit{t}}_{\mathrm{1}}\mathrm{\right)}\\ \begin{array}{c}\mathit{.}\\ \mathit{.}\\ \mathit{.}\end{array}& {\begin{array}{c}\mathit{.}\end{array}}^{\mathit{.}}\mathit{.}& \begin{array}{c}\mathit{.}\\ \mathit{.}\\ \mathit{.}\end{array}\\ {\mathit{w}}_{\mathrm{1}}\mathrm{\left(}{\mathit{t}}_{\mathit{M}}\mathrm{\right)}& \mathrm{·}\mathrm{·}\mathrm{·}\hspace{0.17em}& {\mathit{w}}_{\mathit{N}}\mathrm{\left(}{\mathit{t}}_{\mathit{M}}\mathrm{\right)}\end{array}\right]\mathrm{·}\left[\begin{array}{c}\\ {\mathit{F}}_{\mathrm{1}}\\ \begin{array}{c}\mathit{.}\\ \mathit{.}\\ \mathit{.}\end{array}\\ {\mathit{F}}_{\mathit{N}}\end{array}\right]\mathrm{=}\left[\begin{array}{c}\\ {\mathit{F}}_{\mathrm{occ}}\mathrm{\left(}{\mathit{t}}_{\mathrm{1}}\mathrm{\right)}\\ \begin{array}{c}\mathit{.}\\ \mathit{.}\\ \mathit{.}\end{array}\\ {\mathit{F}}_{\mathrm{occ}}\mathrm{\left(}{\mathit{t}}_{\mathit{M}}\mathrm{\right)}\end{array}\right]$(9)where t₁...t_M are the times when the M photometric points were collected.

In principle, there is no restriction on the maximum number of cells. Nonetheless, large values of N are discouraged for two reasons. First of all, if N>M then the system in Eq. (9) is overdetermined and the data are “overfit”. This corresponds to extracting from the data an amount of information larger than the information carried by the data themselves.

The second reason is given by the Nyquist sampling theorem. Assuming that the transit is evenly sampled at positions separated by the distance Δx₀ on the stellar disk, the planet samples the transit chord with the spatial frequency ${\mathit{f}}_{\mathrm{0}}\mathrm{=}\frac{\mathrm{1}}{\mathrm{\Delta }{\mathit{x}}_{\mathrm{0}}}$. The corresponding Nyquist frequency is f = f₀/ 2, which is the maximum frequency that allows the reconstruction of the transit chord without aliases. The Nyquist theorem thus suggests to divide the chord into cells whose minimum width is Δx_min = 2Δx₀. Plugging this into Eq. (5), the requirement on the spatial frequency translates into N<M/ 2.

We let the user of TOSC decide the optimal N for their purposes, cautioning that in our tests values of N larger than M/ 2 have proven to introduce aliases in the reconstruction.

Noise in the data, numerical precision and the fact that the system is overdetermined if N<M generally lead to the lack of an exact solution for Eq. (9). For these reasons, we solve Eq. (9) through a χ² minimization procedure rather than by using canonical linear methods, aiming at the vector $\begin{array}{c}{\mathrm{􏽥}}^{}\\ \mathit{x}\end{array}$ which best approximates Eq. (9) via the condition $\begin{array}{c}\mathrm{􏽥}\\ \mathit{x}\end{array}\mathrm{=}\underset{\mathit{x}}{\min }\mathrm{(}\mathrm{\parallel }\mathit{A}\mathrm{·}\mathit{x}\mathrm{-}{\mathit{F}}_{\mathrm{occ}}{\mathrm{\parallel }}^{\mathrm{2}}\mathrm{)}\mathit{,}$(10)where A is the M × N matrix of the coefficients in Eq. (9)$\mathit{A}\mathrm{=}\left[\begin{array}{c}\\ {\mathit{w}}_{\mathrm{1}}\mathrm{\left(}{\mathit{t}}_{\mathrm{1}}\mathrm{\right)}& & & & & \mathit{...}\hspace{0.17em}& & & & & {\mathit{w}}_{\mathit{N}}\mathrm{\left(}{\mathit{t}}_{\mathrm{1}}\mathrm{\right)}\\ \begin{array}{c}\mathit{.}\\ \mathit{.}\\ \mathit{.}\end{array}& & & & & & & & & & \begin{array}{c}\mathit{.}\\ \mathit{.}\\ \mathit{.}\end{array}\\ {\mathit{w}}_{\mathrm{1}}\mathrm{\left(}{\mathit{t}}_{\mathit{M}}\mathrm{\right)}& & & & & \mathit{...}\hspace{0.17em}& & & & & {\mathit{w}}_{\mathit{N}}\mathrm{\left(}{\mathit{t}}_{\mathit{M}}\mathrm{\right)}\end{array}\right]\mathrm{;}$(11)and ${\mathit{F}}_{\mathrm{occ}}\mathrm{=}\left[\begin{array}{c}\\ {\mathit{F}}_{\mathrm{occ}}\mathrm{\left(}{\mathit{t}}_{\mathrm{1}}\mathrm{\right)}\\ \begin{array}{c}\mathit{.}\\ \mathit{.}\\ \mathit{.}\end{array}\\ {\mathit{F}}_{\mathrm{occ}}\mathrm{\left(}{\mathit{t}}_{\mathit{M}}\mathrm{\right)}\end{array}\right]\mathit{.}$(12)To improve the capability of TOSC to converge to a physically acceptable solution, in TOSC there is the possibility to constrain the fit cell flux¹. In the assumption that the transit LCs are not affected by bright faculae, the maximum flux expected from the cells cannot exceed the photospheric flux. The user can thus constrain the maximum flux F_i,max corresponding to T_eff,max = T_eff + nσ_Teff, where σ_Teff is the measurement uncertainty on the stellar T_eff, while n is a free parameter set by the user. Setting n = 0, the maximum allowed temperature equals the nominal stellar T_eff; setting n to a very large number corresponds to putting no constraint on T_eff,max. Moreover, the user can also constrain the minimum flux F_i,min emerging from the grid cells by setting the minimum temperature T_eff,min for the reconstructed transit chord. To this purpose, the study of Berdyugina (2005) and Andersen & Korhonen (2015) provide some useful indications on the spots’ temperature contrast as a function of the stellar effective temperature.

Mathematically speaking, solving Eq. (10) is an ill-posed problem. Hence, we regularize the inversion problem using Tikhonov’s method, that is, we insert an additive term into Eq. (10) such that the best solution $\begin{array}{c}{\mathrm{􏽥}}^{}\\ \mathit{x}\end{array}$ is found as: $\begin{array}{c}\mathrm{􏽥}\\ \mathit{x}\end{array}\mathrm{=}\underset{\mathit{x}}{\min }\mathrm{(}\mathrm{\parallel }\mathit{A}\mathrm{·}\mathit{x}\mathrm{-}{\mathit{F}}_{\mathrm{occ}}{\mathrm{\parallel }}^{\mathrm{2}}\mathrm{+}\mathit{\lambda }\mathrm{\parallel }\mathit{T}\mathrm{·}\mathit{x}{\mathrm{\parallel }}^{\mathrm{2}}\mathrm{)}\mathit{,}$(13)where λ is a non-negative Lagrangian multiplier and T is, in the Tikhonov regularization framework, an ad-hoc matrix which gives preference to solutions with acceptable characteristics. We give preference to smooth solutions of Eq. (13) by minimizing the differences between neighboring elements of x, with the aim of providing solutions with minimum jumps between neighboring cells along the transit chord². $\mathit{T}\mathrm{=}\left[\begin{array}{c}\\ \mathrm{1}& -1& \mathrm{0}& \mathrm{0}& \mathit{...}\hspace{0.17em}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& \mathrm{1}& -1& \mathrm{0}& \mathit{...}\hspace{0.17em}& \mathrm{0}& \mathrm{0}\\ \begin{array}{c}\mathit{.}\\ \mathit{.}\\ \mathit{.}\end{array}& & & & & & \begin{array}{c}\mathit{.}\\ \mathit{.}\\ \mathit{.}\end{array}\\ \mathrm{0}& & & & \mathit{...}\hspace{0.17em}& \mathrm{1}& -1\end{array}\right]\mathit{.}$(14)The constraints and the regularization through the Lagrangian multiplier λ are useful to reduce the propagation of errors and numerical instabilities introduced by approximations and noisy data. A posteriori we also find that the regularization in the regression algorithm is able to reduce the autocorrelation in the residuals introduced by numerical approximations in the code (see Sect. 3). The side effect of the regularization is to smooth the solution of the reconstruction problem, lowering the spatial resolution of the reconstruction. In particular, at large λ the solution of Eq. (13) is over-smoothed, leading again to autocorrelated first neighbors. We will further discuss the impact of λ on the solution in Sect. 4.2, where we apply TOSC to a clear transit feature detected in real data.

In TOSC we implement an automatic search of the optimal multiplier λ. This search consists in increasing iteratively the value of λ until the absolute value of the autocorrelation between first neighbors decreases below the threshold $\mathrm{1.96}\mathit{/}\sqrt{\mathit{N}}$. This threshold corresponds to the 95% confidence interval for the hypothesis that residuals are randomly distributed around zero (Chatfield 1980). We remark that this is an automatic optimization criterion, which can be overridden by the user in case of unsatisfactory results.

Once the optimal λ multiplier is fixed, TOSC computes the array of fluxes $\begin{array}{c}{\mathrm{􏽥}}^{}\\ \mathit{x}\end{array}$ in Eq. (13), which is then converted into temperatures by means of the grid of fluxes predicted by the synthetic spectral models (Sect. 2.1). The temperature contrast ΔT discussed in the following sections is the difference between the reconstructed and the effective temperature, and is thus negative in the presence of cool spots.

3. Simulated transits

To test TOSC, we check its robustness and accuracy in reconstructing a simulated transit chord. In the following subsections, we simulate a planetary system with the parameters of the planet system hosted by HAT-P-11 (listed in Table 1), and we analyze the output of the reconstruction algorithm when different configurations of spots are simulated.

3.1. Transit curve simulation

There are several codes capable to simulate the LCs of planetary transits. These codes generally simulate a clean stellar disk, with no signature of stellar activity. Since in this work we are interested in analyzing the effect of photospheric spots on transit LCs, we write a simple code to simulate the planetary transits over spotted photospheres, to provide to TOSC some test cases. In brief, this code pixelates the stellar disk along two perpendicular directions, and for each pixel the emerging flux is computed convolving the photometric passband with the same limb-darkened stellar spectrum used in TOSC. The number of pixels along each direction is chosen one order of magnitude larger than the number of the simulated photometric points (see Sect. 3.2). A posteriori, we find that increasing the number of pixels does not have any effect on the reconstruction of the transit chord.

Then, the flux dimming for any position of the planetary disk over the stellar disk is computed. The integration of the occulted stellar flux is performed using the IDL routine aper.pro, which computes aperture photometry with the capability of counting subpixels. A posteriori, we find that this feature provides a major improvement in the computation of the LC, which in turn allows a more precise reconstruction of the transit chord.

3.2. Transit over a spot-free stellar disk

The simplest case to test TOSC is the transit of the planet over a spot-free stellar disk. To this purpose, we have simulated the transit of the planet uniformly sampling the transit chord with 100 points between $\mathrm{-}{\mathit{R}}_{\mathrm{\ast }}\mathrm{(}\sqrt{\mathrm{1}\mathrm{-}\mathrm{(}\mathit{b}\mathrm{+}{\mathit{R}}_{\mathrm{p}}{\mathrm{)}}^{\mathrm{2}}}\mathrm{-}{\mathit{R}}_{\mathrm{p}}\mathrm{)}$ and ${\mathit{R}}_{\mathrm{\ast }}\mathrm{(}\sqrt{\mathrm{1}\mathrm{-}\mathrm{(}\mathit{b}\mathrm{+}{\mathit{R}}_{\mathrm{p}}{\mathrm{)}}^{\mathrm{2}}}\mathrm{-}{\mathit{R}}_{\mathrm{p}}\mathrm{)}$ (see Fig. 1). This is the typical number of photometric points in transit LCs provided by dedicated space-borne observatories. Following the Nyquist sampling theorem (Sect. 2.3), this gives a number of cells equal to 50.

First, we neglected the LD effect and did not include measurement errors in the simulations. In the absence of LD, the occulted flux is a constant fraction ${\mathit{R}}_{\mathrm{p}}^{\mathrm{2}}$ of the stellar flux between the second and third contact. Nonetheless, we used the transit simulator discussed in Sect. 3.1 to check the inaccuracies introduced by numerical approximations in the transit simulation.

We ran TOSC on the simulated LC, setting λ = 0 for the regularization discussed in Sect. 2.3. We find that TOSC is not capable of returning exactly ΔT = 0 as expected. The best fit is systematically ~10^-3 K off the true value, and it shows a regular autocorrelated pattern with amplitude of the order of ~10^-4 K (Fig. 3, top panel). These features are due to the propagation of numerical approximations both in the LC simulation and in the transit chord reconstruction. We remark that these inaccuracies are negligible compared with the uncertainty usually associated with the measurement of stellar effective temperatures (≲100 K), and are thus irrelevant in any practical purpose.

We then simulated the more realistic case of a planetary transit over a limb-darkened stellar disk, using the quadratic LD law with the coefficients in Table 1. In this case, the unregularized algorithm was able to reconstruct the transit chord temperature profile with an uncertainty of ≲1 K (see Fig. 3, second panel from top). We thus find that the presence of LD increases the offset of the best fit solution with respect to the true value. This is due to the fact that the flux assigned to each grid cell is an approximation of the true one (see Eq. (6)).

Finally, we added some randomly-generated noise to the transit LC. First of all we simulated a 1 mmag photometric noise, typical of ground-based high-precision photometric observations. A posteriori, we find that the output of the regression is heavily affected by the propagation of the photometric uncertainties in the matrix inversion such that, in the most extreme cases, the automatic optimization of λ does not converge to a smooth solution. When it does, we find that the λ multiplier is of the order of ≲10 and the accuracy in the reconstruction of the transit chord is σ_T ~ 100 K (Fig. 3, second panel from bottom).

If we decrease the photometric noise down to 0.1 mmag, similar to Kepler’s precision (Borucki et al. 2010) and to what is expected from future photometers (e.g., CHEOPS, Scandariato et al. 2016), we find that the convergence of TOSC is more robust: the accuracy in the reconstruction is σ_T ≲ 30 K (Fig. 3, bottom panel).

Fig. 3 Examples of reconstruction of the thermal profile of the transit chord for the case of an unspotted stellar disk. The orbital parameters used in the simulation are listed in Table 1. From top to bottom, we show the case of (i) a stellar disk with no LD, noise-free LC and λ = 0; (ii) a stellar disk with LD, noise-free LC and λ = 0; (iii) a stellar disk with LD and LC with a photometric precision of 1 mmag and λ = 23; (iv) a stellar disk with LD and LC with a photometric precision of 0.1 mmag and λ = 2.9.

To better understand how λ and σ_T depend on the photometric uncertainty, we ran the following test: we simulated the five different values of the photometric uncertainty, from 0.1 mmag to 1 mmag, shown in Fig. 4. For each of these values, we generated 10⁵ random samples of photometric noise to add to the model LC. For each noisy LC we ran TOSC and we stored the corresponding λ and σ_T. For each value of the photometric noise we can thus compute the median λ and σ_T, which are plotted in Fig. 4. We find that λ is generally of the order of a few units or less and does not depend on the photometric accuracy, while σ_T increases with the photometric uncertainty of the LC, consistently with the test cases discussed above.

Fig. 4 λ (red squares) and σT (black circles) as functions of the photometric uncertainty of the LC.

3.3. Transit over a spotted stellar disk

Photospheric spots modify the transit LC in two ways. On the one hand, if spots are not crossed by the transiting planet, then the transit LC is deeper, because the planet occults a larger fraction of the stellar light (see e.g., Czesla et al. 2009; Ballerini et al. 2012). The effect is thus that the planet has an apparent radius larger than the true one. On the other hand, as discussed in the introduction of this paper, if spots are crossed during the planetary transit then the transit LC shows localized anomalies.

To test the robustness of TOSC in the presence of spots, we first simulated the case of a stellar disk with a spot not occulted by the planet during its transit over the disk. For this purpose, we used the same model of the HAT-P-11 system used above, and we simulated a circular spot 1000 K cooler than the photosphere and radius equal to 0.15 R_∗ (i.e., approximately three times larger than the planet). These parameters depict one of the largest and coolest spots found by Béky et al. (2014b) for HAT-P-11. The transit LC was simulated as above, and no random noise was added in order to test the net effect of unocculted starspots on the transit chord reconstruction.

Given this test case, the reconstructed transit chord computed by TOSC with λ = 0 is systematically 20 K hotter than the expectations. This is due to the fact that, assuming a spot-free photosphere, TOSC overestimates the fraction of stellar light occulted by the planet, and thus the temperature of the occulted stellar photosphere. To fix this bias, we fit the simulated LC using the JKTEBOP code version 34 (Southworth et al. 2004), obtaining an apparent planetary radius of 0.058965 R_∗, that is, larger than the true radius. We ran TOSC with this apparent radius, and we find that the simulated transit chord is reconstructed with an accuracy of ≲3 K, comparable with the accuracy of ~1 K discussed in Sect. 3.2. This indicates that the assumption of a slightly larger planet counteracts the effect of the neglected unocculted photospheric spots, thus improving the accuracy of the algorithm.

This test favors the use of the apparent radius of the planet rather than the true one, thus we advise the user of TOSC to fit the transit LC under examination and derive the apparent planetary radius for that given transit. This represents an advantage for TOSC, because the apparent planetary radius is easily measurable from the transit LC, while the true radius is not precisely measurable in case of active photospheres. Nonetheless, we remark that in the most optimistic scenario the LC is affected by a photometric uncertainty as small as 0.1 mmag, which introduces an uncertainty of ≲30 K in the transit chord reconstruction (Sect. 3.2). Thus, if the true planetary radius is used, then the corresponding bias (≲20 K in the simulated case) turns out to be negligible, unless more extremely spotted photospheres have to be taken into account.

Next, we tested the robustness of TOSC in the presence of transit anomalies given by spot-crossing events, which is the main goal of our algorithm. We first simulated the same photospheric spot as above placed at the center of the transit chord (Fig. 5, top panel), and we simulated the corresponding light curve, adding a random photometric noise of 1 mmag (Fig. 5, middle panel). We fit the LC assuming the apparent planetary radius obtained with JKTEBOP (0.058965 R_∗, see above), and show the result in the bottom panel of Fig. 5: the reconstructed chord obtained with λ = 3 (as returned by the automatic optimization) clearly shows a cooler interval corresponding to the anomaly in the LC, with width and temperature contrast consistent with the simulated spot within the uncertainties. For the sake of comparison, the bottom panel of Fig. 5 shows the average of the temperature inside the cells drawn over the simulated photosphere.

To compute the uncertainty σ_T on the temperature along the transit chord we proceeded in the following way. First of all, we simulated the transit LC with the same parameters of the given star and planet, assuming a clear photosphere: the fit of this noise-free LC would provide a noise-free temperature profile of the unspotted transit chord, with the exception of small numerical instabilities. Then, we computed the residuals of the best fit of the observed LC and we add them to the simulated noise-free LC. This gives us a simulated LC with the same photometric accuracy of the observed one. Finally, we ran TOSC on this noisy simulated LC using the same λ that was optimized for the observed LC, and we compute the standard deviation of the temperature of the grid cells as an estimate of σ_T. This is an average uncertainty on the reconstructed temperature profile, obtained with the same noise model of the data.

Fig. 5 Left panel: sketch of the simulated transiting system with a photospheric spot along the transit chord. The horizontal dashed lines delimit the transit chord. Stellar flux increases from orange to white. Right panel – top panel: simulated LC with the unperturbed model and the transit fit represented by the solid black line and the solid green line respectively. Middle: residuals of the transit fit. Bottom: the reconstructed transit chord. In the bottom panel, the uncertainty on the reconstruction is shown in dashes, while the solid green line represents the simulated temperature profile of the chord (see text for better details).

Fig. 6 As in Fig. 5, for the test case with a smaller and hotter spot as discussed in the text.

Finally, we considered a scenario in which the spot is smaller and hotter than that described above, in order to reduce the amplitude of the anomaly with respect to the transit profile depth. In particular, we simulate a circular spot with half the planetary radius and 500 K cooler than the photosphere, placed along the transit chord but off-center, in order to introduce also some asymmetry in the system. Given that this spot has less contrast with respect to the previous test case, its signal is below the 1 mmag photometric uncertainty simulated above. We thus added a smaller photometric noise of 0.1 mmag, which is slightly smaller than the signal expected from the spot-crossing event (see right panel in Fig. 6). We find that, even if the anomaly is at the limit of the photometric noise, the spot can be detected with a regularization parameter as small as λ = 1. Nonetheless, we remark that the spot-to-photosphere contrast (ΔT ≳ − 200 K) underestimates the simulated one (− 500 K). This is for two reasons. First of all, the regularization operated by the algorithm tends to smooth out the discontinuities in temperature profile of the chord. Secondly, and most importantly, in this case the spot does not cover the full height of the reconstructed chord, which is equal to the diameter of the planet (see Fig. 1). The algorithm is not capable of reconstructing the flux distribution perpendicular to the transit chord, then the net result is that along any cell in the transit chord (Fig. 1) the flux missing due to the cooler temperature of the spot is redistributed over the entire surface of the cell. As a consequence, the temperature of the cells returned by the algorithm is an upper limit of the spot temperature.

4. Comparison with previous results on real transits

In Sect. 3 we tested TOSC against a few simulated cases. In this section we compare it with other models that have been already discussed in the recent literature. To this purpose, we have applied TOSC to a few test cases that have been analyzed with similar reconstruction algorithms.

4.1. CoRoT-2

CoRoT-2 (GSC 00465-01282) is an active young G dwarf (T_eff = 5630 K) hosting the Hot Jupiter CoRoT-2b (R_p = 1.466R_J) which transits the stellar disk every ~1.74 days (Alonso et al. 2008). The planetary transit LCs observed with the CoRoT satellite (Baglin et al. 2006) show both deformations with respect to the unperturbed transit profile and changes in the transit depth. These features have been attributed to the ubiquitous presence of photospheric dark spots (Lanza et al. 2009; Silva-Valio et al. 2010; Silva-Valio & Lanza 2011; Bruno et al. 2016).

Wolter et al. (2009) analyze the most prominent transit anomaly in the whole photometric monitoring of CoRoT-2, found in the transit close to JD = 2 454 335.0, with the aim of modeling the spot crossed by the planet. The authors claim that the LC does not constrain the spot’s size, shape and contrast with respect to the photosphere. For this reason they limit their models to circular spots or circle segments enclosed in the transit chord. They also adopt two fixed values for the spot’s contrast (i.e., its temperature), consistently with the analysis of Lanza et al. (2009). The authors conclude that the spot is either ~350 K cooler than the photosphere and is of the same size of the planet, or it is ~1200 K cooler and is as small as half the planet size. They do not discuss the scenario of a large spot centered out of the chord.

We analyzed the same transit LC as that discussed by Wolter et al. (2009). We downloaded the latest reduction of the data from the CoRoT archive, we extracted the time interval around the analyzed transit, we kept the photometric points with good quality as flagged by the reduction pipeline, and we normalized the transit profile with a straight line fitting the LC right before the ingress and after the egress of the planet. The final LC is shown in the top panel of Fig. 7.

We fit the LC with TOSC assuming the parameters in Table 2 and using the photometric passband of CoRoT. We also adopt the apparent planetary radius (see Sect. 3.3) R_p = 0.1617 R_∗ as obtained fitting the LC with JKTEBOP after masking out the transit anomaly, in order to fit only the data points not affected by the intersected spot. A posteriori, we find that the Lagrangian multiplier selected by TOSC (λ = 1.08) does not give a smooth reconstruction of the transit chord. This is likely to be an effect of the photometric uncertainty of 1.5 mmag of the LC (Sect. 3.2). We thus artificially set λ = 6.25 and find clear evidence of an area ~300 K cooler than the photosphere slightly after the mid transit and extending over ≃0.3 R_∗ in the projected plane.

As discussed in Sect. 3.3, we remark that TOSC does not provide any information on the flux distribution inside each individual grid cell and that, as a consequence, it provides only an average of the temperature contrast in each cell. Thus we cannot exclude that the cells are only partly covered by the spot, and that the spot’s temperature is lower than the one returned by TOSC. Thus we conclude that the spot crossed by CoRoT-2b is ≳300 K cooler than the photosphere, consistently with the results of Wolter et al. (2009).

Fig. 7 Top panel: CoRoT-2 LC discussed in the text. The unperturbed model and the transit fit are represented by the solid black line and the solid green line respectively. Middle: residuals of the transit fit. Bottom: the reconstructed transit chord with the uncertainty on the reconstruction shown in dashes.

Fig. 8 As in Fig. 7 for the HAT-P-11 system discussed in the text. The left and right panels show the result for the transits occurring on BJD = 2 454 967.6 and BJD = 2 455 671.4 respectively.

4.2. HAT-P-11

HAT-P-11 is a V = 9.6 K4 dwarf (T_eff = 4780 K) in the Kepler field (Borucki et al. 2010), orbited by a hot Neptune every 4.9 days (Bakos et al. 2010). Béky et al. (2014b) aim at analyzing the transits observed by the Kepler telescope. To this purpose, they develop SPOTROD, a Monte Carlo algorithm which models the transit anomalies assuming a number of circular photospheric spots. In particular, in their paper they show extensively the results for the transits occurring at BJD = 2 454 967.6 and BJD = 2 455 671.4. In the first case, they find evidence of a spot entirely lying under the transit chord located ~0.31 R_∗ ahead the mid-transit, with radius ≃0.05 R_∗ and spot-to-photosphere flux ratio of ≃0.4. For the second transit, the authors find two degenerate models for the crossed spot, both with radius ≃0.23 R_∗ and flux ratio of ~0.8. The two degenerate models are centered off the transit chord, ~0.4 R_∗ after the mid-transit, and they basically differ by the fact that the two circular spots are centered either above or below the chord. Based on the geometry of their model, we derive that the length of the intersection between the spot and the chord in the projected plane is 0.08 R_∗ and 0.3 R_∗ in the two transits respectively.

For the sake of comparison, we ran TOSC on the same two transit LCs discussed above, using Kepler’s passband. The parameters of the planetary system needed as input are listed in Table 1. The apparent planetary radius obtained with JKTEBOP (see Sect. 3.3) is 0.06006 R_∗ and 0.06063 R_∗ for the two transits respectively. The output of the algorithm is shown for both cases in Fig. 8.

For the first transit TOSC optimizes λ = 0.34. The reconstructed anomaly is centered at ≃− 0.3 R_∗ and is d ≃ 0.15 R_∗ wide. The fit spot-to-photosphere temperature contrast is ΔT ≃ − 900 K, which corresponds to a flux ratio of 0.53.

For the second transit TOSC returns λ = 0.89. The reconstructed anomaly is centered at ≃0.4 R_∗ and its width is ≃0.3 R_∗. The fit spot-to-photosphere temperature contrast is ΔT ≃ − 350 K, which leads to a flux ratio of 0.8.

In both cases our results are thus generally consistent with what is obtained using SPOTROD. Only for the first transit we find some discrepancy in the spot-to-photosphere flux ratio. A closer look at the best fit of the LC returned by TOSC reveals that the algorithm systematically underestimates the observed flux at the center of the transit anomaly by more than the photometric uncertainty. This is due to the smoothing of the best-fit model introduced by the regularizing constraints in the algorithm. As a matter of fact, the best fit returned by TOSC follows more closely the data if we artificially reduce λ down to 0.01. Correspondingly, the temperature contrast decreases to ΔT ≃ − 1200 K, which corresponds to a flux ratio of ~0.38, consistent with the results of Béky et al. (2014b). Nonetheless, the residuals show some degree of autocorrelation, introduced by the inversion algorithm (see Sect. 3.2). The purpose of the automatic optimization discussed in Sect. 2.3 is to find the smallest Lagrangian multiplier which returns non autocorrelated residuals. In this particular case, the best Lagrangian multiplier is λ = 0.34. If, on the other hand, we artificially use large multipliers, the best fit tend to be over-smoothed and, correspondingly, the residuals get more autocorrelated as λ increases. In the left panel of Fig. 9 we plot the best fit obtained with three different values of λ, while in the right panel we plot the autocorrelation of the residuals computed at different values of λ.

Fig. 9 Left panel: zoom in of the left panel in Fig. 8, showing the effect of different values of λ in the reconstruction algorithm. Right panel: autocorrelation ACF1 between neighboring cells in the residuals vs. λ. The three examples shown in the left panel are marked in colors, while the gray dashes represent the $\mathrm{1.96}\mathit{/}\sqrt{\mathit{N}}$ threshold adopted in Sect. 2.3 (see text for discussion).

4.3. HAT-P-36

HAT-P-36 is a V = 12.5 G5 dwarf (T_eff = 5620 K) orbited by a hot Jupiter (R_p = 1.3R_J) approximately every 1.3 days (Bakos et al. 2012). In their detailed study of the system, Mancini et al. (2015) analyzed the anomalies in the transit LCs of the system. They find the most evident anomalies in the transit LC in the Gunn r band, collected on April, 15th 2013 with the BFOSC imager mounted on the 1.52 m Cassini Telescope at the Astronomical Observatory of Bologna in Loiano, Italy. They use the GEMC+PRISM codes (Tregloan-Reed et al. 2013, 2015), which basically combines Monte Carlo and genetic algorithms to simulate circular spots on the stellar surface. In their analysis of the aforementioned LC, they find evidence of two spots (one after the second contact and the other before the third contact) with T ≃ 5000 K and radius ≃0.2 R_∗. In particular, the planet disk does not enter completely into the spots or, in other words, the spots do not completely cover the occulted belt along the vertical direction. Given the geometrical model used by the GEMC+PRISM algorithm, we derive that the length of the intersection of the chord with each of the stellar spots is ≲0.16 R_∗ in the projected plane.

Fig. 10 As in Fig. 7 for the transit LCs of the HAT-P-36 system discussed in the text.

We analyzed the normalized LC provided by Mancini et al. (2015), using TOSC with the parameters listed in Table 3. We do not list the argument of periastron ω as it is meaningless in case of circular orbits. We also used the apparent radius R_p = 0.1273 R_∗ obtained fitting the LC with JKTEBOP. The output of the transit chord reconstruction, obtained with regularization λ = 1.08, is shown in Fig. 10. We find that the two spot-crossing events occur at the borders of the reconstructed chord. Nonetheless, our fit suggests that the two sections of the projected transit chord intersected by the spots are ≳0.12 R_∗, while the temperature contrast is between ΔT = − 100 K and ΔT = − 200 K, consistently different from what Mancini et al. (2015) found (ΔT ≃ − 600 K). This may be due to the fact that, according to Mancini et al. (2015), the two spots do not intersect the whole height of the occulted belt. As we discussed in Sect. 3.3, when the grid cells are only partly covered by the spot, the temperature returned by TOSC is an intermediate value between the T_eff and the true temperature of the spot T_s. Let us assume that the fit flux F is given by the weighted mean $\mathit{F}\mathrm{=}\mathrm{(}\mathrm{1}\mathrm{-}{\mathit{\alpha }}^{\mathrm{)}}{\mathit{F}}_{\mathrm{\ast }}\mathrm{+}\mathit{\alpha }{\mathit{F}}_{\mathrm{s}}\mathit{,}$(15)where α is the fraction of the grid cells covered by the spot, while F_∗ and F_s are the temperature of the quiet photosphere and the spot respectively. After performing some algebra, we obtain $\frac{{\mathit{F}}_{\mathrm{s}}}{{\mathit{F}}_{\mathrm{\ast }}}\mathrm{=}\mathrm{1}\mathrm{-}\frac{\mathrm{1}\mathrm{-}\mathit{F}\mathit{/}{\mathit{F}}_{\mathrm{\ast }}}{\mathit{\alpha }}\mathit{,}$(16)where $\frac{{\mathit{F}}_{\mathrm{s}}}{{\mathit{F}}_{\mathrm{\ast }}}$ is the spot-to-photosphere flux contrast and $\frac{\mathit{F}}{{\mathit{F}}_{\mathrm{\ast }}}\mathrm{\simeq }\mathrm{0.8}$ is the contrast corresponding to the ΔT ≃ 150 K returned by TOSC. Following the scenario depicted by Mancini et al. (2015) in their Fig. 8, if we assume α ≲ 0.5 then we obtain $\frac{{\mathit{F}}_{\mathrm{s}}}{{\mathit{F}}_{\mathrm{\ast }}}\lesssim \mathrm{0.6}$, which corresponds to a temperature contrast of ≲− 400 K, closer to their result.

4.4. GJ1214

GJ1214 is a M4.5 red dwarf orbited by a super-Earth (R_p ≃ 2.6 R_⊕) every 1.58 days (Charbonneau et al. 2009). On May,17th 2012, Nascimbeni et al. (2015) observed a planetary transit simultaneously in the B and R photometric bands using the LBC camera mounted at the double 8.4 m Large Binocular Telescope. The photometric uncertainty of the observed LCs is ≲1 mmag. Both transit LCs show a clear signature of the crossing of a spot near the egress of the planet. Analyzing the shape (width and height) of the transit anomaly in the two bands, the authors conclude that the projected spot is as large as the planet, and is ~110 K cooler than the photosphere.

Fig. 11 As in Fig. 7 for the transit LCs of the GJ1214 system discussed in the text. The left panel shows the fit of the B band LC, obtained with regularization parameter λ = 0.42. The right panel shows the fit of the R band LC, obtained with regularization parameter λ = 0.19.

Nascimbeni et al. (2015) also analyze the residuals with respect to the best fit, dividing the LCs in three non overlapping sections: the off-transit, the spot (corresponding to the transit anomaly) and the in-transit (which does not include the anomaly). They find no correlation between the B and R off-transit residuals, which is indicative of the fact that no correlated noise (due to e.g., instrumental or atmospheric effects) affects the data. Conversely, in the transit anomaly the B and R residuals show a clear correlation, as expected when a large spot is occulted. Most interestingly, the same correlation is present in the in-transit residuals, suggesting that the whole transit chord is characterized by inhomogeneities.

We have used TOSC to analyze in further detail the same B and R LCs discussed above. In Table 4 we list the parameters injected in the algorithm. As in Sect. 4.3, we do not list the argument of periastron ω as it is meaningless in case of circular orbits. We also use the apparent planetary radius derived by Nascimbeni et al. (2015), that is 0.1186 R_∗ and 0.1178 R_∗ in the B and R bands respectively. The output of the reconstruction is shown in Fig. 11. Unfortunately, the spot-crossing event occurs too close to the egress and TOSC is not able to entirely model it. Nonetheless, both reconstructed transit chords clearly show a photospheric spot ≲80 K cooler than the photosphere near the egress.

The transit LCs also show that after the mid-transit and before the anomaly the transit depth increases by ~0.2% in the B band and ~0.1% in the R band. This suggests that the part of the stellar disk occulted by the planet after mid-transit is brighter than the average stellar disk. Accordingly, the temperature contrast returned by TOSC suggests that the second half of the transit chord is hotter by ~50 K than the quiet photosphere. The increase in brightness may be due to a large group of photospheric faculae, extending for ~0.5 R_∗ and neighboring the spot detected near the stellar limb. This scenario is consistent with the solar case, where spots are often accompanied by faculae (Wiegelmann et al. 2014). Moreover, we remark that the brightening effect seems to increase approaching the stellar limb, as indicated by Foukal et al. (1991), Solanki et al. (2006), and that it is more pronounced in the B band than in the R band, consistently with the fact that the facular contrast increases toward shorter wavelengths (Chapman & McGuire 1977; Ermolli et al. 2007).

We remark that TOSC is only suited for modeling cold spots on the stellar surface, that is, it is not capable to model warm regions such as faculae or plages, for which a proper spectrum should be taken into account. Even though the numerical value of the temperature contrast is unreliable, still the best fit of the LCs clearly indicates an extended warm section of the transit chord, consistently between the B and R data sets. The presence of this region supports the hypothesis of Nascimbeni et al. (2015), which postulated the presence of additional inhomogeneities on the stellar surface to explain the correlated noise in the in-transit residuals. The presence of faculae in the photosphere of GJ1214 has also been recently postulated by Rackham et al. (2017).

Another interpretation of the analyzed transit LC is that cool spots are ubiquitous along the transit chord. In this scenario, even masking out the most evident transit anomaly, the lower envelope of the LC is shallower than the LC expected in the case of a quiet transit chord, that is, the apparent planetary radius is underestimated (see Sect. 3.3). If we arbitrarily increase the apparent R_p used to reconstruct the transit chord by 5%, then the ΔT is smaller (i.e., the transit chord is cooler) by ~20 K, i.e., there are no faculae in the transit chord and the temperature contrast is of the order of ~− 100 K, as found by Nascimbeni et al. (2015).

5. Conclusions

In this paper we present TOSC, a new fast algorithm which analyzes the anomalies in the transit LCs with the aim of reconstructing the flux distribution, and the temperature profile, of the transit chord. The only inputs needed by the code are the parameters of the planetary system, the LD coefficients and the model spectrum for unspotted and spotted photosphere. The presence of unocculted spots, not reconstructed by TOSC, can be accounted for using the apparent radius of the planet in place of the true radius.

We test the performances of TOSC using simulated transits. We show that TOSC is generally able to reconstruct the crossed spots for a photometric accuracy better than 1 mmag, and that the temperature contrast is returned with an uncertainty of ≲100 K. We also test TOSC with real transits available in literature which have been already analyzed with alternative codes. These comparisons show that TOSC is consistent with previous approaches, even in the case of more sophisticated algorithms.

TOSC is available as a web interface³, where the user can run the algorithm feeding a few input files and retrieving the output. The web page contains also extensive instructions and some example data. The source code is also available for download.

Acknowledgments

G.S., V.N. and I.P. acknowledge financial support from “Accordo ASI–INAF” No. 2013-016-R.0 July 9, 2013. G.S. dedicates this paper to his alive and kicking father, after whom TOSC is named.

References

All Tables

All Figures

	Fig. 1 Schematic representation of the geometry discussed in the text. The circle represents the transiting planet, while the grid rectangles show the cells used in the chord reconstruction.
In the text

	Fig. 2 Graphical representation of the computation of the weights discussed in the text. The circle is the planetary disk centered in O, while the gray grid shows the cells drawn on the stellar disk partially covered by the planetary disk.
In the text

Fig. 3 Examples of reconstruction of the thermal profile of the transit chord for the case of an unspotted stellar disk. The orbital parameters used in the simulation are listed in Table 1. From top to bottom, we show the case of (i) a stellar disk with no LD, noise-free LC and λ = 0; (ii) a stellar disk with LD, noise-free LC and λ = 0; (iii) a stellar disk with LD and LC with a photometric precision of 1 mmag and λ = 23; (iv) a stellar disk with LD and LC with a photometric precision of 0.1 mmag and λ = 2.9.

In the text

	Fig. 4 λ (red squares) and σT (black circles) as functions of the photometric uncertainty of the LC.
In the text

Fig. 5 Left panel: sketch of the simulated transiting system with a photospheric spot along the transit chord. The horizontal dashed lines delimit the transit chord. Stellar flux increases from orange to white. Right panel – top panel: simulated LC with the unperturbed model and the transit fit represented by the solid black line and the solid green line respectively. Middle: residuals of the transit fit. Bottom: the reconstructed transit chord. In the bottom panel, the uncertainty on the reconstruction is shown in dashes, while the solid green line represents the simulated temperature profile of the chord (see text for better details).

In the text

	Fig. 6 As in Fig. 5, for the test case with a smaller and hotter spot as discussed in the text.
In the text

	Fig. 7 Top panel: CoRoT-2 LC discussed in the text. The unperturbed model and the transit fit are represented by the solid black line and the solid green line respectively. Middle: residuals of the transit fit. Bottom: the reconstructed transit chord with the uncertainty on the reconstruction shown in dashes.
In the text

	Fig. 8 As in Fig. 7 for the HAT-P-11 system discussed in the text. The left and right panels show the result for the transits occurring on BJD = 2 454 967.6 and BJD = 2 455 671.4 respectively.
In the text

Fig. 9 Left panel: zoom in of the left panel in Fig. 8, showing the effect of different values of λ in the reconstruction algorithm. Right panel: autocorrelation ACF1 between neighboring cells in the residuals vs. λ. The three examples shown in the left panel are marked in colors, while the gray dashes represent the $\mathrm{1.96}\mathit{/}\sqrt{\mathit{N}}$ threshold adopted in Sect. 2.3 (see text for discussion).

In the text

	Fig. 10 As in Fig. 7 for the transit LCs of the HAT-P-36 system discussed in the text.
In the text

	Fig. 11 As in Fig. 7 for the transit LCs of the GJ1214 system discussed in the text. The left panel shows the fit of the B band LC, obtained with regularization parameter λ = 0.42. The right panel shows the fit of the R band LC, obtained with regularization parameter λ = 0.19.
In the text