Quantitative lateral and axial flow imaging with optical coherence microscopy and tomography

Arno Bouwens; Daniel Szlag; Maciej Szkulmowski; Tristan Bolmont; Maciej Wojtkowski; Theo Lasser

doi:10.1364/OE.21.017711

1. Introduction

Optical coherence tomography (OCT) is an interferometric imaging technique that allows three-dimensional label-free imaging of biological tissues and cells [1]. Besides structural imaging, OCT has been used extensively for functional imaging, where blood flow imaging in living tissue is arguably one of the most important applications. Blood flow imaging methods based on OCT can be divided into two categories: angiographic methods that image the structure of the vasculature in tissue [2–8], and velocimetric methods measuring the velocity of red blood cells [9–18]. While the former techniques use some form of motion contrast to yield a qualitative image of blood vessels, the latter allow a quantitative estimation of the blood flow speed. Most quantitative methods are ultimately based on measuring the Doppler frequency shift of light scattered by moving particles, but use different approaches to estimate the speed of the scatterers from the acquired data. All these methods require the acquisition of a time-series of A-scans or B-scans. Hence, the highest Doppler frequency that can be measured without aliasing is given by half the acquisition rate. The most straightforward approach is to evaluate the local phase difference between time samples, from which the Doppler frequency is then obtained. In joint spectral and time domain OCT, Fourier transformations are performed along both time and wavenumber dimensions, yielding the spatially resolved Doppler spectrum [14, 19, 20]. At each location, the Doppler shift is estimated as the frequency corresponding to the maximum magnitude position or to the center of mass using circular statistics [21]. Autocorrelation-based methods involve a similar acquisition scheme but calculate the Doppler shift by searching the maximum in the power spectral density calculated as the Fourier transform of the autocorrelation [16].

In all these quantitative methods, the estimated Doppler frequency is eventually used to calculate the axial velocity component of moving scatterers. However, in a number of biological applications (e.g. neurology and ophthalmology), many capillaries are oriented nearly perpendicular to the optical axis, making the axial flow component too small to be measured in small vessels such as capillaries. It is therefore crucial to provide a reliable measurement of the lateral component.

The lateral velocity component can be estimated from the width of the Doppler frequency spectrum [22–24, 24–27]. Indeed, because wave vectors in the beam illuminating the sample and in the collected backscattered light have an angular spread with respect to the optical axis, they will experience different Doppler shifts, causing broadening of the Doppler spectrum. This effect is called spectral broadening due to beam geometry [28]. Proskurin et al. argued that the finite transit time of scatterers through the focal volume causes spectral broadening [25]. This effect can be shown to be equivalent to broadening due to beam geometry [28]. The equivalence can be understood intuitively by considering the effects of increasing the numerical aperture (NA) of an optical system. A higher NA results in a larger angular spread of wave vectors, increasing broadening of the Doppler spectrum due to the beam geometry argument. Equivalently, the increased angular spread of wave vectors yields an improved lateral resolution, thereby reducing transit time through the focus, which results in broadening of the Doppler spectrum.

A conceptually different method to measure axial and lateral flow velocity components as well as diffusion was presented by Lee et al. [44]. Their method is based on dynamic light scattering measurements with OCT and involves fitting the autocorrelation of the OCT signal to a model. While this method allows to extract a large amount of information about the local flow properties, the non-linear fitting of the model’s many parameters requires an intricate algorithm resulting in long computation times.

Phase noise in the interferometer creates Doppler broadening even in the absence of flow, thus imposing a limit on the minimum detectable lateral flow velocity [14, 29]. Therefore, the larger the Doppler frequency broadening for a given flow velocity vector, the better the sensitivity to small flow changes. Higher NA increases Doppler broadening due to the lateral flow component [23]. Hence, to image slow lateral flow or to detect small flow changes, a higher NA is advantageous. In this case, optical coherence microscopy (OCM) can have an advantage over OCT. For OCM systems employing Gaussian beams however, this results in a reduced depth of focus. In this work, we use extended focus optical coherence microscopy (xfOCM) to counter this problem. xfOCM uses a Bessel mode to illuminate the sample and a Gaussian mode to collect scattered light [30]. The Bessel-Gauss configuration achieves a high lateral resolution of 1.3 μm, maintained over an extended depth of focus of 400 μm.

To date, the models developed to analyze Doppler broadening in OCT have all assumed Gaussian modes and light source spectra. In this work, we develop a model that allows arbitrary light source spectra, and arbitrary illumination and detection modes. This allows modeling of arbitrary axial and lateral parts of the point spread function, respectively. Given the optical system’s parameters, our model can be used to correctly combine the measured central Doppler frequency and Doppler broadening to provide quantitative images of axial and lateral flow velocities. To validate the model, we perform quantitative imaging of lateral and axial flow in a capillary channel using both an OCT system with Gaussian modes and an xfOCM instrument.

2. Theory

2.1. Fourier domain coherent imaging for static scattering

Typical Fourier domain implementations of OCT record an interferogram as a function of wavenumber at each position in a lateral raster scan. After background subtraction and assuming the backscattered light intensity is small compared to the reference light intensity, the signal of interest in each wavenumber channel is

I (k) \propto ℜ [U_{r}^{*} (k) U_{m} (k)],

where U_r is the reference field amplitude and U_m the scattered field amplitude measured by the system. Each wavenumber channel can be treated as being monochromatic. The optical system is characterized by its (achromatic) illumination mode m_i(a, b) and detection mode m_d(a, b), parametrized in the entrance principal plane with coordinates p = (a, b), as in Fig. 1(a).

Fig. 1 Schematic representations of wave vectors and the principal planes of an objective. (a) The axes and definitions used in the model. P and Q are the entrance and exit principal planes of the objective, respectively. (b) A low NA system has a small angular spread of wave vectors; lateral flow is therefore not detected because the flow velocity vector v_n is then perpendicular to K = k_i(p_i) − k_d(p_d). Equal Gaussian illumination and detection modes are shown in the schematic. (c) Higher NA systems have a larger angular spread of wave vectors and the scalar product is no longer always zero, allowing the detection of lateral flow.

Download Full Size | PDF

We will assume that the objective lens obeys Abbe’s sine condition [31]. For scalar waves and under the Debye approximation [32], the incident light field in the focal volume at a position r measured from the focal point, can be written as a superposition of plane waves with different angular directions [33]:

U_{in} (r; k) = \frac{j k}{2 π f} S (k) \int_{P} m_{i} (p_{i}) exp [- j r \cdot k_{i} (p_{i})] d^{2} p_{i},

where f is the focal length of the objective, k_i(p_i) describes the illumination wave vectors at the exit principal plane and S(k) is the spectral envelope of the field. In the notation U_in(r;k), r is a variable for the function U_in whereas k is a parameter indicating the (monochromatic) wavenumber channel. The direction cosines of k_i are functions of p_i = (a_i, b_i):

α_{i} = \frac{a_{i}}{f}, β_{i} = \frac{b_{i}}{f}, γ_{i} = \sqrt{1 - α_{i}^{2} - β_{i}^{2}} .

The plane waves are scattered upon interaction with the sample. Using the first order Born approximation [34], the scattered field in the principal plane is

U_{s} (p_{d}; k) = \int_{V} U_{in} (r; k) F (r; k) \frac{exp (- j k | r_{d} (p_{d}) - r |)}{| r_{d} (p_{d}) - r |} d^{3} r,

where r_d(p_d) describes the exit principal plane of the objective (as shown in Fig. 1(a)) and F(r;k) is the scattering potential. The scattering potential can be expressed in terms of the sample susceptibility (assumed to be non-dispersive):

F (r; k) = \frac{k^{2}}{4 π} χ (r) = \frac{k^{2}}{4 π} (n^{2} (r) - 1)

Because |r_d| = f ≫ |r|, we can write

U_{s} (p_{d}; k) \approx \frac{e^{- j k f}}{f} \int_{V} U_{in} (r; k) F (r; k) exp [j r \cdot k_{d} (p_{d})] d^{3} r

Finally, the coupling into the detection mode needs to be considered:

U_{m} (k) = \int_{P} m_{d} (p_{d}) U_{s} (p_{d}; k) d^{2} p_{d} .

For a detailed analysis of the tomogram reconstruction from Eq. (7), we refer to the works of Villiger et al. [35] and Sheppard et al. [36]. In the case of NA → 0, both k_d and k_i are parallel to the optical axis, as shown in Fig. 1(b). Then, k_d = −k_i = k1_z and hence the measured scattered field would be:

U_{m} (k) \propto k S (k) \int_{V} F (r; k) exp (- 2 j k z) d^{3} r .

If the scattering potential is taken to be a set of scattering layers, the classical layer model for OCT image formation is retrieved. In this work however, we will start from Eq. (7) without further approximation because the NA is crucial to measure lateral and axial flow components.

2.2. Fourier domain coherent imaging for dynamic scattering

When scatterers in the sample are moving, the scattering potential becomes a function of time and can be separated in a static and dynamic part:

F (r; t, k) = F_{s} (r; k) + F_{d} (r; t, k) .

The dependence of the dynamic term F_d on the time parameter t is assumed to be on a much slower time-scale than the period of the light wave, such that the monochromatic approximation leading to Eq. (4) remains valid. A sample containing different ensembles of uniformly moving rigid scatterers can be modeled by:

F_{d} (r; t, k) = \sum_{n = 0}^{N} F_{d, n} (r - v_{n} t; k),

where the v_n are the velocity vectors of the ensembles and n runs over the contributing ensembles. To evaluate the volume integral in Eq. (6) for each of the ensembles of moving scatterers, we can adopt a coordinate system r₀ = r − v_nt that is moving along with the ensemble. The time-dependent part U_m(t;k) of the measured field is then:

U_{m} (t, k) \propto \sum_{n = 0}^{N} k S (k) \int_{V} d^{3} r_{0} [F_{d, n} (r_{0}; k) \int_{P} m_{d} (p_{d}) \int_{P} m_{i} (p_{i}) exp [- j (r_{0} + v_{n} t) \cdot K] d^{2} p_{i} d^{2} p_{d}],

where K ≡ K(p_d, p_i) = k_i(p_i)−k_d(p_d). In Eq. (11), the Fourier transformed scattering potential can be identified:

{\tilde{F}}_{d, n} (K) = \int F_{d, n} (r_{0}; k) e^{- j r_{0} \cdot K} d^{3} r_{0} .

From the measurement of the field U_m, OCT systems attempt to reconstruct the depth resolved scattering potential of the sample at each position in a lateral raster scan. Typically, an inverse Fourier transform of Eq. (1) is calculated over k to this end. For this reconstruction to be accurate at some distance from the focal point, either the range of K(p_d, p_i) must be limited by m_d(p_d) and m_i(p_i) such that r₀ · K effectively measures the optical path length along the depth direction (i.e. r₀ · K ≈ 2zk), or advanced reconstruction algorithms need to be applied. The former solution is valid within good approximation in the focal volume of low NA or extended focus systems [35]. The latter solution has been investigated under the name of interferometric synthetic aperture microscopy (ISAM) [37–39]. In this work, we will assume either one of the solutions is adopted to ensure correct reconstruction. Defocussing effects will therefore not be considered.

Finally, the time-dependent part of the measured field can be rewritten as:

U_{m} (t, k) \propto \sum_{n = 0}^{N} k S (k) \int_{P} m_{d} (p_{d}) \int_{P} m_{i} (p_{i}) {\tilde{F}}_{d, n} (K) exp (- j v_{n} t \cdot K) d^{2} p_{i} d^{2} p_{d} .

In this expression, the functions m_d(p_d) and m_i(p_i) can be interpreted as distribution functions for the variables p_d and p_i. Integration over these variables creates the Doppler spectrum as follows. For each realization of p_d and p_i, the corresponding illumination and detection wave vectors can be found, as shown in Fig. 1(a). Their direction cosines are given by Eq. (3). Every such pair of illumination and detection wave vectors determines a Doppler frequency shift per wavenumber channel:

f_{D} (v_{n}, p_{d}, p_{i}) = \frac{1}{2 π} v_{n} \cdot K (p_{d}, p_{i}) .

The weight of each frequency in the Doppler spectrum is then given by:

m_{d} (p_{d}) m_{i} (p_{i}) {\tilde{F}}_{d, n} (K) .

In the extreme case of zero NA shown in Fig. 1(b), we find that the distribution reduces to a single Doppler frequency, as expected. Moreover, lateral flow will not be detected in that case, since v_n · K = 0. Conversely, when the NA is increased, K spans a range of angles due to the spread of wave vectors in the illumination and detection modes, as shown in Fig. 1(c). Therefore, a broadened Doppler frequency distribution will be seen and lateral flow can be measured. An analysis of the Doppler spectrum allows the measurement of both flow components, as will now be shown.

2.3. Doppler spectrum per wavenumber channel: mean and standard deviation

To analyze the Doppler spectrum per wavenumber channel, the weighted mean μ[f_D;k] and standard deviation σ[f_D;k] of the Doppler frequency are evaluated:

μ [f_{D}; k] = \frac{\int_{P} \int_{P} f_{D} (v_{n}, p_{d}, p_{i}) m_{d} (p_{d}) m_{i} (p_{i}) d^{2} p_{i} d^{2} p_{d}}{\int_{P} \int_{P} m_{d} (p_{d}) m_{i} (p_{i}) d^{2} p_{i} d^{2} p_{d}},

σ^{2} [f_{D}; k] = μ [f_{D}^{2}] - μ^{2} [f_{D}] .

For the sake of simplicity, a single point scatterer has been taken into account (F̃_d,n(K) ∝ k² exp[− jr_n · K]). In other words, back-scattering is approximated to be uniform within the system’s aperture. The model in principle also allows different types of scatterers by incorporating F̃_d,n(K) in the Doppler frequency weight, but such an analysis is beyond the scope of this work. Since most optical systems exhibit rotationally symmetric illumination and detection modes, polar coordinates (ρ, θ) will be used to calculate the integrals in the entrance principal plane. Considering one ensemble of point scatterers with velocity vector v = (v_x, v_y, v_z), we find for the mean Doppler frequency for one wavenumber channel:

μ [f_{D}; k] = \frac{n k}{2 π f} v_{z} (\frac{M_{f, i}}{M_{1, i}} + \frac{M_{f, d}}{M_{1, d}}),

where n is the refractive index of the medium and the following notation for integrals of the illumination mode has been adopted:

M_{1, i} = \int_{0}^{R} ρ m_{i} (ρ) d ρ,

M_{3, i} = \int_{0}^{R} ρ^{3} m_{i} (ρ) d ρ,

M_{f, i} = \int_{0}^{R} \sqrt{f^{2} - ρ^{2}} ρ m_{i} (ρ) d ρ,

with R the radius of the entrance pupil. Analogous notation is used for integrals of the detection mode. Eq. (18) confirms that the mean Doppler frequency is determined solely by the axial flow component. Moreover, when a moderate NA objective is used, m(ρ) is significantly different from zero only for ρ ≪ f and therefore M_f ≈ fM₁. This approximation is further verified in Fig. 2(b). The mean Doppler frequency then reduces to the known formula:

μ [f_{D}; k] \approx n v_{z} k / π .

In general, the variance of the Doppler frequency is dependent on all flow components:

σ^{2} [f_{D}; k] = \frac{n^{2} k^{2}}{4 π^{2} f^{2}} \frac{A}{2} (v_{x}^{2} + v_{y}^{2}) + \frac{n^{2} k^{2}}{4 π^{2} f^{2}} C v_{z}^{2},

with

A = \frac{M_{3, i}}{M_{1, i}} + \frac{M_{3, d}}{M_{1, d}}, B = \frac{M_{f, i}^{2}}{M_{1, i}^{2}} + \frac{M_{f, d}^{2}}{M_{1, d}^{2}}, C = 2 f^{2} - A - B .

The magnitudes of the ratios M_f / M₁ and M₃ / M₁, are shown in Fig. 2 for Gaussian and Bessel modes. The dependence of the Doppler broadening on A, f and B can be interpreted as follows. As mentioned in the introduction, it has been shown that broadening of the Doppler spectrum due to beam geometry is equivalent to broadening due to the finite transit time of scatterers moving through the focal volume [28]. Hence, Doppler broadening due to lateral flow is higher when the lateral resolution is higher. Higher lateral resolution requires increasing the NA, which in turn will increase A, as shown in Fig. 2(c). This explains the presence of A in the first term of Eq. (23). The second term measures broadening due to pure axial flow, which is caused by the angular spread of wave vectors with respect to the axial component of the velocity vector. Wave vectors with different angles measure different Doppler frequencies. Equivalently, this can be interpreted by the finite depth of focus of the system, which limits the focal volume axially. The depth of focus therefore creates Doppler broadening due to the finite axial transit time of scatterers moving through the focal volume. For high NA Gaussian systems, the depth of focus will be reduced and the factor C increases, as shown in Fig. 2(d). However, this effect is only of importance when the flow is mostly axial, since C ≪ A/2, as can be seen by comparing Fig. 2(c) with (d). Note that for extended focus systems with ideal Bessel illumination and detection modes (m_i,d(ρ) = δ (ρ − ρ₀)), we find C = 0 for all NA.

Fig. 2 (a) Illumination or detection mode in the principal plane for Gaussian mode (blue) and Bessel mode (green) in a system with NA = 0.3. (b) The ratio M_f / M₁ normalized to the focal length f for Gaussian and Bessel modes as a function of NA. For lower NA, M_f / M₁ approach f. (c) M₃ / M₁ normalized to f² as a function of NA. (d) C = 2f² − A − B normalized to f² as a function of NA for Gaussian illumination and detection modes. Note that C ≪ M3/M1.

Download Full Size | PDF

Figure 3 shows the mean and standard deviation of the Doppler frequency calculated using Eq. (18) and Eq. (23) for different angles of the velocity vector, and for different NAs. The mean Doppler frequency shown in Fig. 3(c) is seen to be almost independent of the NA, as predicted by the approximation leading to Eq. (22). For pure axial flow, however, we see that the mean Doppler frequency is slightly lower for higher NA. This is due to the smaller mean z-component of K for higher NA. On the other hand, a higher NA increases the Doppler broadening, as evidenced in Fig. 3(d). When a typical OCT system of NA = 0.05 is compared to a system with NA = 0.3, the Doppler broadening for lateral flow is increased about sixfold. Fig. 3(d) confirms that even for pure axial flow, some Doppler frequency broadening remains, especially for higher NA. This is caused by the second term of Eq. (23), as explained above.

Fig. 3 (a) The velocity vector v of a point scatterer makes an angle α with the optical axis. (b) Equal Gaussian illumination and detection modes for an optical system with f = 16.4 mm. NA = 0.05, 0.3 and 0.5 for the blue, green and red curves respectively. (c) Mean Doppler frequency as a function of α. Parameters: λ = 780 nm, |v| = 1 mm/s, n = 1.33. (d) Standard deviation of the Doppler frequency.

Download Full Size | PDF

2.4. Spatially resolved Doppler spectrum

To yield depth-resolved imaging of the Doppler spectrum, we build on the method of joint spectral and time domain OCT (jstdOCT) [14]. After subtraction of the background spectrum, the interferograms recorded by the spectrometer are first mapped from the wavelength to the wavenumber space. Then, the inverse Fourier transform is calculated along the wavenumber axis:

{\tilde{U}}_{m} (t, l) = \frac{1}{2 π} \int U_{m} (t, k) e^{j 2 k l} d k,

where 2l represents optical path length, and the factor 2 takes into account the double path due to the reflection configuration. As discussed above, we will assume that Ũ_m(t, l) correctly represents the depth resolved sample structure. The Doppler spectrum can be estimated from the Fourier transform of Eq. (25) along t. After this Fourier transform, we have a measurement of the Doppler spectrum for each depth resolution. Indeed, the integral over k in Eq. (25) selects one resolution volume at depth l. More details on the full algorithm are given in section 3.4. The integration over the k-spectrum will cause the Doppler frequency distribution to broaden. Indeed, f_D in Eq. (14) depends linearly on k, the scattering potential in the Doppler frequency weight in Eq. (15) contains the factor k² and upon inspection of Eq. (13), it is seen that the additional weighting factor kS(k) must be taken into account. The effect on the weighted mean and standard deviation of the Doppler frequency for a point scatterer as calculated in section 2.3, can now be evaluated.

The mean Doppler frequency can be calculated using the law of total expectation:

μ [f_{D}] = \frac{\int k^{3} S (k) μ [f_{D}; k] d k}{\int k^{3} S (k) d k} .

In the case of point scatterers, we can substitute Eq. (18) for μ[f_D;k],

μ [f_{D}] = \frac{n}{2 π f} v_{z} (\frac{M_{f, i}}{M_{1, i}} + \frac{M_{f, d}}{M_{1, d}}) \frac{\int k^{4} S (k) d k}{\int k^{3} S (k) d k} .

As before, for moderate NA, M_f ≈ fM₁. When S(k) is a Gaussian function with mean k₀ and standard deviation k_σ, the moments of S(k) can easily be calculated. Moreover,

k_{σ}^{2} ≪ k_{0}^{2}

for typical OCT/OCM spectra, which allows to simplify the ratio of the fourth to the third moment of S(k) to k₀. With these two approximations, Eq. (27) reduces again to nv_zk₀/π.

Similarly, the variance is found using the law of total variance:

σ^{2} [f_{D}] = \frac{\int k^{3} S (k) σ^{2} [f_{D}; k] d k}{\int k^{3} S (k) d k} + σ_{k}^{2} [μ [f_{D}; k]] .

Substituting Eq. (23) for σ²[f_D;k],

\begin{array}{l} σ^{2} [f_{D}] & = & \frac{n^{2}}{4 π^{2} f^{2}} (\frac{A}{2} (v_{x}^{2} + v_{y}^{2}) + C v_{z}^{2}) \frac{\int k^{5} S (k) d k}{\int k^{3} S (k) d k} \\ + & (\frac{\int k^{5} S (k) d k}{\int k^{3} S (k) d k} - {(\frac{\int k^{4} S (k) d k}{\int k^{3} S (k) d k})}^{2}) {(\frac{n}{2 π f} v_{z} (\frac{M_{f, i}}{M_{1, i}} + \frac{M_{f, d}}{M_{1, d}}))}^{2} . \end{array}

To interpret the different terms of Eq. (29), we first assume S(k) to be Gaussian with central wavenumber k₀, standard deviation k_σ and

k_{σ}^{2} ≪ k_{0}^{2}

, which is valid for typical OCT/OCM spectra. Finally, with M_f ≈ fM₁ for moderate NA we find:

σ^{2} [f_{D}] \approx \underset{T_{1}}{\underset{︸}{\frac{n^{2} k_{0}^{2}}{4 π^{2} f^{2}} (\frac{A}{2} (v_{x}^{2} + v_{y}^{2}) + C v_{z}^{2})}} + \underset{T_{2}}{\underset{︸}{k_{σ}^{2} \frac{n^{2} v_{z}^{2}}{π^{2}}}} .

The two terms of Eq. (30) can be interpreted as follows. The first term T₁ is the standard deviation for a single wavenumber channel at the central wavenumber as given by Eq. (23); it accounts for broadening due to the angular spread of incident and scattered wave vectors, as discussed above. Note that T₁ is independent of spectral width. The second term T₂ measures an extra broadening due to the axial flow component, and depends on the spectral width k_σ. It can be interpreted using the transit time argument: a larger spectral bandwidth creates a higher axial resolution, resulting in a shorter axial transit time, and therefore a larger Doppler broadening. Indeed, the axial resolution is given by the coherence length l_c = 1/(nk_σ).

Figure 4 shows the mean and standard deviation (from Eq. (29)) of the Doppler frequency for different NA and light source bandwidth. The mean Doppler frequency is independent of these two parameters, as seen in Fig. 4(a) and (c), and as predicted by the approximations. The behavior of the standard deviation in Fig. 4(b) and (d) can be understood in terms of spatial resolution and the transit time argument. When axial resolution (determined by the source spectrum) and lateral resolution (determined by the NA) are equal, the Doppler frequency broadening is equal for axial and lateral flow components, since the transit time through the focal volume is then independent of flow direction. Conversely, an asymmetric focal volume will have higher Doppler broadening for either axial flow (better axial than lateral resolution) or for lateral flow (better lateral than axial resolution).

Fig. 4 (a) Mean Doppler frequency for a system with equal Gaussian illumination and detection modes (Gauss-Gauss). Parameters: NA = 0.05, |v| = 1 mm/s, λ₀ = 780 nm, n = 1.33. Color code for the graphs is given in panel b. (b) Standard deviation of the Doppler frequency for different spectral bandwidths Δλ (full-width at half-maximum of Gaussian S(λ)). NA = 0.05. (c) Mean Doppler frequency for a Gauss-Gauss system with NA = 0.3. (d) Standard deviation of the Doppler frequency. NA = 0.3.

Download Full Size | PDF

2.5. Approximations for Gaussian modes and low NA

To compare our formula in Eq. (30) to previous work [23,24,27], we need to set the illumination and detection modes to narrow (low NA) and equal Gaussians with waist w₀ ≪ R at the principal plane. Under this condition, $A = 2 w_{0}^{2}$ and C ≈ 0. The beam waist after focussing by the objective is then w′₀ = 2 f / (nk₀w₀), yielding for the Doppler broadening:

σ^{2} [f_{D}] \approx \frac{v_{t}^{2}}{π^{2} {w_{0}^{'}}^{2}} + \frac{v_{z}^{2}}{π^{2} l_{c}^{2}},

where

v_{t} = {(v_{x}^{2} + v_{y}^{2})}^{1 / 2}

is the transverse flow component. The first term in Eq. (31) corresponds to the formula for Doppler broadening as given by Ren et al. [23] and Piao etal. [24]. Their works do not include broadening due to the axial flow component. Srinivasan et al. calculate the broadening from the field autocorrelation for equal axial and lateral resolution [27]; their result is equivalent to Eq. (31). In conclusion, with the right approximations for low NA OCT, our model reduces to the results in previously published work. Table 1 summarizes the different approximations for the mean and standard deviation of the Doppler spectrum.

Table 1. Approximations and corresponding equations for μ[f_D] and σ²[f_D].

View Table

2.6. Doppler spectrum: numerical simulation

To investigate the shape of the Doppler spectrum, a numerical simulation is performed. In the entrance principal plane, two sets of uniformly distributed random coordinates are chosen: p_i for the illumination mode and p_d for the detection mode. For each combination of p_i and p_d, the Doppler frequency and its weight are calculated using Eq. (14) and Eq. (15), including the spectral weighting factor for point scatterers k³S(k) for a Gaussian spectrum. This calculation is repeated for each wavenumber channel. The resulting data are then binned according to their frequency and represented in a weighted histogram. Figure 5 reveals how the Doppler spectrum changes as a function of flow angle for a system with equal Gaussian modes and a system with a Bessel illumination mode and Gaussian detection mode. It is seen that the Bessel-Gauss configuration slightly increases the Doppler broadening for the same NA.

Fig. 5 (a) Mean Doppler frequency as a function of flow velocity angle α. Blue and green curves represent a system with equal Gaussian modes (GG) and a system with Bessel illumination mode and Gaussian detection mode (BG) respectively. Asterisks correspond to values obtained by the numerical simulations shown in panels c and d, and solid lines are calculated using Eq. (27) and Eq. (29). Parameters: NA = 0.3, |v| = 1 mm/s, λ₀ = 780 nm, Δλ = 120 nm, n = 1.33. (b) Standard deviation of the Doppler frequency. (c) Doppler spectra for different flow velocity angles α, calculated by numerical simulation for a system with equal Gaussian modes. (d) Doppler spectra for a system with a Bessel illumination mode and Gaussian detection mode, NA = 0.3.

Download Full Size | PDF

3. Methods

3.1. xfOCM

The xfOCM set-up is based on a Mach Zehnder interferometer and has decoupled illumination and detection paths. Light from a broad-bandwidth source (Ti-Sa laser, Femtolasers, λ₀ = 780 nm, Δλ = 120 nm) is collimated and split into reference and sample beams. In the sample arm, the beam passes through an axicon lens and a telescope. The resulting Bessel-like illumination beam is guided through a relayed beam scanning system. The scanned beam enters a microscope stage, consisting of a tube lens and a 10× Zeiss Neofluar objective with a NA of 0.3. The Bessel-like illumination beam provides a uniform, high lateral resolution (1.3 μm) over an extended depth of field (400 μm). The axial resolution in tissue is 2.5 μm. The system operates with a Bessel illumination mode and a Gaussian detection mode, as shown in Fig. 2(a). A custom-built spectrometer with a high-speed line detector (Basler Sprint spL4096-140km) records the interference spectrum. For the flow measurements presented here, an A-scan rate of 20 kHz was used.

3.2. OCT

The OCT system is based on a Michelson interferometer with equal Gaussian illumination and detection modes and employs a broad-bandwidth source (Ti-Sa laser, Femtolasers, λ₀ = 780 nm, Δλ = 120 nm). The sample objective lens has a focal length of f = 25 mm. The Gaussian mode has a full width at half maximum of 2.4 mm in the back focal plane of this lens, creating a lateral resolution of 10 μm. A custom-built spectrometer with a high-speed line detector (Basler Sprint spL4096-140km) records the interference spectrum. For the flow measurements, an A-scan rate of 100 kHz was used. More details about this set-up can be found in a previous publication [40].

3.3. Flow system

The flow system consists of a glass capillary tube with a measured cross-section of 4.95 × 10⁴ μm² for xfOCM and 3.14 × 10⁴ μm² for OCT experiments. A syringe pump allows to set a constant flow speed through the capillary. The system is filled with a solution of polystyrene beads (diameter 0.2 μm) in water for the xfOCM-based measurements and with an intralipid solution for the OCT-based measurements.

3.4. Data acquisition & algorithm

A lateral raster scan is performed with oversampling (xfOCM: 16× Nyquist frequency; OCT: 125× Nyquist frequency) along the fast scan axis. Due to the oversampling, repeated A-scans are taken within the same lateral resolution volume such that they can be seen as separated in time, not space. The acquired interference spectra I (λ, x(t)) from each lateral resolution volume are first windowed along the time-dimension using a Hanning window, yielding I_w(λ, x, t). Then, following the joint spectral and time-domain method [14], the spectral background is subtracted, the spectra are mapped to the wavenumber domain and the absolute value of the two-dimensional fast Fourier transform (2D-FFT) is calculated, yielding the Doppler frequency distribution as a function of depth. Figure 6 summarizes the data acquisition protocol and algorithm. For the xfOCM measurements a time-window size of 32 points is used, for OCT 256 points. To estimate the mean Doppler frequency, the center of mass of the Doppler spectrum is calculated using circular statistics [21]. The Doppler frequency standard deviation can be estimated by calculating the weighted standard deviation directly, but we have found the following fitting procedure to be more reliable in the presence of noise. First, to prevent aliased Doppler frequencies from biassing the fit, each Doppler spectrum is circularly shifted to have its center of mass at 0 Hz. Then, the Doppler spectrum is fitted with a non-linear least squares algorithm (Matlab, The Mathworks, Inc.) to

p (f_{D}) = a exp (- f_{D}^{2} / b^{2} - f_{D}^{4} / c^{4}) + d .

This empirical function was found to produce a good fit for both OCT and xfOCM Doppler frequency distributions, as shown in Fig. 6(e) for typical xfOCM data. The fitting parameter c is included to fit the (broader) xfOCM distributions; for OCT we use c → ∞. Indeed, the numerical simulations in Fig. 5(c) and Fig. 5(d) show that Doppler spectra in xfOCM are broader than in OCT. The Doppler frequency standard deviation is then calculated using p(f_D) as the weighting function. Lateral and axial flow components can then be found from the mean and standard deviation of the Doppler frequency using Eq. (27) and Eq. (29). To calculate the system-related parameters (i.e. the spectrum, illumination and detection modes) needed in these equations, we used the design parameters of the systems. A Matlab implementation of the above algorithm is provided, see the supporting information for details. Note that the same principles can be applied along the slow scan axis, by acquiring repeated B-scans; see Grulkowski et al. for a detailed description of scan protocols [20].

Fig. 6 (a) Overview of the algorithm used to obtain functional tomograms containing spatially resolved Doppler frequency spectra and flow velocity components. (b) The interference spectra I(λ, x) are acquired in an over-sampled lateral (x) scan. δx represents the lateral resolution. (c) A window is applied around each x-position in the lateral scan such that the windowed data represent time samples I_w(λ, t) at that position. (d) Spatially resolved Doppler spectra Ĩ(z, f_D) are obtained after spectral background subtraction, λ to k mapping and 2D FFT. (e) The center of mass μ of the measured Doppler spectrum (blue graph) is calculated and the re-centered data are fitted to Eq. (32) (green graph). The standard deviation σ of the Doppler spectrum is calculated using the fitted curve as weighting function. (f) μ and σ are converted to lateral and axial flow components v_t and v_z using Eq. (27) and Eq. (29), or their appropriate approximations.

Download Full Size | PDF

4. Experimental results

To validate the model and algorithm, we performed flow measurements at different flow rates set by a syringe pump through a glass capillary using two optical set-ups: xfOCM and OCT. The lateral scanning was perpendicular to the flow direction, such that its influence on Doppler broadening could be removed using $v_{t}^{2} = {\hat{v}}_{t}^{2} - v_{s}^{2}$ , where v_t is the true transverse flow velocity, v̂_t the transverse velocity measured from the Doppler broadening and v_s the scan speed. Another solution would be to make v_s ≪ v_t, or to stop the scan at each lateral position.

4.1. xfOCM

Figure 7(a) shows a structural tomogram of the capillary, and Fig. 7(b) and (c) show tomograms of the mean μ[f_D] and standard deviation σ[f_D] of the Doppler frequency for a flow rate of F = 0.45 ml/h and angle α = 81°. As can be seen in 7(a), the capillary is somewhat asymmetric, with an elliptic cross-section. The axial diameter is 200 μm, but the lateral diameter is larger, accounting for a larger cross-section. The intensity gradient over the depth of the flow channel seen in 7(a) is due to three reasons. Firstly, although the extended focus maintains the lateral resolution over 400 μm, it does not have constant intensity over the whole depth, a well-known property of Bessel beams [42]. Secondly, the Gaussian detection mode employed by xfOCM creates and additional intensity gradient over depth. Thirdly, scattering of the illumination beam by the sample reduces the intensity with depth.

Fig. 7 (a) xfOCM tomogram of the capillary used for flow measurements, dynamic range is 20dB. The channel has an angle α = 81° and the flow rate is F = 0.45 ml/h. The lateral sampling frequency is 16 times the Nyquist frequency. (b) xfOCM tomogram of the mean Doppler frequencies. (c) xfOCM tomogram of the standard deviation of the Doppler frequencies. 10 B-scans were averaged for the images in (b) and (c).

Download Full Size | PDF

For each lateral scan position, μ[f_D] and σ[f_D] are converted to lateral and axial flow profiles as a function of depth, as shown in Fig. 8(a) and (b) for the center of the capillary. By assuming laminar flow, the volume flow rate F set by the syringe pump is converted to flow velocities in the capillary. The maximum velocity v_max of the parabolic flow profile in a cylindrical capillary is then found from v_max = 2F/S, with S the measured cross-section of the capillary. The flow velocities set by the syringe pump can now be compared to the velocities measured using our model and algorithm, as shown in Fig. 8(c) and (d) for an angle α = 81°. The measurements of both axial and lateral flow components are seen to be consistent with the values expected from the flow system parameters. The measured flow speeds in Fig. 8(a) and (b) show some deviation from the parabolic flow profiles most likely due to experimental noise and/or deviations from the supposed laminar flow profile. Also, note that when a particular axial velocity is higher than the expected value, the corresponding transverse velocity component is higher as well, indicating that these errors are most probably due to an inaccuracy in the actual flow rates set by the syringe pump.

Fig. 8 Depth profiles at the center of the capillary (α = 81°) measured with xfOCM of (a) the axial flow component and (b) the lateral flow component, for different flow rates. Parabolae are fitted to the measurements, assuming v = 0 at the capillary wall. The maxima of the parabolae are extracted and compared with the expected velocity set by the syringe pump under laminar flow conditions for (c) the axial flow velocity component and (d) the lateral component. Error bars represent ± standard deviation over 10 measurements.

Download Full Size | PDF

4.2. OCT

Figure 9(a) and (b) shows the depth profiles measured with the OCT system on a capillary with angle α = 87°. Due to the lower lateral resolution of this system (i.e. smaller Doppler broadening), it is not as sensitive to slow lateral flow velocities as xfOCM. Hence, the flow rates used for these measurements were about ten times higher and the flow angle smaller to create larger lateral flow velocities. Nevertheless, the measurements of both axial and lateral flow components in Fig. 9(c) and (d) are again seen to match the values expected from the flow system parameters.

Fig. 9 Depth profiles at the center of the capillary (α = 87°) measured with OCT of (a) the axial flow component and (b) the lateral flow component, for different flow rates. The maxima of the parabolae are compared with the expected velocity set by the syringe pump for (c) the axial flow velocity component and (d) the lateral component. Error bars represent ± standard deviation over 10 measurements.

Download Full Size | PDF

5. Discussion & conclusion

We have introduced a general model for the Doppler frequency spectrum caused by moving scatterers that is applicable to any OCT or OCM system, given its parameters. Once all system related factors have been calculated (this must be done only once per set-up), Eq. (27) and Eq. (29) (or their appropriate approximations, see table 1) yield a straightforward relation between mean and standard deviation of the Doppler spectrum and the flow velocity components.

For the specific case of low NA OCT and Gaussian source spectra, our model yields the same equations as previously published models [23,24,27]. The model presented here therefore provides a generalization of these earlier models to arbitrary point spread functions, through the possibility of arbitrary illumination and detection modes, and light source spectra, respectively. This generalization is important for high NA OCM systems as well as xfOCM and dark-field OCM [45]. Additionally, the model does not make assumptions about the flow angle.

While point scatterers were investigated here, the effect of different types of scatterers could be investigated by including the scattering potential F̃_d,n(K) as a weight in the calculation of the mean and standard deviation of the Doppler spectrum. In this way, the Doppler spectrum created by moving erythrocytes in blood could be investigated.

The model assumes that the scattering structures are rigid bodies, moving with a constant velocity during the acquisition time. This assumption excludes effects such as Brownian motion or significant velocity gradients within a resolution volume. For the model to be valid, these effects must be small in comparison to the average velocity.

To measure the Doppler spectrum, we presented an algorithm building on the joint spectral and time domain method, thereby adopting the advantages of this method in terms of precision and signal-to-noise ratio over phase-resolved methods [14]. To estimate the mean frequency of the Doppler spectrum, either circular statistics or the standard definitions can be used [21]. A fitting procedure was used to allow a reliable estimation of the Doppler frequency standard deviation. For an interesting comparison of Doppler frequency estimators in the presence of noise, we refer to the recent work by Chan et al. [41]. An interesting subject for future research would be to evaluate and compare the performance of OCT/OCM systems in terms of quantitative flow imaging. To this end, the influence of noise should be included in the model.

Our model and algorithm were successfully applied to image axial and lateral flow in a capillary using two different system configurations. No system modifications (digital or optical) to the point spread function are needed to use our method, only the scan protocol needs to be adapted for over-sampling of the lateral resolution. As mentioned before, the scan protocol can be modified to perform Doppler spectral analysis along the slow scan axis, allowing to concentrate on a lower set of velocities [20, 27]. By allowing arbitrary pupil functions, this work now enables high resolution extended-focus systems to be used for quantitative imaging of both axial and lateral flow components. Especially the ability to quantify the lateral flow component will be instrumental in a number of biological applications, such as retina or brain imaging where blood vessels are often oriented perpendicular to the optical axis. Furthermore, xfOCM has been shown to achieve a resolution and penetration depth comparable to multiphoton microscopy but with significantly faster acquisition rates [43]. The additional capability to quantify blood flow velocities will allow the investigation of the three-dimensional vascular functionality at speeds currently unattainable by multiphoton microscopy. Our future work will therefore be focused on the quantitative imaging of lateral and axial flow in biological tissue using the concepts developed here.

Supporting information

A Matlab implementation of the algorithm developed here together with example data can be downloaded from our website at http://lob.epfl.ch/ under “Research/Quantitative flow imaging”.

Acknowledgments

This work was supported by the Commission for Technology and Innovation (CTI, Bern, grant 13964.1 PFLS-LS) and the Swiss National Science Foundation (grant 205321L-135353/1). D. Szlag acknowledges support from the SCIEX-NMSch grant. M. Szkulmowski acknowledges the National Centre for Research and Development Grant No. PBS1A9/20/2013 and from the Polish Ministry of Science and Higher Education (years 2011–2015). M. Wojtkowski and M. Szkulmowski acknowledge the TEAM project financed by the European Union in the frame of Innovative Economy coordinated by the Foundation for Polish Science. We thank Corinne Berclaz and Christophe Pache for their help with the xfOCM set-up, Iwona Gorczynska for making the OCT setup available for measurements and Daniel Ruminski for his help during the OCT measurements.

References and links

1. A. F. Fercher, W. Drexler, C. K. Hitzenberger, and T. Lasser, “Optical coherence tomography - principles and applications,” Rep. Prog. Phys. 66, 239–303 (2003) [CrossRef] .

2. S. Makita, Y. Hong, M. Yamanari, T. Yatagai, and Y. Yasuno, “Optical coherence angiography,” Opt. Express 14, 7821–40 (2006) [CrossRef] [PubMed] .

3. R. K. Wang, S. L. Jacques, Z. Ma, S. Hurst, S. R. Hanson, and A. Gruber, “Three dimensional optical angiography,” Opt. Express 15, 4083–97 (2007) [CrossRef] [PubMed] .

4. J. Fingler, D. Schwartz, C. Yang, and S. E. Fraser, “Mobility and transverse flow visualization using phase variance contrast with spectral domain optical coherence tomography,” Opt. Express 15, 12636–53 (2007) [CrossRef] [PubMed] .

5. A. Mariampillai, B. A. Standish, E. H. Moriyama, M. Khurana, N. R. Munce, M. K. Leung, J. Jiang, A. Cable, B. C. Wilson, I. A. Vitkin, and V. X. D. Yang, “Speckle variance detection of microvasculature using swept-source optical coherence tomography,” Opt. Lett. 33, 1530–2 (2008) [CrossRef] [PubMed] .

6. L. An and R. K. Wang, “In vivo volumetric imaging of vascular perfusion within human retina and choroids with optical micro-angiography,” Opt. Express 16, 11438–11452 (2008) [CrossRef] [PubMed] .

7. V. J. Srinivasan, J. Jiang, M. Yaseen, H. Radhakrishnan, W. Wu, S. Barry, A. Cable, and D. Boas, “Rapid volumetric angiography of cortical microvasculature with optical coherence tomography,” Opt. Lett. 35, 43–45 (2010) [CrossRef] [PubMed] .

8. C. Blatter, B. Grajciar, C. M. Eigenwillig, W. Wieser, B. R. Biedermann, R. Huber, and R. A. Leitgeb, “Extended focus high-speed swept source OCT with self-reconstructive illumination,” Opt. Express 19, 12141–55 (2011) [CrossRef] [PubMed] .

9. J. A. Izatt, M. D. Kulkarni, S. Yazdanfar, J. K. Barton, and A. J. Welch, “In vivo bidirectional color Doppler flow imaging of picoliter blood volumes using optical coherence tomography,” Opt. Lett. 22, 1439–41 (1997) [CrossRef] .

10. Z. Chen, T. E. Milner, D. Dave, and J. S. Nelson, “Optical Doppler tomographic imaging of fluid flow velocity in highly scattering media,” Opt. Lett. 22, 64–66 (1997) [CrossRef] [PubMed] .

11. R. A. Leitgeb, L. Schmetterer, C. K. Hitzenberger, A. F. Fercher, F. Berisha, M. Wojtkowski, and T. Bajraszewski, “Real-time measurement of in vitro flow by Fourier-domain color Doppler optical coherence tomography,” Opt. Lett. 29, 171–3 (2004) [CrossRef] [PubMed] .

12. B. Vakoc, S. Yun, J. de Boer, G. Tearney, and B. Bouma, “Phase-resolved optical frequency domain imaging,” Opt. Express 13, 5483–93 (2005) [CrossRef] [PubMed] .

13. A. H. Bachmann, M. L. Villiger, C. Blatter, T. Lasser, and R. A. Leitgeb, “Resonant Doppler flow imaging and optical vivisection of retinal blood vessels,” Opt. Express 15, 408–22 (2007) [CrossRef] [PubMed] .

14. M. Szkulmowski, A. Szkulmowska, T. Bajraszewski, A. Kowalczyk, and M. Wojtkowski, “Flow velocity estimation using joint spectral and time-domain optical coherence tomography,” Opt. Express 16, 6008–6025 (2008) [CrossRef] [PubMed] .

15. R. K. Wang and L. An, “Doppler optical micro-angiography for volumetric imaging of vascular perfusion in vivo,” Opt. Express 17, 8926–40 (2009) [CrossRef] [PubMed] .

16. V. J. Srinivasan, S. Sakadzić, I. Gorczynska, S. Ruvinskaya, W. Wu, J. G. Fujimoto, and D. A. Boas, “Quantitative cerebral blood flow with optical coherence tomography,” Opt. Express 18, 2477–94 (2010) [CrossRef] [PubMed] .

17. B. Baumann, B. Potsaid, M. Kraus, J. Liu, D. Huang, J. Hornegger, A. Cable, J. Duker, and J. Fujimoto, “Total retinal blood flow measurement with ultrahigh speed swept source/Fourier domain OCT,” Biomed. Opt. Express 2, 1539–1552 (2011) [CrossRef] [PubMed] .

18. S. Zotter, M. Pircher, T. Torzicky, M. Bonesi, E. Götzinger, R. A. Leitgeb, and C. K. Hitzenberger, “Visualization of microvasculature by dual-beam phase-resolved Doppler optical coherence tomography,” Opt. Express 19, 1217–1227 (2011) [CrossRef] [PubMed] .

19. A. Szkulmowska, M. Szkulmowski, D. Szlag, A. Kowalczyk, and M. Wojtkowski, “Three-dimensional quantitative imaging of retinal and choroidal blood flow velocity using joint spectral and time-domain optical coherence tomography,” Opt. Express 17, 10584–10598 (2009) [CrossRef] [PubMed] .

20. I. Grulkowski, I. Gorczynska, M. Szkulmowski, D. Szlag, A. Szkulmowska, R. A. Leitgeb, A. Kowalczyk, and M. Wojtkowski, “Scanning protocols dedicated to smart velocity ranging in spectral OCT,” Opt. Express 17, 23736–54 (2009) [CrossRef] .

21. J. Walther and E. Koch, “Enhanced joint spectral and time domain optical coherence tomography for quantitative flow velocity measurement,” in “Optical Coherence Tomography and Coherence Techniques V,”, R. A. Leitgeb and B. E. Bouma, eds. (Proc. of SPIE, 2011), p. 80910L.

22. Y. Zhao, Z. Chen, C. Saxer, Q. Shen, S. Xiang, J. de Boer, and J. Nelson, “Doppler standard deviation imaging for clinical monitoring of in vivo human skin blood flow,” Opt. Lett. 25, 1358–1360 (2000) [CrossRef] .

23. H. Ren, K. M. Brecke, Z. Ding, Y. Zhao, J. S. Nelson, and Z. Chen, “Imaging and quantifying transverse flow velocity with the Doppler bandwidth in a phase-resolved functional optical coherence tomography,” Opt. Lett. 27, 409–11 (2002) [CrossRef] .

24. D. Piao, L. L. Otis, and Q. Zhu, “Doppler angle and flow velocity mapping by combined Doppler shift and Doppler bandwidth measurements in optical Doppler tomography,” Opt. Lett. 28, 1120–1122 (2003) [CrossRef] [PubMed] .

25. S. G. Proskurin, Y. He, and R. K. Wang, “Determination of flow velocity vector based on Doppler shift and spectrum broadening with optical coherence tomography,” Opt. Lett. 28, 1227–1229 (2003) [CrossRef] [PubMed] .

26. L. Yu and Z. Chen, “Doppler variance imaging for three-dimensional retina and choroid angiography,” J. Biomed. Opt. 15, 016029 (2011) [CrossRef] .

27. V. J. Srinivasan, H. Radhakrishnan, E. H. Lo, E. T. Mandeville, J. Y. Jiang, S. Barry, and A. E. Cable, “OCT methods for capillary velocimetry,” Biomed. Opt. Express 3, 612–29 (2012) [CrossRef] [PubMed] .

28. R. V. Edwards, “Spectral analysis of the signal from the laser Ddoppler flowmeter: Time-independent systems,” J. Appl. Phys. 42, 837 (1971) [CrossRef] .

29. B. Park, M. Pierce, B. Cense, and S. Yun, “Real-time fiber-based multi-functional spectral-domain optical coherence tomography at 1.3 m,” Opt. Express 13(2005) [CrossRef] .

30. R. A. Leitgeb, M. Villiger, A. H. Bachmann, L. Steinmann, and T. Lasser, “Extended focus depth for Fourier domain optical coherence microscopy,” Opt. Lett. 31, 2450–2452 (2006) [CrossRef] [PubMed] .

31. H. Gross, Handbook of Optical Systems, vol. 1 (Wiley-VCH, Weinheim, 2005) [CrossRef] .

32. M. Gu, Advanced Optical Imaging Theory (Springer, Berlin, 1999).

33. J. W. Goodman, Introduction to Fourier Optics (Roberts & Company Publishers, Greenwood Village, 2005), 3rd ed.

34. M. Born and E. Wolf, Principles of Optics : Electromagnetic Theory of Propagation, Interference and Diffraction of Light (Cambridge University Press, Cambridge, 2009), 7th ed.

35. M. Villiger and T. Lasser, “Image formation and tomogram reconstruction in optical coherence microscopy,” J. Opt. Soc. Am. A 27, 2216–2228 (2010) [CrossRef] .

36. C. J. R. Sheppard, M. Roy, and M. D. Sharma, “Image formation in low-coherence and confocal interference microscopes,” Appl. Opt. 43, 1493–1502 (2004) [CrossRef] [PubMed] .

37. T. S. Ralston, D. L. Marks, P. Scott Carney, and S. A. Boppart, “Interferometric synthetic aperture microscopy,” Nature Phys. 3, 129–134 (2007) [CrossRef] .

38. S. G. Adie, B. W. Graf, A. Ahmad, P. S. Carney, and S. A. Boppart, “Computational adaptive optics for broadband optical interferometric tomography of biological tissue,” PNAS 109, 7175–7180 (2012) [CrossRef] [PubMed] .

39. C. J. R. Sheppard, S. Kou, and C. Depeursinge, “Reconstruction in interferometric synthetic aperture microscopy: comparison with optical coherence tomography and digital holographic microscopy,” J. Opt. Soc. Am. A 29, 244–250 (2012) [CrossRef] .

40. M. Szkulmowski, I. Gorczynska, D. Szlag, M. Sylwestrzak, and M. Wojtkowski, “Efficient reduction of speckle noise in Optical Coherence Tomography,” Opt. Express 20, 15149–15169 (2012) [CrossRef] .

41. A. Chan, E. Lam, and V. J. Srinivasan, “Comparison of kasai autocorrelation and maximum likelihood estimators for Doppler optical coherence tomography,” IEEE T. Med. Imaging (2013) [CrossRef] .

42. Z. Ding, H. Ren, Y. Zhao, J. S. Nelson, and Z. Chen, “High-resolution optical coherence tomography over a large depth range with an axicon lens,” Opt. Lett. 27, 243–245 (2002) [CrossRef] .

43. T. Bolmont, A. Bouwens, C. Pache, M. Dimitrov, C. Berclaz, M. Villiger, B. M. Wegenast-Braun, T. Lasser, and P. C. Fraering, “Label-free imaging of cerebral beta-amyloidosis with extended-focus optical coherence microscopy,” J. Neurosci. 32, 14548–14556 (2012) [CrossRef] [PubMed] .

44. J. Lee, W. Wu, J. Y. Jiang, B. Zhu, and D. A. Boas, “Dynamic light scattering optical coherence tomography,” Opt. Express 20, 22262–77 (2012) [CrossRef] [PubMed] .

45. M. Villiger, C. Pache, and T. Lasser, “Dark-Field optical coherence microscopy,” Opt. Lett. 35, 3489–3491 (2010) [CrossRef] [PubMed] .

Approximations	Equations
Scalar fields, point scatterers	Eq. (27), Eq. (29)
+ moderate NA (≤∼ 0.5)	Eq. (22), Eq. (29)
+ Gaussian spectrum S(k), $k_{σ}^{2} ≪ k_{0}^{2}$	Eq. (22), Eq. (30)
+ low NA (≤∼ 0.05), equal Gaussian modes	Eq. (22), Eq. (31)

Quantitative lateral and axial flow imaging with optical coherence microscopy and tomography

Abstract

1. Introduction

2. Theory

2.1. Fourier domain coherent imaging for static scattering

2.2. Fourier domain coherent imaging for dynamic scattering

2.3. Doppler spectrum per wavenumber channel: mean and standard deviation

2.4. Spatially resolved Doppler spectrum

2.5. Approximations for Gaussian modes and low NA

2.6. Doppler spectrum: numerical simulation

3. Methods

3.1. xfOCM

3.2. OCT

3.3. Flow system

3.4. Data acquisition & algorithm

4. Experimental results

4.1. xfOCM

4.2. OCT

5. Discussion & conclusion

Supporting information

Acknowledgments

References and links

Cited By

Figures (9)

Tables (1)

Equations (32)

Optics Express