Burg algorithm for enhancing measurement performance in wavelength scanning interferometry

Wavelength scanning interferometry (WSI) is a technique for measuring surface topography that is capable of resolving step discontinuities and does not require any mechanical movement of the apparatus or measurand, allowing measurement times to be reduced substantially in comparison to related techniques. The axial (height) resolution and measurement range in WSI depends in part on the algorithm used to evaluate the spectral interferograms. Previously reported Fourier transform based methods have a number of limitations which is in part due to the short data lengths obtained. This paper compares the performance auto-regressive model based techniques for frequency estimation in WSI. Specifically, the Burg method is compared with established Fourier transform based approaches using both simulation and experimental data taken from a WSI measurement of a step-height sample.


Introduction
Wavelength scanning interferometry (WSI) can be used as a non-contact profilometry technique that can measure surface topography without the need for the mechanical movement of either the instrument or the measurand. This brings advantages over comparable techniques such as coherence scanning interferometry by allowing significant reductions in measurement time. Kikuta et al first demonstrated distance measurement using injection current driven wavelength scanning in diode lasers [1]. A Fourier transform based method for determining distance using wavelength scanning was then introduced by Suematsu and Takeda [2]. Kuwamura and Yamaguchi demonstrated WSI for surface profilometry [3], using a tuneable dye laser to extend the wavelength scanning range thus achieving axial resolution of a few tens of nanometres. Ruiz et al reported a version of WSI implemented using a tuneable external-cavity diode laser [4], although the limited tuning range limited resolution to several microns albeit over a relatively large range. The authors have previously reported on a variation of WSI using an acousto-optic tuneable filter (AOTF) combined with a halogen source to implement the required wavelength scanning [5]. This method has the advantage of providing a large wavelength scanning range, and thus good axial resolution, at relatively low cost. This paper reports on an investigation into applying autoregressive (AR) model based frequency estimation techniques to WSI and evaluates the resulting performance enhancements.
WSI is implemented using a swept wavelength light source to generate a set of spectral interferograms, detected by a CCD area image sensor, where each pixel corresponds to a point of the surface being measured. Each pixel records a sinusoidal intensity variation as the wavelength is scanned during the measurement. By calculating the phase shift that occurs during the wavelength scan, the height of the surface at each point may be determined.
The achievable axial (height) resolution and range for surface topographs obtained using WSI is dependent on the algorithm used to evaluate the phase shift obtained from the wavelength scanning process. For instance, the more accurate the determination of the phase shift, the higher the resolution that can be obtained. At present, Fourier transform based methods are used as the basis for phase shift estimation in WSI.
A simple way to determine the phase shift is by using a fast Fourier transform (FFT) to find the power spectral density (PSD) of the sinusoidal interference pattern [1]. The fundamental frequency as determined by the peak of the PSD function corresponds to the number of the cycles across the wavelength scanning range. However the resolution using this method is limited to 2, s l where s l is the synthetic wavelength of the light source. For operation using wavelength scan ranges in the visible region, this limits the axial resolution at the micron scale [2]. The resolution of this method can be improved to some degree by applying curve fitting techniques to better estimate the peak of the PSD. Takeda et al [3] developed another technique based on determining the rate of change of instantaneous phase during wavelength scanning. By sweeping the wavenumber of the light source linearly, the rate of change of instantaneous phase will also be linear. The optical path difference at each pixel can be calculated by finding the rate of change of phase as the wavenumber is swept. This yields better axial height resolution than is possible using the FFT based PSD method. However, since the method relies on phase unwrapping, the resolution decreases as the number of sampled points per cycle of the interference pattern decreases. Thus the axial resolution is not constant across the axial measurement range.
WSI is effectively a frequency estimation problem, because the rate of change of phase in the spectral interferogram determines the axial position of a given point the surface under test. The simple FFT and Takeda's method both have limitations in terms of axial resolution and resolution consistency respectively [6]. An alternative approach to the problem of frequency estimation is the use parametric methods which are based on the modelling a signal as a stationary process, the spectral density of which can be described by a specific set of parameters [7]. For the signal produced during the wavelength scanning process in WSI, if the parameters of the model can be estimated, the peak of the PSD function can be determined. A key benefit of parametric estimation of PSD is that it is the spectral resolution does not degrade as a function of data length unlike those methods based on Fourier analysis. This is advantageous for WSI as the data length, and thus axial resolution, is effectively limited by the sweep range and/or step resolution available from the light source. Additionally, if a given axial resolution wavelength can be maintained through a shorter swept wavelength range, using a parametric method, this can have benefit of reducing measurement time.
This paper discusses the application of the Burg method which estimates the parameters, and thus the PSD, when the signal is represented using an AR model. The signal in this case is the modulated interferogram intensity acquired during wavelength scanning. Both simulated data and experimental data from the measurement of a step-height sample are analysed using the Burg method and a comparison is made to the Fourier based approaches.

WSI principle of operation
The authors' implementation of WSI achieves rapid and absolute measurement of surface topography, with no phase ambiguity limitations [5]. Combining WSI and implementing parallelised signal processing with multi-core GPU can produce areal topographs in under 1.5 s, providing significant advantage for applications in embedded measurement and those requiring large numbers of repeated measurements. Figure 1 shows the structure of the WSI apparatus.
The interferometer is a Linnik configuration comprising a reference and measurement arm, formed by a beamsplitter and two identical objective lenses, OL1 and OL2 respectively. Light is coupled into the interferometer from the output of an AOTF via an optical fibre and collimating lens. The AOTF acts as a narrowband filter on the light provided by a broadband halogen white light source. The resulting interferogram produced is imaged onto an areal CCD sensor by the tube lens. The Linnik configuration is unbalanced in terms of optical path length; a virtual reference plane (VRP) designates the axial position in the measurement arm at which there is zero optical path difference. A PC controls the wavelength scanning via the acousto-optic tunable filter and acquires the image data from the CCD camera via a framegrabber.
The wavelength scanning process is achieved by tuning the driving frequency of the AOTF such that the optical wavenumber of the AOTF passband changes linearly. During the scanning process, 256 spectral interferograms are captured by the CCD camera.  Figure 2 shows an idealised sinusoidal intensity response captured by two pixels from on the camera, representing two points at differing axial positions (heights) on the measurand surface. In practice, this sinusoidal response is modulated by an envelope function resulting from the combined spectral characteristics of the various optical elements, in particular the halogen source and the AOTF. The generalised interference pattern can be described by the following equation: where a xy (λ) and b xy (λ) are the interference intensity bias and the amplitude of the fringe contrast respectively, λ is the scanned wavelength, and j xy (λ) is the phase. The x and y indices are the pixel number in the horizontal and vertical directions of the CCD detector respectively. The phase shift of the interference pattern due to the wavelength shift, Δj xy (λ), is given by where h xy is the axial distance between the VRP (point of zero optical path difference) and the measurand.

The Burg method
The Burg method is one way of estimating the AR coefficients and estimating the PSD. It works by minimising sums of squares of forward and backward linear prediction errors with respect to the coefficients [4]. From the definition of an AR model in (3), the residual t e can be calculated from the signal x t by

The tapered Burg method
The Burg method suffers a number of limitations, especially when applied to short data records. For instance, there are frequency errors that depend on the phase of sinusoidal components of the signal and spectral line splitting can occur, particularly when there is a high signal-to-noise ratio (SNR) [8]. These unwanted effects can be reduced by replacing the rectangular taper used in the original Burg method with a taper that is calculated to minimise the frequency error across all phases [9]. The modified algorithm works by 'tapering' the original Burg algorithm using a calculation for the frequency error the original algorithm produces for a sinusoid. Swingler [10] showed that the estimated frequency error for a sinusoid, x k cos , If the phase is considered to be a random uniformly distributed variable then the mean frequency error is zero and the variance is The taper is the chosen to minimise the average frequency error variance with respect to the taper.

Comparing Burg and tapered Burg methods
A numerical simulation was carried out in order to compare the efficacy of the Burg and tapered Burg methods for determining a fundamental sinusoidal frequency from a spectral interferogram generated by WSI. A discretised sinusoid having 256 samples in length and having a normalised frequency of 0.1083 was used as the input data set. The former parameter being the number of interferograms captured by our WSI apparatus, the latter corresponding to an axial position from the virtual reference of 30 μm. Assessing performance in terms of normalised frequency allows the assessment of algorithm in general terms which may then be applied to any WSI apparatus. Additive Gaussian white noise (AGWN) was added to the simulated sinusoid. The SNR, as defined by the ratio of the signal mean and the standard deviation of the applied noise, was 58.2 dB which is representative of typical operating conditions. Figure 3(a) shows the PSDs calculated using the Burg algorithm for three independent runs. The peak of the PSD estimates for each run appear to be effected by the presence of the applied noise which is clearly problematic if the method is going to be applied to WSI. Figure 3(b) shows the same datasets processed using the tapered Burg algorithm, and across all three runs this provides a stable result.
With the tapered Burg method established as the most stable in terms of being able to resolve spectral interferogram frequency, the study was then extended to look at the performance of the algorithm across the range of axial positions with respect to the VRP. Figure 4 shows the error in frequency estimation, now converted to the corresponding axial position error, as the axial position with respect to the VRP is moved in 100 increments between 20 and 50 μm (equivalent to a normalised frequency range of 0.072-0.180 in out apparatus). The study was run three times, with a doubling in the amplitude of AGWN applied each time, to realise an SNR of 31.6 dB, 45.3 dB and 59.7 dB for each run respectively. It can be seen that the overall trend of the error in frequency estimation increase as the axial position moves away from the VRP (increasing normalised frequency). Reduced the SNR appears to have little effect on the ability of the tapered Burg method to resolve the interferogram frequency.

Comparing the tapered Burg method with fourier transform based techniques
In order to evaluate any performance benefits associated with using auto-regressive model based frequency estimation techniques for WSI it is necessary to compare the performance with established Fourier based methods. A numerical simulation was performed where the axial position with respect to the VPR was moved over 100 increments between 20 and 50 μm, the results are shown in figure 5. The simple FFT shows a large cyclic error resulting from the discretisation of the frequency estimation; this may be alleviated to some degree applying interpolation techniques. The Takeda method, based on phase gradient estimation, shows the best performance in terms of error in frequency estimation, but it begins diverging and becomes effectively unusable beyond 41 μm. The tapered Burg shows a slightly higher estimation error than the Takeda but crucially it remains operational across the whole of the investigated range. This result points to the potential of the tapered Burg method of frequency estimation to increase the measurement range of WSI systems within the constraints of the physical apparatus i.e. depth-of-field and coherence length.
In order to validate the potential increase in measurement range using the tapered Burg method, the WSI apparatus shown in figure 1 was used to measure a small step height at two axial positions: 30 and 50 μm with respect to the VRP. The artefact measured was an etched waffle plate sample having a nominal well depth of 100 nm. An areal measurement of a single well was carried out for this study, the stated axial position corresponds to the position of the upper face of the plate. The sample was placed on a precision translation stage which allowed the axial positioning to be  easily adjusted. The measurement protocol was 256 interferograms taken between a source wavelength range of 590.98 and 683.42 nm using a x5 magnification objective. Figures 6(a) and (b) shows the resulting areal measurement obtained by processing the acquired spectral interferograms using the Takeda method. Where the upper face is located 30 μm from the VRP, the measurement appears to have been well resolved. But when the upper face is placed at 50 μm with respect to the virtual reference the measurement result clearly invalid as the algorithm becomes unstable.
The same measurement datasets were then processed using the tapered Burg algorithm. Figure 7 shows that both measurements with the upper face of the artefact at axial positions of 30 and 50 μm are successfully resolved. This validates the use of AR model based methods for extending the measurement range.