Increasing the raw contrast of VLT/SPHERE with the dark-hole technique II. On-sky wavefront correction and coherent differential imaging

Context. Direct imaging of exoplanets takes advantage of state-of-the-art adaptive optics (AO) systems, coronagraphy, and post-processing techniques. Coronagraphs attenuate starlight to mitigate the unfavorable ﬂux ratio between an exoplanet and its host star. AO systems provide di ﬀ raction-limited images of point sources and minimize optical aberrations that would cause starlight to leak through coronagraphs. Post-processing techniques then estimate and remove residual stellar speckles due to hardware limitations, such as noncommon path aberrations (NCPAs) and di ﬀ raction from telescope obscurations, and identify potential companions. Aims. We aim to demonstrate an e ﬃ cient method to minimize the speckle intensity due to NCPAs and the underlying stellar di ﬀ raction pattern during an observing night on the Spectro-Polarimetric High-contrast Expolanet REsearch (SPHERE) instrument at the Very Large Telescope (VLT) instrument without any hardware modiﬁcations. Methods. We implement an iterative dark-hole (DH) algorithm to remove stellar speckles on-sky before a science observation. It uses a pair-wise probing estimator and a controller based on electric ﬁeld conjugation, originally developed for space-based application. This work presents the ﬁrst such on-sky minimization of speckles with a DH technique on SPHERE. Results. We show the standard deviation of the normalized intensity in the raw images is reduced by a factor of up to ﬁve in the corrected region with respect to the current calibration strategy under median conditions for VLT. This level of contrast performance obtained with only 1 min of exposure time reaches median performances on SPHERE that use post-processing methods requiring ∼ 1h-long sequences of observations. The resulting raw contrast improvement provides access to potentially fainter and lower-mass exoplanets closer to their host stars. We also present an alternative a posteriori calibration method that takes advantage of the starlight coherence and improves the post-processed contrast levels rms by a factor of about three with respect to the raw images. Conclusions. This on-sky demonstration represents a decisive milestone for the future design, development, and observing strategy of the next generation of ground-based exoplanet imagers for 10-m to 40-m telescopes.


Introduction
The imaging of faint stellar companions, exoplanets, and circumstellar disks with ground-based high-contrast imaging instruments, such as the Spectro-Polarimetric High-contrast Expolanet REsearch (SPHERE) instrument at the Very Large Telescope (VLT; Beuzit et al. 2019) and the Gemini Planet Imager (GPI; Macintosh et al. 2014), is currently limited by wavefront aberrations and diffraction effects. Extreme adaptive optics (XAO) systems control dynamic optical perturbations caused by atmospheric turbulence and compensate for static phase aberrations introduced by manufacturing defects and misalignments of optical components from the telescope primary mirror to the XAO wavefront sensor (WFS). Classical AO systems split the beam and analyze the wavefront aberrations in a dedicated optical channel. This two-channel configuration (one for the ⋆ Based on observations collected at the European Southern Observatory, Chile, 108.22NJ. measurement of the aberrations and one for the science detector) cannot constrain the aberrations introduced by the optical components in the science channel. These noncommon path aberrations (NCPAs, Fusco et al. 2006) induce additional starlight leakage that interferes with the underlying stellar diffraction pattern and takes the form of static and quasi-static speckles on the science detector, currently preventing the detection of sources fainter than one millionth of their host star's flux.
The mitigation of static speckles caused by NCPAs has been extensively studied by the community. First, postprocessing techniques such as angular (ADI, Marois et al. 2006), spectral (SDI, Racine et al. 1999), polarimetric (PDI, Kuhn et al. 2001), and reference-star (RDI, Lafrenière et al. 2009) differential imaging are broadly employed to enhance detection levels in the images. Both ADI and SDI, routinely used for the detection of exoplanets, suffer from self-subtraction of the off-axis sources at small angular separations (Esposito et al. 2014). ADI requires additional observing time to increase the total parallactic angle rotation, limiting the total observing efficiency. To enhance the A&A 665, A136 (2022) contrast in the images before post-processing, the NCPAs are also currently compensated, prior to the observing night, with phase retrieval techniques (Gonsalves 1982;Paxman et al. 1992;Sauvage et al. 2007). Methods to compensate for the NCPAs with the XAO high-order deformable mirror (HODM) in parallel with the science observations is an active field of research. Science cameras running at high rates allow for techniques that compare WFS telemetry with focal plane images to characterize the optical system and measure the NCPAs (Rodack et al. 2021;Skaf et al. 2022). But, the long read-out time of science detectors prevents the use of such methods with most of the high-contrast imagers. Other techniques based on focal plane WFSs, initially developed for space-based applications, have also been tested in dynamic conditions and implemented on different groundbased instruments. For instance, the self-coherent camera (SCC, Baudoz et al. 2006) has been validated on-sky to correct for loworder aberrations . Higher-order corrections have been demonstrated on-sky using speckle nulling techniques (Bordé & Traub 2006;Martinache et al. 2014;Bottom et al. 2016), as well as the Zernike sensor for Extremely Low-level Differential Aberration (ZELDA, N'Diaye et al. 2016;Vigan et al. 2019). The coronagraphic phase diversity (COFFEE), as well as pair-wise probing (PWP), associated with the electric field conjugation (EFC) algorithm (Give'On et al. 2007) were previously only tested in the laboratory with a rotating phase plate simulating turbulence residuals (Potier et al. 2019;Herscovici-Schiller et al. 2019) and to calibrate state-of-the-art instruments using their internal sources (Paul et al. 2014;Matthews et al. 2017;Potier et al. 2020b).
The combination of PWP with EFC aims to minimize the speckle intensity in parts of the science detector called dark holes (DHs). PWP estimates the complex focal-plane electric field (E-field) through a temporal modulation of the intensity by introducing phase-shifted surface shapes to the HODM, resulting in an interferometric measurement in the image plane. Once estimated, the static E-field is minimized by the EFC algorithm, which determines the HODM settings that destructively interfere with the measured E-field in the DH region. As shown in Potier et al. (2020b), the technique applied to half of the field of view (commonly called a half DH) enables the compensation of both the static phase and amplitude aberrations, as well as the underlying stellar diffraction pattern, thus improving the contrast in the raw images (i.e., the intensity of the stellar residuals normalized by the off-axis point spread function of the same star). While this method has been shown to be suitable for correcting aberrations with a stable wavefront at extremely high-contrast ratios (Trauger & Traub 2007;Seo et al. 2019;Potier et al. 2020a), its application is more challenging during on-sky observations of stars with ground-based telescopes due to the varying conditions between the different images required by PWP. This is likely why the loop of the iterative PWP+EFC technique has never been successfully closed on-sky previously. However, we recently demonstrated in the lab that leaving the dynamic speckles to average in a halo through long-exposure imaging should overcome previous limitations and enable the correction of speckles that evolve slower than the servo-loop speed .
In this paper, we describe the on-sky results obtained with SPHERE where we successfully minimized the speckle intensity with PWP+EFC. In Sect. 2, we describe the direct on-sky compensation of the NCPAs in closed loop. In Sect. 3, we demonstrate that the PWP estimations can be used to implement an alternative coherent differential imaging (CDI) post-processing method that improves the exoplanet detection limit in the images, even on the uncorrected parts of the science detector.

Method
The observations were performed during the second half of the night of February 15, 2022 starting at 5 a.m. UTC. The seeing ranged from 0.56 ′′ to 1.2 ′′ at 500 nm according to the Paranal DIMM-Seeing monitor, slowly degrading over time with an average of ∼0.8 ′′ . This value is typical of median conditions at Paranal. The coherence time remained stable, varying between 8 and 10 ms during the full observing sequence. Wind originated from between the north-east and the north-west, with a speed below 5 m s −1 throughout the observations. These good conditions allowed us to use the small pinhole for spatial filtering in front on the Shack-Hartmann (SH) WFS (Poyneer & Macintosh 2004). We observed HIP 57013 (RA = 11 41 19.79, δ = −43 05 44.40, R = 5.5, H = 5.5), a bright early-type A0V member of the Sco-Cen association and part of the BEAST survey (Janson et al. 2021), using the H3 filter (λ 0 = 1667 nm, ∆λ = 54 nm). We attenuated the starlight with an Apodized Pupil Lyot Coronagraph (APLC, Soummer 2005), whose focal plane occultor is 185 mas in diameter.
Starting from the original SPHERE wavefront calibration based on phase retrieval, we ran the PWP+EFC focal plane wavefront sensing and control algorithms, following the method described by Potier et al. (2020b). The pair-wise probes and EFC solutions were applied on the HODM by changing the SH WFS reference slopes while the AO system ran in closed loop. The pair-wise probes consisted in the vertical displacement of one individual actuator located within the Lyot stop opening, and poked with a peak-to-valley amplitude of ±400 nm. We used two probes per iteration, which were distributed vertically to create a top DH as described in Potier et al. (2020b). We used an exposure time of 64 s, which averages the intensity variations due to atmospheric turbulence. Assuming the atmospheric turbulence was in steady-state during one iteration and small AO residuals with respect to the wavelength , the turbulent halo was ignored by PWP when the difference between each pair of diversity images was calculated.
The corrected region (DH) is an annulus with an inner and outer edge of 5 λ 0 /D and 15 λ 0 /D, respectively, which corresponds to angular separations ranging from 220 mas to 650 mas. This choice aims to minimize the risk of the loop being perturbed by the starlight leakage due to AO aliasing ( Fig. 6 of Cantalloube et al. 2019) or by the wind driven halo (Cantalloube et al. 2020). At each iteration, the EFC solution is filtered through a Tikhonov regularization (Pueyo et al. 2009), which suppresses modes higher than the 400th mode out of 1377 total modes controlled by SPHERE. A low servo-loop integral gain of 0.1 was used to account for biases in the PWP, especially those induced by small unexpected fluctuations in the intensity (mainly originating from the wind driven halo) during the acquisition of the probe images and errors in the model of light propagation. Given the apparent stability of the loop, we increased the gain to 0.2 from the third iteration to expedite the correction. Each probe image acquisition required 80 s (64 s exposure and ∼16 s overhead, including readout time) for a total of ∼400 s per iteration (four images required for PWP and one image to measure current raw contrast in an unprobed image). The total time required for the on-sky NCPA calibration could be reduced in the future, either by increasing the servo-loop gain, by decreasing the exposure time for the probed images, or by skipping the unprobed image that we used to verify loop convergence. This optimization will be performed before the official commissioning of this technique on SPHERE. Raw image with bad pixels removed. Right column: Same as for the left column, but high-pass filtered using a Gaussian kernel with a standard deviation of 0.57λ 0 /D. The images are normalized by the maximum of the off-axis PSF, and multiplied by 10 5 and 10 6 for the left and right images, respectively. The corrected region is encircled in green. Figure 1 shows the on-sky coronagraphic images (left) using the current SPHERE calibration procedure (top) and the closedloop DH correction inside the green line after five iterations (bottom). Images are normalized by the maximum of the offaxis point spread function (PSF) measured on the same star (i.e., the no-coronagraph PSF) to study the results in the "normalized intensity" dimensionless metric. Normalized intensity differs from raw contrast in that it does not account for potential spatial variations in the off-axis PSF transmission and morphology. The quasi-static speckles that set the detection limit for current SPHERE observations are effectively removed from the "DH correction" image, leaving only the fainter, smooth AO halo. The final raw image is similar to the simulated bottom left image of Fig. 4 in Potier et al. (2020b) where the halo of the SPHERE AO system (SAXO) limits the overall performance after compensating for the static speckles in the DH. The raw image in Fig. 1 is further limited by the fast-evolving low-order aberrations, whose level was underestimated in the former simulation.

Results
To reduce the effect of this halo, and highlight the impact of static and quasi-static speckles in the SPHERE performance, the images are high-pass filtered using a Gaussian kernel with a standard deviation of 2 pixels (0.57λ 0 /D, right column of Fig. 1). In the DH, the visible static speckles have been corrected, minimizing the estimated mean coherent intensity from ∼ 6 × 10 −6 to ∼ 1 × 10 −6 . The remaining photo-electrons in the smooth regions in the high-pass-filtered DH after an exposure time of 64 s with the science detector mostly come from the photon noise induced by the turbulent halo. We show the mean modulated and unmodulated intensities in the DH at each iteration in Fig. 2. The modulated component is the absolute value of the E-field squared estimated by PWP, while the unmodulated component represents the difference between the total recorded intensity and the modulated component. The unmodulated component therefore includes: (1) astrophysical sources located in the field of view; (2) calibration errors in PWP, such as instrument model inaccuracies; (3) the field of speckles whose life time is smaller than the measurement rate of PWP (∼5 min), such as the turbulent halo; and (4) chromatic residuals. Figure 2 demonstrates the convergence of the PWP+EFC closed-loop algorithm because the averaged modulated component in the DH is minimized. However, it also shows that the averaged level of starlight residual in the total image remains almost constant at ∼3 × 10 −5 throughout the process. Indeed, even though the processed high-contrast imaging data are limited by the effects of static and quasi-static aberrations at small separations, the smooth halo of turbulence is often the dominating term in the raw images. This halo is removed in post-processing through the filtering of the low spatial frequency residuals or with differential imaging, leaving behind only its associated photon noise. We also show the radial profile of the filtered images (calculated as the standard deviation in azimuthal rings of size λ/2D in the DH regions versus the angular separation) throughout the different PWP+EFC iterations in Fig. 3. Detection performance in the DH region has been improved from a factor of two to a factor of five in that area. We obtain the best improvement factor (more than three) between 245 and 500 mas. The normalized intensity reaches a level below 10 −6 between 410 and 635 mas. Although we did not try to minimize the starlight below 220 mas, we noticed an improvement between 150 and 220 mas, which we attribute to changing observing conditions (turbulence, wind driven halo, quasi-statics) as iterations progress.

Method
PWP can also be used alone to further suppress stellar speckles in the images in post-processing. Indeed, the stellar E-field estimated by PWP can serve as a reference image to be subtracted from the remaining total intensity, thereby improving the S/N of the astrophysical signal that is not coherent with the starlight. This technique called coherent differential imaging (CDI) has been demonstrated with many focal plane wavefront sensors (Baudoz et al. 2006;Bottom et al. 2017;Jovanovic et al. 2018). CDI can be used to further remove residual speckles inside the DH. In our case, the loop has converged and no more speckles are to be calibrated in post-processing.  In this section, we demonstrate the ability of CDI to enhance the contrast in regions that are not corrected during the PWP+EFC process. We also propose an optimization of the CDI technique to account for the changing atmospheric transmission over time between the recording of the total intensity images and the probe images. This simple algorithm will be investigated further in the near future to quantify potential over-subtraction of the astrophysical signal. At each iteration, the PWP estimation of the focal plane E-field (or reference image) is numerically recentered with respect to the total intensity signal to account for any drift on the science detector during the probing steps. It is then rescaled to minimize the mean intensity of the high-pass-filtered CDI result, in the region symmetrically opposed to the top DH described in Sect. 2. This area is arbitrarily chosen for the sake of illustration, while any region of the science detector well sensed by PWP could be selected. This optimization assumes that most of the total signal in that area is caused by stellar residuals and is, therefore, coherent. This optimization should not be strongly affected by a planet or other astrophysical object present in the field of view but injection of fake signal in the data is out of the scope of this paper. The estimated mean scaling factor is 1.0, which demonstrates that the calibration of the probes was accurate. The standard deviation of the scaling factor over the iterations is 0.1, which we attribute to varying atmospheric transmissions. The scaling method, reprocessed at each iteration, will be used in the PWP+EFC closed loop in the future for a faster convergence of the DH correction. The DH region is encircled in green. All the images were then filtered using a high-pass Gaussian function whose standard deviation is 0.57λ/D, to remove the halo of turbulence, and multiplied by 10 6 .

Results
The high-pass-filtered total intensity, the optimized reference image, and the CDI result at each iteration are shown in Fig. 4. Qualitatively, the static speckle intensity in the CDI results is minimized with respect to the total intensity in the region sensed by PWP. The DH geometry is defined either by the regions of the science detector reachable by PWP with the chosen probes, or by the spatial frequencies properly corrected by the AO system. For instance, the CDI performance is limited in the central region along the horizontal axis for the former reason (see Fig. 3 in Potier et al. 2020b) and in the regions above the AO cutoff frequency or where aliasing caused by the SH WFS dominates for the latter. Also, the result of CDI is limited by bright residual speckles at small separations. We conclude that PWP is not as sensitive to these residuals, either because they are caused by large wavefront errors that violate the linear approximation made by PWP (Give'On et al. 2007), or because they are too close to the fast-evolving low-order aberrations that violate the steady-state assumption required by PWP. For the former case, in principle the iterative PWP+EFC algorithm should be able to correct such aberrations during the observation. Therefore, better performance of the CDI algorithm is expected after a few iterations of DH correction in the targeted region. This approach is well demonstrated by the CDI result at iterations three and four, where the speckles in the top DH are almost perfectly removed. Figure 5 shows the radial profile of both the total intensity and the result of CDI in the optimized region (encircled in blue in Fig. 4). The CDI technique enhances the contrast level in the post-processed images between 150 and 640 mas for each estimation of the coherent signal with PWP. With this promising result, we envision observing strategies where PWP is applied at a regular cadence and the resulting CDI results are used as inputs for post-processing methods, or for a direct analysis. However, the contrast after applying CDI degrades with PWP+EFC iterations. Indeed, the minimization of the speckle intensity due to amplitude aberrations in the top DH with one HODM is achieved by degrading the contrast in the bottom part of the field of view. CDI then starts from a higher speckle intensity, with more noise, in the bottom part of the image. Hence, applying CDI before any half DH correction would be more efficient on the "degrade" part of the image. The implementation of CDI should, therefore, be defined depending on the science objective and telescope time. Ideally, either CDI would be applied on top of a DH activewavefront correction for science analysis in that region, requiring the implementation of PWP+EFC successively in different DHs to cover a broader region; or CDI would be directly applied in a larger area by applying PWP alone, without EFC correction, at a regular cadence. The former case would provide the best contrast performance, while the latter would reduce the required telescope time.

Conclusions
We present two promising implementations of a dark hole (DH) technique on VLT/SPHERE. First, relying on an efficient adaptive optics (AO) system combined with good observing conditions, we demonstrated a pair-wise probing (PWP) + electric field conjugation (EFC) algorithm directly on-sky that minimized the static speckle field intensity in a defined region using VLT/SPHERE. The standard deviation of the normalized intensity was improved by a factor of up to five in the DH. Therefore, 1min-exposure images, after the PWP+EFC correction is applied (30 min here), achieve the same level of performance as ADI  with 1h-sequences of observations. This promising technique is expected to be even more robust with the improvement of the AO systems planned for the next generation of high-contrast imaging instruments (Chilcote et al. 2018;Boccaletti et al. 2020). Second, we showed that estimation of the coherent light in the field of view with PWP could be used to further calibrate the static speckles in post-processing through coherent differential imaging (CDI). PWP could also be used at a regular cadence during a high-contrast observation to estimate the spatial distribution of the coherent speckles and then subtract them from the raw images. The residuals, containing the incoherent signal of interest, could then be further analyzed with conventional methods such as angular (ADI), spectral (SDI), or reference (RDI) differential imaging.
Here, the performance was quantified via a simple highpass-filtered image. Therefore, the CDI method also solves the self-subtraction issue that emerges when performing ADI at small angular separations with a coronagraph instrument, such as those that will equip the upcoming extremely large telescopes (Baudoz et al. 2014;Fitzgerald et al. 2019;Kasper et al. 2021). We will assess this capability in future experiments, where the technique will target the inner edges of known circumstellar disks or face-on disks whose morphologies would be altered by ADI and SDI post-processing techniques. We will also combine this method with star-hopping algorithms (Wahhaj et al. 2021) in the near future to understand the speckle lifetime and how the A136, page 6 of 7 A. Potier et al.: Increasing the raw contrast of VLT/SPHERE with the dark hole technique. II. DH correction is impacted by the telescope alignment from one target to another. All of these experiments aim to determine the most optimized strategy of observation for the search for and the characterization of exoplanets in direct imaging.