© C. Theys et al. 2022

1 Introduction

In this paper, we study the reconstruction of an object based on an observation that does not result from a relationship of mere convolution between the object and the response of the instrument but, rather, requires a complete Fredholm integral of the first kind for its description, Groetsch (2007). In the present application, the object is the solar corona and the instrument is the coronagraph ASPIICS, using an external occulter in front of a Lyot (1933) coronagraph, as per the configuration first proposed by Evans (1948).

The external occulter creates an artificial eclipse for the coronagraph set in its shadow and allows for a better reduction of the straylight, as analyzed by Newkirk & Bohlin (1963) and Koomen et al. (1975), for example. The drawback of such system is the difficulty in carrying out an observation very close to the solar limb due to vignetting. The aperture of the telescope becomes progressively reduced to a crescent shape by the umbra of the external occulter located in front of it. This was well illustrated in Llebaria et al. (2004, 2006) for the LASCO-C2 coronagraph of the SOHO mission, Brueckner et al. (1995). In spite of these vignetting effects, which do not exist in a pure Lyot (1933) coronagraph, the strong reduction of the stray light prevails; and an Evans configuration has been adopted in almost all space experiments, as noted by Koutchmy (1988).

As a first approximation, the loss of the inner observing angle close to the limb is the ratio of the diameter of the telescope to the diameter of the occulter, in units of solar radius. A small vignetting zone requires a very large occulter that is set far away from the telescope. Such a result cannot be obtained with a single spacecraft and the formation flying project ESA’s Proba-3 ASPIICS is an excellent way to solve this problem, Lamy et al. (2010); Renotte et al. (2014); Galano et al. (2018).

The first satellite holds a 1.42-m diameter occulting disk and the second satellite, positioned at 144-m away in its umbra, carries a 5-cm Lyot coronagraph, where the final observation is made, Galano et al. (2019). ASPIICS aims to observe the solar corona very close to the solar limb, as has never been done before beyond eclipse conditions. The first proposals for ASPIICS (Lamy et al. 2007, 2010, 2017) carried the ambition of an observation of the solar corona very close to the limb, as low as 1.01 R_⊙ (solar radius), and even the goal of including chromospheric and prominence physics. The current project, for reasons of reliability in the event of positioning errors, will place the telescope pupil in the shadow cone at a position where the shadow disk is significantly larger than the pupil. The drawback is that this increases the blind area near the solar limb. Recent studies have envisaged observations starting at no sooner than 1.05 R_⊙ or even 1.1 R_⊙, depending on the experimental parameters used, Shestov et al. (2021) and Rougeot (2020).

In the region of vignetting, the point spread function (PSF) degrades rapidly before disappearing when the shadow of the occulter covers the entire entrance pupil. While the main effect is geometrical, diffraction makes the response more complex and the frontier between light and darkness less well defined. The spatial variation of the PSF in the region of vignetting is very fast, and solutions considered in the field of adaptive optics, such as in La Camera et al. (2015), Denis et al. (2015), Thiébaut et al. (2016), cannot apply to the present work.

The present paper follows and complements the study of Aime et al. (2019), based on the assumption of a Richardson-Lucy type algorithm, (Richardson 1972; Lucy 1974), applied and adapted to the case of a non-constant PSF. This keeps the complete matrix formalism and upgrades the analysis with regard to several fundamental aspects. A realistic solar object is used instead of simple point sources. The stray light present in the observation is taken into account, and not only the direct Fredholm transformation as in the former study. This stray light increases since the internal occulter used is smaller and a compromise between resolution near the solar limb and reduced stray light is to be considered.

This work has a signal processing aspect and a physical aspect that are considered in that order in this paper. In Sect. 2, we describe the Fredholm inversion in the case of a space variant PSF. The matrix formalism and the reconstruction algorithm are given there. Different statistical models adapted to the raw data (Poisson model) and to the de-biased data (Gauss-Poisson and Gaussian models) are presented in Sect. 3. The PSFs and biases calculated for ASPIICS with two values of the internal occulter are given in Sect. 4. Numerical results performed on an eclipse image are shown in Sect. 5. A robustness study in the case of instrumental defects is presented in Sect. 6 for formation flying errors that are present but ignored in the inversion procedure and for the surface roughness of the optics. Prospects for the extension of the matrix formalism to other instruments are proposed in Sect 7. Our conclusions are given in Sect. 8. Peculiarities of the PSFs in the region dominated by diffraction are detailed in the Appendix.

2 Fredholm inversion principle

As stated in the introduction, the objective of this work is to reconstruct the object X(α, ß) in the vignetting area of an externally occulted coronagraph. The intensity Y(α, ß) that is formed in the final observing plane of the coronagraph can be written as: (1)

where h(σ,ξ,β,η) is the space-variant PSF in the region of vignetting, the integral a two-dimensional (2D) Fredholm integral of the first kind and B(α, ß) is the stray light that appears here as a bias. This equation corresponds to Eq. (3) of Aime et al. (2019), with the addition of stray light. We note that dummy variables (ξ, η) and coordinates (α, ß) are associated with the object and image planes, respectively. In the following, we give the elements of the matrix approach to discretize Eq. (1).

2.1 Matrix formalism

We note that in the case of an invariant PSF, the Fredholm relation in Eq. (1) can be simplified to a convolution relation, allowing for a fast resolution in Fourier space via FFTs. This is no longer the case for the variant PSF considered here and matrix manipulations must be carried out. The discretized form of Eq. (1) can be written as: (2)

where y_mn, x_i,j and b_m,n are the M × N discretized values of Y(α, ß), Χ(ξ, η), and B(ξ, η), respectively, while h_m,i,n,j are the M² × N² discretized values of h(α, ξ,β, η).

To use a matrix formalism, we arrange the arrays in vectors. The order chosen is lexicographic, that are the values of each array of dimension M × N (read from left to right and from top to bottom) are arranged in this order in a vector. An image of size M × N then yields a vector y of length M × N: (3)

In doing so, Eq. (2) can be written in the following matrix form: (4)

where y, x and b are vectors of size M × N containing the discretized values of the image and the object. H is the matrix containing the elements h_m,i,n,j, of dimension MN × MN: (5)

Each column of H corresponds to one PSF ordered lexicographically. The first column corresponds to the PSF of the point source at the upper left corner of the object (i = j = 1), while the last column corresponds to the PSF of the point source at the bottom-right corner of the object (i = M, j = N). It should be noted that inevitable boundary problems exist in the discrete representation of the Fredholm integral, as they exist for the simple convolution relation.

2.2 Reconstruction algorithm

The estimation of x given the data, y, matrix, H, and bias, b, namely, Eq. (4), is a so-called ill-posed problem. In fact, the direct inversion of this equation provides unacceptable solution due to very bad conditioning of H. The classical method to solve such an inverse problem is to use an iterative algorithm to minimize a criterion which measure a discrepancy between the observation y and the linear model Hx. Iterations have to be stopped before instabilities in the solution appear unless an explicit regularization is used, (Demoment 1989; Titterington 1985). Since we are dealing with intensities, a non-negativity constraint has to be imposed on x. The problem is then to solve: (6)

with J a criterion depending on x, and, of course, on y as well. As is classically done, we decompose J into two terms: (7)

where J₁(x) and J₂(x) are, respectively, the data consistency and the regularization terms. The coefficient γ makes it possible to adjust the weight of J₁ with respect to J₂, (see Titterington 1985). Depending on the statistical model for the noise, the data consistency term J₁ (x) takes different forms, as we show in Sect. 3.

The regularization term J₂ is independent from the data consistency term, that is, it is a regularization term and must, in the context of an ill-posed problem, impose a constraint of smoothness on the solution x, Titterington (1985). The choice of J₂ is further discussed in Sect. 5.2, once the particular form of the ASPIICS PSF has been presented (in Sect. 4).

J(x) is assumed to be convex relative to x and, therefore, the simplest way to minimize it is to use a gradient descent algorithm with non-negativity constraint. The procedure is described in Lantéri et al. (2002) and summarized here.

The Karush Kuhn Tucker (KKT) conditions, (Karush 1939; Kuhn & Tucker 1951), once resolved, can be reduced to an unique equation that the optimal value of the solution, x*, must satisfy, that is, for a component r: (8)

This condition will then serve as a descent direction in the iterative algorithm: (9)

with , the descent step-size chosen to have at each iteration k. The search procedure is carried out in an interval computed at each k to ensure ≥ 0, starting from a positive initial value, making it an interior point algorithm, (see Armijo 1966; Lantéri et al. 2002).

Making use of Eq. (7), the gradient ∇_x J necessary for the implementation of the algorithm of Eq. (9), can be split into the two terms: (10)

The gradient of the data consistency term ∇_x J₁ depends on the given model and is given in Sect. 3, ∇_x J₂ will be given in the Sect. 5.2.

3 Data modeling

This section is dedicated to the development of statistical models for raw data and for debiased data further on. We show in what follows that the drawback of subtracting the bias is that it makes it more difficult to set up of the statistical model needed for the Fredholm inversion.

3.1 Raw data (raw) and the Kullback-Leibler divergence

We may consider that the discretized values of y_m,n in Eq. (2) correspond to a classical integrated intensity in the pixels of the detector. These intensities, converted into units of photoelectrons, undergo a photodetection process, which is well represented, for an ideal detector by a Poisson distribution P. We remark that in the case of this model, the y_raw data are non-negative.

The corresponding noisy data are denoted as y_raw: (11)

Under the assumption of Eq. (11), the raw data follow a Poisson statistics. Assuming moreover that all pixels of the image y_i are independent, the likelihood or conditional probability p(y_raw/x) may be written as: (12)

In order to achieve a minimization algorithm, we take the opposite of the log likelihood to obtain J₁, in Eq. (7): (13) (14)

and the gradient w.r.t x, to get ∇_x J₁, Eq. (10) is: (15)

where ./ indicates element-wise division.

This expression of the gradient of J₁ , ∇_x J₁ associated with the gradient of J₂, Eq. (29), will be used in the iterative reconstruction algorithm, Eq. (9), to process raw data.

We note that in the non-regularized case, namely, γ = 0 and with a descent step-size, the algorithm obtained from Eq. (9) in the case of a Poisson statistic is equivalent – although not strictly equal – to the so-called Richardson Lucy algorithm.

3.2 Debiased data (deb)

Since we assume a knowledge of the bias, that is a numerical computation of the straylight, it seems an a priori natural choice to subtract it as part of the pre-processing of the data. Let us denote the new bias-removed data, y_deb: (16)

where is the estimate of the straylight. Contrary to the Poisson case, the y_deb data can be negative, that is the case once the bias is removed from the data. We can approximate P(b) – with a Gaussian law of zero mean and mean , because of the large number of photons in this solar experiment. Two statistics are derived from there, presented below.

3.2.1 Gauss-Poisson model

We considered the data after removing the bias, Eq. (16). Due to the subtraction of an estimate of the bias, some pixels became negative and the observed quantity can no longer be written as a pure Poisson process. The model of Eq. (16) then becomes a Gauss-Poisson model, such that: (17)

with (0, diag(b)), where diag(b) is the diagonal covariance matrix with diagonal b.

This Gauss-Poisson model has been the subject of several astrophysical publications in e.g Lantéri & Theys (2005); Benvenuto et al. (2008). The probability function in this case is the convolution of the Poisson and Gauss probability laws for a pixel, as per Eq. (3) of Benvenuto et al. (2008): (18)

and the gradient for x is, for a pixel, i: (19)

with (20) (21)

and (22)

It is important to note that the terms to be calculated in the gradient expression above can lead to numerical instability for large values of n_i. This leads us to consider a simplified model.

3.2.2 Gaussian model

A simplified and more tractable model for debiased data is a full Gaussian model with a covariance matrix depending on the data: (23)

where diag(y_deb) is the covariance matrix with y_deb on the diagonal. Then, the conditional probability, assuming that all pixels of the image y_i are independent, is: (24)

In order to achieve the minimization algorithm, we take the opposite of the log likelihood: (25) (26)

and the gradient for x: (27)

In the same way as for the KL divergence, this expression of ∇_xJ₁, associated with the gradient of J₂, will be used in the iterative reconstruction algorithm, namely, Eq. (9).

4 Application to ASPIICS

In order to derive the values of the matrix H and the straylight b for ASPIICS, we made use of the numerical model developed by Rougeot et al. (2017). This approach was further used as a starting point for studies including optical defects of the coronagraph (Rougeot et al. 2017, 2019; Shestov & Zhukov 2018; Shestov et al. 2021).

The program, written in MATLAB, computes the propagation of a wave originating from a point source at the level of the Sun through the whole ASPIICS system until the final observing plane, with a spatial sampling of 0.7497 × 0.7497 arcsec. This model takes into account the Fresnel diffraction of the external occulter mounted on the Occulter Spacecraft and its propagation across the Lyot-style coronagraph carried by the Coronagraph Spacecraft. The numerical model makes it possible to compute the resulting PSF for each point source, inside the solar disc or in the corona, and by integrating the stray light attributed to the entire solar disc. In this section, we assume that the complete coronagraph system is perfect. The effects for an imperfect system is analyzed in Sect. 6.

In the present study, we perform the calculations for two values for the radius of the IO, corresponding to 1.0206 and 1.0759 times the size of the EO image (physical values of 1.662 mm and 1.752 mm inside the experiment). The larger the diameter of the IO, the better the straylight rejection is, with the drawback of preventing the observation of the corona very close to the solar limb. In the following, we refer to these IOs as the Small IO (SIO) and the Large IO (LIO). The inversions of the Fredholm transform are made on rectangular arrays of 128 (transverse) by 236 (radial) points (96 × 177) arcsec, both for the object and the evanescent PSFs. The radial dimension is so to cover the two different vignetting areas of LIO and SIO. The matrix H are of size 30 208 × 30 208, since an image of dimension N × M requires the use of a matrix of dimension (N × M)². The results for SIO given in the present paper have never been published before.

Fig. 1 Expected level of the corona (model of Cox 2015) and straylights for the small (SIO) and the large (LIO) internal occulter

4.1 Straylight, b

Stray lights corresponding to the SIO or LIO internal occulter for the perfect coronagraph are given in Fig. 1. Stray light curves exhibit the characteristic two-bump shape that is observed in all studies for ASPIICS. As expected, the stray light is lower for LIO than for SIO. The maximal value is 1.3210⁻⁶ at ρ = 1.0542 (in units of solar radius) for SIO, and 1.112710⁻⁷ at ρ = 1.1112 for LIO. We note that in a small region, the second bump of straylight of LIO overcomes the stray light of SIO.

To compare the stray light with respect to the corona, the expected theoretical level of the K-corona in white light is plotted in the same figure. It is drawn following the model of Cox (2015).

4.2 Evanescent PSF and the matrix H

The combined effect of the external and internal occulters produces a vignetting on the telescope aperture of the Lyot coronagraph. Since the internal occulter (IO) is placed on the image of the external occulter (EO) in the coronagraph and is larger than it, we can consider, as a first geometric approximation, that the IO defines the vignetting domain. In that region, the PSF evolves both in shape and integrated intensity, ultimately disappearing close to the limb. Figure 2 represents radial cuts of a few PSFs for the two cases SIO and LIO. For SIO, the PSF vanishes below roughly 1.04 solar radius, while the same phenomenon occurs already at 1.1 solar radius for LIO. In this representation, only a few PSF are drawn, for the sake of clarity. A 2D image of a PSF will be shown in Appendix A for comparison between perfect and imperfect optics.

Figure 3 gives integrated quantities for the small (SIO) and the large (LIO) internal occulter. The integrated flux correspond to the integral of the PSF, variable from zero to 1 in the vignetting zone, this integral being normalized to 1 outside this region. In the same graph pseudo Strehl ratios are given, also variable from zero to 1 in the vignetting zone. This Strehl ratio was already drawn in Fig. 2 as the maximal value of the PSF, being normalized to 1 outside the vignetting zone. We note that these two representations corresponds to different normalizations of the PSF (integral or maximal point normalized to 1). The exact curves obtained after propagation of waves through the entire coronagraph system are denoted with an index D, Φ_D, and S_D, with D for Diffraction. Also, Φ_D and S_G correspond to the simple geometrical model. At this scale, geometric and diffraction curves are almost identical. A more detailed analysis of the particular behavior of PSFs in the vignetting domain is given in Appendix A.

Figure 4 gives the level of the corona from the model of Cox (2015), namely, the level of this signal transmitted by the coronagraph after the Fredholm transformation and the final observed intensity after the addition of straylight, which appears as a bias in the observation. These curves are drawn in a linear scale. Denoting by C the true level of the solar corona, this operation is summarized as C × Φ + S L, where Φ and S L are the values for SIO and LIO (shown in Figs. 1 and 3).

These curves highlight the loss of signal near the solar limb caused by the use of LIO in comparison with SIO. On the contrary, it also highlight the much greater bias due to straylight when using the SIO occulter. This figure clearly demonstrates the fact that a compromise has to be found for the size of the internal occulter between straylight and closeness of observation near the solar limb. A zoom on the central part of matrix H is given for LIO and SIO (in Fig. 5, top and bottom panels, respectively).

Fig. 2 Radial cuts of a few PSFs computed for SIO (blue) and LIO (red). Point sources are distant one another of Δρ = 7.810−3 Rsun, or about 7.5 arcsec. Envelopes (black curves) correspond to Strehl ratios shown in Fig. 3.

Fig. 3 Integrated flux ΦG (geometrical model) and ΦD (end to end diffraction computation) for the small (SIO) and the large (LIO) internal occulter. Corresponding Strehl ratios SD and SG are also given. At this scale, geometric and diffraction curves are almost identical.

Fig. 4 Expected level of the corona (denoted by C) and observed intensities obtained as the product C × ΦD + S L, with ΦD the integrated flux (end-to-end diffraction computation) and SL is the solar straylight, for the small (SIO) and the large (LIO) internal occulter.

4.3 Object x

The inversion of the Fredholm transformation has to be verified on a realistic representation of the solar corona near the limb. In their recent study, Shestov et al. (2021) based their simulation on a direct image obtained during a real solar eclipse and synthetic object that present finer details than the eclipse image, Shestov et al. (2019). The object used in the present study lies halfway between these two objects.

We used a solar corona picture obtained by D. Sabers, R. Royer and M. Druckmüller during the total solar eclipse of March 9, 2016 in Indonesia, where the fine details of the coronal structure have been enhanced by image processing. We impose the synthetic object to follow the dynamic range and radial decrease of intensity of the true corona, as done by Shestov et al. (2021) on their synthetic image. An illustration of this object is given in Fig. 6 represented directly and in polar coordinates, as in DeForest et al. (2018), for example. Because of the size of the matrix H, a small part of this image can be used as a reference object in the present study. An example is given in Fig. 7. We took this region because it is radially homogeneous and allows a fair comparison between both internal occulters. It corresponds to a small region of 96 arcsec × 177 arcsec of the solar corona near the solar limb (from 1.0206 to 1.2041 in solar radius). Another part of the image is used further on in the paper. For the purposes of our study, these images are re-sampled with the sampling of 0.7497 arcsec per pixel used for the PSFs. The intensity levels have been modified to recover the dynamic of the corona given in Fig. 1. To do so, we forced the radial mean image to follow the model of Cox (2015) by multiplying the image by a factor dependent on ρ. Of course, other images of the corona could be used, such the ones shown in Habbal et al. (2010), Boe et al. (2020) for example, or other pictures from Druckmüller’s eclipse home page.

Fig. 5 Region of 1000 × 1000 points of the 30 208 × 30 208 matrix for LIO (top panel) and SIO (bottom panel) internal occulters.

4.4 Simulated observation of coronal features y

Figures 7 and 8 show the data material used for the simulations, in accordance with Eqs. (4) and (11). Figure 7 shows the images corresponding to the noiseless classical intensity for the SIO and LIO cases. Top panel of the figure corresponds to the effects of the vanishing PSF alone, while bottom panel of the figure includes the straylight, which appears very strong in the SIO case.

Figure 8 (top panel) gives an example of the observed images of Fig. 7 after photodetection, for an average of 10⁴ photons per pixel. In the bottom panel, the figure shows the same image with a subtraction of the bias. It should be noted that in this procedure, the calculation corresponds to performing a Poisson transformation of the bias (photonisation) and subtraction of its mean value. The top panel of Fig. 8 is the biased observation, bottom panel is the de-biased one. We note that in these two figures, the upper third of the images reproduce the same object in order to facilitate the comparison of the successive effects.

The subtraction of an estimate of the bias get rid of the vertical structures due to the scattered light, at the expense of a few negative values. At this high level of light, the debiased data gives an image much closer to the reality of the solar corona than the raw data.

Fig. 6 Solar eclipse observation of March 9, 2016 (thanks to Pr Miloslav Druckmüller) from which the data material is chosen. Top panel: direct image. Bottom panel: representation in polar coordinates of the highlighted region around the limb. The object x is a small fragment drawn from this image.

5 Numerical results

This section is dedicated to the results of numerical simulations. The data material used, y, H, and b are those mentioned more precisely in the previous section.

5.1 Correction of data by simple flux normalization

A very simple compensation of the variation of the flux of PSFs due to vignetting can be applied to the debiased noisy observation of Fig. 8 (bottom panel). The result is shown in Fig. 9. Here, again, the upper third of the image corresponds to the object. The observed images are simply divided by the flux Φ_D of Fig. 3 for SIO and LIO. To prevent divisions by zero or by very small terms, a constant value (here: 0.01) is added to Φ_D.

The result is that, in the vignetting region, the variation of the integrated flux of the PSF is compensated to a large extent, but not its variable blurring effect. The image quality sharply deteriorates with increasing vignetting effects. Moreover, the bias subtraction has made some of the data to become negative, as it can be seen in the range of the color bar (about 2000 pixels of the 30 000 pixels of the image are negative).

What we seek to obtain in the present study is to also compensate for this lack of image quality, or, in other words, to perform a reconstruction of the observed image inside the vignetting zone and recover a positive image. This fully justifies the choice of using a more sophisticated reconstruction algorithm, as presented in Sect. 2.2. The application to ASPIICS makes it possible to choose the term of regularization adapted to the object that is to be reconstructed.

Fig. 7 Triple representation of the same region of the solar corona: the object, its Fredholm transformation for SIO and LIO. The solar limb is on the left. Top panel: without the stray light. Bottom panel: with the stray light (from Fig 1). Color bar is in solar units (1 at disk center).

5.2 Choice of the regularization term

The regularization term J₂, as its name suggests, must regularize, namely, it must prevent the appearance of high frequencies in the reconstruction during the iterative process. In addition, it must take into account the specificities of the object to be reconstructed. The solar corona is an object with smooth variations and this justifies the choice of the distance, L2. This choice of distance also guarantees the convexity of the functional.

More specifically, we chose for J₂(x) a discrepancy function between the object, x, to a smoothed version of it called A_Rx, which corresponds to the effect of a perfect circular aperture of radius, R. At each iteration, the matrix, A_R, operates the convolution of x with an Airy pattern corresponding to the PSF of a telescope of radius R. Here, R becomes a parameter that has to be defined. After a few tests, we found that the value of R = 2.5 cm was suitable. This corresponds to the entrance aperture of the ASPIICS coronagraph. However, we did not carry out an exhaustive study on the value of this parameter. This regularization term prevents the amplification of high spatial frequencies beyond the cut-off of the perfect instrument.

We further give a particular and complex structure to J₂(x) by including the Strehl, S_i, in the regularization function, as follows: (28)

and the gradient is given by: (29)

Here, S_i plays a complementary role to the regularization factor γ of Eq. (7). The product γS_i weights the strength of the regularization according to the image quality. The regularization coefficient then decreases from γS_i = γ in the external region where S_i = 1 toward γS_i = 0 in the shadowed region close to the solar limb where S_i = 0. There, the algorithm becomes unregularized.

In practice, the maximum value of each column of H (hence, of each PSF) gives the quantity associated with the generalized Strehl ratio described in Sect. 4.4, which is represented in Fig. 3. The Strehl ratio, as with the PSF, is the same for a given radial distance to the solar edge. Therefore, in the matrix formulation, we denote using S_i the value of the Strehl ratio in the column i. We note that due to the lexicographic approach, S_i is periodic of the period, N.

We can see that in a region characterized by very strong vignetting, the PSFs are dominated by diffraction effects and present the surprising behavior shown in Fig. A.3 (as already discussed in Sect. 4.2). These PSFs are no longer centered on their corresponding point sources and their responses are shifted towards areas of less vignetting (see Fig. A.3). This aspect and related consequences are further discussed in Sects. 5.4 and 5.5. Simulations that are not presented here indeed confirm the efficiency of such a regularization by including the variable term S_i in comparison with a fixed term S_i = 1, for example.

Fig. 8 Triple representation of the same region of the solar corona: the object and the observations for SIO and LIO. Top panel: raw photodetected image for an average of 104 photons per pixel. Bottom panel: De-biased data (subtraction of an estimate of the straylight, inducing a few negative values). Color bar is given as the number of photons.

Fig. 9 Correction of data by simple flux normalization: debiased observations of Fig. 8 are divided by the flux ΦD for SIO and LIO (given in Fig. 3). Correction factors are limited to 100 to prevent division by zero. Some values are negative. From top to bottom: object, the results for SIO and LIO.

5.3 Implementation of the algorithmic reconstruction

The reconstruction is then performed according to the algorithm of Eq. (9). The regularization used is that described in the Eq. (28). Thus, to speed up the calculation at each iteration, k, the convolution of x by an Airy function (the term A_Rx) is performed by spatial filtering using matlab 2D discrete cosine transform dct2.

The quality criterion chosen is the reconstruction error, ϵ^k, calculated by: (30)

where x_t is the true object. The calculation of the error is restricted to areas where Φ_D is not too small, leading to two different areas for LIO (from 1.1050 to 1.2034 solar radius) and SIO (from 1.0472 to 1.2034 solar radius).

In the case of early stopping, that is, without regularization, γ = 0, the chosen reconstruction object is such that: (31) (32)

In the case of a regularized procedure, γ ≠ 0, the reconstructed object is given by: (33) (34)

We draw attention to the fact that early stopping requires some a priori knowledge of x_t, and is therefore difficult to implement in real cases. Conversely, the regularization procedure in its principle stabilizes the error and the algorithm can be stopped at any time after a large number of iterations is made, without allowing the errors to dominate.

A key problem is thus a correct choice of γ. In the simulations, a numerical set of regularization coefficient is used in the algorithm and the value chosen is the one giving the minimum of the error as a function of γ for large k, as shown in Eq. (34).

The regularization coefficient will depend on both the data and the functions J₁ and J₂. In the present case, the value γ = 10⁻⁶ is found to be the best for the Kullback functional, y_raw and the value γ = 10 is found to be the best for the Gaussian functional, y_deb. The order of magnitude difference in γ for the two models comes from the different expressions of J₁, Eqs. (14) and (26).

5.4 Object reconstruction using raw data

The results obtained in this section are those obtained with the raw noisy data, that is, without the subtraction of the estimate of the bias. Furthermore, Eq. (14) and the corresponding Poisson model described in Sect. 3.1 are used for processing the data.

Figure 10 gives the result of the direct inversion by minimisation of J₁(x), that is, using Eq. (7) with γ = 0, and stopping the iterations when the minimum of the criterion ϵ^k of Eq. (30) is obtained. Early stopping allows to obtain what is the best reconstruction that can be achieved in this un-regularized inversion.

In addition to a noisy granular structure, there is a surprising quasi-deterministic phenomenon of pseudo-fringes, which is especially clearly visible for the LIO case. The interpretation of the presence of these pseudo-fringes was not immediate. After several trials and errors based on simulations, we came to the conclusion that these pseudo fringes arose mainly due to the very particular shape of the PSF in the area of strong vignetting, especially, with the variation of flux of the PSF. Indeed, when a simulation is carried out by imposing a normalization of the PSF, that is to say a constant integral also in the vignetting zone, this phenomenon disappears. The phenomenon of parasitic fringes also disappears if the images are processed without the presence of bias. These simulation tests are not transcribed here for the sake of brevity. It can therefore be deduced from these observations that the phenomenon of parasitic fringes is linked both to the evanescent nature of the PSFs and to the presence of additive bias. A more in-depth study of these aspects would be very interesting from a data processing point of view, but this obviously goes beyond the objectives of this work.

The resolution for most of these disturbances was to implement the regularization, namely, to use the full Eq. (7) with γ = 10⁻⁶. The result is shown in Fig. 11. For SIO, the fringes have completely disappeared, while for LIO they remain partially present, but for a zone outside the domain of exploitation of these observations, namely, below the geometric cutoff due to vignetting.

Fig. 10 Reconstructed images from the raw data given in Fig. 8 (top panel), Poisson model (KL divergence of Sect. 3.1, γ = 0, early stopping). From top to bottom: object and the best reconstructions for SIO and LIO.

Fig. 11 Regularized reconstructed images from the raw data given in Fig. 8 (top panel), Poisson model (KL divergence of Sect. 3.1, γ = 10−6). From top to bottom: object and the best reconstructions for SIO and LIO.

5.5 Object reconstruction using debiased data

We have proposed two models for the debiased data: a Gaussian Poissonian model (Sect. 3.2.1) and a pure Gaussian model (Sect. 3.2.2). In practice, the implementation of the first model leads to a digital dead end, due to the very large number of photons; in particular, (Hx)_i and n_i can take values of several thousand photons, so the terms ! cannot be evaluated numerically. We therefore had to abandon this for the purposes of reconstructing the solar corona. Consequently, the results obtained in this section are those obtained with the Gaussian model described in Sect. 3.2.2.

As for the raw data, y_raw, Fig. 12 gives the result of the direct inversion by minimization of J₁(x), that is using Eq. (7) with γ = 0, and stopping the iterations when the minimum of the criterion ϵ^k of Eq. (30) is obtained. As already indicated for the Poisson case, this simulation assumes the knowledge of the object x_t, which makes it possible to obtain the best reconstruction in this non regularized inversion. Unlike what is obtained in the previous processing based on raw data, no parasitic fringes are present here and we find that these non-regularized results are already quite promising.

Figure 13 gives the result obtained by applying a regularization with a weight γ = 10. As noted previously, the optimal numeric values of γ are those which give the smallest error, as the number of iterations increases, (k = 1000 here). These results are also excellent. We can note a surprising appearance of very weak fringes for LIO in the completely obscured zone where no measurement can be expected.

In Fig. 14, another part of the solar corona is processed with the debiased data (see Sect. 3.2.2). Reconstructions are given for LIO (Fig. 14, top panel) and SIO (Fig. 14 bottom). Again, the use of SIO allows for a reconstruction of the corona closer to the solar limb.

Fig. 12 Reconstructed images from debiased data of Fig. 8 (bottom panel) based on the Gaussian model of Sect. 3.2.2, with γ = 0, early stopping. From top to bottom: object and the best reconstructions for SIO and LIO.

5.6 Error analysis as a function of ρ

A quantitative estimate of the reconstruction error is given in Fig. 15. For each point of the reconstructed image, an error taking the form of Eq. (30) is computed and then its average across the image is drawn as a function of the radial position, ρ.

Figure 15 gives the errors corresponding to the LIO (top panel) and SIO (medium panel) internal occulters. In each figure, the PSF geometrical limits of vignetting are drawn as vertical dashed lines (for LIO from 1.10093 to 1.17312 solar radius, and for SIO from 1.04249 to 1.11468 solar radius). The two regimes, below and above the lower geometrical limits, clearly appear.

Both algorithms, based on Poisson and Gauss statistics, give better results than the simple correction, which is only a satisfactory minimum. The Gaussian model is an approximate model, but it performs a bias subtraction before processing, giving, better results than the exact Poisson statistic. In particular, it does not present the problem of fringes which seems, at least in part, to be linked to the presence of the bias. Despite a strong stray light, the SIO gives a much better result than the LIO, in terms of the proximity of observation to the solar limb. However, if we analyze the error very closely (as in Fig. 15, bottom panel) the LIO shows an improvement over the SIO outside ρ > 1.14 solar units, but it is very slight.

Fig. 13 Regularized reconstructed images from debiased data of Fig. 8 (bottom panel) based on the Gaussian model of Sect. 3.2.2, with γ = 10. From top to bottom: object and the best reconstructions for SIO and LIO.

Fig. 14 Result obtained in another part of Fig. 6. Regularized reconstructed images from debiased data (Gaussian model of Sect. 3.2.2, γ = 10). Top panel: for LIO; top to bottom: object, data, and reconstructed object. Bottom panel: for SIO, with the same representations.

Fig. 15 Error as a function of ρ, for LIO (top panel) and SIO (medium panel). Here, SC: simple correction; REG and NR: regularized and non-regularized algorithms; raw (red curves) and deb (blue curves): raw data and data after bias subtraction, respectively. Vertical dashed lines represent the geometrical limits of vignetting. Bottom panel: zoom on the error of the above figures. Raw data: red and magenta curves. Bias subtracted data: Blue and cyan curves. See more details in the text.

6 Inverse Fredholm reconstruction for the imperfect coronagraph

In the previous section, the ASPIICS coronagraph was assumed to be perfect for the Fredholm inversion. To discuss further, we now consider some “defaults” from the optics or the satellite flying formation. These aspects have been analysed by Shestov et al. (2021) and Rougeot (2020). We limit here our analysis to formation flying errors and the effect of inevitable residual roughness of the optical surface of the telescope of the coronagraph.

6.1 Effects of formation flying errors

The novelty of the Proba-3 formation flying mission is that the external occulter is not mechanically attached to the coronagraph. While it offers of a clear advantage in terms of the stray light reduction and vignetting (Rougeot et al. 2017), it also comes with small and varying misalignment and off-pointing errors between the occulter and the instrument, due to the imperfect in-flight formation of the two satellites. The detailed description of the different errors and their respective impact on the diffraction is presented in Shestov & Zhukov (2018). As a reminder, the ideal configuration is when the entrance aperture is well centred in the umbra and the coronagraph boresight perfectly points toward the center of the sun and, thus, the center of the external occulter. Here, we consider the lateral displacements and off-pointing of the coronagraph.

Such formation errors make both the response of the coronagraph in the vignetting region and the straylight be asymmetric. On one side, the occulter masks more the inner field of view. The stray light is thus weaker and the vignetting more pronounced with respect to the ideal configuration. On the other side, it is the opposite, as the occulter edge is located closer to the rim of the solar disk. The direct consequence is that the bias b and the matrix H of the PSFs responses are naturally different from the ideal case that was previously considered. In the following, we denote, with H_FF and b_FF, respectively, the new matrix and stray light bias in the case of formation flying errors. Thus, Eq. (4) becomes: (35)

Moreover, the misalignment and coronagraph pointing attitude are time-varying and not steady, due to the periodic actuation of both satellites to maintain and control the formation (with a characteristic time of ten seconds). It is thus incorrect to assume that the new matrix, H_FF, and the estimated bias are known a priori at any given time.

To simulate such a situation, which is likely to occur in reality, we consider the case where the observations are made with positioning errors (i.e., with H_FF and b_FF) and the Fredholm inversion is performed as if the observations were perfect (i.e., with H and b). This approach corresponds to refraining from the “inverse crime”, as it is called in the field of signal processing (see, e.g., Colton et al. 1998 and Wirgin 2004). Therefore, the de-biased data of Eq. (16) now appear as follows: (36)

where we synthesize the data y_deb with H_FF and b_FF which contain formation flying errors effects, from which we subtract the stray light b expected for the perfect coronagraph.

In order to perform the Fredholm inversion in the case of such formation flying errors, we need to compute H_FF and b_FF. Here, we considered an extreme but realistic case of a 5mm lateral displacement between the entrance aperture and the occulter, combined to an off-pointing of the coronagraph boresight by 5 arcsec (Galano et al. 2019). The directions of both errors have been chosen such that their effects are cumulative. The data material used here were computed with the MATLAB code of Rougeot (2020; see Chapter 5 in particular). We note that the straylight is no more axis-symmetric in the case of a formation flying error: it is brighter in one part of the field of view and weaker in the other part. For the purposes of illustration, we only considered the latter case. Figure 16 shows cuts of the stray-light, namely, b_FF, which is compared to the ideal case, for both occulters SIO and LIO. The impact is more pronounced in the LIO case, as the relative difference reaches 17% at the stray light peak.

Figure 17 shows the best results that can be obtained when inverting the Fredholm transform using H and b instead of the unknown H_FF and b_FF true values. Obviously, the results are not as good as in the case of the perfect coronagraph, but we do not see an exponential increase of errors. The effect appears more visible for SIO, in particular in the region around 1.06 solar radius. This is probably due to the more important error in the bias subtraction there. Pseudo fringes are slightly visible in y_deb, and are amplified in the final reconstruction despite the use of a strong regularization term. The larger internal occulter appears more favourable in that case.

Fig. 16 Expected level of straylight for SIO in blue, and LIO in red. The plain curves show the nominal case without formation flying errors, as in Fig.1. The dotted curves correspond to the case with the formation flying error, in the region where the straylight is weaker due to the displacement of the external occulter.

6.2 Effects of surface roughness of the telescope

The surface roughness of optical surfaces produces micro phase defaults whose diffraction degrades the PSF structures. The modeling of such random micro-structures is non-trivial, as described in Rougeot et al. (2019). The same approach is used in this paper.

The surface roughness is modeled by a 2D random function which represents micro-deviations of glass thickness with respect to the ideal surface. Such defects are well known in optics to produce speckles. A rough approximation is to consider that the final degraded PSF is the sum of the ideal PSF plus a speckle figure, the two figures interfering in terms of complex amplitude.

For the almost perfectly polished surfaces of the coronagraph, the considered defaults are very small. The root mean square (rms) amplitude is on the order to λ/300, but with an extremely small correlation length of λ/2. The resulting speckle pattern is of a very small amplitude, while extending very far away. The speckle structure, which is very weak compared to the Airy spot, is only revealed very far from the center of the response, as was first clearly illustrated in Fig. A1 in the appendix of Rougeot et al. (2019). Indeed, the degraded PSF due to surface roughness was found to be almost identical to the perfect PSF up to an angular distance of 2 arcmin and the speckle regime becomes clearly dominant outside of 4 arcmin.

The same behavior is observed for PSFs in the vignetting region. We give an example in Fig. 18 for the LIO case and a point source situated at 1.1011 in units of solar disk. Pinned speckles are visible at the location of the airy rings, which appear to be broken; see Bloemhof et al. (2001) and Soummer et al. (2007) for more details. These speckles appear very far from the center of the response in Fig. 18 (bottom panel). It is necessary to make a comparison with the perfect case shown in the same Fig. 18 (top panel) to see them with some difficulty. The region shown in the figure is much larger (384 x 290 arcsecs) than the region used for the computation of PSFs to establish the matrix H.

For the order of magnitude of the microstructures considered in this analysis, these speckle effects do not modify the matrix, H, which is quite the same as for the perfect instrument. There is a noticeable effect for the straylight, b, which is very sensitive to the long-range behavior of diffraction originating from all point sources from the solar disk, as shown in Rougeot (2020). Here, this effect also remains quite weak in the region concerned by the Fredholm inversion, as shown in Fig. 19. It can be seen there that the straylight curve is slightly above that of the perfect instrument while integrated fluxes and Strehl for the perfect or speckle case remain strictly identical. Because of these small differences, we have not offered the same study here as for the formation flying errors.

Fig. 17 Processing under unknown formation flying errors. Top panel: for LIO (from top to bottom), the object, the image ydeb of Eq. (36), and the reconstructed image using H and b instead of the unknown HFF and bFF true values. Bottom panel: same details but for SIO. The wrong subtraction of the bias slightly visible in ydeb appears to be largely responsible for the rise in errors in the final result.

Fig. 18 PSF from the point source at 1.1011 Rsun obtained for perfect optics (top panel) and optics with surface roughness expected for ASPIICS (bottom panel). The pinned speckles are only visible with difficulty far away from the pic of the response, in the region x > 1.35 in the image. The field is of about 384 × 290 arcsec.

Fig. 19 In black, left y-axis, stray lights for the perfect and surface roughness (speckle) cases. Right y-axis, integrated flux Φ and Strehl values, S, for the two experiments are found to be identical.

7 Perspectives

This section considers the application of Fredholm inversion to other solar coronagraphs and then to exoplanet imaging using starshades.

7.1 Applicability to other externally occulted solar coronagraphs

In the present work, we deal with the case of imagery in the vignetting area of ASPIICS. The C2 and C3 coronagraphs from SOHO-LASCO, which are externally occulted coronagraph with strong vignetting, Llebaria et al. (2006), could be treated in a similar way. The problem will be to obtain a satisfactory estimate of H and b. Their computation will be more delicate than for ASPIICS for two reasons. One is the presence of the support of the external occulter in the field. The other is the short distance between the occulters and the telescope aperture, leading to very large Fresnel numbers, a difficulty analysed from a theoretical and numerical point of view in Aime (2020) and in Rougeot & Aime (2018). A study could also be done for the new concept of coronagraph METIS aboard the Solar Orbiter ESANASA observatory, Antonucci et al. (2020), but a complete study of the propagation would need to be carried out.

It is also possible to envisage a similar study for a coronagraph using an external serrated occulter. Such an occulter have long been considered for its qualities of straylight rejection, Newkirk & Bohlin (1965). According to Rougeot & Aime (2018), to be efficient, an external serrated occulter must have a very large number of teeth; typically, a 1000-tooth occulter is expected to provide a gain of 100 in the darkness of the shadow. Such a free-flying occulter would produce an angular periodicity of the impulse response of the order of a thousandth of a circumference, that is a periodicity of barely 6 arcseconds around the solar limb, a structure that the matrix, H, considered in the present work could handle numerically without difficulty.

7.2 Inverse Fredholm reconstruction for exoplanet imaging using starshades

Another totally different field of application may be the detection of exoplanets and the study of circumstellar environments. Indeed, if we consider that circumstellar fields may be very complex and structured, the problem of imaging of exoplanets does not simplify to imaging an ensemble of point sources, and the inversion of the Fredholm transformation may be of interest there for observations using starshades, Cash (2006). The same problem of an evanescent PSF arises, in a rather complex manner. In this case, the shaped occulters have a low number of petals (for example 20 in the computation of Cady 2012) with a particularly studied profile to obtain a very dark shadow, Vanderbei et al. (2007), Flamary & Aime (2014). These occulters are typically of tens of meters to tens of thousands of kilometers in front of a telescope of a few meters in diameters (see Arenberg et al. 2021 for a very pedagogic discussion of this issue).

Although the physical characteristics are very far from the problem of ASPIICS studied in this paper, it can be very similar for the numerical aspects. The use of occulter with a small number of petals will lead to a very structured vignetting. The calculus has to be extremely precise for the straylight (the term b), but probably less demanding for the vanishing PSF and the matrix, H. If the external occulter rotates, as in the Starshade Rendezvous Mission project, Seager et al. (2019), and if the exposure time is so to include an integer number of rotations, the observation recovers a radial symmetry similar to that of ASPIICS. For a non-spinning occulter, there will be an angular periodicity depending on the number of petals. In both cases the main problem will remain a precise estimation of H and b.

A rough illustration of the principle is made in Fig. 20, showing the results obtained by a simple Fresnel filtering. We present a shaped external occulter, with its Fresnel diffraction at the level of the telescope aperture, the light that this aperture would receive for three directions of source points around the star, and corresponding PSFs. Reversing the Fredholm transformation as we simulated here for ASPIICS could presumably improve near-star imagery close to the inner working angle of the system. It is likely that more elaborate techniques using the Maggi-Rubinowicz approach (Born & Wolf 1999, Cady 2012, Rougeot & Aime 2018) will need to be used to obtain precise results in a realistic case.

Fig. 20 Rough illustration of vignetting for exoplanet imaging using starshade. Top-left panel: example of an external occulter with 12 petals and a shaped profile following a Sonine curve of the form (1 − ρ2)2. The background image gives its Fresnel diffraction pattern at the level of the pupil of the telescope. The discs give the illumination on the pupil for 3 directions of point sources in the sky (in fact the pupil is fixed and the shadow moves according to the direction of observation, but the relative effect is the same). Bottom panel: PSFs corresponding to the three above directions on the sky.

8 Conclusion

In this paper, we develop a general matrix method of image reconstruction for data that have undergone a Fredholm transformation. The application is made for ASPIICS, a test solar object is drawn from an eclipse observation and the effect of vignetting is made using the same procedure as in Rougeot et al. (2017). The stray light is taken into account and simulations are carried out for two sizes of the internal occulter of Lyot’s coronagraph.

Signal processing and modeling aspects are deeply analyzed, in particular, the very important question of whether it is better to reconstruct the observation from raw data or from data from which an estimate of the straylight has been suppressed. We used specific algorithms for these two types of data, with a subtraction (or without) of an estimate of the stray light. A physical regularization was used in both cases. According to the study of errors we can derive a notable advantage from the procedure where a subtraction of an estimate of the straylight is done before processing the data; and this is regardless of the fact that we had to use a mere Gaussian model instead of the fundamental Poisson model. In summary, the results are excellent. Despite the fact that it produces a greater amount of stray light, the preference goes to the small internal occulter. However, this is a solution that will not be chosen in the final project, for reasons of experimental safety and tolerance to positioning errors.

The main study was carried out for the purposes of a perfect experiment, for which we know exactly the PSFs and the level of straylight. We have also processed the case were the observation and inversion were not realized using the same PSF and bias.This was done in Sect 6.1 for formation flying errors, using H_FF and bias, b_FF, to synthesize the data and just H and bias, b, used to reconstruct the object. Of course the results are not as good as if every thing is known, but the image reconstruction remains interesting. A study made for optical quality defects, in the form of microstructures producing speckles in the PSF, showed that these effects were rather negligible in the region of vignetting. This can be explained by the fact that the main effect, in the vignetting zone, is the geometric variation in the shape of the pupil while the wave defects intervene as second-order effects.

This work led us to define a procedure for processing future ASPIICS data. However, this study remains theoretical, as in practice, such a fine knowledge of the PSFs may not be achievable. The next treatments will have to be done based on actual observations.

We considered generalizations of this work for other solar coronagraphs using external occulters, such as SOHO-LASCO C2 or even the new coronagraph METIS. Because of the difficulties of calculation for high Fresnel numbers involved in these experiments, as well as because of a more complex geometry, the calculation of H and b will be difficult to obtain there. It would be simpler to consider a system similar to ASPIICS equipped with a serrated occulter as computed by Rougeot & Aime (2018), or even with two circular occulters, as studied by Aime (2020), but such projects are not envisaged at the present time.

A very promising perspective on this work might be to apply the inverse Fredholm technique to imaging circumstellar fields and exoplanets, Arenberg et al. (2021). This point is presented briefly in Sect 7.2 and this is a study that we intend to develop in a future work.

Acknowledgements

We thank Prof. RNDr. Miloslav Druckmüller, CSc. of the Institute of Mathematics, Faculty of Mechanical Engineering Brno University of Technology, Czech Republic for the use of eclipse images.

See: www.zam.fme.vutbr.cz/~druck/. Thanks are due to Dr. Mamadou N’Diaye for interesting discussions.

Appendix A PSFs in the extreme zone of vignetting

In this appendix, we describe the behavior of PSF in the extreme vignetting zone and, in particular, the existence of a response where simple geometry would give total darkness. Figure 3 of Sect. 4.2 already shows that the flux Φ and the Strehl for diffraction computations extended over a larger zone than expected from a pure geometrical point of view.

Fig. A.1 Zoom on the very low values of the flux for the small internal occulter (SIO), shown in Fig. 3. Below ρ = 1.0425 Rsun (dark point) the light is geometrically blocked by the occulter (red points) but faint diffracted light enters the entrance aperture (blue points).

The term “geometric” corresponds to the area of the pupil free of occultation, as shown in Bayanna et al. (2011) and Shestov et al. (2021). It is indeed a problem of intersecting circles (the telescope aperture and the external occulter). According to Weisstein (2003), we have: (A.1)

where R is the radius of the external occulter, increased by internal occulters to 1.02060 (SIO) or 1.0759 (LIO) times its size, r is the radius of the Lyot coronagraph aperture, reduced to 0.97 times its size by the Lyot stop and d as the angular radial position, p, converted in meters by a factor which corresponds, for a point at the aperture center, to the radius of the projected solar image onto the external occulter. It depends, of course, on the (variable) distance to the Sun. For the parameters of our study, this radius was 0.6718-m.

The difference between Φ_G and Φ_D is clearly visible in Fig. A.1 where PSFs vanish. The example given in the figure is for SIO, but a similar behavior is observed for LIO. As a result, diffraction makes a PSF response to exist where the geometrical calculation gives a total darkness. For example, Fig. A.1 shows that the geometric transmission becomes equal to zero below ρ = 1.0425 Rsun. Fig. A.2 shows the response to a point source at that exact position. In the figure, the position of the point source is marked with a white circle, and the response is surprisingly shifted to the right of that point.

This behavior is general for the whole region of strong vignetting. An illustration of that is given in Fig. A.3 that shows five positions of point sources (vertical lines) and cuts of corresponding PSFs for SIO (x-axis between 1.03 to 1.08) and five for LIO (x-axis between 1.08 and 1.14). To easily match point sources with corresponding PSFs, vertical lines are drawn with heights corresponding of the maximum value of corresponding responses and drawn in the same color. For both SIO and LIO cases, PSFs are shifted to the right (toward larger radial positions) and appear at the same angular position (green dashed lines). For example, for LIO, the five point sources are situated at positions: 1.0933, 1.0948, 1.0964, 1.0979, and 1.0995, while the responses all have their maximum at about 1.11, in units of solar radius. For SIO, the same behavior is present, shifted by about 0.0156 solar unit.

Fig. A.2 Point source at the zero geometrical transmission position (ρ = 1.0425 of Fig. A.1) and corresponding shifted response.

Fig. A.3 Positions of point sources (vertical lines) and corresponding PSFs drawn with the same color. Responses are shifted to almost the same position (vertical green lines) for SIO and LIO, respectively.

These PSFs, that exist only because of diffraction effects of the external occulter, are very weak, thus they have little impact on the final image. However these weak responses have a noticeable effect for the inversion of the Fredholm transform. Indeed, since there is no more bijection between responses and point sources, it is difficult to recover the direction from which an intensity has come.

A probable explanation for this behavior may be that The Fresnel diffraction by the occulter induces effects on both the modulus and the phase of the diffracted wave at the entrance aperture of the coronagraph. In the bright region of the Fresnel diffraction, the phase variation is very small, while it rotates many times 2π in the shadowed region. This phase effect produces a tilt of the wave front. The light is directed toward the centre of the pattern. This effect is at the origin of the well known spot of Arago behind the occulting disk. As the pupil of the coronagraph enters into the shadowed region, the diffracted wave gets more and more affected by this phase variation, which amplifies the effect of the tilt.

References

All Figures

	Fig. 1 Expected level of the corona (model of Cox 2015) and straylights for the small (SIO) and the large (LIO) internal occulter
In the text

	Fig. 2 Radial cuts of a few PSFs computed for SIO (blue) and LIO (red). Point sources are distant one another of Δρ = 7.810−3 Rsun, or about 7.5 arcsec. Envelopes (black curves) correspond to Strehl ratios shown in Fig. 3.
In the text

	Fig. 3 Integrated flux ΦG (geometrical model) and ΦD (end to end diffraction computation) for the small (SIO) and the large (LIO) internal occulter. Corresponding Strehl ratios SD and SG are also given. At this scale, geometric and diffraction curves are almost identical.
In the text

	Fig. 4 Expected level of the corona (denoted by C) and observed intensities obtained as the product C × ΦD + S L, with ΦD the integrated flux (end-to-end diffraction computation) and SL is the solar straylight, for the small (SIO) and the large (LIO) internal occulter.
In the text

	Fig. 5 Region of 1000 × 1000 points of the 30 208 × 30 208 matrix for LIO (top panel) and SIO (bottom panel) internal occulters.
In the text

	Fig. 6 Solar eclipse observation of March 9, 2016 (thanks to Pr Miloslav Druckmüller) from which the data material is chosen. Top panel: direct image. Bottom panel: representation in polar coordinates of the highlighted region around the limb. The object x is a small fragment drawn from this image.
In the text

	Fig. 7 Triple representation of the same region of the solar corona: the object, its Fredholm transformation for SIO and LIO. The solar limb is on the left. Top panel: without the stray light. Bottom panel: with the stray light (from Fig 1). Color bar is in solar units (1 at disk center).
In the text

	Fig. 8 Triple representation of the same region of the solar corona: the object and the observations for SIO and LIO. Top panel: raw photodetected image for an average of 104 photons per pixel. Bottom panel: De-biased data (subtraction of an estimate of the straylight, inducing a few negative values). Color bar is given as the number of photons.
In the text

	Fig. 9 Correction of data by simple flux normalization: debiased observations of Fig. 8 are divided by the flux ΦD for SIO and LIO (given in Fig. 3). Correction factors are limited to 100 to prevent division by zero. Some values are negative. From top to bottom: object, the results for SIO and LIO.
In the text

	Fig. 10 Reconstructed images from the raw data given in Fig. 8 (top panel), Poisson model (KL divergence of Sect. 3.1, γ = 0, early stopping). From top to bottom: object and the best reconstructions for SIO and LIO.
In the text

	Fig. 11 Regularized reconstructed images from the raw data given in Fig. 8 (top panel), Poisson model (KL divergence of Sect. 3.1, γ = 10−6). From top to bottom: object and the best reconstructions for SIO and LIO.
In the text

	Fig. 12 Reconstructed images from debiased data of Fig. 8 (bottom panel) based on the Gaussian model of Sect. 3.2.2, with γ = 0, early stopping. From top to bottom: object and the best reconstructions for SIO and LIO.
In the text

	Fig. 13 Regularized reconstructed images from debiased data of Fig. 8 (bottom panel) based on the Gaussian model of Sect. 3.2.2, with γ = 10. From top to bottom: object and the best reconstructions for SIO and LIO.
In the text

	Fig. 14 Result obtained in another part of Fig. 6. Regularized reconstructed images from debiased data (Gaussian model of Sect. 3.2.2, γ = 10). Top panel: for LIO; top to bottom: object, data, and reconstructed object. Bottom panel: for SIO, with the same representations.
In the text

Fig. 15 Error as a function of ρ, for LIO (top panel) and SIO (medium panel). Here, SC: simple correction; REG and NR: regularized and non-regularized algorithms; raw (red curves) and deb (blue curves): raw data and data after bias subtraction, respectively. Vertical dashed lines represent the geometrical limits of vignetting. Bottom panel: zoom on the error of the above figures. Raw data: red and magenta curves. Bias subtracted data: Blue and cyan curves. See more details in the text.

In the text

	Fig. 16 Expected level of straylight for SIO in blue, and LIO in red. The plain curves show the nominal case without formation flying errors, as in Fig.1. The dotted curves correspond to the case with the formation flying error, in the region where the straylight is weaker due to the displacement of the external occulter.
In the text

Fig. 17 Processing under unknown formation flying errors. Top panel: for LIO (from top to bottom), the object, the image ydeb of Eq. (36), and the reconstructed image using H and b instead of the unknown HFF and bFF true values. Bottom panel: same details but for SIO. The wrong subtraction of the bias slightly visible in ydeb appears to be largely responsible for the rise in errors in the final result.

In the text

	Fig. 18 PSF from the point source at 1.1011 Rsun obtained for perfect optics (top panel) and optics with surface roughness expected for ASPIICS (bottom panel). The pinned speckles are only visible with difficulty far away from the pic of the response, in the region x > 1.35 in the image. The field is of about 384 × 290 arcsec.
In the text

	Fig. 19 In black, left y-axis, stray lights for the perfect and surface roughness (speckle) cases. Right y-axis, integrated flux Φ and Strehl values, S, for the two experiments are found to be identical.
In the text

Fig. 20 Rough illustration of vignetting for exoplanet imaging using starshade. Top-left panel: example of an external occulter with 12 petals and a shaped profile following a Sonine curve of the form (1 − ρ2)2. The background image gives its Fresnel diffraction pattern at the level of the pupil of the telescope. The discs give the illumination on the pupil for 3 directions of point sources in the sky (in fact the pupil is fixed and the shadow moves according to the direction of observation, but the relative effect is the same). Bottom panel: PSFs corresponding to the three above directions on the sky.

In the text

	Fig. A.1 Zoom on the very low values of the flux for the small internal occulter (SIO), shown in Fig. 3. Below ρ = 1.0425 Rsun (dark point) the light is geometrically blocked by the occulter (red points) but faint diffracted light enters the entrance aperture (blue points).
In the text

	Fig. A.2 Point source at the zero geometrical transmission position (ρ = 1.0425 of Fig. A.1) and corresponding shifted response.
In the text

	Fig. A.3 Positions of point sources (vertical lines) and corresponding PSFs drawn with the same color. Responses are shifted to almost the same position (vertical green lines) for SIO and LIO, respectively.
In the text