The cross-correlation in spectral domain based Doppler optical coherence tomography

Optics Communications 338 (2015) 433–437

Contents lists available at ScienceDirect

Optics Communications journal homepage: www.elsevier.com/locate/optcom

The cross-correlation in spectral domain based Doppler optical coherence tomography Chaoliang Chen, Jiuling Liao, Wanrong Gao n Nanjing University of Science and Technology, Department of Optical Engineering, 200 Xiao Ling Wei, Nanjing 210094, Jiangsu, China

art ic l e i nf o

a b s t r a c t s

Article history: Received 1 June 2014 Received in revised form 16 October 2014 Accepted 18 October 2014 Available online 12 November 2014

We propose a cross-correlation method which is based on the cross-correlation of two adjacent A-scans of interferogram fringe for imaging the velocities of the blood flowing in vessels. The method was tested by measurements of the velocities of flowing particles within a glass capillary with known mean velocities. Mean standard deviations of flow velocities of the particles determined through the proposed method were compared to those by the conventional phase resolved method. In vivo imaging of a mouse ear was performed and the Doppler flow velocity maps were reconstructed by both methods. The experimental results demonstrated that the proposed method can significantly suppress the phase noises caused by phase instabilities and improve the signal to noise ratio of blood flow map. & Elsevier B.V. All rights reserved.

Keywords: Optical coherence tomography Blood flow Velocity measurement Signal to noise ratio

1. Introduction Spectral domain optical coherence tomography (SDOCT) [1,2] is a non-invasive imaging technology, capable of providing high resolution, depth resolved cross-sectional images of highly scattering sample, such as biological, with high speed. In order to obtain additional biological information of tissue, many functional SDOCTs have been developed. Thanks to the advantages of the improved imaging speed and sensitivity in SDOCT, a functional extension of SDOCT, spectral domain phase resolved Doppler OCT (PRDOCT) [3,4], which relies on the phase difference between successive A-lines at the same depth, was developed to extract velocity information of blood flow in functional vessels within the scanned tissue beds. To optimize the performance of PRDOCT, Ren et al. proposed a moving-scatterer- sensitivity optical Doppler tomography (MSS-ODT) technique to improve sensitivity for imaging blood flow in vivo [5], and Wang and An also proposed a Doppler optical micro-angiography (DOMAG) method to remove the characteristic texture pattern noise caused by the heterogeneous property of sample [6]. A modified phase-resolved method has also been reported to reduce the underestimation of the velocity by an extended averaging before calculating the phase difference [7]. Other approaches have also been tried to eliminate the deterioration by evaluation of the flow velocity from fluctuations [8,9] or variance [10] of intensity instead of phase shift. However, two main factors, low signal-to-noise ratio (SNR) [11,12] n

Corresponding author. E-mail address: [email protected] (W. Gao).

http://dx.doi.org/10.1016/j.optcom.2014.10.073 0030-4018/& Elsevier B.V. All rights reserved.

and motion artifacts causing severe phase instabilities [13], deteriorate and even preclude accurate determination of velocity using the conventional phase resolved method. The SNR of phase difference is influenced by the overlap degree of light spots on the sample and time interval between adjacent A-lines [14]. As the time interval decreases, SNRs of both phase difference and intensity signal detected by CCD decrease. In this manuscript, we propose a cross-correlation method which is based on the cross-correlation of two adjacent A-scans of interferogram fringes for imaging of velocity information in a SDOCT system. The proposed method can significantly decrease the influence of noise caused by phase instabilities on the images of velocity information without need of increasing the overlap degree of light spots and time interval between adjacent A-lines.

2. Theorem of the proposed method Fig. 1 is a block diagram for computing structural and Doppler flow images from the spectral fringe intensity profiles I j (λ) using cross-correlation method, where j is the index for A-scans and λ is wavelength. First, the spectral fringe intensity profiles are converted from wavelength domain to wavenumber domain (I ′j (k)) by a standard spline interpolation algorithm to keep the axial resolution at different depth the same, where k represents wavenumber. The results I ′j (k) for the jth axial scan are correlated with the high frequency component I ″j + 1 (k) for the (j + 1)th axial scan which is extracted by a high pass filter to yield a new spectral fringe intensity profiles Fj(k). with a double number of pixels. Here,

434

C. Chen et al. / Optics Communications 338 (2015) 433–437

Fig. 1. Block diagram of signal processing for cross-correlation method.

the high pass filter is used to remove the influence of DC component in the correlation progress. In this work, the cut off frequency of high pass filter is the frequency of spectral fringes corresponding to the depth of 2 mm, which means that the spectral fringes corresponding to the depth 4 2 mm can be remained. It should be noted that the maximum of the cut off frequency is the frequency of spectral fringes corresponding to the depth of sample surface. Then the DC component Fref (k) of new spectral fringe profiles obtained by calculating the median of all new A-lines is subtracted from each new A-line profiles F j (k). The (z ) depth-dependent complex analytic signal Γ˜ (j, i) for the jth axial scan was obtained by taking the inverse Fourier transform of F ′j (k) and then removing the redundant mirror signal for z < 0, where i is the index for axial pixels. A structural image can be produced by taking the magnitude of Γ˜ (j, i) with a procedure of calibration for all A-lines within a frame. The calibration procedure is to combine two adjacent pixels as one point through averaging, which makes the number of pixels in axial direction as same as that of the conventional phase resolved method. The Doppler frequency shift fD (p, q) at any given pixel (p, q) is then calculated from Γ˜ (j, i) by [7]

⎡ p+N q+M ˜ ⁎ ˜ ⎢ Im ∑ j = p ∑i = q (Γj, i • Γj + 1, i ) fD (p, q) = tan−1⎢ p+N q+M ⁎ ⎢⎣ Re ∑ j = p ∑i = q (Γ˜j, i • Γ˜j + 1, i )

( (

= 1, 2, …,

) )

⎤ ⎥ ⎥/(2π⋅dt) ⎥⎦

N (1)

where N is the horizontal window size (determined by the sampling frequency of A-lines in the transverse direction), M is the axial window size, Γ˜j + 1, i⁎ is the conjugate of Γ˜j + 1, i , and dt is the time interval of two adjacent A-lines. The two adjacent A-scans of spectral interference fringe can be approximately regarded as from the same lateral position because the movement between two adjacent A-scans is much smaller than the lateral resolution. In the cross-correlation method, the new A-lines are obtained by cross correlation of two adjacent A-scans of spectral interference fringes, which has two advantages: first the random noise in the spectral fringe induced in the signal detection process can be effectively suppressed; second, the pixel number of new A-lines is doubled without changing the frequency of spectral fringe intensity profiles, leading to doubling the axial size of window (M ) over that in conventional phase resolved method and keeping the smallest distance between two adjacent blood vessels unchanged. Because the number of phase difference vectors used in calculating Doppler frequency shift ( fD (p, q)) with Eq. (1) is doubled, the cross-correlation method can further suppress the noise caused by phase instabilities and improve the SNR of Doppler frequency shift maps. In SDOCT system, it is clear that the depth interval ( Δh) between two adjacent pixels has to be smaller than axial resolution (Δz ). The intensity of M pixels in the range of Δz can then be approximated to represent the information of sample at the same depth. Due to the maximum frequency of spectral fringes detected

by SDOCT system is determined by CCD, according to Nyquist sampling theorem, the depth (zmax ) corresponding to maximum frequency of spectral fringes can be expressed by

z max =

Lπ 2Δk′

(2)

where L is the pixel number of CCD, Δk′ is the range of wavenumber from first pixel to last pixel of CCD. Then the depth interval (Δh) between two adjacent pixels can be expressed by

Δh =

z max L/2

(3)

Because the wavenumber Δk′ can be expressed as 2π 2π 2π ⋅ Δλ ′ , where Δλ′ is the wavelength range reΔk′ = λ − λ ≈ 2 min

max

λ0

ceived by the CCD from first pixel to last pixel, λ0 is the central wavelength of light source. The depth resolution can be expressed 2 2 ln 2 λ

as Δz = π Δλ0 , where Δλ is the full width at half maximum (FWHM) of light source. So the value of M can be determined by the following equation:

M=

Δz +1 = Δh

0. 88Δλ′ +1 Δλ

(4)

where 〈 〉 represents an operation which returns the nearest bigger integer. As pointed out in Ref. [7], the distribution of phase difference between adjacent A-lines can be expressed as a constant phasor plus a random phasor, and the probability density function can be expressed as a Gaussian function. As the mean velocity (v¯ ) of the moving particles increases, the phase difference (Δφ ) distribution moves closer to the maximum −π or π , the probability of phase wrap increases, which decreases the accuracy of velocity measurement and increase the measuring error (standard deviation, STD). The growing progress of STD can be described by the accumulated probabilities integral of Gaussian distribution, which can be expressed as Gauss error function, so the STD of calculated Doppler frequency shift versus mean velocity (v¯ ) can be given by ∞

STD =

(v¯ − vmax )2m + 1 1 1 + ⋅ ∑ ( − 1)m , (0 ≤ v¯ ≤ vmax ) π m=0 m! (2m + 1) 2

(5)

where vmax is the maximum velocity that can be measured by SDOCT which is determined by the time interval of two adjacent A-scans, m is an integer.

3. Experimental results and discussions 3.1. System setup The schematic of the fiber optics based SDOCT system used in our study is shown in Fig. 2. The broadband light source is a super luminescent diode (SLD) with a center wavelength of 830 nm and

C. Chen et al. / Optics Communications 338 (2015) 433–437

435

3.2. Experiments in vitro

SLD

λ0 = 830nm Δλ = 50nm

PC 2™2 Coupler

L3 L5

PC

L1

L2

P M Scanner

FG CCD

Sync

L4 PC

DG

Sample

Fig. 2. Schematic of the fiber optics based SDOCT system. PC, polarization controller; L1, L2, L3, L4, L5, lenses; P, prism pair; DG, transmission diffraction grating; FG, function generator; Sync, galvanometer scanner drive synchronization.

bandwidth of  50 nm. The reference arm contains a prism pair for dispersion compensation. In the sample arm, the beam is transversely scanned by using a galvanometer. The galvanometer is driven by a ramp function. Light back reflected from the sample and reference arms is combined by the Michelson interferometer and sent to an imaging spectrometer consisting of a collimating lens, a transmission diffraction grating, an achromatic focusing lens, and a fast-line-scan CCD camera (2048 pixels, 14  14 mm2 pixel size). The spectral interference fringes are digitized at 8 bit resolution. Digitized spectral fringe profiles from the camera are acquired by a frame grabber card and transferred to a computer for further signal processing. The data acquisition and transfer are triggered by a signal generated by computer synchronized with the ramp function that drives the lateral scanning galvanometer. In our experiments, the frequency of A-lines is set to be 50 kHz, and the imaging frame rate is about 13 frames/s with each frame of 2000 A-lines.

The improved accuracy in velocity measurement afforded by cross-correlation method was demonstrated by imaging a scattering sample. The sample is a glass capillary tube with an inner diameter of 300 mm, and a syringe pump controlled the flow of 2% Intralipid solution in the tube. The glass capillary was mounted at an angle of 78° to the direction of the incident light so that the Doppler frequency shift caused by moving particles can be detected in the spectral interference fringe. The SDOCT data are analyzed by using both methods with the same horizontal averaging window of 2 A-lines. Fig. 3(a) and (b) are Doppler frequency shift maps obtained by conventional phase resolved method and cross-correlation method, respectively. On comparison of the two maps, we can see that the number of noise points indicated as discontinuous changes on intensity caused by phase instabilities in Fig. 3(b) is much smaller than that in Fig. 3(a). In our SDOCT system, the value of M calculated by Eq. (4) is 3, so the axial window sizes are 3 and 6 in conventional phase resolved method and in cross-correlation method, respectively. Fig. 3(c) shows the structural image. To show the differences between results obtained by each method more intuitively, one column data corresponding to the dark line in each Doppler frequency shift map is displayed in Fig. 3(d) and (e), the horizontal axis represents depth. The intensity of velocity in Fig. 3(e) varies much more smoothly than that in Fig. 3(d) and the SNRs of Fig. 3 (d) and (e) are 53 dB and 56 dB, respectively, which shows that cross-correlation method can suppress the instable phase noises more effectively than the conventional phase resolved method. The Doppler frequency shift fD corresponding to velocity V is given by nV cos (θ) = fD λ0/2, where n is the refractive index of sample, θ is the Doppler angle, and λ0 is the central wavelength of light source. The scale bars in Fig. 3 are 100 mm. The velocity profile in capillary tube [15] can be expressed as V (z) = Vmax [1 − (2z /d)2], where z representing depth, d is the inner diameter of capillary tube and Vmax is a constant. The velocity profile in the tube can be found with known mean velocity controlled by syringe pump, and the theoretical predicted data is shown as red lines in Fig. 3(d) and (e). The STD in one column of

a

b

c

d

e

f

Fig. 3. The Doppler frequency shift maps of Intralipid in a capillary tube of a 300 mm inner diameter obtained by (a) conventional phase resolved method and (b) crosscorrelation method. (c) Structural image of scattering sample. The black dotted curves in (d) and (e) represent calculated axial Doppler frequency shift distribution at the lateral position of the black line in (a) and (b), respectively. The red curves are the theoretical predicted data. (f) Mean STDs of the middle 50 column data within Doppler frequency shift maps calculated by conventional phase resolved method (asterisks) and cross-correlation method (circles) versus mean velocity, red line and green line are fitted data of asterisks and circles respectively. Scale bars in (a)–(c) represent 100 mm. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

436

C. Chen et al. / Optics Communications 338 (2015) 433–437

a

b

c

d

Fig. 4. (a) and (b) The Doppler frequency shift maps of blood flow calculated by conventional phase resolved method and cross-correlation method, respectively. The regions marked by red arrows and white arrows are blood vessels. (c) Structural image of a mouse ear obtained by SDOCT. (d) The flow profile of blood vessels in the regions marked by the white arrows shown in (a) and (b), blue curves and red curves represent the data calculated by conventional phase resolved method and cross-correlation method, respectively. The scale bars represent 300 mm. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Doppler frequency shift maps calculated by each method can be employed to demonstrate the improved performance of crosscorrelation method quantitatively. Doppler frequency shifts corresponding to the different velocities were measured to find mean STDs by averaging the STDs of middle 50 column data within Doppler frequency shift maps. The mean STDs of the Doppler frequency shifts versus the mean velocities changed by syringe pump are given in Fig. 3(f), the asterisks and the circles correspond to the mean STDs of Fig. 3(a) and (b), respectively. In experiments, due to the frequency of A-lines is 50 kHz, and the range of mean velocity of flowing particles is from 10 mm/s to 25 mm/s which is a small part of 0–vmax , so the value of m in Eq. (5) can be approximated to be zero and the model of STD can be simplified in a linear form: STD = Av¯ + B , where A and B are constants. The red line and green line in Fig. 3(f) are fitted data of asterisks and circles, and the fitting model is the linear form mentioned above. It can be seen that the values of mean STDs grow with the increment of mean velocity, which is consistent with the theoretical results. Furthermore, the mean STDs calculated by the proposed method decrease by about 2  103, compared with those calculated by conventional phase resolved method, demonstrating that the cross-correlation method can suppress the noise caused by phase instabilities more effectively and obtain more accurate velocities than the conventional phase resolved method.

curve and red curve represent the data calculated by conventional phase resolved method and cross-correlation method, respectively. It can be seen that cross-correlation method can suppress the phase instability noise more effectively than conventional phase resolved method and improve the SNRs of Doppler frequency shift maps from 45 dB to 49 dB. When calculate Doppler frequency shift using conventional phase resolved method, the phase difference vectors over M adjacent pixels are averaged to decrease the influence of phase noise, so the distance between two adjacent blood vessels which can be separated by SDOCT system is limited. The averaging progress can be equivalent to a low pass filter, and the highest frequency in blood flow map after the filter is half of that in original map with M = 3, and the smallest distance between two adjacent blood vessels has to be bigger than 2Δh, so that they can be separated. In our SDOCT system, the smallest distance between two adjacent blood vessels is 5.5 mm. Compared with conventional phase resolved method, the phase difference vectors over 6 adjacent pixels are averaged in the proposed method, but the depth interval is only Δh/2, so the smallest distance is the same as the conventional phase resolved method.

3.3. Experiments in vivo

4. Conclusion

The improved performance of cross-correlation method is also demonstrated by imaging blood vessels in a mouse ear in vivo. The structural image, Doppler frequency shift maps, and flow profile of blood vessels in a mouse ear are shown in Fig. 4. The scale bars in Fig. 4 represent 300 mm. The range of Doppler velocity in Fig. 4(a) (conventional phase resolved method) and Fig. 4(b) (cross-correlation method) are consistent with that in Ref. [5]. The regions marked by red arrows and white arrows in Fig. 4(a) and (b) are blood vessels. Fig. 4(c) is the structural image, and the imaging area is 2.2 mm wide and 750 mm deep (optical depth in air). Fig. 4(d) shows the flow profile of blood vessels in the regions marked by the white arrows shown in Fig. 4(a) and (b), where blue

In this manuscript, we have proposed a cross-correlation method based on the cross-correlation of two adjacent A-scans of interferogram fringe to decrease the influence of noises caused by phase instabilities. The capability of the method of greatly suppressing the noise caused by phase instabilities and improving the SNR of Doppler frequency shift maps were demonstrated by the experimental results. Because the SNR of velocity maps can be improved without increasing the sampling frequency of A-lines in transverse direction in our method, resulting in an increase in imaging speed. This advantage is especially appealing for the real time imaging in DOCT system.

C. Chen et al. / Optics Communications 338 (2015) 433–437

Acknowledgments This research was supported by the National Natural Science Foundation of China (61275198 and 60978069).

References [1] A.F. Fercher, C.K. Hitzenberger, G. Kamp, S.Y. Elzaiat, Measurement of intraocular distances by backscattering spectral interferometry, Opt. Commun. 117 (1995) 43–48. [2] R. Leitgeb, C.K. Hitzenberger, A.F. Fercher, Performance of Fourier domain vs. time domain optical coherence tomography, Opt. Express 11 (2003) 889–894. [3] Z. Chen, T.E. Milner, S. Srinivas, X. Wang, A. Malekafzali, M. Gemert, S. Nelson, Noninvasive imaging of in vivo blood flow velocity using optical Doppler tomography, Opt. Lett. 22 (1997) 1119–1121. [4] Y. Zhao, Z. Chen, Z. Ding, H. Ren, S. Nelson, Real-time phase-resolved functional optical coherence tomography by use of optical Hilbert transformation, Opt. Lett. 27 (2002) 98–100. [5] H. Ren, T. Sun, D.J. MacDonald, M. Cobb, X. Li, Real-time in vivo blood-flow imaging by moving-scatterer-sensitive spectral-domain optical Doppler tomography, Opt. Lett. 31 (2006) 927–929. [6] R.K. Wang, L. An, Doppler optical micro-angiography for volumetric imaging of vascular perfusion in vivo, Opt. Express 17 (2009) 8926.

437

[7] A. Szkulmowska, M. Szkulmowski, A. Kowalczyk, M. Wojtkowski, Phase-resolved Doppler optical coherence tomography—limitations and improvements, Opt. Lett. 33 (2008) 1425–1427. [8] Y. Wang, R. Wang, Autocorrelation optical coherence tomography for mapping transverse particle-flow velocity, Opt. Lett. 35 (2010) 3538–3540. [9] M. Szkulmowski, A. Szkulmowska, T. Bajraszewski, A. Kowalczyk, M. Wojtkowski, Flow velocity estimation using joint spectral and time domain optical coherence tomography, Opt. Express 16 (2008) 6008–6025. [10] A. Mariampillai, B.A. Standish, E.H. Moriyama, M. Khurana, N.R. Munce, M. Leung, J. Jiang, A. Cable, B. Wolson, I. Vitkin, V. Yang, Speckle variance detection of microvasculature using swept-source optical coherence tomography, Opt. Lett. 33 (2008) 1530–1532. [11] B. Park, M.C. Pierce, B. Cense, S. Yun, M. Mujat, G.J. Tearney, B. Bouma, J. Boer, Real-time fiber-based multi-functional spectral-domain optical coherence tomography at 1.3 mm, Opt. Express 13 (2005) 3931–3944. [12] L. Huang, Z. Ding, W. Hong, C. Wang, T. Wu, Higher-order cross-correlationbased Doppler optical coherence tomography, Opt. Lett. 36 (2011) 4314–4316. [13] S.H. Yun, G.J. Tearney, J.F. De Boer, B.E. Bouma, Motion artifacts in optical coherence tomography with frequency-domain ranging, Opt. Express 12 (2004) 2977. [14] S. Makita, F. Jaillon, I. Jahan, Y. Yasuno, Noise statistics of phase-resolved optical coherence tomography imaging: single-and dual-beam-scan Doppler optical coherence tomography, Opt. Express 22 (2014) 4830–4848. [15] S. Yazdanfar, Noninvasive microstructural and velocity imaging in humans by color Doppler optical coherence tomography, Case Western Reserve University, Cleveland, Ohio, 2003 (Ph.D. dissertation).