Noise-tolerant wavefront shaping in a Hadamard basis

Light scattering is the main limitation for optical imaging. However, light can be focused through or inside turbid media by spatially shaping the incident wavefront. Wavefront shaping is ultimately limited by the available photon budget. We developed a new `dual reference' wavefront shaping algorithm that optimally uses the available light. Our method allows for multi-target wavefront shaping, making it suitable for transmission matrix measurements or transmitting images. We experimentally confirmed the improvement of the focus intensity compared to existing methods.


Introduction
Scattering and diffusion of light are the main limitations for optical imaging. However, wavefront shaping techniques make it possible to shape the incident wavefront so that it exactly matches the scattering properties of the medium, and form a focus through or inside it [1]. Feedback-based wavefront shaping algorithms search for the optimum incident field that maximizes some feedback signal, such as transmission into a specific output mode [2]. The task of finding the optimum wavefronts for multiple targets (output modes) is equivalent to measuring the transmission matrix (TM) of the scattering medium, up to an unknown phase factor for each row of the matrix [3]. The TM describes the input-output response of the scattering medium. Knowledge of this matrix makes it possible to focus at any arbitrary position behind the medium [4,5], to transmit images through it [3,6], or to study fundamental aspects of light scattering, such as the distribution of transmission channels [7,8].
A common procedure to find the optimum incident field with a phase-only spatial light modulator (SLM) is using the stepwise sequential algorithm (SSA) [2], as illustrated in Fig. 1(a). In this algorithm, the SLM is divided into segments, and the phase of each individual segment is varied between 0 and 2π consecutively while keeping the phase of all other segments fixed. Due to interference between light originating from the controlled segment and light originating from all other segments, the feedback signal will respond sinusoidally. By fitting these sinusoids, the TM is reconstructed [2].
The drawback of the SSA algorithm becomes apparent when many segments are used. In this case the contribution of each segment is small compared to the reference field coming from the rest of the segments, giving a low interferometric visibility of the feedback signal. This results in a low signal to noise ratio (SNR) in the measured TM elements. This low SNR is especially noticeable when the photon budget for the measurements is limited, such as in high-speed wavefront shaping [9], and in microscopy.
A method that was introduced by Popoff et al. [10] organizes the segments in a Hadamard basis. The feedback signal is measured for each Hadamard vector, allowing the TM to be reconstructed. In this method, always a large fraction of the segments is modulated, increasing the visibility of the interference signal, and consequently improving the SNR [10]. However, this technique requires part of the incident light to remain unmodulated serving as a reference field for the measurements. Since part of the light is not modulated, the contrast of the final focus is not optimal. This method is arXiv:2103.12637v1 [physics.ins-det] 15 Mar 2021 illustrated in Fig. 1(b), showing 65% of the SLM are always modulated (in red) and the other 35% of the segments are the reference part [10].
A different way to increase the SNR for a given photon budget is by performing a pre-optimization [11,12]. In pre-optimization schemes, one first performs an optimization with part of the SLM segments and uses that solution as a starting point for a second optimization step. With pre-optimization, the feedback signal in the second optimization step will be higher, causing the second step to be more robust to noise in most cases [11].
A special case of pre-optimization was developed by Tao et al. [12]. In the first step, half of the segments is modulated using a Hadamard-based algorithm while the other half is used as a reference. Since the number of controlled segments is equal to the number of reference segments, this approach leads to an optimally balanced interference of modulated and reference field. In the second step, the optimized half of the wavefront is displayed on the SLM. The roles of modulated and reference segment are now switched (see Fig. 1(c)), and the other half of SLM segments is optimized.
The method of Popoff et al. can be used to simultaneously find the wavefronts to focus light to multiple targets, or even to project arbitrary images through the scattering medium. In effect, this method measures multiple rows of the TM simultaneously. Unfortunately, pre-optimization approaches can not be used for this purpose. After the pre-optimization step, one has to choose what wavefront to use as a starting point for the next step, so the procedure can only be used on one point at a time.
In this work, we present a new algorithm, called dual reference, which is a combination of the methods used by Popoff et al. [3], and Tao et al. [12]. We take advantage of the optimal balanced interference of the modulated and reference field as in Tao et al., while keeping the ability to measure all rows of the TM simultaneously, as in Popoff et al.
We will first describe our dual reference method, and explain why it outperforms other methods for a wide range of photon budgets. Afterward, we compare the experimental results of the performance of the proposed method with the existing feedback-based wavefront shaping methods for various photon budgets.

Dual reference algorithm
In our dual reference algorithm, we split the SLM segments into two groups of the same size with a small number O of overlapping segments, as illustrated in Fig. 1(d). The groups have size N 1 = (N + O)/2, where N 1 is a power of two. The field in the b-th output mode is now given by where E a is the field in the a-th input mode, N is the total number of SLM segments, and t ba are the elements of the transmission matrix. In the first step, the algorithm of Popoff et al. [10] is used: a Hadamard pattern is displayed on the segments of group 1 (segments 1 to N 1 ) and the remaining segments are used as a reference. The relative phase between the two groups is changed from 0 to 2π in P steps. The measured feedback signal for Hadamard vector j ∈ [1, N 1 ] and phase step k ∈ [1, P ] is given by where h ja ∈ {−1, 1} is the a-th element of the (non-normalized) j-th Hadamard basis vector, E 0 is a constant to normalize the incident field, and the reference field from group 2 is given by In order to retrieve the optimized wavefront, we first perform a Fourier transform of the measured intensity Finally, a Hadamard transform is used to convert S j to S a1 = E * ref,2 t ba E a . Assuming that we have an even illumination on the SLM, this gives us the transmission matrix elements t ba for the elements of group 1 (elements 1 to N 1 ), up to an unknown overall factor E * ref,2 . Our modification to Popoff 's algorithm is to add a second step, where group 1 is now used as a reference, and the segments of group 2 are modulated. These two steps are shown in Fig. 1 (d) with the modulated SLM segments in red and the overlap segments as crossed segments. Since the number of controlled segments is close to the number of reference segments for all measurements, the modulated and reference fields have almost equal average contribution to the field at point b. This 50/50 choice for the reference and modulated part guarantees the maximal interferometric visibility and minimal sensitivity to noise for a given photon count.
After completing the two steps, we have the corrected wavefronts for group 1 and group 2. To find the final optimal wavefront, these two partial corrections need to be combined. However, since the two groups used different reference fields in the measurements, the two partial wavefronts will have a different phase offset. Since the overlap segments are present in both groups, these segments can be used for determining the relative phase between the two corrections. Denoting the overlap segments as a = 1..O, we compute the relative phase from the dot product of the two measurements for these overlapping segments and the relative amplitude follows from After compensating the relative phase and amplitude of the two partial wavefronts, these two wavefronts can be combined into a single wavefront.

Experiment
Now, we will experimentally compare the performance of the dual reference algorithm to the original SSA algorithm and to the methods of Popoff and Tao. Since Tao's method can only measure one row of the TM at the time (hence only focus on one target), the comparison with that method is not completely fair. For completeness, we do include the method in the comparison, using only a single target in the measurements.
The experimental setup is illustrated in Fig. 2. In this setup, light from a HeNe laser (Thorlabs, 2mW, λ = 632.8 nm) is expanded and modulated by a phase-only spatial light modulator (Hamamatsu X13138-07). A 4f system images the SLM onto the back focal plane of a microscope objective (A-Plan 100x/0.8, Zeiss), which focuses the light onto the sample. After transmitting through the sample, the light is collected by an identical objective lens and recorded by a CMOS camera (Basler acA640-750um). The camera records the intensity distribution at the back focal plane of the second objective. A set of 100 independent targets on the camera is used as the feedback for the wavefront shaping algorithms.
The sample is an 11±3 µm thick layer of zinc-oxide (Sigma Aldrich, average grain size 200 nm) on a coverslip with a thickness of 170 µm. The transport mean free path of similar zinc oxide samples was measured to be around 0.6 µm at a wavelength of λ = 632.8 nm [13]. Consequently, the sample is optically thick so that there is no transmitted ballistic light. The sample is mounted on a translation stage (Zaber T-LSM050A) to scan the sample laterally for ensemble averaging.
The beam illuminating the SLM has a Gaussian intensity profile. It is well known that such a non-uniform illumination causes a sub-optimal enhancement [14]. To overcome this effect, we choose a circular geometry instead of the common rectangular one, as is shown in Fig. 2(inset). The sizes of the rings are adjusted such that each segment, on average, contributes equally to the target signal. This way, the illumination effectively is uniform, and the enhancement is maximized.
To control the photon budget used in the experiments, we used various neutral density filters with different transmission ratios in front of the camera. The power of detected light was measured by a calibrated power meter (S121C, Thorlabs). The measured power is converted to the photon count collected by the camera. We kept increasing the total photon budget as long as the camera did not reach saturation during the measurement.
In the wavefront shaping experiments, we measured the enhancement of the focus, defined as the ratio of the optimized intensity at the focus location to the reference intensity which is the averaged intensity over 100 positions for the sample [1].
In our experiments, measuring the corrected wavefront for one target took 50 seconds. For more than one target, this time interval does not change for all the algorithms except for Tao's method, which is linearly increased with the number of targets.

Optimal number of phase steps
All four algorithms require the choice of two parameters: the number of SLM segments (N ) and the number of phase steps per segment (P ). Together, these parameters determine the total number of measurements N P . It is well known that in absence of noise the expected enhancement scales linearly with N [15], so N should be chosen as high as possible, resulting in a minimal number of phase steps P = 3. In presence of noise, however, the relation is more complicated, and different types of noise (such as read-out noise, shot noise and laser excess noise) affect the results in different ways [11]. When the number of measurements N P is fixed, a higher value of P will give a more accurate measurement of each transmission matrix element, at the expense of decreasing N .
To find the optimum value for P for the different algorithms, we measured the enhancement for different combinations of P and N . To enable a fair comparison, we choose the number of SLM segments in a way that all algorithms have the same number of measurements. For the lowest phase steps number (P =3), we set 128 segments for Tao's method [12] and SSA. For Popoff's method [3] the SLM is divided into 197 segments, out of which 128 segments are modulated, and the remaining 69 segments (35%) serves as the reference. For our dual reference algorithm, we had 123 modulated segments with 5 segments mutual in two 64-segment groups. As a result, the number of measurements for all the algorithms is the same, and equal to 384. We kept this number fixed by reducing the number of segments by two or four when P is respectively 6 and 12. Figure 3 shows the averaged enhancement over 100 targets and 3 trials versus total photon budget when P is 3 (solid curve), 6 (dash curve), and 12 (dot curve). For the dual reference method Fig. 3(a), Popoff's method Fig. 4(b), and Tao's method Fig. 4(d), the highest enhancement for every photon budget is obtained for P = 3. However, Fig. 4(c) shows that for a photon budget lower than 5 · 10 4 , SSA shows a slightly better performance for larger values of P .
The results from Fig. 3 show that a lower P with higher N is better than performing accurate measurements for fewer segments. Therefore, we set P = 3 for all algorithms and compare the performance of the dual reference algorithm with the other algorithms. Figure 4 illustrates the measured enhancement for the different methods, averaged over 100 targets and 3 trials. For all algorithms, the enhancement increases as the total photon budget increases since the signal to noise ratio increases.

Results and discussion
The measured enhancement for Tao's method, shown in the red dot curve in Fig. 4, is higher than the other algorithms due to the pre-optimization step. However, finding the corrected wavefronts for all 100 targets took 100 times longer than the other three algorithms. This is because the corrected wavefront from step 1 is different for every target and it prevents us from simultaneous multi-target optimization. During step 2, the corrected wavefront from step 1 is displayed on the SLM, causing the camera to be overexposed for higher photon budgets. For this reason, there are fewer reported measurements for Tao's method than for the other methods.
The measured enhancement for the dual reference algorithm, shown as the blue solid curve in Fig. 4, confirms that using this method gives us higher enhancement for multi-target optimization compared to SSA (black dash-dot curve) and Popoff's method (green dash curve).
For a high photon budget, we theoretically expect to have an enhancement of N π/4 for SSA, dual reference, and Tao's methods, and 0.65N π/4 for Popoff's method. In our experiments, the measured enhancement for the highest photon budget (1.8 · 10 5 photons) is 36 % lower than the expected enhancement from the theory. This difference can be explained by practical imperfections, such as the static (non-modulated) field caused by the reflection of the SLM front surface. In our case, this static field has an amplitude of 10 % of the input field amplitude.  [10], and (d) Tao et al. [12]. Bars represent the standard error of the measurement set.

Conclusion
With the advance of high-speed SLMs, the photon budget will ultimately be the limiting factor for wavefront shaping. We have presented a new feedback-based wavefront shaping algorithm to optimally use the available light.
We demonstrated that for a fixed number of measurements the best performance of all algorithms is reached by minimizing the number of phase steps.
Our algorithm achieves the maximal interferometric visibility during the measurements, resulting in an optimal SNR. Moreover, there is no need to reserve segments for the reference beam. We have experimentally demonstrated that this method enables us to perform simultaneous multi-target optimization with a higher enhancement than the popular methods used by  and Tao et al. [2017] (for multi-target) for a wide range of photon budgets.
We envision that the presented method could be beneficial to the experiments having the transmission matrix measurement involved, especially when the photon budget is limited.

Funding
This work was financially supported by the Nederlandse Organisatie voor Wetenschappelijk Onderzoek (TTW-NWO, Vidi grant 14879).

Disclosures
The authors declare no conflicts of interest.