Robust Phase Retrieval with Green Noise Binary Masks

Many measurement systems can only detect the intensity of the optical wave field. Phase retrieval, as the name suggests, aims to reconstruct a complex-valued signal from intensity-only measurements [1, 2, 3, 4, 5, 6, 7]. The phase retrieval problem is encountered in several disciplines, including crystallography, optical imaging, astronomy, X-ray, and electronic imaging. There has been considerable recent progress in phase retrieval algorithms due to the advent of modern optimization theories [8]. In particular, the so-called coded diffraction pattern (CDP) phase retrieval algorithm adopts pre-defined optical masks to improve the reconstruction performance with extra constraints [9, 10]. Specifically, the algorithm reconstructs a signal $\mathbf{U}\in\mathbb{C}^{n}$ (complex-valued) based on the Fourier intensity measurements $y_{i}=|\mathcal{F}(\mathbf{I}_{i}\circ\mathbf{U})|^{2},i=1,\dots,M$, where $M$, $\mathbf{I}$ and $\circ$ denote the number of measurements, pre-defined optical masks and elementwise multiplication, respectively; and $\mathcal{F}$ is the Fourier transform. In practice, the pre-defined masks can be implemented using a spatial light modulator (SLM) or digital micromirror device (DMD) [11, 12, 13, 14]. The number of measurements is related to the complexity of the signal and the masks. For instance, the single-shot lensless phase retrieval methods in [12, 13] require only a single measurement, but it can reconstruct only sparse objects. On the other hand, using multiple-level phase masks can also reduce the number of measurements. However, they require optical devices that can perform accurate multi-level phase modulation for their implementation. The cost of these devices is one concern; the error of these devices due to the global grayscale-phase mismatch and spatial non-uniformity is another. If binary masks are used, it is empirically shown in [9, 10] that at most $6$ measurements are required for exact recovery via convex programming. Thus, the number of measurements is not excessive. Binary masks can provide accurate reconstruction as the two levels can be reliably reproduced; these can be realized with both DMDs and SLMs. For this reason, binary masks are used in this work.

The importance of the randomness of the masks is emphasized in some research papers to ensure good phase retrieval performance [9, 10]. Traditionally, white noise masks are commonly used since randomness is guaranteed. However, the Fourier spectrum of white noise masks is non-bandlimited in theory. When using them for phase retrieval, the intensity measurements can have many significant high-frequency components located beyond the $0$th diffraction order (an example is shown in Figure 1(a), the yellow circle denotes the boundary of the 0th diffraction order). In general, we always focus on the $0$th diffraction order such that the data outside the $0$th diffraction order are not detected. Even if we also detect the data outside the $0$th diffraction order, they will still be ignored in the phase retrieval optimization process. This is because traditional CDP phase retrieval algorithms assume the masked data as a grid of sampled data so that its spectrum will be periodic based on the sampling theory. The algorithm can only handle the frequency components within the Nyquist frequency (i.e., the red square in Figure 1(a)). Therefore, the high-frequency components located beyond the Nyquist frequency of the system will be ignored by the phase retrieval algorithm. It is equivalent to filtering out the high-frequency components of the intensity measurements, and this introduces errors to the optimization process. In this paper, we base on our previous work [15] on digital halftoning to propose a green noise optical masking scheme for phase retrieval. The proposed scheme allows the optical masks to have the green noise property, which means that the energy of the masks is concentrated in a ring shape mid-frequency region. It allows the intensity measurements to also have the data concentrated in the same region (plus those around the zero frequency, an example is shown in Figure 1(b)). Thus, the amount of the filtered high-frequency components can be reduced. Halftoning is a process commonly used in digital printers to simulate shades of gray by varying the size of tiny black dots. Recently, we apply the multiscale error diffusion (MED) technique to the halftoning process so that we can control the black dot patterns to have the green noise property [15]. Such a technique can also be used in optical mask generation. One advantage of the scheme is that we can generate the binary masks of different densities (the ratio of the number of ones to zeros in the masks) while maintaining the green noise property. The density of the masks needs to be carefully chosen to avoid the saturation problem when taking the measurements. By using the proposed green noise masking scheme, the error of the reconstructed phase images can be reduced compared with using the traditional white noise masks. It has been verified with our experimental results in Section 4.

The objective of the CDP phase retrieval algorithm [9, 10] can be described by Eq. (1).

where $\mathbf{y}_{i}$ is a Fourier intensity measurement and $\circ$ denotes the elementwise multiplication. The use of the optical masks $\mathbf{I}$ provides the constraints to reduce the ill-posedness of the problem. It has been shown in [9, 10] that $\mathbf{U}$ can be perfectly reconstructed from noiseless measurements $\mathbf{y}$ with high probability given $\mathbf{I}$ are random and $M$ is sufficiently large. However, in a coherent imaging system, speckle noise inherent in the laser and shot noise due to the limited photon measurements are prevalent. They bring errors to the Fourier intensity measurements and corrupt the reconstructed phase. The maximum a-posterior (MAP) method is often used to ensure the performance of the estimation. It is applied for estimating the maximum point of an unknown quantity with the help of a posterior distribution. To further improve the optimization performance, different regularization terms are added to the optimization cost function. In this work, we adopt the total variation (TV) regularization term to smooth out the noise while preserving the image edges [16]. Combining with the MAP formulation, we adopt the following TV-MAP approach [16] for the implementation of the CDP phase retrieval algorithm in this paper:

In Eq.(2), $\nabla$ is the gradient operator, $\|\cdot\|_{1}$ refers to the $\ell_{1}$ norm and DFT refers to the discrete Fourier transform; and $\dot{\mathbf{y}}$ represents the discrete Fourier intensity measurements which are the sampled version of $\mathbf{y}$. Similarly, $\dot{\mathbf{I}}$ and $\dot{\mathbf{U}}$ are the samples of $\mathbf{I}$ and $\mathbf{U}$, respectively. The first term of Eq. (2) is the TV regularization term that has been popularly used in noise removal. It is based on the principle that images with spurious details have higher total variation. Reducing the total variation of the measurement can remove the unwanted details. The second term of Eq. (2) is the MAP data fidelity term that evaluates the performance of the estimation. Minimizing these two terms can obtain an estimate that is close to the original $\dot{\mathbf{U}}$ while reducing the spurious detail due to noise. Furthermore, the optimal solution allows discontinuities along the curves; therefore, edges can be preserved in the restored image. The regularization strength is controlled by the constant $\alpha$. It is usually selected empirically according to the complexity of the signal. Similar to [9, 10], Eq. (2) can be solved via the Alternating Direction Method of Multipliers (ADMM) algorithm [17].

As mentioned in Section 1, binary masks are used in this work. Both SLM and DMD can be used to implement $\mathbf{I}$. When SLM is used, $\mathbf{I}$ can be implemented as a set of binary phase masks such that $\mathbf{I}\in\{1,-1\}$. If DMD is used, $\mathbf{I}$ becomes a set of binary amplitude masks such that $\mathbf{I}\in\{0,1\}$. Figure 2 shows a typical experimental setup for implementing the CDP phase retrieval system. It is also adopted in our experiments. In this optical system, the light from the HeNe laser is projected onto a SLM or DMD board which implements the random masks. The masked light is projected onto an object placed on the image plane. The light then goes through a lens so that the image captured by the CMOS camera (placed on the back focal plane of the lens) is the Fourier transform intensity $\mathbf{y}_{i}$ of the transmittance function of the object $\mathbf{U}$ multiplied with the $i$-th mask $\mathbf{I}_{i}$.

As is indicated in [18, 19], $\mathbf{y}_{i}$ appear as a kind of diffraction pattern as shown in Figure 1. It should be noted that $\mathbf{y}_{i}$ is obtained from the Fourier transform of the multiplication of $\mathbf{I}_{i}$ and $\mathbf{U}$ in spatial domain. This is equivalent to the convolution of $\mathcal{F}(\mathbf{I}_{i})$ and $\mathcal{F}(\mathbf{U})$ in frequency domain. However, $\mathcal{F}(\mathbf{I}_{i})$ is non-bandlimited if $\mathbf{I}_{i}$ is white. Thus, the convolution of $\mathcal{F}(\mathbf{I}_{i})$ and $\mathcal{F}(\mathbf{U})$ is also non-bandlimited. It is the reason why we can find many high-frequency data in Figure 1(a). In this paper, we propose to use the green noise masks that have a ring-shaped frequency spectrum concentrated in the mid-frequency region shown in Figure 3. If the green noise masks are used, the convolution of $\mathcal{F}\left(\mathbf{I}_{i}\right)$ and $\mathcal{F}\left(\mathbf{U}\right)$ will also have the data concentrated in the mid-frequency region presented in Figure 1(b) (since $\mathcal{F}\left(\mathbf{U}\right)$ has large magnitude around zero frequency). It can be seen that the high-frequency components of the $0$th diffraction order are much weaker than those using the white noise masks.

The above property of the green noise masks is important for solving Eq. (2). In general, we always focus on the $0$th diffraction order such that the data outside the $0$th diffraction order are not detected. Even if we also detect the data outside the $0$th diffraction order, they will still be ignored in the phase retrieval optimization process. As mentioned above, the reconstruction algorithm is a discrete process. In each iteration, the samples of $\mathbf{U}$, i.e. $\dot{\mathbf{U}}$, are used to compute ${{\dot{\mathbf{g}}}_{i}=\left|DFT\left({\dot{\mathbf{I}}}_{i}\circ\dot{\mathbf{U}}\right)\right|}^{2}$. Then the MAP loss function evaluates the difference between ${\dot{\mathbf{g}}}_{i}$ and the discrete intensity measurement $\dot{\mathbf{y}}$. The difference is then used to update the estimation of $\dot{\mathbf{U}}$. The sampling period of ${\dot{\mathbf{I}}}_{i}$ and $\dot{\mathbf{U}}$ is $d$, i.e. the distance of the SLM pixels or the micro-display units of the DMD. Therefore, the Nyquist frequency of this discrete system is $\frac{\lambda f}{2d}$ denoted by the red square in Fig 1 (coincides with the boundaries of the $0$th diffraction order circled in yellow with radius $\frac{\lambda f}{2d}$). Based on the sampling theory, ${\dot{\mathbf{g}}}_{i}$ is periodic, and only one period of the replicated spectrum, i.e. from $-\frac{\lambda f}{2d}$ to $\frac{\lambda f}{2d}$ will be used by the algorithm to compare with the same region of ${\dot{\mathbf{y}}}$. Therefore, the high-frequency data of ${\dot{\mathbf{y}}}$ beyond the red square will not be considered in the phase retrieval process. It is equivalent to filtering out these high-frequency data.

As mentioned above, the intensity measurements obtained using the traditional white noise masks will contain many significant high-frequency components. Some of them can even distribute outside the red square shown in Figure 1(a). These data will not be considered in the phase retrieval process and thus degrades the performance. In Figure 1(b), the data of the intensity measurements with the green noise masks concentrates in the mid-frequency region (plus around the zero frequency) and are well within the red square. They will be fully utilized in the phase retrieval process. We will verify the advantages of the green noise masks by the experimental results in Section 4.

The aim of the experiments is to demonstrate the performance achieved by using the proposed green noise-masking scheme in an actual CDP phase retrieval system. It is compared with the traditional white noise masks [9]. The experimental setup comprises a Thorlabs 10mW HeNe laser with wavelength $\lambda=632.8$ nm; and a $12$-bit $1920\times 1200$ Kiralux CMOS camera with pixel pitch $5.86\mu$ m. For the implementation of the noise masks, a $1920\times 1080$ Holoeye Pluto phase-only SLM with pixel pitch $\delta_{SLM}=8\mu m$ is used. Only the central $256\times 256$ SLM pixels are utilized in all experiments. The focal length of the lens behind the SLM is $75mm$. As explained in Section 2, the Nyquist frequency of the system is $\frac{\lambda f}{2d}$ corresponding to the red square in Figure 1. For all experiments, we use the same set of optical devices hence the Nyquist frequency is the same as $\frac{\lambda f}{2d}$. The data of the same rectangular areas are captured for the phase retrieval algorithm using different kinds of random masks. The whole experimental setup has already been shown in Section 2 (Figure 2). Figure 7 shows a photograph of the actual hardware setup.

In the experiments, we also use the SLM to generate the testing objects. Specifically, three images namely, Cameraman, Vortex and Cell (see Figure 8) are loaded to the SLM to generate the phase testing objects that can impose multiple-level phase changes. In this experiment, the binary phase masks $\mathbf{I}\in\{1,-1\}$ are used. The images are pre-multiplied with the random masks before sending to the SLM. Therefore, the output of the SLM has already implemented the function $\mathbf{I}_{i}\circ\mathbf{U}$ in Eq. (1). The green noise masks are generated using the FMEDg method described above. Based on the simulation result described in Section 3.3 and further experimental work, we choose $R_{1}=1.5$ and $\sigma=0.5$ when generating the green noise masks. Only $3$ measurements are taken for each kind of binary phase masks. The number of measurements $M$ does not affect the conclusions of the experiment as long as it is not too small. If it is very small, the phase retrieval algorithm sometimes cannot converge. It is worth noting that the fixed constant $\alpha$ in the optimization algorithm controls the strength of regularization. For different experimental objects, the optimal $\alpha$ are different because of their different complexities (curves, sizes, edges, etc.). In our experiments, we fine-tune the values of $\alpha$ to obtain the minimum error in each case. The experimental results and the corresponding Fourier intensity measurements are presented in Figure 8. It can be clearly seen that the performance of the reconstructed images via the green noise masks are much better than those reconstructed via the white noise masks, both quantitativel y and qualitatively. The improvement in MSE can be up to $47.57\%$ . The improvement on the vortex images is particularly significant. A zoom-in region of the reconstructed cell image is presented to better visualize the performance of different masking schemes. It can be seen that textures and edges of the reconstructed cell image via the green noise masks are much clearer and sharper than those via the white noise masks.

We also use a real object which is a periodic pillar optical grating (shown in Figure 9) fabricated with SU-8 to demonstrate the effectiveness of the proposed green noise masks. The SU-8 is transparent under $632.8$ nm. The period of the grating is $200\ \mu m$ and the height is $d=600$ nm. Hence the phase difference between air and SU-8 area is $\Delta\phi=\frac{2\pi d(n_{1}-n_{2})}{\lambda}\approx 0.98\pi$, where $n_{1}\approx 1.52$ and $n_{2}=1$ represent the refractive index of the SU-8 and air, respectively. We use the SLM in the last experiment to generate the same set of the phase-only random masks. The central $250\times 250$ pixels of the SLM are used in the experiment. Only $1$ measurement is used for each kind of binary masks. Using more masks would not change the conclusion of the experiment.

The experimental results are shown in Figure 9. For the amplitude part, periodic rings appear in the reconstructed images because of the sharp change along the edges of the pillars. As can be seen, the amplitude image via the green noise mask is smoother and more uniform than that of the white noise mask. For the phase part, the phase difference between the pillar areas and the air areas is approximately $0.98\pi$ as stated in the above paragraph. It is clearly shown that the reconstructed image with the green noise mask can better reconstruct the phase part while there exists phase errors via the white noise mask. The phase errors mainly locate close to the edges of the pillar, which means that the lost high-frequency information degrades the performance. From the above experimental results, we have shown that the proposed green noise mask outperforms the traditional white noise mask both quantitatively and qualitatively when using in SLM-based phase retrieval systems.

To demonstrate the generality of the proposed green noise masking scheme, we also apply the scheme to a DMD-based CDP phase retrieval system. Different from the experimental setup in Figure 7, a $1024\times 768$ pixels DMD with pixel pitch $13.68\mu m$ installed on a Discover $1100$ DLP board is used to generate the random masks. The other parts of the setup can be found in Figure 2. Figure 10 shows a photograph of the actual hardware setup. The optical path includes a reflective mirror (placed in the bottom-right corner) is used to match the blazed grating angle of the DMD board. Several reflective mirrors are also used in the laser beam magnification path because of the limited size of the optical table. In this experiment, the binary amplitude masks $\mathbf{I}\in\{0,1\}$ are loaded to the DMD board. A Newport USAF chart is placed on the object plane as the testing object. The projected light can only pass through the chart through its holes thus introduces sharp amplitude changes along the boundaries of the holes. The light that can pass through the chart should have the same phase, and we do not expect there is any variation in phase for the regions with no light. For the experiment, $3$ measurements are captured by the camera for each kind of binary masks. The green noise masks are generated using the same $R_{1}$ and $\sigma$ as in the SLM experiments. Similar to the SLM experiment, the number of measurements chosen does not affect the conclusions of the experiment as long as it is not too small. They are fed to the same TV-MAP algorithm to retrieve the amplitude and phase of the USAF chart. And we fine-tune the values of $\alpha$ to achieve the minimum error in each case. Figure 11 shows the experimental results and the corresponding Fourier intensity measurements with different binary masks. It can be seen that the amplitude and phase images retrieved through the use of different kinds of binary masks have clear differences. Even though the background areas of the amplitude images are smooth for all kinds of binary masks, the quality of the images in the regions where light can pass through has a considerable difference. The amplitude image resulting from using the green noise masks remains sharp and relatively uniform. Comparing with the ground truth, the $SSE(dB)$ of the amplitude images with the green and white noise masks are $-9.400$ dB and $-3.722$ dB, respectively. As to the phase images, the one reconstructed with the green noise masks has a much smoother background and fewer errors compared with the white noise masks.

This paper proposes a green noise binary masking scheme for the coded diffraction pattern (CDP) phase retrieval systems. We have explained why the high-frequency data in the Fourier intensity measurements will be ignored in the CDP phase retrieval process. To reduce the information loss, we suggested using green noise masks that have very low high-frequency concentrations. The green noise masks are designed using a digital halftoning method we developed earlier. It has the advantage that the principal frequency and the density of the noise mask can be controlled independently. Also, we have demonstrated how the high-frequency data in the Fourier intensity measurements can affect the quality of the reconstruction images. It is verified by the experimental results, which show that the quality of the amplitude and phase images retrieved using the proposed green noise binary masking scheme outperforms the commonly used white noise binary masks when applying to two CDP phase retrieval systems.

This work was supported by the Hong Kong Research Grant Council under General Research Fund no. PolyU 152478/16E.

The authors declare no conflicts of interest.

Qiuliang Ye: Conceptualization, Methodology, Simulation, Experiment, Writing - original draft, Writing - review & editing. Yuk-Hee Chan: Conceptualization, Methodology. Michael G. Somekh: Supervision, Writing - review & editing. Daniel P.K. Lun: Supervision, Project administration, Funding acquisition, Conceptualization, Methodology, Writing - review & editing.