Calibration of residual aberrations in exoplanet imagers with large numbers of degrees of freedom

Imaging faint objects, such as exoplanets or disks, around nearby stars is extremely challenging because host star images are dominated by the telescope diffraction pattern. Using a coronagraph is an efficient solution for removing diffraction but requires an incoming wavefront with good quality to maximize starlight rejection. On the ground, the most advanced exoplanet imagers use extreme adaptive optics (ExAO) systems that are based on a deformable mirror (DM) with a large number of actuators to efficiently compensate for high-order aberrations and provide diffraction-limited images. While several exoplanet imagers with DMs using around 1500 actuators are now routinely operating on large telescopes to observe gas giant planets, future systems may require a tenfold increase in the number of degrees of freedom to look for rocky planets. In this paper, we explore wavefront correction with a secondary adaptive optics system that controls a very large number of degrees of freedom that are not corrected by the primary ExAO system. Using Marseille Imaging Testbed for High Contrast (MITHiC), we implement a second stage of adaptive optics with ZELDA, a Zernike wavefront sensor, and a spatial light modulator (SLM) to compensate for the phase aberrations of the bench downstream residual aberrations from adaptive optics. We demonstrate that their correction up to 137 cycles per pupil with nanometric accuracy is possible, provided there is a simple distortion calibration of the pupil and a moderate wavefront low-pass filtering. We also use ZELDA for a fast compensation of ExAO residuals, showing its promising implementation as a second-stage correction for the observation of rocky planets around nearby stars. Finally, we present images with a classical Lyot coronagraph on MITHiC and validate our ability to reach its theoretical performance with our calibration.


Introduction
Imaging extrasolar planets is one of the most demanding endeavors in astronomy today. Collecting photons from these planets will yield precious astrophysical information on their nature, leading to knowledge of their orbital and rotational parameters (Galicher et al. 2014;Vigan et al. 2016a;Chilcote et al. 2015;Maire et al. 2016;Zurlo et al. 2016), clues on their location and formation history (Crepp et al. 2012;Öberg et al. 2011;Piso et al. 2016), and constraints on the chemistry at stake in their atmosphere (Knutson 2007;Phillips et al. 2020). However, imaging an exoplanet, and in particular an Earth analog, is extremely challenging because of the flux ratio (also referred to as the contrast) of at least 10 −8 with its host star at angular separations shorter than 50 mas (Traub & Oppenheimer 2010). Using a coronagraph on a large telescope is a compelling approach for attenuating the light of an observed star image and providing the sensitivity to collect sufficient photons from the planet at high angular resolution.
For large telescopes on the ground, angular resolution does not scale with the telescope's primary mirror diameter as Earth's atmospheric turbulence blurs the observations. Adaptive optics (AO) systems have been developed with a wavefront sensor (WFS) and a deformable mirror (DM) to measure and compensate for the aberrations introduced by the atmosphere in real time, leading to images at nearly the theoretical resolution limit (e.g., Beuzit et al. 1997;Herriot et al. 1998;Rousset et al. 2003;Roddier 2004). The accurate calibration of AO systems relies on the computation of the interaction matrix, a numerical tool relating the response of WFS subapertures to DM actuations.
The latest generation of high-contrast imaging instruments, such as the Very Large Telescope (VLT) Spectro-Polarimetric High-contrast Exoplanet REsearch (SPHERE) (Beuzit et al. 2019), Gemini Planet Image (GPI) (Macintosh et al. 2014), and the Subaru Coronagraphic Extreme Adaptive Optics (SCExAO) (Jovanovic et al. 2015), include novel extreme adaptive optics (ExAO). These systems are characterized by a wavefront correction with a high temporal frequency, up to 3.5 kHz, and a DM with a large number of actuators, up to 40 across the telescope pupil diameter. While allowing for correction in a field of view up to ∼1 ×1 , ExAO enables the production of images with Strehl ratios up to 95% in the near infrared (e.g., Sauvage et al. 2016). Such image quality translates into an increase in the amount of light gathered in the core of the point spread function (PSF) of the planet and, therefore, in the signal-to-noise ratio (S/N). In addition, the high image quality allows for effi-Article number, page 1 of 14 arXiv:2103.15844v1 [astro-ph.IM] 29 Mar 2021 A&A proofs: manuscript no. paper cient starlight suppression at large separations to observe faint structures, such as protoplanetary and debris disks, around an observed bright star.
These instruments have proven to be remarkably efficient, with contrasts down to 10 −6 at 0.3 angular separations, leading to several new young substellar mass companion discoveries (Macintosh et al. 2015;Perrin et al. 2015;Konopacky et al. 2016;Chauvin et al. 2017;Keppler et al. 2018;Cheetham et al. 2018), many newly imaged disks (e.g., Avenhaus et al. 2018;Ginski et al. 2016;Stolker et al. 2016;de Boer et al. 2016), and studies of other astrophysical objects at unprecedented resolutions (e.g., Kervella et al. 2015;Schmid et al. 2017). Extensive surveys such as the Gemini Planet Imager Exoplanet Survey (GPIES) with Gemini/GPI (Nielsen et al. 2019) or the SpHere INfrared survey for Exoplanets (SHINE) with VLT/SPHERE  have brought new constraints on the population of young giant planets.
Their success has also shed light on the limiting factors that prevent observations at even deeper contrasts and smaller separations. One of the main limitations for coronagraphy are noncommon path aberrations (NCPAs), which are optical distortions on the instrument science path that are unseen by the AO system on its WFS path. As they slowly vary with time (Soummer et al. 2007), these quasi-static aberrations do not average over time and cannot simply be calibrated as a bias. Over the past few years, several strategies have been proposed to calibrate these residual errors (see, for example, the review by Jovanovic et al. 2018).
Among these methods, the Zernike wavefront sensor (ZWFS; Zernike 1934;Bloemhof & Wallace 2003;Dohlen 2004;Wallace et al. 2011) is one of the leading approaches for NCPA correction in exoplanet imagers. Based on this approach, the Zernike sensor for Extremely Low-level Differential Aberration (ZELDA; N'Diaye et al. 2013) was developed, installed, and tested on VLT/SPHERE, showing promise as a calibration tool (N'Diaye et al. 2016;Vigan et al. 2019). This device allowed for diagnosing the low-wind effect, which has been crippling early SPHERE observations (Sauvage et al. 2015). It has also been explored in simulations for alternative applications, such as fine telescope segment co-phasing (Janin-Potiron et al. 2017a). The ZWFS is set to be implemented on the Nancy Grace Roman Space Telescope, a future NASA flagship mission (Shi et al. 2016) with high-contrast capabilities for exoplanet observations, as well as in future ground-based extremely large telescope (ELT) instruments (e.g., Carlotti et al. 2018).
When looking at the instrumentation requirements for the ELTs currently being built, there is a clear trend for larger numbers of actuators. For example, the M4 mirror in the European Southern Observatory (ESO) ELT 1 will have approximately 80 actuators across the pupil diameter (Vernet et al. 2012). This DM will provide an actuator pitch of ∼50 cm and run at a maximum rate of 1 kHz, which will be sufficient to provide images with a ∼70% Strehl ratio in the K band for the first light instruments (Thatte et al. 2016;Davies et al. 2016) in single-conjugated AO mode. While this will enable detection capabilities for studying warm or massive gaseous planets Houllé et al. submitted), it will not be sufficient to provide ExAO performance for observing cold and light rocky planets (Kasper et al. 2010). For such a challenging task, an actuator pitch of ∼20 cm and an AO loop running at 2-3 kHz is required (Kasper et al. 2011(Kasper et al. , 2013, resulting in nearly 200 actuators across the pupil diameter. Similarly, and although it is composed of several inde-1 https://www.eso.org/sci/facilities/eelt/ pendent primary mirrors and associated deformable secondary mirrors, the Giant Magellan Telescope (Johns et al. 2012) will provide a pitch of ∼20 cm at first light with close to 120 actuators across the pupil diameter. The number of actuators also tends to increase for space applications. For example, Ruane et al. (2020) performs phase conjugation with two DMs from Boston Micromachines, yielding 50 actuators across the diameter of the active area for each device. It is therefore crucial to investigate the wavefront correction possibilities that are offered with a very large number of degrees of freedom.
To handle this new range of high-order spatial frequencies, we explore the use of a second-stage AO system, including a ZWFS and a spatial light modulator (SLM) implemented downstream from the first ExAO system. To make the most of the ZWFS formalism developed by N' Diaye et al. (2013), we propose the use of a phase conjugation technique with a calibration procedure that uses probes that are performed in a two-frame acquisition sequence. We thus do not compute any standard interaction matrix: Each iteration implies solving a second-order polynomial equation per phase measurement bin.
In Sect. 2, we present the Marseille Imaging Testbed for High Contrast (MITHiC) with its 2019 upgrade and a ZELDA prototype to test the measurement and correction on a very large number of degrees of freedom. We also describe the closed-loop tests using ZELDA as a WFS. In Sect. 3, the experimental procedure is detailed with the ZELDA formalism, the Fourier filtering process, and the two-frame calibration of pupil distortion. In Sect. 4, we present the test results for static phase error correction and for an on-bench simulation of real-time SPHERE residual turbulence. Under certain assumptions, we show that it is possible to measure and compensate for aberrations down to a few nanometers in standard deviation in the optical path difference (OPD) with various input phase errors.

The MITHiC test bed
MITHiC was initiated in 2010 to develop and test technologies for exoplanet high-contrast imaging (N'Diaye et al. 2012b(N'Diaye et al. , 2014. It was used as a test bed for wavefront sensing techniques such as ZELDA (N'Diaye et al. 2013) or COronagraphic Focalplane wave-Front Estimation for Exoplanet detection (COFFEE) , coronagraphic concepts such as the apodized pupil Roddier & Roddier or the dual-zone phase mask coronagraphs (Soummer et al. 2003a,b;N'Diaye et al. 2010N'Diaye et al. , 2012b, and dark hole algorithms (Herscovici-Schiller 2018; Leboulleux 2018), as well as other science, such as detector characterization (Gach et al. 2019). Figure 1 presents the current optical setup of the bench. It is divided into three major parts, which are briefly described here and in much greater detail in Appendix A.
The telescope part goes from the light source to a dedicated phase screen and simulates a telescope observing an unresolved star. The light source is a fiber-coupled superluminescent diode, which has a central wavelength λ = 670.7 nm and a spectral bandwidth of 10 nm. The source is linearly polarized to match the SLM requirements precisely. A phase screen manufactured by SILIOS Technologies is then used to introduce various types of aberration patterns that are typically encountered in large optical telescopes. It includes both static patterns and two continuous rings of simulated residual turbulence for VLT/SPHERE and the ELT High Angular Resolution -Monolithic -Optical and Nearinfrared -Integral field spectrograph (HARMONI) (Vigan et al. 2016b).
The following part is the second-stage AO system that is based on an SLM for the phase correction and a ZELDA WFS. The pupil of the system is reimaged on an SLM that mimics a DM to apply phase corrections. Our implementation allows for a particularly large number of subapertures on the WFS and a very high density of correction elements on the SLM, with, respectively, 418 measurement points and 274 pixels across the pupil diameter, D. As a downside, using an SLM as a wavefront corrector imposes the use of linearly polarized monochromatic light. The ZELDA WFS is implemented using a phase mask with a depth of 350 nm, corresponding to nearly λ/4 in optical depth at 670.7 nm, and a diameter of 64 µm, corresponding to a relative size of 1λ/D (N'Diaye et al. 2013). It is mounted on a three-axis mount that enables its positioning with 1 µm accuracy. The pupil is then reimaged with a diameter of 418 pixels on the detector of ZELDA, which provides more than 130,000 single-resolution elements for wavefront sensing. The acquisition procedure and the associated formalism for wavefront measurements are detailed in Sect. 3. The WFS can be operated for single wavefront measurements, but it is best used in a closed loop with the SLM, as presented in the following sections. Without any correction from the SLM, the total amount of aberrations is approximately equal to 35.5 nm root mean square (RMS) on the wavefront and is mainly due to the spherical aberration from cumulative wave-front errors from the lenses. Using the phase error compensation described in Sect. 3, the level of residual phase aberrations can ultimately be reduced down to ∼2 nm RMS, which is close to the theoretical quantization limit of the SLM at 0.76 nm RMS (see details in Sect. A).
Finally, the last part of the bench corresponds to the coronagraph and associated science imaging channel. The test bed currently implements a classical Lyot coronagraph (CLC) based on an opaque focal plane mask (FPM) and an associated Lyot stop. In the final focal plane, the science imaging camera provides a fast readout with low noise to obtain high-contrast images. The camera can be used to image the coronagraphic signal either in the focal plane or in the pupil plane with an optional lens.

Principle of wavefront correction in a closed loop
We used the bench to perform closed-loop corrections of various types of phase errors. The phase measurements were performed with ZELDA, whose principle is recalled in Sect. 3.1. As this sensor naturally provides phase error maps, we worked with an interaction-matrix-free correction by applying the opposite measured phase map to the SLM. The algorithm is described in Sect. 3.2. To overcome distortion in pupil conjugation, which is automatically dealt with in a classical system relying on an interaction matrix, we performed an additional calibration step, detailed in Sect. 3.3. We finally performed an additional lowpass spatial filtering to avoid aliasing propagation and calibration errors (see Sect. 3.4).

Phase computation with ZELDA
ZELDA uses a focal plane phase mask on the point source image to produce interferences between the light going through and surrounding the mask in a downstream pupil plane. This results in pupil intensity variations, I c , that are related to the phase aberrations, ϕ, upstream from the phase mask. In the regime of small aberrations (ϕ 1 rad) and assuming a second-order Taylor expansion for ϕ, the intensity-phase relation from N' Diaye et al. (2013) yields where θ is the phase shift introduced by the mask, P denotes the amplitude profile in the entrance pupil (a top-hat function in our case), and b represents the real amplitude profile, which only depends on the size of the ZELDA mask and can be calculated once and for all. Retrieving the phase map is obtained by solving a second-degree polynomial equation for each pixel in the detector pupil plane. The OPD map is simply derived by inverting the relation ϕ = 2πOPD/λ. This formalism is cost-efficient in terms of computation.

Aberration compensation on MITHiC
On MITHiC, the pupil is sampled with 145,000 subapertures in the pupil on the WFS camera and more than 60,000 DM-like actuators on the SLM. With such large numbers of degrees of freedom, we implemented the direct correction of the aberrations without the requirement for a full-fledged interaction matrix; it would be extremely time-consuming to compute on MITHiC, which is not optimized for high-speed operations.
The first step is to acquire an initial OPD map, which requires two acquisitions on the WFS camera: one with the ZELDA mask and one without it for flux normalization, as described in N'Diaye et al. (2013). Both images are dark-corrected and re-centered. The OPD map is then computed using Eq. (1). In practice, we use the pyZELDA 2 package for the phase reconstruction . With 418 pixels across the pupil diameter, the wavefront measurement with ZELDA allows us to retrieve extremely fine phase information with spatial frequencies up to 209 cycles per pupil (c/p).
Since the pupil is sampled with 1.5 times more pixels by the ZELDA camera than the SLM, we first tried to directly apply the negative of a resampled phase map onto the SLM to compensate for the measured aberrations. However, as we will see in Sect. 4, this leads to a poor wavefront correction in a closed loop due to the geometrical distortion of the pupil. We therefore perform an additional calibration, as well as a spatial filtering of the computed OPD. Once these calibrations are applied, the processed OPD map is downscaled using spline interpolation to get a pupil diameter of 274 pixels to match the pupil size on the SLM, converted into phase, and displayed on the SLM.
For phase errors larger than 40 nm RMS, the ZELDA wavefront reconstruction is outside the sensor linearity range and the error is underestimated (N'Diaye et al. 2016). Therefore, a single correction step will, in general, not be sufficient to fully compensate for the phase errors of the bench. Repeating the previous procedure in a closed-loop fashion allows for convergence in a small number of iterations.

Compensation for the geometrical distortions
In the absence of a classical interaction matrix, we calibrated the system to compensate for geometric pupil distortion between the WFS camera and the SLM. The calibration step was performed by introducing a grid of N ×N Gaussian phase spots on the SLM, for which we know the exact theoretical position. The input coordinates, x i , y i , of each spot, i, can be represented by a column vector, V i . In the OPD map measured with ZELDA, we then retrieve the output position vector of the spots, V i . If we consider an affine distortion, there exists a 2 × 2 matrix, M, and a 2 × 1 column vector, C, such that, for all spots i, the V i coordinates satisfy: We fit this transformation using a least-square algorithm to determine the M and C coefficients to solve the following optimization problem: Finally, using the estimatedM andĈ from the fit, we applied the inverse transform to each spot of our OPD maps to derive V i : (4) Figure 2 illustrates this process with a 7 × 7 grid. This square grid does not probe all of the surface of the pupil disk but only an inscribed square centered in the pupil. Nevertheless, the 49 introduced dots already ensure the convergence of the determination of the six parameters of M and C. A more complete sampling would have to be tested even if it might not yield increased precision. The computing cost of this calibration sums up to: two frame acquisitions, 7×7 Gaussian fits on 32 ×32 pixel subarrays, and a linear regression to fit a six-parameter function on 49 points, which is performed within a few tens of seconds. Despite the fact that the bench is not stabilized in temperature or humidity, we found that performing the calibration once per day is sufficient for our usage. Figure 2 (left) shows a clear difference between the introduced calibration spots (arrow ends) and the measured calibration spots (arrow tips). The error shows not only a shift to the left, which can be seen as a misalignment, but also an additional, more complex distortion pattern. When removing the pure shift component, as in Fig. 2 (center), the maximum distance between the introduced and measured spot positions is 6.96 pixels on the ZELDA camera, with an average of 4.27 pixels. In terms of pupil percentage, this average displacement writes as 0.01D, which is a sub-actuator displacement for a VLT/SPHERE-like DM. The affine approximation to compensate for this error provides accurate results as the position of the spots after correction closely matches the expected position. The maximum distance between the introduced and corrected spots is shorter than 0.6 pixels, which corresponds to less than 0.002D.

Initial displacement -4x
Shift component removed -4x Residuals -40x The final displacements are presented in Fig. 2 (right). The average distance between an introduced point and a corrected point is 0.15 pixels, which shows an efficient correction. Different grid sizes were tested, but N larger than 7 do not provide any quantitative improvements of the distortion correction; as such, we use this value from here on out. This calibration method can be easily implemented on any instrument already hosting a Zernike WFS, such as VLT/SPHERE (Beuzit et al. 2019) and the Roman Space Telescope (RST) Coronagraphic Instrument (CGI) (Kasdin et al. 2020), or on test beds, such as the High Contrast Imaging Testbed (HCIT) ). This calibration is, in fact, applicable to any wavefront sensor that provides measurements in the pupil plane, such as the pyramid WFS (Ragazzoni 1996), or more generally with most of the class of Fourier-based WFSs (Fauvarque et al. 2016).

Low-pass filtering
Aliasing effects will appear when a wavefront measurement at very high spatial frequency is applied on a much lower-order DM (or equivalent). This was previously highlighted for NCPA correction with ZELDA in VLT/SPHERE, where the OPD maps had to be either low-pass filtered in the Fourier space (N'Diaye et al. 2016) or alternatively projected onto the controllable modes of the ExAO system ) before being applied on the DM. In the present case, we implemented an optional low-pass filtering procedure by multiplying the Fourier transform of the OPD map by a Hann window of a specific width, as illustrated in Fig. 3. We adapted the window width, w, to match the desired cutoff frequency. As a definition of the filter value, we set w to be the spatial frequency value at which the Hann window reaches zero. Therefore, a low-pass filter at w c/p yields a Hann window diameter of 2wλ/D. Using this method of filtering is analogous to choosing the number of actuators and the influence function of the equivalent DM. With a filter cutoff frequency at w c/p, we modeled a DM with 2w actuators. With a Hann window, the full width at half maximum of the corresponding actuators' influence function is D/2w. Figure 4 highlights the radial profile of the corresponding influence function, computed by a Fourier transform of the filtering window. The influence function is Gaussian-like and is therefore a good approximation of the influence function

Results
In this section we present our results in different configurations. We demonstrate first a static configuration that corresponds to the compensation for various static phase maps, and second, a dynamic configuration where we introduce residual atmospheric aberrations representative of the on-sky performance of the SPHERE ExAO (SAXO) system. In this configuration we tried to correct for the ExAO residual aberrations in pseudo-real time by using the tools developed in the static case. Finally, we looked at the improvement in MITHiC coronagraphic images following a fine residual aberration compensation. For all the results hereafter, the tip and tilt aberrations have been removed. While ZELDA is highly sensitive to these modes, in most AO systems they are managed by a dedicated tip and tilt sensor and/or a dedicated mirror that runs at a very high frequency. Furthermore, MITHiC is affected by a small amount of thermal turbulence inside its enclosure, which induces tip and tilt variations of typically 2-3 nm RMS in the 5-10 Hz frequency range. These variations appear at a higher frequency than our correction loop (currently ∼1 Hz), so we do not address them directly in the present work.

Static aberrations
We addressed the correction of the quasi-static aberrations on MITHiC, which are illustrated in Fig. 5 (top left). The standard deviation in the OPD inside the pupil, σ OPD , is typically 35 nm RMS at the beginning of the correction step. With the application of the closed-loop correction procedure, including the geometrical distortion correction and for various values of w, we obtained the results shown in Fig. 5. The top plot shows the OPD maps from ZELDA at three representative iterations of the procedure with w =115 c/p. The bottom plot displays the evolution of σ OPD at different iterations.
The general behavior is identical for most values of w, with a steady decrease in σ OPD during the first two or three iterations followed by a stabilization during the following iterations. This is due to the limited linearity range of ZELDA, which is reduced compared to other classical WFSs, such as the Shack-Hartmann WFS. As highlighted in N' Diaye et al. (2013), the linear behavior of ZELDA is typically limited to ±0.03λ (±20 nm here), and larger absolute values of phase error will always be underestimated by the second-order reconstruction. As long as there is no phase wrapping, the estimation will always favor a decrease in the wavefront error, which means that a correction in a closed loop will eventually converge. The same behavior is used for co-phasing techniques proposed in multi-wavelength approaches with the ZWFS (Vigan et al. 2011;Cheffot et al. 2020). Ultimately, the usable capture range of ZELDA in a closed loop depends on the phase shift imposed by the mask and ranges between −0.14 λ and +0.36 λ for the MITHiC prototype (N'Diaye et al. 2016).
For w between 75 and 137 c/p, σ OPD stabilizes around ∼1.9 nm RMS, which is above the theoretical limit of 0.76 nm RMS. Most of the residual aberrations can be at-tributed to low-order aberrations from the internal turbulence on MITHiC, vibrations of the bench, and dust particles on optical surfaces that create small errors in the ZELDA reconstruction. For w =20 c/p, the convergence is slower than for higher filtering values and it does not quite reach this low value of 1.9 nm RMS, which indicates that there is a non-negligible amount of aberrations with spatial frequencies between 20 and 75 c/p. The w =5 c/p case is even more extreme, with an even slower decrease in σ OPD and a convergence toward a larger value of ∼14 nm RMS after 15 iterations.
The w =200 c/p and the no-filter curves show the same steady decrease during the first iterations, showing that the filter can be optional for up to five iterations, for example to record a flat position to perform other science later. However, the loop is no longer stable and diverges after about seven iterations when no filter is applied and after 15-20 iterations for w =200 c/p. We note that, even though the w =200 c/p filter maintains all the spatial frequencies measured by ZELDA, the attenuation of the highest frequencies by the apodization of the filtering window allows for a significant stabilization of the loop. We attribute this effect to a misestimation or mis-correction of the highest spatial frequencies -which originate from pupil edge effects aliasing since the correction maximum frequency is at 137 c/p -as well as higher-order distortion residuals since we only calibrate for linear transformations.
To emphasize the need for the pupil distortion calibration detailed in Sect. 3.3, we present the convergence results without the distortion correction in Fig. 6. In these tests, the measured OPD map is simply filtered, scaled down, and applied to the SLM. For the two smallest w (5 and 20 c/p), the results are strictly identical to the previous test. For larger w, however, the behavior is very different: After reaching a minimum around 5 nm RMS, the standard deviation quickly diverges and converges to values larger than 20 nm RMS after ten iterations.
The divergence is caused by the apparition of high-spatialfrequency wrinkles that slowly degrade the wavefront quality. The issue remains insignificant for w =20 c/p, corresponding to a VLT/SPHERE-like case, but becomes critical for w larger than 75 c/p. If the computation of an interaction matrix can be performed to get rid of these ripples, new adaptive ELTs might struggle to construct them due to the absence of an appropriate calibration source (Heritier et al. 2018). Our calibration therefore shows another way of sensing the distortion modes and improving the correction with on-sky measurements with a simple implementation as long as the phase errors are not too large. In systems with large phase errors, such as a system without a first step of AO correction, further strategies would have to be explored. For example, to measure the distortion, the position of the introduced Gaussian spots could be estimated through the turbulence by averaging turbulence with long exposures. Further studies would have to be performed to assess the precision of such measurements.
Based on the successful implementation of our wavefront correction procedure, we highlight one emblematic case of quasi-static aberrations where a very large number of degrees of freedom for wavefront correction would be particularly interesting for future applications: the compensation for residual co-phasing errors for telescopes with segmented primary mirrors. To test ZELDA in this practical case, we selected a pattern that combines typical NCPAs and a segmented pupil with random tip, tilt, and piston on individual segments on the MITHiC SILIOS phase screen. The results of the correction procedure with w =115 c/p are illustrated in Fig. 7  ZELDA to measure fine differential piston errors. The high 115 c/p cutoff value for w ensures a good correction for the segment discontinuities.
The use of ZELDA for fine segment co-phasing was proposed by Janin-Potiron et al. (2017b); in their approach, they measure the residual piston, tip, and tilt of individuals segments and then apply the correction directly on the segments using their control actuators. In our case, as we benefit from a large number of subapertures and actuators, we go one step further and correct not only for the piston, tip, and tilt, but also for higher-order phase errors such as low-order aberrations on individual segments. On a pupil with 30 segments across the diameter, such as the ELT, 120 actuators across the pupil diameter will yield roughly 4x4 actuators per segment, which would allow for subsegment aberration correction. The curve in Fig. 7 that is the closest to this case is the orange curve, simulating a case with 40 actuators for 11 segments. This sort of possibility could prove extremely valuable for exoplanet imagers on future ELTs.

Temporal error analysis with residual ExAO aberrations
When using ZELDA, it is also worth considering the compensation for turbulence residuals in real time, for example downstream an ExAO system that provides first-stage correction. Sim-ilarly to the co-phasing case, this approach could potentially provide a second-stage WFS to increase the ExAO correction level in extremely demanding cases.
The phase screen installed on MITHiC provides an external ring of engraved turbulence that simulates residuals from the VLT/SAXO system (Fusco et al. 2006) in median seeing conditions at the Paranal Observatory. This ring can be rotated at a controlled speed to emulate the wind speed. To overcome the bench limitation in terms of temporal frequency (∼1 Hz), we discretized the steps to simulate a VLT/SPHERE residual turbulence for an ExAO system, which, in practice, means that the screen rotates very slowly and in discrete steps that are synchronized with the ZELDA measurements. The delay associated with one-frame processing is therefore 0.67 ms.
To explore several simple, yet realistic, AO systems, we included a delay, d, expressed in numbers of frames. The turbulence starting point on the phase screen was identical for each experiment. At each time step, k, the following steps were performed: (i) acquire a phase measurement, ϕ k ; (ii) compute the corresponding correction, Ψ k , using the procedure from Sect. 3.2; (iii) move the phase screen by a step of δ, simulating the wind-induced motion of the atmosphere in a 1/1500 s interval on an eight-meter-size telescope (a lateral shift of a×δ is also applied, a denoting a random number in the range [−1 : 1]); and (iv) apply the correction, Ψ k−d , on the SLM with a gain of 0.5. Article number, page 7 of 14 A&A proofs: manuscript no. paper   8 presents σ OPD for the standard deviation as seen by ZELDA in the OPD as a function of iterations for three different cases. While the starting standard deviation in the OPD is between 21 and 23.5 nm RMS for all the configurations, there is a clear reduction in the wavefront error after fewer than ten iterations. However, the mean level reached for each curve that represents a filtered OPD is different. The highest threshold at 16.9 nm RMS is reached with the two-frame delay configuration with the lowest filter at 44 c/p (in red): The controllable modes are restricted by the low-pass filter. As a result, the contribution to wavefront error coming from uncontrolled high spatial frequencies, or fitting error, is higher. By increasing the filter frequency cutoff to 55 c/p (in blue) while maintaining the twoframe delay, we reduced this fitting error and therefore improved the final correction quality to 15.6 nm RMS. We finally considered the hypothetic case where the lag error could have been reduced to one frame (in green) with the filter at 44 c/p. This lag reduction provides a clear improvement in the correction quality, which reduces to 14.6 nm RMS, yielding a gain of 2.3 nm RMS, on average, with respect to its two-frame lag counterpart. The final threshold in all three configurations is noticeably higher than what was achievable in Sect. 4.1 due, in large part, to the fitting error of approximately 10 nm. Furthermore, we used a loop gain of 0.5, which will undeniably limit the ultimate performances in a dynamic case. We did not identify the exact origin of the peri-odic pattern shown in Fig. 8. Nevertheless, the fact it is identical for each realization when the same turbulence is used leads us to conclude that it is inherent to the phase screen rather than a result of our analysis. Figure 9 provides an additional description of the distribution of aberrations with spatial frequency in these three different cases after 40 ms of simulated time, which is equivalent to 60 iterations. The plateaus reached in Fig. 8 show that other iteration values would yield the same power density spectrum (PSD) behavior. The gain with respect to the initial phase step (in blue) is clear for the three PSD curves in the low-and mid-spatial frequencies, between 0 and 20 c/p. Between 10 and 40 c/p are the differences in terms of the OPD of Fig. 8: For the two-frame delay curves, the stronger apodization of the filter at 44 c/p translates into a higher PSD around the 20-40 c/p area. Except for the low frequencies that are crippled by the bench turbulence, the best correction is provided by the reduced lag simulation (in green). The peak close to 20 c/p corresponds to the cutoff frequency of the VLT/SAXO system that was simulated by the rotating phase screen.
In this experiment we performed two camera acquisitions for each phase measurement: one without the ZELDA mask to obtain a clear pupil image and one with the mask. This approach is not really suitable for an on-sky system as moving a mask at this cadence would be difficult. Furthermore, it would need to  work solely based on the ZELDA image to maximize efficiency. Relying on a clear pupil measurement made at the beginning for the rest of the sequence has been investigated on MITHiC; unfortunately, it causes the correction to diverge, due to Fresnel propagation effects that modify the pupil intensity distribution, and biases the phase reconstruction. Based on the work of Vigan et al. (2019), we expect this effect to not be present in VLT/SPHERE data. In Fig. C.1 of that paper, the authors show the on-sky clear pupil image that they use for their correction sequence, which shows no sign of the turbulent phase. Although it is not directly reported in their work, they confirm that the clear pupil image is highly stable over time (private communication). However, MITHiC has not been designed to minimize these amplitude effects, and the test bed is very sensitive to misalignment or to phase errors along the optical axis because we scale down the eight-meter-class telescope turbulence on the 5 mm test bed pupil.

Impact on coronagraphic images
The contrast performance of coronagraphs is crucially dependent on the quality of the incoming wavefront at the level of the FPM. Our bench calibration procedure allows for an accurate compensation for the MITHiC aberrations down to a few nanometer RMS. With such a high level of wavefront correction, MITHiC becomes diffraction-limited and should offer a significant attenuation of the diffraction with the CLC. While the CLC provides a moderate contrast (compared with the apodized-pupil Lyot coronagraph ;Soummer 2005), and a moderate inner working angle (compared with the vortex coronagraph; Mawet et al. 2005), it remains a user-friendly coronagraph for simple coronagraphic imaging in a monochromatic light system that works with a circular pupil, such as MITHiC.
We performed coronagraphic measurements in the initial configuration of MITHiC (i.e., without aberration compensation and after five iterations of the compensation procedure with w =125 c/p), which means that most of the controllable spatial frequencies with the SLM are corrected. Figure 10 shows the coronagraphic images before and after the aberration compensation. Due to the use of a beam splitter to send part of the beam to ZELDA, the left part of the images is dominated by relatively bright ghosts, which are circled in blue. The two symmetric ghosts circled in white on either side of the PSF are induced by the other beam splitter located in the collimated beam that hits the SLM. Ignoring the ghosts, the residual aberration correction results in a substantial improvement of the diffraction attenuation of the coronagraph, with a visible decrease in the amount of light scattered in speckles all over the field of view. Before calibration, the focal plane shows diffraction rings, Power spectrum (nm/cycles/pupil) 2 Iteration 0 1-frame lag; 44 c/p 2-frames lag; 44 c/p 2-frames lag; 55 c/p Fig. 9: PSD of the OPD maps at simulated times of 0 ms (cyan) and 40 ms for the two-frame delay with w=55 c/p (blue) and w=44 c/p (red) and the one-frame delay with w=44 c/p (green). The vertical black line represents the cutoff frequency at 55 c/p. speckles, and pinned speckles (Sivaramakrishnan et al. 2002) up to 15 λ/D. After calibration, the brightness of the speckles and residual diffraction appears largely attenuated, indicating an improvement in the CLC performance. Apart from the four ghost images, some irregular speckle patterns are visible at larger separations in the coronagraphic images. This indicates the presence of residual high-spatial-frequency aberrations that could originate from NCPAs between the ZELDA arm and the coronagraph arm.
The gain in contrast performance of the CLC is confirmed by the PSF and coronagraphic image profiles in Fig. 11, which show a clear distinction before and after residual aberration compensation. This is particularly visible for the coronagraphic profiles, with a mean gain of a factor of five between 3 and 15 λ/D.
We compare the aberration-corrected coronagraphic profile to a CLC simulation that takes a residual OPD map measured with ZELDA prior to the coronagraphic measurement as input. Between 3 and 7 λ/D, our measurements are clearly at the theoretical limit of the CLC implemented on MITHiC. At larger separations, there is a slow departure between the theoretical and experimental profiles, which can be explained by some residual aberrations and by the readout noise floor of the camera above ∼12 λ/D.

Conclusion
The wavefront quality is a driving parameter in the performance of high-contrast imaging systems that use coronagraphy to attenuate telescope diffraction on an observed star image. ZELDA has already been demonstrated to be a powerful secondary WFS that can calibrate aberrations unseen by an ExAO system.
In the present work, we demonstrated the ability of ZELDA to work in a closed loop with an SLM to perform a second-stage AO correction, controlling a very large number of degrees of freedom. Thanks to the use of a high-density SLM, the MITHiC test bed offers tens of thousands of DM-like actuators for wavefront correction. While these numbers appear large by today's standards, they are in line with the expectations for high-contrast imaging of Earth-like planets with future ground-based ELTs. Rather than building a large interaction matrix, we instead implemented a direct wavefront correction based on the OPD maps obtained with ZELDA and complemented by a simple and fast distortion calibration procedure. Our results demonstrate that this procedure is very efficient at compensating for the static phase patterns of the bench and remaining stable once a minimum has been reached. With the current procedure, the correction for the wavefront errors in a static case reaches a level of ∼2.6 nm RMS.
The pupil distortion calibration is a crucial step. While such effects are naturally taken into account in systems relying on a classical interaction matrix, this is not the case here. We instead  Fig. 10: Coronagraphic images obtained before (left) and after (right) compensation of the bench aberrations. The compensation is performed with a low-pass filtering value of 125 c/p. The CLC FPM has a radius of 2.7 λ/D, the footprint of which is circled in black. The left part of the images is contaminated by several ghosts, circled in blue, which were introduced by the beam splitter located between the SLM and the coronagraphic FPM, which is used to send part of the beam to ZELDA. The speckles circled in white are created by ghost reflections in the cubic beam splitter in front of the SLM.  : PSF and coronagraphic profiles before (black) and after (red) compensation of the bench aberrations. The profiles for the PSF and coronagraphic data are obtained with the azimuthal average and standard deviation of the images, respectively. The profiles are compared to the theoretical simulated profile (blue) obtained using as input an OPD map measurement with the bench aberration compensation prior to the coronagraphic measurement. The ghosts visible in Fig. 10 are masked when computing the profiles. The grayed out part of the plot corresponds to the area masked by the FPM in the coronagraphic images.
used a two-frame calibration where Gaussian reference spots are introduced on the SLM. Based on their known input position, we were able to easily calibrate the distortions in the system and subsequently correct for them when applying the phase corrections on the SLM. Without this calibration, the wavefront errors quickly diverge in the phase compensation procedure due to the apparition of high-spatial-frequency ripples induced by the distortion.
We also demonstrate a dynamic "fast-ZELDA" mode with a more dynamic case where the sensor is used to compensate for the residuals of a VLT/SPHERE-like ExAO system. This setup acts as a second-stage correction that could possibly be used in very demanding science cases. For this demonstration, we used the MITHiC phase screen, which allows for the simulation of VLT/SPHERE-like turbulence residuals. We simulated a system running at 1.5 kHz with a delay of one or two frames in the correction of the residuals. The results show an efficient wavefront correction that is stable over a thousand iterations. While the simulated residual aberrations present a wavefront error of around 20-25 nm RMS, the closed-loop correction, dominated by the uncontrollable fitting errors and the bench internal turbulence, reduces this to approximately 15 nm RMS.
Our setup on MITHiC presents some shortcomings. Firstly, we work with an SLM in monochromatic light. From the WFS point of view, similar results could certainly be obtained with a reasonable broadband since ZELDA has been shown to accommodate bandwidths of 20% well (N'Diaye et al. 2013(N'Diaye et al. , 2016. However, the SLM is chromatic and cannot be used in a broadband setup. Broadband demonstration would require the implementation of a classical high-density DM. Secondly, our measurement procedure requires two camera acquisitions for one OPD measurement: with and without the ZELDA mask. This would be an issue when running a system at several kilohertz. Further studies need to be performed on the stability of the amplitude in the pupil, as well as on whether or not taking the acquisition without the mask is mandatory. Nevertheless, we could foresee a fast synchronized tip and tilt mirror that would temporally switch the beam from one configuration to the other or a setup with a beam splitter to perform the two simultaneous