## Abstract

We present a new approach for constructing optical phase-space-time-frequency tomography (OPSTFT) of an optical wave field. This tomography can be measured by using a novel four-window optical imaging system based on two local oscillator fields balanced heterodyne detection. The OPSTFT is a Wigner distribution function of two independent Fourier Transform pairs, i.e., phase-space and time-frequency. From its theoretical and experimental aspects, it can provide information of position, momentum, time and frequency of a spatial light field with precision beyond the uncertainty principle. Besides the distributions of *x* – *p* and *t* – *ω*, the OPSTFT can provide four other distributions such as *x* – *t*, *p* – *t*, *x* – *ω* and *p* – *ω*. We simulate the OPSTFT for a light field obscured by a wire and a single-line absorption filter. We believe that the four-window system can provide spatial and temporal properties of a wave field for quantum image processing and biophotonics.

© 2011 OSA

## 1. Introduction

Most of optical imaging methods are limited by the temporal and spatial resolutions because of the unwanted scattering light coupled into their detection systems. This limitation is due to the fundamental concept in the process of measurement, that is, the uncertainty principle. Better resolution in position will reduce the resolution in momentum (angle) of an optical imaging system. Similarly, better resolution in time domain will reduce the flexibilities of spectroscopic analysis on a physical object. According to the uncertainty principle, the position *x* and momentum *p*, time *t* and frequency *ω* of a spatial light field cannot be measured simultaneously with high resolution. However, the distribution of *x* and *p*, *t* and *ω* of the spatial light field can be measured simultaneously with high resolution by using two local oscillator fields. The use of two local oscillator fields in a balanced heterodyne detection scheme is also called two-window technique [1, 2]. The optical phase-space tomography (OPST) [1, 3] is Wigner distribution [4, 5] associated with a Fourier Transform pair of position and momentum (angle) coordinates. Wigner distribution associated with two independent Fourier Transform pairs such as position-momentum and time-frequency is called optical phase-space-time-frequency tomography (OPSTFT) or Wigner phase-space-time-frequency distribution. In this paper, we develop a four-window heterodyne detection scheme based on two local oscillator (LO) fields for measuring the OPSTFT. The OPSTFT, *𝒲* (*x*, *p*, *ω*, *t*), will offer the correlation information of *x*, *p*, *t*, and *ω* of a wave field through the distributions of *x* – *p*, *ω* – *t*, *x* –*t*, *p*–*t*, *x* – *ω*, and *p* – *ω*, where the two variables are plotted by fixing the other two variables.

In quantum optics, Wigner function is usually used to represent the quantum mechanical wave function or quantum state of a physical system because there is no quantum device can directly measure the wave function. The Wigner function can have negative values [6–8]. The negative value is a unique property of quantum mechanics because its implies that a particle cannot have definite position and momentum values at the same time or [*x̂, p̂*] = *i*. Existing approaches of reconstructing Wigner function [9, 10] for optical fields usually involve tomographic inversion (Radon transform) of the rotated form of Wigner distribution. Raymer [11,12] has pioneered the measurement of Wigner distributions for quadrature-field amplitude of non-classical state by using optical homodyne detection. The method involves tomographic inversion (Radon transform) of a set of measured probability distributions of quadrature amplitudes. The measurement method developed by Raymer has also been used to measure the Wigner distribution for transverse spatial state (spatial mode) [13,14] and time-frequency domain [15] of an optical electromagnetic field. In addition, Wigner function for spatial properties of a quantum state can be measured by using a parity-inverting Sagnac interferometer at the single-photon level [16].

Recently, residual spatial fluctuation in terms of small displacements and tilts (momentum) of a whole optical beam has been observed and used to exhibit EPR entanglement [17]. Quantum imaging [18] has played a central role in understanding spatial fluctuations in quantum regions. Controlling these quantum fluctuations can improve image resolution and beam positioning for targeting technology. Entanglement with a large number of modes such as images has been accomplished by using four-wave mixing in an atomic vapor [19]. The property of multimode squeezed light [20] which allows us to increase the sensitivity beyond the shot-noise limit, could create many interesting new applications in optical imaging, high-precision optical measurements, optical communications and optical information processing. The Wigner function contains sub-Planck phase-space structures [21,22], which can be used to detect small movement of a large object.

The main advantage of applying Wigner distribution in optical imaging and sensing is certain wave-particle features of Wigner distributions, related to optical coherence, can survive over distances that are large compared to the transport mean free path [23]. Another advantage is that Wigner phase-space distributions could bridge the gap between phenomenological transport equations and rigorous wave equation treatments. Since rigorous transport equations can be derived for Wigner distributions, they are essential for obtaining fundamental new insights into the nature of light propagation in multiple scattering media. Evolution equations for Wigner distributions, which include optical coherence scatterings, are generally non-local and are relatively unexplored. With suitable approximations, these non-local equations reduce to the usual radiative transport equations [24]. Establishing the physical relationship between Wigner distributions and the phenomenological specific intensity will impact most existing methods of imaging in multiple scattering media.

For a wave field varying in one spatial dimension and one spectra domain, *ℰ* (*x, ω*), the Wigner phase-space-time-frequency distribution is given by,

*dp*∫

*dt𝒲*(

*x*,

*p*,

*ω*,

*t*) = |

*ℰ*(

*x, ω*)|

^{2}, ∫

*dx*∫

*dω 𝒲*(

*x, p, ω, t*) = |

*ℰ*(

*p,t*)|

^{2}, and etc. Most important, Eq. (1) shows that the Wigner distribution is Fourier Transform related to the two-point mutual coherence function in position and spectra. Therefore, it is sensitive to the spatially and spectrally varying phase and amplitude of the field.

## 2. Theoretical approach: measurement of OPSTFT

In our previous work [2], we use two-window technique based on two local oscillator fields for directly measuring Kirkwood-Rihaczek (KR) position-momentum distribution. Then, the Wigner function is obtained through linear transformation of the KR distribution.

In this work, we develop a four-window heterodyne detection scheme based on two local oscillator (LO) fields for measuring the OPSTFT. As shown in Figs. 1(a) and (b), a local oscillator (LO) field is a phase coherent superposition of two fields i.e a focused and a collimated Gaussian fields. The focused LO Gaussian beam has spatial width of Δ*x _{a}* =

*a*and broad optical spectrum with bandwidth of Δ

*ω*=

_{a}*α*. The collimated (large) LO Gaussian beam has spatial width of Δ

*x*=

_{A}*A*and narrow optical spectrum with bandwidth of Δ

*ω*=

_{A}*β*. The position resolution

*a*is provided by the focused LO Gaussian beam. The momentum resolution 1/

*A*is provided by the collimated LO Gaussian beam. The purpose of this arrangement is to obtain independent control of position and momentum (angle) resolution such that the product of (Δ

*x*· Δ

_{a}*p*) =

_{A}*a*· 1/

*A*≤ 1 clearly surpasses the uncertainty principle limit as shown in the shaped area in Fig. 1(a). Simultaneously, the method can provide independent control of frequency (spectra) and time (path-length) resolution such that (Δ

*ω*·Δ

_{A}*t*) =

_{a}*β*·1/

*α*≤ 1 to surpass the uncertainty principle limit as shown in the shaped area in Fig. 1(b). The spectra-resolved resolution

*β*is provided by the narrowband collimated LO Gaussian beam. The path-resolved resolution 1/

*α*is provided by the broadband focused LO Gaussian beam. In other words, we use position (time) window to obtain position (temporal) information of an optical field and simultaneously use momentum (frequency) window to reject noises due to other propagating angles (frequencies) from coupling into the detection system.

A schematic setup of the four-window technique is shown in Fig. 2. In this scheme, optical phase-space-time-frequency tomography (OPSTFT) is measured by scanning transverse position (*d _{x}*) and transverse momentum (

*d*), and by tuning wavelength (

_{p}*d*) and path delay (

_{ω}*τ*). Our system employs balanced heterodyne detection of the probe signal field, which is overlapped with two strong local oscillator fields (LO). The beat amplitude

*V*is determined by the spatial and spectra overlapping of two local oscillator (LO) and signal (S) fields at the plane of the detector (

_{B}*Z*=

*Z*) as,

_{D}*E*

_{LO}_{(}

_{S}_{)}(

*x′, ω, z*) is used to represent the two local oscillator fields (signal field), respectively. The

_{D}*x′*denotes the transverse position at the detector plane. As shown in Fig. 2, when the LO fields are translated off-axis for a distance of

*d*and the tunable filter is tuned

_{x}*d*away from optical center frequency, the LO field will have its spatial and spectra arguments shifted in Eq. (2) as given by,

_{ω}*E*(

*x′, ω, z*), to the fields at the source planes,

_{D}*E*(

*x, ω, z*= 0), of lenses L1 and L2 as follow; (i) The LO and signal fields each will experience a spatially varying phase of ${e}^{-\frac{i{\mathit{\text{kx}}}^{2}}{2f}}$ after passing through the lens. (ii) Since the lens L2 in the signal beam is translated off-axis by

*d*and the path-length of the signal beam is delayed by

_{p}*τ*, the signal field will experience additional phase-shift of ${e}^{-\frac{ik{(x-{d}_{p})}^{2}}{2f}}$ and delay of

*e*

^{−}

*, respectively. (iii) The LO and signal fields propagate to the detector plane through a distance of*

^{iωτ}*d*=

*f*, which is the focus-length of lenses L1 and L2, and will experience the phase-shift of ${e}^{-\frac{ik{(x-{x}^{\prime})}^{2}}{2f}}$. From (i), (ii) and (iii), the LO and signal fields in Eq. (3) can be rewritten in terms of the input field

*z*= 0 as described in [2],

Then, the mean square heterodyne beat signal is finally obtained as given by,

*x*

_{∘}→

*x*and

*ω*

_{∘}→

*ω*. The |

*V*|

_{B}^{2}is proportional to the phase-space-time-frequency convolution integral of the Wigner distributions for the two local oscillator and the signal fields in the planes of the input lenses L1 and L2, respectively.

We will use a broadband light source for this experiment. The LO fields are engineered in the form as given by,

*γ*and

*ϕ*are the relative amplitude and phase of the two LO fields, respectively. The spectrum of LO fields are assumed to be Gaussian function. This can be accomplished by using a single mode fiber with a tunable bandpass Gaussian filter from Newport. The phase-dependent part of Wigner function for the LO takes the form

*A*≫

*a*and

*α*≫

*β*. Then, the range of integration for the momentum, position, frequency and time coordinates in Eq. (11) is limited by the signal field. In this scheme (similar electronic components as in Refs. [1,2]), the two LO fields which differ in frequency by

*ϕ*kHz are phase-locked. The signal field is modulated such that the heterodyne beat signals with the focused LO field and the collimated LO field are about Ω

*MHz and Ω*

_{ω}*MHz +*

_{ω}*ϕ*kHz, respectively. The root mean square beat amplitude is measured with an analog spectrum analyzer with a bandwidth of 100 kHz (

*> ϕ*kHz) centered at Ω

*MHz. The output of the spectrum analyzer is squared in real time with a low noise amplifier and an analog multiplier, then the amplified signal is sent to the lock-in-amplifier. Substituting Eq. (13) into Eq. (11) and changing notations (*

_{ω}*d*→

_{x}*x*, $\frac{k{d}_{p}}{f}\to {p}_{o}$,

_{o}*d*→

_{ω}*ω*,

_{o}*τ*→

*t*), we find that the in- and out-of phase quadrature amplitudes in the lock-in-amplifier are directly corresponding to the real and imaginary parts of the quantity,

_{o}*𝒦*(x, p,

*ω*, t) is the Kirkwood-Rihaczek phase-space-time-frequency distribution. Eq. (14) is readily inverted to yield the Wigner phase-space and time-frequency function or the OPSTFT of the signal field by a linear transformation. We obtain,

*S*and

_{R}*S*are the real and imaginary parts of Eq. (14), i.e. the in- and out-of-phase quadrature amplitudes, which are simultaneously measured.

_{I}## 3. Simulation of OPSTFT

#### 3.1. A thin wire

The measurement of OPSTFT based on two local oscillator fields is called four-window technique. In practice, the four independent variables *a*, 1/*A*, 1/*α* and *β* are chosen to be small to resolve scales of interest on the spatial and temporal properties of a wave field. Physical properties of an object can be extracted through measuring the OPSTFT of the scattered light field through the object. The OPSTFT can provide *x* – *p* and *ω* – *t* distributions, and other four distributions such as *ω* – *x*, *t* – *x*, *ω* – *p* and *t* – *p* distributions. These six distributions can provide new types of information for the light field under study. We numerically simulate the measurement of OPSTFT for the light field scattered through a thin wire with the diameter of 0.6 mm. We use a Gaussian beam with a wave field as given by,

*σ*and

_{x}*σ*are the spatial and temporal bandwidths. R and

_{t}*ω*are the radius of curvature and center/carrier frequency of the light field. We are interested in looking at the light scattered right after the wire, where the scattered light field,

_{c}*ℰ*(

_{wire}*x,t*), can be written as the product of Eq. (16) and a wire function (slitfun[x]=If[–0.3mm ≤

*x*≤ 0.3mm, 0.0, 1.0] in Mathematica program). We have used this field for exploring phase-space interference analog to superposition of two spatially separated coherent states [2]. However, the phase-space distribution associated with time-frequency distribution incorporated into Wigner distributions haven’t been explored. We use

*σ*= 0.85 mm,

_{x}*σ*=200 fs, and

_{t}*R*=−10000 mm and

*ω*= 0 for the simulation. We first obtain

_{c}*ℰ*(

_{wire}*p, ω*) by numerically Fourier transforming the

*ℰ*(

_{wire}*x,t*). This can be accomplished by numerically generating 20 points in position co-ordinate and 20 points in time coordinate from the product field of Eq. (16) and the wire function. Then, we generate the Kirkwood-Richaczek phase-space-time-frequency distribution, $\mathcal{K}(x,\hspace{0.17em}p,\hspace{0.17em}\omega ,\hspace{0.17em}t)={\mathcal{E}}_{\mathit{\text{wire}}}^{*}(x,\hspace{0.17em}t){\mathcal{E}}_{\mathit{\text{wire}}}(p,\hspace{0.17em}\omega )\text{exp}(ixp-i\omega t)$, as in Eq. (14) for the scattered field

*ℰ*(

_{wire}*x,t*). We first plot the real and imaginary parts of

*𝒦*(

*x, p,*0

*,*0) as shown in Fig. 3(a) and (b). We adopt the units for position (x) in mm, momentum (p) in rad/mm, frequency (

*ω*) in 10

^{13}Hz, and time in 10

^{−13}

*s*. The frequency (

*ω*) is plotted as frequency difference from the center frequency

*ω*. The time (t) is plotted as time difference from the zero delay (

_{c}*τ*= 0). The

*𝒦*(0, 0,

*ω, t*) is zero because there is zero beat signal (|

*V*|

_{B}^{2}) at the position x=0 and p =0. To explore the time-frequency distribution associated with nonzero beat signal (|

*V*|

_{B}^{2}) at the phase-space point(x=0.4, p=0), we plot the real and imaginary parts of

*𝒦*(0.4, 0,

*ω, t*) as shown in Figs. 3(c) and (d). From the

*𝒦*(

*x, p, ω, t*) distribution, we can obtain Wigner phase-space-time-frequency distribution by using linear transformation as in Eq. (15). We use Mathematica program for performing the numerical integration. Figure 4(a) is the plot of

*𝒲*(

*x, p*, 0, 0) where the tunable filter in the collimated LO beam is set at center frequency and the delay (

*τ*) in the signal beam is set to zero. The phase-space oscillation along the x=0 is due to the coherent phase-space interference of two spatially separated wave packets after the wire. This oscillation exhibits sub-Planck phase-space structure [21, 22] and negative values of Wigner function [6–8]. Now, we plot the

*𝒲*(0, 0,

*ω, t*) as shown in Fig. 4(b), where the phase-space point is at (0, 0). Note that the

*𝒦*(0, 0,

*ω, t*) is zero everywhere but not the Wigner function of

*𝒲*(0, 0,

*ω, t*). The reason is KR distribution is suitable for describing local properties of a wave function, while the Wigner function is suitable for describing wave properties of a particle wave function [2]. Figure 4(c) shows an interesting result of

*𝒲*(0, 2,

*ω, t*), which is negative, i.e, the inverse of Fig. 4(b). This is because the Wigner distribution of

*𝒲*(0, 2, 0, 0) has negative value. These properties are important to explore hyper-entanglement of intrinsic properties of single-photon spatial qubit states. Other four distributions such as

*𝒲*(0,

*p*, 0,

*t*),

*𝒲*(0,

*p*,

*ω,*0),

*𝒲*(

*x,*0, 0,

*t*), and

*𝒲*(

*x,*0,

*ω,*0) are plotted as shown in Figs. 4(d), (e), (f) and (g), respectively. The

*𝒲*(0,

*p*, 0,

*t*)) and

*𝒲*(0,

*p*,

*ω*, 0) exhibit oscillation/interference behavior along momentum coordinate. The phase-space interference due to the two spatially separated wave packets influenced the distribution of momentum (angle) in the time and spectra domains.

#### 3.2. An absorption filter

We numerically simulate the measurement of OPSTFT for a broadband light field passing through an absorption filter with bandwidth of 2 THz. We use a Gaussian beam with a Gaussian linewidth as given by,

*σ*is the spectra bandwidth. It is much easier to work on the spectra domain of the field so that the light field passing through the filter can be written as the product of Eq. (17) and an absorption filter function (slitfun[x]=If[–0.1 ≤

_{ω}*x*≤ 0.1, 0.0, 1.0] in Mathematica program). We use

*σ*= 0.85 mm,

_{x}*σ*=5.0 THz,

_{t}*R*=−10000 mm and

*ω*= 0 for the simulation. First, we generate the Kirkwood-Richaczek phase-space-time-frequency distribution, $\mathcal{K}(x,\hspace{0.17em}p,\hspace{0.17em}\omega ,\hspace{0.17em}t)={\mathcal{E}}_{\mathit{\text{filter}}}^{*}(x,\hspace{0.17em}\omega ){\mathcal{E}}_{\mathit{\text{filter}}}(p,\hspace{0.17em}t)\text{exp}(ixp+i\omega t)$, where the

_{c}*ℰ*(

_{filter}*p,t*) is obtained by numerically Fourier transformed the

*ℰ*(

_{filter}*x, ω*). The

*𝒦*(

*x*,

*p*, 0, 0) is zero because there is zero beat signal (|

*V*|

_{B}^{2}) at

*ω*= 0 and

*t*= 0. Figures 5(a) and (b) show the real and imaginary parts of

*𝒦*(0, 0,

*ω*,

*t*). To explore the position-momentum distribution associated with nonzero beat signal (|

*V*|

_{B}^{2}) at the time-frequency point(

*ω*=0.2, t=0), we plot the real and imaginary parts of

*𝒦*(0.2, 0,

*ω*,

*t*) as shown in Figs. 5(c) and (d). Since we have the

*𝒦*(

*x*,

*p*,

*ω*,

*t*) distribution, then we can obtain Wigner phase-space-time-frequency distribution by using linear transformation as in Eq. (15), by means of numerical integration.

Figure 6(a) is the plot of *𝒲* (0, 0, *ω*, *t*) where the position *d _{x}* of the mirror for the LO beam is set to zero and the momentum

*d*of the lens L1 in the signal beam is set to zero. The time-frequency oscillation along the

_{p}*ω*= 0 is due to the coherent time-frequency interference of two spectrally separated wave packets after the filter. This observation has been observed in Wigner time-frequency distribution [15]. We plot the

*𝒲*(

*x*,

*p*, 0, 0) as shown in Fig. 6(b), where the time-frequency point is at (0, 0). Figure 6(c) shows an interesting result of

*𝒲*(

*x*,

*p*, 0, 3), which is negative, i.e, the inverse of Fig. 6(b). This is because the Wigner distribution of

*𝒲*(0, 0, 0, 3) has negative value. Other four distributions such as

*𝒲*(0,

*p*, 0,

*t*),

*𝒲*(0,

*p*,

*ω*, 0),

*𝒲*(

*x*, 0, 0,

*t*), and

*𝒲*(

*x*, 0,

*ω*, 0) are plotted as shown in Figs. 6(d), (e), (f) and (g), respectively. The

*𝒲*(0,

*p*, 0,

*t*)) and

*𝒲*(

*x*, 0, 0,

*t*) exhibit oscillation/interference behavior along p =0 and x=0, respectively.

In real experiment, if we want to measure 10 points in position, 10 points in momentum, 10 points in frequency and 10 points in time, then we need to conduct 10 *×* 10 *×* 10 *×* 10 =10000 measurement points for constructing the real and imaginary parts of *𝒦* (*x, p, ω, t*) distribution as given in Eq. (14). Then, the Wigner distribution *𝒲* (*x*, *p*, *ω*, *t*) is obtained through linear transformation of KR distribution as given in Eq. (15), where the numerical integration is used to transform 10000 measurement values of the real and imaginary parts of KR distribution.

## 4. Conclusion

In conclusion, we have developed four-window optical heterodyne imaging technique based on two local-oscillator (LO) beams for measuring OPSTFT of a wave field. The four-window technique can independently control and simultaneously provide high resolution in position, momentum (angle), time, and spectra for characterizing spatial-temporal properties of a wave function. This newly developed optical phase-space-time-frequency tomography can be used for exploring hyper-entanglement (spatial-temporal) qubit states of photon wave function and for early disease detection in biophotonics.

## References and links

**1. **K. F. Lee, F. Reil, S. Bali, A. Wax, and J. E. Thomas, “Heterodyne measurement of Wigner distributions for classical optical fields,” Opt. Lett. **24**, 1370–1372 (1999). [CrossRef]

**2. **V. Bollen, Y. M. Sua, and K. F. Lee, “Direct measurement of the Kirkwood-Rihaczek distribution for the spatial properties of a coherent light beam,” Phys. Rev. A **81**, 063826 (2010). [CrossRef]

**3. **F. Reil and J. E. Thomas, “Observation of phase conjugation of light arising from enhanced backscattering in a random medium,” Phys. Rev. Lett. **95**, 143903 (2005). [CrossRef] [PubMed]

**4. **E. P. Wigner, “On the quantum correction for thermodynamic equilibrium,” Phys. Rev. **40**, 749–759 (1932). [CrossRef]

**5. **M. Hillery, R. F. O’Connell, M. O. Scully, and E. P. Wigner, “Distribution functions in physics: fundamentals,” Phys. Rep. **106**, 121–167 (1984). [CrossRef]

**6. **Ch. Kurtsiefer, T. Pfau, and J. Mlynek, “Measurement of the Wigner function of an ensemble of Helium atoms,” Nature **386**, 150–153 (1997). [CrossRef]

**7. **D. Leibfried, T. Pfau, and C. Monroe, “Shadows and mirrors: reconstructing quantum states of atom motion,” Phys. Today **51**(4), 22–28 (April 1998). [CrossRef]

**8. **K. Wodkiewicz and G. H. Herling, “Classical and non-classical interference,” Phys. Rev. A **57**, 815–821 (1998).

**9. **K. Vogel and H. Risken, “Determination of quasiprobability distributions in terms of probability distributions for the rotated quadrature phase,” Phys. Rev. A **40**, 2847–2849 (1989). [CrossRef] [PubMed]

**10. **J. Bertrand and P. Bertrand, “A tomographic approach to Wigner function,” Found. Phys. **17**, 397–405 (1987). [CrossRef]

**11. **D. T. Smithey, M. Beck, M. G. Raymer, and A. Faridani, “Measurement of the Wigner distribution and the density matrix of a light mode using optical homodyne tomography: application to squeezed states and the vacuum,” Phys. Rev. Lett. **70**, 1244–1247 (1993). [CrossRef] [PubMed]

**12. **A. I. Lvovsky and M. G. Raymer, “Continuous-variable optical quantum-state tomography,” Rev. Mod. Phys. **81**, 299–332 (2009). [CrossRef]

**13. **M. G. Raymer, M. Beck, and D. F. McAlister, “Complex wave-field reconstruction using phase-space tomography,” Phys. Rev. Lett. **72**, 1137–1140 (1994). [CrossRef] [PubMed]

**14. **D. F. McAlister, M. Beck, L. Clarke, A. Meyer, and M. G. Raymer, “Optical phase retrieval by phase-space tomography and fractional-order Fourier transforms,” Opt. Lett. **20**, 1181–1183 (1995). [CrossRef] [PubMed]

**15. **M. Beck, M. G. Raymer, I. A. Walmsley, and V. Wong, “Chronocyclic tomography for measuring the amplitude and phase structure of optical pulses,” Opt. Lett. **18**, 2041–2043 (1993). [CrossRef] [PubMed]

**16. **B. J. Smith, B. Killett, M. G. Raymer, I. A. Walmsley, and K. Banaszek, “Measurement of the transverse spatial quantum state of light at the single-photon level,” Opt. Lett. **30**, 3365–3367 (2005). [CrossRef]

**17. **K. Wagner, J. Janousek, V. Delaubert, H. Zou, C. Harb, N. Treps, J. F. Morizur, P. K. Lam, and H. A. Bachor, “Entangling the Spatial Properties of Laser Beams,” Science **321**, 541–543 (2008). [CrossRef] [PubMed]

**18. **M. I. Kolobiv, *Quantum Imaging*, 1st ed. (Springer, 2006).

**19. **V. Boyer, A. M. Marino, R. C. Pooser, and P. D. Lett, “Entangled Images from Four-Wave Mixing,” Science **321**, 544–547 (2008). [CrossRef] [PubMed]

**20. **M. I. Kolobov, “The spatial behavior of nonclassical light,” Rev. Mod. Phys. **71**, 1539–1589 (1999). [CrossRef]

**21. **W. H. Zurek, “Sub-Planck structure in phase space and its relevance for quantum decoherence,” Nature **412**, 712–717 (2001). [CrossRef] [PubMed]

**22. **F. Toscano, D. A. R. Dalvit, L. Davidovich, and W. H. Zurek, “Sub-Planck phase-space structures and Heisenberg-limited measurements,” Phys. Rev. A **73**, 023803 (2006). [CrossRef]

**23. **S. John, G. Pang, and Y. Yang, “Optical Coherence Propagation and Imaging in a Multiple Scattering Medium,” J. Biomed. Opt. **1**, 180–191 (1996). [CrossRef]

**24. **A. Ishimaru, *Wave Propagation and Scattering in Random Media* (Academic, 1978), Vols. I and II.