Expand this Topic clickable element to expand a topic
Skip to content
Optica Publishing Group

Model of dynamic speckle evolution for evaluating laser speckle contrast measurements of tissue dynamics

Open Access Open Access

Abstract

We introduce a dynamic speckle model (DSM) to simulate the temporal evolution of fully developed speckle patterns arising from the interference of scattered light reemitted from dynamic tissue. Using this numerical tool, the performance of laser speckle contrast imaging (LSCI) or speckle contrast optical spectroscopy (SCOS) systems which quantify tissue dynamics using the spatial contrast of the speckle patterns with a certain camera exposure time is evaluated. We have investigated noise sources arising from the fundamental speckle statistics due to the finite sampling of the speckle patterns as well as those induced by experimental measurement conditions including shot noise, camera dark and read noise, and calibrated the parameters of an analytical noise model initially developed in the fundamental or shot noise regime that quantifies the performance of SCOS systems using the number of independent observables (NIO). Our analysis is particularly focused on the low photon flux regime relevant for human brain measurements, where the impact of shot noise and camera read noise can become significant. Our numerical model is also validated experimentally using a novel fiber based SCOS (fb-SCOS) system for a dynamic sample. We have found that the signal-to-noise ratio (SNR) of fb-SCOS measurements plateaus at a camera exposure time, which marks the regime where shot and fundamental noise dominates over camera read noise. For a fixed total measurement time, there exists an optimized camera exposure time if temporal averaging is utilized to improve SNR. For a certain camera exposure time, photon flux value, and camera noise properties, there exists an optimized speckle-to-pixel size ratio (s/p) at which SNR is maximized. Our work provides the design principles for any LSCI or SCOS systems given the detected photon flux and properties of the instruments, which will guide the experimental development of a high-quality, low-cost fb-SCOS system that monitors human brain blood flow and functions.

© 2022 Optica Publishing Group under the terms of the Optica Open Access Publishing Agreement

1. Introduction

Diffuse correlation spectroscopy (DCS) [13] and laser speckle contrast imaging (LSCI) [48] are two non-invasive optical techniques for measuring cerebral blood flow (CBF) based on speckle fluctuations of coherent light. DCS quantifies a blood flow index by measuring the temporal intensity autocorrelation function of the speckle intensity fluctuations induced by tissue dynamics. This technique is routinely used for bedside monitoring of CBF changes in humans but requires expensive single photon counting devices, which is not feasible for some measurement geometries including high-density approaches. LSCI uses relatively low-cost complementary metal–oxide–semiconductor (CMOS) cameras to measure blood flow changes based on the spatial speckle contrast reduction arising from temporal integration of the speckle fluctuations over the camera exposure time [9,10]. Traditional LSCI is utilized to obtain two dimensional images of superficial cerebral blood flow of the tissue with wide-field illumination using laser light, mostly for small animals such as mice with thinned skulls or cranial window preparations, which is not feasible for human brain measurements. Recently, a new technique named speckle contrast optical spectroscopy (SCOS) and its tomographic extension speckle contrast optical tomography (SCOT) have demonstrated speckle contrast measurements with free-space imaging of larger source-detector separations, enabling non-invasive measurements of blood flow in deeper tissues in phantoms, human forearm and small animal brains [1113]. But generalizing free space techniques to human brain function measurements in a large area is challenging due to hair, the limited area of the focus, and the limited field of view of the camera [14]. Fiber based systems have the potential to expand upon the capabilities of free-space SCOT for non-invasive monitoring of CBF for large brain areas. Predictions of the performance of a fiber based SCOS (fb-SCOS) system with respect to the measurement parameters will guide the experimental development of low-cost, high-density, wearable systems for everyday monitoring of CBF. Thus, a proper noise model for LSCI or SCOS is desired. Recently, Xu et al. have established an analytical model [15] that calculates the signal to noise ratio (SNR) defined as $K^2/\sigma (K^2)$ described by the number of independent observables (NIO), taking into account the effects of fundamental noise due to finite speckle sampling and shot noise, where $K^2$ is the spatial speckle contrast squared calculated over window size squared $w^2$, and the noise or the standard deviation of $K^2$ is denoted as $\sigma (K^2)$. But this model does not specify how NIO is related to experimentally controllable parameters including the speckle to pixel size ratio $s/p$, camera exposure time $T_{exp}$, camera frame rate, photon flux per speckle $c_s$, and the number of camera pixels etc.. Instrumental noise such as read noise is also not included in the analytical model.

In this paper, we present a complete noise model for LSCI or SCOS measurements using the dynamic speckle model (DSM) that simulates the temporal evolution of speckle patterns, to which shot noise, camera dark and read noise can be added. This numerical tool is an extension of a previous model that generates static speckle patterns to include tissue dynamics and the temporal evolution of speckle patterns [16,17]. We have validated our modeling results with the existing analytical expression [15], and established the link between NIO and the experimentally controllable parameters for arbitrary noise levels, and included camera read noise into the analytical model that was initially established in the fundamental and shot noise regime. We have used the numerical model as well as the calibrated analytical model to identify the impact of $T_{exp}$, $s/p$, and $c_s$ on the performance of SCOS, particularly for two commercially available cameras, to guide the future development of fb-SCOS systems for human brain blood flow and functional measurements.

2. Methods

2.1 Analytical model

For LSCI or SCOS, the speckle contrast, $K$, is calculated as the standard deviation of the speckle intensity measured across the pixels divided by the mean of the intensity values, $K=\sigma (I)/\langle I \rangle$. When there is no contribution from noise sources, the fundamental contrast $K_f$ and the field autocorrelation function $g_1(\tau )=\langle E(t)E^*(t+\tau )\rangle / \langle |E(t)|^2\rangle$ typically used in DCS are related via [18]

$$K_{f}^2 =\frac{2\beta}{T_{exp}}\int_0^{T_{exp}} g_1^2(\tau)\Big(1-\frac{\tau}{T_{exp}}\Big)d\tau,$$
when $g_1(\tau )$ follows the analytical expression of $g_1(\tau )=exp(-\tau /\tau _c)$ as for the case with typical source detector separation in most human brain measurements for deep blood flow measurements when fitting for the early delay time $\tau$ [19,20], where $\tau _c$ is the decay constant, the fundamental speckle contrast is expressed in terms of $\tau _c$ and $T_{exp}$ as [7,8]
$$K_{f}^2 = \beta \frac{\tau_c}{T_{exp}} \Big[1 + \frac{\tau_c}{2T_{exp}} \Big(exp\Big(\frac{-2T_{exp}}{\tau_c}\Big)-1\Big)\Big] .$$
where $\beta$ takes into account the loss of coherence due to a finite coherence length of the light source and the detection of multiple modes such as the two polarization states and multiple speckles detected by a single detector or pixel. We have found numerically that the coherence parameter has a dependence on $s/p$ given by $\beta = \beta _l(s/p)^2 / (1 + (s/p)^2 )$ (Supplementary Material Fig. S1) where $\beta _l$ is the value of $\beta$ when $s/p\gg 1$, which can be less than 1 due to the effect of the laser coherence length and the polarization states, for instance. Here, $s/p$ is the ratio of the speckle size $s$ given by the full-width half maximum (FWHM) of the spatial autocorrelation map to the pixel size ($p$) of the detector that images the speckle patterns. The values of $K_f^2$ as a function of $s/p$ is consistent with the results obtained from another numerical tool [21], but that tool is not capable of modelling the full temporal statistics of the dynamic speckle pattern [22] with experimentally controllable parameters.

The noise sources in LSCI or SCOS measurements include the fundamental noise arising from the finite sampling of the speckle statistics, and instrument noise including shot, dark, and read noise. Note that in the literature "dark noise" often includes read noise but we have separated the dark noise and read noise contributions in this manuscript for quantitative analysis of each parameter. Both shot noise and dark noise increase with exposure time $T_{exp}$ while read noise stays constant. Since the sum of the variances arising from all the noise sources add up, the measured speckle contrast is given by

$$K_{all}^2 = K_f^2 + K_s^2 + K_d^2 + K_r^2$$
Here $K_{all}$ and $K_f$ represent the speckle contrast that takes into account all noise sources and the fundamental contrast respectively, while $K_s$, $K_d$ and $K_r$ represent the bias terms in the contrast measurements that arise from shot, dark and read noise respectively. These bias terms cannot be ignored at low photon flux values typically observed for human brain measurements, and needs to be corrected for to obtain $K_f^2$ which is related to the blood flow as described by Eqs. (12). Note that for each of these $K$’s the mean intensity is given by the fundamental signal, i.e. we assume that the mean intensity signal offset from the dark and read noise has been corrected for. The bias term arising from shot noise depends on the average photon flux per pixel $c_p$ and $T_{exp}$ [15]:
$$K_{s}^2 = \frac{1} {Q_e c_p T_{exp}} .$$
Here $c_p$ is related to the average photon flux per speckle $c_s$ as $c_p = c_s / (s/p)^2$. $Q_e$ is the quantum efficiency for converting photons to electrons. We use $c_s$ of $10,000$ photons/s and $100,000$ photons/s in this manuscript, which are within the range of typical photon flux values in human measurements for a source-detector separation in the range of 10-30 mm [23,24].

Here we neglect the dark noise contribution $K_d=0$, which is generally small for the exposure time range we are considering and can be incorporated into shot noise by adding the dark electrons into $Q_ec_p$ if needed. The read noise $\sigma _r$ contribution to the bias is

$$K_{r}^2 = \frac{\sigma_r^2}{(Q_ec_p T_{exp})^2} .$$
The noise in LSCI or SCOS measurements is defined as the standard deviation of $K^2$, i.e. $\sigma (K^2)$. Recently, Xu et al. [15] analytically derived the relation between $K^2$ and $\sigma (K^2)$ in terms of the number of independent observables ($NIO$) when the measurement noise arises only from shot-noise and fundamental speckle sampling statistics, and found $\sigma (K_{all}^2) = K_{all}^2\sqrt {2/NIO}$ in the limit of large $T_{exp}$, $c_p$, and $NIO$, i.e., $T_{exp}\gg \tau _c$, $c_p T_{exp} \gg 1$ to ensure the system is in the shot-noise limit, and $NIO\gg 1$ [15]. However, a link between $NIO$ and experimental parameters including $s/p$ and the total number of measured speckles is desired. Also, when shot noise is present, the definition of $NIO$ is further complicated since shot noise contribution does not depend on the speckle intensity correlation between adjacent pixels explicitly while fundamental noise does. Thus, we have utilized the method that Xu et al. has developed to obtain the fundamental noise [15] for all $T_{exp}$, and included the shot and read noise contribution to $\sigma (K_{all}^2)$ as
$$\sigma(K^2_f)=K_f^2\sqrt{\frac{2(1+2K_f^2)} {NIO}}, $$
$$\qquad\sigma(K^2_s)=c_{K_s}K^2_s\sqrt{\frac{1}{w^2}}, $$
$$\qquad\sigma(K^2_r)=c_{K_r}K^2_r\sqrt{\frac{1}{w^2}}, $$
where $\sigma (K^2_f)$, $\sigma (K^2_s)$, $\sigma (K^2_r)$ are the noise in $K^2_f$, $K^2_s$, $K^2_r$ respectively, and $w^2$ is the window size squared or the number of pixels used for calculating the speckle contrast. Here $c_{K_s}$ and $c_{K_r}$ are constants which we have calibrated to be $c_{K_s}\approx 1.90$, $c_{K_r}\approx 1.47$ (Supplementary Material Fig. S2) for the measurement conditions we have explored. The total noise $\sigma (K^2_{all})$ can be calculated as
$$\sigma(K^2_{all})=\sqrt{\sigma^2(K^2_f)+\sigma^2(K^2_s)+\sigma^2(K^2_r)}.$$
For the case of fb-SCOS, when a measurement is made through a fiber with a fixed number of speckles $m$, we define $NIO$ as
$$NIO = m(s/p)^2 f(s/p) = w^2 f(s/p),$$
where $f(s/p)$ accounts for the reduction in $NIO$ due to the spatial correlation of measurements across pixels and $w^2=m(s/p)^2$ is the number of pixels used to image the intensity of light from the fiber. We have fixed the number of fiber modes in this study since modern CMOS cameras have enough pixels to measure all the modes from a fiber. The effect of the number of pixels can also be studied using our model for particular systems where the number of pixels is limited. If $f(s/p)=1$, $NIO=m(s/p)^2=w^2$. However, due to the intensity correlation between pixels that arise from the correlation between speckles, $f(s/p)\le 1$, and thus $NIO\le w^2$ as will be shown in Section 2.2. When $s/p\gg 1$, the speckle intensities are sufficiently sampled and $NIO$ approaches the total number of independent fiber modes which is proportional to the number of speckle modes $m$ obtained as the total area of the speckle image divided by the area of a single speckle $s^2$. Here the shape of the speckles is considered rectangular instead of circular for simplicity, and any effect of the difference in the speckle area due to the speckle shapes is included in the constant $c_{s/p}$ below. When $s/p\ll 1$, the correlation between neighboring pixels is minimized thus $NIO$ approaches $w^2$. With these asymptotic behaviors, we have constructed an empirical model for $f(s/p)$ as
$$f(s/p)=1/(1+c_{s/p}(s/p)^2),$$
where the constant $c_{s/p}$ is calibrated to be $0.292$ from our numerical model as described in Section 2.2. Here $c_{s/p}$ is the ratio of the number of mathematically well-defined orthogonal fiber modes to the number of speckles that are correlated with long tails overlapping with each other resulting in their boundaries not well defined. We then define the signal to noise ratio (SNR) as
$$SNR = \frac{K_f^2}{\sigma(K_{all}^2)}.$$
Here, we are using $K_f^2$ instead of $K_{all}^2$ in the numerator of the SNR calculation, since $K_f^2$ is the signal that contains blood flow information while $K_{all}^2$ contains bias terms from the noise sources as shown in Eq. (3). Experimentally, $K_f^2$ needs to be estimated from the measurements by subtracting the noise bias terms using Eqs. (3)–5. Experimentally, temporal averaging over multiple frames is sometimes utilized to improve SNR. When temporal averaging of $N_{fr}$ camera frames is used, Eqs. (6a), (6b), (6c) is revised as
$$\sigma_{fa}(K^2_f)=K_f^2\sqrt{\frac{2(1+2K_f^2)} {NIO_{fa}}}, $$
$$\quad\sigma_{fa}(K^2_s)=c_{K_s} K^2_s\sqrt{\frac{1}{w^2N_{fr}}}, $$
$$\quad\sigma_{fa}(K^2_r)=c_{K_r} K^2_r\sqrt{\frac{1}{w^2N_{fr}}}, $$
where the subscript $fa$ denotes the corresponding variable obtained for frame averaging i.e. the average of multiple frames. Here $NIO_{fa}$ can be expressed as
$$NIO_{fa} = m(s/p)^2 f(s/p)N_{fr}h(T_{fr}/\tau_c)$$
where $T_{fr}$ is the inverse of the camera frame rate or the time between adjacent camera frames, and $h(T_{fr}/\tau _c)$ takes into account possible temporal correlations between adjacent frames. For most cameras and human measurements, $T_{fr}\gg \tau _c$ thus there is no correlation between speckle images of adjacent camera frames and $h(T_{fr}/\tau _c)=1$ (Supplementary Material Fig. S3). The total noise for frame averaging $\sigma _{fa}(K_{all}^2)$ can be simplified as
$$\sigma_{fa}(K_{all}^2)=\sigma(K_{all}^2)/\sqrt{N_{fr}}.$$

2.2 Dynamic speckle model

Here we describe our physics-based model to simulate the temporal evolution of a dynamic speckle pattern from first principles. The Matlab code of our numerical model is publicly available [25]. A schematic representation of the model is shown in Fig. 1(a). The system consists of $N_p$ particles acting as point sources that are moving in time ($N_p=5000$ for our simulations). The initial position of the $n^{th}$ particle ($x_{n,i}, y_{n,i}, z_{n,i}$) is randomly distributed within a volume of size [1000 $\mu m$ 1000 $\mu m$ 1000 $\mu m$] which represents the scattering medium. We use a standard Monte-Carlo method of Gaussian random walks to simulate the Brownian motion of the particles by updating the location of the $n^{th}$ particle at each time step as

$$\begin{aligned}&x_n(t+\Delta t) = x_n(t)+ Norm(0,2D\Delta t), \\ &y_n(t+\Delta t) = y_n(t)+ Norm(0,2D\Delta t), \\ &z_n(t+\Delta t) = z_n(t)+ Norm(0,2D\Delta t), \end{aligned}$$
where $Norm(0,2D\Delta t)$ is a random sample drawn from a normal distribution with mean 0 and variance $2D\Delta t$, where $D$ is the diffusion coefficient which is related to $\tau _c$ by $D = 1/\sqrt {\tau _c k^2}$ ($k$ is the wavenumber of light (i.e. $2\pi$/$\lambda$) and we consider a wavelength of light of 0.785 $\mu m$), and $\Delta t$ is the simulation time step. A periodic boundary condition is used for the situation when a particle moves outside of the volume, such that the particle re-enters the volume on the other side. A camera is fixed at a distance $z_{cam}$ (3000 $\mu m$ in this manuscript) away from the center of the volume to record the speckle patterns. The speckle patterns are sampled at discrete sampling locations. At a particular time $t$, the field at a sampling location ($x_{samp}, y_{samp}, z_{cam}$) on the camera is calculated as
$$\begin{array}{r}E(x_{samp}, y_{samp},z_{cam}, t) = \sum\limits_{n=1}^N \frac{exp(ikr_n)}{r_n} \\ r_n^2 = (x_n-x_{samp})^2+(y_n-y_{samp})^2+(z_n-z_{cam})^2, \end{array}$$
where $r_n$ is the distance between the $n^{th}$ particle and the sampling location on the camera. The sensor’s size is $25\;\mu m\;\times \;25\;\mu m$ with a pixel size of 0.5 $\mu m$, both of which can be adjusted. In principle, the number of scattering particles $N_p$ needs to be larger than the total number of speckles on the sphere surface of the detector that encloses the scattering medium to ensure enough $NIO$ and prevent any long range correlations among the simulated speckles. But using a large number of scattering particles ($N_p \sim 10^7$ for the parameters in this manuscript) can be computationally expensive. In the simulation, we need to choose the number of particles used to generate the speckle patterns. We observe that $N_p=5000$ is sufficient for the detector size and $s/p$ used and introduces $<4\%$ long range correlations (Supplementary Material Fig. S4). This ensures long range correlations have a small impact on our fundamental speckle contrast estimates. Additionally, the speckle intensity probability distribution function obtained with $N_p = 5000$ particles for large $s/p$ is a single exponential decay function (Supplementary Material Fig. S5) consistent with [21]. Also, we have observed that $f(s/p)$ approaches $1$ at small $s/p$ (Fig. 2(d)), indicating the effect of this long range correlation is not significant. Experimentally, the intensity at a pixel is proportional to the integration of all the photons that fall within the area of the pixel. Numerically, the intensity at each pixel is calculated as the mean of the speckle intensity at a few discrete sampling points within the pixel size $p$,
$$I_p(x,y,t) = \frac{1}{N_{SampPerPixel}^2 } \quad \sum_{|x_{samp}-x| < p/2, |y_{samp}-y|< p/2}^{} |E(x_{samp}, y_{samp},z_{cam}, t)|^2$$
where $N_{sampPerPixel}^2$ are the total number of sampling points within each pixel. When $s/p<1$, the number of sampling points per speckle, $N_{SampPerSpeckle}$, needs to be large enough to properly account for the intensity variations with a speckle size. These two parameters are simply related via $N_{SampPerSpeckle} = (s/p) * N_{SampPerPixel}$. We have used $N_{SampPerSpeckle}= 4$ to satisfy the Nyquist theorem. The speckle size $s$ in the simulation increases linearly with increasing $z_{cam}$. To calculate the temporal blurring of the speckle pattern obtained for a specified $T_{exp}$, the intensity at a pixel $p$ is integrated from $t = 0$ to $T_{exp}$ as
$$I_p(x,y,T_{exp}) = \int_{0}^{T_{exp}} I_p(x, y,t) dt = \sum_{k=1}^{T_{exp}/\Delta t} I_p(x, y, k \Delta t).$$
Here the time step $\Delta t$ has to be smaller than $\tau _c$ to obtain enough temporal sampling points, and we have used $\Delta t$ = 0.1$\tau _c$. To account for the noise terms, we need to scale the continuous speckle intensity obtained from the DSM simulation to the discrete camera counts on the detector. The camera count $Z_p(T_{exp})$, for the pixel $p$ (in photo-electrons $e^-$) integrated over time $T_{exp}$ is [26,27]
$$Z_p(T_{exp}) = N_{pe}(T_{exp}) + N_{de}(T_{exp}) + N_{re}$$
where $N_{pe}$ is the number of photo-electrons from the fundamental signal and shot noise. Here we consider a camera gain value of 1 ADU/$e^-$. If the average number of photons within a bin time $\Delta t$ is $N_{bin}$, we normalize the intensity $I_p$ obtained from DSM such that the average value of $\langle I_p(t_n)\rangle _n=N_{bin}$, where $I_p(t_n)$ is the intensity within the $n^{th}$ bin time $t_n$. Then we perform a Poisson draw for each exposure time as, $N_{pe}(T_{exp})=poissrnd(Q_eI_p(T_{exp}))$, and $I_p(T_{exp})$ is the sum of $I_{p}(t_n)$ within the exposure time $T_{exp}$. Here $\textit {poissrnd}$ is the Poisson draw function in MATLAB. $N_{de}$ is the number of dark electrons that depends on $T_{exp}$ and follows a Poisson distribution. $N_{re}$ is the number of electrons generated by the read noise that follows a Gaussian distribution with mean $B$ and standard deviation $\sigma _r$. In our simulations, the mean $B$ is set to zero since it is a constant offset subtracted from the measured intensity. The values of $N_{re}$ drawn from a Gaussian distribution are rounded to the nearest integer value. The time-dependence of the camera count, $Z_p$, is implicit. The variance of $Z_p$ can be expressed as the sum of variances of the individual noise terms.
$$var(Z_p) = var(N_{pe}) + var(N_{de}) + var(N_{re}) .$$
The first term in Eq. (19) can be written as:
$$var(N_{pe}) = var(Q_e I_p (T_{exp})) +{<}Q_e I_p(T_{exp}) > \\$$
where the first term on the right hand side comes from the fundamental signal and the second term arises from shot noise [27]. The variance in the number of electrons from dark noise and read noise is given by:
$$var(N_{de}) ={<} I_d T_{exp}> $$
$$var(N_{re}) = \sigma_r^2 $$
where $I_d$ is the dark count rate which is set to be zero in this manuscript, since it is in general small for the $T_{exp}$ we encounter for human brain measurements. Here $\sigma _r$ is the read noise. In our model, the raw speckle contrast squared $K_{all}^2=var(Z_p)/<Z_p>^2$ and $\sigma (K_{all}^2)$ are calculated over $w^2$ pixels (see Fig. 1(c)). Note that $\sigma (K_{all}^2)$ decreases as the square-root of the number of pixels used to estimate $K_{all}^2$. In the DSM, we calculate $K_{all}^2$ and $\sigma (K_{all}^2)$ over a small number of pixels $w^2$ for computational efficiency. DSM runs a configuration of 10 ms in 1 $\mu$s time steps (i.e. 10,000 time points), $\tau _c$= 10 $\mu s$, $40 \times 40$ sampling points in $\sim$ 10 min on a computer equipped with 2.4 GHz Intel Xeon E5-2680v4 processor using 8-core and 8 GB of RAM. To obtain the results for a larger number of pixels on the order of $10^7$, we scale $\sigma (K_{all}^2)$ as $\sim \sqrt {1/w^2}$.

 figure: Fig. 1.

Fig. 1. (a) Schematic representation of the dynamic speckle model (DSM). The particles undergo Brownian motion in a box with dimensions [$x_{max}\;y_{max}\;z_{max}$]. The speckle pattern is obtained as a superposition of the spherical waves emitted from all the moving particles detected by a camera at a distance $z_{cam}$ from the center of the scattering volume that contains all the particles. The camera has $N_{pixels}$ along both x and y axes. (b) The intensity measured at a pixel $p$ is obtained by numerically integrating the intensity over the area of the pixel, where $x_{samp}$, $y_{samp}$ are the discrete sampling points. (c) Illustration of normalized speckle patterns obtained at different $T_{exp}$ using DSM. $K^2(T_{exp})$ is calculated from the speckle patterns for all of the $N$ unique realizations of the simulations. Then $\sigma (K^2)$ is obtained as the standard deviation of the $N$ instances of $K^2$ values at each $T_{exp}$.

Download Full Size | PDF

 figure: Fig. 2.

Fig. 2. (a) Intensity autocorrelation function $g_2(\tau )$ obtained from DSM (solid line) for varying lag times $\tau$ at $s/p = 1$ and theoretical expected behavior ($g_2(\tau )$= 1+$\beta$ exp(−2$\tau /\tau _c$)) (dashed line): $\tau _{c,expected}$= 10$\mu s$, $\tau _{c,fitted}$ =10.04 $\mu s$. (b) The fundamental contrast squared $K_f^2$ versus camera exposure time for speckle dynamics numerically calculated at $s/p = 1$ with the DSM (solid line) is in good agreement with the theoretical prediction (dashed line) (see text for details). (c) The SNR of the fundamental signal $K_f^2/\sigma (K_f^2)$ for a fixed number of speckles (m = $1.74 \times 10^7$) as a function of $T_{exp}$ obtained from DSM (solid line) for $s/p=0.5$, $1$ as compared with the analytical expression (dashed line) from Eq. (6a) to obtain NIO and thus $f(s/p)$ from fitting. The results of $f(s/p)$ are $1.02$ and $0.78$ for $s/p=0.5$ and $1$ respectively. (d) $f(s/p)$ obtained from DSM obtained from Eqs. (6a), (6b), (6c) and 8 as compared with the empirical expression in Eq. (9). The coefficient obtained from fitting is $c_{s/p}=0.292$.

Download Full Size | PDF

Our DSM is able to correctly model the full temporal statistics as predicted by the analytical expression $g_2(\tau ) = \langle I(t)I(t+\tau )\rangle / \langle |I(t)|\rangle ^2= 1+\beta |g_1(\tau )|^2=1 + \beta$exp$(-2\tau /\tau _c)$ [28] shown in Fig. 2(a). In DSM, the intensity autocorrelation function $g_2(\tau )$ is calculated for each pixel and averaged over the entire window with $w^2$ pixels. We also compared the DSM result of $K_f^2$ with the analytical model in Eq. (2) as shown in Fig. 2(b), with $K_f^2=var(I_p)/\langle I_p\rangle ^2$ using the intensity $I_p$ simulated from DSM before converting it to the camera counts $Z_p$. To compare DSM and the analytical noise model, we have obtained results of $K_f^2/\sigma (K^2_f)$ as shown in Fig. 2(c) to obtain NIO from Eq. (6a) for a particular $s/p$ value and the corresponding $f(s/p)$ from Eq. (8). The number of speckles $m$ is fixed at $m$ = $1.74 \times 10^7$, which is on the same order of the number of fiber modes estimated for a 2.5 mm diameter fiber bundle comprised of $\sim$2000 fibers with a 50 $\mu m$ core diameter and a 0.66 numerical aperture (NA) with the wavelength of light $\lambda$ = 0.785 $\mu m$ typically used in functional near infrared spectroscopy. In Fig. 2(d), $f(s/p)$ as a function of $s/p$ from DSM is compared with the empirical model (Eq. (9)) and the constant $c=0.292$ is calibrated from fitting. We can use both the DSM and the analytical noise model with the calibrated $f(s/p)$, $h(T_{fr}/\tau _c)$ parameters to quantify the performance of a particular SCOS system.

2.3 Experimental setup

We have developed a prototype fb-SCOS system to validate the predictions of DSM simulations. The system schematic is shown in Fig. 3. A narrow linewidth laser source (CrystaLaser, $\lambda$ = $785$ nm) [17] is coupled to a multimode fiber that delivers light to the dynamic phantom sample. The source fiber and the detector fiber are always in contact with the surface of the phantom. The other end of the multimode detection fiber (1000 $\mu m$ core diameter, NA = 0.39) is imaged onto a sCMOS camera (Orca Fusion BT, Hamamatsu, pixel size 6.5 $\mu m$) by an objective (20x 0.4 NA, Olympus) and zoom lens (12.5 - 75 mm FL 6X Manual Zoom Video Lens, Edmund Optics). The effective magnification and thus $s/p$ of this system can be adjusted by changing the focal length of the zoom lens. Here $s/p$ is estimated using a static phantom from the empirical relationship of $s/p$ to $K^2$ given by $K^2=\beta _{l}(s/p)^2 / (1 + (s/p)^2 )$ (Fig. S1), where $\beta _l$ = 0.5 takes into account the two polarizations. The dynamic sample used for imaging is fabricated by using 8% v/v concentration of Intralipid 20% (batch 10LI4282, Fresenius Kabi, Germany) in deionized water at room temperature. The decay time of the sample is $\sim 13\mu$s at a source-detector separation of $\sim$ 10 mm.

 figure: Fig. 3.

Fig. 3. Schematic of the fb-SCOS experimental setup to measure speckle contrast for varying $s/p$ by changing the focal length of the zoom lens. MMF = Multi-Mode Fiber.

Download Full Size | PDF

3. Results

3.1 Experimental validations

We have measured $K_{all}^2$ and $\sigma (K_{all}^2)$ for different zoom lens focal lengths (FL) resulting in different $s/p$ ratios experimentally using the setup described in the Section 2.3, and compared with the numerical results obtained using the DSM. Specifically, we have used a FL = 50 mm and 75 mm with the $s/p$ ratios calibrated using a static phantom to be $s/p=$ 0.6, 0.9 respectively. Here $K_{all}$ is calculated as the average of the contrast values obtained for each 7 by 7 sliding window within the region of interest (ROI) to minimize the impact of potential non-uniformity in illumination. Then $\sigma (K_{all}^2)$ are obtained from 500 speckle patterns. The results of $K_{all}^2$, $K_f^2$, and $\sigma (K_{all}^2)$ versus $T_{exp}$ measured from the dynamic phantom is shown in Fig. 4, where $K_f^2$ is calculated as $K_{all}^2-K_s^2-K_r^2$ for each 7 by 7 window and then averaged over all the windows, with $K_s^2$, and $K_r^2$ obtained from Eqs. (4), (5), respectively. The decay time $\tau _c$ of the sample is $\sim 13\;\mu$s which is obtained by fitting Eq. (2) to the calculated $K_f^2$ with $\beta$ obtained from the static phantom results for a given $s/p$. We then input all the parameters from the experiments including $\tau _c$, $s/p$, $Q_e$, $T_{exp}$, and $c_s$ to DSM to simulate speckle patterns and obtain $K_{all}^2$, $K_f^2$ and $\sigma (K_{all}^2)$ numerically. We see in Fig. 4 that the DSM results of $K_{all}^2$, $K_f^2$, and $\sigma (K_{all}^2)$ as functions of the $T_{exp}$ are in good agreement with experimental results.

 figure: Fig. 4.

Fig. 4. (a)-(b) Experimental results of $K^2$ as a function of exposure time at 50 mm zoom lens calibrated to (a) $s/p$ = $0.6$ and (b) 75 mm zoom lens calibrated to $s/p$ = $0.9$. Experimental noise subtracted $K_f^2$ (green line) is fitted to $K_f^2$ from DSM (magenta line) which follows Eq. (2) to extract $\tau _c$. Experimentally measured $K_{all}^2$ (blue line) is consistent with Eq. (3) (red line) considering all noise sources. (c)-(d) $\sigma (K_{all}^2)$ as a function of exposure time $T_{exp}$ (blue line) are in good agreement with DSM results (red line) for (c) $s/p=0.6$ and (d) for $s/p=0.9$. The parameters in the simulation matching the experimental conditions are: $\tau _c\sim$ 13 $\mu s$, $w^2$ = 25 $\times$ 25 pixels, $c_s \;\sim$ 250,000 photons/s. Here each data point of $K_{all}^2$, $K_{f}^2$ and $\sigma (K_{all}^2)$ are obtained by averaging over 500 experimental data points repeated 5 times and by running 25 unique instances of DSM simulations repeated 4 times .

Download Full Size | PDF

3.2 Modeling results

We have used DSM as well as the calibrated analytical model to evaluate the impact of each of the noise sources on $K^2$, $\sigma (K_{all}^2)$, and the signal to noise ratio $K_f^2/\sigma (K_{all}^2)$ with various $T_{exp}$ and $s/p$ values. Specifically, we have considered two commercially available cameras: Orca Fusion BT (Hamamatsu) with 2304 $\times$ 2304 pixels, 89 Hz frame rate, $\sigma _r$ = 1.6 $e^-$ [29], and the acA2000-340kmNIR (Basler) with 2048 $\times$ 1088 pixels, 340 Hz frame rate, $\sigma _r$ = 14.2 $e^-$ [30]. These cameras differ in acquisition frame rates, read noise value, and costs (Hamamatsu: $\sim{\$} $25k, Basler: $\sim{\$}$4k). In Fig. 5, we first show the results of the speckle contrast $K_f^2$, $K_f^2+K_s^2$, $K_f^2+K_{r}^2$, $K_{all}^2$ as functions of $T_{exp}$ for two typical $c_s$ values of 10,000 and 100,000 photons/s, $\tau _c$ = 10 $\mu s$, and $s/p$ = 1. We see that the DSM and analytical results perfectly align with each other. For the case of shot noise, Eq. (18) is set as $Z_p=N_{pe}$ to obtain $K_f^2+K_s^2$. To evaluate the impact of read noise, we obtain $K_f^2+K_r^2=K_{all}^2-K_s^2$. Here $K_s^2$ is obtained by subtracting $K_f^2+K_s^2$ (calculated from $Z_p=N_{pe}$) and $K_f^2$ (calculated from $Z_p= Q_eI_p(T_{exp})$). For a dynamic sample, $K_f^2$ and all the noise contributions to $K_{all}^2$ decreases with $T_{exp}$. As shown in Fig. 5, the dominant contribution to $K_{all}^2$ (red line) is read noise (blue line) at small $T_{exp}$ and switches to shot noise (black line) at longer $T_{exp}$. The cross-over happens at a smaller $T_{exp}$ for a larger $c_s$. The higher read noise value for the Basler camera moves the cross-over point to a longer $T_{exp}$ as compared to the Hamamatsu camera, and we see that at $c_s=10,000$ photons/s, read noise is the dominate noise source for all the $T_{exp}$ values we have evaluated. To minimize the contribution of instrument noise sources to $K_{all}^2$, $T_{exp}$ in principle needs to be set at a value larger than the cross over $T_{exp}$ experimentally. However, since the signal $K_f^2$ itself also decreases with $T_{exp}$, it has to be compared with the noise $\sigma (K^2)$ to determine the optimized $T_{exp}$ for best performance, as will be demonstrated later in Figs. 7 and 8. In Fig. 6, $K_f^2$, $K_f^2+K_s^2$, $K_f^2+K_{r}^2$, $K_{all}^2$ are calculated with both the DSM and analytical models for varying $s/p$ at a $T_{exp}$ for the maximum frame rate of these cameras ($T_{exp} = 11.2$ ms for the Hamamatsu and $T_{exp} = 2.9$ ms for the Basler camera) at $c_s$ = 10,000 and 100,000 photons/s. The DSM and the analytical models are in perfect agreement with each other. With an increasing $s/p$ value, the signal $K_f^2$ increases until $s/p\sim 2$, which provides enough sampling of the speckle pattern from Nyquist theorem [21]. Since the intensity per pixel decreases with increasing $s/p$ for a fixed $c_s$, the contribution of the noise sources to the bias in $K^2$ also increase, as seen in Fig. 6, $K_{all}^2$ (red line) deviates more from $K_{f}^2$ (magenta line).

 figure: Fig. 5.

Fig. 5. The square of the speckle contrast $K_f^2$ (magenta), $K_f^2+K_s^2$ (black), $K_f^2+K_r^2$ (blue), $K_{all}^2$ (red) corresponding to fundamental signal, fundamental + shot noise, fundamental + read noise and all noise sources respectively as a function of $T_{exp}$ for (a) the Hamamatsu camera with $\sigma _r$ = 1.6 $e^-$ and $c_s=10,000$ photons/s, (b) the Hamamatsu camera with $\sigma _r$ = 1.6 $e^-$ and $c_s=100,000$ photons/s, (c) the Basler camera with $\sigma _r$ = 14.2 $e^-$ and $c_s=10,000$ photons/s and (d) the Basler camera with $\sigma _r$ = 14.2 $e^-$ and $c_s=100,000$ photons/s. The DSM results are represented by cross markers and the analytical results are given by solid lines. Other parameters: number of speckles $m= 1.74 \times 10^7$, $s/p$ = 1, $\tau _c$ = 10 $\mu s$, $Q_{e,Hamamatsu}$ = 0.55, $Q_{e,Basler}$ = 0.45. Read noise values are provided by the manufacturers

Download Full Size | PDF

 figure: Fig. 6.

Fig. 6. The square of the speckle contrast $K_f^2$ (magenta), $K_f^2+K_s^2$ (black), $K_f^2+K_r^2$ (blue), $K_{all}^2$ (red) corresponding to fundamental signal, fundamental + shot noise, fundamental + read noise and all noise sources respectively as functions of $s/p$ at long exposure time for (a) the Hamamatsu camera with $\sigma _r$ = 1.6 $e^-$, $T_{exp} = 11.2$ ms and $c_s=10,000$ photons/s, (b) the Hamamatsu camera with $\sigma _r$ = 1.6 $e^-$, $T_{exp} = 11.2$ ms and $c_s=100,000$ photons/s, (c) the Basler camera with $\sigma _r$ = 14.2 $e^-$, $T_{exp} = 2.9$ ms and $c_s=10,000$ photons/s and(d) the Basler camera with $\sigma _r$ = 14.2 $e^-$, $T_{exp} = 2.9$ ms and $c_s=100,000$ photons/s. The DSM results are represented by cross markers and the analytical results are given by solid lines Other parameters: $m = 1.74 \times 10^7$, $\tau _c$ = 10 $\mu s$, $Q_{e,Hamamatsu}$ = 0.55, $Q_{e,Basler}$ = 0.45. Note that we have simulated more data points for smaller $s/p$ values where the variation of $K^2$ versus $s/p$ is relatively larger.

Download Full Size | PDF

 figure: Fig. 7.

Fig. 7. The signal to noise ratio ($K_f^2$/$\sigma (K_{all}^2)$) as a function of exposure time for a single camera frame and frame averaging (a) over 89 Hz for the Hamamatsu camera with $\sigma _r$ = 1.6 $e^-$ and $c_s$ = 10,000 photons/s, (b) over 89 Hz for the Hamamatsu camera with $\sigma _r$ = 1.6 $e^-$ and $c_s$ = 100,000 photons/s, (c) over 340 Hz for the Basler camera with $\sigma _r$ = 14.2 $e^-$ and $c_s$ = 10,000 photons/s and (d) over 340 Hz for the Basler camera with $\sigma _r$ = 14.2 $e^-$ and $c_s$ = 100,000 photons/s. The analytical noise model results are indicated by blue dotted line. Parameters: $m = 1.74 \times 10^7$, $s/p$ = 1, $T_{meas}$ = 100 ms, Hamamatsu frame rate = 89 Hz, Basler frame rate = 340 Hz, $Q_{e,Hamamatsu}$ = 0.55, $Q_{e,Basler}$ = 0.45, 300 unique configurations of speckle patterns are simulated to estimate $K_f^2$ and $\sigma (K_{all}^2)$.

Download Full Size | PDF

 figure: Fig. 8.

Fig. 8. The signal to noise ratio ($K_f^2$/$\sigma (K_{all}^2)$) as a function of $s/p$ for single frame and frame averaging (a) over 89 Hz for the Hamamatsu camera with $\sigma _r$ = 1.6 $e^-$, $T_{exp} = 11.2$ ms and $c_s$ = 10,000 photons/s, (b) over 89 Hz for the Hamamatsu camera with $\sigma _r$ = 1.6 $e^-$, $T_{exp} = 11.2$ ms and $c_s$ = 100,000 photons/s, (c) over 340 Hz for the Basler camera with $\sigma _r$ = 14.2 $e^-$ and $T_{exp} = 2.9$ ms and $c_s$ = 10,000 photons/s, (d) over 340 Hz for the Basler camera with $\sigma _r$ = 14.2 $e^-$ and $T_{exp} = 2.9$ ms and $c_s$ = 100,000 photons/s, over a 100 ms measurement time. The analytical noise model results are indicated by blue dotted line. The peak in SNR is at $s/p$= 1.1 and 1.2 respectively for 10,000 and 100,000 photons/s for the Hamamatsu camera and at $s/p$= $<0.2$ and 0.6 for 10,000 and 100,000 photons/s for the Basler camera. Parameters: $m = 1.74 \times 10^7$, $T_{meas}$ = 100 ms, Hamamatsu frame rate = 89 Hz, Basler frame rate = 340 Hz, $Q_{e,Hamamatsu}$ = 0.55, $Q_{e,Basler}$ = 0.45, 300 unique configurations of speckle patterns are simulated to estimate $K_f^2$ and $\sigma (K_{all}^2)$.

Download Full Size | PDF

We have then calculated SNR defined as $K_f^2/\sigma (K_{all}^2)$ using DSM and the calibrated analytical model for the two cameras described above with $c_s$ = 10,000 photons/s and 100,000 photons/s in Figs. 7 and 8 respectively. Experimentally, we can also perform temporal averaging of frames within a total measurement time $T_{meas}$ to reduce the contribution of noise and the noise follows Eq. (13). Results are shown for a fixed number of speckles of $m = 1.74 \times 10^7$ at $\lambda$ = 0.785 $\mu m$ which is of the order of the fiber modes estimated for a fiber bundle typically used for functional near infrared spectroscopy measurements [31], with a diameter of 2.5 mm consisting of multiple strands of 50 $\mu m$ core diameter fibers with an NA = 0.66. The calculated SNR is shown for a single camera frame, i.e., $T_{exp} = T_{meas}$ (black line) and for the average of multiple frames over a $T_{meas}=100$ ms measurement time when the exposure time is less than the measurement time (red line). In Fig. 7, we see that increasing $T_{exp}$ for a single frame measurement (black line) increases the SNR and reaches a plateau at large $T_{exp}$, which marks the regime where camera read noise is negligible and the system is in the shot and fundamental noise limit (Fig. 7(a,b,d)). For the Basler camera at low photon flux (Fig. 7(c)), the read noise dominates till $T_{exp}$= 20 ms and thus SNR does not plateau. When averaging over multiple camera frames (red line) is considered within a fixed $T_{meas}=100$ ms, SNR peaks at a $T_{exp}$ determined by the inverse of the maximum camera frame rate, i.e. $T_{exp} = 11.2$ ms for the Hamamatsu and $T_{exp} = 2.9$ ms for the Basler camera if the SNR does not plateau before these exposure times.

In Fig. 8, SNR is calculated for varying $s/p$ at the $T_{exp}$ given by the inverse of the camera frame rate, i.e., $T_{exp} = 11.2$ ms for the Hamamatsu and $T_{exp} = 2.9$ ms for the Basler camera for a single camera frame (black line) and average of multiple frames over $T_{meas}=100$ ms based on the frame rate of the camera (red line). There exists an optimized $s/p$ value at which SNR is the highest which marks the transitional point where camera read noise starts to dominate due to decreased $c_p$ with increasing $s/p$ for a fixed $c_s$. As $c_s$ increases, the transitional point moves to larger $s/p$ as expected. The optimized $s/p$ value is $s/p=1.1$ at $c_s=10,000$ (Fig. 8(a)) and shifts to $s/p=1.2$ at $c_s=100,000$ (Fig. 8(b)) for the Hamamatsu camera and $s/p < 0.2$ (Fig. 8(c)) at $c_s=10,000$ and shifts to $s/p=0.6$ at $c_s=100,000$ (Fig. 8(d)) for the Basler camera. From Figs. 7 and 8, we observe that there exists an optimized $s/p$, $T_{exp}$ value for a particular $c_s$ used for each camera.

4. Discussion

We have introduced the numerical tool of DSM, which simulates the evolution of the speckle patterns by calculating the superposition of the light fields scattered from multiple source locations. Using this model, we have evaluated the impact of fundamental and instrumental noise sources in LSCI and SCOS measurements that quantify the blood flow by the speckle contrast within a certain $T_{exp}$. We have identified the optimal parameter space for $s/p$ and $T_{exp}$ for two commercial CMOS cameras: Hamamatsu (Orca Fusion BT) and Basler (acA2000-340kmNIR) for fb-SCOS studies relevant for human brain measurements.

The DSM generates dynamic speckle patterns by diffusing the locations from which the multiple light fields are emitted. This produces a temporally correlated evolution of the speckle patterns that mimic those typically measured on the surface of the tissue with large source-detector separations of several millimeters or more. Note that DSM should be considered as a fast-computing tool to generate speckle pattern evolution using the summation of random spherical wavelets in the single scattering regime. This is able to generate the same universal speckle statistics as expected in the human brain measurement by mimicking the multiple scattering in the tissue which is generally computationally expensive. Thus, the diffusion coefficient of sources in Eq. (14) is not the same as the diffusion coefficient of the tissue, but the two diffusion coefficients are related to a single parameter of $\tau _c$. An additional important feature is that the speckle to pixel size ratio $s/p$ can be precisely tuned by changing the distance between the scattering centers and the simulated camera and/or by changing the simulated pixel size on the camera over which the speckle intensity is integrated. This makes it easy to use the DSM to investigate the effect of the $s/p$ ratio on the speckle intensity statistics. The DSM does not presently include the effects of static scattering that can arise for free space LSCI measurements of mouse brain imaging, but this is straightforward to be included in the model with the addition of scattering centers that do not move in time. But static scattering is generally not relevant for measurements made on humans with source-detector separations greater than a few millimeters.

In the noise model and the DSM, we neglected dark noise, fixed pattern noise, camera saturation effects and the contributions of ambient light. We did not consider these factors since they are negligible or they can be corrected experimentally to have minimal impact on the performance of LSCI and SCOS. These parameters can be incorporated for experimental conditions when they have significant impacts on the noise statistics. Our results indicate that increasing $T_{exp}$ decreases the relative contribution of read and dark noise making the signal shot-noise dominant. and the SNR plateaus at long $T_{exp}$. In reality, we cannot go to a $T_{exp}$ that is too long such that $K^2$ is too small to be accurately measured due to digitization errors of the detector, which is not discussed in this manuscript. We have compared $K^2$ obtained from the digitized and non-digitized intensities and found that for the parameters including the photon flux values we have utilized and the properties of the two cameras, digitization error is negligible (Supplementary material Fig. S6). In reality, increasing $T_{exp}$ is advantageous till the plateau region is achieved and going beyond that $T_{exp}$ does not give any advantage for SNR. We have also demonstrated that when utilize temporal averaging to obtain $K^2$ within a fixed measurement time $T_{meas}$, there exist an optimized $T_{exp}$ that provides the best performance for a particular LSCI or SCOS system.

In summary, we have developed a noise model for a fb-SCOS system considering fundamental and instrumental noise sources and evaluated the performance of two commercially available high-speed CMOS cameras. Measurement conditions other than those been analyzed in this work can be achieved by adjusting parameters of the model which has been made publicly available [25]. We have provided guidelines for choosing experimental parameters relevant for human brain activation measurements using fb-SCOS. Ultimately, the number of parallel measurements that can be made is limited by the number of speckle modes available and it is interesting to note that CMOS cameras already have enough pixels to measure all available speckle modes for the typical detector areas utilized for measurements on humans. In the future, we will compare the potential of fb-SCOS with DCS [19,32] as well as other human brain imaging techniques including the traditional functional near infrared spectroscopy [33] and the recently developed interferometric techniques [34,35].

Funding

Meta Platforms Inc.

Acknowledgments

This research was supported by funding provided by Meta Platforms Inc.

Disclosures

Dr. Boas consulted with Meta Platforms Inc. on topics related to speckle contrast optical spectroscopy and diffuse correlation spectroscopy that ultimately led to the ideas presented in this paper. Dr. Boas’s interests were reviewed and are managed by Boston University in accordance with their conflict of interest policies.

Data Availability

The numerical model can be found in [25]. The data that support the results reported in the paper can be requested by contacting the corresponding author.

Supplemental document

See Supplement 1 for supporting content.

References

1. D. A. Boas, L. Campbell, and A. G. Yodh, “Scattering and imaging with diffusing temporal field correlations,” Phys. Rev. Lett. 75(9), 1855–1858 (1995). [CrossRef]  

2. T. Durduran and A. G. Yodh, “Diffuse correlation spectroscopy for non-invasive, micro-vascular cerebral blood flow measurement,” NeuroImage 85, 51–63 (2014). [CrossRef]  

3. E. M. Buckley, A. B. Parthasarathy, P. E. Grant, A. G. Yodh, and M. A. Franceschini, “Diffuse correlation spectroscopy for measurement of cerebral blood flow: future prospects,” Neurophotonics 1(1), 011009 (2014). [CrossRef]  

4. D. A. Boas and A. K. Dunn, “Laser speckle contrast imaging in biomedical optics,” J. Biomed. Opt. 15(1), 011109 (2010). [CrossRef]  

5. J. Senarathna, A. Rege, N. Li, and N. V. Thakor, “Laser speckle contrast imaging: theory, instrumentation and applications,” IEEE Rev. Biomed. Eng. 6, 99–110 (2013). [CrossRef]  

6. A. K. Dunn, “Laser speckle contrast imaging of cerebral blood flow,” Ann. Biomed. Eng. 40(2), 367–377 (2012). [CrossRef]  

7. D. Briers, D. D. Duncan, E. R. Hirst, S. J. Kirkpatrick, M. Larsson, W. Steenbergen, T. Stromberg, and O. B. Thompson, “Laser speckle contrast imaging: theoretical and practical limitations,” J. Biomed. Opt. 18(6), 066018 (2013). [CrossRef]  

8. D. D. Duncan and S. J. Kirkpatrick, “Can laser speckle flowmetry be made a quantitative tool?” J. Opt. Soc. Am. A 25(8), 2088–2094 (2008). [CrossRef]  

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

10. H. Cheng, Q. Luo, S. Zeng, S. Chen, J. Cen, and H. Gong, “Modified laser speckle imaging method with improved spatial resolution,” J. Biomed. Opt. 8(3), 559–564 (2003). [CrossRef]  

11. H. M. Varma, C. P. Valdes, A. K. Kristoffersen, J. P. Culver, and T. Durduran, “Speckle contrast optical tomography: A new method for deep tissue three-dimensional tomography of blood flow,” Biomed. Opt. Express 5(4), 1275–1289 (2014). [CrossRef]  

12. C. P. Valdes, H. M. Varma, A. K. Kristoffersen, T. Dragojevic, J. P. Culver, and T. Durduran, “Speckle contrast optical spectroscopy, a non-invasive, diffuse optical method for measuring microvascular blood flow in tissue,” Biomed. Opt. Express 5(8), 2769–2784 (2014). [CrossRef]  

13. T. Dragojević, H. M. Varma, J. L. Hollmann, C. P. Valdes, J. P. Culver, C. Justicia, and T. Durduran, “High-density speckle contrast optical tomography (scot) for three dimensional tomographic imaging of the small animal brain,” NeuroImage 153, 283–292 (2017). [CrossRef]  

14. C. P. Lin, I. E. Orukari, C. Tracy, T. Durduran, and J. P. Culver, “Exploring the feasibility of fiber-based speckle contrast optical spectroscopy,” in Microscopy Histopathology and Analytics, (Optical Society of America, 2020), pp. JW3A–37.

15. J. Xu, A. K. Jahromi, and C. Yang, “Diffusing wave spectroscopy: A unified treatment on temporal sampling and speckle ensemble methods,” APL Photonics 6(1), 016105 (2021). [CrossRef]  

16. X. Cheng, Y. Lockerman, and A. Z. Genack, “Phase singularity diffusion,” Opt. Lett. 39(11), 3348–3351 (2014). [CrossRef]  

17. D. D. Postnov, X. Cheng, S. E. Erdener, and D. A. Boas, “Choosing a laser for laser speckle contrast imaging,” Sci. Rep. 9(1), 2542 (2019). [CrossRef]  

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

19. C. Zhou, G. Yu, D. Furuya, J. H. Greenberg, A. G. Yodh, and T. Durduran, “Diffuse optical correlation tomography of cerebral blood flow during cortical spreading depression in rat brain,” Opt. Express 14(3), 1125–1144 (2006). [CrossRef]  

20. X. Cheng, E. J. Sie, S. Naufel, D. A. Boas, and F. Marsili, “Measuring neuronal activity with diffuse correlation spectroscopy: a theoretical investigation,” Neurophotonics 8(03), 035004 (2021). [CrossRef]  

21. S. J. Kirkpatrick, D. D. Duncan, and E. M. Wells-Gray, “Detrimental effects of speckle-pixel size matching in laser speckle contrast imaging,” Opt. Lett. 33(24), 2886–2888 (2008). [CrossRef]  

22. L. K. Frisk, M. Verma, C.-H. P. Lin, S. Chetia, F. Beslija, J. Trobaugh, J. P. Culver, and T. Durduran, “Modelling the experimental and fundamental signal-to-noise in speckle contrast optical tomography,” in Dynamics and Fluctuations in Biomedical Photonics XIX, vol. PC11959 International Society for Optics and Photonics (SPIE, 2022).

23. E. J. Sie, H. Chen, E.-F. Saung, R. Catoen, T. Tiecke, M. A. Chevillet, and F. Marsili, “High-sensitivity multispeckle diffuse correlation spectroscopy,” Neurophotonics 7(03), 035010 (2020). [CrossRef]  

24. W. B. Baker, A. B. Parthasarathy, T. S. Ko, D. R. Busch, K. Abramson, S.-Y. Tzeng, R. C. Mesquita, T. Durduran, J. H. Greenberg, D. K. Kung, and A. G. Yodh, “Pressure modulation algorithm to separate cerebral hemodynamic signals from extracerebral artifacts,” Neurophotonics 2(3), 035004 (2015). [CrossRef]  

25. S. Zilpelwar, D. Postnov, A. Chen, et al., “Dynamic speckle model,” Github, 2022, https://github.com/BUNPC/DynamicSpeckleModel.

26. M. Konnik and J. Welsh, “High-level numerical simulations of noise in ccd and cmos photosensors: review and tutorial,” arXiv, arXiv preprint arXiv:1412.4031 (2014). [CrossRef]  

27. B. Mandracchia, X. Hua, C. Guo, J. Son, T. Urner, and S. Jia, “Fast and accurate scmos noise correction for fluorescence microscopy,” Nat. Commun. 11(1), 94 (2020). [CrossRef]  

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

29. https://www.hamamatsu.com/us/en/product/type/C15440-20UP/index.html.

30. https://www.baslerweb.com/en/products/cameras/area-scan-cameras/ace/aca2000-340kmnir/.

31. F. Scholkmann, S. Kleiser, A. J. Metz, R. Zimmermann, J. M. Pavia, U. Wolf, and M. Wolf, “A review on continuous wave functional near-infrared spectroscopy and imaging instrumentation and methodology,” NeuroImage 85, 6–27 (2014). [CrossRef]  

32. X. Cheng, H. Chen, E. J. Sie, F. Marsili, and D. A. Boas, “Development of a monte carlo-wave model to simulate time domain diffuse correlation spectroscopy measurements from first principles,” J. Biomed. Opt. 27(08), 083009 (2022). [CrossRef]  

33. J. J. Selb, D. A. Boas, S.-T. Chan, K. C. Evans, E. M. Buckley, and S. A. Carp, “Sensitivity of near-infrared spectroscopy and diffuse correlation spectroscopy to brain hemodynamics: simulations and experimental findings during hypercapnia,” Neurophotonics 1(1), 015005 (2014). [CrossRef]  

34. W. Zhou, O. Kholiqov, J. Zhu, M. Zhao, L. L. Zimmermann, R. M. Martin, B. G. Lyeth, and V. J. Srinivasan, “Functional interferometric diffusing wave spectroscopy of the human brain,” Sci. Adv. 7(20), eabe0150 (2021). [CrossRef]  

35. J. Xu, A. K. Jahromi, J. Brake, J. E. Robinson, and C. Yang, “Interferometric speckle visibility spectroscopy (ISVS) for human cerebral blood flow monitoring,” APL Photonics 5(12), 126102 (2020). [CrossRef]  

Supplementary Material (1)

NameDescription
Supplement 1       Supplementary Material

Data Availability

The numerical model can be found in [25]. The data that support the results reported in the paper can be requested by contacting the corresponding author.

25. S. Zilpelwar, D. Postnov, A. Chen, et al., “Dynamic speckle model,” Github, 2022, https://github.com/BUNPC/DynamicSpeckleModel.

Cited By

Optica participates in Crossref's Cited-By Linking service. Citing articles from Optica Publishing Group journals and other participating publishers are listed here.

Alert me when this article is cited.


Figures (8)

Fig. 1.
Fig. 1. (a) Schematic representation of the dynamic speckle model (DSM). The particles undergo Brownian motion in a box with dimensions [$x_{max}\;y_{max}\;z_{max}$]. The speckle pattern is obtained as a superposition of the spherical waves emitted from all the moving particles detected by a camera at a distance $z_{cam}$ from the center of the scattering volume that contains all the particles. The camera has $N_{pixels}$ along both x and y axes. (b) The intensity measured at a pixel $p$ is obtained by numerically integrating the intensity over the area of the pixel, where $x_{samp}$, $y_{samp}$ are the discrete sampling points. (c) Illustration of normalized speckle patterns obtained at different $T_{exp}$ using DSM. $K^2(T_{exp})$ is calculated from the speckle patterns for all of the $N$ unique realizations of the simulations. Then $\sigma (K^2)$ is obtained as the standard deviation of the $N$ instances of $K^2$ values at each $T_{exp}$.
Fig. 2.
Fig. 2. (a) Intensity autocorrelation function $g_2(\tau )$ obtained from DSM (solid line) for varying lag times $\tau$ at $s/p = 1$ and theoretical expected behavior ($g_2(\tau )$= 1+$\beta$ exp(−2$\tau /\tau _c$)) (dashed line): $\tau _{c,expected}$= 10$\mu s$, $\tau _{c,fitted}$ =10.04 $\mu s$. (b) The fundamental contrast squared $K_f^2$ versus camera exposure time for speckle dynamics numerically calculated at $s/p = 1$ with the DSM (solid line) is in good agreement with the theoretical prediction (dashed line) (see text for details). (c) The SNR of the fundamental signal $K_f^2/\sigma (K_f^2)$ for a fixed number of speckles (m = $1.74 \times 10^7$) as a function of $T_{exp}$ obtained from DSM (solid line) for $s/p=0.5$, $1$ as compared with the analytical expression (dashed line) from Eq. (6a) to obtain NIO and thus $f(s/p)$ from fitting. The results of $f(s/p)$ are $1.02$ and $0.78$ for $s/p=0.5$ and $1$ respectively. (d) $f(s/p)$ obtained from DSM obtained from Eqs. (6a), (6b), (6c) and 8 as compared with the empirical expression in Eq. (9). The coefficient obtained from fitting is $c_{s/p}=0.292$.
Fig. 3.
Fig. 3. Schematic of the fb-SCOS experimental setup to measure speckle contrast for varying $s/p$ by changing the focal length of the zoom lens. MMF = Multi-Mode Fiber.
Fig. 4.
Fig. 4. (a)-(b) Experimental results of $K^2$ as a function of exposure time at 50 mm zoom lens calibrated to (a) $s/p$ = $0.6$ and (b) 75 mm zoom lens calibrated to $s/p$ = $0.9$. Experimental noise subtracted $K_f^2$ (green line) is fitted to $K_f^2$ from DSM (magenta line) which follows Eq. (2) to extract $\tau _c$. Experimentally measured $K_{all}^2$ (blue line) is consistent with Eq. (3) (red line) considering all noise sources. (c)-(d) $\sigma (K_{all}^2)$ as a function of exposure time $T_{exp}$ (blue line) are in good agreement with DSM results (red line) for (c) $s/p=0.6$ and (d) for $s/p=0.9$. The parameters in the simulation matching the experimental conditions are: $\tau _c\sim$ 13 $\mu s$, $w^2$ = 25 $\times$ 25 pixels, $c_s \;\sim$ 250,000 photons/s. Here each data point of $K_{all}^2$, $K_{f}^2$ and $\sigma (K_{all}^2)$ are obtained by averaging over 500 experimental data points repeated 5 times and by running 25 unique instances of DSM simulations repeated 4 times .
Fig. 5.
Fig. 5. The square of the speckle contrast $K_f^2$ (magenta), $K_f^2+K_s^2$ (black), $K_f^2+K_r^2$ (blue), $K_{all}^2$ (red) corresponding to fundamental signal, fundamental + shot noise, fundamental + read noise and all noise sources respectively as a function of $T_{exp}$ for (a) the Hamamatsu camera with $\sigma _r$ = 1.6 $e^-$ and $c_s=10,000$ photons/s, (b) the Hamamatsu camera with $\sigma _r$ = 1.6 $e^-$ and $c_s=100,000$ photons/s, (c) the Basler camera with $\sigma _r$ = 14.2 $e^-$ and $c_s=10,000$ photons/s and (d) the Basler camera with $\sigma _r$ = 14.2 $e^-$ and $c_s=100,000$ photons/s. The DSM results are represented by cross markers and the analytical results are given by solid lines. Other parameters: number of speckles $m= 1.74 \times 10^7$, $s/p$ = 1, $\tau _c$ = 10 $\mu s$, $Q_{e,Hamamatsu}$ = 0.55, $Q_{e,Basler}$ = 0.45. Read noise values are provided by the manufacturers
Fig. 6.
Fig. 6. The square of the speckle contrast $K_f^2$ (magenta), $K_f^2+K_s^2$ (black), $K_f^2+K_r^2$ (blue), $K_{all}^2$ (red) corresponding to fundamental signal, fundamental + shot noise, fundamental + read noise and all noise sources respectively as functions of $s/p$ at long exposure time for (a) the Hamamatsu camera with $\sigma _r$ = 1.6 $e^-$, $T_{exp} = 11.2$ ms and $c_s=10,000$ photons/s, (b) the Hamamatsu camera with $\sigma _r$ = 1.6 $e^-$, $T_{exp} = 11.2$ ms and $c_s=100,000$ photons/s, (c) the Basler camera with $\sigma _r$ = 14.2 $e^-$, $T_{exp} = 2.9$ ms and $c_s=10,000$ photons/s and(d) the Basler camera with $\sigma _r$ = 14.2 $e^-$, $T_{exp} = 2.9$ ms and $c_s=100,000$ photons/s. The DSM results are represented by cross markers and the analytical results are given by solid lines Other parameters: $m = 1.74 \times 10^7$, $\tau _c$ = 10 $\mu s$, $Q_{e,Hamamatsu}$ = 0.55, $Q_{e,Basler}$ = 0.45. Note that we have simulated more data points for smaller $s/p$ values where the variation of $K^2$ versus $s/p$ is relatively larger.
Fig. 7.
Fig. 7. The signal to noise ratio ($K_f^2$/$\sigma (K_{all}^2)$) as a function of exposure time for a single camera frame and frame averaging (a) over 89 Hz for the Hamamatsu camera with $\sigma _r$ = 1.6 $e^-$ and $c_s$ = 10,000 photons/s, (b) over 89 Hz for the Hamamatsu camera with $\sigma _r$ = 1.6 $e^-$ and $c_s$ = 100,000 photons/s, (c) over 340 Hz for the Basler camera with $\sigma _r$ = 14.2 $e^-$ and $c_s$ = 10,000 photons/s and (d) over 340 Hz for the Basler camera with $\sigma _r$ = 14.2 $e^-$ and $c_s$ = 100,000 photons/s. The analytical noise model results are indicated by blue dotted line. Parameters: $m = 1.74 \times 10^7$, $s/p$ = 1, $T_{meas}$ = 100 ms, Hamamatsu frame rate = 89 Hz, Basler frame rate = 340 Hz, $Q_{e,Hamamatsu}$ = 0.55, $Q_{e,Basler}$ = 0.45, 300 unique configurations of speckle patterns are simulated to estimate $K_f^2$ and $\sigma (K_{all}^2)$.
Fig. 8.
Fig. 8. The signal to noise ratio ($K_f^2$/$\sigma (K_{all}^2)$) as a function of $s/p$ for single frame and frame averaging (a) over 89 Hz for the Hamamatsu camera with $\sigma _r$ = 1.6 $e^-$, $T_{exp} = 11.2$ ms and $c_s$ = 10,000 photons/s, (b) over 89 Hz for the Hamamatsu camera with $\sigma _r$ = 1.6 $e^-$, $T_{exp} = 11.2$ ms and $c_s$ = 100,000 photons/s, (c) over 340 Hz for the Basler camera with $\sigma _r$ = 14.2 $e^-$ and $T_{exp} = 2.9$ ms and $c_s$ = 10,000 photons/s, (d) over 340 Hz for the Basler camera with $\sigma _r$ = 14.2 $e^-$ and $T_{exp} = 2.9$ ms and $c_s$ = 100,000 photons/s, over a 100 ms measurement time. The analytical noise model results are indicated by blue dotted line. The peak in SNR is at $s/p$= 1.1 and 1.2 respectively for 10,000 and 100,000 photons/s for the Hamamatsu camera and at $s/p$= $<0.2$ and 0.6 for 10,000 and 100,000 photons/s for the Basler camera. Parameters: $m = 1.74 \times 10^7$, $T_{meas}$ = 100 ms, Hamamatsu frame rate = 89 Hz, Basler frame rate = 340 Hz, $Q_{e,Hamamatsu}$ = 0.55, $Q_{e,Basler}$ = 0.45, 300 unique configurations of speckle patterns are simulated to estimate $K_f^2$ and $\sigma (K_{all}^2)$.

Equations (26)

Equations on this page are rendered with MathJax. Learn more.

K f 2 = 2 β T e x p 0 T e x p g 1 2 ( τ ) ( 1 τ T e x p ) d τ ,
K f 2 = β τ c T e x p [ 1 + τ c 2 T e x p ( e x p ( 2 T e x p τ c ) 1 ) ] .
K a l l 2 = K f 2 + K s 2 + K d 2 + K r 2
K s 2 = 1 Q e c p T e x p .
K r 2 = σ r 2 ( Q e c p T e x p ) 2 .
σ ( K f 2 ) = K f 2 2 ( 1 + 2 K f 2 ) N I O ,
σ ( K s 2 ) = c K s K s 2 1 w 2 ,
σ ( K r 2 ) = c K r K r 2 1 w 2 ,
σ ( K a l l 2 ) = σ 2 ( K f 2 ) + σ 2 ( K s 2 ) + σ 2 ( K r 2 ) .
N I O = m ( s / p ) 2 f ( s / p ) = w 2 f ( s / p ) ,
f ( s / p ) = 1 / ( 1 + c s / p ( s / p ) 2 ) ,
S N R = K f 2 σ ( K a l l 2 ) .
σ f a ( K f 2 ) = K f 2 2 ( 1 + 2 K f 2 ) N I O f a ,
σ f a ( K s 2 ) = c K s K s 2 1 w 2 N f r ,
σ f a ( K r 2 ) = c K r K r 2 1 w 2 N f r ,
N I O f a = m ( s / p ) 2 f ( s / p ) N f r h ( T f r / τ c )
σ f a ( K a l l 2 ) = σ ( K a l l 2 ) / N f r .
x n ( t + Δ t ) = x n ( t ) + N o r m ( 0 , 2 D Δ t ) , y n ( t + Δ t ) = y n ( t ) + N o r m ( 0 , 2 D Δ t ) , z n ( t + Δ t ) = z n ( t ) + N o r m ( 0 , 2 D Δ t ) ,
E ( x s a m p , y s a m p , z c a m , t ) = n = 1 N e x p ( i k r n ) r n r n 2 = ( x n x s a m p ) 2 + ( y n y s a m p ) 2 + ( z n z c a m ) 2 ,
I p ( x , y , t ) = 1 N S a m p P e r P i x e l 2 | x s a m p x | < p / 2 , | y s a m p y | < p / 2 | E ( x s a m p , y s a m p , z c a m , t ) | 2
I p ( x , y , T e x p ) = 0 T e x p I p ( x , y , t ) d t = k = 1 T e x p / Δ t I p ( x , y , k Δ t ) .
Z p ( T e x p ) = N p e ( T e x p ) + N d e ( T e x p ) + N r e
v a r ( Z p ) = v a r ( N p e ) + v a r ( N d e ) + v a r ( N r e ) .
v a r ( N p e ) = v a r ( Q e I p ( T e x p ) ) + < Q e I p ( T e x p ) >
v a r ( N d e ) = < I d T e x p >
v a r ( N r e ) = σ r 2
Select as filters


Select Topics Cancel
© Copyright 2024 | Optica Publishing Group. All rights reserved, including rights for text and data mining and training of artificial technologies or similar technologies.