RESOLUTION ENHANCEMENT FROM SCATTERING IN PASSIVE SENSOR IMAGING WITH CROSS CORRELATIONS

It was shown in [Garnier et al., SIAM J. Imaging Sciences 2 (2009), 
396] that it is possible to image re ectors by backpropagating cross correlations 
of signals generated by ambient noise sources and recorded at passive sensor 
arrays. The resolution of the image depends on the directional diversity of 
the noise signals relative to the locations of the sensor array and the reflector. 
When directional diversity is limited it is possible to enhance it by exploiting 
the scattering properties of the medium since scatterers will act as secondary 
noise sources. However, scattering increases the fluctuation level of the cross 
correlations and therefore tends to destabilize the image by reducing its signal- 
to-noise ratio. In this paper we study the trade-o in passive, correlation-based 
imaging between resolution enhancement and signal-to-noise ratio reduction 
that is due to scattering.

1. Introduction. It has been shown recently that the Green's function of the wave equation in an inhomogeneous medium can be estimated by cross correlating signals emitted by ambient noise sources and recorded by a passive sensor array [1,7,22,23]. The cross correlation C(τ, x 1 , x 2 ) of the signals recorded at two sensors x 1 and x 2 : is related to the Green's function between x 1 and x 2 . In a homogeneous medium and when the source of the waves is a space-time stationary random field that is also delta correlated in space and time, it has been shown [19,16] that the τ -derivative of the cross correlation of the recorded signals is proportional to the symmetrized Green's function between the sensors: where G is the Green's function. In an inhomogeneous medium and when the sources completely surround the region of the sensors the approximate identity (1.1) is still valid and it can be shown using the Helmholtz-Kirchhoff theorem [22,9]. This is true even with spatially localized noise source distributions provided the waves propagate within an ergodic cavity [1]. More generally, the cross correlation as a function of the lag time τ can have a peak at plus or minus the inter-sensor travel time T (x 1 , x 2 ), provided the ambient noise sources are well distributed around the sensors. The inter-sensor travel times obtained from peaks of cross correlations can then be used tomographically for background velocity estimation [4,14,18,24]. However, when the ambient noise sources have spatially limited support, the signals recorded by the sensors are generated by wave energy flux coming from the direction of the noise sources, which results in an azimuthal dependence of the quality of travel time estimation [21]. In general, travel time estimation by cross correlation of noise signals is possible when the line between the sensors is along the direction of the energy flux and difficult or impossible when it is perpendicular to it. This can be explained using a stationary phase analysis [9,13]. It is, however, possible to enhance the quality of travel time estimates by exploiting the enhanced directional diversity provided by the scattering of waves in a randomly scattering medium [9,20]. We consider this briefly in Section 6.
Cross correlations of signals emitted by ambient noise sources and recorded by a passive sensor array can also be used for imaging of reflectors [14,9,10]. The data to be used for imaging is the matrix of cross correlations (C(τ, x j , x l )) j,l=1,...,N between pairs of sensors (x j ) j=1,...,N . The objective is to image reflectors in the medium. In this paper: -We consider imaging functionals that use the cross correlation matrix of the signals recorded by a passive sensor array and emitted by ambient noise sources. The form of these functionals is motivated by the analysis of the cross correlation matrix in asymptotic regime in which the coherence time of the sources is small compared to typical travel times. The imaging functionals can then be applied in a general context. The analysis of Section 4, in particular, shows the limitations of the imaging functionals proposed in [9,10], which were derived for a homogeneous background. We introduce in Section 6 some new imaging functionals that work well in randomly scattering medium but we do analyze them here.
-We analyze the resolution and the signal-to-noise ratio (SNR) of the imaging functionals. We consider simple models for the reflector and the background medium that allow a rather complete analysis and explain the main phenomena that we want to study. In particular we show that resolution and SNR are not independent notions in a randomly scattering medium. In the analysis weuse the Born approximation for the reflectors and a single-and double-scattering approximations for the random medium. We also use extensively the stationary phase method.
-We carry out numerical simulations to illustrate the results. This confirms that the theoretical predictions obtained in asymptotic regimes can be observed in realistic situations.
Homogeneous background medium. Here we summarize the main findings of [9,10] when the background medium is homogeneous. When there is only one reflector located at z r , the cross correlation C(τ, x j , x l ) between the two sensors at x j and x l may have additional peaks at lag times different from the inter-sensor travel times. Stationary phase analysis [9] shows that these additional peaks appear at lag times equal to the sum or the difference between the travel times from the sensors to the reflector. This suggests that the reflector can be imaged by backpropagating or migrating the cross correlation matrix of the recorded signals. This form of passive sensor array imaging depends in an essential way on the illumination configuration, that is, the relative position of the sensors, the noise sources and the reflector as follows: -If the noise sources are spatially localized and the sensors are between the sources and the reflector, then additional peaks in the cross correlations occur at lag times equal to plus or minus the sum of travel times T (x l , z r ) + T (x j , z r ). We call this the daylight configuration (see Figure 1.1a).
-If the noise sources are spatially localized and the reflector is between the sources and the sensors, then an additional peak in the cross correlations occurs at lag time equal to the difference of travel times T (x l , z r ) − T (x j , z r ). We call this the backlight configuration (see Figure 1.1b). The imaging functional that we propose therefore depends on the illumination configuration. To form an image in a daylight illumination configuration, each element of the cross correlation matrix is evaluated at lag time equal to the sum of the travel times T (x l , z S ) + T (z S , x j ) between the sensor x l and a search point z S in the image domain, and between the search point z S and the sensor x j . The Kirchhoff migration imaging functional is the sum of the migrated matrix elements over all pairs of sensors, as discussed in Subsection 3.3. To form an image in a backlight illumination configuration, each element of the cross correlation matrix is evaluated at lag time equal to the difference of the travel times T (x l , z S ) − T (z S , x j ) and the migrated matrix elements are summed over all pairs of sensors. The backlight imaging functional has poor range resolution because it exploits only differences of travel times, which are much less sensitive to the range than the sums of travel times [10].
Randomly scattering medium. It is known that scattering by random inhomogeneities in the medium can enhance the directional diversity of waves [17]. This was studied in the context of time-reversal experiments and it has been shown that time-reversal refocusing can be enhanced in a randomly scattering medium [8,15]. This was also considered in the context of travel time estimation by cross correlation of noisy signals [9,20]. In this paper we show how wave scattering affects both the resolution and the SNR of passive sensor imaging when using cross correlations of noise signals. We show in Section 4 that the scattering medium enhances the directional diversity of waves, which in turn improves the resolution of the imaging functional. We show also that the random fluctuations in the cross correlations due to the scattering reduce the SNR of the imaging functional. These are the two competing phenomena that we analyze here in detail in a regime of weak scattering.
We consider passive sensor imaging with cross correlations in configurations such as the one shown in Figure 1.1c. Now the ambient noise sources provide only backlight illumination but the scatterers in the medium provide a secondary daylight illumination. The question we address here is how to exploit this secondary daylight illumination due to scattering. We compute both the heights and widths of the peaks of the cross correlation at special lag times associated with the peaks influencing the migration functionals. We also calculate the standard deviations of the fluctuations of the cross correlations due to scattered waves. The main results of the SNR analysis of the peak in the cross correlations generated by the secondary illumination from the scatterers are summarized in Subsection 4.6. When the scattering region is far enough from the sensor array then migrating directly the cross correlations works well, as the theory of Section 4 predicts and as shown in Figure 4.3. However, in unfavorable scattering situations the peaks generated by the secondary illumination are weak compared to the random fluctuations of the cross correlations. An example is given in Section 6 where the scattering region is too close to the array. As a consequence migration of the cross correlation matrix gives the blurred and speckled image shown in Figure 6.1. This is a typical example where the potential resolution enhancement is not achieved because of low signal-to-noise ratio. In order to enhance the signal-tonoise ratio we time window the cross correlations, select the tails (or coda) and cross correlate them. By migrating this special fourth-order cross correlation matrix it is possible to get a much better image than by migrating directly the cross correlation matrix. This is what we show in Section 6.
When scattering is very strong then imaging with cross correlations cannot be done with migration or related coherent imaging methods. Fluctuations in the cross correlations from wave scattering are then large and statistical stability becomes a critical issue [12].
The paper is organized as follows. In Section 2 we introduce the framework for the analysis of cross correlation of noise signals. In Section 3 imaging by cross correlation of noise signals is presented when the medium is homogeneous or has a slowly varying background. In Section 4 we address the case when the medium contains random inhomogeneities responsible for scattering. In Section 5 we consider the case in which scattering is not generated by random inhomogeneities but by a reflecting interface in the medium. In Section 6 we show how the use of iterated cross correlations and appropriate imaging functionals can improve the resolution and SNR of passive sensor imaging in a scattering medium.

2.
The cross correlation of ambient noise signals.
2.1. The wave equation with noise sources. We consider the solution u of the wave equation in a three-dimensional inhomogeneous medium with propagation speed c(x): The term n(t, x) models a random field of noise sources. It is a zero-mean stationary (in time) random process with autocorrelation function Here · stands for statistical average with respect to the distribution of the noise sources. For simplicity we will consider that the process n has Gaussian statistics. The time distribution of the noise sources is characterized by the correlation function F (t 2 − t 1 ), which is a function of t 2 − t 1 only by time stationarity. The function F is normalized so that F (0) = 1. The Fourier transformF (ω) of the time correlation function F (t) is a nonnegative, even, real-valued function proportional to the power spectral density of the sources: (2. 3) The spatial distribution of the noise sources is characterized by the autocovariance function δ(y 1 − y 2 )K(y 1 ). The process n is delta-correlated in space and K characterizes the spatial support of the sources.
2.2. Statistical stability of the cross correlation function. The stationary solution of the wave equation has the integral representation where G(t, x, y) is the time-dependent outgoing Green's function. It is the fundamental solution of the wave equation starting from G(0, x, y) = ∂ t G(0, x, y) = 0 (and continued on the negative time axis by G(t, x, y) = 0, ∀t ≤ 0). The empirical cross correlation of the signals recorded at x 1 and x 2 for an integration time T is It is a statistically stable quantity, in the sense that for a large integration time T , the empirical cross correlation C T is independent of the realization of the noise sources. This is stated in the following proposition proved in [9]. Proposition 2.1. 1. The expectation of the empirical cross correlation C T (with respect to the distribution of the sources) is independent of T : where the statistical cross correlation C (1) is given by andĜ(ω, x, y) is the time-harmonic Green's function (i.e. the Fourier transform of G(t, x, y).
2. The empirical cross correlation C T is a self-averaging quantity: (2.9) in probability with respect to the distribution of the sources.

Small decoherence time hypothesis.
We assume from now on that the decoherence time of the noise sources is much smaller than the typical travel time that we want to estimate, that is, the travel time from the reflector(s) to the sensor array. If we denote by ε the (small) ratio of these two time scales, then we can write the time correlation function F ε of the noise sources in the form where t 1 and t 2 are scaled relative to typical travel times. The hypothesis ε 1 is both natural and useful: 1) In experimental situations (surface wave tomography) noise records are first bandpass-filtered and then cross correlated [18]. If the central frequency ω 0 of the filter is high enough so that the corresponding wavelength λ 0 is much smaller than the distance d from the sensor array to the reflector, then we have ε = λ 0 /d 1. As we will see below, the range resolution of migration imaging by cross correlation is inversely proportional to the bandwidth, so the hypothesis ε 1 turns out to be natural in order to get some resolution.
2) The Fourier transform of the time correlation function of the sources has the formF ε (ω) = εF (εω), so that the statistical cross correlation (2.8) involves a product of Green's functions evaluated at high frequencies: This expression shows that a geometric asymptotics approach and a stationary phase analysis are appropriate tools to study the cross correlations in the regime ε 1.
3. Passive sensor imaging with cross correlations in a homogeneous background medium. In this section we revisit the results of [9,10] and show that it is possible to image reflectors by cross correlations of signals generated by ambient noise sources and recorded by passive sensors.
We want to image reflectors in the medium from signals recorded by the array of sensors located at (x j ) j=1,...,N . In active array imaging the sensors can be used as emitters as well as receivers. The data set collected is the impulse response matrix (P (t, x j , x l )) j,l=1,...,N,t∈R , where the (j, l)-entry of this matrix is the signal (P (t, x j , x l )) t∈R recorded by the jth sensor when the lth sensor emits a Dirac impulse. The usual migration techniques [2,6] that backpropagate the impulse responses numerically in a fictitious medium produce images of the reflectors. However, in passive sensor imaging the sensors of the array do not have emission capacity and they can only act as receivers. The data collected in passive array imaging are the signals (u(t, x j )) t∈R generated by ambient noise sources and recorded by the jth sensor j = 1, . . . , N . It turns out that the matrix of cross correlations of the recorded signals (C(τ, x j , x l )) j,l=1,...,N,τ ∈R (with C defined as in (2.6)) behaves like the impulse response matrix of an active sensor array. It is therefore possible to image the reflectors by backpropagating the cross correlations.
3.1. Differential cross correlations. In order to image reflectors with cross correlations it is often necessary to have available data sets {C(τ, x j , x l )} j,l=1,...,N and {C 0 (τ, x j , x l )} j,l=1,...,N , with and without the reflectors, respectively, so that we can compute the differential cross correlations {C − C 0 } and migrate them. In general, the primary data set {C} cannot be used directly for imaging because peaks in the cross correlations due to the reflectors may be very weak compared both to the peaks of the direct waves, at lag times equal to the inter-sensor travel times, as well as to the non-singular components due to the directionality of the energy flux [9]. This is so both in homogeneous and in scattering media. By singular component of a cross correlation we mean the leading contribution for ε small coming from the stationary phase approximation. By non-singular component we mean the remainder after the stationary phase approximation.
In many applications where we want to image localized changes in the environment, such as in reservoir or volcano monitoring, both data sets are usually available.
When the peaks in the cross correlations due to the reflectors to be imaged are well separated from peaks at inter-sensor travel times and from effects from noise energy flux directionality, we may be able to migrate directly the cross correlations of the primary data set {C(τ, x j , x l )} j,l=1,...,N . This issue is discussed in detail in Section 6 in [9].
3.2. Stationary phase analysis. In this section we carry out the analysis when the background medium is homogeneous with background speed c 0 and there is a point reflector at z r . Since we assume that the reflector is weak and small, we can use the point interaction approximation for the Green's function [10]: HereĜ 0 is the Green's function of the background medium, that is, in the absence of reflector,Ĝ T (x, y) is the travel time from x to y, σ r is the reflectivity of the point reflector, and l 3 r is the effective scattering volume. The statistical differential cross correlation is given by where C (1) 0,r is the statistical cross correlation in the presence of the reflector, that is, equation (2.11) with the full Green's function (3.1). The cross correlation C (1) 0 is obtained in the absence of the reflector, that is, equation (2.11) with the background Green's function (3.2). We collect the terms with the same power in σ r l 3 r . The terms of order O(1) cancel and we retain only the terms of order O(σ r l 3 r ) consistently with the Born or lowest order scattering approximation: In the differential cross correlations the contributions of the direct waves are removed so that the small singular components of the reflected waves can be observed. The next proposition 3.1 was proved in [10]. It is a quantitative form of the analysis first carried out in [9] where only the stationary paths were identified, while Proposition 3.1 requires a full stationary phase analysis. By identifying the singular components of the differential cross correlations we determine the appropriate imaging functional that should be used to migrate the cross correlations as we discuss in Subsection 3.3. This proposition is also the first step in the SNR analysis that we carry out in Section 4.  The proof of the proposition uses stationary phase analysis for (3.4) in order to show that the singular contributions of the differential cross correlation correspond to (stationary) pairs of ray segments joining a source point y, the reflector z r , and the sensors x 1 and x 2 . The proposition shows that there are two main types of configurations of sources, sensors, and reflectors. 1) The noise sources are spatially localized and the sensors are between the sources and the reflectors (cases (a) and (b) in Figure 3.1). We call this the daylight illumination configuration. In this configuration the singular components of the differential cross correlation are concentrated at lag times equal to (plus or minus) the sum of travel times T (x 2 , z r ) + T (x 1 , z r ).
2) The noise sources are spatially localized and the reflectors are between the sources and the sensors (cases (c) and (d) in Figure 3.1). We call this the backlight illumination configuration. In this configuration the singular components of the differential cross correlation are concentrated at lag times equal to the difference of travel times Proposition 3.1. In the backlight imaging configuration, in the regime ε → 0, the differential cross correlation ∆C (1) 0 has a unique singular contribution at lag time equal to the difference of travel times T (x 2 , z r ) − T (x 1 , z r ) and it has the form: where K z,x is the energy released by the noise sources along the ray joining x and z: In the daylight illumination configuration, the differential cross correlation ∆C (1) 0 has two singular contributions at lag times equal to plus or minus the sum of travel times T (x 2 , z r ) + T (x 1 , z r ). The peak centered at plus the sum of travel times has the form: The peak centered at minus the sum of travel times has the form: (3.8) Note that K zr,xj is not zero only if the ray going from x j to z r extends into the source region, which is the backlight illumination configuration, while K xj ,zr is not zero only if the ray going from z r to x j extends into the source region, which is the daylight illumination configuration. Proposition 3.1 shows that the heights of the peaks are inversely proportional to the square of the distance from the reflector to the sensor array, that they do not depend on the distance from the sources to the reflector or the sensor array, and that they are proportional to the thickness of the source region. The second property was already noted in [10]. It is a consequence of two competing phenomena that cancel each other: on the one hand the geometric decay of the amplitude of the Green's function as a function of the distance from the sources to the reflector, and on the other hand the decay of the Hessian determinant of the rapid phase at the stationary points, which by the stationary phase theorem, produces a multiplicative factor in the evaluation of the integral (3.4).
More quantitatively, if we denote by d(A, R) the distance from the source array to the reflector, by W K the width of the source region, and by K 0 the typical value of K, then we find that the typical height of the singular peaks is If we denote by ω 0 the central frequency of the sources and by B their bandwidth, then and the time width of the peak is B −1 .

Migration imaging of cross correlations.
We consider first migration imaging with daylight illumination. The imaging functional at a search point z S is the daylight imaging functional defined by It is a consequence of Proposition 3.1 that the migration functional should be evaluated at lag time equal to (plus or minus) the sum of travel times Since the cross correlation satisfies C(−τ, x j , x l ) = C(τ, x l , x j ), we obtain the form of the functional (3.11) The resolution analysis of the daylight imaging functional is carried out in [10]. The cross range resolution of the daylight imaging functional for a linear sensor array with aperture a is given by λ 0 a/d(A, R). Here d(A, R) is the distance between the sensor array and the reflector and λ 0 is the central wavelength. The range resolution for broadband noise sources is equal to c 0 B −1 where c 0 is the background velocity and B is the bandwidth. The range resolution for narrowband noise sources is λ 0 a 2 /d(A, R) 2 . We consider next migration imaging with backlight illumination. The imaging functional at a search point z S is the backlight imaging functional defined by (3.12) The sign of the travel time in the argument of the imaging functional is determined by Proposition 3.1. It is shown there that the singular component of ∆C(τ, The resolution analysis of the backlight imaging functional is carried out in [10]. The cross range resolution is λ 0 a/d(A, R) while the range resolution is λ 0 a 2 /d(A, R) 2 whatever the bandwidth, which means that the range resolution is very poor compared to the daylight imaging functional.

Numerical simulations.
In the numerical simulations presented here we evaluate the statistical cross correlations C for different configurations of noise sources, reflector, and sensors. The statistical cross correlation is what is obtained with the empirical cross correlation C T for an infinitely large integration time T . The statistical stability (i.e. the fluctuations of C T with respect to its statistical average for finite T ) has been studied in detail theoretically and numerically in [9]. It is not a limiting factor in this type of problems as long as the recording time window can be taken arbitrarily large.
We consider a three-dimensional homogeneous background medium with velocity c 0 = 1. We compute the image in the plane (x, z) and use the homogeneous background Green's function (3.2). The random sources are a collection of 100 randomly located point sources in a layer of size 100 × 15 with power spectral densitŷ F (ω) = exp(−ω 2 ). We consider a point reflector at position (−5, 60) with reflectivity σ r l 3 r = 0.01 and 5 sensors located at (−37.5 + 7.5j, 100), j = 1, . . . , 5.
In Figure 3.2 we consider a backlight illumination configuration. We apply both the backlight imaging functional (3.12) and the daylight imaging functional (3.11). As predicted by the theory, the backlight imaging functional has good cross-range resolution but very poor range resolution. The daylight imaging functional is not appropriate.
In Figure 3.3 we consider a daylight illumination configuration. We apply both the backlight imaging functional (3.12) and the daylight imaging functional (3.11). As predicted by the theory, the daylight imaging functional has good cross-range and range resolutions. The backlight imaging functional is not appropriate.
4. Passive sensor imaging in a randomly scattering medium. In this section we consider imaging by cross correlations, as in the previous section, but in addition to the homogeneous background there are random inhomogeneities that cause scattering. In the first three subsections we formulate the multiple scattering problem and introduce the setup for the analysis of the differential cross correlations. In the following subsections we focus our attention on the backlight illumination configuration with scatterers behind the sensor array (see Figure 4.3a). This is an interesting configuration because the sources provide a backlight illumination and only the backlight imaging functional gives an image of the reflector in a homogeneous medium. In a scattering medium we anticipate that the scatterers have the advantage that they can enhance the directional diversity and provide a secondary daylight illumination, but they also have the drawback that they introduce fluctuations in the cross correlations. These are the phenomena that we want to analyze in detail. Subsections 4.4-4.6 give a quantitative analysis of the first two moments (mean and variance) of the differential cross correlation at appropriate lag times. Subsection 4.7 summarizes this quantitative analysis and applications to migration imaging are discussed in Subsection 4.8.

4.1.
A model for the scattering medium. In order to analyze the cross correlation technique in a scattering medium, we first introduce a model for the inhomogeneous medium. We assume that the speed of propagation of the medium has a homogeneous background speed value c 0 and small and weak fluctuations responsible for scattering: where V (x) is a random process with mean zero and covariance function of the form Here E stands for the expectation with respect to the distribution of the randomly scattering medium, σ s is the standard deviation of the scattering amplitude, l s is the correlation length of the fluctuations of the speed of propagation, and the function ρ(x) characterizes the spatial support of the scatterers. Note that we have assumed that the correlation length is small enough so that we can consider that the process V is delta-correlated. The cluttered Green's functionĜ clu , that is to say, the Green's function of the medium with clutter in the absence of reflector, is the fundamental solution of with c clu (x) given by (4.1) (and with radiation conditions). We now assume that a reflector is embedded at z r in the cluttered medium. We model the reflector by a local variation V r (x) of the speed of propagation: Here Ω r is a compactly supported domain with volume l 3 r . The full Green's function G clu,r , that is to say, the Green's function of the medium with clutter in the presence of the reflector at z r , is solution of It can be written as the sum where the cluttered Green's functionĜ clu is solution of (4.2). The correctionĜ cor is solution of the Helmholtz equation with the velocity c clu and the source term so that it can be represented aŝ This expression is exact and it is called the Lipmann-Schwinger equation. The Born approximation (or single-scattering approximation) consists in replacingĜ clu,r on the right side of (4.3) byĜ clu , which gives [3, section 13.1.2]: This approximation is valid if the correctionĜ cor is small compared toĜ clu , i.e., in the regime in which σ r 1. We also assume that the diameter l r of the scattering region Ω r is small compared to typical wavelength. We can then model the reflector by a point reflector: , and we can write the correction in the form and the full Green's function iŝ 4.2. The differential cross correlation. We consider the statistical differential cross correlation, that is the difference between the statistical cross correlations in the presence and in the absence of reflector. It is given by whereĜ clu is the cluttered Green's function (4.9) in the absence of the reflector and G clu,r is the full Green's function (4.4) in the presence of the reflector. We substitute the approximation (4.4) into (4.5-4.6). Consistent with the Born approximation, we neglect the terms of order O(σ 2 r ) in ∆C (1) (τ, x 1 , x 2 ) and we find When the medium is homogeneous we recover the expression (3.4). When the medium is cluttered the expression (4.8) involves the cluttered Green's function solution of (4.2) that can have a complicated structure. We simplify the expression of the cluttered Green's function in the next subsection in order to get a model that is mathematically tractable, and complex enough to explain the main phenomena we want to study, in terms of resolution enhancement and SNR reduction.

Expansion of the cluttered Green's function.
We use the so-called second Born approximation for the cluttered Green's function solution of (4.2). This approximation takes into account single scattering and double scattering events for the interaction of the waves with the cluttered medium: whereĜ 0 is the Green's function (3.2) of the homogenous background medium and G 1 andĜ 2 are given bŷ We take into account single and double scattering events because we will collect and study all terms of order O(σ r ), O(σ r σ s ), and O(σ r σ 2 s ) in the differential cross correlation, while we neglect all terms of order O(σ r σ 3 s ) and O(σ 2 r σ s ). This means that we consider a regime in which σ r σ s 1, in words, scattering due to the random medium is weak, and the scattering amplitude of the reflector is even weaker.
That is why we have first expanded the differential cross correlation with respect to σ r , and then with respect to σ s .
In fact, we will see thatĜ 2 does not contribute to the differential cross correlation to this order, but this is not clear in advance. It is a consequence of the highfrequency (stationary phase) analysis for the particular geometric configuration of the scattering region, the array and the reflector (see Figure 4.3a) that only the cross correlation terms that involve two (single-scattering) componentsĜ 1 give a significant contribution to the mean and the variance of the overall cross correlation, while the cross correlation terms that involve one (double-scattering) componentĜ 2 do not make any significant contribution. This is discussed further in the second half of Appendix B.
By substituting the expansion (4.9) into (4.8) we can distinguish three contributions in the differential cross correlation: is the contribution of the direct waves (i.e. those which have not interacted with the random scatterers) which is given by (3.4), -∆C (1) 1 is the contributions of the waves that have been scattered once by the random scatterers, which is given by a sum of terms of the form (3.4) in which one of the factorŝ G 0 is replaced byĜ 1 , -∆C (1) 2 is the contributions of the waves that have been scattered twice by the random scatterers, which is given by a sum of terms of the form (3.4) in which two of the factorŝ G 0 are replaced byĜ 1 or one of the factorsĜ 0 is replaced byĜ 2 , -we neglect higher-order terms.

4.4.
Statistical analysis of the differential cross correlation. The contributions of the direct waves to the differential cross correlation are described in Proposition 3.1. In particular, the singular contribution (3.5) of the direct waves gives a peak to the backlight imaging functional (3.12) at the target location, but it does not give any contribution to the daylight imaging functional (3.11) at the target location.
We now analyze the contributions ∆C of the scattered waves to the differential cross correlation in a backlight illumination configuration. Proposition 4.1 shows that scattering generates zero-mean random fluctuations in the cross correlations. This will reduce the SNR as we will discuss in the next subsections. Proposition 4.2 shows that scattering can increase the resolution by enhancing the directional diversity, which generates additional singular peaks that were not created by the direct illumination and that can be used in migration imaging.
The following proposition is proved in Appendix A. It describes the first two moments of cross correlation ∆C   1 (x 1 , x 2 , τ ) have mean zero and variance whose leading-order terms are as follows.
At lag time equal to the difference of travel times τ = T (x 2 , z r ) − T (x 1 , z r ) the variance of the fluctuations is Var ∆C At lag time equal to the sum of travel times τ = T (x 2 , z r )+T (x 1 , z r ) the variance of the fluctuations is At lag time equal to minus the sum of travel times τ = −T (x 2 , z r ) − T (x 1 , z r ) the variance of the fluctuations is Var ∆C In order to complete the picture, we can add that: 1) The coherence time of the fluctuations is the inverse of the bandwidth of the noise sources.
2) The time-integrated variance is: 3) The quantity ρ K x1,x2,zr is the integral of the function z → ρ(z) K 2 zr,z + K 2 zr,x1 over the surface S x1,x2,zr that has the form of an ellipsoid whose main axis is (z r , x 2 ). The support of this function is the one of the scattering region. Therefore, the quantity If the sensor array is centered at x 0 and the diameter of the array is smaller than the distance from the array to the reflector, then the ellipsoid is an oblate spheroid centered at x c = (z r + x 0 )/2, its main axis is the line from z r to x 0 , the polar radius is (3/2)|z r − x 0 | and the equatorial radius is √ 2|z r − x 0 | (see Figure 4.1): The following proposition is proved in Appendix B. It describes the first moment of cross correlation ∆C  has a singular contribution at lag time equal to the sum of travel times T (x 2 , z r ) + T (x 1 , z r ) provided that the ray going from z r to x 1 reaches into the scattering region behind the sensors (see the left picture in Figure 4.2). This singular contribution has the form (4.20) The mean of the cross correlation ∆C (1) 2 has a singular contribution at lag time equal to minus the sum of travel times T (x 2 , z r )+T (x 1 , z r ) provided that the ray going from z r to x 2 reaches into the scattering region behind the sensors (see the right picture in Figure 4.2):   The variance (4.13) of the fluctuations of the cross correlation at lag time equal to the difference of travel times depends only on the scatterers localized along the rays that participate to the singular peak at the difference of travel times and plotted in Figure 3.1c-d to leading order. The variance is especially sensitive to the scatterers close to the reflector at z r and close to the sensor array, as can be seen from the denominators in (4.14).
We consider two typical situations: -if there are scatterers only beyond the sensor array (as in Figure 4.3a), then these scatterers do not produce fluctuations in the cross correlation at the difference of travel times, so that the singular contribution at this lag time is not affected and the backlight imaging functional will exhibit a clear but elongated peak at the target location.
-if there are scatterers only between the sources and the target, then there are fluctuations in the cross correlation at lag time equal to the difference of travel times. If we denote by d(A, R) the distance from the sensor array to the reflector, by d − (A, S) the distance from the sensor array to the scattering region that is between the sources and the reflector, by W K the width of the source region, by W − s the width of the scattering region that is between the sources and the reflector, then we have .  This formula allows us to discuss the performance of the cross correlation technique in a scattering medium from the point of view of the SNR of the cross correlation at lag time equal to the difference of travel times, which exhibits a singular contribution used in the backlight imaging functional.
1) The SNR does not depend on the distance from the sensor array to the sources.
2) The SNR depends only on the scatterers localized along the rays that participate to the singular peak at lag time equal to the difference of travel times and plotted in Figure 3.1c-d.
3) The SNR depends on the distance from the sensor array to the reflector and it decays with this distance.
4) The SNR depends on the distance from the sensor array to the scattering region and it increases with this distance.
The noise sources have usually very large bandwidth, so it is possible to filter the recorded signals to select a given frequency band of the form [ω 0 − B/2, ω 0 + B/2] (with B < ω 0 ). We can see that the resolution is proportional to B −1 while the SNR is proportional to B 1/2 /ω 2 0 . By increasing the bandwidth, it is therefore possible to increase the SNR and to enhance the resolution.
When forming the imaging functional the cross correlations are summed as in (3.11). If we assume that the distances between the sensors are larger than c 0 /B, then the fluctuations ∆C (1) 1 (τ, x j , x l ) j,l=1,...,N are uncorrelated from each other for all distinct pairs. Therefore the noise level is reduced by a factor N/ √ 2: This means that, even if the singular peaks of the cross correlations are buried into the fluctuations due to the clutter, the imaging functional helps in building up the SNR.
4.6. The cross correlation at the sum of travel times. At a lag time equal to the sum of travel times T (x 1 , z r ) + T (x 2 , z r ), the cross correlation exhibits a singular peak and random fluctuations. The width of the peak and the coherence time of the fluctuations are equal to B −1 . It is possible to give the order of magnitude of the height of the peak and of the standard deviation of the fluctuations of the cross correlation at lag time equal to the sum of travel times from (4.19) and (4.15). The important remark is the following one: The mean (4.19) of the cross correlation has a singular contribution at lag time equal to the sum of travel times provided that the rays going from the reflector to the sensors intersect the scattering region. The variance (4.15) of the fluctuations of the cross correlation at lag time equal to the sum of travel times depends only on the scatterers localized along the ellipsoid (4.18) to leading order. If the scattering region intersects the rays going from the reflector to the sensors but is outside the volume delimited by the ellipsoid (4.18), then the cross correlation exhibits a clear peak at the sum of travel times.
If the scattering region intersects the rays going from the reflector to the sensors and the ellipsoid (4.18), then we have the following. If we denote by d(A, R) the distance from the sensor array to the reflector, by d(A, S) the distance from the sensor array to the scattering region, by d(K, S) the distance from the source region to the scattering region, by W K the width of the source region, and by W + s the width of the scattering region that is behind the sensor array, then we have (4.26) The SNR of the main peak at the sum of travel times T (x 1 , z r ) + T (z 1 , z r ) over the standard deviation of the fluctuations of the cross correlations at this lag time is Let us discuss the main properties of the SNR of the cross correlation at lag time equal to the sum of travel times: 1) The SNR is high if the scattering region is far enough from the sensor array (more exactly, if it is outside the volume delimited by the ellipsoid (4.18)).
2) It is sensitive to the scattering region behind the sensor array when not too far.
3) It increases with the source volume but decays with the square distance from the source region to the scattering region. 4) It increases with the bandwidth. This shows that both the resolution and the SNR increases with the bandwidth. Note that the third point simply means that the SNR is proportional to the energy emitted by the noise sources and received by the scattering region, which is then scattered to provide the secondary daylight illumination.

4.7.
On the trade-off between resolution enhancement and signal-tonoise ratio reduction. We can now discuss the trade-off between resolution enhancement and SNR reduction due to scattering. 1) Proposition 4.2 shows that scattering leads to the appearance of a singular contribution in the cross correlation at lag time equal to (plus or minus) the sum of travel times. This happens provided there are rays issuing from the scattering region and going to the sensors and then to the reflector. The scatterers can be thus seen as secondary sources that can provide a daylight illumination for the reflector, as shown in Figure 4.2. This ensures that the daylight imaging functional (3.11) will have a peak at the reflector location with very good range resolution. This is very interesting result, the one we wanted to get, because in the backlight illumination configuration considered here, the direct waves do not generate such a peak and only the backlight imaging functional (3.12) can be used in absence of scattering, which produces an elongated peak with poor range resolution.
2) Proposition 4.1 shows that the scattered waves also involve fluctuations in the cross correlations that can be larger than the additional peaks exhibited in Proposition 4.2. As a consequence, the peak in the daylight imaging functional (3.11) at the reflector location that we have just mentioned can be buried in the contributions of the non-singular random components. This happens in the setup of Figure 6.1, which we discuss in Subsection 6.2.
3) There exist special configurations which are favorable for imaging with secondary illumination from scatterers. If the scatterers are far enough from the sensor array, in the sense that the scattering region is outside the volume delimited by the ellipsoid (4.18), then the fluctuations of the cross correlation at lag time equal to the sum of travel times are small. This situation does not prevent the existence of a ray going through the reflector, a sensor, and a scatterer, which ensures the existence of a singular peak in the cross correlation at lag time equal to the sum of travel times.
4) It may happen that these favorable conditions are not met in practice. We will see in Section 6 that an iterated cross correlation technique can enhance the contributions of the scattered waves in order to strengthen the singular components of the scattered waves.

Migration imaging of cross correlations.
In the previous subsection we saw that the scatterers can play the role of secondary sources. Here we illustrate this fact by showing that it is possible to use the daylight imaging functional (3.11) in configurations that are not a priori in a daylight illumination configuration.
In Figure 4.3 we consider a configuration with a backlight illumination of the region of interest and with a layer of scattering medium behind it. We plot the images obtained with the backlight imaging functional (3.12) and with the daylight imaging functional (3.11). The backlight imaging functional gives a good image (Figure 4.3b), as in Figure 3.2b, which is expected since it is the result of the contributions of the direct waves described in Proposition 3.1. More striking is the daylight imaging functional (Figure 4.3c) that also gives a good image, in contrast to Figure 3.2c. This shows that the layer of scattering medium succeeds in redirecting part of the flux of energy so that the region of interest experiences a secondary daylight illumination and the daylight imaging functional is effective. By comparing with the primary daylight configuration of Figure 3.3, one can see that the daylight imaging functional in Figure 4.3c has range and cross-range resolutions that are comparable to the ones obtained in a primary daylight configuration. Here the situation is favorable in that the scattering region does not intersect the ellipsoid (4.18), so that the fluctuations of the cross correlation at lag time equal to the sum of travel times are small and the SNR of the image in Figure 4.3c is high. The reduction of the SNR can become problematic if the scattering region intersects the ellipsoid (4.18). We will see in Section 6 how to enhance the contributions of the scattered waves by an iterative cross correlation technique.

5.
Passive sensor imaging with a reflecting interface. In this section we consider the case in which scattering is not generated by random inhomogeneities but by a partially or totally reflecting interface in the medium. This configuration is of interest for instance in geophysics where the earth surface plays the role of a reflecting interface for the earth's crust. The Green's function in the presence of the reflecting interface and in the absence of the reflector is then simplŷ with R m the reflection coefficient of the interface andĜ 0 is the free space Green's function (3.2). In the case of a perfect mirror we have R m = −1. In the presence of the reflector at z r and in the point interaction approximation for the reflector the Green's function is of the form G m,r ω, x, y =Ĝ m ω, x, y + ω 2 c 2 0 σ r l 3 rĜm ω, x, z r Ĝ m ω, z r , y , and it contains the contributions of rays reflected by the interface. In the next proposition 5.1, proved in Appendix C, the virtual source region denotes the image of the support of the function K through the symmetry with respect to the plane I m . Proposition 5.1 shows that the mirror generates virtual sources which play the role of secondary noise sources and can provide a daylight illumination of the region of interest (the reflector and the sensors). Beyond the interesting additional singular components described below, at lag times equal to plus or minus the sum of travel times that contribute to the daylight imaging functional, there are additional components at lag times which are at least of the order of twice the travel time from the sensors to the reflecting interface, so that they do not play any role in the daylight imaging functional provided the reflecting interface is far enough.
Proposition 5.1. Let us consider the backlight illumination configuration with a mirror behind the sensor array ( Figure 5.1a). The differential cross correlation has several singular components.
There is a singular contribution at lag time equal to the difference of travel times as described by (3.5).
There is a singular contribution at lag time equal to the sum of travel times T (x 2 , z r ) + T (x 1 , z r ) which has the form: there is a singular contribution at lag time equal to minus the sum of travel times −[T (x 2 , z r ) + T (x 1 , z r )] which has the form:

2)
where K x,z is the energy released by the virtual sources along the ray joining x and z: Note that K xj ,zr is not zero only if the ray going from z r to x j extends into the virtual source region, which is the secondary daylight illumination configuration mentioned above. Figure 5.1 we plot the images obtained in a backlight illumination configuration with a mirror behind the region of interest. The image obtained with the backlight imaging functional (3.12) is expected since it is the result of the contributions of the direct waves described in Proposition 3.1. The good image obtained with the daylight imaging functional (3.11) illustrates the theoretical predictions of Proposition 5.1: the mirror generates virtual sources that illuminate the region of interest with a secondary daylight illumination.

Migration imaging of cross correlations. In
The configuration with a mirror and the configuration with a randomly scattering medium studied in the previous section have therefore common points: They both provide a secondary daylight illumination. As a result, the daylight migration functional can use the scattered or reflected waves and gives a good range resolution of the reflector, while the backlight migration functional uses the direct waves and gives a good cross range resolution of the target. By multiplying the two functionals, one obtains the location of the reflector with a good accuracy (see Figure 5.2).
There are two major differences between the configuration with a mirror and the configuration with a randomly scattering medium: 1) A very significant part of the energy flux is reflected in the case of the mirror, while only a small part is scattered in the case of a randomly scattering medium (in the weakly scattering regime that we address). 2) The randomly scattering medium redirects the energy flux in all directions, while the reflecting interface redirects the energy flux only in the specular direction. This means that the mirror induces an enhancement of the directional diversity that is strong but only in a special direction. Therefore the enhancement of the directional diversity by a mirror is not as effective as the one provided by a randomly scattering medium. If the orientation of the interface is not adjusted properly, in the sense that the condition "there is a ray issued from a virtual source point and going to a sensor and the reflector" is not fulfilled, then the daylight imaging functional does not perform well as predicted by Proposition 5.1 and as shown in Figure 5.3. In contrast to this result, the orientation of the randomly scattering layer plays no role, as predicted by Proposition 4.1 and as shown in Figure 5.4. 6. Iterated cross correlations for passive imaging in a randomly scattering medium. We have noticed in Section 4 that the peaks in the differential cross correlation ∆C (1) (τ, x j , x l ) that are relevant to imaging of reflectors can be buried in fluctuations. This happens when the SNR of the peaks at lag time equal to plus or minus the sum of travel times ±[T (x j , z r ) + T (x l , z r )] is low. In this section we propose to use an iterated cross correlation technique that masks the contributions of the direct waves and increases the effective SNR of the peaks. This technique was shown to be efficient for inter-sensor travel time estimation in [9,20] and is briefly reviewed in Appendix D. In the next subsections we describe this technique for reflector imaging with secondary daylight illumination from scattering and present the results of numerical simulations.
6.1. The coda cross correlation. It is possible to form a special fourth-order differential cross correlation matrix ∆C (3) T (τ, x j , x l ) between sensors (x j ) j=1,...,N from the differential cross correlations ∆C T (τ, x j , x l ) obtained from the recorded data. This is done as follows.
2) Cross correlate the tails of the differential cross correlations and sum them over all complementary sensors in the array to form the coda cross correlation between x j and x l : (6.1) This algorithm is parameterized by three important times: 1) The time T is the integration time and it should be large so as to ensure statistical stability with respect to the distribution of the noise sources.
2) The time T c1 is chosen so that the parts of the Green's functions t → G(t, x k , x j ) and t → G(t, x k , x l ) limited to [T c1 , T c2 ] do not contain the contributions of the direct waves. This means that T c1 depends on the index of the sensors j, l, k and should be a little bit larger than max(T (x k , x j ), T (x k , x l )).
3) The time T c2 should be large enough so that the parts of the Green's functions t → G(t, x k , x j ) and t → G(t, x k , x l ) limited to [T c1 , T c2 ] contain the contributions of the incoherent scattered waves. This means that T c2 should be of the order of the power delay spread.
As shown in [9], it follows from the statistical stability of the cross correlation C T that the differential coda cross correlation ∆C (3) T is a self-averaging quantity and it is equal to the statistical differential coda cross correlation ∆C (3) as T → ∞: The time windowing is very important because it selects the contributions that we want to use for imaging the reflector at z r . The asymptotic analysis of the functional ∆C (3) can be carried out under the same conditions and with the same tools as in Section 4. It involves long and tedious calculations whose detailed results do not add to the qualitative understanding. We can summarize them by stating that the differential coda cross correlation ∆C (3) has singular components at lag time equal to (plus or minus) the sum of travel times ±[T (x 1 , z r ) + T (x 2 , z r )], and has less additional terms than the usual differential cross correlation studied in Subsection 4.2. As a consequence we anticipate that migration of the differential coda cross correlation using the daylight migration functional will produce an image of the reflector with a higher SNR.
An alternative version of the differential coda cross correlation can be computed by the following algorithm: 1) Calculate the tails of the cross correlations C T,coda using the data C T collected in the presence of the reflector and the tails of the cross correlations C T,coda,0 using the data C T,0 collected in its absence: 2) Compute the coda cross correlations C T using C T,coda and C T,0 using C T,coda,0 : T,0 (τ, x j , x l ) = N k=1,k ∈{j,l} C T,coda,0 (τ , x k , x j )C T,coda,0 (τ + τ, x k , x l )dτ .
This algorithm produces a version ∆C (4) T of the differential coda cross correlation that is different from ∆C T . It is a statistically stable quantity. Its statistical average ∆C (4) contains the singular components at lag time equal to (plus or minus) the sum of travel times ±[T (x 1 , z r ) + T (x 2 , z r )] that we want to use for imaging and it has less unwanted contributions compared to ∆C (1) , similarly as ∆C 3) . As a result we anticipate that both versions ∆C (3) T and ∆C (4) T of the differential coda cross correlation will produce images by migration.
6.2. Passive imaging of a reflector with coda cross correlation. In Figure  6.1 we consider a situation in which the scattering region intersects the ellipsoid (4.18),  Figure a plots the image obtained with the backlight imaging functional (3.12) applied with the differential cross correlation ∆C (1) . Figure b plots the image obtained with the daylight imaging functional (3.11) applied with the differential cross correlation ∆C (1) . Figure d plots the image obtained with the backlight imaging functional (3.12) applied with the differential coda cross correlation ∆C (3) . Figure e plots the image obtained with the daylight imaging functional (3.11) applied with the differential coda cross correlation ∆C (3) . Figure f plots the image obtained with the daylight imaging functional (3.11) applied with the differential coda cross correlation ∆C (3) multiplied by the backlight imaging functional (3.12) applied with the differential cross correlation ∆C (1) . Figure g plots the image obtained with the backlight imaging functional (3.12) applied with the differential coda cross correlation ∆C (4) . Figure h plots the image obtained with the daylight imaging functional (3.11) applied with the differential coda cross correlation ∆C (4) . Figure i plots the image obtained with the daylight imaging functional (3.11) applied with the differential coda cross correlation ∆C (4) multiplied by the backlight imaging functional (3.12) applied with the differential cross correlation ∆C (1) .
which involves fluctuations of the cross correlation at lag time equal to the sum of travel times. As a result the daylight imaging functional with the usual differential cross correlation technique does not exhibit a peak at the reflector location ( Figure  6.1b), because the peaks at lag times equal to the sums of travel times ±[T (x 1 , z r ) + T (x 2 , z r )] are buried in the other components, contrarily to what happens in Figure  4.3 when the scattering region does not intersect the ellipsoid. However, the daylight imaging functional with the differential coda cross correlation technique ∆C (3) (figure 6.1e) or ∆C (4) (figure 6.1h) gives a much better image. The overall result is that the backlight imaging functional I B used with ∆C (1) has good cross-range resolution and SNR while the daylight imaging functional I D used with ∆C (3) or ∆C (4) has good range resolution and SNR. The multiplication of these two functions gives the location of the reflector with greater accuracy (Figure 6.1f or i).

Conclusion.
In this paper we have analyzed the role of wave scattering in passive sensor imaging. The main result is that a randomly scattering medium can increase the directional diversity of the energy flux that illuminates the region of interest consisting of the sensor array and the reflectors to be imaged. The range resolution can be enhanced significantly when the scattering medium provides a secondary daylight illumination of the region of interest, which generates a peak at the reflector location in the daylight imaging functional (3.11). We have shown that a randomly scattering layer is more efficient for directional diversity enhancement than a reflecting surface.
However, a randomly scattering medium also introduces random fluctuations in the cross correlations which reduce the signal-to-noise ratio of the peaks in the cross correlations that contribute to the formation of the image. When the scattering region is far enough from the sensor array, these fluctuations are small. When the scattering region is close, we have remarked that the contributions of the reflectors in the cross correlation data can be weak compared to the standard deviation of the fluctuations. Migration (backpropagation) of the cross correlation then gives blurred and speckled images. We have seen that, in such configurations, it is possible to migrate a special fourth-order cross correlation that evaluates the correlations of the tails (coda) of the second-order cross correlations. This technique aims at eliminating the contributions of the direct waves and enhancing the contributions of the scattered waves. As a result, passive imaging using differential fourth-order coda cross correlations and a product of two different (backlight and daylight) migration functionals turns out to be quite effective, in that it gives an image with good SNR and good resolution.
This trade-off between enhancement of the directional diversity of the illumination along with reduction of the SNR of the image makes sense only when the scattering is weak. As scattering strength increases migration of cross correlations, including fourth order cross correlations, will eventually fail. It may be possible to estimate inter-sensor travel times, as was done in a recent study [5] in the microwave regime, in an indoor environment, using special statistical techniques, but it will be difficult to image in such strongly scattering media. Therefore, correlation based methods for passive sensor imaging in scattering media are likely to be effective in regimes of weak to intermediate scattering.
In the backlight illumination configuration this means that the four points should be along the same ray with y 1 , y 2 → z r → z or y 1 , y 2 → z → z r .
We first study the contribution (i) of the scatterers z such that y 1 , y 2 → z r → z. We introduce the unit vectorê 1 = z r − z |z r − z| and complete it with two other unit vectors (ê 2 ,ê 3 ) so that (ê 1 ,ê 2 ,ê 3 ) is an orthonormal basis. We make the change of variables y j → (a j , b j , c j ), for j = 1, 2: with a j > 0 since we consider here that y j → z r → z. The Jacobian of the change of variables is ε|z r − z| 3 . This gives a parameterization of the variable y j along the ray joining z and z r . We obtain Var ∆C We compute the integrals in h, b 1 , c 1 , b 2 , c 2 by making use of the identity and we make the changes of variables a j = l j /|z r − z|, j = 1, 2, so that we obtain Var ∆C 1,1 (τ, x 1 , x 2 ) (i) = σ 2 r l 6 r σ 2 s l 3 s 2 15 π 7 c 5 0 ε 5 We second study the contribution (ii) of the scatterers z such that y 1 , y 2 → z → z r . We introduce the unit vectorê where the rapid phase is This is the cross correlation of waves that have interacted once with the scattering medium and once with the reflector with waves that have interacted once with the scattering medium. That is why the average (with respect to the scattering medium) is not zero. The stationary points satisfy the three conditions The rapid phase does not depend on y, so the stationary points need to fulfill the two conditions The first condition means that x 2 and z r must be on the same ray issued from z. If z → z r → x 2 , then the second condition means τ = −T (x 1 , z r ) + T (x 2 , z r ). This corresponds to a backlight illumination from the secondary source z.
If z → x 2 → z r , then the second condition means τ = −T (x 1 , z r ) − T (x 2 , z r ). This is the contribution corresponding to a daylight illumination from the secondary source z and plotted in the right Figure 4.2. We focus our attention to the backlight illumination configuration. We introduce the unit vector and we complete it with two other unit vectors (ê 2 ,ê 3 ) so that (ê 1 ,ê 2 ,ê 3 ) is an orthonormal basis. We make the change of variables z → (a, b, c) with This gives a parameterization of the variable z along the ray joining x 2 and z r . Here the scatterers are behind the sensor array which means that we restrict ourselves to a > 0. We also parameterize the lag time τ around the difference of travel times: We perform a Taylor series expansion of the rapid phase The second condition (in ∇ z ) cannot be fulfilled in our geometry (Figure 4.3a), in which the scatterers z are on the opposite side of the sensor x 2 , the reflector z r , and the sources y. Therefore the term E[∆C (1) 2,2 (τ, x 1 , x 2 )] gives a negligible contribution. The cross correlation terms that involve one wave with one componentĜ 2 and one wave with onlyĜ 0 have similar properties as the term ∆C (1) 2,2 studied here above and it is not possible to find stationary points or maps of points for their expectations. This is again due to our geometry in which the scatterers are on the opposite side of the sensor array, the reflector, and the sources. Therefore they do not contribute to the mean of the cross correlation.
Finally, none of the terms considered in this subsection can contribute to the variance at order σ 2 s (the variances of these terms are of order σ 4 s ), so the variance of the cross correlation is determined by the terms ∆C (1) 1 studied in Appendix A (whose variances are of order σ 2 s ). Appendix C. Proof of Proposition 5.1. The proof of Proposition 5.1 is not completely obvious because the virtual sources are fully correlated with the noise sources. By expanding the expression of the cross correlation in terms of the background Green's functionĜ 0 , we obtain that the differential cross correlation is the contribution of sixteen terms, which correspond to take into account the virtual sources and the virtual reflector. We will study in detail one of them: where the rapid phase is T 0 (y) = −T (ỹ, x 1 ) + T (ỹ, z r ) + T (z r , x 2 ) − τ.
The analysis using the stationary phase theorem then goes along the same lines as in the previous appendices. The analysis of the other terms is similar and completes the proof.
Appendix D. Inter-sensor travel time estimation with coda cross correlations. It is possible to estimate the travel time between two sensors x 1 and x 2 in a scattering medium by looking at the main peaks of a special fourth-order cross correlation C T (τ, x 1 , x 2 ). This fourth-order cross correlation uses the data recorded by an array of auxiliary sensors x a,k , k = 1, . . . , N a , and it is evaluated as follows.
1) Calculate the cross correlations between x 1 and x a,k and between x 2 and x a,k for each auxiliary sensor x a,k : x a,k )u(t + τ, x l )dt, l = 1, 2, k = 1, . . . , N a .
3) Cross correlate the tails of the cross correlations and sum them over all auxiliary sensors to form the coda cross correlation between x 1 and x 2 : T (τ, x 1 , x 2 ) = Na k=1 C T,coda (τ , x a,k , x 1 )C T,coda (τ + τ, x a,k , x 2 )dτ .
The statistical coda cross correlation C (3) was studied in [9] by a stationary phase analysis. It differs from the statistical cross correlation C (1) in that the contributions of the direct waves are eliminated and only the contributions of the scattered waves are taken into account (note that some of the contributions of scattered waves are also eliminated, but only those which correspond to small additional travel times, which are also those which induce small directional diversity). Since scattered waves have more directional diversity than the direct waves when the noise sources are spatially localized, the coda cross correlation C (3) (τ, x 1 , x 2 ) usually exhibits a stronger peak at lag time equal to the inter-sensor travel time T (x 1 , x 2 ). In particular, in contrast with the cross correlation C (1) , the existence of a singular component at lag time equal to the travel time T (x 1 , x 2 ) does not require that the ray joining x 1 and x 2 reaches into the source region, but only into the scattering region. We illustrate these results in Figures D.1-D.2 in which the five sensors are aligned perpendicularly to the energy flux coming from the noise sources and the cross correlation C (1) (τ, x 1 , x j ), j = 1, . . . , 5 does not have a peak at lag time equal to the travel time between the sensors x 1 and x j . In Figure D has a peak at lag time equal to the travel time between the sensors. In Figure D.2 the ray going through the sensors x 1 and x j does not intersect the scattering layer and the coda cross correlation C (3) (τ, x 1 , x j ) does not have a peak at lag time equal to the travel time between the sensors. Note the presence of the auxiliary sensors that are necessary for the evaluation of the coda cross correlation.