Correlations between reflected and transmitted intensity patterns emerging from opaque disordered media

The propagation of monochromatic light through a scattering medium produces speckle patterns in reflection and transmission, and the apparent randomness of these patterns prevents direct imaging through thick turbid media. Yet, since elastic multiple scattering is fundamentally a linear and deterministic process, information is not lost but distributed among many degrees of freedom that can be resolved and manipulated. Here we demonstrate experimentally that the reflected and transmitted speckle patterns are correlated, even for opaque media with thickness much larger than the transport mean free path, proving that information survives the multiple scattering process and can be recovered. The existence of mutual information between the two sides of a scattering medium opens up new possibilities for the control of transmitted light without any feedback from the target side, but using only information gathered from the reflected speckle.

The propagation of monochromatic light through a scattering medium produces speckle patterns in reflection and transmission, and the apparent randomness of these patterns prevents direct imaging through thick turbid media. Yet, since elastic multiple scattering is fundamentally a linear and deterministic process, information is not lost but distributed among many degrees of freedom that can be resolved and manipulated.
Here we demonstrate experimentally that the reflected and transmitted speckle patterns are correlated, even for opaque media with thickness much larger than the transport mean free path, proving that information survives the multiple scattering process and can be recovered. The existence of mutual information between the two sides of a scattering medium opens up new possibilities for the control of transmitted light without any feedback from the target side, but using only information gathered from the reflected speckle.
In multiply scattering materials, the random inhomogeneities in the refractive index scramble the incident wavefront, mixing colors and spatial degrees of freedom, resulting in a white and opaque appearance [1]. Under illumination with coherent light and for elastic scattering, interference produces large intensity fluctuations that is not averaged out by a single realization of the disorder, resulting in a seemingly random speckle pattern [2]. In principle the speckle pattern encodes all the information on the sample and the incident light [3]. A complete knowledge of the scattering matrix allows one to reverse the multiple scattering process and to recover the initial wavefront, thus permitting imaging through turbid materials [4,5]. Conversely, if the scattering matrix is not known, a multiply scattering material effectively behaves as an opaque screen.
Speckle patterns are not as random as they appear at first sight. Interference between the possible scattering paths in the medium are known to produce spatial correlations between the intensity measured at different positions [6][7][8], and correlations of different ranges have been identified [9]. Short-range correlations determine the size of a speckle spot. Long-range correlations emerge as a consequence of constraints such as energy conservation or reciprocity [10][11][12]. Spatial correlations have not been used for imaging, a notable exception being the optical memory effect [13], a correlation of purely geometrical origin that has been exploited for non-invasive imaging through an opaque scattering layer [14,15].
At first glance, as transmitted and reflected waves are expected to undergo very different multiple scattering sequences, correlations between transmitted and reflected wavefronts are expected to quickly average to zero. Very little attention has been given to cross-correlations between transmitted and reflected speckles, their existence being only mentionned in passing [16,17]. However, a recent theoretical study suggested that a long-range correlation should survive even for thick (opaque) scattering media [18]. The existence of this reflection-transmission correlation suggests that one could non-invasively extract information on the transmitted speckle from a measurement restricted to the reflection half-space.
Here we report the first measurement of the intensity correlation between transmitted and reflected speckle patterns, for scattering materials with thickness L and scattering mean free path covering all the range from single scattering (L ) to diffusive transport (L ). The data are supported by 3D numerical simulations, and by a theoretical analysis of the lineshape of the correlation function, and its dependence on the experimental parameters. The experiments and the theory embrace the complexity and the richness of the phenomenon, thus opening the way to its use as a basic ingredient in the design of new approaches for sensing, imaging or communicating through opaque scattering media.
The experimental apparatus is shown in Fig. 1a. A monochromatic wave is incident at an angle ∼ 45 • on a suspension of TiO 2 particles in glycerol, held between two microscope slides to form a scattering slab. The slab thickness L is controlled using calibrated spacers, and the mean free path is controlled by varying the TiO 2 concentration. Typical samples with different optical thickness b = L/ , from semitransparent to fully opaque are shown in Fig. 1b. For a set of given L and we record the intensity patterns R(r) and T (r) on the surface of the sample in reflection and transmission, respectively, with two identical imaging systems (see Methods). As the samples are liquid the resulting speckle patterns change in time due to Brownian motion of the scatterers, with a correlation time τ that depends strongly on the sample thickness. Choosing an integration time < τ , and a time arXiv:1707.03622v1 [physics.optics] 12 Jul 2017 The speckle patterns on the two surfaces, T (x, y) and R(x, y) respectively, are recorded with two identical imaging systems. (b) Examples of samples with thickness L = 20µm but different TiO2 concentrations: from left to right 5 g/l, 10 g/l and 40 g/l, which correspond to a mean free path of (60, 20.4 and 9.8) ± 2.5 µm, respectively.
interval between successive measurements > τ , allows us to measure speckle images R(r) and T (r) for a large ensemble of configurations of the disordered medium. An example pair of images measured for a given realization of disorder is shown in Fig. 2a,b. For each pair of R and T we calculate the correlation function C RT defined as (1) where ∆r = (∆x, ∆y) is a transverse shift between the images, and the overline denotes the spatial average over the coordinates r. We have introduced the notation δf = f − f for the statistical fluctuations of a random variable f , with · representing the ensemble average over disorder. Plotted as a 2D map the correlation function C RT (∆x, ∆y) appears random, with a granularity similar to that of a speckle image (Fig. 2c). After ensemble averaging over the realizations of the disorder a clear pattern emerges in C RT (∆x, ∆y) (Fig. 2d), showing that the transmitted and reflected speckle patterns are indeed correlated.
Speckle correlations are commonly divided into three categories: Short-range correlations (C 1 ) that decay with the separation between the observation points on the scale of the wavelength, long-range correlations (C 2 ) that have a polynomial decay, and infinite-range correlations (C 3 ) [8,9]. The short-range correlation C 1 corresponds to the approximation of a field obeying Gaussian statistics [2], while C 2 and C 3 are non-Gaussian corrections. An additionnal infinite-range correlation (C 0 ) appears under illumination by a point source located inside the medium [19]. One can see in Fig. 2d that the lineshape of C RT (∆r) is much wider than a speckle spot, indicating that the dominant contribution to this correlation is long-range in nature.
In order to characterize the lineshape of the correlation function, and to probe its dependence on the sample parameters, we measured C RT (∆r) for different values of and L, covering the full range from the single scattering (L ) to the diffusive (L ) regime. The results are summarized in Fig. 3 (center and right columns), where both 2D maps C RT (∆x, ∆y) and cross-sections along the line ∆y = 0 (indicated as a dotted line in the 2D maps) are displayed. It is interesting to note that both the shape and the sign of the reflection/transmission correlation substantially depend on L and . In the single scattering regime (optical thickness b 1), C RT is dominated by a narrow peak (still much larger than a single speckle spot) with a negative side lobe. In the multiple scattering regime (b 1), C RT is dominated by a wide negative dip.
The short-range contributions to C RT (C 1 ) decay on the scale of the wavelength [2], and are thus negligible in all measurements since in the reflection-transmission geometry the observation points are separated by a distance L 2 + |∆r| 2 which is much larger than λ for even the thinniest sample (see SI section 2). Hence, C RT is necessarily a long-range correlation of the C 2 type. It is interesting to note that the reflection/transmission geometry naturally favors the observation of non-Gaussian long-range correaltions, without requiring any post-processing to remove the C 1 contribution that dominates in the pure transmission geometry [10,11,20,21]. Another feature of our experiment is the illumination/detection geometry that excludes any contribution from specularly reflected and transmitted averaged fields (see Methods). Indeed, in the geometry in Fig. 1a, the averaged field is not collected by the detectors, and we directly correlate T (r) = |δE T (r)| 2 and R(r) = |δE R (r)| 2 , thus avoiding spurious interference terms in the correlation function for the single scattering regime (see SI section 3).
To support the experimental data, we have performed full numerical simulations of wave propagation in threedimensional disordered media. In the simulations, the samples consist of slabs of dipole scatterers with random positions. The scalar wave equation is solved numerically using the coupled-dipole method (see Methods). From the computation of the complex amplitude of the scattered field in reflection and transmission, we deduce C RT performing the ensemble average over a large set of samples. The results of the simulations are displayed in Figs. 3 (left column), and are in very good agreement with the experimental data. The general shape of the correlation in the regimes b < 1, b 1 and b > 1 is well reproduced in the simulations.
In order to refine the analysis, and to get more physical insight, we have also used a formal pertubation theory to calculate the shape of the correlation function C RT . In this formal approach, the correlation function C RT (∆r) = δR(r)δT (r + ∆r) / R(r) T (r) is directly computed from a statistical ensemble averaging, without going through the intermediate spatial average in Eq. (1) used for the experimental data. Both averaging processes coincide provided that λ, a condition that is always satisfied in our experiments (see Methods and SI section 1). Formal pertubation theory uses 1/(k ) as a small parameter, with k = 2π/λ, and relies on a diagrammatic formalism that allows one to derive explicit expressions of intensity correlation functions [6][7][8]. In the reflection/transmission geometry, care must be taken to properly account for leading contributions [16,22].
Let us first discuss the regime of large optical thickness L corresponding to Fig. 3f-i. Strikingly, we observe in this regime that C RT (∆r) is negative, in agreement with the prediction in Ref. [18]. This means that for every bright spot in reflection (transmission) the corresponding area in transmission (reflection) is more likely to be darker, and vice versa. This feature can be inferred from flux conservation arguments. Indeed, defining T ∝ T (r)dr and R ∝ R(r)dr, energy conservation imposes T + R = 1 for a non-absorbing medium, from which we can show that C RT (∆r) d∆r ∝ δT δR = − δT 2 < 0 (see SI section 4). Note that the existence of negative long-range C 2 correlations has been previously pointed out in Refs. [17,23]. Refining the analysis performed in Ref. [18] (see SI section 5), we have studied theoretically the regime L λ, and found that both the amplitude and the width of the correlation function depend on L and , as in the experimental data in Figs. 3f-i. For L , the dominant diagrams belong to the class represented in Fig. 4a, that are typical of long-range C 2 correlation functions. They predict a correlation function that is isotropic, independent of the angle of incidence, and scales as C RT (∆r) = C RT 2 (∆r) = −f (|∆r|)/(kL) 2 , where f is a dimensionless function that decays on a range |∆r| L [18]. This long-range character of the correlation function originates from the crossing of two diffusive paths that probe a transverse distance L, as represented in Fig. 4a. Moreover, the correlation function in this regime is independent of the disorder strength k , which makes it strikingly different from that observed in a pure transmission geometry, for which C T T 2 ∼ 1/[(k )(kL)] ∝ 1/g, where g is the dimensionless conductance of the sample [24]. Another important difference between C RT 2 and C T T 2 is the evolution of their information content with respect to the detection scheme. Although C T T 2 contains the same information whether it is measured on the sample surface or in the far field, this is not the case for C RT 2 . Indeed, in the far field, we have C RT (k b , k b ) ∼ C RT (∆r) d∆r = const. for any pair of observation directions k b , k b , as the information content is spread uniformly over all degrees of freedom. Therefore we will focus our discussion on the correlations measured on the sample's surface.
In the regime of moderate optical thickness ∼ L λ, where single scattering is expected to dominate, an intensity correlation extending far beyond the size of a single speckle spot is still observed (see Fig. 3c-e), but with a positive peak appearing in the vicinity of the negative contribution. The apparent relative position and amplitude between the peak and the dip depends on the angle of incidence of the illumination (see SI section 6). Contrary to the negative dip in the correlation function observed at large optical thickness, the lineshape is anisotropic, with negative side lobes (hardly visible in Fig. 3c-e but visible in a calculation shown in SI section 6) that are more pronounced along the direction of the incidence plane. Moreover, the amplitudes of both the positive peak and the side lobes substantially depend on the incidence angle. These two features of the correlation function in the regime ∼ L λ (long-range extent and dependence on the angle of incidence) suggest a qualitative description based on diagrams of the class represented in Fig. 4b, that satisfy both properties simultaneously (see SI section 7). Computing these diagrams lead to a contribution to the correlation function scaling as 1/(kL) 4 for b 1. This is consistent with the fact that, according to the measurements and the numerical simulations, this contribution has to be negligible at large optical thickness, where the C RT 2 correlation function discussed previously scaling as 1/(kL) 2 dominates. The observed anisotropy in the correlation function is also well reproduced, supporting the relevance of the analysis based on the diagrams in Fig. 4b. Interestingly, these diagrams are formally similar to those leading to the infinite-range correlations C 0 observed when the sample is excited by a point source [19,25]. In a non-absorbing medium, as a consequence of energy conservation, the C 0 contribution is related to the fluctuations of the local density of states at the source position [26,27]. In the present context, where a plane wave excitation is used, the C 0 -type contribution to the reflection-transmission correlation function is long-range and satisfies C RT 0 (∆r) d∆r = 0 (see SI section 7). This property also leads to the conclusion that the C 0type contribution observed is specific to speckle patterns measured on the surface of the sample, and vanishes in the case of far-field angular measurements.
Note finally that in the quasi-ballistic regime L λ, which is not the focus of our experiment, we expect the correlation C RT to contain additional contributions to C 2 and C 0 (see SI section 7), that still result in an overall positive peak.
In summary, we have demonstrated experimentally the existence of a cross-correlation between the speckle patterns measured in reflection and transmission on the surface of a disordered medium. The correlation persists in the regime of large optical thickness L in which the sample is opaque due to multiple scattering. The measurements are supported by 3D numerical simulations, and have been analysed using a perturbative theory (valid when λ). We have found that the reflectiontransmission correlation has two contributions: a positive peak dominant at moderate optical thicknesses L , and a negative dip dominant in the diffusive regime L . The existence of this transmission-reflection correlation proves that mutual information between the two sides of a strongly scattering medium can in principle be exploited, opening new possibilities for the non-invasive control of the transmitted light (or other kind of coherent waves) from measurements restricted to the reflection half-space. This offers new strategies for the detection of objects hidden behind opaque scattering media, including ghost imaging schemes, and the control of wave propagation by wavefront shaping techniques [28,29].

METHODS
Experimental set-up. The experiments are performed using a 2 mW He-Ne laser (632.8 nm model HNLS008L-EC, Thorlabs) incident on the sample at 45 • . The imaging system sketched in Fig. 1a consists of two identical microscope objectives (10X Olympus Plan Achromat Objective, 0.25 NA) and two plano-convex 150 mm lenses. The intensity speckle patterns on the two surfaces of the sample are imaged with two identical CCD cameras (Allied Vision Manta G-146) (see SI section 9). The integration time of the cameras is set to 1 ms, chosen to be shorter than the decorrelation time of the speckle patterns due to Brownian motion of the scatterers. The samples are a mixture of TiO 2 particles and glycerol. We prepared the samples with three different TiO 2 concentrations: 50, 150 and 400 mg of TiO 2 in 10 ml of glycerol, which resulted in a mean free path of (60, 20.4 and 9.8) ± 2.5 µm, respectively (see SI section 8). The thickness of each sample is fixed using two feeler gauges of the required thickness. To insure a proper mixing and absence of clustered grains, the TiO 2 powder is mixed with glyc-erol using a magnetic stirrer for two hours, and after that was sonicated for 30 minutes.
Definition of the correlation function. The definition of the correlation function used to analyze the experimental data is given by Eq. (1), while in the theoretical analysis we used C RT (∆r) = δR(r)δT (r + ∆r) / R(r) T (r) . (2) The two definitions appear different, but they are rigorously equivalent under two conditions: ergodicity and Gaussian fields.
Ergodicity allows us to transform Eq. (1) into C RT (∆r) = δR(r)δT (r + ∆r) / δR(r) 2 δT (r) 2 1/2 , which is identical to Eq. (2) when δR(r) 2 = R(r) 2 and δT (r) 2 = T (r) 2 , i.e. when R(r) and T (r) obey the Rayleigh statistics (Gaussian fields). For systems close to the Anderson localization (low dimensionless conductance g) deviations from the Rayleigh distribution are expected [11,21], but for the diffusive media under examination here the two formulas are identical to a very high accuracy (see SI section 1 for a detailed analysis).
Numerical simulations. The 3D numerical simulations are performed using the coupled-dipole method [30]. Scatterers represented by electric point dipoles are randomy distributed in a rectangular box, with longitudinal thickness L and transverse size ten times larger to mimic the slab geometry used in the experiments. Since the measurements are not resolved in polarization, and since the input light is expected to depolarize on a length scale on the order of [31], we neglect polarization and numerically solve the scalar wave equation. To limit the number of scatterers and save computational time, the polarizability α of each scatterer has been chosen to maximize the scattering cross-section σ s = k 4 |α| 2 /(4π) leading to α = 4iπ/k 3 . Adjusting the number density of scatterers ρ, we can vary the scattering mean-free path = 1/(ρσ s ) and simulate different kinds of samples. Solving numerically the coupled-dipole equations, we compute the scattered field at any point on the input and exit surfaces of the slab, and deduce the correlation function C RT (∆r) . The ensemble averaging is performed by computing the field for many realizations of the positions of the scatterers. As an example, in the regime k = 10 and b = 1.5, we have used N = 2685 dipoles and 26 millions of configurations.