R Rougeot

R Flamary

D Galano

C Aime

Performance of the hybrid externally occulted Lyot solar coronagraph Application to ASPIICS

Keywords: corona -Instrumentation, high angular resolution -Method, analytic

Context. High-contrast hybrid coronagraphs, which combine an external occulter and a Lyot-style coronagraph became a reality in recent years, despite the lack of analytic and numerical end-to-end performance studies. The solar coronagraph ASPIICS which will fly on the future ESA Formation Flying mission Proba-3 is a good example of such a hybrid coronograph. Aims. We aim to provide a numerical model to compute theoretical performance of the hybrid externally occulted Lyot-style coronagraph, which we then aim to compare to the performance of the classical Lyot coronagraph and the externally occulted solar coronagraph. We will provide the level and intensity distribution of the stray light, when the Sun is considered as an extended source. We also investigate the effect of different sizes for the internal occulter and Lyot stop. Methods. First, we have built on a recently published approach, to express the diffracted wave front from Fresnel diffraction produced by an external occulter at the entrance aperture of the coronagraph. Second, we computed the coherent propagation of the wave front coming from a given point of the Sun through the instrument. This is performed in three steps: from the aperture to the image of the external occulter, where the internal occulter is set, from this plane to the image of the entrance aperture, where the Lyot stop is set, and from there to the final image plane. Making use of the axis-symmetry, we considered wave fronts originating from one radius of the Sun and we circularly average the intensities. Our numerical computation used the parameters of ASPIICS. Results. The hybrid externally occulted Lyot coronagraph rejects sunlight below 10 -8 B from 1.3R -in the particular configuration of ASPIICS. The Lyot coronagraph effectively complements the external occultation. We show that reducing the Lyot stop allows a clear gain in rejection, being even better than oversizing the internal occulter, that tends to exclude observations very close to the solar limb. As an illustration, we provide a graph that allows us to estimate performance as a function of the internal occulter and Lyot stop sizes. Conclusions. Our work consists of a methodological approach to compute the end-to-end performance for solar coronagraph.

Introduction

The Sun's corona consists of a fully ionized plasma, with a strong magnetic field. Its physical structure and dynamics are governed by multiple processes, theoretical models of which still need to be assessed and investigated. Indeed, its properties, such as density of the plasma, temperature, and magnetic field structures, are far more complex than any other planet magnetosphere, as described by [START_REF] Cox | Allen's astrophysical quantities[END_REF]. Coronal mass ejections (CMEs), heating processes operating in the corona, and even solar wind interaction and acceleration are still not perfectly understood. The active study of the corona of the Sun needs both simultaneous and complementary multi-wavelength observations. A very high angular resolution, of the order of the arcsecond, is required to constrain the finest coronal structures [START_REF] Zhukov | [END_REF]Peter et al. 1965), as well as sporadic events such as CMEs.

Observing the solar corona in white light requires perfect eclipse conditions, because the coronal brightness in this spectral band is much fainter than the halo of diffraction produced by the Sun, typically from 10 -6 B to 10 -10 B , where B is the mean solar brightness (Allen 2005). By creating artificial eclipses, the first Lyot solar coronagraph was a breakthrough for the study of the solar corona (Lyot 1939;Dollfus 1983). The development of the external occultation technique (Evans 1948) coupled with advanced stray light rejection concepts, such as toothed or multiple discs (Newkirk & Bohlin 1965;[START_REF] Purcell | Report of NRL Progress[END_REF], and the advent of space-borne instruments considerably improved the performance of solar coronagraphs, as described in the review paper by Koutchmy (1988). The solar coronagraph LASCO C2 (Large Angle Spectroscopic Coronagraph) of the Solar and Heliospheric Observatory mission can be held as a representative and successful example. This instrument combined an external occulter made of multiple discs and a Lyot-style solar coronagraph to achieve a sufficient rejection of residual diffracted sunlight, and managed to observe the corona in white light beyond 1.5R , where R is the radius of the Sun, with a resolution of 11.4arcsec per pixel (Brueckner et al. 1995). However, historically, observing the solar corona very close to the solar limb, where instrumentally scattered sunlight usually predominates, has never been successful without an eclipse of the Sun by the Moon.

The development of Formation Flying space missions will pave the way for new advanced concepts of space instrumentation by virtually enlarging instruments to unprecedented size, especially in coronagraphy. The solar coronagraph ASPIICS (Association de Satellites Pour l'Imagerie et l'Interférométrie de la Couronne Solaire) described in Lamy et al. (2010) and [START_REF] Renotte | SPIE[END_REF] takes advantage of the future ESA Formation Flying mission Proba-3, and is split between two spacecraft. The concept of the optical design, detailed in Galy (2015), is similar to SOHO LASCO C2. Its 1.42m diameter external occulting disc is mounted on the Occulter Spacecraft while the Lyot-style solar coronagraph of 50mm diameter aperture is carried by the Coronagraph Spacecraft positioned 144m behind. Its theoretical angular resolution is 2.77arcsec -with Rayleigh criterion at 550nm. Such a large size is new in the domain of solar coronagraphy and is the main feature of ASPIICS. The hybrid coronagraph is expected to be able to observe the solar corona in white light from 1.08R (1036arcsec) to 3R (2880arcsec) (Lamy et al. 2010), revealing fine scale structures very close to the solar limb.

As already mentioned, solar coronagraphy is mainly constrained by the halo of diffraction from the direct sunlight, which limits any observations and drives the performance. Sensitive analysis of this particular stray light is rather complex as it includes an extended light source -the Sun -which makes such a work much more difficult than considering one single point source, as is done in the particular domain of exoplanet coronagraphy where numerous analytic studies have been donesee for instance (Cash 2006;Vanderbei et al. 2007;Flamary et al. 2014). In solar coronagraphy, only a few numerical and analytical studies of rejection performance have been published (Lenskii 1981;Aime et al. 2002;Aime 2007Aime , 2013)). We also note that the purely analytical studies of Ferrari (2007) and Ferrari et al. (2010) come with the drawback that the Lyot stop must equal the entrance pupil and cannot be reduced. There are also some experimental approaches, such as the works described in Fort et al. (1978), [START_REF] Bout | [END_REF], Venêt et al. (2010) and Landini et al. (2010). In contrast, extensive and complete analytic and numerical analysis appears nowadays mandatory, for modern advanced high-contrast instrumentation. We address this need by presenting here a general study on the performance of solar coronagraphic systems. We compute the global response of the hybrid externally occulted Lyot coronagraph, that we compare to the classical Lyot coronagraph and the external solar coronagraph. We also investigate the impact of the size of the Lyot mask and the stop on stray light rejection.

The paper is organized as follows. The model and the framework adopted for this study are given in Sect. 2. The mathematical wave propagation into the Lyot-style coronagraph is derived in Sect. 3, standing as a new computation. The comparison of the response of the different coronagraphic systems and further analysis on sizing both Lyot mask and stop are discussed in Sect. 4. Conclusions are given in Sect. 5.

Model of the coronagraph

Presentation of the model

The classical Lyot coronagraph is made of four key planes representing the instrument. In a previous theoretical study of this system by Aime et al. (2002), these planes are denoted as A (entrance aperture), B (focal plane), C (image of the entrance . We note that we distinguish the Lyot mask, denoting the occulting disc set in plane B, and the internal occulter, being in plane O' of the hybrid coronagraphic system S EL . The four systems S Ø , S L , S E and S EL include the same circular entrance aperture in plane A. Table 2 summarizes the description of the four systems.

Our model is generic, and provides a methodological study on performance of such coronagraphic systems when considered as perfect. In order to have a realistic set-up, we have used the parameters of ASPIICS coronagraph (Lamy et al. 2010;[START_REF] Renotte | SPIE[END_REF] for the numerical computation. The external occulting disc of radius R = 710mm is located at z 0 = 144.348m before the 50mm diameter entrance aperture. The telescope consists of a converging lens of focal length f = 330.385mm. The Sun is assumed to be at infinity. Its angular radius is R = 0.0046542rad as seen from the centre of the aperture, so ∼ 960arcsec. The angular radius of the external occulter is 1.0568R , as viewed from the centre of the entrance aperture. Table 3 summarizes the numerical parameters.

The radius of the Lyot mask set in plane B will be given in solar units R , since this plane is the conjugate of the solar disc. However, the internal occulter is set in plane O' which is the conjugate image of plane O. We will thus speak in terms of units of external occulter image. A simple proportional relationship applies here to convert this particular unit system to solar units, or metric units if needed. In plane O', the image of the external occulter radius R corresponds to 1.0568R , so to 1.629mm. We emphasize that using solar unit has no real meaning in plane O', since it is not conjugated with the Sun. Finally, as plane C is the conjugate image of plane A, the dimension of the Lyot stop will be given in units of the image of the entrance pupil, meaning that a Lyot stop of 1.00 has the exact same size as the image of the radius of the pupil in plane C, so it corresponds to R p = 25mm.

Analytic and numerical framework

All the planes previously defined are assumed to be perfectly parallel and perpendicular to the optical axis, so that the geometry is axis-symmetric. We note that our model is general enough to cover transverse off-sets, but this would require further computations that are left to future works. To each plane we set a (r, θ, z) cylindrical coordinate system. The z-axis refers to the optical axis, oriented positively towards the detection plane.

The corresponding Cartesian coordinate system (x, y, z) is defined by x = r cos θ and y = r sin θ. In the remainder of the article, we will sometimes use both coordinates simultaneously, because this slight abuse of notation allows more compact and readable equations. To provide a better understanding, we will use as subscript the letter O, A, B, O', C or D referring to the corresponding plane for every quantity.

To model the perfect sharp-edged disc, the transmission in plane O is a radial gate function τ(r) = 0 if r ≤ R and τ(r) = 1 else. The Lyot mask (internal occulter) in plane B (plane O') is similarly modelled. The entrance pupil in plane A is a perfect circular aperture of radius R p = 25mm.

Our study uses monochromatic light, here λ = 550nm. We have adopted a Fresnel regime to describe diffraction induced by the external occulter, as suggested by the large value of Fresnel number N f = R 2 /λz 0 = 6350 [START_REF] Born | Principles of Optics[END_REF]. The analytic propagation of wave front is based on paraxial Fourier optics formalism [START_REF] Goodman | Introduction to Fourier Optics[END_REF]. Under this assumption, Fresnel freespace propagation of a wave front Ψ 0 (x, y) over a distance z is written as convolution product. The complex amplitude Ψ z (x, y)

Ψ z (x, y) = Ψ 0 (x, y) 1 iλz exp iπ x 2 + y 2 λz = ϕ z (r) iλz × F λz Ψ 0 (x, y) × ϕ z (r) , (1)
where ϕ z (r) = exp(iπr 2 /λz), and r = x 2 + y 2 is the transverse radius in cylindrical coordinates. F λz denotes the 2D Fourier transformation with spatial frequencies u = x/λz and v = y/λz. Equation ( 1) is the so-called Fourier-Fresnel transformation of the function Ψ 0 (x, y), where the phase term exp(2iπz/λ) for the longitudinal propagation has voluntary been omitted. Moreover, in the Fourier formalism, a converging lens of focal length f is modelled by the quadratic phase factor ϕ -f (r) = exp(-iπr 2 /λ f ).

Propagating through a lens consists of multiplying the complex amplitude of the incoming wave front by ϕ -f (r). A well known result is the propagation to the focal plane of a lens, that is, z = f in Eq. (1). In this case, the two quadratic phase factors ϕ -f (r) and ϕ + f (r) cancel each other out. Consequently, the wave in the focal plane is directly proportional to the Fourier transformation of the incoming wave at its entrance, to a scale factor λ f and a quadratic term ϕ f (r) that is canceled when computing the intensity.

The Sun is modelled by a collection of incoherent point sources. The global response of any system is given by the incoherent sum of their respective elementary intensities. Every point source is identified by a set of angular coordinates (α, β) on the sky, with α 2 + β 2 ≤ R . As it will be discussed in Sec. 4.1, the number of point sources needs to be carefully fixed to meet Shannon criteria for the numerical sampling of the Sun. To provide a better understanding, we will also use as subscripts the coordinates (α, β) to refer to a precise point source for every quantity. As the Sun is at infinity, the light coming from every point source is modelled by tilted planar wave, whose unitary complex amplitude is written as Ψ ,α,β (x, y) = exp(-2iπ/λ(αx + βy)).

We use the center-to-limb variation of the Sun B(α, β) from [START_REF] Hamme | [END_REF], given in Eq. ( 2), to model the non-uniformity of the viewed brightness of the solar disc. This choice has been driven by the need to have a representative limb-darkening function for the specific wavelength λ = 550nm.

B(ρ) = 1 -0.762 1-1 -ρ 2 -0.232 1-ρ 2 log 1-ρ 2 , (2)
where ρ = α 2 + β 2 /R is the angular radial coordinate on the solar disc, expressed in solar units.

Fresnel diffraction by the external occulter

Fresnel diffraction produced by a sharp-edged disc has already been described by Aime (2013) (his Eq. ( 5)). We will briefly recapitulate the main results of interest. The tilted planar wave front Ψ ,α,β (x, y) coming from the point source at (α, β) arrives onto the occulter in plane O, and then propagates in free-space. Related complex amplitude Ψ A,α,β of the wave front arriving on plane A is thus written as a Fourier-Fresnel transformation. 0,0 (x + αz 0 ,y + βz 0 ),(3) where

Ψ A,α,β (x, y) = Ψ ,α,β (x, y) × τ 1 iλz 0 exp iπ x 2 + y 2 λz 0 = T α,β (x, y) × Γ α,β × Ψ A,
T α,β (x, y) = exp -2iπ αx + βy λ (Tilt) Γ α,β = exp -iπ (α 2 + β 2 )z 0 λ (Offset) Ψ A,0,0 (x, y) = 1- 1 iλz 0 η 2 +ξ 2 ≤R 2 exp iπ (x -ξ) 2 +(y -η) 2 λz 0 dξdη, (4)
with ξ, η the Cartesian variables for integration over the occulting disc. As a result, an off-axis point source produces the same complex amplitude as the on-axis point source -Ψ A,0,0 -but shifted of the quantity (z 0 × α, z 0 × β) towards negative (x, y) directions. The constant phase term Γ α,β accounts for the offset of position, and the original tilt T α,β of the wave is conserved. Let us now consider the particular case of the on-axis point source. Taking advantage of the cylindrical symmetry, we naturally change for polar coordinates (r, θ). Equation ( 4) is then written as a radial Hankel transformation:

Ψ A,0,0 (r) = 1 - ϕ z 0 (r) iλz 0 R 0 2πρ exp iπ ρ 2 λz 0 J 0 2π rρ λz 0 dρ, ( 5
)
where ρ is the radial variable for integration over the disc, r is the transverse radial coordinate on plane A, ϕ z 0 (r) = exp(iπr 2 /λz 0 ) and J 0 (r) is the Bessel function of the first kind. Equation ( 5) is the exact analytic expression of Fresnel diffraction in the particular case of the on-axis point source. The computation of the Hankel transformation Ψ A,0,0 (r) (Eq. ( 4)) remains a delicate operation, as described by Lemoine (1994). We chose to use NIntegrate in Mathematica (Wolfram 2012), since it has been proved to give sufficient numerical precision (Aime 2013). An analytic expression using the Lommel series can alternatively be used. In Fig. 2 (curve (a)), we plot Fresnel diffraction pattern |Ψ A,0,0 (r)| 2 , known as the bright spot of Arago, for the 710mm diameter disc at the distance z 0 = 144.348m, in a logarithmic scale. The distance between the first zeroes of the Arago spot is approximately 1.53λz 0 /2R, being 171µm in our case. Using an analogy to Shannon's criteria, the radial sampling must be much tighter than half of this value. We voluntary chose to over-sample at 0.1µm, for better 2D-interpolation. In Fig. 2, we also show the central spot in linear scale (curve (b)), whose peak intensity is 1 for r = 0 as expected, and the transition zone between shadow and light (curve (c)). We note that the intensity decreases below 10 -4 for larger values of r.

Propagation through the coronagraph

The propagation of every wave front consists in a coherent process through each successive plane of the coronagraph. We successively derive Ψ B , Ψ O , Ψ C , and finally Ψ D , for each system S Ø , S L , S E , or S EL . However, the observed response is the result of the incoherent summation of every elementary intensity. So, let us consider one point source located at (α, β) on the solar disc.

In this section, we will voluntary omit the subscripts (α, β) for a better readability, and we will preferably use polar coordinates (r, θ) rather than the Cartesian coordinates (x, y).

Classical Lyot coronagraph

As already described, the classical Lyot coronagraph S L is modelled by planes A, B, C and D (Aime et al. 2002), and does not include the external occulter. At the entrance aperture, Ψ A corresponds then to a simple tilted planar wave front. Here, the coronagraph acts as a mere imaging system, adding the Lyot mask in the focal plane. Using the approach of Fourier formalism described in Sec. 2.2, the propagation process through the whole instrument consists in scaled Fourier transformations between each of the successive planes, that is, from A to B, from B to C and from C to D. The images are of different sizes, depending on the lenses used for imaging, but these variations in size do not affect the result. In terms of Fourier analysis, the Lyot mask in B behaves as a high-pass filter and the Lyot stop in C behaves as a low-pass filter. It is the conjugate effect of these two masks that makes the Lyot coronagraph efficient for the rejection of the direct sunlight where we want to observe the corona. We name P(r), M(r) and L(r) the radial transmission functions of the entrance pupil in A, the Lyot mask in B and the Lyot stop in C respectively. The wave front Ψ B in the focal plane of the objective L1 results of a Fourier transformation, as described in Sec. 2.2.

Ψ B (r, θ) = ϕ f (r) iλ f × Ψ B (r, θ), (6)
where

Ψ B (r, θ) = F λ f [Ψ A (r, θ) × P(r)].
Therefore, the intensity in plane B is merely proportional to the Fourier transformation of Ψ A (r, θ) × P(r). The wave front in plane B encounters the Lyot mask M(r) and the second objective L2 of focal f 2 . Plane C is the image of the entrance aperture that is located at a distance

d = f × f 2 /( f -f 2 ) from plane B,
Ψ C (r, θ) = ϕ d (r) iλd × F λd Ψ B (r, θ) × M(r) × ϕ -f 2 (r) × ϕ d (r) = -ϕ d (r) λ 2 f d × Ψ C (r, θ), (7)
where

Ψ C (r, θ) = F λd Ψ B (r, θ) × M(r)
. Again, the intensity in plane C is proportional to the Fourier transformation of Ψ B (r, θ)×M(r) (Aime et al. 2002). Based on the same principle, the wave front in plane D is obtained by performing a Fourier transformation of Ψ C (r, θ) × L(r), corresponding to the image on the focal plane of the whole imaging system.

Ψ D (r, θ) = ϕ f (r) iλ f × F λ f Ψ C (r, θ) × L(r) . (8)
The wave propagation for the reference telescope S Ø is limited to the first propagation to plane B as described in Eq. ( 6). We also applied this analytic formulation to the external solar coronagraph S E by considering Fresnel diffraction for the complex amplitude Ψ A at the entrance aperture, given by Eq. (3).

Hybrid externally occulted Lyot solar coronagraph

Let us now consider the hybrid coronagraph S EL . This system varies from the classical Lyot coronagraph S L , since it has the internal occulter in plane O'. Plane B has no more actual interest in this particular case and shall be skipped. Moreover, Ψ A consists now of Fresnel diffracted wave front as given in Eq. ( 3), because of the external occulter in plane O. We directly write the Fourier-Fresnel propagation over the distance z 1 = z 0 f /(z 0f ) between planes A and O'. The wave front Ψ O in plane O' is then expressed as

Ψ O (r, θ) = ϕ z 1 (r) iλz 1 × F λz 1 Ψ A (r, θ) × P(r) × ϕ -f (r) × ϕ z 1 (r) Ψ O (r, θ) = ϕ z 1 (r) iλz 1 × Ψ O (r, θ), (9)
where

Ψ O (r, θ) = F λz 1 Ψ A (r, θ) × P(r) × ϕ -z 0 (r) .
The main difference between Ψ B in Eq. ( 6) and Ψ O in Eq. ( 9) is the quadratic phase factor ϕ -z 0 (r) = exp(-iπr 2 /λz 0 ). It can be interpreted as a virtual converging lens of focal length z 0 which rejects the external occulter at infinity. Consequently, the image of plane O made by the primary objective is now moved into the focal plane, and so it is computed as a simple Fourier transformation, as previously. This reasoning makes the computations much more convenient than first considering the wave in plane B and then propagating it to O' using Fresnel propagation over the distance f 2 /(z 0f ).

Then, the wave front in plane O' encounters the internal occulter M(r) and the second objective L2. Here, we can directly apply Eq. ( 7) to derive the complex amplitude of the wave front in plane C, where the distance d becomes now to d = z 1 × f 2 /(z 1f 2 ). However, the quadratic phase factor ϕ -z 0 (r) remains. Since we want to obtain in plane C the exact image of the pupil, we have to get rid of this unwanted factor. This is simply obtained by multiplying the complex amplitude in plane C by ϕ +z 0 (r), which corresponds to a diverging lens of focal z 0 that compensates the first virtual converging lens.

Ψ C (r, θ) = -ϕ d (r) λ 2 z 1 d × Ψ C (r, θ) × ϕ +z 0 (r), (10)
where

Ψ C (r, θ) = F λd 2 Ψ O (r, θ) × M(r)
. Finally, the wave front in plane D is given by Equation (8).

Observed intensities

The total intensity I K on plane K ∈ {A, B, O , C, D} is the incoherent sum of the elementary intensities due to every points source describing the whole solar disc. From here, the complex amplitude Ψ K,α,β (r, θ) will be written as a function of four variables Ψ K (α, β, r, θ), to clarify the integration process. Taking into account the center-to-limb darkening function B(α, β), the integrated intensity due to the whole solar disc is

I K (r, θ) = B(α, β) × |Ψ K (α, β, r, θ)| 2 dαdβ, (11)
where K ∈ {A, B, O , C, D}. This Fredholm integral of the first kind cannot be computed as a mere convolution since Ψ K (α, β, r, θ) is not shift invariant with respect to (α, β). So a 2D numerical summation must be performed, and the integral in Eq. ( 11) shall be transformed into a finite sum, using discrete values α k and β l for α and β. We emphasize that for each α k and β l we obtain a elementary 2D image. The required number of sampling points on the solar disc, that is, k and l, can be derived using Shannon criteria of the interpolation formula. The sampling must be tighter than 0.5λ/D p radian, with D p = 2R p the diameter of the entrance pupil. This corresponds to an upper limit of 1.13arcsec in our case, being a minimum of 1692 points in a solar diameter, or about 2.25 millions points on the whole solar disc. For the general Fredholm integral, the derivation of the required number of samples is not so straightforward, but the result is the same due to the inherent nature of band limited images. The computation of the final observed intensity in plane D therefore requires three times as much 2D Fast Fourier Transformations as the number of sampling points on the solar disc -6.75 million.

At this stage, we can actually take advantage of the axissymmetry of the system. We remind the reader that it assumes that the Sun and every remarkable planes are parallel and aligned to the optical axis. In Eq. ( 11), we now replace solar angular coordinates (α, β) by (ρ, θ s ), with α = ρ cos θ s and β = ρ sin θ s , meaning having ρ ∈ [0, R ] and θ s ∈ [0, 2π]. The assumed symmetry makes the 2D image of one point source |Ψ K (ρ, θ s , r, θ)| 2 in plane K rotate identically with respect to the point source on the solar disc, that is, θ s . In other words, it only depends on the relative angular difference φ = θ-θ s . As a result, integrating over the solar polar angle θ s is equivalent to circularly integrating on the 2D image plane, so over θ. Moreover, the solar brightness is a radial function, so B(α, β) = B(ρ). By substituting θ by φ, the integrated intensity given in Eq. ( 11) becomes the following radial function

I K (r) = 2π 0 R 0 B(ρ) × |Ψ K (ρ, θ s , r, φ + θ s ) | 2 ρdρ dφ. (12)
From a numerical point of view, it is much more convenient to compute this last integral than the rough full two-dimension summation. In Eq. ( 12), θ s can arbitrarily be fixed to 0, since the integration is performed over 2π, so we choose α = ρ and β = 0. This means that we only need to propagate the wave fronts coming from the point sources of one elementary radius of the Sun. The integration is thus performed as a weighted numerical summation of two-dimension images, using discrete values ρ k to sample the solar radius, using the sampling requirement defined above, followed by a circular average of the result.

Analysis and discussion

Numerical implementation

We now present the results of the complete computation of the observed intensities I K (r) on each plane K ∈ {A, B, O , C, D}, accordingly for each configuration S Ø , S L , S E and S EL . Using the complex amplitude Ψ A computed with Mathematica, the wave front is linearly interpolated and the propagation is performed with Matlab 2D Fast Fourier transformation combined with the re-centering routine fftshift when necessary. In addition to this paper, we provide a dedicated Matlab/Octave toolbox, whose content will be continuously updated, for the sake of reproducible research. This package can be found at the following link: https://github.com/rrougeot/FourierOptics.

A difficult and sensitive parameter in this numerical study is the choice of the sampling in each plane. Indeed, as discussed in Soummer et al. (2007), in successive planes (A to B or A to O', B to C or O' to C, and C to D), the sampling requirements are opposite. This problem is known as the two-fold sampling requirement. The point of view which has been adopted in the present work is somewhat empirical. We imposed the same number of points in the occulter image in plane O' and in the telescope aperture in plane A. We note that this a priori is sensible since it provide a similar resolution in all planes, but other sampling strategies might be of interest. The telescope aperture is padded inside an array of N × N points, and n p < N points are used in the radius R p of the aperture. The spatial sampling in plane A is thus s A = R p /n p , for a total spatial field F A = N × R p /n p . Due to the properties of Fourier spatial frequencies, the field F A in plane A produces a sampling s O = λd/F A in plane O', where d = z 0 f /(z 0f ) is the distance between plane A and plane O'. Moreover, the size of the image of the external occulter radius is R × d/z 0 . Therefore, the corresponding number of points n r in plane O' is n r = Rd/z 0 s O . By imposing n r = n p , we obtain:

n r = R × R p × N z 0 × λ ∼ 14.9525 √ N. (13)
For numerical reasons, N should preferably be a power of two, and at least 4096 points are required for a correct sampling of the image to respect Shannon criterion. Of course, the larger N the better the result due to zero-padding effect, being a compromise between computation time and precision.

The results of computation reported in this paper have been made using a machine with two 14 core Intel Xeon processors and 512GB of RAM, using N = 2 13 = 8192, giving n p = 1353, which corresponds to a sampling of 18.5µm and a spatial field of 15.6cm in plane A. We chose to over-sample the solar radius by setting 1000 point sources rather than 846. Each step in the numerical computation has been verified in particular with point sources. We used a Lyot mask or internal occulter of 1.065R radius and a Lyot stop sizing 0.99 times the entrance aperture image, to illustrate the numerical study. Despite the fact that the size of the Lyot stop is only scaled to the size of the entrance aperture, we warn the reader that the results strongly depend on the value of R p .

Impulse response in the Lyot coronagraph

We investigate the impulse response originating from one point source in both plane B and O', while adding the external occulter in plane O. We first analyzed the on-axis point source. In this case, the response is symmetric, and the image is a bright circle which perfectly fits the image of the external occulter edge. In plane O', this circle is very thin, as illustrated in Fig. 3. The intensity is focused as if the edge of the external occulter emits light as a real object. The response in plane B consists in a larger blurred circle, as expected. Second, we analyzed the response from the off-axis point source at α = 768arcsec, i.e. 0.8R .

Response of the different coronagraphic systems

Intensity in plane A

The rejection from the external occulter can first be assessed by computing the intensity I A (r) at the entrance aperture of the telescope. In Fig. 5, we plot the penumbra profile of diffracted light and its corresponding geometrical umbra profile. The horizontal axis represents the radial coordinate in mm, starting at the center of the umbra cone. The occulter disc is 1.0568 larger than the solar stenope image, which corresponds to a geometrical umbra of Rz 0 tan R = 38mm radius, and the intensity is equal to the full solar irradiance beyond R + z 0 tan R = 1382mm. Because of diffraction, the scattered light remains at a level of 10 -4 B at the center of the umbra cone, as a flat plateau. The external occultation has therefore reduced direct sunlight by four orders of magnitude at the entrance aperture, which is a first significant advantage for both externally occulted systems S E and S EL . Here, apodization techniques (Aime 2013) or more complex shapes of occulter [START_REF] Bout | [END_REF] may improve the performance. In Appendix A, variations of the distance z 0 are investigated, and we provide a plot illustrating different penumbra profile.

Intensity in plane B

Figure 6 shows the radial intensities I B (r) in plane B, limited to 3.2R , in logarithmic scale. Here, the image of the Sun (blue curve) is perfectly focused, and is used as a reference for normalization. This consists of the global response of the raw telescope S Ø . The center-to-limb variation is apparent as a slight decrease in the range 0 -1R . Sunlight falls abruptly to 10 -3 B at 1R , then extends as a large tail of residual light brighter than 10 -5 B . This comes from the summation of the Airy rings at large radius. The diffracted light pattern produced by the external occulter consists of a bell-like curve out-of-focus, as expected, since the focal plane is not the conjugate image plane of the external occulter. The width of this peak is function of the size of the entrance aperture, like the Airy radius. We note that the peak is not symmetric, and reaches a maximum of 10 -3 B around 1.05R . This last curve models the response of the externally occulted solar coronagraph S E .

Intensity in plane O'

Similarly, Figure 7 shows the radial intensities I O (r) in plane O' in logarithmic scale, using the same scaled axis as Figure 6 for a purpose of comparison. We remind the reader that using solar units here has no real meaning, since plane O' is not conjugated with the Sun as discussed in Sec. 2.1. The image of the Sun is very similar to the one in the focal plane -Fig. 6, but is slightly out-of-focus here, as it is at df = 0.758mm ahead. The drop to 10 -3 B is consequently smoother. The large tail of residual light is still present. The diffracted light by the external occulter is now focused in a very narrow peak of 10 -2 B amplitude, located at the exact angular position of the edge of the external occulter image, i.e. 1.0568R or 1.629mm. This feature is expected because plane O' is the conjugate image plane of the external occulter.

Residual light

Before inspecting the observed intensity in plane C and plane D, it is interesting to look at the integrated residual light on both planes B and O'. In the classical Lyot coronagraph S L , the Lyot mask is set in the focal plane and blocks the direct focused sunlight (Fig. 6). A relatively large amount of residual light yet propagates further inside the instrument. In the case of the hybrid externally occulted Lyot coronagraph S EL , the internal occulter blocks the diffracted light fringe (Fig. 7), and its dimension governs the rejection. For comparison, we looked at the integrated residual light denoted as

L(r) = 2π 0 dθ r r I K (r)rdr, ( 14
)
where K ∈ B, O , and with the numerical upper limit r = 3.2R . We analyzed the three following cases: residual light that is not blocked for a given Lyot mask of radius r in the classical Lyot coronagraph S L . When adding the external occulter, the residual light is already reduced by a factor 10 -4 ∼ I A (0), and there is an appreciable difference between planes B and O'. The integrated residual light decreases relatively slowly in plane B, from 10 -4 B at r = 1.02R and loosing two orders of magnitude over 0.1R . In plane O', the decrease is very abrupt, from 10 -4 B at r = 1.06R and loosing two orders of magnitude over 0.06R . As a conclusion, this shows that the internal occulter will filter out a larger amount of diffracted sunlight by being set in plane O' rather than plane B, for a given size. We provide a definitive confirmation of this last statement in Appendix B, by looking at the final response in plane D.

L Ø,B (

Intensities in plane C

In Fig. 9 we present the radial intensities I C (r) in plane C in logarithmic scale. The transverse radius is given in units of entrance aperture image, meaning that a radius of one corresponds to the image of R p = 25mm. A 1.065R Lyot mask and an internal occulter of equivalent angular size have been used here, for the systems S L and S EL respectively. We normalized the intensities using the exact image of the entrance pupil. The classical Lyot coronagraph S L shows a narrow peak at r = 1, being the exact position of the image of the pupil edge. This 10 -1 B fringe is produced by the diffraction of sunlight by the entrance aperture. A similar feature is observed in the case of the hybrid externally occulted Lyot coronagraph S EL . The diffraction fringe is however much less bright, reaching about 10 -4 B , due to the external occultation beforehand. The role of the Lyot stop is to block this diffracted light peak.

Intensities in plane D

Finally, in Fig. 10 we give the final response of the four imaging systems, in logarithmic scale, using the same occulting masks of 1.065R in planes B and O', and a Lyot stop of 0.99 times the image of the pupil -24.75mm in 1:1-scale. We superimposed the four observed intensities, meaning I B (r) for the reference telescope S Ø and for the external coronagraph S E , and I D (r) for the classical Lyot coronagraph S L and for the hybrid coronagraphic system S EL . The reference image of the Sun in plane B is used as a reference for normalization. Both systems S L and S E show a relatively bright (10 -3 B ) diffraction fringe located around 1.065R , while the hybrid coronagraphic system S EL already rejects sunlight below 10 -5 B . Outside 1.5R , the three systems S L , S E and S EL reject below 10 -6 B , 10 -7 B and 10 -8 B respectively. As for now, this analysis has proved the efficiency of combining external occultation with an internally occulted Lyot-style coronagraph, compared to the classical Lyot coronagraph, with a gain of at least two orders of magnitude.

Sizing the internal occulter and the Lyot stop

We now consider the hybrid externally occulted Lyot coronagraph S EL only. We investigated the impact of sizing the internal occulter on the observed intensity in plane D, while keeping the external occultation ratio R/z 0 constant. We looked at radii of 1.005, 1.01, 1.02, 1.03 and 1.04 times the external occulter image. They respectively correspond, in angular units (metric units), to 1.0621R (1.637mm), 1.0674R (1.645mm), 1.0779R (1.662mm), 1.0885R (1.678mm) and 1.0991R (1.694mm). We superimpose onto Fig. 11 the radial cuts of intensities I D (r), in logarithmic scale, using a Lyot stop of 0.99 (plot (a)) and one of 0.96 (plot (b)) entrance pupil image. The plot is given in the range 0.5 -2R to zoom in the diffraction fringe area.

In a similar way, we analyzed the effect of sizing the Lyot stop, keeping a fixed internal occulter. In Fig. 12, we compare the radial cuts of final intensities I D (r), in logarithmic scale, using an internal occulter of 1.01 (plot (a)) and one of 1.03 external occulter image (plot (b)). We investigated the following sizes of Lyot stop: 1.00, 0.99, 0.98, 0.96, and 0.92 entrance pupil image. The interested reader will find in Appendix C the same study for the classical Lyot coronagraph.

These plots give an intuition of the behaviour of the result. Indeed, we show here that the Lyot stop mainly acts over the diffracted light in the range 1.2 -3R , and does not significantly impact the main diffraction feature. As shown in Fig. 12, reducing the radius of the Lyot stop from 1.00 to 0.99 already gives appreciable improvement on the rejection, of about one order of magnitude. In parallel, increasing the size of the internal occulter mainly contributes to reducing the level of residual sunlight around the edge of the external occulter image -Fig. 11. We observe an improvement of two orders of magnitude from 1.005 to 1.02 external occulter image. However, it seems that the performance in rejection in the range 1.5 -3.2R is more impacted by the size of the Lyot stop than the internal occulter. It is also interesting to observe that the position of diffraction peak matches the position of the image of the internal occulter in plane D. So, when increasing the size of the occulting mask, the peak of diffraction translates accordingly. As a result, the residual diffracted sunlight keeps contaminating the inner region of the field of view.

It is of course the combined effect of the internal occulter and the Lyot stop that makes the performance of the coronagraph. Figure 13 illustrates this point as we have plotted the residual diffracted sunlight level observed at 1.3R , versus the radius of the internal occulter and versus the radius of the Lyot stop. At that stage, one can already get an idea of possible theoretical performance (for an ideal and perfect instrument), of such a hybrid externally occulted Lyot solar coronagraph -in this particular configuration. At least, rejecting diffracted sunlight below 10 -8 B at 1.3R using a 1.065R internal occulter looks feasible, but this shall be considered as a theoretical lower limit.

Analysis of the vignetting

Finally, we discuss the vignetting induced by both external and internal occulters since it is a characteristic feature in coronagraphy. A dedicated analytic study of issues with external occultation can be found in Raja [START_REF] Raja Bayanna | [END_REF]. This vignetting affects the transition region where the coronograph removes the direct light from the solar disc and transmits that of the solar corona. In our present study, we take advantage of our model of light wave propagation to estimate the vignetting coming from off-axis point sources outside the solar disc, that is, ρ = α 2 + β 2 > R . Again, the complex amplitude Ψ A,α,β incoming into the entrance aperture is given by Eq.(3). It is interesting to understand the consequences of the external occultation. We remind the reader that the figure of the Arago bright spot in plane A, plotted in Fig. 2, is shifted of the quantity (z 0 × α, z 0 × β) towards negative (x, y) directions. The transition between the shadow region and the high-intensity plateau is located around R = 710mm from the centre of the spot. Denoting ρ 0% = (R -R p )/z 0 = 1.0196R = 979arcsec and ρ 100% = (R + R p )/z 0 = 1.094R = 1050arcsec, we can roughly say that ρ ≤ ρ 0% only faint diffracted light from the shadow region will enter the entrance pupil, corresponding to the occulted region of the sky. ρ 0% ≤ α ≤ ρ 100% the high-intensity plateau is partially captured by the entrance pupil. This region of the sky is partially vignetted. α ≥ ρ 100% The part of the wave front that enters into the entrance pupil is poorly affected by the external occulter and can be approximated by a planar tilted wave front. There is no vignetting.

The complete vignetting function of the coronagraph has been inspected by propagating the wave front from unitary amplitude point sources in the range ρ 0% ≤ α ≤ ρ 100% , and β = 0. In Fig. a 1.065R internal occulter and a 0.99 Lyot stop. Plot (a) represents the impulse responses when there is no external nor internal occulters. In that case, we observed perfect Airy patterns. Plot (b) illustrates the actual response of the coronagraphic system, i.e. with both occulters. We can see in Fig. 14 the strong nonlinear perturbation of the inner corona, in addition to the attenuation of the intensity. The image reconstruction in this region will be a delicate problem of inversion in a Fredholm regime. Beyond ρ 100% , we find again a perfect Airy pattern, as expected.

Conclusions

We have presented a dedicated analytic and numerical analysis of the theoretical rejection performance of the classical Lyot coronagraph, the externally occulted solar coronagraph, and the hybrid externally occulted Lyot coronagraph. We first computed Fresnel diffraction produced by the external occulter, and second the coherent propagation of the wave fronts through the instrument. Here, our results applies to the geometry of ASPIICS coronagraph. We provide the observed intensity of residual sunlight in the final focal plane, that may contaminate the observation of the solar corona in white light, especially close to the solar limb. Using a Lyot mask of 1.065R radius and a Lyot stop being 0.99 times the image of the entrance aperture radius, we showed that the perfect classical Lyot coronagraph manages to reject below 10 -6 B from 1.3R . The externally occulted solar coronagraph provides a better performance, with a gain of one order of magnitude outside the diffraction peak intensity. Finally, the hybrid externally occulted Lyot coronagraph improves the global performance by rejecting diffracted sunlight below 10 -8 B from 1.3R . We also refined our study to exhibit the coupled effects of sizing both internal occulter and Lyot stop, in the case of the hybrid coronagraphic system. Oversizing the mask allows us to decrease the intensity of the diffraction peak intensity located around the image of the external occulter, and reducing the radius of the Lyot stop allows to globally reduce the residual sunlight. As a concrete result, we have provided in Fig. 13 a graph estimating the rejection at 1.3R , as a function of both sizes of the internal occulter and Lyot stop. Rather than investigating whether the hybrid externally occulted Lyot solar coronagraph could meet any requirements in stray light rejection, our work claims to be a methodological model in order to estimate the end-to-end performance of such instruments. Solar astronomy will benefit from it to identify any benchmark for ongoing or future activities. However, we would like to emphasize that the given results remain a lower bound. In practice, performance would be degraded by any other sources of stray light and scattering.

As future activities, one can investigate in the nature of the external occulting disc itself. Radially apodized external occulters have already been analytically proved by Aime (2013) to be more efficient than the sharp-edged disc, in the context of solar coronagraphy. In an experimental approach, [START_REF] Bout | [END_REF] and Landini et al. (2010) investigated 3D-shaped occulters which deviate from a simple radial apodization. Model validation against experimental results would be very interesting and instructive. However, the main difficulty stands in the two-dimensional representation of such complex external occulter shapes and computing the Fresnel diffraction, as it has been recently investigated by [START_REF] Sirbu | [END_REF], in addition to the extend source that is represented by the Sun.

Fig. 1 .

 1 The light encounters these successive planes in the order O, A, B, O', C and D. As described previously, our work is focused on three different coronagraphic systems plus the related reference imaging system. Figure 1 presents a schematic illustration of the four systems named as follows. S Ø denotes the raw telescope used as a reference, consisting of plane A and plane B -Fig. 1a. S L is the classical Lyot coronagraph, consisting of plane A, the Lyot mask in plane B, and the Lyot stop in plane C -Fig. 1a. S E is the externally occulted solar coronagraph, consisting of the external occulter in plane O, and ending at the focal plane B -Fig. 1b. S EL denotes the hybrid externally occulted Lyot solar coronagraph composed of the external occulter in plane O, an internal occulter and the second objective (L2) in plane O', and the Lyot stop in plane C -Fig. 1b

Fig. 1 .

 1 Fig. 1. Schematic representation of the four imaging systems. (a) The classical Lyot coronagraph S L made of planes A, B, C and D. The raw telescope S Ø ends at focal plane B. (b) The hybrid externally occulted Lyot solar coronagraph S EL made of planes O, A, O', C and D. The externally occulted solar coronagraph S E ends at focal plane B. Figures not to scale.

R.Fig. 2 .

 2 Fig. 2. Fresnel diffraction pattern Ψ A,0,0 (r) 2 for unitary on-axis point source at infinity, for the 710mm radius occulting disc at z 0 = 144.348m. (a) Full range, in logarithmic scale. (b) Zoom in the Arago bright spot in the central region [0, 500µm], in linear scale. (c) Transition region between shadow and light around 710mm, in linear scale.

  as given by the relation of conjugation for lens. Writing a Fourier-Fresnel transformation of Ψ B over the distance d, the three quadratic phase factors ϕ f (r), ϕ -f 2 (r), and ϕ d (r) cancel each other out.

  Figure 4 shows the two-dimension intensities in both plane B (plots (a) and (b)) and O' (plots (c) and (d)). The shift of the Arago spot produces strongly asymmetric light pattern, while the onaxis case was perfectly symmetric. The sharpness and the fine scale structure of the diffraction features, as shown in plots (b) and (d) in Fig. 4, for plane B and O' respectively, prove the need of a very high sampling on both planes. A very dominant point is that the light pattern in plane B tends to be spread perpendicularly to the local edge of the image of the external occulter. On the contrary, plane O' shows a diffracted light pattern fitting locally the image of the external occulter.

Fig. 3 .

 3 Fig. 3. Two-dimension observed intensity in plane O' of diffracted wave front originating from the on-axis point source at the center of the Sun, in logarithmic scale. (a) Full field. (b) Zoomed in region of interest around 1.08R .

Fig. 4 .

 4 Fig. 4. Two-dimension observed intensities of diffracted wave front originating from the off-axis point source at α = 768arcsec, in logarithmic scale. (a) I B (x, y) in plane B. (b) I B (x, y) enlarged in region of interest around 1.08R . (c) I O (x, y) in plane O'. (d) I O (x, y) zoomed in region of interest around 1.08R .

Fig. 5 .Fig. 6 .

 56 Fig. 5. Intensity I A (r) in plane A, in logarithmic scale. The intensity is normalized to the mean solar brightness. Red: radial profile of diffracted light produced by the external occulter. Black dashed: related geometrical profile. (a) Penumbra profile in full range. (b) Enlargement in [0mm; 100mm] corresponding to the umbra region.

R.Fig. 9 .Fig. 10 .

 910 Fig. 9. Observed intensities I C (r) in plane C, in logarithmic scale. The transverse radius is given in unit of image of the entrance pupil. The intensities are normalized to the mean solar brightness. Blue: system S L . Red: system S EL . Vertical line: Lyot stop of 0.99 times the image of the entrance pupil.

Fig. 11 .Fig. 12 .Fig. 13 .Fig. 14 .

 11121314 Fig. 11. Radial cuts I D (r) for internal occulters of 1.005, 1.01, 1.02, 1.03 and 1.04 external occulter image, in logarithmic scale. (a) Fixed Lyot stop of 0.99 entrance pupil image. (b) Fixed Lyot stop of 0.96 entrance pupil image. The transverse radius is given in solar units, in the range [0.5R , 2R ]. The intensities are normalized to the mean solar brightness.

Table 1 .

 1 Key planes of the coronagraphic systems.

	Plane Description
	O	External occulter plane
	A	Aperture of the telescope
	B	Focal plane of the telescope
	O'	Image of plane O made by the telescope
	C	Image of plane A
	D	Final focal plane

Table 2 .

 2 Definition of the four imaging systems.

	Name of the system	Combination of planes
	S Ø	Raw telescope	A + B
	S L	Classical Lyot coronagraph A + B + C + D
	S E	External solar coronagraph O + A + B
	S

EL Hybrid externally occulted O + A + O + C + D Lyot solar coronagraph aperture) and D (final focal plane)

. By adding an external occulter, two additional planes must be introduced. On one hand, the external occulter is positioned in plane O at a finite distance from the entrance aperture of the telescope. On the other hand, plane O' denotes the image of the external occulter made by the telescope. It is located further behind the focal plane. In our model, we assume that the primary objective (L1) coincides with the pupil in plane A. Table

1

recalls the names and descriptions of all the planes, also illustrated in

Table 3 .

 3 Parameters for the numerical study -from ASPIICS configuration. See text for details.

	Parameter	Value
	Wavelength	λ = 550nm
	Angular radius of the Sun	R = 0.0046542rad
		= 960arcsec
	Distance to the Sun	∞ (1 Astronomical Unit)
	Radius of the external occulter R = 710mm
	Distance plane O -plane A	z 0 = 144.348m
	Radius of the pupil	R p = 25mm
	Focal length of the telescope	f = 330.385mm
	of the propagated wave front at distance z is

  14 we plot the two-dimension observed intensities in plane D for the externally hybrid external Lyot coronagraph S EL , using Observed intensities I O (r) in plane O', in logarithmic scale. The transverse radius is given in solar units. The intensities are normalized to the mean solar brightness. Blue: raw image of the Sun in plane O'. Red: diffracted sunlight by the external occulter. Vertical line: 1.065R radius internal occulter.

		0					S Ø , SL SE, SEL
		-2
	Intensity (log)	-6 -4
		-8
		-10
		0	0.5	1	1.5	2	2.5	3
					Radius (R ⊙ )
	-6 -4 -2 0 Fig. 7. 0.8 Residual intensity (log)	0.9	1	1.1	1.2	1.3	L Ø,B L E,B L E,O'	1.4
					Radius (R ⊙ )

Fig. 8. Integrated residual light L(r), as given in Eq.(14), in logarithmic scale, and normalized. Enlargement in the range [0.8R , 1.4R ]. Black: L Ø,B . Blue: L E,B . Red: L E,O' .

Acknowledgements We would like to thank the European Space Agency and Proba-3 project for having supported this activity.

The shape of the umbra cone is intimately linked to the distance z 0 between the external occulter and the telescope. In the context of Proba-3 Formation Flying mission, this is of particular interest, since the inter spacecraft distance separated the coronagraph and the occulter may vary. We computed the Fresnel diffraction pattern (Eq.( 5)) and the penumbra profile I A (r) in plane A (Eq.( 12)) for different values of z 0 , while keeping a constant radius R for the external occulter. In Fig.

A.1, we give the radial penumbra profile of diffracted sunlight in logarithmic scale at z 0 = 1m, 10m, 50m, 100m, 144.348m and 200m. First of all, we must state that Fresnel diffraction theory usually requires small angles approximation, which may not be the case at small z 0 . The smaller z 0 , the smaller the stenope image of the Sun in plane O, computed as z 0 tan R , as given in Table A.1. From z 0 ≤ 50m, we observe that the penumbra is bell-shaped in the central region. This comes from the two-dimension convolution of this stenope image with the Arago bright spot in plane A. At larger z 0 , this feature vanishes and we obtain a smooth penumbra profile. We notice that the umbra is about 2, 5 darker at z 0 = 100m than 144.348m, as plane A is then closer to the occulter. The geometrical umbra is reduced to a point at z 0 = 152.55m, and the so-called ante-umbra region extends behind, where the external occulter cannot mask the whole solar disc any more. This last case is illustrated by Of course, setting the internal occulter in plane O' instead of plane B is much more efficient about stray light rejection, as already discussed in Sec. 4.3.4. We observe an average gain in performance of one order of magnitude. The diffraction peak intensity around 1.06R is even two orders of magnitude lower.

Appendix C: Reducing the Lyot stop for the classical Lyot coronagraph

We report in this appendix the impact of reducing the Lyot stop for the classical Lyot coronagraph S L . The approach is the same as the one presented in Sec. 4.4. We computed the global response in the final focal plane by varying the radius of the Lyot stop, from 1.00 to 0.96 times the image of the entrance pupil, and using a fixed 1.065R

Appendix D: Two-dimension intensities

In this appendix, we report the two-dimension image of the global response for the four different imaging systems studied. In