Predictive wavefront control on Keck II adaptive optics bench: on-sky coronagraphic results

The behavior of an adaptive optics (AO) system for ground-based high contrast imaging (HCI) dictates the achievable contrast of the instrument. In conditions where the coherence time of the atmosphere is short compared to the speed of the AO system, the servo-lag error can become the dominant error term of the AO system. While the AO system measures the wavefront error and subsequently applies a correction (typically taking a total of one or a few milliseconds), the atmospheric turbulence above the telescope has changed resulting in the servo-lag error. In addition to reducing the Strehl ratio, the servo-lag error causes a build-up of speckles along the direction of the dominant wind vector in the coronagraphic image, severely limiting the contrast at small angular separations. One strategy to mitigate this problem is to predict the evolution of the turbulence over the delay time. Our predictive wavefront control algorithm minimizes, in a mean square sense, the wavefront error over the delay and has been implemented on the Keck II AO bench. In this paper, we report on the latest results of our algorithm and discuss updates to the algorithm itself. We explore how to tune various filter parameters based on both daytime laboratory tests and on-sky tests. We show a reduction in the residual-mean-square wavefront error for the predictor compared to the leaky integrator (the standard controller for Keck) implemented on Keck for three separate nights. Finally, we present contrast improvements for daytime and on-sky tests for the first time. Using the L-band vortex coronagraph for Keck's NIRC2 instrument, we find a contrast gain of up to 2 at a separation of 3 lambda/D and up to 3 for larger separations (3-7 lambda/D).


INTRODUCTION
Single conjugate adaptive optics (AO) systems provide a good correction over a modest field-of-view, thereby enabling a vast array of science. Within a high contrast imaging (HCI) system, this simplest form of AO becomes a powerful tool when observing bright, on-axis guide stars. By delivering an extreme AO correction through high speed operation and maximizing the actuator density of the deformable mirror (DM), HCI systems can search for exoplanets and disks around bright stars. The resulting AO systems provide Strehl ratios (SR) of greater than 60% in near-infrared (NIR) wavelengths. 1,2 With these advancements in AO, HCI systems have reached post-processed contrasts of 10 −6 at spatial separations of 200 milli-arc-seconds. 3 Closer to the host star, a typical AO system is dominated by either noncommon path aberrations or the servo-lag error. The servo-lag error, in the worst-case scenario, can visually be seen as a wind-driven halo in a coronagraphic image for instruments such as the Very Large Telescope's (VLT) SPHERE HCI instrument, contaminating the region near the inner working angle of the coronagraph. 4 For other HCI systems, the servo-lag error is a more subtle build-up of speckles. The servo-lag error is due to the atmospheric turbulence above the telescope evolving in the time between when the wavefront is measured and the new DM commands are applied. The delay time includes half the wavefront sensor (WFS) integration time, the wavefront sensor camera read-out time, computation time, and the DM response time; it ranges from 1-4 wavefront sensor frames, depending on the system. One solution to minimize the servo-lag error is to predict the evolution of the wavefront over the known delay using predictive control.
Over the last five years, renewed efforts to implement predictive control on-sky for HCI have resulted in algorithms such as empirical orthogonal functions (EOF) 5 which is being tested at SCExAO 6, 7 and W. M. Keck Observatory . 8 Van Kooten et al. 2017, 2019 9, 10 proposed a recursive solution to the minimum mean square cost function and within their framework test prediction for non-stationary turbulence as well as VLT/SPHERE telemetry. 11 Most recently, data-driven subspace predictive control 12 has been proposed to operate in closedloop, with plans to test on-sky with the MagAOX system. 13 This work differs from previous work in the early 2000s focused on optimal control algorithms such as Linear-Quadratic-Gaussian control (LQG), [14][15][16] H2 optimal control, 17 and Predictive Fourier control 18 that all rely on a combination of modeling and data-driven methods to identify input disturbances (both from telescope vibrations and turbulence) to predict the disturbances. Some of these methods have been successfully commissioned for tip-tilt control including on SPHERE/Very Large Telescope, 19 and GPI/Gemini Telescope 20 and are being tested for other low order modes on-sky. 21 Finally, efforts to implement machine learning for predictive control are underway various groups. Early results from CANARY showed machine learning algorithms on-sky as a reconstructor. 22 More recent simulation work show promising results when comparing EOF to a neural network under certain assumptions and advocate for operating non-linear wavefront sensors. 23 However, more work into necessary to understand the require training for on-sky deployment and the required hardware for real-time-control. These are also discussed by J. Nousianen et al (2021) 24 and R. Swanson et al. (2021) 25 while showing results for other machine learning algorithms.
In this paper, we focus on W. M. Keck Observatory's HCI system, Keck II's near-infrared camera (NIRC2).
In Xuan et al. 2018, 26 the authors found a correlation between the final post-processed contrast for NIRC2 and the coherence time of the atmosphere divided by the WFS integration time when using the Shack-Hartmann wavefront sensor (SHWFS); both variables influence the impact of the servo-lag error. These results show that predictive control would improve the performance of the system. Predictive control was implemented on the Keck II AO bench and Jensen-Clem et al. 2019 8 show a factor of 2.2 reduction in the residual-mean-square (RMS) error during daytime tests with simulated turbulence projected onto the DM. Since then, updates to the AO system have been made along with the predictive controller. In this paper, we present the latest results from daytime tests and on-sky engineering time at Keck, including a new recursive implementation tested during the day and contrast measurements from both day and night time tests using EOF. We evaluate the performance of the new controller by looking at the residual wavefront error but also the contrast curve for NIRC2 specifically looking for improved contrast at small angular separations where the servo-lag error is expected to have the most impact.
We outline the paper as follows: an overview of Keck II AO and the NIRC2 instrument is given in Sec. 2.
In Sec. 3, we give details of the controller implemented on Keck including the predictive control (Sec. 3.2). We present results from daytime and night time tests in Secs. 4 and 5, respectively. Discussion of the results are given in Sec. 6. We outline our plans for the system in Sec. 7. Finally, we end with our summary and conclusions in Sec. 8.

KECK II ADAPTIVE OPTICS SYSTEM AND NIRC2
NIRC2 (PI: Keith Matthews) contains two vortex coronagraphs to cover the full wavelength range of the instrument (K-, and L'-/M-band vortex, 27 respectively) to enable HCI. NIRC2 is fed by the Keck II AO bench that operates with either a visible SHWFS 28 or a near-infrared pyramid wavefront sensor (PyWFS). 1,29 The PyWFS enables AO correction for much fainter, redder guide stars, supporting unique science capabilities. For the remainder of this paper we will be referring to the PyWFS configured AO system.

Integral control law
Operating at 1kHz with 21 actuator across the aperture in H-band allows for the PyWFS to achieve good performance with the standard AO control law: a leaky integrator. In Eq. 1, the leaky integrator calculates the DM commands, V int (t), at time step t.
The leak, k, is nominally set to 0.99 for standard operation on Keck and for our tests. The gain (g) has a typical value of 0.4 which is adjusted depending on seeing conditions. S(t) and S ref are the slopes measured by the PyWFS at time t and the reference slopes encoding the non-common path aberrations for NIRC2, respectively.

Predictive Control
Using the separation principle, our predictive controller is implemented in two steps: 1. predict the input disturbance (i.e. wavefront error due to atmospheric turbulence above the telescope) and 2. control the system plant (i.e., DM). Within this flexible framework, we can implement different predictive algorithms. We are also able to pick various control laws including open-loop control or the already existing leaky integrator implemented for the PyWFS.
In this work, we predict the DM commands one lag time into the future by finding a filter that uses a subset of previous measurements subject to some cost function. We focus on data-driven approaches where we use recent measurements to determine the future state of the input disturbance. With that we rely on spatial and temporal relationships of the data which not only allow for us to predict but also help reduce other error terms in the system that is present in the data. Since we aim to predict the full input disturbance (e.g., the wavefront due to turbulence), we must first reconstruct the pseudo-open loop (POL) wavefront as the PyWFS is downstream from the DM measuring the residual wavefront error.

Pseudo-open loop reconstruction
Making use of the shared memory structure within the CACAO RTC framework, the POL DM commands are calculated using the current PyWFS residual wavefronts and the previous DM commands from lag frames behind.
Only the modes controlled by the leaky integrator are reconstructed to form the PyWFS residuals in this step and therefore are reconstructed using the CM , the filtered command matrix. At any point, up to 120000 frames (2 minutes at 1 kHz) of filtered POL DM commands are available in the shared memory and can be saved on command. The POL DM commands can be converted to optical-path-difference (OPD) using 0.6 µm V . For calculating the predictive filter we make use of these filtered POL DM commands.

Predictive filter
There are a few different flavors of predictive control implemented on the PyWFS RTC that all aim to minimize the same cost function: the linear minimum mean square error (LMMSE). We write the LMMSE cost function in Eq. 2 where a unique predictive filter, F i is determined for each mode or in our case DM actuator.
The filtered POL DM commands we want to predict (i.e., true wavefronts) are represented by w i at δt steps into the future while w is our POL measured data with dt being the time step between wavefront measurements.
Our regressors are contained in the history vector h(t), which is formed from the POL DM commands (Eq. 3) containing temporal order, n, of previous measurements (i.e., w). Each measurement vector is of length m.
In 2019, EOF was implemented on the PyWFS RTC in python such that the predictive filter was calculated in a batch sense by storing up to l wavefront sensor measurements. 8 In this implementation, the cost function is solved in a regularized least-squares sense as shown in Eq. 4 where α is the regularization parameter and I is the identity matrix. P and D contain l pairs of the true wavefront (w i ) and history vectors, h, respectively.
The matrix multiplication of the predictive filter has now been implemented in the GPU on the PyWFS RTC, allowing us to predict at the speed of the system, an upgrade since 2019, 8 and therefore run on-sky. The EOF filter, using Eq. 4 is first calculated offline in python and then passed to the shared memory once prediction is enabled. At the moment, the filter takes 9 seconds to calculate using 1 minute of telemetry (l=60000 frames).
We have also updated the filter such that the number of temporal orders (how many previous measurements used to form h), n, can vary. Previously, we had been limited to a temporal order of 10 but now have the option to vary this parameter. Currently, however, doing so greatly increase the computation time to minute timescales.
More recently, the recursive LMMSE (RLMMSE) solution 10 has also been implemented on the Keck RTC within python. The RLMMSE allows for the system to continuously update the predictive filter, effectively tracking any changes in the input disturbance (atmospheric turbulence). The solution to the cost function, Eq. 2, can be written in terms of the auto-covariance of h and the cross-covariance between the w i and h (C hh and c hwi , respectively) as shown in Eq. 5. Using this equation, we can calculate the batch solution by once again storing up data and estimating the auto-covariance and the cross-covariance matrices. We can also build-up our covariances using Eqs. 6 and 7. In practice, the filter first runs in the background, once it has initially converged using a few thousand frames (e.g., l=O(1000)), the filter is then sent to the shared memory and prediction is enabled. The filter can then update recursively as fast as possible within the python implementation (each calculation takes roughly 0.02 seconds giving us an update rate of 50 Hz) at the moment.

Control Law
Our current implementation of the predictive filter is a semi-open loop controller. The predicted wavefronts return the DM commands for correcting the full atmospheric turbulence, The predicted commands can be mixed with the integrator to generate the signal sent to the DM, V (t). This configuration is shown in Eq. 8. The integral controller (the second term) is the same as Eq. 1. Now the slopes encode the residuals of the predictor. We tune the mixing parameter m and the gain g to have stable correction where m values between 0.1 and 0.75 provide stability with values closer to 0.5 providing the best correction. In the end, this mimics a leaky integrator operating on the predicted values and residuals that allows for us to control the system plant, combining with our predictive filter to become our predictive controller.

Residual wavefront error
To determine the performance of the AO system we look at the residual-mean-square (RMS) calculated from the residual (close-loop; CL) reconstructed wavefronts. Note we use the RMS as an intermediate performance metric and then use the measured contrast as our final metric. We make use of two different CL wavefront reconstructed values: filtered and full CL wavefront. For both cases, the CL values (i.e., residuals) are calculated in DM units of voltages and then a conversion factor of 0.6 µm V is applied to convert between DM voltages and optical-path-difference/wavefront.
The filtered CL DM commands are those reconstructed in real-time by the RTC from the residual/CL PyWFS slopes. Modal filtering is embedded in the CM and so the resulting residual is filtered depending on the applied modal gains and mode cut-off for the leaky integrator. In post-processing the saved AO telemetry is used to reconstruct the full CL DM commands using the full command matrix without any modal filtering.
The filtered and full wavefront residuals are then computed from their respective CL DM commands using the conversion factor and finally the filtered and full RMS wavefront errors are computed represented by

DAYTIME TESTING
As outlined in Sec. 3 changes were made to the current predictive filter implementation and new features such as RLMMSE were added. In this section, we present results from our daytime tests on the Keck II AO bench including results from these new features that have been tested during the day.  31 Operating at 1kHz, the turbulence moves less than 1 sub-aperture each frame. Hence, it is more relevant to think of each sub-aperture as a time series than to consider the spatial correlations between sub-apertures as done using Frozen Flow Hypothesis. More day time testing is needed to find the optimal parameters for the RLMMSE (such as the initial training time; data was taken with a few thousands frames of training, much less than EOF) as well as to improve the rate at which the RLMMSE updates. Currently, updating the filter is done within python and takes around 0.02 seconds; ideally the algorithm will be implemented in C enabling much faster update rates. With the current implementation, however, we see a stable performance indicating that the statistics of the system are not changing faster than our update rate with the AO loop running at 1kHz.
During the day of July 27, 2021, we tested EOF under a variety of different configurations and focused on varying the number of modes and the HO gain while taking K-band coronagraphic images with NIRC2.
We are limited to K-band by the internal light source, however, for on-sky tests we make use of the L-band vortex. We applied multi-layer turbulence with the parameters outlined in the Keck Adaptive Optics Note No 303 (KAON303) 32 removing piston, tip, and tilt. Setting α = 1, and training on 2 minutes of data, the EOF provided an improvement over the nominal integrator configuration, as described below. We plot the median raw noise curves for 325 controlled modes in Fig. 2 and for 300 controlled modes in Fig. 3 of the NIRC2 images throughout the daytime tests. These images show the contrast out to the DM control radius of 10λ/D. The noise as a function of separation was computed using the Vortex Image Processing (VIP 33 ) python package: a ring of FWHM-diameter circular apertures was constructed at each separation, and the noise was taken to be the standard deviation of the aperture sums at a given separation. Figures. 2 and 3 show that the predictor is able to perform better (by a factor of 3.5 for 300 modes and 4-5 for 325 modes) than the integrator at small separations (4-7 λ/D; where we expect predictive control to provide an improvement) under a variety of different system configurations. Figure 2 shows that varying the HO gain for the predictor does not have a large impact on its performance. We, however, acknowledge that this behavior might be very different on-sky because for this experiment we are introducing turbulence using the DM and therefore we do not have the same amount of HO spatial frequencies of atmospheric turbulence that would cause aliasing on the PyWFS that we would on-sky.
Finally, comparing the daytime integrator's performance, we see that the integrator is more affected by the number of corrected modes than the predictor and that fewer modes provide a better noise curve (comparing the two integrators plotted in Fig 3). This might be due to our Zernike modes being defined on a circular aperture but applied to a segmented pupil or that this is the response of our imaging arm with the vortex coronagraph.

ON-SKY TESTING
In this section, we present the results from our two engineering nights and one science night in 2021 that allowed for on-sky testing of predictive control.

Night of June 20, 2021
We were able to close the AO loop on J21095901+2859392 -a bright star with an H-band magnitude of 4.24.
Issues with the secondary mirror of the telescope and the suboptimal performance of the AO system resulted in low SR and the decision to take non-cornographic images. These issues mean the results from June 20, 2021 cannot be used conclusively to demonstrate the performance of predictive control compared to the integrator.
However, the night allowed us to verify the functionality of our implementation and observe the overall behavior The mixing factor (m) varied slightly while we changed the HO gain. The applied atmospheric turbulence on the DM was a multi-layer atmosphere model 32 with piston, tip, and tilt removed and a r 0 of 20 cm.
of the predictive controller. We note that during our actual tests no seeing data was available from the Canada-France-Hawaii Telescope's (CFHT) MASS and DIMM instruments. 34 In Fig. 4, we show the filtered RMS wavefront error (i.e., CL f ilt ) for the entire night. We also show the SR as measured from NIRC2 images. Studying the time series of the RMS, we can see the improvements provided by prediction as the RMS becomes smaller when predictive control is turned on (shaded orange regions; it is helpful to track the minimum RMS values to see the performance difference and not the maximum RMS values that are contaminated by telescope offloading, DM saturation, and other effects). The time series also reveals the effects of all the changes made to the AO controller while observing. The continued improvement of prediction in a shaded region (i.e., between times 14:40 and 14:45) is a result of the mixing factor and integrator gain being adjusted to find the optimal configuration. During this observation period, we also tested the effects of changing the assumed frame lag in the POL step. Figure 5 shows the results from a lag of 1 and lag of 2. With a lag of 1 we are able to improve the performance, but the histogram is skewed with a tail tending towards larger RMS values. This tail is removed when a lag of 2 is used by the predictor. We also confirm the lag of 1.7 frames as determined by S. Cetre et al. 2018. 30 We will investigate implementing a non-integer lag in the near future for our predictive control filter for further improvement. Furthermore, the lag of 2 histogram shows a greater improvement in RMS wavefront error. Comparing the NIRC2 SR for the integrator and predictor with lag 2 in Fig. 6, we see a clear shift towards larger SR values when predictive control is turned on. While we do find that a lag of 2 is better for the predictor and the predictor performs better than the integrator it is difficult to draw concrete conclusions on the provided improvement due to the lack of seeing data and the secondary mirror issues.

July 27, 2021
On the night of July 27, 2021, we received a half night of engineering time for predictive wavefront control tests during which the performance of the AO system was typical. In this section, we present the results for one of our targets, HIP 117578. With an H-band magnitude of 6.1, HIP 117578 is bright enough to take advantage of the 1kHz frame rate of the PyWFS camera. We took NIRC2 images with the L-band coronagraph in place and made use of QACITS to keep the point-spread-function (PSF) centered on the vortex coronagraph. 35  of time for both the full and filtered RMS values (i.e., CL f ull and CL f ilt ). We also plot the seeing as measured by the DIMM from CFHT. The figure indicates four data sets that were taken for this target. For each data set we configured the AO system as desired (see Tab. 1) and then started the QACITS sequence which first runs an optimization sequence to center the vortex and then takes a pre-determined number of science images. Half way through this set of science images we switch from the integrator controller to the predictive controller. These full QACITS sequences, which include NIRC2 frames taken using the integrator followed by predictive control, constitute a "dataset" (Tab. 1). The NIRC2 images were all taken with an exposure time of 0.3 seconds and 90 coadds at the full frame of 1024-by-1024 pixels. When we swap controllers, one to two NIRC2 frames are contaminated and hence dropped from this analysis. Therefore, we compare the predictor to integrator data that was taken just before the predictor data set. In Fig. 8, we plot the histogram of the RMS wavefront error for the four different data sets once again for the filtered and full RMS values. These histograms only include data that were taken during the NIRC2 exposures that are used for our contrast analysis. The ratio of the medians shows data for 2 frame lag predictor compared to the integrator. The integrator data was taken before and after prediction was turned on and off. The probability is the counts per bin normalized to the total counts in each dataset such that the sum of all the bins in each histogram is 1. for all cases when looking at the filtered wavefront error which makes sense as the LMMSE aims to minimize the residual wavefront error and we trained on the filtered wavefront. For the full wavefront, the predictor's median is worse than that of the integrator for dataset 3. The standard deviation for dataset 2 using both the filtered and full residuals wavefronts is improved by the predictor which indicates that the predictor provided a more stable AO correction than the integrator during this time.
We present the contrast curves for the four different data sets taken while observing HIP 117578 in Fig. 9.
The NIRC2 data was first bad pixel, flat-field, and sky corrected, as well as registered, using the automated pipeline described in Xuan et al. 2018. 26 We used VIP to de-rotate, median subtract, and median combine the data, modifying the number of images in each dataset such that the integrator and predictive components had the same total parallactic rotation. We then used VIP to create student-t corrected, algorithmic throughput corrected contrast curves, as detailed in Gomez Gonzalez et al. 2017. 33 From the contrast curves we find a contrast gain of 2 for dataset 2 at a separation of 3 λ/D as reported in Tab. 1. In all of the datasets, however, we see an improvement in contrast located at different separations, which can be seen in Fig. 9. Figure 15 in   App. A shows the contrast gain (ratio of the integrator and predictor contrast curves) for the four data sets, clearly showing the achieved improvement with predictive control. In data set 2 we see contrast gains close to 3 at separations of 4-6 λ/D while we see more modest gains around 1.5 for data sets 0 and 3. From these contrast curves, we have repeatable improvement in contrast with EOF predictive control take over a period of 1.5 hours while the only time we do not have a contrast improvement (dataset 3) we made a significant change to order of the predictive filter which we had not previously tested. We also plot the post-processed corongraphic PSFs for dataset 4 in Fig 10, further illustrating the spatial improvement provided by the predictor.
Studying all the data, it appears that the ratio of the full std of the RMS wavefront error is the best indicator for contrast improvement as shown in Tab 1. When the full std ratio is less than 1 the contrast is worse, when it is equal to 1 we have slight improvement and for larger than 1 we have the greatest contrast improvement.

October 25 2021
We did more testing with EOF predictive control during science observations on the night of Oct 25, 2021 for which we had typical AO performance. Not only were we able to demonstrate an improvement in residual   wavefront error again but also an improvement in contrast. Figure 11 shows the ADI median subtracted contrast curve for HIP 25645 which has a H-band magnitude of 6.77. The contrast is improved between 5-8 λ/D which is in agreement with the results from July and also close between 2-3λ/D. For these observations we corrected 212 modes with a HO and integrator gain of 0.28 and 0.4, respectively. These parameters were set by Keck personnel following standard AO operations. A mixing factor of 0.5 and gain of 0.23 where used for predictive control.

Summary of controller performances
The RMS wavefront error results from our latter two nights (Jul and Oct for which we had good AO performance) are summarized in Fig 12. We plot the seeing against the filtered residual RMS wavefront error. We take the CFHT seeing value and resample the RMS wavefront error to every 10 seconds to get a mean correction and then take the mean RMS wavefront values within 30 seconds of the CFHT seeing value. We make use of all the available predictor data that was taken. The upper and side subplots show the kernel density functions. The plots show that distribution of seeing values for all the data is roughly the same for both controllers (with the seeing actually being better in some cases for the October integrator than the predictor) on the two nights while the distribution in residual wavefront error is different for the two controllers over both nights. This indicates that the difference in the predictor and integrator is not due to a change in seeing. The predictor on both nights shows a consistent reduction in wavefront error compared to the leaky integrator as illustrated by the scatter plots being shifted to the left for the predictor datasets compared to the two integrator datasets. There is also a clear reduction in the spread of wavefront errors for predictive control showing that for different conditions the predictive controller has similar performance resulting in a more reliable and consistent performance. This is in agreement with what was found by van Kooten et al 2020. 11 6. DISCUSSION

Performance metric
We compare, in this paper, the residual wavefront error and raw contrast of predictive control to the current controller performance at Keck, the leaky integrator. More specifically, we define an improvement compared to the integrator as setup by Keck personnel such as the telescope operator, observing assistant, or staff astronomer who follow guidelines to provide robust, good AO performance depending on the conditions. This protocol was strictly followed for the Oct 25 night while loosely followed for observations in June and July (e.g., we manually increased the number of modes such that the integrator had the best performance). During June and July we altered the leaky integrator to perform even better than what was suggested by the operator by mainly correcting for more modes than recommended. In the end, our implementation of predictive control on Keck aims to produce the best performance possible with the current system through a data-drive approach. Therefore, our chosen metric is the current performance of the system.
Under the assumption that the servo-lag error is a large contributor to the wavefront error and dominating the contrast at small angular separations on NIRC2, 26 we implemented a data-driven predictor, EOF, to reduce the temporal wavefront error. Since the real goal is to improve the contrast at small angular separations, we, therefore, attribute an improvement in contrast to a reduction in the servo-lag error. Our method, however, relies on both spatial and temporal correlations, therefore we are also able to correct for other sources of errors in the system (i.e., calibration and registration alignment), making this an attractive method in general for improving the overall performance of an AO system.

Discussion of results
With the our predictive filter in place and demonstration of contrast improvement over a couple of nights and reduction in RMS wavefront error over several nights, we can transition towards making an comprehensive review of optimum parameters (see Sec. 7).
We began such testing in July 2021 where we aggressively varied the HO gain between 0.3 and 0, where 0 value indicates that we are not controlling any HO modes. It is when we are not controlling the HO gains (dataset 2) that we achieved the largest improvement in contrast. Comparing dataset 1 and 2 where the HO gain goes from 0.3 to 0, the contrast curves in Fig. 9 show the integrator is able to perform slightly better at larger separations than the predictor due to HO modes that are being controlled. Comparing then to the data taken in Oct shown in Fig. 11, where the HO gain value is similar to dataset 1, we are able to perform better with the predictor. We note that the number of total corrected modes are less in Oct as the number of modes needed to be lowered to provide a stable correction for the leaky integrator. Therefore, at larger separations, it is difficult to know what HO and number modes to use for the best improvement with prediction as it is clearly driven by external factors such as the atmospheric conditions. Close-in to the star, however, the predictor performs better regardless of the HO gains in those datasets. The transition separation between these two regimes appears to be around 4 λ/D.
Finally, we point out the behavior of the filtered and full residual wavefront error as shown in Fig. 7. In the bottom panel we plot CL f ilt which is used to generate the POL wavefront error that is used to generate the predictive filter by minimizing the mean square error. As we change various parameter, the CL f ilt varies in clear explainable way where for example by increasing the g from dataset 0 and dataset 1 decreases the residual wavefront error. In the top panel of Fig. 7 we plot CL f ull where we include the HO modes that we are not controlling into the residual wavefront error. Between datasets it varies in the opposite way from CL f ilt . Again taking dataset 0 and dataset 1, we are increasing the gain along with the HO gain such that we are improving the residual wavefront error for the controlled modes. We, however, appear to be increasing the full residual wavefront error. Testing along with end-to-end simulations are needed to better understand this behavior starting for the leaky integrator.

Power spectral density analysis
In Fig. 13 we present the estimated power spectral densities (PSDs) for the on-sky July data. From the POL PSDs we see that temporally the optical turbulence did not vary significantly between datasets and swapping of the controller. In closed-loop, the predictor reduces the wavefront error at small frequencies while also reducing the peak near 60 Hz. While the predictor PSD is more flat compared to the integrator, it also indicates that a greater improvement in performance is still possible by removing the peak around 60 Hz. To further improve the system and the performance of the AO system, we must understand where different features of the PSDs originate from. In the CL PSD the peak at near 60Hz has the same amplitude as the low frequencies suggesting that this peak is now a main source of the residual wavefront error. The frequency of the peak is similar to one seen in a tip/tilt PSD indicating that both peaks might have the same source but further investigation is needed.
Since EOF trains on the POL, the peak is a weak signal compared to the turbulence such that increasing the temporal order does not improve the performance. Therefore applying EOF to the closed-loop coefficients might help to remove this feature. We outline the steps necessary in Sec 7.

Comparing EOF to other controllers
We outlined above our performance metric used in this paper. Understand the performance of predictive control in comparison to other proposed controllers is outside the scope of this paper. We instead outline work already underway and provide comments on EOF and RLMMSE as well as gain optimization.
Comparing EOF to other controllers is underway in a laboratory setting at University of California Santa Cruz on the SEAL test bench. 36 We will compare predictive control not only to other types of predictive controllers but also to other controllers such as modal gain control. In our own simulations and with telemetry data the EOF and the RLMMSE algorithms have the same performance when the appropriate regressors and training data are selected. Even with fewer regressors the RLMMSE can achieve the same performance as shown in Fig. 14 where two different RLMMSE results are shown to have the same performance as EOF. The regressors and amount of training data influence the achievable performance and contribute to over-and under-fitting to the data. In the case where we use the most recent measurement only (i.e., temporal order 1) and a spatial order of one (i.e., temporal regressors only) the RLMMSE solution is equivalent to the optimal gain for the integrator for each actuator. If we change our basis set to modes, we will have the optimal gain for each mode. Adding more regressors, depending on the system setup, might not result in better performance and in the case where the lag is one frame, more regressors might result in over-fitting. Therefore care must be taken when setting up the predictive filter. Similarly modal gain optimization relies on pseudo-open loop reconstruction that incorporates the temporal lag. The optimal gains are found by solving a specific cost function and while the optimal gain is not a predictive algorithm, it improves the temporal performance by maximizing the correction bandwidth. The performance of the optimal gain optimization and a data-driven predictor might be indistinguishable especially when the lag is 1 frame and/or for non-stationary input disturbances. The data, regressors, optical setup of Figure 14: EOF on 30 seconds of Keck AO telemetry compared to two version of RLMMSE that use no spatial regressors and have the same temporal order (s1t10) and another which has 5 temporal orders and spatial order of a 3x3 grid around each actuator so that it correlates to its neighbors (s3t5). With EOF we train the data in advance on an earlier 30 seconds of the data while the RLMMSE trains online. The RLMMSE has many fewer regressors allowing it to converge quickly. We also plot the true open-loop (OL) and the Keck integrator (integrator).
the system, and the performance metric will all contribute to how different modal gain and predictive control methods will behave.

FUTURE WORK
Given the contrast improvement illustrated in Fig. 9, we have a clear path forward to improve the performance of the predictive control including: 1. implement predictive control in a true closed-loop configuration, 2. determine the optimal parameters for both EOF and RLMMSE, 3. match predictive control results to end-to-end simulations expanding on our integrator simulations that match on-sky contrast, 37 4. investigate modal implementation and direct use of slopes for prediction, and 5. continue on-sky verification of predictive control under different conditions and guide star magnitudes.
To address some of these topics above, we are currently developing a high fidelity model of Keck II AO bench for which we can better test and match on-sky data by having the exact infrastructure in simulation. With this simulation infrastructure, we hope to better understand the tunable parameters for future facilitization. Already, we find that the performance is improved and more stable within seconds of swapping to predictive control and little tuning is needed with selecting the nominal parameters presented in this work.
Finally, we outline below the next steps to operate predictive control in a true closed-loop configuration.
Implementing integral action on the predicted POL wavefront will provide a more stable performance while also limiting the number of tunable parameters, further simplifying the use of the algorithm.
In a closed-loop configuration the error measured by the PyWFS, e(t), is the difference between the input disturbance d(t) and the corrected phase given by the control signal at the last time step, V (t − 1), as shown in Eq. 9.
Assuming the PyWFS measures the true error signal perfectly we can write e(t) as a function of the slopes as in Eq. 10.
The predictive filter is predicting the full input disturbance, d(t), that will occur when the commands are sent.
We can plug Eq. 9 into Eq. 1 to get the closed-loop configuration of predictive control with integral action, Eq. 11.

CONCLUSIONS
We present the latest results from predictive control on Keck II AO system in this paper. We show that predictive control: 1. improves the on-sky SR of NIRC2, 2. provides a repeatable contrast gain of a factor of 1.5-3 on-sky at small separations (3-7 λ/D) with the vortex coronagraph, 3. provides a more reliable correction under different seeing conditions, 4. and the improved performance is not driven by seeing.
Close inspection of the RMS wavefront error during the exposure of the NIRC2 images used for contrast analysis shows a reduction in the standard deviations of the distribution of wavefront errors for predictive control, revealing that prediction provides a more consistent correction over these time periods. Studying the noise and contrast

APPENDIX A. CONTRAST GAINS
For further reference, we show the contrast gains achieved for the four different datasets in Fig. 15. These curves were calculated by dividing the contrast curves shown in Fig. 9.

ACKNOWLEDGMENTS
The predictive wavefront control demonstration is funded by the Heising-Simons Foundation.  and 'full', respectively) for each data set along with the ratio of contrast (C; contrast gain) at a separation of 3 λ/D. For each measured performance metric a value greater than 1 indicates the predictor is performing better than the integrator. The controller settings for each of the data sets are shown and the number of modes controlled for all data sets was 300 modes with the LO gain being set to 1. For the NIRC2 images, the exposure was set to 0.3 seconds and the co-adds to 90 for all images. The temporal order (Order) and cut-off-frequency (CoF) are also given for each data set. These values correspond to the histograms shown in Fig. 8 showing data taken the night of July 27, 2021.  14 EOF on 30 seconds of Keck AO telemetry compared to two version of RLMMSE that use no spatial regressors and have the same temporal order (s1t10) and another which has 5 temporal orders and spatial order of a 3x3 grid around each actuator so that it correlates to its neighbors (s3t5). With EOF we train the data in advance on an earlier 30 seconds of the data while the RLMMSE trains online. The RLMMSE has many fewer regressors allowing it to converge quickly.