Wavefront control in adaptive microscopy using Shack-Hartmann sensors with arbitrarily shaped pupils

: In adaptive optical microscopy of thick biological tissue, strong scattering and aberrations can change the effective pupil shape by rendering some Shack-Hartmann spots unusable. The change of pupil shape leads to a change of wavefront reconstruction or control matrix that should be updated accordingly. Modified slope and modal wavefront control methods based on measurements of a Shack-Hartmann wavefront sensor are proposed to accommodate an arbitrarily shaped pupil. Furthermore, we present partial wavefront control methods that remove specific aberration modes like tip, tilt and defocus from the control loop. The proposed control methods were investigated and compared by simulation using experimentally obtained aberration data. The performance was then tested experimentally through closed-loop aberration corrections using an obscured pupil.


Introduction
Adaptive optics (AO) systems generally use a wavefront sensor like a Shack-Hartmann sensor (SHS) to measure aberration and a wavefront corrector like a deformable mirror (DM) to perform phase-conjugate compensation. An AO system running in closed-loop continuously converts the slope measurements of the SHS to driving commands of the DM. In most cases, the pupil shape is circular and unchanged during correction. In other words, the amplitude variations at the pupil plane are typically small and have little effect on slope measurements. However, in some applications with severe aberrations the pupil shape can change during correction and the wavefront slopes in some pupil regions will become unmeasurable. Recently, SHS based AO has been used in confocal and two-photon microscopes to correct specimen induced aberrations to improve image quality [1][2][3][4]. Strong scattering in specimens can create varying dark regions within the pupil and spots in some subapertures can be totally missing (see Fig. 11 in [2] and Supplementary Fig. 1 in [4]). The wavefront slopes in dark regions were replaced by zero value in [2], which would however misestimate the aberration. In essence, the change of pupil shape leads to a change of wavefront reconstruction or control matrix, which should be updated accordingly during correction.
Wavefront control can be a two-step process that first converts the slope data to optical path difference or phase distribution (i.e. wavefront reconstruction) and then calculates corrective commands for the DM. In [5], the pupil shape is assumed circular and unchanged. Wavefront reconstruction algorithms from gradients within arbitrarily shaped pupils has been reported in previous literature for optical metrology or testing [6][7][8][9][10], however, not yet concerning wavefront control. These algorithms are classified as either zonal or modal, depending on whether the representation of reconstructed wavefront is values at sampling grid points or coefficients of orthogonal polynomials defined over the pupil area. Zou & Zhang proposed a zonal least-squares algorithm using a domain extension technique to establish a normal equation set for irregular-shaped pupils [6]. The accuracy of this algorithm can be improved by introducing Gerchberg-type iterations, which however is very time consuming [7,8]. Mochi & Goldberg developed a modal, non-iterative wavefront reconstruction algorithm based on Gram-Schmidt orthonormalization [9]. Ye et al. presented a modal algorithm based on numerical orthogonal gradient polynomials that can be used to directly fit the slope data within arbitrarily shaped pupil [10]. As speed of operation is important, we will not adopt any iterative zonal reconstruction algorithms and will only focus on modal reconstruction algorithm. After wavefront reconstruction, the control signals of the DM are calculated by simple matrix manipulations, provided that the response matrix of DM is pre-calibrated by an SHS sensor or interferometer. We call this modal control hereafter.
If we do not need to know the phase distribution exactly, the intermediate wavefront reconstruction step can be omitted and the slope data are directly translated to DM commands, namely the slope control in [11]. However, in many applications, only a part of phase aberration is assigned to the DM for correction. It is common in adaptive microscopy to remove tip, tilt and defocus modes from the control matrix, as these modes cannot be detected in an epi-configuration microscope. Furthermore, AO systems consisting of two DMs (namely woofer and tweeter) have to correct low and high order aberrations in a decoupled manner. In this case, it is essential to constrain each DM from generating some aberration modes to realize decoupling control. This is achievable by adding proper constraints on the slope control matrix [12,13].
For our current aim of using SHS based AO in biological microscopy, both modal and slope wavefront control approaches compatible with arbitrarily shaped pupils are investigated in this paper in order to find an improved solution to this issue. It is important to discuss both approaches, as the effects of missing SH spots could be different in each case, as are the mathematical methods to deal with this problem. We also investigate partial wavefront control methods that exclude specific aberration modes from the correction loop. The principles and procedures of slope and modal wavefront control for arbitrary pupils are presented in Section 2. Numerical simulations and results are given in Section 3. Experimental demonstrations are shown in Section 4.

Slope Control with arbitrarily shaped pupil
The principle of slope control is expressed as S is a slope vector measured by the SHS in both x and y directions; A is the DM command vector; Rs is the slope response matrix of the DM; M is the number of subapertures of the SHS and N is the number of actuators of the DM. Each column of Rs is composed of the slope variations after sending a unit driving command to an actuator. Rs typically has many more rows than columns and is a sparse matrix because of the finite spread of each actuator's influence function. So Eq. (1) is essentially an overdetermined linear equation and A can be solved by a least-squares method through a pseudo-inverse; ( ) A more robust way to calculate A is using singular value decomposition (SVD). The matrix Rs can be represented as the product of three matrices.
The pseudoinverse of Rs is obtained by Σ is a diagonal matrix consisting the singular values of Rs. + Σ is formed by replacing every non-zero diagonal entry by its reciprocal and transposing the resulting matrix. The ratio of the maximum singular value to the minimum is called the condition number, which is required to be small to achieve stable closed-loop correction. If any diagonal element of Σ is close to zero, the condition number will be very large and Rs is approximately singular. To get around this problem, one can simply zero the elements in + Σ that are related to small singular values. Furthermore, to have a small low condition number, the subaperture number across an actuator should be optimized. As shown in [14], at least 1.5 subapertures per actuator are required. The command vector A is given by + s A = R S (5) In practice, Rs is measured by the same SHS that is also used for aberration detection and + s R should be calculated in advance before correction. The pupil area should be fully illuminated during calibration in order to obtain complete slope data in Rs. However, the pupil shape may change during correction and the wavefront slopes in some subapertures will become unmeasurable. We have two options to deal with this problem. The first one is simply replacing the unmeasurable slopes with zero in the vector S and while still using the precalculated + s R [2]. This is a simple approach but could lead to lower corrective accuracy. The second option is to remove rows related to invalid subapertures from the matrix Rs and recalculate its pseudoinverse. In the second way, Eq. (5) is revised as where s R  and S  are composed of available slope data.

Modal control with arbitrarily shaped pupils
In modal control, the wavefront aberration represented by modal coefficients is reconstructed from discrete slope data measured by the SHS. The wavefront aberration can be decomposed into orthogonal modes defined over the pupil area where ( ) , i Q x y are orthogonal modes in the pupil area and i c are modal coefficients.
Using matrix notation, Eq. (7) can be rewritten as The measured slope vector S is given by where ∇ denotes evaluation of gradient in the x and y directions. The coefficient vector of wavefront aberration is solved from Eq. (9) as ( ) After wavefront reconstruction, the DM command vector A can be solved from Rc is the modal response matrix of the DM, which is converted from the pre-calibrated slope response matrix Rs as From Eqs. (11) and (12), vector A is calculated by where nk D is the conversion coefficient and K is the number of Zernike polynomials used. Equation (14) can be rewritten in matrix form as For arbitrary pupils, Eqs. (11), (12) and (13) can be rewritten as ( ) where S  , s R  and ∇Z  are calculated from valid subapertures. The matrix P should be updated by Cholesky decomposition of T Z Z   if the pupil shape varies, as the samples used to calculate Z will change.

Partial correction with arbitrarily shaped pupils
The wavefront control methods based on Eq. (6) or Eq. (21) try to correct all aberration modes. However, partial correction that keeps some modes uncorrected is required in some cases. This is common, for example, in biological microscopy in an epi-illumination microscope, as tip, tilt and defocus modes cannot be detected, so are usually removed from the control loop. For modal control, partial correction is easy to implement by setting specific modal coefficients in vector C to zero and Eq. (21) can be modified as where where n F is the influence function of actuator n and ( ) (23) can be rewritten as 1 1 , 0 Using the orthogonality , (1) and using matrix notation, we can get  is a row vector consisting of each influence function's modal expansion coefficient of mode j Q . d R can be easily expanded (i.e. to contain more rows) to inhibit the correction of more aberration modes. Please note that d R is actually composed of some selected rows of c R in Eq. (12). γ is a gain factor. Since the row number of Rs is much bigger than that of Rd, γ should be large enough to ensure the suppression of modes in d R . Here we choose 1000 γ = empirically in both simulation and experiment. The command vector A can be solved from Eq. (27) still using the SVD method.

Wavefront control procedures for arbitrarily shaped pupils
The flowcharts of slope control and modal control for arbitrary shape-changing pupils are shown in Fig. 1. The blue boxes in Fig. 1(a) denote extra steps due to requirement of partial correction in slope control. These steps are very similar to those in the blue boxes of Fig. 1(b). So if partial correction for a shape-changing pupil is required, the workload of one cycle in slope control is comparable to that in modal control. Otherwise, slope control is much simpler and faster than modal control.

Simulation Model
To investigate the performance of the proposed wavefront control methods, an adaptive optics system model including a SHS and a DM has been established in MATLAB. The specifications of the SHS and the DM in our simulation were consistent with the devices used for experimental demonstration (detailed in Section 4). The mapping of the SHS and the DM is depicted in Fig. 2. The pupil diameter was set as 85% of the side length of the lenslet array. Some actuators were deliberately arranged outside the illuminated pupil to get better correction results [16]. In each subaperture, the spot intensities were quantized to an 8 bit image. Additive Gaussian read out noise with RMS value of 10 (on 0 to 255 grayscale) was applied to the images. The centroid position of each spot was estimated by thresholded center of gravity method [17].
where ( , ) n n x y are the actuator position coordinates; c is the inter-actuator coupling factor which was set as 0.13; d is the actuator spacing.
The input wavefront errors were specimen-induced aberrations measured by a Mach-Zehnder phase stepping interferometer [18,19]. The measured aberrations were then unwrapped using a minimum L p -norm method [20]. 256 wavefronts were recorded for each specimen at 16 × 16 grid points. As an example, a wavefront aberration of C. elegans sample is shown in Fig. 3. Both aberration distribution map and the corresponding first 65 Zernike coefficients (without piston) are given.

Simulation Results for Circular Pupil
In this section, we will give the wavefront correction results using modal and slope control with a fully illuminated circular pupil. Wavefront correction results using modal control are shown in Fig. 4    The wavefront correction results using slope control are depicted in Fig. 6. For comparison, the wavefront aberration directly fitted by the DM and the residual error are shown in Fig. 7 where the measurement error of the SHS is avoided.  In principle, the aberration modes of tip, tilt and defocus are self-correcting and are not detectable by a SHS sensor in the epi-configuration normally used in biological microscopy [19]. However, these modes were measured by the interferometer in [18]. Here we perform partial wavefront correction that excluded tip, tilt and defocus for aberration shown in Fig. 3 and the results are given in Figs. 8 and 9. The coefficients of tip, tilt and defocus in the residual wavefront were calculated to compare with their initial values given in Fig. 3. From Figs. 8(b) and 9(b), both modal and slope control could ensure that the first three aberration modes are almost uncorrected. However, high order modes in slope control were not corrected as well as in modal control. We used R c in Eq. (29) to evaluate the residual high order aberrations in partial correction.
where C i are modal coefficients except for the first three terms, thus excluding the contributions of tip, tilt and defocus. In modal control R c was only 0.0029 μm RMS but in slope control was 0.017 μm RMS. So the aberration modes were not fully decoupled in partial slope control by adding constraints in control matrix in Eq. (27).

Simulation Results of Arbitrary Shape Pupil
For two-photon microscopy, the ballistic component of the excitation light decreases exponentially with imaging depth [4]. The effect of scattering combined with aberration in specimen leads to corrupted or missing spots on SHS. These missing measurements should not be involved in wavefront control, thus in effect they change the pupil shape over which wavefronts should be fitted. We use the phase measurement results of thicker brain tissue of about 90μm in [18] to emulate this effect. The wavefront aberration of thick brain tissue and the corresponding spots pattern on SHS are shown in Figs. 10(a) and 10(b). Spots in some subapertures became speckles or were even missing because of severe local high order aberrations. Here we used image sharpness (IS in Eq. (30) as a criterion for spot selection. where I(x, y) is the intensity distribution within a subaperture. In our simulation, an unaberrated spot had an IS of 0.103. The spots with IS lower than 0.03 were discarded and were not involved in wavefront control. Figure 10 (c) shows the wavefront with a pupil mask to filter out invalid subapertures. Following the wavefront control procedures in Fig. 1, we performed wavefront correction for all aberration modes and the residual wavefront errors are shown in Figs. 11(a) and 11(b). Residual wavefront error using direct DM fitting is given in Fig. 11(c) for comparison. If we use zero replacing invalid slope measurements and drive the DM using its original slope response matrix, the residual wavefront is given in Fig. 11(d), which shows much lower correction accuracy than our methods. We also applied partial correction excluding first three modes on an irregularly shaped pupil in Fig. 10. The nth term of orthogonal modes defined over the irregularly shaped pupil was formed by a combination of the first n terms of Zernike polynomials. So the first three terms of the new orthogonal modes would still represent image displacement. The coefficients of new orthogonal modes were calculated by Eq. (19) and shown in Fig. 10(d). Partial wavefront correction results are given in Figs. 12 and 13. In modal control, the first three modal coefficients were very close to their initial values and other high order modes were well corrected. In slope control, the first three modes remained uncorrected but with noticeable drift from their initial value and the high order modes were not fully suppressed.

Experimental Demonstration
The experimental system layout is depicted in Fig. 14. A single-mode, fiber-coupled laser (λ = 527nm) was used to mimic a point source. The collimated light was reflected by a 140-actuator DM (Thorlabs) and then divided into two beams by a beam splitter (BS). One beam was imaged by lens L 2 onto a camera to observe the focal image. The other beam was detected by a 25 × 25 SHS (WFS150-7AR, Thorlabs) after passing the relay lenses L 3 and L 4 . Two irises in front of the DM were used to block out beam to create a non-circular shape pupil to enable testing of the methods with obscured subapertures. Specimen-induced wavefront aberrations were simulated by the DM which also acted as the wavefront corrector. The wavefront control software was developed in LABVIEW.
The influence functions of the DM were measured by the SHS with fully illuminated pupil to obtain the slope response matrix Rs. Wavefront correction results with a fully illuminated circular pupil are shown in Fig. 15. Zernike coefficients of wavefront before and after correction are plotted for comparison. A total of 436 lenslets were illuminated and used by the SHS for aberration estimation. For obstructed pupil, 371 lenslets were illuminated and the wavefront correction results are given in Fig. 16. Wavefronts before and after correction are fitted by orthogonal modes.
The experimental demonstrations show that the methods work effectively for the compensation of aberrations when spots are missing from SH sub apertures in a closed loop system. As the DM was used in this experimental demonstration both to induce and correct the aberrations, we have avoided the possibility of aliasing of higher order non-correctable modes. However, in combination with the results of the previous sections, where higher-order aberrations were indeed present, these results indicate that the method should function well in the presence of both complex aberrations and obscured pupils. Fig. 15. Wavefront correction results with fully illuminated circular pupil. (a) Initial wavefront, RMS = 0.145μm, (b) Residual wavefront using modal control, RMS = 0.00917μm (c) Residual wavefront using slope control, RMS = 0.00814μm (d) Residual wavefront using modal control and partial correction, RMS = 0.0906μm (e) Residual wavefront using slope control and partial correction, RMS = 0.0870μm (f) Zernike coefficients of wavefront before and after correction. The signal processing time (excluding exposure time and DM settling time) using slope control in full correction was about 3ms with a general-purpose computer (3.5 GHz CPU). This time increased to 6ms for partial slope control or modal control in full or partial correction. This speed is compatible with changing of imaging position of specimens.

Conclusions
We have investigated the performance of slope control and modal control for both circular and irregularly shaped pupils, including when sub-apertures are obscured by specimen scattering, by simulation and experiments. Although the simulation was based on the correction of specimen-induced aberrations in biological microscopy, the conclusions can be generalized to other AO systems. For full aberration correction on arbitrarily shaped pupils, both slope and modal control can suppress the aberration to a low level depending on the corrective capability of the DM and the sensing accuracy of the SHS. Slope control is recommended in this case because of its simplicity and because there is no need to consider the optimization of modal control. For partial correction, modal control is superior to slope control because the aberration modes are decoupled directly in the wavefront reconstruction process.

Funding
National Natural Science Foundation of China (No. 61505008). European Research Council (ERC) under the Horizon 2020 research and innovation program (AdOMiS, grant agreement no. 695140).