## Abstract

The slope of the inverse square of the contrast values versus camera exposure time at multi-exposure speckle imaging (MESI) can be a new indicator of flow velocity. The slope is linear as the diffuse coefficient in Brownian motion diffusion model and the mean velocity in ballistic motion model. Combining diffuse speckle contrast analysis (DSCA) and MESI, we demonstrate theoretically and experimentally that the flow velocity can be obtained from this slope. The calculation results processes of the slop don’t need tedious Newtonian iterative method and are computationally inexpensive. The new indicator can play an important role in quantitatively assessing tissue blood flow velocity.

© 2014 Optical Society of America

## 1. Introduction

Laser speckle contrast analysis (LASCA) [1, 2], firstly introduced in the 1980s, is a particularly attractive method because of its simple and inexpensive manner. It has been widely used to image biological blood flow [3–5]. The very essence [6] of LASCA is to quantify the visibility of the speckle pattern in terms of moments of the distribution of recorded intensities for a given exposure duration. The visibility of the speckle pattern is defined as the ratio of the standard deviation to the mean intensity,${K}_{s}={\sigma}_{s}/\u3008I\u3009$, where the subscript *s* refers to the spatial or to temporal [3, 7].

LASCA is mostly based on single scattering, which is limited on a superficial layer inside a relatively transparent medium. Diffuse speckle contrast analysis (DSCA) [8, 9] exploits the diffusive nature of the transport of light for the relatively turbid medium which exhibits strong multiple scattering. This new technique is formed by the union of diffusing wave spectroscopy (DWS) [10] and LASCA. The autocorrelation function of DSCA can be obtained from the correlation transport equation, which considers the effects of absorption, scattering and different dynamical motions and the influence of boundaries and heterogeneities in this physical model. Thus by means of this autocorrelation function, we can obtain quantitative determination of dynamical properties of the medium. DSCA has more advantages and can overcome some difficulties emerging from LASCA.

Traditionally, the correlation times (*τ _{c}*) is regarded as the indicator of mean velocity of the scatters in LASCA and DSCA. The correlation time is the decay time of the speckle autocorrelation function and is inversely proportional to the mean velocity of the scatters. If the ratio of the exposure time to the correlation time is much more than 1, the square of the contrast values (

*K*

^{2}) is linearly related to

*τ*[11]. The inverse square of the contrast values (1/

_{c}*K*

^{2}) is similarly regarded as the indicator of the flow index. Nevertheless, this relationship is speculative and there are no first principles arguments to verify its validity [2]. Meanwhile, the calculation processes of

*τ*need to use tedious Newtonian iterative method and to consider the reduction in the measured contrast due to noise, averaging over uncorrelated speckles [12] and polarization [13]. It is computationally expensive.

_{c}In this paper, we present a new indicator for the flow index, i.e. the slope of 1/*K*^{2} vs. camera exposure time (T) at multi-exposure speckle imaging (MESI) [14]. It is theoretically derived that the slope of the function of 1/*K*^{2} vs. T (*k _{slope}*) changing over exposure time depends on source-detector distance, flow information and optical parameters of medium (the absorption coefficient

*μ*, the scattering coefficient

_{a}*μ*and anisotropy parameter

_{s}*g*et al.). When the ratio of the minimum exposure time to the correlation time is much more than 1,

*k*is in linear relationship with the diffuse coefficient (

_{slope}*D*) in Brownian motion diffusion model and the mean velocity (

_{B}*v*) in ballistic model (random flow and shear flow [15]). We perform a tissue phantom experiment in multiple scattering slab medium and Monte Carlo simulation to demonstrate this linear relationship between this new indicator and the mean velocity. The experiment results and MC simulation results fit well.

*k*can become a quantitative indicator for tissue blood flow velocity. The new indicator has more advantages and can play an important role in quantitatively assessing tissue blood flow velocity.

_{slope}## 2. Theory

#### 2.1 Diffuse speckle contrast analysis

In homogeneous multiple scattering slab medium, the electric field autocorrelation function ${G}_{1}\left(\overrightarrow{r},\tau \right)=\u3008E\left(\overrightarrow{r},t\right){E}^{*}\left(\overrightarrow{r},t+\tau \right)\u3009$ can be obtained from the correlation diffusion equation [16]

*k*

_{0}is the wavenumber of the light in medium,

*S*(

*r*) is the light-source distribution and $\u3008\Delta {r}^{2}\left(\tau \right)\u3009$ represents the mean square displacement of the moving scatterers after a delay time

*τ*.

When the slab thickness [17] is $d/{l}_{tr}\ge 3$, the unnormalized electric autocorrelation function has the following expression [18] by standard diffusion approximation,

*d*is the thickness of slab medium,

*l*is the transport mean free path,

_{tr}*R*is effective reflection coefficient to account for the index mismatch between media and surrounding medium,

_{eff}*n*=

*n*/

_{in}*n*is the ratio of the index of refraction inside and outside the diffusing medium,$\alpha $ is the fraction of dynamic photon scattering events in medium,

_{out}*m*is an integer,

*r*is source-detector distance in x-y plane and

*z*is along the thickness direction of slab medium. A scheme of the geometry together with a description of some of the notations is shown in Fig. 1.

The normalized *g*_{1}(*r, z, τ*) is given by

*τ*as correlation time, when the function

_{c}*g*

_{1}(

*r, z, τ*) = 1/

*e*is established. For the strong multiple scattering ${\mu}_{s}^{\text{'}}\gg {\mu}_{a}$, $3{\mu}_{s}^{\text{'}}{\mu}_{a}$ can be omitted in Eq. (2). The speckle contrast is related to the normalized electric autocorrelation function, given by

*T*is exposure time and

*β*is a constant accounting for the reduction in the measured contrast due to averaging over uncorrelated speckles and polarization [13].

Substituting Eqs. (2)–(10) into Eq. (11), the diffuse speckle contrast can be written as

*ξ*

_{±}

_{m}_{, ±}

*is constant at the fixed source-detector distance and depends on*

_{x}*r*

_{±}

*,*

_{m}*r*

_{±}

*and*

_{x}*G*

_{1}(

*r, z*, 0).

For Brownian motion $\u3008\Delta {r}^{2}\left(\tau \right)\u3009=6{D}_{B}\tau $, *${g}_{1}\left(r,z,\tau \right)$*is the sum of the similar function form ${g}_{1}\left(\tau \right)=\mathrm{exp}\left(-\sqrt{\Gamma \tau}\right)$, correspond to contrast function form [6]

For the ballistic motion model $\u3008\Delta {r}^{2}\left(\tau \right)\u3009={v}^{2}{\tau}^{2}$, *${g}_{1}\left(r,z,\tau \right)$* is the sum of the similar function form ${g}_{1}\left(\tau \right)=\mathrm{exp}\left(-\Gamma \tau \right)$, correspond to contrast function form [6]

When the ratio of the exposure time to the correlation time is much more than 1, i.e. ${\Gamma}_{B}T\gg 1$ and ${\Gamma}_{b}T\gg 1$, the numerator function $\left(3+6\sqrt{{\Gamma}_{B}T}+4{\Gamma}_{B}T\right){e}^{-2\sqrt{{\Gamma}_{B}T}}$ and ${e}^{-2\sqrt{{\Gamma}_{b}T}}$ can be omitted in Eq. (13) and Eq. (15). *V*(*T*) of Brownian motion and ballistic motion can be simplified to the same following expression

Exploiting the simplified Eq. (17) and Eq. (18), the inverse square of the contrast values 1/*K*^{2}(*T*) is given by

Because of ${\Gamma}_{1}/T\ll 1$ and ${\Gamma}_{3}/T\ll 1$, Eq. (19) and Eq. (22) can be simplified to the following expression

From the above results, we can deduce that the function of 1/*K*^{2}(*T*) has a linear relationship with $T$ for Brownian motion and ballistic motion over a broad range. The slope of the function of 1/*K*^{2}(*T*) vs T (*k _{slope}*) depends on the source-detector distance, flow information and optical parameters of medium. When the source-detector distance is constant,

*k*has a linear relationship with diffuse coefficient

_{slope}*D*in Brownian motion and mean velocity

_{B}*v*in ballistic motion.

*k*can be used as the new indicator for quantitatively assessing flow index.

_{slope}#### 2.2 Monte Carlo simulation

For checking the accuracy of the theory, we ran simulation for a point source vertically entering the slab medium with different dynamical properties. For our simulation, one hundred million is applied and the code takes ~2 h to run on an Intel i5 3450 3.1 GHz processer. The well-developed MC codes [19, 20] have been widely used to simulate light transport in tissues. We modify these codes to be able to obtain the momentum transfer of each photon-scatter collision and finally success to record the numerical distribution of total dimensionless momentum transfer. The normalized autocorrelation can be calculated as

*P*(

*Y*,

*r*,

*z*) is normalized probability density function of total dimensionless momentum transfer at different distances and $Y={\displaystyle {\sum}_{i=1}^{n}\left(1-\mathrm{cos}{\theta}_{i}\right)}$ is the total dimensionless momentum transfer experienced

*n*uncorrelated scattering events. The contrast can be calculated by Eq. (11). The whole calculation process avoids the standard diffusion approximation.

Figure 2(a) shows MC simulation results and the analytical
results Eq. (26) when the mean velocity is 0.80
mm/s. We try much larger m values to calculate Eq.
(26). When m tends to be infinite, the inverse square of the contrast values change
very slowly and are almost constant. The slope difference between m = 0, ± 1,…,
± 30 and m = 0, ± 1,…, ± 300 is 0.3%. However, a systematic
difference remains between the analytical expressions and MC simulation results. This small
discrepancy is caused by the photons which don’t experience enough scattering events and
are non-diffuse after passing through the slab. As shown in Fig. 2(b) the single scattering photons and the ballistic-diffusive transition photons
account for a certain proportion, which make
1/*K*^{2}(*T*) become smaller than the analytical values
calculated by diffusion approximation. However, we can clearly find that the functions of
1/*K*^{2}(*T*) both have a linear relationship with
*T* for ballistic motion over a broad range.

## 3. Experimental set-up

To validate the role of the new indicator, a flow phantom was built as shown by Fig. 3.The continuous He-Ne laser source (Melles Griot 25-LHP-925-230; *λ*
= 632.8*nm*, 30*mW*, Carlsbad California, USA) with a long
coherence length is vertically entering the rectangular vessel (2*cm* height and
*d* = 2*cm* wide in cross section). The intralipid solution with
a concentration of 0.15% is pumped through the polyethylene tubing (4.8*mm* inner
diameter), at different mean flow velocities and different exposure times. The optical
properties of the flow phantom are *μ ^{’}_{s}* =
2

*cm*

^{−1}and

*μ*= 0.01

_{a}*cm*

^{−1}, which satisfy

*μ*>>

^{’}_{s}*μ*, and the diffusion approximation is valid.

_{a}The mean flow velocities are 0.80 *mm*/*s* to 2.40 *mm*/*s*, drived by a digital peristaltic pump. The exposure time is set at 0.3 *ms*, 0.5 *ms*, 0.8 *ms*, 1 *ms*, 2 *ms* and 3 *ms*. The transmission speckle is imaged using a prime lens, 120*mm* lens tube and a CCD (Microvision MV-VS141FM/C; 12 bits, pixel size 4.65 *μm**4.65 *μm*, resolution ratio 1392*1040). The imaging system has a magnification of 1.7. The detector position is set at the position *z* = 2.0 *cm r* = 0.2 *cm*.

## 4. Results and discussion

20 raw speckle images are recorded for each mean velocity and each exposure time. Because of the influence of the CCD background noise, the raw speckle images need to subtract the average noise intensity, which is averaged over 50 noise images recorded for each exposure time, before the calculation of the spatial contrast value *K _{s}*.

*K*is calculated from every segment of 101 × 101 pixels for reducing the noise in laser speckle correlation [12]. 20 speckle contrast images

_{s}*K*are averaged for each exposure time and each mean velocity.

_{s}Figure 4 shows the spatial contrast of different mean
velocities at different exposure times. The mean velocities are 0.80 *mm*/*s*, 1.60
*mm*/*s* and 2.40 *mm*/*s*
corresponding to the black, the red and the blue curves respectively. If we expect to obtain the
linear relationship in correlation time *τ _{c}*, we need
nonlinearly fit the contrast by MESI equation. However,

*β*is one of the unknown quantities and mainly depends on the experiment conditions. Usually we obtain

*β*by calculating the contrast of the reflectance static speckle. Holding

*β*constant, we make the nonlinear fitting to obtain

*τ*. This calculation processes need tedious Newtonian iterative method and are computationally inexpensive. Meanwhile, to improve the signal to noise ratio (SNR), the number of exposure time must be adequate. This process is cumbersome.

_{c}The new indicator *k _{slope}* can overcome these shortcomings as motioned
above. As shown in Fig. 5, we can clearly see that the
experiment results and the analytical results have the same tendency.Both 1/

*K*

_{s}^{2}(

*T*) have a linear relationship with

*T*at different mean velocity. The difference mainly caused by

*β*is at the slope as shown in Eqs. (25) and (26). The unpolarized recorded speckle images and the reduction in the measured contrast due to averaging over uncorrelated speckles make

*β*< 1. The

*k*of the experiment results is equivalent to enlarge 1/

_{slope}*β*times. Nevertheless, we don’t need consider the influence of

*β*and don’t calculate its value like the calculation process in MESI, because the effect of

*β*is same for different velocities. Meanwhile, the linear fitting is relatively simple and doesn’t need too much exposure times for achieving the same SNR. The new indicator has more advantages than the correlation time

*τ*.

_{c}In order to verify our argument on the analytical expressions Eq. (26), we calculate experiment results *k _{slope}*
at different mean velocities as shown in Table 1and
Fig. 6(a).We can clearly find that

*k*varies linearly with the scattering particle flow velocity. For eliminating the influence of

_{slope}*β*and plotting on the same scale in an image, the normalization is performed with respect to

*k*of the experiment results and the analytical results at their respective median values as shown in Table 1. The normalization process is performed for the different mean velocities in the same way to observe easily this linear relationship between the slope and the mean velocity.

_{slope}Figure 6(b) shows that the normalized *k _{slope}* of experiment results and the analytical results fit well. The median mean velocity can be regarded as the baseline measure value, so we can utilize this linear relationship to obtain other quantitative velocity information.

*k*can be a new indicator to quantitatively assesse flow velocity.

_{slope}## 5. Conclusion

In this work, we have presented a new indicator *k _{slope}* to measure quantitatively mean velocity. The calculation theory is based on diffuse speckle contrast analysis (DSCA) and multi-exposure speckle imaging (MESI). We demonstrate that we are able to use the linear relation to fit the inverse square of the contrast values (1/

*K*

^{2}) at different exposure times.

*k*is in linear relationship with the mean velocity. We provide MC simulation results and experiment results in a liquid phantom with varying flow velocity in the transmission geometry. The results are in support of our theory.

_{slope}We can utilize this linear relationship and the baseline measure value *k _{slope}* of known mean velocity to obtain other quantitative velocity information. Meanwhile the calculation process of

*k*doesn’t need consider the influence of the unpolarized recorded speckle images and the reduction in the measured contrast due to averaging over uncorrelated speckles. The normalization can be performed to overcome these influence as motioned above. The computation can be greatly simplified for multi-exposure speckle imaging. The new indicator can play an important role in quantitatively assessing tissue blood flow velocity.

_{slope}## Acknowledgment

The work is supported by the National Natural Science Foundation of China under Grant No 11174151.

## References and links

**1. **A. Fercher and J. Briers, “Flow visualizaiton by means of single-exposure speckle photography,” Opt. Commun. **37**(5), 326–330 (1981). [CrossRef]

**2. **J. D. Briers and S. Webster, “Laser speckle contrast analysis (LASCA): a nonscanning, full-field technique for monitoring capillary blood flow,” J. Biomed. Opt. **1**(2), 174–179 (1996). [CrossRef] [PubMed]

**3. **P. Li, S. Ni, L. Zhang, S. Zeng, and Q. Luo, “Imaging cerebral blood flow through the intact rat skull with temporal laser speckle imaging,” Opt. Lett. **31**(12), 1824–1826 (2006). [CrossRef] [PubMed]

**4. **K. Murari, N. Li, A. Rege, X. Jia, A. All, and N. Thakor, “Contrast-enhanced imaging of cerebral vasculature with laser speckle,” Appl. Opt. **46**(22), 5340–5346 (2007). [CrossRef] [PubMed]

**5. **Z. Wang, S. Hughes, S. Dayasundara, and R. S. Menon, “Theoretical and experimental optimization of laser speckle contrast imaging for high specificity to brain microcirculation,” J. Cereb. Blood Flow Metab. **27**(2), 258–269 (2007). [CrossRef] [PubMed]

**6. **R. Bandyopadhyay, A. Gittings, S. Suh, P. Dixon, and D. Durian, “Speckle-visibility spectroscopy: A tool to study time-varying dynamics,” Rev. Sci. Instrum. **76**(9), 093110 (2005). [CrossRef]

**7. **H. Cheng, Y. Yan, and T. Q. Duong, “Temporal statistical analysis of laser speckle images and its application to retinal blood-flow imaging,” Opt. Express **16**(14), 10214–10219 (2008). [CrossRef] [PubMed]

**8. **R. Bi, J. Dong, and K. Lee, “Deep tissue flowmetry based on diffuse speckle contrast analysis,” Opt. Lett. **38**(9), 1401–1403 (2013). [CrossRef] [PubMed]

**9. **R. Bi, J. Dong, and K. Lee, “Multi-channel deep tissue flowmetry based on temporal diffuse speckle contrast analysis,” Opt. Express **21**(19), 22854–22861 (2013). [CrossRef] [PubMed]

**10. **D. J. Pine, D. A. Weitz, P. M. Chaikin, and E. Herbolzheimer, “Diffusing wave spectroscopy,” Phys. Rev. Lett. **60**(12), 1134–1137 (1988). [CrossRef] [PubMed]

**11. **H. Cheng and T. Q. Duong, “Simplified laser-speckle-imaging analysis method and its application to retinal blood flow imaging,” Opt. Lett. **32**(15), 2188–2190 (2007). [CrossRef] [PubMed]

**12. **S. E. Skipetrov, J. Peuser, R. Cerbino, P. Zakharov, B. Weber, and F. Scheffold, “Noise in laser speckle correlation and imaging techniques,” Opt. Express **18**(14), 14519–14534 (2010). [CrossRef] [PubMed]

**13. **P.-A. Lemieus and D. J. Durian, “Investigating non-Gaussian scattering processes by using nth-order intensity correlation functions,” J. Opt. Soc. Am. A **16**(7), 1651–1664 (1999). [CrossRef]

**14. **A. B. Parthasarathy, W. J. Tom, A. Gopal, X. Zhang, and A. K. Dunn, “Robust flow measurement with multi-exposure speckle imaging,” Opt. Express **16**(3), 1975–1989 (2008). [CrossRef] [PubMed]

**15. **X. L. Wu, D. J. Pine, P. M. Chaikin, J. S. Huang, and D. A. Weitz, “Diffusing-wave spectroscopy in a shear flow,” J. Opt. Soc. Am. B **7**(1), 15–20 (1990). [CrossRef]

**16. **D. Boas and A. Yodh, “Spatially varying dynamical properties of turbid media probed with diffusing temporal light correlation,” J. Opt. Soc. Am. A **14**(1), 192–215 (1997). [CrossRef]

**17. **Z. Q. Zhang, I. P. Jones, H. P. Schriemer, J. H. Page, D. A. Weitz, and P. Sheng, “Wave transport in random media: The ballistic to diffusive transition,” Phys. Rev. E Stat. Phys. Plasmas Fluids Relat. Interdiscip. Topics **60**(4), 4843–4850 (1999). [CrossRef] [PubMed]

**18. **T. Durduran, R. Choe, W. B. Baker, and A. G. Yodh, “Diffuse optics for tissue monitoring and tomography,” Rep. Prog. Phys. **73**(7), 076701 (2010). [CrossRef]

**19. **L. Wang, S. L. Jacques, and L. Zheng, “MCML--Monte Carlo modeling of light transport in multi-layered tissues,” Comput. Methods Programs Biomed. **47**(2), 131–146 (1995). [CrossRef] [PubMed]

**20. **L. Wang and S. L. Jacques, http://omlc.ogi.edu/software/mc/.