1. Introduction

Fourier ptychography (FP) is an emerging microscopy [1] which borrows ideas from synthetic aperture [2,3] and phase retrieval [4,5] concepts. In the FPM system, the traditional light source is replaced by a programmable LED array to provide angle-varied illumination. A series of low-resolution (LR) images with different incident illumination wave vectors are recorded by an imaging sensor and employed to sequentially update the corresponding sub-aperture information. In this way, oblique illumination realizes aperture scanning in the Fourier domain and expands the spectral range of objective lens acquisition, thus FPM can bypass the trade-off between spatial resolution and FOV in conventional optical systems.

In typical FPM platforms, the systematic errors, including inherent aberrations, inevitable noises, and model misalignment [6], are mixed in LR measurements and affect the quality of reconstructed HR images. For inherent aberrations, the embedded pupil function recovery (EPRY) algorithm has been proposed to recover object spectrum and pupil function simultaneously [7,8]. For inevitable noises, the noise effect can be alleviated via adjusting acquisition exposure time and conducting some noise removal operations [9,10]. Nevertheless, the model misalignment during the reconstruction process has not been properly addressed. The spatial positional deviation of the LED element is converted to the incorrect location of sub-aperture in the Fourier domain, resulting in error accumulation with iterations and corrupting the reconstruction consequently [11].

The model misalignment calibration methods can be divided into two categories: algorithmic self-calibration (ASC) and physics-based methods. These two methods can also integrate with a deep-learning-based FP architecture [12–15] to implement the misalignment error minimization. Algorithmic self-calibration approaches [6,16–20] involve jointly solving two inverse problems: FPM reconstruction and systematic parameter calibration. It generally re-evaluates the difference between the captured and updated LR images at each iteration. In this respect, Sun et al. [17] first introduced the global positional misalignment (GPM) model to characterize the illumination wave vectors. Pan et al. [6] proposed a system calibration procedure termed SC-FPM to iteratively optimize the sub-aperture positions, while the simulated annealing (SA) algorithm, LED intensity correction module, and adaptive step-size strategy are integrated into the calibration procedure. Although the ASC methods eliminate the demand for pre-knowledge and test ground truth, the LR intensity measurements still contain a mixture of information from the sample and various errors. These ASC methods are mostly based on the GPM model, however, their effectiveness is limited by the optimization algorithms, resulting in a tendency to reach locally optimal solutions rather than globally optimal results [20].

Physics-based methods [21–23] intend to separate the systematic parameter calibration from the reconstruction process and localize the precise sub-aperture position through the physical model. As such, Regina et al. [21] localized bright-field (BF) illumination angle through the auto-correlation spectrum of LR intensity image. Zhang et al. [22] achieved BF alignment using bright-dark-field (BDF) transition images acquired without placing the sample. Zheng et al. [23] inferred misalignment parameters by analyzing the shifted BF position of a defocused object. The first method [21] utilizes an additional $2\times$ digital camera adapter to reduce the size of sampling aperture, ensuring all BF Fourier spectra are within the LR image size. The second method [22] requires attention to the adjustment of the distance from the illumination source to the sample to ensure the presence of the BDF transition images [24]. The defocused strategy [23] even necessitates a high-precision focusing knob. In addition, the calibration process primarily takes place in the BF area and does not include the dark-field (DF) area. While extracting information from images in the BF area is easier due to their higher intensity, the positional parameters in BF don’t accurately represent the real situation in DF.

In this paper, we report an independent model misalignment correction strategy based on the transformed pupil function, termed TPF-FPM, achieving global position calibration through the localization of DF illumination angles. First, we deduce the connection between the frequency positional deviation and pupil function for individual DF illumination during the iterative path of EPRY. Second, based on the fact that the frequency positional deviations of sub-apertures form certain distribution characteristics under different misalignment situations, we propose a judgment criterion to determine the search directions of the shifting pupil centers. Then, the offset pupil centers corresponding to corner LED elements can be precisely identified through a boundary-finding scheme and the global positional parameters can be further calculated. In TPF-FPM, using the amplitude of pupil function as the calibration reference can prevent interference from aberrations in the phase. Besides, the pupil function can also serve as an evaluation indicator to visualize whether the correct calibration has been performed in the DF area. Simulations and experiments demonstrate that TPF-FPM can precisely calibrate model misalignment of LED array and the reconstructed results are superior to the ASC method. The proposed pupil-function-based strategy allows for accurate DF localization using a physical model and promises to be applied in the precise positioning correction in conventional ptychography [25,26].

2. Principle

2.1 Global positional misalignment model

The forward imaging model of FPM is shown in Fig. 1, an LED array with quasi-monochromatic wavelength $\lambda$ illumination is placed far enough beneath the sample. The positions of LED elements in the LED array can be characterized using a four-dimensional model $(\Delta x, \Delta y, \theta, h)$ [17]. Here, $\Delta x$ and $\Delta y$ represent the lateral offsets of the LED array along the $x$ and $y$ directions, respectively. $\theta$ denotes the rotation angle at the horizontal plane and $h$ is the vertical distance of the LED array from the sample.

 

Fig. 1. The forward imaging model of FPM with positional misalignment. The spatial positional offset of the LED element leads to the change of wave vector and results in the shifting of sub-aperture position in the Fourier domain. The wave vector of the central LED element (green circle) is labeled with a red arrow.

Download Full Size | PDF

Assuming the global offset is larger than displacement errors between individual LEDs, the spacing between adjacent LEDs can be considered equal, denoted by $d$. For the $l^{th}$ LED element which locates on $m$ row and $n$ column, its spatial position $(x_l,y_l)$ can be expressed as

Suppose the sample is a thin object and far away from the LED array, the illumination beam emitted from the LED element to the sample can be regarded as a plane wave. For the central segment of the sample, the incident wave vector of the $l^{th}$ LED element can be formulated as

In the process of image acquisition, LED elements are turned on sequentially to complete the scanning of the sample’s spectrum. However, if there is any spatial deviation in the LED element’s position, it may cause changes in sub-spectral zones during the aperture scanning, which can lead to undesired reconstruction.

2.2 Transformation between frequency deviation and pupil function

The aperture scanning process can be understood as moving the coherent transfer function(CTF) of the microscope system on the sample spectrum plane. Illuminating the sample by the $l^{th}$ LED element with the wave vector $(k_{xl},k_{yl})$ is equivalent to shifting the center of CTF by $(-k_{xl},-k_{yl})$ in the Fourier domain. The corresponding sub-spectrum information is low-pass filtered by the CTF. As a result, the captured LR intensity measurement $I_{l}^{c}(\boldsymbol {r})$ can be expressed as

The Eq. (3) demonstrates that the function of the illumination wave vector $\boldsymbol {k}_l$ can be included in the expression of pupil function $P(\boldsymbol {k})$ according to the Fourier shift theorem. In the conventional FPM reconstruction algorithm, the pupil function can be recovered using an EPRY iterative path. At the start of the recovery process, an initial HR sample spectrum solution denoted by $O_0(\boldsymbol {k})$ can be generated by interpolating the amplitude of the central LR image and setting the phase to zero. The initial guess of pupil function $P_0(\boldsymbol {k})$ is often set as a circular low-pass filter with all ones inside the objective passband. Corresponding to the particular sub-aperture position, the estimated LR complex amplitude field $\varphi _l^e(\boldsymbol {r})$ is simulated as

During iterations, the amplitude component of $\varphi _l^e(\boldsymbol {r})$ is replaced by the square root of the actual LR measurement $I_l^c(\boldsymbol {r})$ to acquire the updated complex amplitude field $\varphi _l^u(\boldsymbol {r})$

In FPM, the complex-valued reconstruction of HR sample spectrum is iteratively achieved by enforcing two constraints on the updated quantities. The captured LR intensity images are utilized as the amplitude constraint, whilst the support constraint is imposed by the pupil function in the Fourier domain. In cases when these two constraints are defined, the valid solution to Eq. (3) can be enclosed via the EPRY updating route

Due to the model misalignment in practice, varying incident illumination wave vector deviations lead to diverse offsets of sub-aperture positions. For the $l^{th}$ illumination wave, its wave vector deviation can be expressed as

 

Fig. 2. The schematic of the principle of offset pupil. (a) The captured frequency information during FPM forward imaging; (b) Recovered offset pupil due to frequency information correspondence.

Download Full Size | PDF

Meanwhile, the intensity of the deviated pupil is not the same, depending on the corresponding sub-spectral energy. In the field of conventional ptychography, Maiden et al. [27] have analyzed the update function of embedded ptychographical iterative engine (ePIE), which is homologous to EPRY. The object update of ePIE arises from a weighting equal to the normalized intensity of the probe, i.e., the good update in places where the probe is bright should be accepted, whereas the previous object estimate should be retained for the dim probe. In FPM, the probe update of EPRY shares similar regions of the weight function since the object that exists in the Fourier domain has a dramatically high dynamic range, which is analogous to the probe in ptychography. To illustrate the weight function, the Eq. (9) can be formulated as

During the implementation of pupil function recovery, DF images used in reconstruction are essential in addition to ensuring sub-apertures overlap. However, the positional deviations of BF sub-apertures will not result in offset pupils as BF images don’t contribute to the recovery of pupil function in the framework of EPRY [6,13], making the transformation in Eq. (9) only applicable to the DF area. To illustrate this difference, all captured LR images are used for reconstruction and the model misalignment is imposed on different regions of the LED array. In this simulation, the spatial translation offset ($\Delta x=-1$ mm) is progressively added from the central to external LED rings across different situations and the corresponding recovered pupil functions are shown in Fig. 3. The system parameters here are the same as we set in the subsequent relevant simulations and experiments. The light source is a $13\times 13$ LED array (spacing 8.128 mm) with a central wavelength 505 nm and placed at 98 mm beneath the sample plane. The NA and magnification of the objective lens are 0.13 and 4, respectively. The pixel size of the imaging sensor is 6.5 $\mu$m. It’s apparent from Figs. 3(b)-(c) that the pupil function remains almost unchanged compared to the situation without model misalignment in Fig. 3(a), indicating that the sub-aperture positional deviation in BF area does not conform to the above transformation. In contrast, the circular outlines beyond the passband become brighter as incorrect-located DF sub-apertures increase under the same number of iterations, as shown in Figs. 3(d)-(h).

 

Fig. 3. The amplitude of recovered pupil function by EPRY and 10 iterations. (a) and (b)-(h) are results without and with model misalignment ($\Delta x=-1$ mm), respectively; the number of LEDs with translation offset in (b)-(h) are 1, $3\times 3$, $5\times 5$, $7\times 7$, $9\times 9$, $11\times 11$ and $13\times 13$, respectively.

Download Full Size | PDF

In general, areas with more overlapping pupils tend to have higher intensity. However, due to the retention of the CTF constraint during the recovery of the object spectrum, the performance of pupil function differs inside and outside the objective passband, as shown in Figs. 3(d)-(h). The pupil function inside the passband is updated simultaneously with the object spectrum, whereas outside the passband, it is solely superimposed in intensity with iterations. As the number of incorrectly located DF sub-apertures and iterations increases, the contrast inside and outside the passband will gradually become more distinct, resulting in a darker appearance inside the passband. To balance the overall contrast and facilitate the extraction of shifted pupil profiles, it’s not necessary to perform plenty of iterations.

2.3 Misalignment situation determination

The mapping relationship between the transformed pupil function and the DF frequency positional deviation has been deduced and then the distribution of deviated pupil centers is further explored based on the spatial arrangement of the LED array. Despite the degradation of reconstruction quality attributed to the combined effect of multiple misalignment parameters, the misalignment parameters affect wave vectors in different ways. In the following simulations, the main focus will be on the translation and rotation of the LED array. The frequency positional deviation caused by height factor $h$ is very small as it will not significantly change during the focusing process for different samples once the working distance of microscopy is determined.

Figures 4(a1)-(a2) show the ideal and actual frequency positions of sub-apertures when translation ($\Delta x=\Delta y=1$ mm) or rotation ($\theta =2^{\circ }$) appears independently. Figures 4(b1)-(b2) are the distribution maps of the absolute value of frequency positional deviation. The centers of shifted pupils are displayed as solid shapes in Figs. 4(c1)-(c2). Since the corner LED elements are furthest from the center of the LED array and more susceptible to severe deviation, they are marked with diamond shapes for subsequent analysis. As shown in Fig. 4(c1), in the case of translation, pupil centers are located in the same quadrant. The right-down corner LED element acquires the smallest offset, and its corresponding shifted pupil constructs a portion of the entire pupil function boundary in the translation direction. As for the rotation situation in Fig. 4(c2), the corresponding shifted pupils of four corner LED elements compose a part of the entire pupil function boundary in the approximate diagonal direction.

 

Fig. 4. The influence of translation($\Delta x=\Delta y=1$ mm) and rotation($\theta =2^{\circ }$) on incident wave-vectors. (a) The ideal and actual frequency positions; (b) The distribution maps of the absolute value of frequency deviations; (c) The distributions of actual frequency deviations, where solid shapes indicate the centers of DF pupils.

Download Full Size | PDF

According to the above analysis, the boundary of pupil function in certain directions can indicate the centers of deviated pupils corresponding to particular LED elements. It’s a pure circular boundary-finding issue and can be addressed by converting the two-dimensional frequency coordinates $\Delta \boldsymbol {k}_l=(\Delta k_{xl},\Delta k_{yl})$ to polar coordinates $\Delta \boldsymbol {k}_l=(\rho _l,\theta _l)$ and implementing the relevant search algorithm [21]. Since the passband of the objective lens is known and unchanged, a sharp decline will appear at the passband radius $R_{obj}$ along several radial lines out from the correct center $\Delta \boldsymbol {k}_l$. This boundary-finding method is effective for the non-overlapping portion of the circle. Nevertheless, the recovered pupil function is the result of multiple pupils overlapping, so that only those pupils whose centers are at the outermost part of the distribution maps can be localized. Thus an extra search direction constraint termed Pupil search direction criterion is necessary to define the search range and precisely localize the centers of particular shifted pupils. In this criterion, the distribution characteristics of shifted pupil centers are divided into three categories:

i) Rotation only. The impact of rotation on the wave vectors is reversed for centrosymmetric illumination positions. As shown in Fig. 5(a), $R_i$ is approximately equal at the dramatic drop along any two radial lines out from the image’s center in complementary directions. As a result, the four diagonal directions can be used to judge misalignment cases first. If $R_1\approx R_2$ and $R_3\approx R_4$ are satisfied, it indicates that the misalignment parameters only involve $\theta$, and we refer to this case as "Rotation only". In this situation, four diagonal directions are employed as initial directions to search the shifted pupil centers of four corner LED elements. By reducing the number of LR images used for reconstruction, 16 DF illumination angle deviations can be obtained and their corresponding locations are labeled with solid diamonds in Fig. 5(b).

 

Fig. 5. The schematic of search directions determination when "Rotation only". (a) The pupil function recovered by 169 LR images with model misalignment $\theta =2^{\circ }$; (b) Pupil centers in which the deviation directions can be determined.

Download Full Size | PDF

ii) Translation dominates. Taking further observation on Figs. 4(b1)-(b2), the wave vector deviations decrease from inside to outside of the LED array in the case of translation. On the contrary, the DF sub-apertures suffer from relatively bigger frequency positional deviations in the case of rotation. Therefore, it’s a complicated situation when translation and rotation occur simultaneously, which is more in line with practice.

Spatial translation inevitably leads to an asymmetric offset of the entire pupil function. When the judgment condition of the previous case is not met, an interval $\alpha =[0,\pm 180^{\circ }]$ is employed to search the directions $[\alpha _1,\alpha _2]$ corresponding to the $[R_{max}, R_{min}]$ at the sharp fall of the radial lines. A series of judgments around $[\alpha _1,\alpha _2]$ will be further implemented to determine whether the effect of rotation can be ignored. As shown in Fig. 6(e), when the effect of rotation is much smaller than that of translation, the shifted pupils corresponding to the fourth ring LED elements form the profile of pupil function in the $\alpha _1$ direction. The pupil functions reconstructed from different numbers of LR images are shown in Figs. 6(a)-(d). If $R_{max}$ along $\alpha _1$ direction does not satisfy the stepwise decreasing rule, it means "Translation dominates". At this point, we can only consider the effect of translation and utilize $\alpha _2$ as the initial search direction. Because both the distribution characteristics of shifted pupil centers in the cases of "Translation dominates" and "Translation only" follow a gradually decreasing rule in $\alpha _2$ direction, as shown in Figs. 6(e) and (f). Finally, 4 DF illumination angle deviations can be obtained and their corresponding locations are labeled with solid diamonds in Fig. 6(f).

 

Fig. 6. The schematic of search directions determination when "Translation dominates". (a)-(d) The pupil functions recovered by 169, 121, 81, and 49 LR images with model misalignment ($\Delta x=\Delta y=1$ mm, $\theta =0.2^{\circ }$); (e) Distribution characteristics of shifted pupil centers; (f) Pupil centers in which the deviation directions can be determined.

Download Full Size | PDF

iii) Rotation dominates. As the rotation degree increases and its effect exceeds that of translation, the extent of both is worthy of attention. The $R_{max}$ will gradually decrease when using the second judging condition, as can be seen in the distribution characteristics of pupil centers corresponding to the right-top corner LED element in Fig. 7(e). To this end, the radial line direction $\alpha _1$ corresponding to $R_{max}$ can be utilized. Whether $\alpha _2$ is a viable search direction depends on the clarity of the boundary in that direction. Taking Fig. 7(e) as an example, if the deviation of the left-down corner LED element is minimal, the boundary of pupil function may coincide with the passband edge in $\alpha _2$ direction, which is not conducive to segmenting effective information. For pupil functions reconstructed by different numbers of LR images, the respective $\alpha _1$ and $\alpha _2$ need to be acquired through the search process described above. Eventually, $4\sim 8$ DF illumination angle deviations can be obtained.

 

Fig. 7. The schematic of search direction determination when "Rotation dominates". (a)-(d) The pupil functions recovered by 169, 121, 81, and 49 LR images with model misalignment ($\Delta x=\Delta y=1$ mm, $\theta =1^{\circ }$); (e) Distribution characteristics of shifted pupil centers; (f) The diagram of the decomposition of search direction.

Download Full Size | PDF

Despite the locations of the other two corner LED elements (left-top and right-down) are equally prominent compared to the surroundings in Fig. 7(e), their radial line directions cannot be determined through a simple interval search. To provide a denser global position correction constraint, we propose an optional approximate direction determination method. Reordering items in the formula of wave vector deviation and neglecting the correlation coefficient term $(-1/\lambda )$, the formation of wave vector deviation at each illumination position can be regarded as a vector addition operation of the independent effects of translation and rotation. When the ideal location $(d\sqrt {m^2+n^2})$ of LED element is sufficiently large compared to the translation degree $(\sqrt {\Delta x^2+\Delta y^2})$ in the spatial domain, the assumption that the translation effect in the denominator item can be ignored is valid

As shown in Fig. 7(f), for four corner LED elements equidistant from the center, the rotation effect can be approximately described as unit vectors along diagonal directions, marked with orange arrows. The translation vector $\beta =(u,v)$ is also unique to the four LED elements and is indicated by a blue arrow, which characterizes the effect relative to the rotation vector and can be calculated by the connection between $\alpha _1$ and $\alpha _2$ (marked with red arrows)

Two solutions in terms of $u$ could be obtained via Eq. (15). If the solution is smaller than 1, it indicates that the translation effect is not considerable compared to the rotation effect, as it will not cause $R_{min}$ to be less than $R_{obj}$. Another solution greater than 1 describes a larger translation vector that causes the reversal of deviation direction for the left-down LED element in Fig. 7(e). We can choose the appropriate solution based on the comparison between $R_{min}$ and $R_{obj}$. Under the premise that translation vector $\beta =(u,v)$ is obtained in the above analysis, the remaining directions $\alpha _3$ and $\alpha _4$ can be further calculated through the vector addition rule

After obtaining the positional deviations of particular LED elements, the exact sub-aperture positions can be further identified since the locations (row $m$, column $n$) of these elements are known in the LED array. However, determining their spatial locations in rotation cases is challenging because the pupil functions are identical when having the same angular magnitude and opposite direction. The main distinguishing factor is that the distribution characteristics of pupil centers are substantially symmetrical to each other. In this regard, we present a supplementary criterion termed Half-sampling criterion to determine the direction of rotation, as shown in Fig. 8. We can use only half of the LR images for reconstruction and the green circles in Fig. 8(a) represent the lighted LED elements. It’s obvious in Figs. 8(c1)-(c2) that the pupil centers’ distributions of two rotation cases ($\theta =2^{\circ }$ and $\theta =-2^{\circ }$) are completely opposite. When the rotation direction is positive, the pupil function retains only its upper half, otherwise, it’s the lower half. As shown in Figs. 8(d1)-(d2), based on the four initial search directions obtained from previous judgments, the $R_i$ at the sharp drop is identified again. There must be a significant decrease in $R_i$ along two directions, highlighted with orange circles. Observe whether these two directions lie in the interval $[0,180^{\circ }]$ or $[0,-180^{\circ }]$. If they belong to the former interval $[0,180^{\circ }]$, it proves that the rotation direction is positive, and the contrary is negative.

 

Fig. 8. The schematic of the Half-sampling criterion which determines the direction of in-plane rotation. (a) The illumination pattern of LED array; (b) The recovered spectrum; (c) Distribution characteristics of offset pupil centers in the cases of $\theta =2^{\circ }$ and $\theta =-2^{\circ }$, respectively; (d) Recovered pupil function, while the orange circles indicate the decreasing of $R_i$ along the radial line in that direction.

Download Full Size | PDF

2.4 Positional misalignment correction strategy

The process of TPF-FPM is shown in Fig. 9 and can be divided into the following steps:

 

Fig. 9. The flowchart of TPF-FPM.

Download Full Size | PDF

Step 1. Reconstruct pupil function from captured LR images. Note that although the recovered pupil function only characterizes the deviations of DF illumination angles, joint reconstruction with BF images can improve the contrast relative to the background. Besides, the offset pupil is not affected by the number of iterations like the object function. More iterations only increase contrast for easier segmentation. When the number of DF images is greater than that of BF images, within 5 iterations is sufficient to maintaining effective information about the boundary of pupil function, and also saving time. For cases where the number of DF images is less than that of BF images, more iterations are required. Based on the distinct performances inside and outside the passband, we suggest using different thresholds for segmentation separately.

Step 2. Calibrate the DF sub-aperture positions. First, determine the initial search directions based on the Pupil search direction criterion. Second, a boundary-finding scheme is performed to identify shifted pupil centers $\Delta \boldsymbol {k}_{l}$ that occur near the sharp decline at objective radius $R_{obj}$

Step 3. Once the localization of DF sub-apertures is accomplished, a compact solution about global misalignment parameters $[\Delta x^c,\Delta y^c,\theta ^c,h^c]$ can be enclosed by Eq. (2). Arbitrary two sub-aperture positions correspond to a set of approximate solutions. For the case of "Rotation only", the lateral offsets can be considered as 0 and each sub-aperture position is enough to obtain an independent solution about $[\theta ^c,h^c]$. Take the average values of all obtained solutions and the final result will be used to update the HR images in the exact sub-apertures.

3. Simulations

To verify the effectiveness of TPF-FPM and cover the three mentioned misalignment cases, we perform several groups of simulations with different positional parameters. The platform is a laptop with a CPU (Intel i7-6700) and no parallel computing framework is utilized. The systematic parameters stay the same as in Sec. 2.1. Two images Baboon and Aerial are employed as HR intensity and phase samples, respectively, each containing $512\times 512$ pixels. To simulate the actual imaging circumstance, the Gaussian noise with the mean of 0 and standard deviation of $5\times 10^{-5}$ is artificially introduced to each LR image. For comparison, the algorithmic self-calibration approach SC-FPM is chosen as a cross-reference. The recovered positional parameters of the two methods are shown in Table 1.

Table 1. Recovered positional misalignment parameters of SC-FPM and TPF-FPM

View Table

The calibrated parameters of the two methods are highly consistent with the actual positional misalignment parameters. The height factor $h^c$ by SC-FPM is slightly higher than that of TPF-FPM. As an ASC method using simulated annealing algorithm, SC-FPM typically requires about 30 iterations to converge to stable solutions, while TPF-FPM has a shorter processing time. As described in Sec. 2.4, the calibration process of TPF-FPM mainly consists of two parts: pupil function recovery and boundary-finding implementation. The former generally takes about 15 s, while the time of the latter depends on the complexity of the judgment steps required to determine the initial search directions. For the same misalignment case, the time required for correction is roughly equal.

To further evaluate the improvement of reconstruction quality, we employ the calibrated parameters to correct the frequency positions of sub-apertures during reconstruction. The recovered HR images of two algorithms are shown in Fig. 10, where image quality assessment indexes PSNR and SSIM are selected to quantitatively evaluate the reconstructed results. Figures 10(a) and (b) are HR samples. The labels $①\sim ④$ correspond to different misalignment cases in Table 1. The reconstructed results by EPRY, SC-FPM, and TPF-FPM are shown in red, green, and blue boxes, respectively. The quality of recovered HR intensity images using both calibration algorithms is satisfactory, with TPF-FPM slightly outperforming in evaluation values. As for the phase images, it can be intuitively seen from the results of $② \sim ④$ that TPF-FPM performs better than SC-FPM. Some dark artifacts appearing in EPRY are effectively eliminated after the calibration of TPF-FPM. However, SC-FPM exhibits a certain degree of mutual interference between intensity and phase, as shown in Figs. 10(f2)-(f3), demonstrating its performance is relatively poor when "Translation dominates". It’s related to the difference in correction logic between the two algorithms. SC-FPM inevitably falls into the local optimum as the evaluation metric based on the entire image is sensitive to noise. In contrast, TPF-FPM utilizes the local contour features of the recovered pupil function, and a few iterations further avoid excessive accumulation of noise. Despite these two methods involving more or fewer iterations, TPF-FPM achieves more accurate positional parameters and better reconstruction compared to SC-FPM.

 

Fig. 10. Recovered HR images by EPRY, SC-FPM, and TPF-FPM with 10 iterations. (a)-(b) The HR intensity and phase samples; (c1)-(h4) The recovered HR intensity and phase images by three algorithms, labeled with red, green, and blue boxes, respectively. The evaluating values of PSNR and SSIM are marked in the lower right corner.

Download Full Size | PDF

4. Experiment

Based on the analysis of the intensity of shifted pupils in Sec. 2.2, we chose a testicular tissue slice with more high-frequency details as the calibration sample to suppress the noise in the experimental system and make the boundary of the transformed pupil function easier to segment. The model misalignment parameters of both calibration methods were obtained using this sample and employed to reconstruct other samples (USAF-1951 board and blood smear) to demonstrate the generalizability of corrected positions. A $13\times 13$ programmable LED array (8.128 mm spacing) is utilized to provide angle-varied illumination, which is placed arbitrarily in the horizontal plane to produce the model misalignment error. In our experimental system, the time spent in collecting 169 LR images is 51.11 s. Each LED element is embedded with tricolor illumination modules whose central wavelengths are 460 nm, 505 nm, and 629 nm, respectively. A 4 $\times$ /0.13 NA objective lens and a CMOS camera with pixel pitch 6.5 $\mu$m constitute the acquisition part of the FPM system.

The central LR images of the calibration sample captured under the illumination of 629 nm and 505 nm are shown in Figs. 11(a1)-(a2). Based on the correction theory of our method, we attempt to qualitatively analyze the calibrated results through pupil functions. The recovered pupil functions by three algorithms are shown in Figs. 11(b1)-(d2). It’s apparent from Figs. 11(b1)-(b2) that some bright profiles appear in the lower right section outside the passband. We can presumably determine that there is a negative translation in the LED array. In the first group of calibrations, SC-FPM converged on $[\Delta x=-0.896$ mm, $\Delta y=-0.845$ mm, $\theta =-0.582^{\circ }, h=116.10$ mm] at 31 iterations, taking 31.425 s. The results of TPF-FPM are $[\Delta x=-1.177$ mm, $\Delta y=-0.804$ mm, $\theta =-0.701^{\circ }, h=98.42$ mm], taking 23.527 s. In the second group of calibrations, SC-FPM converged on $[\Delta x=-1.082$ mm, $\Delta y=-0.909$ mm, $\theta =-1.234^{\circ }, h=108.63$ mm] at 60 iterations, taking 64.761 s. The results of TPF-FPM are $[\Delta x=-1.638$ mm, $\Delta y=-0.941$ mm, $\theta =-0.785^{\circ }, h=99.47$ mm], taking 23.643 s. The height factor by SC-FPM far deviates from the practical measurement of $h=98$ mm, yet the height factor by TPF-FPM is closer to the actual situation and provides better guidance for aligning the LED array. Reflected in the pupil functions of Figs. 11(c1)-(c2), there are still numerous bright spots outside the passband even after the correction of SC-FPM. While the absence of any noticeable bright spots outside the passband in Figs. 11(d1)-(d2) indicates that TPF-FPM is more effective in eliminating the positional deviations of DF area. In addition, the running time of TPF-FPM is shorter than that of SC-FPM, and in the second group of calibrations, it’s even less than half of SC-FPM.

 

Fig. 11. The recovered pupil functions by three algorithms. (a1)-(a2) The central LR images of the selected calibration sample under the illumination of 629 nm and 505 nm, respectively. (b1)-(d2) The recovered pupil functions by EPRY, SC-FPM, and TPF-FPM, respectively.

Download Full Size | PDF

One segment (128$\times$128 pixels) of USAF-1951 board under 629 nm illumination is chosen as the test sample. Figure 12(a) presents the captured LR image. The comparative effects of three algorithms are shown in Figs. 12(c1)-(c3) and (f1)-(f3), marked with the same color as in the Fig. 11. An area containing two sets of stripes (group 7, elements 5 and 6) is enlarged to illustrate the boost on intensity image. It’s evident from the intensity diagrams that TPF-FPM alleviates the dark-spot effect present within the fringes in EPRY and enhances the integrity of the contour. However, the details of stripes are still unclear and even worse than EPRY after the correction of SC-FPM, demonstrating that its calibrated positions are imprecise. Meanwhile, the cross-section mean intensity distributions of the purple boxes in Figs. 12(b) and (d1)-(d3) are presented to further support the effectiveness of TPF-FPM. As shown in Figs. 12(e1)-(e4), TPF-FPM improves the intensity contrast, while the result of SC-FPM is inferior to that of EPRY. The comparison can also be proved in the recovered phase maps. TPF-FPM produces better performance in removing remarkable artifacts in the background.

 

Fig. 12. Experimental results of one segment (128$\times$128 pixels) of USAF-1951 resolution target recovered by three algorithms under 629 nm illumination. (a) The captured LR intensity image corresponding to central LED illumination; (b) The zoom-in of the orange rectangle in (a); (c1)-(c3) The recovered intensity images by three algorithms, respectively; (d1)-(d3) The zoom-ins of orange rectangles in (c1)-(c3); (e1)-(e4) The cross-section mean intensity distributions of the purple boxes in (b) and (d1)-(d3); (f1)-(f3) The recovered phase images by three algorithms, respectively.

Download Full Size | PDF

Besides, we further show the results of recovering a blood smear with the corrected positions under 505 nm illumination. The HR images recovered by three algorithms are illustrated in Figs. 13(c1)-(c3) and (e1)-(e3). The zooms-ins of orange boxes in recovered intensity and phase images are shown in Figs. 13(d1)-(d3) and (f1)-(f3) to compare the reconstructed details. The results of EPRY show some wrinkle artifacts before misalignment calibration, while both correction algorithms eliminate this issue. Although SC-FPM successfully eliminates most artifacts, it doesn’t perfectly reconstruct the profiles of blood cells, especially those overlapping ones. In contrast, TPF-FPM yields more satisfying results with improved blood cell contour details. Additionally, when it comes to phase maps, SC-FPM may produce dark shadows in the background, whereas TPF-FPM delivers a cleaner and more uniform overall phase.

 

Fig. 13. Experimental results of one segment (40$\times$40 pixels) of blood smear recovered by three algorithms under 505 nm illumination. (a) The segment of central LR intensity image; (b) The zoom-in of the orange rectangle in (a); (c1)-(c3) The recovered HR intensity images; (d1)-(d3) The zoom-ins of orange rectangles in (c1)-(c3); (e1)-(e3) The recovered HR phase images; (f1)-(f3) The zoom-ins of orange rectangles in (e1)-(e3).

Download Full Size | PDF

5. Conclusion

In this paper, we propose a novel model misalignment calibration strategy called TPF-FPM, which is based on the conversion between the spatial frequency deviation and shifting pupil function. The feasibility and effectiveness of TPF-FPM are demonstrated by the great performances achieved in simulations and experiments. Compared with the joint optimization concept in algorithmic self-calibration approaches, we constructed an explicit physical mapping model to separate the misalignment calibration from the reconstruction. In addition, the utilization of the pupil function greatly weakens the interference of other error sources. Unlike calibration methods based on the central region of the BF area, we achieve overall calibration by minimizing the positional deviations of edge corners. The analysis of positional deviation distribution not only enables the determination of the comparative relationship between different positional parameters in an LED array, but also provides the potential for revealing the positional error threshold that seriously degrades reconstruction.

It’s worth discussing there are some points that could be improved for TPF-FPM in future work. Firstly, the pupil-function-based strategy can be extended to the dark-field localization under full-pose-parameter misalignment [28,29] or other illuminators. We only need to predetermine the distribution maps of shifted pupil centers and which LED elements correspond to the offset pupils that form the contour of pupil function in particular directions. Secondly, in addition to the boundary-finding scheme, the estimation of illumination position can also be posed as an object detection problem in neural networks to improve the estimation accuracy [30]. Thirdly, the recovery of pupil function is limited by the ePIE-based algorithm. Due to the direction and speed of error accumulation in the iterative route, the reconstruction may collapse in some special misalignment cases. When the positional deviation does not exceed the fault threshold, this issue can be addressed by appropriately reducing the update step size or adopting more robust iterative algorithms [31]. Fourthly, to improve the robustness to noise, TPF-FPM imposes certain requirements on the selected calibration samples. The intensity of the shifted pupil depends not only on the energy of current sub-spectrum, but also on the update weight function $w$ in Eq. (11). We expect the profile of offset pupils corresponding to the high frequency sub-spectrum could be effectively pronounced during the iterative process. This issue can be related to the construction of weight functions [27].