Accurate and reliable peak extraction of axial response signals plays a critical role in confocal microscopy. For axial response signal processing, nonlinear fitting algorithms, such as parabolic, Gaussian or sinc2 fitting may cause significant systematic peak extraction errors. Also, existing error compensation methods require a priori knowledge of the full-width-at-half-maximum of the axial response signal, which can be difficult to obtain in practice. In this paper, we propose a generalised error compensation method for peak extraction from axial response signals. This full-width-at-half-maximum-independent method is based on a corrected parabolic fitting algorithm. With the corrected parabolic fitting algorithm, the systematic error of a parabolic fitting is characterised using a differential equation, following which, the error is estimated and compensated by solving this equation with a first-order approximation. We demonstrate, by Monte Carlo simulations and experiments with various axial response signals with symmetrical and asymmetrical forms, that the corrected parabolic fitting algorithm has significant improvements over existing algorithms in terms of peak extraction accuracy and precision.
© 2019 Optical Society of America under the terms of the OSA Open Access Publishing Agreement
Compared with conventional optical microscopy, confocal microscopy (CM) is capable of performing optical sectioning of three-dimensional objects with enhanced lateral resolution and high sectioning capability [1,2]. CM is widely applied in fundamental research and industrial applications, and many excellent work has been done in areas such as biomedical science [3,4], material science , and micro and macro-scale surface metrology [6–9]. In CM, a small illumination pinhole is imaged onto a sample surface and reimaged onto a pinhole detector plane. If the illumination pinhole is accurately focused on the sample surface, the reflected intensity will pass through the detector pinhole with little loss in intensity, otherwise the intensity loss will be significant and dependent on the amount of defocus . A sequence of intensity values called the axial response signal (ARS) is recorded at different axial scanning positions. Usually, a sample surface is scanned in the axial direction with a fixed scanning step and the ARSs at each sampled lateral position are recorded. By extracting the peak position of a sampled ARS, the surface height can be estimated.
The accuracy and precision of the peak extraction process have significant influences on the measurement uncertainty in CM. In general, the performance of peak extraction algorithms is limited by random measurement noise, such as detector noise and the positioning error of the vertical scanning device. For example, Ruprecht et al. studied the noise-induced peak extraction for a centroid algorithm via the error propagation law . Rahlves et al. analysed the complex interaction of noise and the normalised intensity threshold for a conventional centroid algorithm . A correction factor-based systematic error reduction method was proposed by Ruprecht et al.  for accurate height extraction in CM. Another systematic error correction method based on the Fourier transform was proposed for star trackers or Shack-Hartmann wavefront sensors , where similar centroid estimation methods were applied. However, these two methods [11,13] require an additional full-width-at-half-maximum (FWHM) pre-calibration process, which is difficult to perform because the FWHM of the ARSs of a CM system varies at different surface gradients , and with respect to different spectral illumination components [15,16]. Therefore, it is necessary to develop a FWHM-independent systematic error compensation method that can enhance the accuracy of surface topography measurement, especially for surfaces that feature complex geometries and high roughness.
In our previous work , we have found that the sample surface height has a significant impact on the peak extraction accuracy. This effect has been underestimated in much research into CM metrology. Assuming a normally distributed measurement uncertainty in a CM system, if the sample surface is located exactly at the discrete position of the axial scanning layer, no systematic peak extraction error would occur. If the sample surface is located at any position between two adjacent scanning layers (this is often the case in practice), a systematic peak extraction error on the scale of tens to hundreds of nanometres can be observed . As an important uncertainty contributor, the sample surface height dependent peak extraction error should be quantified and compensated.
In this paper, we propose a corrected parabolic fitting algorithm (CPFA) that is FWHM-independent and provides a generalised systematic peak extraction error compensation solution for confocal surface topography measurement.
2. Principle of the parabolic fitting algorithm
In CM, an ARS is recorded at different axial positions, i.e. the vertical sampling sequence (VSS). A generalised VSS with equally spaced intervals can be written as3,18,19]10,14,15,17]. In Fig. 1, two ARSs are illustrated at different sample surface heights. Since the presence of optical aberrations may lead to an asymmetrical ARS, an empirical intensity threshold is usually applied before signal processing . Only the intensity values that are above are retained for further signal processing. These retained intensity values comprise the effective ARS and its corresponding dependent variable is defined as the effective VSS, which is treated as a column vector.
Parabolic fitting of the ARS in the effective area, i.e. with the effective VSS as the dependent variable, is implemented using a least-square fitting method, i.e. to minimise the following evaluation function Q
The evaluation function Q can be treated as a function with three variables, and . Minimisation of Eq. (3) leads to a set of matrix equations with these three variables . According to the properties of the parabola, the peak position of a general parabolic fitting algorithm (PFA) can be estimated as
In Fig. 2, the height-dependent systematic error, i.e. is calculated for the PFA in Eq. (4). Considerable systematic errors are found since the parabola is a coarse approximation to the sinc2-like ARS model. Note that the height effect is periodic with a period of one scanning step [11,13,17].
3. Corrected parabolic fitting algorithm
3.1 General analysis of parabolic fitting via a differential equation
The aforementioned two methods [11,13] are FWHM-dependent, which indicates their limited effectiveness for accurate surface topography measurement with varied FWHMs. We propose a FWHM-independent method for general error compensation, with which the influence of surface height-dependent systematic and random errors can be significantly suppressed. Based on Eq. (4), the first derivative of the systematic error of the PFA with respect to sample surface height is given as
The generalised expression in Eq. (6) is also valid for centroid algorithms, since the centroid estimation  has a similar expression to that in Eq. (4). Evidently, Eq. (6) is a first-order linear non-homogenous differential equation. and are height-dependent coefficients, which are related to the first derivative of with respect to. The first derivative can be approximated using the central difference of the sampled ARS.
3.2 Approximate solution of the differential equation
The solution of Eq. (6) can be used as an approximate error compensation for the PFA in Eq. (4). The key to the solution is to approximate coefficients and appropriately. A numerical illustration of and is shown along with sample surface height in Fig. 3. In Figs. 2 and 3, the height variation within one scanning step is divided into three height sections. To some extent, and can be treated as constants in each individual surface height section. Such an assumption can be justified using a Taylor series approximation of the sampled ARS. Thus, the differential Eq. (6) can be solved as17]Fig. 2. Since the systematic error is a piecewise function as shown in Fig. 2, the expression in Eq. (9) is the solution of Eq. (6), when the sample surface is located at the central height section. For the other two height sections, the solution of Eq. (6) is similar to that in Eq. (8), with the exception that the parameter in substituted with other X-correlated parameters.
In Eq. (8), coefficient can be deduced from the effective ARS and VSS, but is unknown before the measurement. In order to estimate for error compensation, we have introduced , which accounts for the summed intensity differences of the effective ARS,
According to Eqs. (8)-(12), the height-dependent error can thus be estimated asEq. (13) is called the corrected parabolic fitting algorithm (CPFA). The correct estimated peak position of the CPFA is given asEq. (13), including , and , could be expressed using the sampled ARS and VSS at sample surface height without prior knowledge of the FWHM. The error compensation of the CPFA in Eq. (14) automatically updates with external changes, such as sample surface height, scanning step, intensity threshold and FWHM.
The explored solution in Eq. (8) is only one part of the solution of Eq. (6). As we have emphasised, there are other X-related parameters for the other two height sections. This is the reason why we do not use the estimation in Eq. (11) as the final peak extraction position. The expression for these two height sections can be transformed into expressions similar to that in Eq. (11) and the X-related parameters coincide well with the solutions for Eq. (6) at the corresponding height sections.
In Fig. 4, the systematic errors of different peak extraction algorithms, including the PFA, CPFA and Gaussian fitting algorithm (GFA) , are illustrated at different sample surface heights. It can be seen that the CPFA achieves much smaller height-dependent systematic errors than the PFA and GFA. Here the sinc2 fitting algorithm(SFA) is excluded since the numerical simulations are performed based on the sinc2-like ARS assumption. Therefore, the systematic error of the sinc2 fitting is zero in this simulation . The proposed error compensation method with CFPA can also be used for other peak extraction algorithms, such as centroid algorithms or Gaussian fitting.
4. Evaluation of CPFA
4.1 Sample surface height-dependent evaluation model
Based on Liu’s research into peak extraction performance evaluation for nonlinear fitting algorithms based on Monte Carlo (MC) simulations [21,22], we have proposed a sample surface height-dependent (SHD) evaluation model to characterise the SHD peak extraction uncertainty: SHD peak systematic extraction error and peak extraction standard deviation via MC simulations . The SHD model is expressed asEq. (15) for peak extraction, which can be seen in Fig. 1. The peak extraction error with random noise is
Since random noise and sample surface height are influencing factors to peak extraction performance, peak extraction uncertainties including the mean and standard deviation are analysed here, at different sample surface heights. In our MC simulations, random noise settings for Eq. (15), specifically the positioning error of the VSS is assumed to be uniformly distributed while the detector noise is assumed to be normally distributed , similar to elsewhere . We perform 20,000 trails at each given sample surface height to ensure high confidence . The simulations are performed in the following manner: first simulating the sampled ARS at a certain surface height , adding random noise, normalising the discrete ARS with maximum intensity of 1, applying a normalised intensity threshold for the recorded ARS to get the effective VSS and ARS; then calculating the peak extraction errors at height based on Eqs. (15) and (16) for different peak extraction algorithms; and finally, extracting the means and standard deviations of the peak extraction errors at height for different algorithms. Our proposed CPFA is compared with conventional fitting algorithms, such as PFA and GFA.
4.2 Performance of CPFA under SHD model with noise
In Figs. 5(a) and 5(b), the mean and standard deviation of the peak extraction errors with simulated random noise are depicted at different sample surface heights, in which the positioning error of the VSS is uniformly distributed within [-5%,5%] of one scanning step, while the detector noise has a standard deviation (zero-mean) of 1% of the maximum intensity. The CPFA shows a simultaneous improved peak extraction performance for systematic error compensation and standard deviation suppression, even when compared with GFA. Note that, the performance of SFA can be treated as a reference since the ARS shape is assumed to be sinc2-like. When the signal shape deviates from the sinc2 function , the SFA does not necessarily outperform the GFA .
Since the measured surface height can be any arbitrary value, the root-mean-square (RMS) values [24,25] of the means and standard deviations are used to denote statistical means of the systematic errors and standard deviations , thus,
In Fig. 5(a) and 5(b), the RMS values of means of the peak extraction errors are 0.13, 0.032, 0.017 and 0.019 [%FWHM] for PFA, GFA, SFA and CPFA, respectively. The RMS values of the standard deviations are 0.335, 0.339,0.310 and 0.311 [%FWHM] for PFA, GFA, SFA and CPFA, respectively. Compared with GFA, the CPFA shows an accuracy enhancement by 54% and standard deviation reduction of 7.3%.
In order to demonstrate the general advantages of CPFA over GFA, three sets of simulated comparisons have been performed with different noise settings. The noise settings for these three sets are listed in Table 1. In each simulation set, the scanning step is changed while the other conditions, such as the intensity threshold and noise settings, are kept constant. The noise of the VSS is within several percent of one scanning step, while the detector noise has a standard deviation (zero-mean) of several percent of the maximum intensity.
Figures 6-8 show the variations of RMS values of means and standard deviations along with scanning steps under different noise levels. When the noise level is relatively high as in Figs. 6(a) and 6(b), the PFA can have better performance in terms of its standard deviations than the GFA, while worse performance in systematic errors than the GFA, which has been reported elsewhere . When the noise level is relatively low, as in Figs. 7(a) and 7(b), the GFA has better performance in terms of both the systematic errors and the standard deviations than the conventional PFA. It is clear that the CPFA has outperform the conventional PFA in terms of both the systematic error and standard deviation, regardless of the noise levels. Even compared with GFA, the CPFA can have a significant performance enhancement. More specifically, the CPFA has a much smaller RMS value for the systematic errors than the GFA, which indicates that the CPFA can have a significant accuracy enhancement over the GFA. As for the RMS value of the standard deviations, the CPFA can improve the standard deviation performance by approximately 10%. If we take the SFA into comparison, the CPFA have basically the same peak extraction performance as the SFA, without the SFA’s rigorous assumption on the ARS model. In other words, other nonlinear fitting algorithms such as the GFA can outperform the SFA, when the ARS model deviates from the sinc2 function  while our CPFA will not be influenced in such cases.
4.3 Performance of CPFA under asymmetrical ARSs
In practice, the CM suffers from inevitable aberrations, which can result in asymmetrical ARS [21,26]. Although the CPFA was developed based on a symmetrical ARS in section 3, its consistent performance over asymmetrical ARSs should be demonstrated. The recorded ARS with optical aberrations can be expressed as Eq. (20) is a variant of the Zernike polynomials based wavefront expression, in terms of the field-independent wavefront aberration coefficients .
As expressed in Eq. (18), the effect of coma is eliminated between the illumination and detection process while the combined effect of spherical aberration and astigmatism causes an asymmetrical ARS shape as in Fig. 9. The choice of the aberration coefficients can be referred to the previous research . It can be seen that the combined aberrations of the spherical aberration and astigmatism broaden the FWHM of the ARS. Figures 10(a) and 10(b) illustrate the peak extraction performance of aforementioned algorithms including the PFA, GFA and CPFA under asymmetrical ARSs as in Fig. 9. Comparing the Figs. 8 and 10, the combined aberrations-induced asymmetry leads to a significant increase in the systematic error level and standard deviation level. Nevertheless, the CPFA can still have a smaller systematic error and standard deviation than the PFA, GFA and SFA, indicating the CPFA’s robustness against asymmetrical ARSs.
One point that need be emphasised is that the CPFA is a general algorithm with a better adaptability to ARSs of different forms. In practical applications, the ARS model deviates from the conventional sinc2 function since an extended pinhole is used to guarantee appropriate light efficiency [15,19]. The Gaussian function is a good approximation to the sinc2 function, but not necessarily an appropriate approximation for the ARS model with an extended pinhole . In Figs. 6-8, a parabola is a less accurate approximation to a sinc2 function than a Gaussian function, but the CPFA can still have better peak extraction performance than the GFA. Accordingly, the CPFA may have more significant advantages for ARS in other forms. Usually, a signal filtering technique is used before the peak extraction operation. The filtered ARS may seriously deviate from a Gaussian function while our CPFA will not be affected.
5. Experimental comparisons
In order to demonstrate the peak extraction advantage of CPFA over the current peak extraction algorithms, including PFA, GFA and sinc2 fitting algorithms (SFA)  on peak extraction accuracy and precision, experimental operations as described elsewhere  are performed. In the experimental comparisons, the SFA is also analysed because real ARSs with extended detector pinhole usually deviate from sinc2 functions .
A simple chromatic confocal microscope (CCM) was set up as illustrated in Fig. 11. In CCM , a spectral ARS is acquired by using a spectrometer, instead of axial scanning in the laser CM. Different sample surface heights correspond to different spectral ARSs with different peak wavelengths. From a mathematical perspective, the spectral ARS in CCM and the ARS in laser CM are the same. Thus, using the CCM to evaluate peak extraction performance is also appropriate . However, there are two dominant factors that influence the measurement results of CCM including the peak extraction performance and the input-output models. More specifically, a specific peak wavelength corresponds to a specific surface height. This λ-height relationship is not linear, thus the choice of modelling the λ-height curve is also important in CCM. In order to eliminate the potential influence of the modelling on the measurement results, a piecewise linear fitting method is used to code and decode the relationship between peak wavelengths of spectral ARSs and sample surface heights.
In our CCM, a commercial chromatic confocal sensor head (CL20 MG210, STIL, France) is utilised. The broadband source (MWWHF2, Thorlabs, USA) covers an optical spectrum of 380 nm to 780 nm. The spectrometer (Ocean Maya Pro2000, USA) has a resolution of 0.46 nm/ pixel. Multi-mode fibres with 50 μm diameter and matched coupler are used as the light transmitters among the chromatic confocal head, light source and the spectrometer. A flat mirror is measured at different axial positions along the optical axis of the CCM and a calibrated PZT (P721.CDQ, Physik Instrumente, Germany) with closed-loop repeatedability of 10 nm, is used as the displacement actuator. As shown in Fig. 11, the same point on the mirror is measured at different axial positions, which are determined by the indicator value from the PZT. The CCM signal is recorded with a 0.2 μm interval. Around the test position, one spectrometer pixel of 0.46 nm in the wavelength domain corresponds to approximately 0.8 μm in displacement. Therefore, a 0.2 μm displacement increment of the PZT is within one scanning step of the CCM. At each axial position, twenty-five frames of spectrums with 10 ms integration time are obtained. The flat mirror is tested at approximately one hundred axial positions. An example measured axial response signal is shown in Fig. 12.
As illustrated in Fig. 11, A, B and C are three axial positions where the same point on the flat mirror is measured. The displacement values from the closed-loop PZT are taken as the true values of the target surface heights. More specifically, if the sample surface height at A is 0 μm, the heights at B and C are 0.2 μm and 0.4 μm, respectively. If there are no surface height-dependent systematic peak extraction errors, then these three points A, B and C (the transverse coordinate is the expectation of repeated peak wavelength at the same surface height while the vertical coordinate is the calibrated height from PZT) fall in a straight line. In our experiments, the λ-height curve is determined by the linear fit through A and C, as in Fig. 10, while the measured height value at axial position B can be decoded from this curve. The local height deviation at axial position B is calculated as the relative distance between the decoded height value and the true height value from the PZT. There is an implicit assumption that the coding curve between A and C is linear, which is valid since the displacement variation between A and C is close to 0.4 μm. A piecewise linear fitting method as described above is used for the input-output modelling of the CCM, thus the nonlinear impact of the dispersive lens in the CCM is eliminated. In our procedures, the local height deviations are treated as the quantitative indicators of systematic peak extraction errors, i.e. peak extraction accuracy while standard deviations of extracted heights under the same axial position are obtained to represent the peak extraction precision. According to the model evaluating indictors in Eq. (17), the RMS values of the local height deviations and height extraction standard deviations are obtained to represent the peak extraction algorithm’s performance.
In Fig. 13(a), the RMS values of local height deviations are 0.0347 μm, 0.0144 μm, 0.0132 μm and 0.0071 μm for PFA, GFA, SFA and CPFA, respectively. In Fig. 13 (b), the RMS values of height extraction standard deviations are 0.0365 μm, 0.0349 μm, 0.0348 μm and 0.0307 μm for PFA, GFA, SFA and CPFA, respectively. From the RMS perspective of local height deviations and height extraction standard deviations, the CPFA behaves better than PFA with about 80% accuracy enhancement and 16% precision enhancement in peak extraction. Even compared to other fitting algorithms, the CPFA has some advantages in terms of extraction accuracy and precision with at least 40% accuracy and 10% precision enhancement. Note that the RMS operations in Figs. 13(a) and 13(b) are necessary for peak extraction performance comparisons since the systematic error and standard deviation are sample surface height-dependent, as is illustrated in Fig. 5. In Fig. 13 (b), the height extraction standard deviation is unstable for these fitting algorithms, which can be attributed to other nonlinear effects such as the modulation effect of spectrometer on the spectral ARS, light power fluctuation or multiplicative noise.
A corrected parabolic fitting algorithm (CPFA) for general error compensation, which is independent of the FWHM of the ARS, is proposed in this paper for confocal ARS processing. Our simulations show that the CPFA has advantages over conventional fitting algorithms in terms of peak extraction accuracy and precision, providing that the ARS is recorded at different sample surface height and suffers from sampling noise and optical aberrations. Our chromatic confocal microscope experiment has also demonstrated the effectiveness of our CPFA and shown that the CPFA can achieve approximately 40% accuracy enhancement and 10% precision enhancement over conventional fitting algorithms. The different ARS forms between the simulations and experiments have indicated a general potential of CPFA for further applications with enhanced adaptability of ARSs in different forms.
National Natural Science Foundation of China (NSFC) (51875227, 51675167, 51705178); National Instrument Development Specific Project of China (2011YQ160013); The Key Grant Project of Science and Technology Program of Hubei Province of PR China (2017AAA001); Shenzhen Basic Scientific Research Project (JCYJ2017030717134710).
The authors thank Dr Rong Su from the University of Nottingham for valuable suggestions.
The authors declare that there are no conflicts of interest related to this article.
1. T. Wilson, Confocal Microscopy (Academic, 1990).
2. R. Leach, Optical Measurement of Surface Topography (Springer, 2011).
3. M. Gu, Principles of three dimensional imaging in confocal microscopes (World Scientific, 1996).
4. H. Cui, W. Zhao, Y. Wang, Y. Fan, L. Qiu, and K. Zhu, “Improving spatial resolution of confocal Raman microscopy by super-resolution image restoration,” Opt. Express 24(10), 10767–10776 (2016). [CrossRef] [PubMed]
5. B. V. R. Tata and B. Raj, “Confocal laser scanning microscopy: Applications in material science and technology,” Bull. Mater. Sci. 21(4), 263–278 (1998). [CrossRef]
7. J. Liu, Y. Wang, K. Gu, X. You, M. Zhang, M. Li, W. Wang, and J. Tan, “Measuring profile of large hybrid aspherical diffractive infrared elements using confocal profilometer,” Meas. Sci. Technol. 27(12), 125011 (2016). [CrossRef]
12. M. Rahlves, B. Roth, and E. Reithmeier, “Confocal signal evaluation algorithms for surface metrology: uncertainty and numerical efficiency,” Appl. Opt. 56(21), 5920–5926 (2017). [CrossRef] [PubMed]
14. F. Mauch, W. Lyda, M. Gronle, and W. Osten, “Object depending artifacts in confocal measurements,” Proc. SPIE 8466, 846609 (2012). [CrossRef]
17. C. Chen, J. Wang, X. Liu, W. Lu, H. Zhu, and X. J. Jiang, “Influence of sample surface height for evaluation of peak extraction algorithms in confocal microscopy,” Appl. Opt. 57(22), 6516–6526 (2018). [CrossRef] [PubMed]
18. J. Liu, J. Tan, T. Liu, S. Wang, and Y. Cao, “Axial performance parameters developed for analytical design of center shaded filter in high aperture focusing system,” Opt. Commun. 283(21), 4190–4193 (2010). [CrossRef]
20. J. Tan, C. Liu, J. Liu, and H. Wang, “Sinc2 fitting for height extraction in confocal scanning,” Meas. Sci. Technol. 27(2), 025006 (2015). [CrossRef]
21. C. Liu, Y. Liu, T. Zheng, J. Tan, and J. Liu, “Monte Carlo based analysis of confocal peak extraction uncertainty,” Meas. Sci. Technol. 28(10), 105016 (2017). [CrossRef]
22. J. Liu and J. Tan, Confocal Microscopy (Morgan & Claypool, 2016).
23. BIPM, IEC, IFCC, ILAC, ISO, IUPAC, IUPAP and OIML. Evaluation of Measurement Data—Supplement 1 to the ‘Guide to the Expression of Uncertainty in Measurement’—Propagation of Distributions using a Monte Carlo Method. (Joint Committee for Guides in Metrology, Bureau International des Poids et Mesures, 2008).
24. J. Wang, L. Pagani, R. K. Leach, W. H. Zeng, B. Colosimo, and L. P. Zhou, “Study of weighted fusion methods for the measurement of surface geometry,” Precis. Eng. 47, 111–121 (2017). [CrossRef]
25. J. Wang, L. Pagani, L. P. Zhou, X. J. Liu, W. L. Lu, R. K. Leach, and X. Q. Jiang, “Uncertainty-guided intelligent sampling strategy for high-efficiency surface measurement via free-knot B-spline regression modelling,” Precis. Eng., In press, [CrossRef]
26. T. Wilson and A. R. Carlini, “The effect of aberrations on the axial response of confocal imaging systems,” J. Microsc. 154(3), 243–256 (1989). [CrossRef]
27. H. H. Hopkins, Wave Theory of Aberrations (Oxford, Clarendon Press, 1950).