1. Introduction

Rough surfaces and scattering media generate a characteristic speckle pattern [1] when illuminated by coherent light, e.g. from a laser. By analyzing a time sequence of speckle patterns, valuable information can be retrieved on the sample evolution. Broadly speaking, one may distinguish between “static” speckle patterns generated by solid objects and “dynamic” speckles formed by soft matter systems (e.g. colloidal suspensions, emulsions, polymer solutions, biological tissues), whose components undergo Brownian motion, thereby continuously reconfiguring the scattered speckle pattern. In the former case, relevant to metrology and interferometry [2], a rigid displacement or a long wavelength deformation is often measured, e.g. in response to vibrations, applied load, or a change of temperature; speckle patterns are recorded onto a 2D detector such as a CCD or CMOS camera, using an imaging optics. By contrast, the microscopic (e.g. Brownian) dynamics of soft systems is quantified by g₂(τ) − 1, the autocorrelation function of the temporal fluctuations of the scattered intensity. In these dynamic light scattering measurements (DLS, a.k.a. photon correlation spectroscopy [3]) a point-like detector (e.g. a phototube) is placed in the far field, where it collects light within a few speckles.

Recent developments have made the distinction between these two research fields increasingly fuzzy. On the one hand, imaging geometries similar to those for static speckles have been used to detect motion, e.g. in vascular flow ([4, 5] and references therein). Motion is typically quantified by inspecting the contrast of the speckle pattern, <I²>/<I>² −1, where I is the local intensity and 〈 ··· 〉 is an average over a small region centered around the point of interest. The method is based on the fact that speckles are blurred and the contrast reduced in those regions where significant motion occurs [6]. Rigid motion has been measured also by laser speckle velocimetry [7, 8], by applying spatial cross-correlation methods to speckle images, an approach similar to particle imaging velocimetry [9] and image correlation velocimetry [10] often used in fluid mechanics. On the other hand, CCD and CMOS cameras are now routinely used as detectors for DLS, especially for samples exhibiting slow dynamics, i.e. speckle fluctuations on time scales from a fraction of second up to several hours. In the first implementations of these so-called multispeckle approaches [11, 12], the far field detection scheme of traditional DLS was used and the intensity correlation function was calculated from

 

Fig. 1 a): schematic representation of an image J that is a shifted version of image I (only four pixels are shown for clarity). The intensity in a given pixel of J may be obtained as a linear combination of (up to) four pixels of I, with weights proportional to the colored areas. b): quadrant detection scheme for locating the direction of the shift. The nine elements closest to the peak of the cross-correlation between I and J are represented here. The shaded areas are (in pixel units) α₁ = 0.485869913, α₂ = 0.545406041, α₃ = 0.25; see Sec. 6 for more details. c): typical speckle image obtained in a PCI experiment.

Download Full Size | PDF

As highlighted by this short overview, previous works have focussed either on the measurement of the rigid displacement of a set of scatterers, regardless of any relative motion between them, or, conversely, on their relative motion due to the microscopic dynamics, regardless of any average drift component. Here, we present a new method that allows one to quantify in a single measurement the contribution of each of these phenomena to the evolution of speckle patterns formed in the imaging geometry [17]. Figure 2, which will be discussed later, shows the essence of the method in a glimpse: the spatial cross-correlation of two speckle images taken at a time lag τ exhibits a peak, whose position and height yield the sample rigid shift and its internal dynamics over the time τ, respectively. Applications of this method to gels submitted to gravitational [18] or internal [19] stress have been presented in previous publications, but the method itself was not discussed there. In this paper, we provide a detailed description of the algorithm used to implement the method, addressing in particular the challenges inherent to PCI experiments, i.e. the reduced size of the speckles and the need to process a very large number of images in a reasonably short time. Finally, we demonstrate our method by testing it on a model system, a suspension of Brownian particles contained in a cell displaced by a motor, and by measuring the velocity profile and the microscopic dynamics of a sheared colloidal glass.

 

Fig. 2 a)–c) spatial cross-correlation between speckle images generated by a diluted Brownian suspension that is translated along the y direction by a motor. The speckle patterns are recorded on a CCD using an imaging collection optics (see text for more details). As the delay τ between pair of images is increased, the peak position shifts to larger y and its height decreases, due to the relative motion of the Brownian particles. d) cut of the cross-correlation along the Δx = 0 line. Curves are labeled by τ. Inset: same data, replotted as a function of spatial shift with respect to the peak position. e) same as d), but for a ground glass, a scatterer whose internal dynamics are frozen. Note that the peak height remains essentially constant as the sample is translated.

Download Full Size | PDF

2. Sub-pixel digital imaging correlation algorithm

The first step in our method is to find the local rigid displacement with sub-pixel resolution. To this end, we use a cross-correlation technique inspired by particle imaging velocimetry [9] and image correlation velocimetry [10]. A time series of images of the sample is taken, using a PCI setup (see Fig. 1(c) for an example). Each image is divided into a grid of regions of interest (ROIs), corresponding to square regions of side L in the sample. Under a rigid drift, the speckle pattern in a ROI at time t will appear in a shifted position in a successive image taken at time t + τ. If the scatterers undergo relative motion, in addition to a rigid shift, the shifted ROI will not be an identical copy of the original one. Still, the displacement field can be estimated by calculating by what amount a ROI of the second image has to be back-shifted in order to maximize its resemblance with the corresponding ROI of the first image. In practice, for a given ROI the shift along the horizontal and vertical directions, Δx and Δy, is determined by searching the maximum of corr[I, J], the spatial cross-correlation of the intensity of the two images, defined by

The position (k̄, l̄) of the global maximum of corr[J, I] yields the desired displacement along the direction of columns (x axis) and rows (y axis), Δx = l̄ and Δy = k̄ respectively, with pixel resolution. Several schemes have been proposed to improve this resolution, e.g. by calculating the centroid of corr[J, I], or by fitting its peak to a 2-dimensional analytical function such as a Gaussian. While both methods work well for broad, circularly symmetric peaks, they tend to be less robust when the peak is sharp or it has an asymmetric shape. The shape of the peak is determined by the spatial autocorrelation of the intensity pattern; for our speckle images, it depends on the shape and size of the speckles, which may not be symmetrical, depending on the shape of the illuminated sample volume and the imaging optics [1]. Moreover, the peak usually extends over just a few pixels, because one minimizes the speckle size in order to maximize the information content in the image. To overcome the limitations inherent to peak-based schemes, we use an alternative approach based on a least-square method that allows us to obtain the displacement field with a typical resolution of a few hundredths of a pixel, with no requirements on the shape or broadness of the peak and without using any fitting function.

The first step of our shift-finding algorithm is the same as in standard PIV methods: corr[J, I] is calculated and the pixel-resolved shift, (k̄, l̄), is determined from the position of its global maximum. The next step consists in the refinement of such displacement with sub-pixel resolution. For the sake of simplicity, let us first assume that J is simply a shifted version of I, with displacement (Δx = l̄ + δx, Δy = k̄ + δy), with |δx| < 1, |δy| < 1. The intensity J_r,c may then be expressed as a weighted average over a suitable set of pixels of the intensity of the image I. In principle, an infinite number of terms are needed to obtain J_r,c, if δx and δy are non-integer [21]. In practice, linear interpolation is usually sufficient to reconstruct J to a good approximation, thereby greatly simplifying the calculation. Using linear interpolation, the intensity J_r,c at a given pixel is expressed as the weighted sum over the (at most) four pixels of I that partially overlap with that pixel, as exemplified in Fig. 1(a). Thus,

The optimum shift is obtained by imposing ∂χ²/∂a_i = 0, i = 1,..,4. By exchanging the order of the sums in Eq. (6) one recognizes that a is the solution of b = M · a, with:

Once a is computed, the shifts along the x and y directions are calculated with subpixel precision according to

3. Dynamic light scattering: corrections tog₂(τ) − 1 for drifting samples

In order to quantify the internal dynamics, one needs to compute the intensity correlation function g₂ − 1 between image J and a shifted version of I, so as to avoid any artifact due to the rigid shift of the speckles. Denoting by I′ the image I shifted by (Δx, Δy), the (un-normalized) intensity correlation function corrected for the shift contribution is

We introduce below an alternative method that leads to a much faster algorithm, where the only computational cost is that of evaluating the kernel for a few points, with no need for the calculation of convolution and correlation functions. It is convenient to consider only nonnegative fractional shifts $\widetilde{\delta x}$, $\widetilde{\delta y}$, given by

Using Eqs. (13) and (14), the convolution product (12) may be rewritten as

Equation (21) is the central result of our method. It shows that the intensity correlation function corrected for the shift contribution can be simply obtained as a linear combination of a few terms of covar[J, I], weighted by the kernel. Since covar[J, I] has already been calculated to determine the shift, the extra cost is essentially limited to the evaluation of M² values of the kernel, which is typically negligible. Finally, we note that covar[J, I] vanishes on the length scale of the speckle size as its argument departs from (i_y, j_x), which is close to the location of the peak of the covariance. Hence, it is sufficient to take M on the order of a few speckle sizes (in units of pixels), because in Eq. (21) the contribution of the kernel for larger lags would be multiplied by a vanishingly small quantity. For example, we find that for images with a speckle size of about 5 pixels, the correction is virtually independent of M for M ≥ 8.

4. Experimental tests

We test our method on two samples: a suspension of Brownian particles loaded in a cell displaced by a motor, and a colloidal glass to which a shear deformation is applied. In both cases, no rigid displacement is detected when the samples are left at rest, indicating that other potential sources of drift, such as convection or sedimentation, are here negligible. The Brownian sample is a dispersion of polystyrene microspheres (radius R = 0.265 μm) in an aqueous solution of fructose at 75% weight fraction. The particle volume fraction is 10⁻⁵ and the sample is kept at a temperature T = 9 °C. The setup is described in [23]; we use the imaging geometry shown in Fig. SM1(a) of [13], where an image of the sample is formed onto a CCD detector using light scattered at θ = 90 deg, corresponding to a scattering vector q = 2.46 μm⁻¹. The field of view is 1820×364 μm² and images are acquired at a rate of 10 Hz, with an exposure time of 5 msec. To avoid any spurious increase of the base line of g₂ − 1, the speckle images are corrected for non-uniform illumination as detailed in [24]. The sample cell is attached to a motor that can impose a drift in the y direction at a controlled speed, v_y = 10 μms⁻¹.

Figures 2(a)–2(c) show the spatial cross-correlation calculated applying Eq. (2) to pairs of speckle images taken while displacing the sample, for three different time lags. For τ = 0 s, Eq. (2) yields the spatial autocorrelation of the speckle pattern: accordingly, a sharp peak of height one and centered at Δx = Δy = 0 is visible, whose FWHM ≈ 2.9 pixels provides the speckle size. As the lag increases, the peak position drifts in the y direction, due to the translation motion imposed by the motor. Additionally, its height decreases, due to the relative motion of the Brownian particles that reconfigure the speckle pattern. Figure 2(d) shows a cut of the cross-correlation peak along the Δx = 0 direction, for four values of τ. From this plot, it is clear that if g₂(τ) − 1 was to be computed from a purely temporal cross-correlation, as in Eq. (1), one would observe a spurious, fast decay, essentially due to the rigid shift only. This would correspond to follow corr[I(t), I(t + τ)] at Δx = Δy = 0, as a function of τ. By contrast, if the relative motion of the Brownian particles is to be obtained, one has to measure the height of the peak as it drifts, as in the method proposed here. The inset of Fig. 2(d) shows the same data, plotted as a function of distance along y with respect to the (subpixel) peak position. It is worth noting that the peak width remains constant, in contrast to what suggested (but not demonstrated, to our knowledge) in patent literature [25], where it was proposed that the peak would broaden with τ as a result of the internal motion of the scatterers. Thus, the relevant parameter for extracting the relative motion is indeed the peak height, not its width.

As a final remark on the temporal evolution of the cross-correlation for our Brownian particles, it should be mentioned that in the Seventies and Eighties of the past century several groups have studied the evolution of the speckle pattern generated by a static sample undergoing translational motion. These works typically focussed on the far field speckle pattern, for which a loss of correlation is observed even in the absence of any internal dynamics (see [27] for a review). This has to be contrasted with our experiments, where an imaging geometry is used. In this case, the loss of correlation is only due to the internal dynamics. To show this, we plot in Fig. 2(e) a cut of the cross-correlation peak measured for a translating ground glass, a perfectly static sample. Similarly to [27], we observe the speckle pattern scattered in the forward direction. At odd with previous works, however, we adopt an imaging geometry. We find that the height of the peak does not decreases significantly as the sample is translated, confirming that in the PCI geometry any loss of correlation originates from the relative motion of the scatterer, not their drift.

In Fig. 3(a), we show for the Brownian particles the displacement of the speckle pattern as a function of τ, obtained from the sub-pixel peak position averaged over 200 pairs of images (i.e. 20 s), taken while translating the sample. The data are very well fitted by a linear law (red line) as expected for motion at constant speed, thus indicating that our algorithm captures very well the drift component of the speckle pattern, from a fraction of a pixel up to tens of pixels. The error bars represent σ_y, the standard deviation of the temporal fluctuations of the peak position detected for a fixed lag time τ.For τ ≤ 4.2s, σ_y/Δy < 4% and the error bars are smaller than the symbol size, indicating that the detection of the peak position is consistently reliable, even when the peak height is as low as 0.1 or the displacement is just a fraction of a pixel. For τ = 5.6 s, the error bar is significantly larger, because the peak height is typically of order 0.05 − 0.06, comparable to the noise level. As a consequence, the detected peak position departs erratically from its true value for a significant fraction of the 200 pairs of images, leading to a large value of σ_y. Beyond τ = 5.6 s, no peak can be reliably found, thus preventing the displacement to be measured. The slope of the linear fit to the data is 3.30 ± 0.03 pixel s⁻¹. Recalling that the nominal speed of the motor is v_y = 10 μm s⁻¹, this implies that 1 pixel corresponds to 3.03 ± 0.03 μm in the sample, in excellent agreement with 3.15 ± 0.15 μm/pixel as obtained from the magnification of the imaging system, evaluated using geometrical optics. Figure 3(b) shows the intensity correlation function g₂(τ)−1, averaged over 20 s. Data are normalized such that g₂(τ → 0)−1 = 1. If the sample is kept at rest during the measurement (black squares), the intensity correlation function exhibits an exponential decay, as expected for diluted Brownian suspensions [3], as better seen in the inset that shows the same data in a semilog plot. When the sample is translated at constant speed, the uncorrected g₂ − 1 decays on a much shorter time scale (blue circles) and its shape departs from a simple exponential. Clearly, no information on the microscopic dynamics can be obtained from the uncorrected data. The red crosses are the data corrected according to Eq. (21): for τ ≤ 4.2 s the corrected g₂ − 1 is very close to that measured for the stationary sample, thereby demonstrating the effectiveness of our correction scheme. For larger lags, the corrected data tend to overestimate g₂ − 1: this is consistent with the fact that the displacement cannot be reliably measured, as discussed in relation to Fig. 3(a). Indeed, in this case the peak-finding algorithm spuriously interprets the highest level in the noisy base line of corr[I(t), I(t +τ)] as the (higher-than-expected) degree of correlation.

 

Fig. 3 a): displacement versus time, measured for a diluted Brownian suspension translated at a constant speed. The line is a linear fit to the data. b) Intensity correlation functions probing the microscopic dynamics. Solid squares: quiescent sample; open circles: raw g₂ − 1 measured while translating the sample with a motor; crosses: same data, corrected for the contribution of the rigid motion of the speckle pattern.

Download Full Size | PDF

To evaluate more quantitatively how the correction performs, we analyze the correlation functions using a cumulant development, where one assumes ln[g₂[τ) − 1] ≈ Γ₀ − 2Γ₁τ + Γ₂τ² in the limit τ → 0 [26]. The first cumulant is related to the average diffusion coefficient D by Γ₁ = Dq², while the ratio ${\mathrm{\Gamma }}_2/{\mathrm{\Gamma }}_1^2$ is often referred to as to the polydispersity index (PI). For sufficiently narrow particle size distributions, the PI reduces to the variance of the size distribution normalized by the squared mean size. For the sample at rest we find D = 5.4 × 10⁻¹⁶ m²s⁻¹ and a PI of 0.156, typical of a reasonably narrow size distribution, as expected for our polystyrene particles. For the raw g₂ − 1 measured while translating the sample, we find very different values: D = 3.1 × 10⁻¹⁶ m²s⁻¹ and PI = 197, more than three orders of magnitude larger than the expected value. After correcting the data, the cumulant analysis yields D = 6.4 × 10⁻¹⁶ m²s⁻¹ and PI = 0.25, in fair agreement with the data for the unperturbed sample.

Having validated our method on a model sample whose displacement is well controlled, we test it on a more realistic experimental situation, a sheared colloidal glass for which the velocity field is not uniform over the field of view. The sample is a dense suspension of hard-sphere-like colloidal particles, a widely-studied model system for the glass transition [29]. The particles have radius ≈ 100 nm (as in [28]) and volume fraction ≈ 0.6. The setup is similar to that for the Brownian sample, but here the sample is kept in a square cell of section 10 × 10 mm², in which a glass bead of diameter D_b = 5 mm is inserted. The bead is attached to a motor that displaces it in the y (vertical) direction, parallel to the cell wall, at a speed v_y = 0.1 μms⁻¹. The minimum gap e between the wall and the bead surface is 1280 ± 90 μm. For the small displacements studied here (≤ 850 μm) and given that e ≪ D_b, the deformation is close to a simple shear. The sample is illuminated by a laser sheet in the vertical (x, y) plane, of thickness ≈ 100 μm. We image a region of size 710 × 530 μm² using light scattered at θ = 90 deg, corresponding to q = 20.6 μm⁻¹. The CCD acquisition rate and exposure time are 1 Hz and 10 msec, respectively. To obtain space-resolved information on the mesoscopic displacement and the microscopic dynamics, we run our algorithm on ten ROIs of size 40 × 340 μm² regularly spaced at a growing distance x from the moving wall. As for the previous experiment, speckle images are corrected for non-uniform illumination.

Figure 4(a) shows the velocity profiles close to the bead wall, for two times t after starting shearing the sample. In this representation, the slope of the data is the local shear rate γ̇. The solid line shows the velocity profile expected for homogeneous shear, corresponding to an average shear rate across the whole gap of γ̇ = 7.9 × 10⁻⁵ s⁻¹. It is clear that already at t = 120 s γ̇ is non-uniform across the gap, with a highly-sheared band close to the moving surface (x ≥−164 μm, γ̇ ≈ 1.4 × 10⁻⁴ s⁻¹), followed by a low shear region (x ≤−220 μm, γ̇ ≈ 6.0 × 10⁻⁵ s⁻¹). Similar shear banding has been reported for other colloidal glasses [30]. Interestingly, shear banding is seen to evolve with time. At t = 820 s, the shear rate for the high- and low-shear bands is comparable to that at t = 120 s (γ̇ ≈ 1.3 × 10⁻⁴ s⁻¹ and γ̇ ≈ 6.1 × 10⁻⁵ s⁻¹, respectively), but the boundary between the two zones has moved from x = −190 μm to x = −240 μm. Additionally, the occurrence of a marked drop of v_y close to the moving surface suggests slipping. This behavior is reminiscent of that reported for a variety of jammed or glassy soft materials, see e.g. [31], which exhibit complex spatio-temporal shear patterns. Figure 4(b) shows the intensity correlation function measured for the ROI at the position indicated by the arrow in a) [data are normalized to unity for the smallest available delay]. For the unsheared sample (crosses) no dynamics is observed on time scales up to 20s, about 2000 times the Brownian time for the same particles in the diluted regime. This is consistent with the notion that the microscopic dynamics of a sample at rest is slowed down by orders of magnitude on approaching the glass transition. The open symbols show the uncorrected g₂ − 1: a fast decay is observed, essentially due to the translation of the speckle pattern due to the imposed shear. Once corrected, the data still show a decay of g₂ − 1 (albeit a slower one), thus indicating that particles move with respect to each other, in addition to be advected by the shear deformation. We emphasize that the corrected g₂ − 1 is sensitive to the component of the particle displacement along the direction of q, which lays in the horizontal plane, perpendicular to the shear direction. Therefore, the decay of g₂ − 1 is not due to the affine component of the particle displacement along y, but rather to irreversible rearrangements associated with flow in glassy systems [32]. Interestingly, we find that the decay of g₂ − 1 is faster at t = 820 s, when both the local γ̇ and its gradient are larger. This suggest a direct relation between (local) shear rate and plastic rearrangements, as proposed for granular materials [16] and emulsions [33].

 

Fig. 4 a): velocity profiles for a sheared colloidal glass (symbols). Data are labeled by the time after initiating the shear. The dotted line indicates the position of the mobile wall, the solid line is the the velocity profile for uniform shear. The arrow indicates the location for which the data shown in b) have been measured. b) intensity correlation functions for the quiescent glass (crosses) and after applying a constant shear. Open (solid) symbols indicate the raw (corrected for drift) intensity correlation functions.

Download Full Size | PDF

5. Conclusions

We have introduced a method to obtain the mesoscopic displacement field and the microscopic dynamics in soft materials where the constituents undergo both a drift motion and a relative displacement. The algorithm proposed here is optimized for the typical features of speckle images in PCI experiments, where a small speckle size is highly desirable to maximize the spatial resolution and the statistics of the measurement. The algorithm is highly efficient in that the correction of g₂ − 1 does not requires any significant computational extra cost, besides that necessary to determine the displacement field. The method has been successfully tested on a Brownian suspension and a colloidal glass. Although similar information may be in principle obtained using confocal or optical microscopy, our method allows one to investigate samples that are difficult or impossible to visualize in real space, such as the very small particles of our colloidal glass. A generalization to speckle patterns obtained under multiple scattering conditions or for partially coherent illumination, such as in recent microscopy developments [34, 35] is also possible [36]. The method presented here should be particularly valuable for soft materials where slow dynamics is coupled to the effects of an external stress, as in rheological experiments or in samples submitted to an external field such as gravity [18], or in disordered jammed materials, where internal stress is known to play a major role [19, 37] in the sample dynamics.