CC BY
29
An improved algorithm of structural image reconstruction with rapid scanning optical delay line for Optical
Coherence Tomography
*
Denis A. Petrov , Saif N. Abdulkareem, Kamal E. S. Ghaleb, and Sergey G. Proskurin
Biomedical Engineering department, Tambov State Technical University, Tambov, Russia
*
e-mail: [email protected]
Abstract. A new algorithm of structural image reconstruction in Optical Coherence Tomography is described. The modified rapid scanning optical delay (RSOD) line, low numerical aperture, small angle raster scanning with consecutive averaging and multilevel digital filtering have been used to obtain high quality structural images of an onion and the nail bed of a human thumb. The proposed method significantly improves image contrast and allows visualization of small blood capillaries under the nail plate. © 2016 Journal of Biomedical Photonics & Engineering.
Keywords: optical coherence tomography (OCT), coherence probing depth (CPD), rapid scanning optical delay (RSOD), small angle raster averaging, nail bed, blood capillaries.
Paper #2363 received 2015.03.17; revised manuscript received 2016.03.01; accepted for publication 2016.03.01; published online 2016.03.20. doi: 10.18287/JBPE16.02.020302
References
1. A. F. Fercher, C. K. Hitzanberger, W. Drexler, G. Kamp, and H. Sattman, "In Vivo Optical Coherence Tomography," Am. J. Ophtalmol. 116(1), 113-114 (1993).
2. D. Huang, E. A. Swanson, C. P. Lin, J. S. Schuman, W. G. Stinson, W. Chang, M. R. Hee, T. Flotte, K. Gregory, C. A. Puliafito, and J. G. Fujimoto, "Optical coherence tomography," Science 254(5035), 1178-1181 (1991).
3. G. Liu, W. Qi, L. Yu, and Z. Chen, "Real-time bulk-motion-correction free Doppler variance optical coherence tomography for choroidal capillary vasculature imaging," Opt. Express 19(4), 3657-3666 (2011).
4. S. G. Proskurin, A. Y. Potlov, K. E. S. Ghaleb, and S. N. Abdulkareem, "Structural Image Reconstruction of Biological Object by Using Raster Averaging in Optical Coherence Tomography," Bulletin of YuFU, Technical Sciences 9, 129-134 (2012).
5. S. G. Proskurin, "Raster scan and averaging for speckle reduction in optical coherence tomography," Quantum Electronics 42(6), 495-499 (2012).
6. R. Vinayakrishnan, Speckles in laser Doppler perfusion imaging, Ph.D Thesis, University of Twente, Enschede, The Netherlands (2007).
7. M. Y. Kirillin, G. Farhat, E. A. Sergeeva, M. C. Kolios, and A. Vitkin, "Speckle statistics in OCT images: Monte Carlo simulations and experimental studies," Optics Letters 39(12), 3472-3475 (2014).
8. S. G. Proskurin, and S. V. Frolov, "Visualization of blood vessels using optical coherence tomography," Biomedical Engineering 46(3), 9-13 (2012).
9. B. Marco, A. J. Kennerley, I. Meglinski, and S. Matcher, "Application of Doppler optical coherence tomography in rheological studies: blood flow and vessels mechanical properties evaluation," Journal of Innovative Optical Health Sciences 2(4), 431-440 (2009).
10. S. G. Proskurin, A. Y. Potlov, and S. V. Frolov, "One specific velocity color mapping using optical coherence tomography," Journal of Biomedical Optics 20(5), 051034 (2015).
11. T. Q. Xie, Z. G. Wang, and Y. T. Pan, "High-speed optical coherence tomography using fiberoptic acousto-optic phase modulation," Opt. Express 11(24), 3210-3219 (2003).
12. S. G. Proskurin, "Optical Coherence-domain Imaging of Subcutaneous Human Blood Vessels in vivo," Advances in Life Sciences 1(2), 40-44 (2012).
13. A. Rollins, M. Kulkarni, S. Yazdanfar, R. Ung-arunyawee, and J. Izatt, "In vivo video rate optical coherence tomography," Opt. Exp. 3(6), 219-229 (1998).
14. X. Liu, M. J. Cobb, Y. Chen, M. B. Kimmey, and X. D. Li, "Rapid-scanning forward-imaging miniature endoscope for real-time optical coherence tomography," Opt. Lett. 29(15), 1763-1765 (2004).
15. A. F. Fercher, W. Drexler, C. K. Hitzenberger, and T. Lasser, "Optical coherence tomography - principles and applications," Rep. Prog. Phys. 66(2), 239-303 (2003).
16. R. Koprowski, and Z. Wrobel, "Image Processing in Optical Coherence Tomography: Using Matlab.Katowice," Poland: University of Silesia, 2011.
17. Y. Jian, K. Wong, and M. V. Sarunic, "Graphics processing unit accelerated optical coherence tomography processing at megahertz axial scan rate and high resolution video rate volumetric rendering," J. Biomed. Opt. 18(2), 026002 (2013).
18. J. M. Schmitt, S. H. Xiang, and K. M. Yung, "Speckle in Optical Coherence Tomography," Journal of Biomedical Optics 4(1), 95-105 (1999).
19. X. Zha, Image processing of optical coherence tomography for image guided dental drilling, Master thesis, Delft university of technology (2012).
1 Introduction
The signal expression measured by an OCT system is derived, which reveals the possibility of tissue dispersion compensation by introducing the required amount of dispersion in the reference arm and may be implemented by incorporating the grating- based rapid scanning optical delay (RSOD) lines in the reference arm of OCT.
Since the early nineties when this technique became commercially available its characteristics such as sensitivity, field of view, coherence probing depth (CPD), resolution and speed of image acquisition are continuously improved [1].
Biological tissues characterized by a high scattering coefficient and therefore on the great depth most part of the light are scattered and only a little part of it escapes object without being scattered, so OCT can be efficient only at small depth at one or two millimeters. But instead optical coherence tomography can provide much better (less than 10 micrometers) resolution, so that it could be used to visualize even tiny biological structures.
Typical OCT system (Figure 1) consists of Michelson interferometer with low coherence light source (supercontinuum lasers or superluminescent diodes). Light emitted by the source divided by the beam splitter and directed to the reference and the sample arm where the investigated object is located. Light reflects from optical heterogeneities inside sample arm and from mirror inside reference arm and returns back to the beam splitter which directs it to the photodetector. After that the acquired data is processed and displayed as a 2D image of the inner structures of the object under investigation [2].
OCT systems could provide high resolution images where even smallest blood vessels could be visualized. But due to imperfections of the setup components and mistakes in the process of studying objects like bulk motions of the patient acquired images could have low contrast so it could be very hard to distinguish necessary details of the object, which could be very important in the biomedical diagnostics [3]. Therefore, increasing the image contrast is an important procedure that should be applied in the data processing algorithm. One way to do it is to apply small angle raster scanning and averaging over several adjacent A-scans [4,5]. Also if we want to
Fig. 1 Schematic of the OCT system.
get as much useful information as possible, processing steps should be flexible so we can use it for studying variety of biomedical objects; and fast, so we can use it to process images in real-time and at the video rate. This technique along with digital filtering is described in this paper.
Reference and sample arm electric waves are given
as:
E = E. + E cos(kz + m,
r in ro \ r Tl/
E = E + E cos (kz + m ),
s in so \ s T2/'
Interference signal appeared in the detector is: I =< EE >= I. + E2 + E2 +
r s in ro so
+2EroEso cos[2/ + (^ - = = I. +1 +1 + 2J1 I cos(2nA/Am(i)).
in ro so y ro so v J T \ yy
(1)
(2)
Incoherent component could appear as the phase noise and reduce quality of the acquired images [6,7]. Currently it is not possible to completely reduce its value to zero, but reducing radiation power of the low coherent source would sufficiently decrease its ratio. Application of the small angle raster scanning with consecutive averaging leads to additional phase noise reduction [8]. At the same time this increases signal-to-noise ratio (SNR), reduces speckle noise and leads to image contras enhancement. After eliminating direct current and filtering the low frequency 1/f noise the electric signal from (1) and (2) is represented as follows:
Id(t) = Iin +10 cos[(œ0 - mD) t + A9(i)].
(3)
In a linear optical delay line (ODL) the reference arm mirror scans with constant velocity, ±V, in the fiber optic Michelson interferometer giving depth-resolved signal and, at the same time, can introduce the carrier frequency fc = 2!/! determined by the classical Doppler shift of the reflected light [9, 10].
For rapid scanning optical delay line (RSOD) bandwidth is expressed as:
Af = K
AXda X2dt '
(4)
where da / dt is the angular velocity of the scanning mirror, a is the instantaneous scanning angle and K is a coefficient which depends on RSOD implementation. Application of the rapid scanning grating based ODL makes it possible to separate group and phase delay, keeping respectively high scanning speed (1-100 kHz) and tuning the carrier to the frequency range with minimal noise (fc ~25 kHz) enabling considerable increase of SNR.
The purpose of this paper is to describe efficient real time image processing algorithm using RSOD based OCT system based on multilevel analog and digital data filtration.
2 Materials and methods
In our work we use the modified RSOD [11]. Decoupling of the group and phase delay enables to use it together with ultrafast methods (swept source and spectral domain systems). That enables to make image acquisition at video rate in real time mode while keeping carrier frequency respectively low, within 2030 kHz. The modified RSOD can reduce reference arm radiation power by one-tenth of a milliwatt [12].
In this work we present algorithm of image processing for time domain optical coherence tomography (TDOCT). TDOCT couldn't provide fast acquisition speed as Fourier domain optical coherence tomography but it still can make it fast enough to make data acquisition at video rate [13]. Also images acquired with TDOCT using improved RSOD could provide information of even smallest structures which sometimes can't be visualized with some modern OCT techniques [5,8,14].
Other important feature is using relatively small numerical aperture (NA). The diameter of the beam waist expressed as:
ro = -
1.22 A, 2[ N .A.]'
(5)
Decreasing NA and increasing confocal parameter to a size comparable with coherence probing depth (to 1.5-2.5 mm), we increase the sample beam waist up to 50-100 microns. This could affect spatial resolution and reduce it by more than 3 times, but for multi layered objects this resolution is quite enough to obtain good and clear image [15].
The applied algorithm presented in the Fig. 2. It is necessary to distinguish important stages: "bandpass filtering", "splitting the signal", "threshold filtering", "Fourier transform", "raster averaging" and "median filtering".
Butterwoth Sisnal Threshold
filtering splitting filtering STFT
High quality image
Median Raster
filtering averaging
Fig. 2 Improved algorithm of structural image reconstruction.
The first step is to filter the signal by a bandpass Butterworth filter. The purpose of this step is to remove the high- and low-frequency noise to smooth the noise of the interference pattern. The next processing step is to divide the interference signal into equal duration
segments. Note that, the data have a discrete spectrum due to the use of digital methods of its registration that leads to a number of restrictions specific to the discrete representation of continuous data.
After that threshold filtering is performed for each segment of the interference signal. This filtering has a high efficiency because usually interference has much lower amplitude than the amplitude of the desired signal. The next stage of the data processing is the application of the short time Fourier transform (STFT) to each segment. Since the path difference of the interferometer arms changes continuously, STFT window must also be equal to the coherent envelope and shifted continuously, for one digital point. It has been shown empirically, however, that the signal processing with the window shifted by 15-35% gives the same image quality as a continuous window shifting.
Performing averaging of neighboring rows and columns of the image pixels gives additional advantage in the algorithm. Number of columns that are unified into one is defined as a function of the number of A-scans, and for 900 A-scans it is 3-10. Number of rows that are unified into one is determined depending on the probing depth and often in between 3 and 7. This considerably increases image contrast and significantly reduces the noise. The final stage of the processing algorithm is a graphical and median filtering which is used to suppress the additive and impulse noise of the image.
One of the key features of the algorithm is use of the parallel computing using central processing unit (CPU). Modern CPU contains 2-128 cores which allow parallelization of the image processing algorithm. Most of the image processing algorithms for optical coherence tomography described in numerous papers are using only single core while the other available cores are inactive [16]. This occurs most often in case
than MATLAB or Lab VIEW computing environment used because most of the inner functions in these packages are single threaded and using all CPU cores means use of additional Toolboxes. Some other algorithms implement parallelization using Graphics Processing Unit (GPU) which allows sufficiently increasing speed of image processing but it means additional money cost for buying GPGPU compatible video card [17]. Programming languages which allow multi-threaded calculation such as C# allows making image processing at video rate speed (24 frames per seconds) using even inexpensive CPU with 2-4 cores onboard.
3 Results and discussion
The described algorithm was applied for two-dimensional mapping of an Allium cepa onion bulb (Fig.3) and human nail bed (Fig.4). There are several advantages we can get with the modified RSOD: transmission bandwidth sufficiently improved its value fourfold, incident radiation power was decreased to 0.4 mW and it enables to shift carrier frequency to region of several tens of kilohertz what helps to avoid noise increasing at the higher frequencies. Decreasing radiation power sufficiently reduces value of the inherent component in acquired interference signal. Applying low numerical aperture helps to register mainly photons which are reflected from the layered structures and to avoid multiply scattered in the tissues photons. Photons which exit from the biological tissue after multiple scattering do not contain information about structure of the object that leads to additional increasing of SNR of the demonstrated images. [18]. Nail bed capillaries what can be seen with a naked eye appear on the images after application of the described image processing only.
Fig. 3 Structural images of an onion bulb (Allium cepa) 180 A-scans (a) and 900 A-scans before (b) and after applying the proposed algorithm (c). The image size is 2x2 mm.
a b
Fig. 4 Structural images of a human thumb nail in vivo. 180 A-scans before (a) and after the averaging (b). Blue circle marks the nail bed capillaries. The image size is 1x2 mm.
Table 1 Comparison of performance of OCT image processing modalities.
Technique Proposed algorithm CPU Koprowski R., Wröbel Z. CPU [16] Enhanced Lee filter CPU [19] Gamma Map Filter CPU [19] Proposed algorithm GPU Yifan Jian etc GPU [17]
Processing time per one frame (seconds) 0.027 0.27-1.2 0.13 0.138 27-29 volumes/s (size 1024x256x200) 26 volumes/s (size 1024x256x200)
Signal-to-noise ratio was calculated using the ratio of the mean value of maximal signal intensity to standard deviation as described earlier [5]. Value of SNR for the presented images increased by 20-30% after the application of the described algorithm.
Acquired structural images of in vivo human nail of a thumb before and after using the proposed algorithm are presented in Fig.4. Fig.3 shows 900 A-scans of the structural image of the onion bulb (Allium cepa). After data processing with the improved algorithm the contrast of the image has been considerably improved. Averaging over 5-7 A-scans sufficiently improves image quality in both cases. In case of 180 A-scans image resolution of the picture after processing is good enough. Averaging over 10 and more A-scans blurs the image and leads to the loss of useful information especially on the bigger depth.
The processing has been performed using program written on the C# programming language using a computer with a standard computational power (four-core processor 3.4 GHz, 2 GB RAM). Computation performed using each core of the CPU with Threading
namespace which included in C# standard. Processing time per image is 0.027 seconds, which suggests the possibility of using presented algorithm for image processing and representation of the acquired images in video mode.
Performance of the algorithm compared to some other image processing modalities is presented in the Table 1. Speed of the algorithm using CPU compared to other techniques has highest performance and allows to process images at the rate higher than 24 frames per second. However Koprowski et al. [16] approach includes layer recognition modality which is not implemented in our algorithm, but can be added to it after the digital filtering. In addition Koprowski et. al. program is written for MATLAB package which uses single threaded functions mostly, this can significantly reduce calculation speed. Our algorithm can also be implemented for GPU. In this case, a slight increase in performance to compare with other techniques based on the use of GPU can be seen. However time domain OCT systems cannot provide such data acquisition speed for three dimensional real time data
representation so our image processing algorithm without technical improvements of TDOCT system described above can be used for swept source and Fourier domain systems what will be the subject of our future research. All algorithms were applied using Intel Core i5-4670 processor and GTX 780 video controller using MATLAB, C# and CUDA C software packages.
4 Conclusion
In this paper we have described high quality OCT image processing algorithm. Multilevel filtering, low numerical aperture, small angle raster scanning enables visualization even small blood capillaries, which is hard to distinguish without the described procedures. New processing algorithm with GPU parallelization allows representing images in real time mode. Opportunity of localization of such small nail bed capillaries allows studying blood flow in them thus extracting velocity information in vivo. Detection of the phase shift and
reconstruction functional images instead of structural ones will enable to represent the direction of blood flow in them. Application of the multilevel filtration algorithm for carrier frequency Doppler shift registration and interferogram phase shift information could also give advantages to study blood flow in the lager vessels.
By using this technique, the frequency change in the reference arm is not required. Consequently, the technique reduces the requirements for the detection bandwidth of the photodetector and simplifies data acquisition and signal processing during real-time imaging. Presented in vivo results demonstrate that this technique might be useful in the clinical management of patients when blood-flow monitoring is essential.
