1. Introduction

Fourier-domain schemes of optical coherence tomography (OCT) provide various attractive features over time-domain OCT (TD-OCT) [1–4]. Because of the higher imaging speed and better sensitivity, Fourier-domain OCT (FD-OCT) has become dominant in most biomedical applications. Two implementation methods have been developed to acquire OCT image signals in the Fourier domain. Spectral-domain OCT (SD-OCT) utilizes a spectrometer that consists of a dispersive element and a line camera to acquire a spectral fringe directly in the spectral domain [3]. On the other hand, swept-source OCT (SS-OCT) takes advantage of a frequency-encoded light source for frequency-domain signal acquisitions [4]. It has been widely acknowledged that the FD-OCT schemes commonly suffer from some types of image artifacts that are not observed in TD-OCT systems [3–6]. Effective elimination or suppression of those artifacts has drawn a great deal of attention either by using optical techniques or signal processing methods.

Among the characteristic artifacts of the FD-OCT, the presence of a fixed pattern results in strong erroneous horizontal lines laid over an image [3–6]. This fixed-pattern noise artifact usually comes from the non-flatness of the reference spectrum that is to be interfered with the signal light. The spurious spectral pattern of the reference can be caused by various mechanisms for a practical system implementation. The source spectrum itself can exhibit a fine spectral structure that generates a considerable Fourier-transform conjugate in the OCT image [7]. A static noise pattern may also be produced by optical interferences that occur in the middle of the optical pathway of the reference light. A photodetector may exhibit an intrinsically fringed spectral responsivity, particularly for a spectrometer based on a dispersive grating element. Besides, an OCT system may have static reflection surfaces at the end of the sample arm. This artifact is especially severe for an OCT endoscope catheter imaging since some reflection points are embodied by the fiber end, surfaces of the micro-lens, etc. Those reflections may contribute to generation of fixed horizontal lines which are not related with the actual sample under imaging. Even when those reflection points are out of the effective imaging range, they can still be visible in the final image due to the mutual interference of those reflected fields (autocorrelation artifact). Due to its static nature, this artifact produces straight lines running horizontally in an OCT image. Of course, the fixed-pattern artifact is suppressed by a balanced detector for the case of SS-OCT since such a noise pattern is usually given as a common-mode noise to the balanced detection ports [4]. However, the effect of the balanced detection is not always satisfactory but vulnerable to a spectral imbalance or limited performance of the electronics [8]. Furthermore, some fixed patterns are associated with fixed real reflections from the end of the sample arm as in the case of the OCT endoscope catheter. This type of fixed patterns cannot be rejected by the balanced detection that only suppresses common-mode signal components.

The fixed-pattern artifact can be removed if one has a good knowledge of the reference spectrum which is roughly defined by the detected signal in the absence of a sample. However, measuring a reference spectrum is not always convenient and effective. Some applications do not allow a simple calibration process done in seconds. For example, an OCT endoscope system cannot be calibrated with ease if the system needs to be calibrated in connection with the catheter concerning its internal reflections. Besides, the reference spectrum may vary not only because of the change of the source spectrum but also due to the possible variation of the optical phase or the polarization of light involved with the optical interference that contributes to the noise generation. The demands for frequent calibrations may make the system very inconvenient to use for practical applications. To avoid those issues, one can take advantage of a signal processing method that extracts the reference information to remove the fixed-pattern artifact. It is widely acknowledged that a reference spectrum can be extracted by averaging a plurality of A-line spectra that are acquired with a sample in place [3]. Each interferometric spectrum is subtracted by the mean spectrum to remove the static component without a separate step of reference spectrum measurement. This mean-spectrum subtraction method is popular for both SD-OCT and SS-OCT due to the practical implementation simplicity and inherent immunity to the variation of the reference. One can utilize a slightly modified method of numerically dividing a detected spectrum by a reference spectrum of the mean spectrum instead of the more common subtraction processing [3,6]. The mean-spectrum subtraction, in this report, collectively includes those methods that are based on a reference spectrum obtained by averaging A-line data.

In this study, we reveal a drawback of the conventional mean-spectrum subtraction method involved with its imperfect nature of the fixed-pattern elimination. We propose alternative signal processing methods to overcome this problem and further enhance the image quality. It has been found that a statistical error in averaging is caused by a small number of high-amplitude points that are present mainly due to the tissue-air interface of a high index contrast. This results in residual horizontal lines that significantly degrade the image quality in the final OCT image. Through an experimental demonstration, we found that this problem can be alleviated by using advanced statistical treatments of the raw data. Complete elimination of a static pattern was accomplished by our new processing methods. Therefore, an exact knowledge of the reference spectrum can be extracted from the OCT image data without using other means. It is expected that our methods can be used not only for suppressing the fixed-pattern noise but also for additional signal processing schemes that need the accurate information of the reference spectrum under a possible variation of the spectrum.

2. Fixed-pattern noise and mean-spectrum subtraction

The fixed-pattern artifact can be described by a non-DC signal component that statically resides in a signal A-line. In an FD-OCT system, the signal spectrum is given by the square-law detection of the interfered light components taken as a raw data set of an image A-line. When a complex reference field is denoted by E_r, and the desired reflection signal from a sample is expressed by a sum of complex sample fields E_n’s, the measured spectrum, G( k ), is given as a function of frequency, k, by

The spectral interferogram of the third term or the fourth term in the right-hand side of Eq. (1) appears in the form of a sinusoidal fringe. Its phase is very sensitive on the relative position of the reflection point. For a biological sample, the probability distribution of the phase is presumed to be random for most of the reflection points in a tissue. Thus, summation of the measured spectra effectively washes the fringe out so that the averaged interferogram converges to zero by the number of spectra increasing in the averaging operation. This is the principle of the mean-spectrum subtraction method. By averaging a large number (N) of A-line spectra, the mean spectrum [denoted by 〈G〉, here] effectively becomes the reference spectrum of the first two terms in Eq. (1) as

However, the assumption of random phase distribution is not always effective and may not allow a complete elimination of the interference fringe. A relatively small number of high-amplitude reflection points are often observed in an OCT image, and they may produce a statistical error that is not effectively removed from the averaged A-line spectrum. In typical OCT imaging, the interface of a sample tissue and the surrounding medium (air or water in many cases) forms a high amplitude signal due to the large index contrast. As it is usually smooth in its morphology, a signal from this interface well correlates to signals of the neighboring A-lines. This characteristic produces a high statistical bias so that the mean-spectrum subtraction cannot eliminate the static pattern completely.

The failure of the mean-spectrum subtraction was experimentally confirmed in this study. We used an SD-OCT system operating at a long-wavelength band of 1.3 μm. Our SD-OCT system was constructed by using a dispersive grating and a high-speed InGaAs line camera. The resolution of the system was 11 μm for the axial and 16 μm for the transverse, respectively. The acquisition speed was 41,000 A-lines per second. A set of OCT data was processed by different methods for comparison. Figure 1 shows the raw OCT image (a), the image processed by the physically measured reference (b), and the OCT image processed by the mean-spectrum subtraction (c) for a human finger tip. The sample was slightly tilted on purpose in order to better visualize the fixed-pattern artifact. Our system exhibited a strong fixed pattern in the raw OCT image. A lot of horizontal lines are clearly visible in the raw OCT image of Fig. 1(a). This noisy pattern was removed by subtracting the reference spectrum from each A-line spectrum as seen in Fig. 1(b). This reference spectrum was obtained by an average of 1,000 spectra acquired with no sample in place. The mean-spectrum subtraction was also demonstrated without using the physically measured reference. The mean spectrum was obtained from the whole spectrum data of the image for N = 1,024. As shown in Fig. 1(c), most of the fixed pattern was removed by this method. But some horizontal lines were still visible as a noisy pattern. The axial positions of those lines coincided with the positions of the horizontal tissue-air interfaces of strong reflectance.

 

Fig. 1 Raw OCT image (a), image processed by the physically measured reference (b), and the OCT image processed by the conventional mean-spectrum subtraction method (c) for the same OCT data of a human finger tip.

Download Full Size | PDF

3. Statistical analysis of the fixed-pattern noise

It is advantageous to perform a noise analysis of the fixed pattern after Fourier transform since each noise source is located at an axial point while its effect is spread over the entire spectral range in the OCT spectrum. Due to the linear property of Fourier transform (denoted by F{⋅⋅⋅} in this report), the Fourier transform of an averaged spectrum is equal to the average of the transformed A-line spectra:

In this study, it was found that a statistical characteristic of the OCT image can result in residual noise lines as a consequence of averaging complex reflectance signals of a non-monotonic probability distribution. In addition, the noise source can be differentiated in the depth-resolved domain of OCT with ease. As a case study, the 162nd horizontal line which is highlighted by a blue arrow in Fig. 1(c) was investigated in detail. Since the raw A-line data were given as a set of complex numbers, a statistical analysis could be performed with a complex plane that represents a distribution of complex data points (amplitude and phase) in a Cartesian-coordinate space. Figure 2 shows the complex-plane representation of the A-line data (a), their histogram on the real axis (b), histogram on the imaginary axis (c), and the magnified version ( × 20) of the complex plane for the central area (d). The blue dots in the complex-plane representation stand for complex reflectance data (x + iy) obtained from the horizontal line of Fig. 1(a) at the axial position of observation.

 

Fig. 2 Complex-plane representation of the A-line data (a), their histogram on the real axis (b), their histogram on the imaginary axis (c), and the magnified version ( × 20) of the complex plane for the central area (d). The data are from the 162nd line of the raw OCT image in Fig. 1 (a).

Download Full Size | PDF

The overall distribution of the data points resembled a two-dimensional Gaussian distribution centered at a certain complex number as depicted by Figs. 2(b) and (c). This distribution was shifted by a systematic bias that a reference signal gives to the measured interference spectra. This shift resulted in the fixed horizontal pattern visible in the final image. The mean-spectrum subtraction method, or the mean-line subtraction, is understood as a simple way to find the center of the data point distribution by taking the mean vector (i.e., the mean complex number). In this method, each dot is shifted back by this mean vector for the final distribution to be relocated around the origin in an even and symmetric distribution manner. However, the mean vector may not coincide with the effective center of the distribution because of the statistical oddness. As observed in Figs. 2(a) and (d), some distinguished points have distinctly high amplitudes while the others are densely located around the origin with low amplitudes. These high-amplitude points are distributed in an odd and asymmetric manner so that their distribution does not look purely random. These were generated by a different reflectance source from those of the central area and obeyed a distinguished distribution function as a consequence. The high-amplitude points came from the high-contrast interface of the tissue and the air shown by a bright horizontal curve in Fig. 1. Their amplitudes were far higher than the others in one or two orders of magnitude, i.e., 20 to 40 dB above the reflectivities of the other ordinary points. Thus, the average of the high-amplitude points has a chance to produce a statistical fluctuation that is stronger than the low amplitudes of the ordinary data points.

Such a residual noise can be produced by two mechanisms. First, the phase of such a high-amplitude point may not be random. As seen in Fig. 2(a), the high-amplitude points did not exhibit completely random phases but rather was concentrated to the lower right corner (x>0, y<0) with a certain average phase. It was obvious that the smooth tissue-air interface of the sample could make a certain correlation of those phases. In Fig. 2(d), the black dot represents the mean complex number which was averaged over the whole data for the 162nd line. It was offset to the lower right with respect to the effective center of the distribution as anticipated by the odd phase distribution of the high-amplitude points.

Second, even under an assumption of nearly random phases without an apparent correlation, a mean value calculated from finite sampled information deviates as a result of sampling fluctuation. A morphological roughness of the tissue-air interface determines the phase of the reflectance. As the deterministic motion of a roulette wheel looks random in its result, the phase of the sample surface can be understood as a sampled value of an apparently random variable. However, a statistical sampling does not exactly conserve the statistical moment in a stochastic process. For example, an expectation value of a numbered cubic die after a finite number of measurement samples is not exactly 3.5 (1/6 + 2/6 + 3/6 + ⋅⋅⋅ + 6/6) but fluctuates around 3.5 because of the sampling fluctuation. Thus, even if the surface morphology provides no apparent phase correlation, a mean of apparently random complex reflectance cannot be zero. The standard deviation of such a statistical fluctuation decreases in an inversely proportional rate to the square root of the number of the statistical samples which are involved with the mean. Hence, for the number of high-amplitude points, M, their mean produces a sampling error of a standard deviation to be (1/M)^1/2 of the initial deviation of the samples which roughly equals the mean amplitude. In a rough estimation for the given case of Fig. 2, the sampling error in the mean of the high-amplitude points must have been larger than 1/10 of their average amplitude for M<100. This error is still higher than the low amplitudes of a majority of data points located around the center. Therefore, even under the most favorable assumption of the random phase distribution, a mean can produce an error in finding the distribution center when a part of the data points have extraordinarily high amplitudes beyond the magnitude expected by a monotonic normal distribution.

4. Median-line subtraction

The disadvantage of a mean value discussed in the previous section is well recognized in various applications of the statistical analysis. A median estimator has been suggested as an alternative to alleviate this problem in those cases. It has a decreased sensitivity on high-amplitude data points and can give more tolerant accuracy in determining the center of the statistical distribution. We propose a new signal processing method of median-line subtraction to avoid the observed drawback of the conventional mean-spectrum subtraction method. A median is defined by the halfway value in a set of data when they are sorted in an ascending or descending order. In this research, a median for complex numbers is defined by a sum of the median of the real terms and that of the imaginary terms for simplicity. For a complex-numbered A-line data g _l, the median line, g′ is obtained by

A reference A-line was extracted by the complex median operation given by Eq. (5) from the OCT data. Each A-line was subtracted by the reference A-line of the median line to get rid of the fixed pattern. Figure 3 shows the OCT image obtained by the physical reference spectrum measurement (a), the image processed by the median-line subtraction (b), and the image processed by the mean-spectrum subtraction (c), displayed together for comparison. Our method of the median-line subtraction successfully removed the fixed-pattern noise without any residual lines remaining in the image. No significant difference was found between Figs. 3 (a) and (b), which are the images obtained by the physical calibration and the median-line subtraction, respectively. This feature was also confirmed by the complex-plane analysis depicted in Fig. 2. The median-line value for the 162nd axial position was marked by a red dot in Fig. 2 (d). The position of the complex median matched the effective center of the data point distribution. In contrast, the complex mean marked by a black dot in Fig. 2 (d) deviated from the center. This characteristic is not expected, however, for a median operation in the spectral domain (k-domain) before Fourier transform. A median of a Fourier transform is not equal to a Fourier transform of a median unlike the linear property of a mean operation suggested by Eq. (3). The effect of a high-amplitude reflection point spreads over the full range of the spectrum before the transform, not to be removed in the spectral domain with ease. Thus, our median-line subtraction should be processed with the A-line data after the domain conversion.

 

Fig. 3 OCT image obtained by the physical reference spectrum measurement (a), image processed by the median-line subtraction (b), and the image processed by the conventional mean-spectrum subtraction (c), displayed together for comparison.

Download Full Size | PDF

Due to the successful features found in the OCT images processed by the median-line subtraction method, the result of those signal processing methods could be utilized to measure the reference spectrum with a sample in place. An inverse Fourier transform of the reference A-line directly produces the reference spectrum described by the first two terms of Eq. (1). Figure 4 shows the reference spectra obtained by a physical calibration (blue dots), transform of a median line (black solid line), and the difference between the two spectra. The difference was less than 1% for the full band. It was believed to come from the sum of sample reflection power as described by the second term of Eq. (1). This term is not visible for the calibration spectrum since it was measured without a sample in place. A slight difference in the spectral shape of the difference spectrum to the reference spectra may have reflected the non-flat and imbalanced spectral transmission of the sample power path with respect to the reference power transmission.

 

Fig. 4 Reference spectra obtained by a physical calibration (blue dots), transform of a median line (black solid line), and the difference between the two spectra.

Download Full Size | PDF

The reference spectrum directly extracted from image data can be used for additional signal processing methods that require an exact knowledge of the reference spectrum. Since this extraction does not need an actual measurement of the reference spectrum, it provides a way of using real-time spectrum to monitor the light source of an OCT system. In addition, the obtained information of the reference spectrum can be fed to signal processing methods such as a digital spectral shaping [7].

5. Minimum-variance mean-line subtraction

The median-line subtraction has a problem of heavy computation. A median operation is implemented by using a sorting algorithm in its computation. The sorting process is known to consume a large amount of computation power and time. For N data points, an efficient sorting algorithm requires an order of N⋅log₂ N computational complexity like a fast Fourier transform (FFT). We observed that a sorting task for 1,024 complex numbers (two separate sorting tasks in real numbers) takes eight times longer than a computation time used in an FFT task for 1,024 real numbers in Matlab^TM 7.0. This comparison was made by using the internal functions of sort() and fft() for randomly generated numbers. In most cases, it is unnecessary to extract a reference for every OCT frame. The variation of a reference spectrum is usually so slow that it just needs to be refreshed with a long interval such as tens of seconds. However, some applications, such as an OCT endoscope, may demand a more frequent reference extraction because of the relatively unstable imaging conditions.

In this section, we propose and demonstrate a faster algorithm that is almost equivalent to the median-line subtraction. Since only a part of a horizontal line is contaminated by strong reflections in most cases, an adaptive algorithm can avoid their effect on the mean line by selectively taking a mean value. After dividing a horizontal line into a plurality of segments, multiple means are calculated for the horizontal line. An optimal mean can be selected by estimating the presence of extraordinarily high amplitudes in the segments. A variance obtained from a segment is a good measure to quantify the presence of such strong reflections. The mean is taken by comparing the variances and selecting the minimum one. In our minimum-variance mean-line subtraction scheme, a segment of the minimum variance is selected for its mean value to be assigned as the mean-line value of that axial position. In other words, a mean-line value for z is determined to be one of the segmental means, ${\mathit{g}}_{\mathrm{\Omega }}^{″}$, calculated by

The minimum-variance mean-line subtraction was compared with the other processing methods. A high-contrast sample of an infrared (IR) detection card was laid over a sheet of paper and was imaged by the SD-OCT system. Figure 5 shows the OCT image processed by mean-spectrum subtraction (a), image obtained by median-line subtraction (b), and the image processed by minimum-variance mean-line subtraction (c), respectively. The IR card that had been sharply cut on one side was shown on the right-hand side of each image with a plastic film layer on the top. The flat layer located in the lower left corner is the paper laid under the IR card. Owing to the high index contrasts, a number of residual lines were produced by the mean-spectrum subtraction. They are clearly visible in Fig. 5 (a), especially in the empty space in the left of the IR card. Those lines were not visible for the cases of the median-line subtraction and the minimum-variance mean-line subtraction methods as in Figs. 5 (b) and (c).

 

Fig. 5 OCT image processed by mean-spectrum subtraction (a), image obtained by median-line subtraction (b), and the image processed by minimum-variance mean-line subtraction (c) for an IR detection card.

Download Full Size | PDF

In the extraction of the minimum-variance mean line, each horizontal line was divided into 10 segments of 100 data points. Variances were calculated for them to find the minimum variance by which a final mean was determined for an axial position. The signal processing of this method had a low computation complexity, just around three times as heavy as that of the mean-spectrum subtraction method. Due to the low computation load, a reference spectrum can be extracted from every frame to track the variation in real time without a significant computational load added. It may be natural to think that this method might produce a slightly higher error in theory compared to the median-line method because a statistical center was estimated by using a smaller number of data points (100 vs. 1,024 for our demonstrations). But we have not observed any considerable error of the minimum-variance mean-line subtraction method in this study. On the other hand, note that both of the median-line subtraction and the minimum-variance mean-line subtraction do not give much better results for a horizontal line that is fully occupied by high-amplitude signals. This may happen when an air-tissue interface is very flat and aligned perfectly horizontally without a noticeable surface roughness and tilt for the most part.

6. Conclusion

In this report, we have presented new OCT image processing methods that remove the fixed-pattern noise from an OCT image. The reference spectrum contains the full information of the fixed pattern in an FD-OCT system. Because of its variation in time, a separate measurement of the reference spectrum does not always give sufficient information to suppress the fixed pattern, especially for an OCT system working with an unstable light source or operated in a harsh environment that brings the change of a spurious interference pattern. We found that the conventional approach of the mean-spectrum subtraction method which has been utilized for real-time reference extraction has a drawback of poor accuracy associated with excessive sensitivity to high-amplitude reflection points. In this method, a residual fixed pattern can be produced by a processing error that is influenced by high-amplitude reflection points in an OCT image. Through a theoretical analysis of the image data, it was revealed that a reference spectrum extraction can be understood as a search for the center of data point distribution in the complex plane after the domain conversion. Taking advantage of this noise analysis, we proposed two new image processing schemes that exhibit higher tolerance against the high-amplitude points. In both the median-line subtraction and the minimum-variance mean-line subtraction methods, the reference A-line was extracted successfully without any remaining fixed pattern in the final images. The reference spectrum was also obtained successfully by an inverse Fourier transform of the reference A-line. Since the reference information is directly extracted from image data in our methods, the fixed-pattern noise can be removed easily without any separate calibration. Also, the obtained information on the reference spectrum can be used for other signal processing methods with better accuracy.