Imaging moving targets through scattering media

: Optical microscopy in complex, inhomogeneous media is challenging due to the presence of multiply scattered light that limits the depths at which diﬀraction-limited resolution can be achieved. One way to circumvent the degradation in resolution is to use speckle- correlation-based imaging (SCI) techniques, which permit imaging of objects inside scattering media at diﬀraction-limited resolution. However, SCI methods are currently limited to imaging sparsely tagged objects in a dark-ﬁeld scenario. In this work, we demonstrate the ability to image hidden, moving objects in a bright-ﬁeld scenario. By using a deterministic phase modulator to generate a spatially incoherent light source, the background contribution can be kept constant between acquisitions and subtracted out. In this way, the signal arising from the object can be isolated, and the object can be reconstructed with high ﬁdelity. With the ability to eﬀectively isolate the object signal, our work is not limited to imaging bright objects in the dark-ﬁeld case, but also works in bright-ﬁeld scenarios, with non-emitting objects.


Introduction
Optical imaging is challenging in turbid media, where multiple scattering of light causes a degradation of resolution and limits the depths at which we can reliably image (< 1mm in biological tissue) without having to resort to destructive optical clearing or sectioning techniques [1]. Many approaches currently exist to filter out the multiply scattered light and detect only the unscattered (ballistic) or minimally scattered photons. These include methods such as time and coherence gating, which separate the ballistic photons from the scattered photons based on their transit time to the detector [2,3]; methods that rely on preserving the initial angular momentum or polarization modulation [4][5][6][7]; and methods that rely on spatial confinement, such as confocal and multi-photon microscopy [1,8]. An issue with methods that rely on detecting only the minimally scattered photons is the maximum achievable depth of penetration, since the chance of detecting a quasi-ballistic photon decreases exponentially with increasing depth.
Instead of rejecting the scattered photons, other approaches have aimed to take advantage of the information inherent within the detected speckle field that arises from multiply scattered light. Wavefront shaping (WFS) techniques exploit the principles of time reversal to undo the effects of scattering and enable focusing of light in thick, scattering media [9][10][11][12]. However, WFS usually requires long acquisition times to measure the transmission matrix, and/or the presence of a guide star. On the other hand, speckle-correlation-based imaging (SCI) approaches exploit the angular correlations inherent within the scattering process to reconstruct the hidden object and do not need long acquisition times or a guide star [13,14]. However, SCI methods are limited to working in dark-field scenarios, with sparsely tagged objects [14], since the detected light must consist solely of light arising from the object.
In this work, we demonstrate imaging of hidden moving objects in a bright-field scenario by leveraging the temporal correlations inherent in the scattering process to separate and remove the dominating contribution from the background [15,16]. To create a spatially incoherent light source, a spatial light modulator (SLM) was used to apply the same set of random phase patterns during different acquisitions. The use of a deterministic phase modulator ensured that the background contribution remained constant across the detected images. By removing the background component, the speckle pattern from the object was isolated, and the object was reconstructed with high fidelity. Using this technique, we experimentally demonstrate successful recovery of moving objects that would otherwise be obscured by scattering media. Figure 1 presents an overview of our system. A moving object, hidden at a distance u behind a scattering media, is illuminated using a spatially incoherent, narrow-band light source. The scattered light is detected by a high-resolution camera that is placed at a distance v from the scattering media.

A. System Setup
Object at 1 , 2 Acquired Images B. Captured Images with and without scatterer

C. Hidden Object Retrieval
Camera u v Fig. 1. Principle behind non-invasive imaging of obscured moving objects. A) A spatially incoherent light source illuminates a moving object hidden behind a visually opaque turbid media. The resultant speckle field is captured by a camera sensor. B) Speckle images are acquired by the camera sensor at different times, with the object moving between the captures. The scattering media prevents us from resolving the object. C) The hidden object can be retrieved from the seemingly random speckle images by taking advantage of inherent angular correlations in the scattering pattern. i) Each captured image I n consists of a background, B, subtracted by the imaged object, where the imaged object is the convolution of the PSF of the scattering media, S, and the object pattern, O. ii) Although the background signal dominates over the object, it can be subtracted out by taking the difference between the two captured images ∆I. iii) The object autocorrelation O O is approximated by autocorrelating the difference image ∆I. iv) The hidden object can be reconstructed from the object autocorrelation by using phase retrieval techniques.
In the absence of any correlations in the scattering pattern, the detected image is merely a speckle intensity field. However, by exploiting the deterministic nature of scattering, the hidden object can be recovered [ Fig. 1(c)]. Let us first consider the case where light is confined to emit solely within an isoplanatic range, as defined by the angular memory effect (ME). In this case, the detected light can be mathematically represented as where S is the point spread function (PSF) of the light scattering process, or equivalently the speckle intensity distribution at the camera arising from a single point source at the object plane; and O is the object, defined as the collection of points through which light can be transmitted [14].
For this paper, we use the operator * to denote convolution. The memory effect region can be approximated as δx = uλ π L , where L is the thickness of the scattering media, λ is the wavelength of light, and u is the distance between the scattering media and the object.
If we now consider the case of an absorptive object in a bright-field scenario, then the majority of the detected light arises from the background. Using superposition, the detected intensity image I can be mathematically described as where B is the speckle intensity image arising from the scattered light transmitted through the medium, and S * O is the portion that the object would have contributed if it were transmitting, as opposed to blocking, light [ Fig. 1(c,i)]. Due to the dominating contribution from the background B, we cannot retrieve O from I alone. By acquiring multiple intensity images with the background, but not the object, constant between acquisitions, we can remove the background signal and thereby retrieve the object. One strategy to achieve this is to use a moving object. If the object dimensions falls within the ME region, the contribution of the object in each image can be represented as the convolution of the object pattern with an acquisition-dependent PSF. As long as the rest of the sample is static, the speckle field arising from the background will remain unchanged and can be subtracted out by taking the difference between captures. That is, where I n denotes the n th captured image. Since the scattering PSF is a delta-correlated process (S n (x) S n (x) ≈ δ(x)), taking the autocorrelation (AC) of the image ∆I yields the object autocorrelation (OAC), plus additional noise terms [ Fig. 1(c,iii)]. That is, where denotes autocorrelation. We shall refer to ∆I n ∆I n as the speckle autocorrelation (SAC). The object can be recovered from the SAC by using phase retrieval techniques, such as the Fienup iterative phase retrieval methods, to recover the Fourier phase [ Fig. 1(c,iv)] [17]. The resultant object will have an image size dictated by the magnification of the system, M = − v u .

Effect of travel distance
Depending on the distance traveled by the object, the PSFs S n , n = 1, 2, ... may or may not be correlated. Figure 2 illustrates the effect of travel distance, relative to the ME range, on the SAC. The speckle intensity images I 1 , I 2 were determined using simulation. For comparison, the autocorrelation of the objects has also been provided [ Fig. 2(A, "Object AC")]. For simplicity, only the case of two image captures (n = 1, 2) has been considered. For a moving object, the associated PSFs S 1 , S 2 will have a degree of correlation C(∆x) based on the object travel distance ∆x. For scattering media with thicknesses L greater than the mean free path, the degree of correlation can be approximated using the angular correlation function where k = 2π λ , L is the thickness of the scattering medium, and Θ ≈ ∆x u [18][19][20]. When C(∆x) > 0.5, the object is considered to have traveled within the ME field of view. The following sections describe three possible cases in more detail: iii.

Speckle Correlation
Distance Object Travelled (A.U.)  A) The scattering PSFs experienced by an object have a degree of correlation C(∆x) that depends on the distance the object traveled. When C(∆x) ≥ 0.5 (shown in red), the object is considered to have traveled within the memory effect (ME) region. For comparison, the object and its autocorrelation (AC) are displayed. B) When the object travels inside the ME region, the SAC contains three copies of the object autocorrelation (OAC): a centered, positive copy and two negative copies shifted by an amount proportional to the object travel distance. The OAC can be determined by either deconvolving the SAC or by thresholding out the negative portions (negative with reference to the mean, background level). The object can be reconstructed from the estimated OAC using phase retrieval techniques. C) When the object travels a distance where C(∆x) ≈ 0, only a single copy of the OAC is seen, with additional noise from the cross-correlation between uncorrelated PSFs. The normalized colormap used to display the AC and reconstructed object, with 0 corresponding to the mean background level.

C. Outside Memory Effect Region
In the case where the object travels a small distance (such that C(∆x) ≈ 1), we have where x = (x, y), x i = (x i , y i ) are coordinates in the object plane and image plane respectively, ∆x is the distance the object traveled in the object plane, and ∆x i = M∆x. We can equivalently consider the PSF to be the same in both captures and have the object travel between captures. That is, where A = O O is the object autocorrelation (OAC). The SAC contains three copies of the OAC: a positive copy centered at x = (0, 0), and two negative copies shifted by an amount commensurate with the object travel distance [ Fig. 2(B, "Speckle AC")].
Since C(∆x) ≈ 1 when ∆x ≈ 0, the object may travel a distance shorter than the extent of its autocorrelation. In this case, the SAC will yield positive and negative copies of the OAC that overlap [ Fig. 2(i)]. The OAC can be recovered using deconvolution [ Fig. 2(i, "Deconv. SAC.")]. Using thresholding to remove the negative portions will adversely impact the positive copy and result in an incomplete estimation of the OAC [ Fig. 2(i,"SAC>0")]. For the results presented in Fig. 2, the objects were reconstructed by applying an iterative phase retrieval algorithm on the deconvolved SAC ( [13,14,17]).
In the regime where the object travels within the angular ME range (C(∆x) > 0.5), S 1 and S 2 are correlated. To highlight the impact of the degree of correlation C(∆x) on the SAC, we can mathematically represent S 2 as: where S is a speckle intensity pattern that is uncorrelated with S 1 .The scatter PSFs in the equation above are mean-subtracted speckle intensities. Representing S 2 in the form above allows us to preserve speckle intensity statistics (that is, the speckle intensity variance and mean satisfy respectively.) Using Eq. (11), Eqs. (4) and (5) become and where the last equation follows from noting that the speckle fields are a delta-correlated process and that the cross-correlation of two uncorrelated speckle intensities yields noise. The SAC still contains three copies of the OAC. However, the ratio of the intensity of the positive and negative OAC copies is determined by the ME correlation function C(∆x). Moreover, since S 2 S 1 , there is an additional noise term that increases with decreasing C(∆x). Since there is no overlap between the positive and negative OAC copies, the OAC can be retrieved by either thresholding out the portions of the SAC that are smaller than the background value [ Fig. 2  In the case where the object travels outside the memory effect region between captures, S 1 and S 2 are uncorrelated, and Eq. (13) can be simplified as Eq. (5). Comparing the SAC in Fig. 2(iii) with those in Fig. 2(i-ii), we see that the SAC in the case where the object travels farther than the ME region exhibits more noise. This is expected due to the additional noise term caused by S 1 S 2 that is not present in Case 1.
From above, in all cases (for C(∆x) ∈ [0, 1)), we can successfully retrieve the object autocorrelation from the acquired speckle images, S 1 , S 2 . From the estimated OAC, phase retrieval techniques can then be applied to reconstruct the object at diffraction-limited resolution.
An SLM was used in place of a rotating diffuser in order to generate a deterministic, temporally varying set of 50 to 100 random phase patterns. This set of patterns was used for all the acquisitions to ensure that the background light captured remained constant. The object and camera (pco.edge 5.5, PCO-Tech, USA) were placed at a distance u = 20 − 30 cm and v = 10 − 15 cm from the scattering media (DG10-120 diffuser; Thorlabs, USA) respectively (Fig. 3).  ii. iii.

iv.
A. Object B. Camera Image C. Object AC D. Speckle AC E. Reconstruction Fig. 4. Experimental imaging of moving targets hidden behind a diffuser. A) The "object" is hidden behind a scattering medium and attenuates light transmission. The object was moved 1.5 mm between acquisitions. B) Due to the presence of the scattering medium, the object is obscured, and the camera image I 1 is dominated by the scattered light from the background. C) The ideal object autocorrelation (AC). D) The speckle autocorrelation ∆I ∆I ≈ O O. E) By applying phase retrieval on the speckle autocorrelation, the hidden object was reconstructed with high fidelity. Scale bar = 500 µm.
To ensure that only the object moved between successive image captures, a transmissive SLM (tSLM; Holoeye LC2002) coupled with a polarizer (Thorlabs, LPVISE200-A) was used for amplitude modulation, and served as the object (Fig. 4). For each object, a set of n=4 images, I 1 , ....I 4 were acquired, with the object moving 1.5mm between each acquisition. The raw camera images [ Fig. 4(b)] display a seemingly random light pattern that is similar for different objects. This is due to the dominant contribution of the background.
From each successive pair of acquired images, the OAC [ Fig. 4(d)] was estimated by deconvolving the SAC. The deconvolved SAC images were then averaged to reduce noise and yield a better estimate of the OAC. A Fienup-type iterative phase retrieval method was applied to reconstruct the hidden object with high fidelity [Fig. 4(e)] [13,14,17]. One modification that was made to the algorithm was to add an object support to the object constraints; this object support was determined from the OAC support [21,22]. In all cases, the obscured object was successfully reconstructed [ Fig. 4(e)].
To experimentally demonstrate the effect of object travel distance, we moved an object a distance of 0.5, 1, and 3 mm between image acquisitions, and looked at the corresponding SAC and reconstructed object (Fig. 5). As expected, the SAC contained three copies of the OAC. We also compared the effect of processing the SAC using deconvolution [ Fig. 5(b)] vs. thresholding [ Fig. 5(c)].
For Case i, the object traveled a distance ∆x < δx, and both the object and SAC overlapped in space between successive acquisitions. In the case of object overlap, only the non-overlapping portion of the object can be retrieved [ Fig. 5(i)]. Comparing the result of deconvolution vs thresholding, the reconstructed image from the deconvolved SAC more closely resembles the original object [ Fig. 5(i,b)]. However, in both cases, what we are left with is an incomplete OAC and reonstructed object. ii.
iii. For Case ii, the object traveled a distance δx < ∆x ≤ 2δx. Since the OAC support is approximately twice the object support [21], the positive and negative copies of the OAC overlapped [ Fig. 5(ii)]. Due to the overlap, thresholding resulted in an imperfect object reconstruction [ Fig.  5(ii,c)]. In contrast, by deconvolving, the signal from the negative copies can be used to gain a better estimate of the OAC, from which the object can be reconstructed [ Fig. 5(ii,b)].
For Case iii, the object traveled a distance ∆x >> 2δx, and there was no overlap in the SAC. Due to the large ∆x, C(∆x) decreased, and correspondingly, the noise increased. Since the signal-to-noise ratio (SNR) of the negative copies decreased, the entire OAC cannot be seen in the negative copies [ Fig. 5(iii,a)]; thus, performing a deconvolution results in a noisy, imperfect OAC [ Fig. 5(iii,b)], and it is more advisable to use thresholding to retain only the positive portion of the SAC [ Fig. 5(iii,c)]. If we compare the reconstructed objects in both cases, we see that the object from the thresholded result more closely resembles the original object. Schematic of the experimental setup. A spatially incoherent light source is generated by reflecting an expanded laser beam off a spatial light modulator (SLM) that applied a temporally variying random phase pattern. The partially developed speckle field component is blocked, and only the fully-developed speckle field transmits through the moving object and two scattering layers. The emitted scattered light is collected by a camera. An aperture controls the resolution and the speckle size at the camera. B) Experimental result of imaging a moving target. Two speckle intensity images , I 1 , I 2 , were captured, with the target present for the first capture and absent for the second. The background halo from I 1 and I 2 were removed prior to computing the difference ∆I = I 2 − I 1 ≈ S 1 * O. The speckle autocorrelation yielded an estimate of the object autocorrelation, from which the target was retrieved by applying Fienup phase retrieval. Lens focal length = 400 mm.

Imaging moving objects hidden between scattering media
To further demonstrate our imaging technique, we placed a moving object between two diffusers (Newport 10 o Light Shaping Diffuser, Thorlabs DG10-220-MD) [ Fig. 6(A)]. A moving object (a bent black wire) was flipped in and out of the light path between image captures, such that I 2 = B.
We blocked the partially-developed speckle field (from the propagation of the SLM phase pattern) and used only the fully-developed speckle pattern [23]. This fully-developed speckled pattern was transmitted through both scattering media and the moving object. The emitted scattered light was detected by a camera. The background halo from each detected speckle intensity image was estimated and removed by performing Gaussian filtering (500x500 kernel, σ = 100), and then dividing each image by the background halo [14]. The SAC was then computed to estimate the OAC, from which phase retrieval was applied to reconstruct the hidden object. Although the object is fully obscured from both sides by scattering media and cannot be resolved from the camera image alone, using our technique, we were able to successfully reconstruct the hidden object with high fidelity [ Fig.  6(B)].

Discussion and conclusion
In this paper, we demonstrated successful reconstruction of moving targets that were hidden behind an optically turbid media. Although the angular memory effect has already been used to demonstrate imaging of hidden targets, to the best of our knowledge, these prior systems were limited to imaging dark-field, sparsely tagged objects [13,14,24]. We extended this work to imaging in the bright-field scenario by exploiting the temporal correlations inherent in the scattering process to remove the dominating contribution from the background and isolate the signal arising from the object [15,16]. Although we demonstrated our results on non-emitting objects in the bright-field scenario, our technique works equally well with transmissive or reflective objects. A cursory examination reveals that, when I n = B + S n * O and ∆I = I n − I n+1 , the speckle autocorrelation is still given by Eq. (5), similar to imaging absorptive objects in the bright-field scenario. In the remainder of this section, we discuss some of the factors that impact system performance.
Firstly, our method depends on the angular correlations inherent in the scattering process. Thus, the object dimension should fall within the angular memory effect field of view (FOV), approximated using the full-width-half-maximum (FWHM) of the correlation function, uλ π L . The axial extent of the object, δz, should also fall within the axial decorrelation length 2λ π u D 2 [25].
Since the ME FOV is inversely proportional to L, our technique works best with thin scattering media, or through more anisotropically scattering media, since anisotropy enhances the angular memory effect range [20]. Strongly anisotropic media, such as biological tissue, also exhibit the translational memory effect, which may be exploited to further the fidelity of imaging through scattering layers [26].
Secondly, to maximize SNR and minimize overlap, the object travel distance should be such that δx < ∆x and C(∆x) ≥ 0.5, since smaller values of C(∆x) results in higher levels of noise. However, if the object moves such a large distance as to not fall within the laser light beam, then I 2 = B, and ∆I = S 1 * O, and we can also retrieve the object with high fidelity. In all these cases, successful retrieval of the object is dependent on the background light pattern remaining constant between successive image captures. Thus, the illuminated portion of the tissue should remain constant between image captures, and the time between image captures should fall well within the temporal decorrelation time of the scattering sample. For biological samples, the temporal decorrelation time is related to the motion of scatterers embedded within [27].
Imaging through biological samples can be achieved using a faster system. The imaging speed in our current design was limited by the refresh rate of the SLM(≈ 8 Hz) and by the exposure time required to capture an image (50-200ms). With a more powerful laser, or a faster deterministic random phase modulator, it would be possible to shorten our imaging time and extend our work to imaging within non-static samples, such as biological tissue.
A third factor in the fidelity of the reconstruction is the complexity of the object and the size of the background relative to the object. The dynamic range of the camera should be large enough to resolve the equivalent speckle signal from the object. Since the signal contrast is inversely related to the object complexity [14], the dynamic range of the camera limits the maximum object complexity. To maximize the SNR, the camera exposure and laser power should be adjusted such that the full well depth of the camera is utilized. A camera with a larger well depth and dynamic range would provide higher SNR and the capability to image more complex objects. The diameter of the aperture in the system can be adjusted to fine-tune the image resolution and control the object complexity.
Lastly, each speckle grain at the camera should satisfy the Nyquist sampling criterion and be easily resolvable. At the same time, the number of speckle grains that are captured in each image should also be maximized in order to maximize SNR. Although the scattering PSFs are ideally a delta-correlated process, in practice, we are only sampling a finite extent of the PSF. Thus, the PSF autocorrelation yields a delta function plus some background noise which can be minimized by increasing the number of captured speckle grains [14]. Due to Nyquist requirements, the maximum number of speckle grains is a function of the camera resolution; thus, a high resolution camera would provide lower noise. Another method to reduce this speckle noise is to take multiple acquisitions and compute the average of the speckle autocorrelation images.
In conclusion, we demonstrated successful imaging of hidden moving targets through scattering samples. The temporal and angular correlations inherent in the scattered light pattern allowed us to reconstruct the hidden object in cases where multiply scattered light dominates over ballistic light. This paper presented a first proof of concept. Although we demonstrated imaging of binary-amplitude targets, our system can also be extended to imaging gray-scale targets [28]. Since our imaging technique utilizes the angular memory effect, it is scalable. Moreover, our method does not require access inside the scattering media and can therefore be used as a black box imaging system. With appropriate optimization, this opens up potential for use in applications involving the tracking of moving object in turbulent atmospheres, such as fog or underwater.

Funding
National Institutes of Health (NIH 1U01NS090577); GIST-Caltech Collaborative Research (CG2016); Natural Sciences and Engineering Research Council of Canada (NSERC PGSD3).

Appendix 1 -Deconvolving the speckle autocorrelation
To deconvolve the speckle autocorrelation (SAC), ∆I ∆I, Weiner deconvolution was applied to reduce the deconvolution noise. We briefly describe the process here. We can rewrite Eq. (13) as g = ∆I ∆I ≈ A * h + n = y + n where h(x i ) = 2δ(x i ) − C(∆x)δ(x i ± ∆x i ), A = O O, and n is the noise term. In this case, Weiner deconvolution estimates A by applying where F is the Fourier transform operator, and k = F(n) F(g) ≈ 1 S N R estimates the SNR level of your signal [29]. Since all object ACs have a peak value of A(x i = (0, 0)) = x O 2 , to determine h from the SAC, we estimated the value of C(∆x) by taking the negative/positive peak values in the SAC. The locations of the negative peaks, with respect to the centered, positive peak, provided the value of the shift ∆x i .