1.

Introduction

Direction sensitive investigation of reversible flows using phase sensitive technique was applied in acoustic reflectometry¹ and optical coherence tomography (OCT).². Velocity image reconstruction based on estimation of phase change between sequential A-scans is represented as the color or Doppler image in both. For that purpose, cross-³ and auto-correlation^4–9 functions of the reflected signal are used. For example, the OCT system presented in Ref. 3 is based on the calculation of cross-correlation function for a certain number, $N$, of adjacent A-scans. The mean phase shift variation is used as

The OCT system presented in Ref. 5 operates in the same way. The difference is that the cross-correlation function is used^5,6 in place of the autocorrelation Kasai function.¹ Flow velocity is then calculated from the phase shift of the consecutive A-scans. This allows obtaining high sensitivity to the motion of sharp surfaces or sparse particles in a flow, but leads to the $2\pi$-ambiguity.

To reduce sensitivity to the changes of the incidence angle and to the velocity pulsations of the moving blood, Doppler spectra variance, ${\sigma }^2$, can be used as^8,9

In Refs. 10 and 11, improved Doppler optical microangiography is described. In this method, the back reflected signal only from the moving particles is registered. This leads to a considerable reduction of the noise to compare with the standard backscattered signal processing. Blood flow velocity is determined from the phase shift of the signal carrier frequency.

A double correlation technique is notable for motion artifacts’ reduction and background noise removal,¹² for which adaptive Wiener filtering is applied. It is also used to increase the structural image quality.

High sensitivity biological fluids’ velocity measurement^13,14 is based on the phase synchronization of the reference signal, $\mathrm{\Delta }{\mathrm{\Phi }}_R$, and in the sample, $\mathrm{\Delta }{\mathrm{\Phi }}_S$. An electro-optic phase modulator with periodic phase tuning is used in the reference arm.¹⁴ In practice, this time domain system allows the reduction of motion artifacts which causes interference pattern blurring, that, in turn, causes the sensitivity reduction of the flow registration. This approach is supposed to quantitatively register velocity change, but the Doppler frequency is, again, represented using the phase shift

Another approach to obtain velocity information is based on the speckle flickering frequency in the structural OCT images.¹⁵ Speckle fluctuations depend on the average flow velocity. Therefore, the ratio of the high and low components of the intensity fluctuation power spectrum could show the average velocity of a flow. The Doppler shift is then calculated using the following equation:^16,17

The description method for one specific velocity (OSV) mapping is based on the registration of the Doppler shift of the carrier, $f_c$, itself. This gives the average integrated information over the probe volume and does not depend on Doppler spectrum broadening. The latter does not give the absolute value of the flow velocity, but primarily has a statistical relationship with it rather than a functional one. In addition, OSV mapping does not have the $2\pi$ ambiguity disadvantage. The signal processing of Doppler spectra themselves is performed in the frequency domain after a short time Fourier transform (STFT).

Since it is difficult to register the Doppler shift of the carrier itself when applying ultrafast scanning techniques,^3–14 $f_D=\pm 2V\cos (\alpha )/{\lambda }_0$, we use the modified double pass rapid scanning optical delay (RSOD)^18,19 at lower scanning frequencies. RSOD allows the decoupling phase and group delay, as well as a shift of $f_c$ to the values of 20 to 30 kHz. The maximum signal-to-noise ratio of the detected photocurrent appears to be in the frequency range between the decreasing $1/f$ noise and the increasing white noise. Using a time-domain optical delay and processing the signal in the frequency domain, it is possible to decompose a standard motionless structural image, and OSV images of a chosen predetermined velocity, $V$, with a certain predetermined accuracy, $\mathrm{\Delta }V$.

In the current work, the method of one predetermined velocity mapping is applied to study a flow with complex geometry in a tilted capillary entry (tilted die entry)²⁰ at a steeper angle, $\alpha =77\mathrm{deg}$. The initial OCT image is decomposed into a structural image and two velocity images corresponding to the positive and negative directions of the flow. The final image is obtained after color coding and complexation of the three.

2.

Materials and Methods

The experimental setup (Fig. 1) is a scanning fiber-based Michelson interferometer. Superluminescent diode (SLD) irradiation (${\lambda }_{\mathrm{o}}=1300\mathrm{nm}$, $\mathrm{\Delta }\lambda =70\mathrm{nm}$) after the fiber coupler goes to the reference arm—an optical delay line (modified RSOD¹⁹). Another part of the irradiation goes to the sample arm for the investigated flow. The axial resolution of the system is about $11\mu \mathrm{m}$. Tilted capillary entry is used to model a flow with complex geometry (Fig. 2). On the way back, irradiation is focused on the balanced detectors. After an amplifier and a tunable band-pass filter (2 to 200 kHz), the electric signal goes to a 10-bit analog-to-digital converter. Final processing is performed in the frequency domain using digital STFT.

Fig. 1

Schematic representation of the experimental setup—scanning Michelson interferometer. SLD—broadband source, superluminescent diode; D—detectors in the balance scheme; ODL—modified rapid scanning optical delay line.

Fig. 2

Tilted capillary entry, $\alpha =77\mathrm{deg}$.

For $\alpha =77\mathrm{deg}$, the speed range in the hydrodynamic phantom is equal to 0 to $28\mathrm{mm}/\mathrm{s}$. A 7.5-kHz Doppler shift corresponds to a velocity of $V=17\mathrm{mm}/\sec$. The dynamic range of the system is about 100 dB. The system has been tested on the parabolic velocity profile in a cylindrical capillary with Doppler spectrum broadening $\mathrm{\Delta }f\sim 0$, 3–1, 6-kHz depending on the frequency of the centroid, $f_{\mathrm{o}}$

The standard deviation of the parabolic profile $\mathrm{\Delta }V/V$ is about 5% to 7%. A hydrodynamic phantom has been made of a transparent plastic material with two different inner diameters $d_1=2.2\mathrm{mm}$, $d_2=0.55\mathrm{mm}$. This gives a complex flow with a converging ratio of $4:1$. The outer diameter of the phantom is 5 mm. This allows the study of the flow upstream (before the entry) and downstream (after the entry), with a resolution which cannot be detected by Doppler ultrasound systems. A 1% intralipid solution in distilled water was used for obtaining the backscattered signal and velocity images of the flow. The flow rate was equal to $0.11\mathrm{ml}/\min$, the average velocity of the fluid in the large lumen ($d_1$) was about $0.48\mathrm{mm}/\mathrm{s}$ and, assuming a parabolic velocity profile, the axial flow velocity was $V_{1\max }=0.97\mathrm{mm}/\mathrm{s}$. In the lumen with the smaller diameter ($d_2$), the average flow rate was about $7.7\mathrm{ml}/\mathrm{s}$, and $V_{2\max }=15.4\mathrm{mm}/\mathrm{s}$. Since the percentage of the added intralipid was low, it was assumed that the density $\rho$ and the dynamic viscosity of the solution were similar to those in water at 20°C.

The applied optical delay allows shifting the carrier to the region with minimal noise with a sufficiently fast repetition rate of A-scans (several kilohertz). This allows obtaining good enough structural images and registers the Doppler shift of the carrier corresponding to the indicated velocities. Odd A-scans always have a positive Doppler shift, while in the even A-scans, it is always negative. This approach gives the possibility of splitting the initial image into the three images corresponding to the steady object and images showing the positive and negative directions of the flow (positive and negative OSV images).

3.

Velocity Mapping Algorithm

Using the aforementioned approach, one specific velocity mapping algorithm was suggested. It includes: (1) separation of the input data into two parts which correspond to the positive and negative shifts of the carrier frequency; (2) the independent construction of a two-dimensional (2-D) structural OCT image and (3) 2-D OSV images with (4) consecutive color coding.

$x[s]$ is the input set of data, where $s=\mathrm{1,2},\dots ,S$ and $S$ is the number of samples in the signal. It is received from the raw data file and represents the amplified and digitized mean intensity of the radiation detected by the system.

Then it is separated to the functions $x_{\text{odd}}[p,k]$ and $x_{\text{even}}[p,k]$ corresponding to the positive and the negative shifts of the carrier frequency as follows:

The resulting signal $x_{\text{sum}}[p,k]$ is the algebraic sum of $x_{\text{odd}}[p,k]$ and $x_{\text{even}}[p,k]$:

Note that from the $x_{\text{sum}}[p,k]$ set the analog band-pass filter of the system separates the high- and low-parts of the spectrum.

To control the mapping process, two parameters have been chosen: the upper $V_{\text{upper}}$ and lower $V_{\text{lower}}$ limits of the mapped velocity

On that basis, the upper ${\omega }_{\text{upper}}$ and the lower ${\omega }_{\text{lower}}$ cut-off frequencies of the band-pass filters were determined:

They are necessary for determination of the frequency band from the functions $x_{\text{odd}}[p,k]$ and $x_{\text{even}}[p,k]$.

To compare with the previously described algorithm,²⁰ this algorithm has better tuning of the chosen velocity diapason. It is determined not using the central velocity and its accuracy, but directly by independently choosing the upper and the lower velocity values. In addition, the ${\omega }_{\text{upper}}$ and ${\omega }_{\text{lower}}$ determinations have been also digitally improved. Currently, they are directly determined from the Doppler angle and the source central wavelength.

Functions $x_{\text{odd}}[p,k]$, $x_{\text{even}}[p,k]$, and $x_{\text{sum}}[p,k]$ are processed in three stages: (1) short-time Fourier transform; (2) Hilbert transform, and (3) taking the logarithm of the signal envelope. After that functions $x_{\text{odd}}[m,k]$, $x_{\text{even}}[m,k]$, and $x_{\text{sum}}[m,k]$ are obtained and considered. They represent three B-scans and correspond to the three parts of the image, Figs. 3(a) to 3(c) and 4(a) to 4(c), respectively. Note that $m=\mathrm{1,2},\dots ,M$, where $M$ is the number of samples per one A-scan after the above processing.

Fig. 3

Structural and one specific velocity (OSV) images of chosen velocities in a flow with a parabolic profile. Initial image (a); structural image of the tube (b); set aside velocity image (c); image of the parabolic profile with band-pass filters 6 to 9 kHz (d); 9 to 12 kHz (e); 12 to 15 kHz (f); 14 to 17 kHz (g); complexed image for 9 to 12 kHz filter (h). Dimensions of the pictures are $0.7\times 0.7\mathrm{mm}$.

Fig. 4

Standard structural image (a); image of odd A-scans (b); image of even A-scans (c); OSV image of positively directed flow (d); OSV image of negatively directed flow (e); OSV image of bidirectional flow (f); the resulting image—(g). Dimensions of the pictures are $0.7\mathrm{mm}\times 0.7\mathrm{mm}$.

Automatic identification of the directions of the flow $x_{\text{trends}}[m,k]$ is performed as the element-wise subtraction of $x_{\text{even}}[m,k]$ from $x_{\text{odd}}[m,k]$ with a shift by one column:

Then B-scans corresponding to the flows $x_{\mathrm{pos}}[m,k]$ and $x_{\mathrm{neg}}[m,k]$ are color coded considering the experimentally determined minimal level of the constant signal, $L$, which is separated from the noise (Figs. 3 and 4):

After that functions $x_{m\mathrm{pos}}[m,k]$ and $x_{m\mathrm{neg}}[m,k]$ are complexed into a common array (OSV map) representing the flow direction and location $x_{\mathrm{cart}}[m,k]$ [Fig. 4(f)]:

Function $x_{\text{sum}}[m,k]$ is encoded to be at least 20 times higher than the functions $x_{\mathrm{cart}}[m,k]$ and then is complexed with them:

The resulting function $x_{\text{last}}[m,k]$ is normalized to the color scale, which includes shades of green to visualize the image of the phantom structure; red and blue colors are used to visualize the direction of the flow; the white color is used as a reference point (no movement, no reflection). To minimize complexation distortions, a proportionality factor between the $x_{\text{sum}}[m,k]$ and $x_{\mathrm{cart}}[m,k]$ values equal to 20 was estimated as the best one. This allows simultaneously increasing the OSV accuracy values by 18% and the computation speed by 13%. The resulting image is presented in Figs. 3(h) and 4(g).

4.

Results and Discussion

The described algorithm is realized in a form of a software package in graphical programming language «G» of LabVIEW. A simple model of a blood vessel was used for software testing and debugging. It is a transparent tube having an inner diameter of 0.55 mm. Figure 3 shows a structural image obtained by the system—(a); structural image of the tube—(b); structural image of a selected speed in the flow—(c) and final complexed image—(h).

The OSV image of the parabolic profile with other band-pass filters: 6 to 9 kHz (3.7 to $5.5\mathrm{mm}/\mathrm{s}$)—Fig. 3(d), 9 to 12 kHz (5.5 to $7.3\mathrm{mm}/\mathrm{s}$)—Fig. 3(e), 12 to 15 kHz (7.3 to $9.1\mathrm{mm}/\mathrm{s}$)—Fig. 3(f), 14 to 17 kHz (8.5 to $10.4\mathrm{mm}/\mathrm{s}$)—Fig. 3(g), demonstrate the capabilities of the proposed algorithm.

Figure 4(a) shows the results obtained when scanning was performed at a distance of 0.5 mm from the entry. The original structural image contains information about the steady parts as well as the positive and negative directions of the flow (bright stripes). But it is not possible to visually distinguish one from the others, even by changing the analog and digital filters. In the automated mode, the proposed algorithm allows the separation of initial image into three parts for about 1 to 2 s. The proposed method of OSV Doppler mapping splits the original B-scan into three parts and independently reconstructs the structural image and 2-D color-coded OSV images [Figs. 4(b) to 4(e)] with their subsequent complexation [Figs. 4(f) and 4(g)]. By choosing the upper, $V_{\text{upper}}$, and the lower, $V_{\text{lower}}$, limits of the mapped velocity, it is possible to change the predetermined specific velocity and the accuracy of the equivelocity map.

5.

Conclusion

The method of Doppler mapping of a chosen velocity is applied on the basis of low coherence depth discrimination of the reflected signal and detection of the introduced carrier frequency changes. It gives 2-D maps with a quantitative functional (rather than statistical) relation between the Doppler shifts of the carrier and flow velocity without a $2\pi$ ambiguity. Additional color coding enables automated separation of the standard OCT structural image in to three parts: steady image and two velocity images, even in the more complex flow used to compare to the previously described image.²⁰ Note that velocity images are predetermined and quantitatively mapped and are proportional to the Doppler shift of the carrier itself. Colors are used as an additional option of the algorithm to separate a stationary object, as well as the positive and negative directions of the flow. It is possible to apply this approach in the study of small animals’ hearts where similar flows with complex geometry appear in nature. Since color (Doppler) ultrasound systems have similar disadvantages as the low coherence optical ones, we suggest applying this algorithm in Doppler ultrasound diagnostics to study bigger biomedical objects.