Digital holographic phase imaging with aberrations totally compensated

: Digital holography is a well-accepted method for phase imaging. However, the phase of the object is always embedded in aberrations. Here, we present a digital holographic phase imaging with the aberrations fully compensated, including the high order aberrations. Instead of using pre-defined aberration models or 2D fitting, we used the simpler and more flexible 1D fitting. Although it is 1D fitting, data across the whole plane could be used. Theoretically, all types of aberrations can be compensated with this method. Experimental results show that the aberrations have been fully compensated and the pure object phase can be obtained for further studies.


Introduction
Digital holographic phase imaging (DHPI) is widely used to quantitatively detect phase information with axial sensitivity in nanoscale [1,2]. With this technique, transparent samples such as cells can be visualized [1,3], without contrast agents (stains or fluorescence dyes) and in a completely noninvasive way. Thus, it has promising applications in the imaging of living cells, such as the monitoring of cell growth or the membrane fluctuations [4]. However, the phase of object is usually embedded in aberrations, which distort the object phase and prohibit accurate measurements unless they are fully compensated.
Methods for aberration compensation can be categorized as physical ones [5] and numerical ones [6][7][8][9]. Physical methods include, for example, using the exact same MO in both reference and object arm [10], using an telecentric imaging system [11][12][13], or taking extra object-free holograms [5,14]. Numerical methods, on the other hand, are typically implemented by quantifying and removing aberrations during the digital reconstruction process, and thus offer more flexibility.
In terms of numerical methods, compared to the first and second order aberrations (tilt and defocusing), higher order aberrations are much more difficult to compensate. Thus, some methods only compensate for the low order aberrations, leaving the higher order aberrations unresolved [15][16][17]. This will be not ideal if the system has high order aberrations.
One solution is two-dimensional (2D) surface fitting using pre-defined (high order aberrations included) aberration models [18,19], with the most acceptable technique being the 2D Zernike Polynomials fitting (ZPF) [18]. Many of these techniques, however, perform 2D fitting across the entire field of view (FOV) without considering the presence of object phase [18,19], leading to somewhat suboptimal results.
Compared to 2D surface fitting, one-dimensional standard polynomials fitting (1DSPF) is more flexible, powerful and computation efficient, and theoretically, any type of line profile can be represented with 1D standard polynomial. Nevertheless, as pointed out by Colomb, T. et al. [20], traditional 1DSPF cannot characterize aberrations with cross terms, such as xy and x 2 y.
In summary, although many efforts have been devoted to compensate the phase aberrations in DHPI systems, most of them request prior knowledge on the aberration models, or are incapable of correcting high order aberrations especially those with cross terms, or negatively affected by the presence of object phase. DHPI system with aberrations totally compensated without the influence of object phase has not yet been demonstrated to our knowledge. Here we report a DHPI system that permits the wide-field phase image of object, with aberrations be totally compensated. Neither information about optical system nor any aberration models need to be known in advance. Experimental results show that the phase aberrations can indeed be fully compensated and the pure phase of object can be obtained for further analysis.

Principles
In this section, we first give a brief introduction on why aberrations always exist in the DHPI systems, and how they affect the phase measurement. Then, we give a detailed description on the aberration compensation procedure. Finally, its ability of aberration compensation is analyzed.

Interference and aberrations
In holographic systems, the interference pattern between object wave O and reference wave R can be described as where (·)* denotes the complex conjugation. The first two terms in Eq. (1) are zero terms and the others are cross terms. Either of the cross terms can be used to reconstruct the image of object. Here we use OR*. The zero and twin terms can be removed with phase shifting [21] or off-axis techniques [11]. After extracting OR* from the holograms, it is numerically propagated to the image plane using angular spectrum method [22]. In the image plane, the object wave O i and reference wave R i can be defined as O i = Aexp(iφ O ) and R i = Bexp(iφ R ) where A and B are the amplitude while φ O and φ R are the phase of the object and reference wave, respectively. The object phase φ O is the quantity of interest.
If there exists optical defects, the object wave and the reference wave will be altered as O ia = A′exp(iφ O )exp(iφ Oa ), R ia = B′exp(iφ R )exp(iφ Ra ) where φ Oa and φ Ra are the optical defects induced phase aberration. Thus the reconstructed wavefront becomes: where φ ta is the total phase aberration and φ ta = φ Oa -φ R -φ Ra . It can be seen that the phase of object (φ O ) is buried in φ ta , which is consisted of the phase of reference wave (φ R ) and the phase aberrations caused by the imperfections of system (φ Oa and φ Ra ). The phase of object will be unaffected if the total phase aberration φ ta is constant, which is impractical in real case. Because even if the optical system has no defects, which means φ Oa and φ Ra are all constant, the object phase is still affected by φ R . For example, in the widely used off-axis holographic system, the angle between the reference wave and the object wave induces tilt into φ R , which causes tilt aberration. Moreover, if a MO is used to magnify the object, quadratic phase is introduced into φ Oa , causing defocusing aberration. Aberrations will be more sophisticated if the system is imperfect, i.e. it has optical defects.
It can be seen that in real case, it is almost unavoidable that the object phase is affected by various kinds of aberrations. Therefore, for an accurate measurement of the object phase, the total phase aberration φ ta must be fully quantified and compensated.

Fitting procedure and aberration compensation ability
It is easy to find out the amount of aberrations in the background and then compensate them. However, the work is difficult for the object areas because in these areas the object phase are messed up with aberrations. Thus, we should use the aberration data in the background area to estimate that in the object area. When estimating aberrations, the object areas should be removed manually [20] or automatically with background detection methods [12], to get rid of their negative influence. Therefore, here, after unwrapping the reconstructed phase image [23], the background without object area is extracted to obtain a 2D phase distribution φ that will be used to quantify the phase aberrations. If aberrations have been fully compensated, the background will be flat or its standard deviation (STD) will be very small [20]. The process of aberration compensation is given in Fig. 1. Here we use φ re and φ q to represent the "residual phase aberration" and the "quantified phase aberration", respectively. Initially, φ re equals to φ. We calculate the STD of φ re . If STD is larger than a threshold value t, we find the direction u along which the STD is maximum using grid search. The grid search is realized by calculating the STDs of line profiles separated by 2°. After this, extract the profiles along direction u and its normal direction v. Then using 1DSPF, two extracted profiles are fitted to functions f u and f v , respectively, as: where a n and b n (n = 0, 1, 2, 3…) are the polynomial coefficients of order n. The order of polynomial function n can be set as needed. Here, we limit the correction to fourth order (n = 4) as it already includes nearly all of the ordinary aberrations encountered in an optical system. Assume the angle between uv and the principle axes xy is θ. According to the coordinates transform, the quantified 2D phase aberration (φ q ) can be calculated as: 4 4 where x, y are the coordinates of the entire FOV (including the object areas), P k,l is the coefficient of aberration x k y l . The relationship between aberration coefficients P k,l and polynomial coefficients (a n and b n ) is given in Table 1. Table 1. Relationship between the aberration coefficients P k,l and the polynomial coefficients (a n and b n )  Table 1, it should be noted that, the polynomial coefficient a n only influences the coefficients P k,l of aberrations x k y l where k + l = n. For example, a 2 only influences the coefficients of aberrations x 2 , y 2 and xy while the other aberrations are not related. On the other hand, to quantify coefficient P k,l , polynomial coefficients a n and b n where n = k + l are both needed. For example, to obtain the coefficient of aberration xy 2 , polynomial coefficients a 3 and b 3 are both needed.

Coefficients Aberrations Description Expression
It can be seen that we use the fitted polynomial coefficients, a n and b n , to calculate the aberration coefficients P k,l according to Table 1. Then use them to calculate the quantified 2D phase aberration (φ q ) according to Eq. (4). Using Eq. (4), aberrations across the whole plane, including that in the object areas, are obtained and compensated as O ia R ia *exp(-iφ q ). After this, the residual phase aberration φ re changes to φ re = φ re -φ q . The φ re will be fitted again and the whole procedure will be repeated, until its STD is smaller than the threshold t.
From Eq. (4), we can see that all kinds of aberration (limited to n ≤ 4) have been quantified including those high order aberrations and those with cross terms. Higher order aberrations can also be included with larger n. As we mentioned before, the order of polynomial fu order aberrati highest order

Experime
Experiments w the schematic 2 (a)) and cs, GCGge (a), the umber "3" even after phase unwrapping. That is because the phase of object was embedded in large amount of aberrations. To obtain a pure phase image of object, aberrations have to be removed. We used the procedure described in Section 2.2 to quantify and remove phase aberrations present in Fig. 3(c). The object area was manually excluded based on the intensity image in Fig. 3(a). We set the STD threshold t = 0.46 and the algorithm terminated after 10 iterations. The residual aberrations φ re (left) and the quantified phase aberrations φ q (right) for iteration 1 to 7 are shown in Fig. 4. The residual aberrations are presented with object areas un-removed and in wrapped status for clear observation. In Fig. 4, we can see that the aberrations have been gradually reduced with increasing number of iterations. For example, there are ~5 periods of wrapped aberration (-π to π) in Fig. 4(d), whereas in Fig. 4(g), it is reduced to less than 1 period. Meanwhile from the color bars of the ri smaller with the STDs of Generally spe All of these p We note 226.39 to 7.9 defocusing) t sophisticated compensation Iteration 0 Directions (degrees) STD 22 Table 2 al are also indic searches the Theoretically, directions ma directions.

Transmis
Here we show imaging of bi breast cancer Stem Cell Ba Eagle's mediu (FBS, Hyclon to trypsin (0. ight figure, it iteration, also f residual pha eaking, in Tabl rove that the ab that, in Table  5. This is beca that are easy and hard to n. the phase nd human ovided by modified ine serum e exposed n DMEM medium to final density of 5 × 10 4 cells/mL, and seeded on cover slips. After 6h incubation at 37°C, the samples were fixed in 4% paraformaldehyde for 4h at 4°C。 The intensity image of the C2C12 cells is shown in Fig. 6(a). Since these cells are mostly transparent, they are hardly visible in the intensity image as indicated by the yellow arrow (except when out-of-focus as marked by the red arrow). Phase imaging is an ideal way to visualize such samples due to index differences between the cells and their surroundings. However, because of the aberrations, the phase of the cells is difficult to observe as in Fig.  6(b).
We then ran our algorithm for 2 iterations, and the compensation result is shown in Fig.  6(c). Meanwhile, the area near the yellow arrow in Fig. 6(c) is enlarged and shown in Fig.  6(d). We observe that after compensation, the cell structures are much more visible with a spindle-like morphology. Their structures can be clearly observed including the nucleus and their inner structure. The flat background suggests that the aberrations have been compensated. Quantitatively, the STD of the residual phase was reduced from 298.85 to 0.36. Figure 7 shows the phase image of MCF-7 cells after aberration compensation with conventional 2D ZPF (a) and our method (b). The STDs of residual aberrations on the background of Fig. 7(a) and Fig. 7(b) are 0.36 and 0.33, respectively. The area near the red arrow in Fig. 7 (b) is enlarged and showed in Fig. 7(c) with pseudo three-dimensional display. Figure 7(d) shows the reversed phase data of (a) and (b) along the yellow line profile in Fig.  7(a).
In Fig. 7(a) and (b), we observe that both methods can remove the aberration and the MCF-7 cells show a typical epithelial like shape. However, the 2D ZPF presents a small amount white region around each cell (indicated by green arrows in Fig. 7(a)). This is more obvious in Fig. 7(d), in which as the green arrows indicated, there are singularities around cells. It is possibly caused by the negative effect of object phase, so the quantified phase aberration is not entirely accurate. Such artifact is not present with our method. From the experiments in Fig. 6 and Fig. 7, it is can be seen that our method can compensate the total aberrations presented in the DHPI of cells. The nucleus and inner structures can be clearly observed for further studies.

Summary
In conclusion, we use 1DSPF to compensate the aberrations existed in DHPI system. The 1DSPF is adaptively implemented along direction with the maximum aberration and its vertical direction. Then using the relationship between the 1DSPF coefficients and the coefficients of aberrations, the quantified 2D aberration is obtained and removed. The residual aberrations will be fitted again. This process is executed iteratively until the residual aberration is small enough. Additionally, the influence of object phase can be removed. Experimental results show that this method can compensate the total phase aberration, including those high order aberrations with cross terms. Neither pre-defined aberration models nor priori information about optical system are required.
This method solves the tough and unavoidable phase aberration problem existed in DHPI systems. It is powerful and aberration can indeed be fully and effectively removed, advancing the development and applications of DHPI.
We should point out that here we only demonstrated the aberration compensation on a single plane. In the case of 3D objects where phase aberrations across different objet layers are different, for aberrations on each object layer, they can be individually compensated with this method. Moreover, this method requires that samples should have some flat background.