Open Access
Add to BibSonomyAbstract
We present a robust simulation technique to model the time-reversed ultrasonically encoded (TRUE) technique for deep-tissue imaging. The pseudospectral time-domain (PSTD) algorithm is employed to rigorously model the electromagnetic wave interaction of light propagating through a macroscopic scattering medium. Based upon numerical solutions of Maxwell’s equations, the amplitude and phase are accurately accounted for to analyze factors that affect the TRUE propagation of light through scattering media. More generally, we demonstrate the feasibility of modeling light propagation through a virtual tissue model of macroscopic dimensions with numerical solutions of Maxwell’s equations.
© 2014 Optical Society of America
1. Introduction
The penetration depth of fluorescence techniques have been limited to a depth of one transport mean free path, or approximately 1 mm in most biological samples [1]. Irregular scattering effect in biological tissues randomizes light amplitude and phase, resulting in a very limited penetration depth of light. Increasing the penetration depth of light through turbid media essentially increases the applicability of optical techniques for biomedical problems.
To achieve light delivery through biological tissues, Xu et al. proposed a technique, time-reversal of ultrasonically encoded light (TRUE) [2], which combines optical phase conjugation (OPC) with ultrasound encoding. An ultrasound pulse that is much less scattered is employed to modulate the illuminated light. The key idea is to overlap focal points of the ultrasound and light; by crossing the incident ultrasound and incident light beam, light within the overlapping focal region is modulated by the ultrasound. This overlapping region serves as a virtual light source that is frequency-shifted by the acousto-optic effect. This virtual light source plays the role of a guide star: after phase-conjugation, light propagates and scatters in reversed directions through the scattering medium towards the virtual light source. The combination of ultrasound-modulation and optical phase conjugation enables directing light through scattering medium to a specific position. As shown by [1, 3], the TRUE technique can achieve imaging through tissue with a penetration depth of 1 mm.
2. Methods
The OPC phenomenon is the key mechanism of the TRUE technique that directs light delivery to a target position through multiple scattering within a turbid medium. By yielding precise amplitude and phase information, a rigorous simulation of the OPC phenomenon of light propagation through a macroscopic scattering medium can shed light on the analysis of the TRUE technique. Recently, a two-dimensional model of the TRUE technique with dimensions 160-μm-by-100-μm using the finite-difference time-domain (FDTD) algorithm has been reported [4]. To adequately model the optical characteristics of biological tissue structures consisting of hundreds or thousands of biological cells, a simulation that can handle the TRUE technique for a problem of macroscopic dimensions is desired.
To study the TRUE phenomenon for light scattering through macroscopic turbid medium such as biological tissue structures, we report an innovative simulation technique based on the pseudospectral time-domain (PSTD) algorithm. The PSTD algorithm is computationally and memory-wise economic, hence advantageous to model large-scale light scattering problems. Simulation of the OPC simulation technique for a pulse illumination has been reported [5, 6]. However, most realistic biomedical optics involve a continuous-wave (CW) illumination; to date, a rigorous simulation that can handle the OPC phenomenon of a CW illumination upon scattering media with macroscopic dimensions has not been reported; previously reported technique [5, 6] falls short due to the tremendous computations involved.
In this manuscript, we present an innovative simulation technique to model the macroscopic problem of OPC. The simulation technique reported in this manuscript can model a macroscopic CW OPC phenomenon that was not possible before, enabling a robust simulation analysis of the TRUE technique. The simulation itself is computationally economic with minimal storage memory required. The numerical derivatives are calculated using the PSTD algorithm [7] which is advantageous for handling large-scale light scattering problems. In the PSTD algorithm, temporal derivatives of Maxwell’s equations are calculated using a 2nd-order central difference scheme, whereas the spatial derivatives are evaluated using Fourier transforms:
where E is the electric field, kx is the wave number, and F represents the Fourier transform. Since computer calculations are based on discrete numbers, discretization of the continuous electromagnetic fields is necessary. According to the Nyquist sampling theorem, the spatial derivatives calculated in Eq. (1), with a coarse grid of 2 spatial samples per wavelength, allows the PSTD technique to achieve similar accuracy as the FDTD technique, which requires 20 spatial samples per wavelength. The coarse spatial grid points makes it possible to simulate large-scale optical phenomena with economic computational memory. In order to model a light scattering system isolated in free space, an anisotropic perfectly matched layer (APML) absorbing boundary condition [8] is implemented to absorb all outgoing waves. If all outgoing waves never re-interact with the medium, the optical simulation is considered an isolated system. With the reported innovative simulation technique, the macroscopic, multiple light scattering phenomenon of the TRUE technique can accurately modeled for the first time. Furthermore, we analyze the effectiveness of light delivery through scattering media, which is the foundation of the TRUE technique.In the TRUE technique setup, a virtual light source is created by an ultrasound frequency-shifting technique [2] to differentiate light from a specific location and background noise. This can be modeled using a novel simulation method that we recently developed [9], which enables modeling a finite-width, CW light source in a PSTD simulation that was not possible before. By modeling the overlapping region of the ultrasound pulse and incident light beam as a localized light source embedded within the scattering medium, we can trace and analyze the light from the virtual source and simulate the phase-conjugated light focusing at the target position.
A two-dimensional PSTD simulation (Fig. 1) is implemented to model the TRUE technique. The virtual light source is surrounded by randomly positioned dielectric cylinders to model light propagation through scattering medium. We use a simulation grid with a uniform spatial resolution of 0.33 μm. The scattering medium is modeled by a rectangular region inside vacuum (n = 1.0) packed with randomly positioned, r-μm-diameter dielectric cylinders with a refractive index of n = 1.2. A standard APML absorbing boundary condition [8] is implemented to absorb outgoing waves, simulating an isolated light-scattering experiment. The TRUE virtual source can be modeled as a piecewise continuous electromagnetic field (Ei) generated from a localized region [9] by adding into the simulated space a finite-width, sinusoidal CW plane-wave at each time-step:
where (x0, y0) is the center of the light source; k and ω is the wavenumber and angular frequency, respectively. Edges of the added electromagnetic field are smoothened by a Gaussian profile in the y-direction to minimize numerical artifacts. A CW light of arbitrary wavelength can be generated from this rectangular region to model light from the virtual light source. To relate with experiments, a realistic wavelength such as 532 nm, 1 μm can be modeled. The optical field is added into space as a “soft source” [10] to prevent numerical artifacts from the source region. This allows modeling the TRUE virtual light source as a finite-width, CW light source where the incident light beam and ultrasound overlap in the scattering medium.
Fig. 1 PSTD simulation of the TRUE technique. In the TRUE technique, a virtual light source (yellow arrow) is formed inside the scattering medium where the incident light (yellow dashed lines) and incident ultrasound (purple dashed lines) overlaps and light is frequency-shifted by the ultrasound pulse. Light reaching the OPC region (blue shaded-area) is phase-conjugated and back-propagates towards the virtual source. The PSTD simulation region is surrounded by an APML absorbing boundary condition (pink shaded area) to model an isolated system.
In Fig. 2, the propagation of phase-conjugated light originated from a virtual light source with dimensions 5-μm-by-10- μm embedded inside a scattering medium is simulated. The virtual light source is modeled as a coherent, finite-width, CW plane wave propagating through a scattering medium. The scattering medium consists of 676 dielectric cylinders randomly positioned in space; the incident wavelength λ = 532 nm with a beam width of 10 μm. Light is incoherently scattered by the irregular geometry, resulting in a network-like sinuous optical paths. Light originated from the virtual light source is phase-conjugated and back-propagates through the scattering medium. Snapshots of light propagation are shown in Fig. 2: (a) 5000 time steps (0.25 ps), (b) 20000 time steps (1.00 ps), and (c) 40000 time steps (2.00 ps). As shown in the Fig. 2(d), by phase-conjugation, light can be delivered through the scattering medium to the location of the virtual light source.
Fig. 2 (Media 1) Simulation of the back-propagation of phase-conjugated light originating from a finite-width, CW plane wave virtual light source. A CW Gaussian beam (wavelength λ = 532 nm, beam width = 10 μm) impinging a scattering medium is simulated. The 200-μm-by-300-μm scattering medium consists of 676 randomly positioned, 5-μm-diameter dielectric (n = 1.2) cylinders. Snapshots at various instances depicting light propagation through the scattering medium are shown: (a) 5000 time steps (0.25 ps), (b) 20000 time steps (1.00 ps), and (c) 40000 time steps (2.00 ps). (d): a 30um-by-30um zoomed-in view of (c). (e) The intensity profile of the back-propagated light reaching the virtual light source (green line) is compared to light just emitted from the virtual light source (black thin line) with a root-mean-square error of 7.8%.
However, the light profile of the virtual light source may not be as simple as a coherent, uniform, finite-width plane wave. As light trespasses through a macroscopic scattering medium, the coherence of light gradually decreases as light is irregularly scattered into all directions. For weakly scattering media, the virtual light source may be suitably modeled by a coherent plane wave of finite extent. For thick scattering media, as light is multiply scattered, the direction of light propagation becomes oriented randomly; thus, the virtual light source formed by acousto-optic effect within the overlapped region may consist of speckle-like point light sources rather than a coherent, finite-width plane wave. Instead of a coherent plane wave, we model the virtual light source as a collection of 4 randomly positioned, speckle-like point light sources (each with a diameter of 0.66 μm) within a 10-μm-by-10- μm region. As shown in Fig. 3, light originated from the speckle-like virtual light source is phase-conjugated, back-propagates through the scattering medium back to the virtual light source. Snapshots of light propagation are shown in Fig. 3: (a) 5000 time steps (0.25 ps), (b) 20000 time steps (1.00 ps), and (c) 40000 time steps (2.00 ps). A zoomed-in view of light successfully delivered to the virtual light source is shown in 3(d).
Fig. 3 (Media 2) Simulation of the propagation of phase-conjugated light originating from the speckle-like virtual light source. The speckle-like virtual source is modeled as an aggregate of 4 randomly positioned 0.66-μm-diameter point light sources. The 200-μm-by-300-μm scattering medium consists of 676 randomly positioned, 5-μm-diameter dielectric (n = 1.2) cylinders. Snapshots at various instances depicting light propagation through the scattering medium are shown: (a) 5000 time steps (0.25 ps), (b) 20000 time steps (1.00 ps), and (c) 40000 time steps (2.00 ps). (d): a 30um-by-30um zoomed-in view of (c). (e) The intensity profile of the back-propagated light reaching the virtual light source (blue line) is compared to light just emitted from the virtual light source (black thin line) with a root-mean-square error of 9.0%.
Next, the effectiveness of light delivery of the speckle-like virtual light source is compared to a single, coherent plane wave virtual light source for different wavelengths. The coherent plane wave virtual light source of various cross-sectional width (w) is modeled: w = 0.67 μm, 2 μm, 5, μm, 10 μm, and 20 μm. The back-propagation of light is simulated at various wavelengths: 532 nm, 694 nm, 1064 nm, up to 3 μm (frequency ranging from 563 THz to 100 THz). The intensity profile of the back-propagated light reaching the virtual light source is compared to the light intensity profile of just emitted from the virtual light source; the root-mean-square (RMS) error is calculated and shown in Fig. 4. Simulation results show that the back-propagation of light through a scattering medium is insensitive to the dimensions of the virtual light source—the TRUE virtual light source can be modeled by a coherent finite-width plane wave or a speckle-like virtual light source with no significant difference in performance.
Fig. 4 Light delivery of various virtual light sources at various frequencies (100 THz to 563 THz) is compared. We simulate TRUE light propagation through a 200-μm-by-300-μm scattering medium consisting of 676 randomly positioned, 5-μm-diameter dielectric (n = 1.2) cylinders. Light delivery from a speckle-like virtual light source is compared to a single, coherent plane wave light source of various widths: w = 0.67 μm, 2 μm, 5, μm, 10 μm, and 20 μm. The root-mean-square (RMS) error of back-propagated light is compared to the original emitted light. Simulation results show insignificant differences for the effectiveness of the speckle-like virtual light source and a coherent, finite-width plane wave virtual light source.
3. Analyzing the density of scattering media
In the following sections, we model the virtual light source as a CW, finite-width, plane wave to facilitate quantitative analysis. Back-propagation of light towards the target position within a scattering medium is simulated; factors upon which the effectiveness of TRUE light delivery depends are analyzed. As shown in Fig. 5, light propagating through a scattering medium is randomized by the irregular geometry. To analyze this problem, we simulate a localized light source embedded inside a 200-μm-by-300-μm scattering medium consisting of N randomly positioned, 5-μm-diameter, dielectric (n = 1.2) cylinders. Scattering characteristic is proportional to the difference of the refractive index between the scatterer and the surrounding medium; for biological structures, the refractive index difference is typically less than 1.2. The scattering medium is modeled by an aggregate of dielectric (n = 1.2) cylinders positioned in vacuum (n = 1); the reduced scattering coefficient μs’ = 0.0113 (1/ μm). Arbitrary wavelength of light can be modeled with the reported simulation technique; here we model an incident wavelength of 1 μm. By phase-conjugating the outgoing light recorded in the OPC region at the right edge of the simulation grid (x = 200 μm, y = 0 to 300 μm), light back-propagates towards the virtual light source. For denser scattering media, scattered light spreads into a wider range of directions. Intensity pattern of light propagating through the scattering medium of various densities are shown in Fig. 5(a), 5(b), 5(c).
Fig. 5 Simulation of a CW, phase-conjugated light (wavelength λ = 1 μm) back-propagating (from right to left) through a scattering medium of various densities. The 200-μm-by-300-μm scattering medium consists of N randomly positioned, 5-μm-diameter dielectric (n = 1.2) cylinders: (a) N = 200, (b) N = 400, (c) N = 600. The phase-conjugated light back-propagates through the scattering medium and converges at the virtual light source. Intensity of light back-propagating through the scattering medium is depicted.
To quantitatively analyze the delivery of light through scattering medium, the cross-sectional amplitude profile recorded at the right edge of the virtual light source is shown in Fig. 6. From bottom to top are five profiles each corresponding to a 200-μm-by-300-μm scattering medium consisting of N = 0, 100, 200, 600, and 676 randomly-positioned, 5-μm-diameter, dielectric (n = 1.2) cylinders. For the case of no scatterers (N = 0), the amplitude profile exhibits a high-frequency wriggle resulting from near field diffraction of light just emerging from the virtual light source. With scatterers randomly positioned in space (N ≠ 0), this high-frequency wriggle is dominated by the interference of light from the irregular geometry and becomes less pronounced. As the scattering medium becomes denser with more cylinders, light is scattered into a wider range of directions. Furthermore, an increased percentage of light is backscattered away, resulting in less light delivered into the scattering medium [6]. Nevertheless, width of the central peak remains unchanged regardless of the density of the scattering medium.
Fig. 6 Varying the density of the scattering medium. The amplitude profile of the phase-conjugated light reaching the target position (as shown in Fig. 5) is compared for various N (number of dielectric cylinders within the scattering medium). Each profile is normalized by its maximum value; from bottom to top: N = 0, 100, 200, 600, and 676. For a denser scattering medium, the amount of light reaching the target position decreases.
4. Radius of constituent scatterers
Next, we vary the radius of the constituent scatterers of the scattering media. If the length-scale of the constituent scatterers affects the TRUE light propagation, varying radius would show a change in the delivery of light back to the target position. We model TRUE light propagation through various scattering media, each consisting of uniform r-μm-radius dielectric cylinders (refractive index n = 1.2). Cross-sectional intensity profiles at the target position are compared (Fig. 7. inset-figure); the root-mean-square error is calculated and compared for various radius r. Simulation results indicate that the effectiveness of the delivery of light through a scattering medium is insensitive to the size constituent scatterers.
Fig. 7 Varying radius of constituent scatterers of a scattering medium. We simulate TRUE light propagation through 600-μm-by-600-μm scattering media consisting of uniform scatterers with refractive index n = 1.2 and radius r, r ranging from 1.2 to 18 μm. The amplitude profile of the phase-conjugated light reaching the target position is compared to the original intensity profile of light emitted from the virtual light source. The root-mean-square error for various radius r is calculated. The error is insensitive to the size of the scatterers.
5. Fraction of phase-conjugated light
In theory, optimal performance can be achieved if all scattered light is time-reversed by optical phase conjugation. However, in practice it is impractical to phase-conjugate all scattered light; typically a fraction of the scattered light is phase-conjugated. By means of PSTD simulation, we vary the amount of back-propagated light and analyze light delivery to the target position of the TRUE technique. The amount of light undergoing phase-conjugation is determined by light reaching the OPC region (as depicted in Fig. 1). With a larger OPC region, more light is phase-conjugated and propagates towards the target position. We investigate the effect of different cross-sectional widths of the OPC region. As shown in Fig. 8 from bottom to top, each profile corresponds to different cross-sectional width of the OPC region: 50 μm, 100 μm, 200 μm, and 300 μm; as more light is phase-conjugated, the focal point width doesn’t change but the signal-to-background ratio increases.
Fig. 8 Varying the width of the OPC region. A CW Gaussian beam (wavelength λ = 1 μm, beam width = 10 μm) impinging a scattering medium is simulated. The percentage of scattered light being phase-conjugated is increased by increasing the span of the OPC region. The time-reversed light profile reaching the target position is compared. To compare the width, each profile is normalized by its maximum value. Width of the OPC region from bottom to top is: 50 μm, 100 μm, 200 μm, and 300 μm. With increased amount of phase-conjugated light, the signal-to-background ratio increases.
6. Modeling TRUE delivery of light through a virtual tissue model
Next we simulate light propagation through a virtual tissue model. Constructed based on the refractive index of a HT29 cell [11], the 2-D virtual tissue model consists of 670 biological cells randomly oriented within a 1000-μm-by-1000-μm region; each cell is approximately 30 μm in diameter. By modeling a virtual light source embedded at the center of the tissue model, we can analyze delivery of light through biological turbid media. CW light with a wavelength of 1 μm is emitted from a 5- μm-by-10-μm virtual light source. In this simulation, light reaching the OPC region at the right edge of the simulation is phase-conjugated and back-propagates towards the virtual light source. Maximum amplitude of the phase-conjugated light emerges at the center of the tissue model. A zoomed-in view of the center of the tissue model is shown in the inset-figure where light converges at the virtual light source and continues to propagate outwards. As shown in Fig. 9, via multiple scattering, the phase-conjugated light propagates through the virtual tissue model forming a pattern of optical paths similar to a river network where smaller creeks flow into larger rivers. Figure 9 demonstrates that the TRUE phenomenon for biological tissue structure can be analyzed with the reported simulation technique. Furthermore, arbitrary irregular geometry, such as a denser or sparser distribution of biological cells can be readily modeled. In addition, the reported simulation enables placing the CW virtual light source at arbitrary position (inside or outside) the scattering medium that was not possible before [4]; this can be very useful for modeling realistic problems. To date, the reported simulation method is the only method that can rigorously model the TRUE technique on a macroscopic scale without heuristic approximations.
Fig. 9 (Media 3) (Top): Simulation of TRUE light propagating through a 1000-μm-by-1000-μm virtual tissue model consisting of ~670 cells. The 2-D virtual tissue model is constructed based upon a HT29 cell. (Bottom): The time-reversed light profile (playback scenario) is compared with the forward scenario light profile. Both field profiles were measured at the center of the simulation space. (Inset-figure): A 80-μm-by-80-μm zoomed-in view of back-propagated light converging at the virtual light source.
With currently available computational resources, a 1-mm-by-1-mm simulation with a resolution of 0.333 μm can be completed within a reasonable timeframe. This simulation takes ~48 hours to run on a 12-core Intel x5650 processor. To model optical properties of finer structures, higher resolution (finer grid) is required, which involves more computations and longer simulation time. For example, to model a 1-mm-by-1-mm scattering medium with a resolution of 0.15 μm, a run time of 24*8 = 192 hours is required. The reported simulation technique can model the TRUE phenomenon for arbitrary, irregular biological tissue structures consisting of hundreds or even thousands of cells with the desired accuracy.
7. Summary
In this manuscript, we report a rigorous PSTD simulation analysis of the TRUE technique, where the frequency-shifted light is phase-conjugated and back-propagates through the scattering medium towards the target position (whereabouts of the virtual light source). The reported simulation technique enables modeling the CW virtual light source at arbitrary position (inside or outside) the scattering medium. Our simulation is based upon numerical solutions of Maxwell’s equations; by accounting for the phase and amplitude, we accurately simulate the light propagation and interference of the TRUE technique. Specifically, we analyze possible factors that affect light delivery of the TRUE technique, including the virtual light source, the density of the scattering media, the size of the constituent scatterers, and the fraction of phase-conjugated light. Simulation results show that the effectiveness of TRUE light delivery is insensitive to the morphology of the virtual light source (Fig. 3), the density of scattering media (Fig. 6) or the size of the constituent scatterers (Fig. 7). However, the amount of phase-conjugated light plays a critical role of light delivery—with increased phase-conjugated light, the signal-to-background ratio of the TRUE technique increases.
Furthermore, we demonstrate a robust modeling the TRUE light propagation through a virtual tissue model with macroscopic dimensions that has not been achieved before. With the specific amplitude and phase, light can be directed through a 1-mm-thick biological tissue and focus at the target position (Fig. 9). More generally, our results show that such macroscopic light scattering problem can be rigorously analyzed based upon numerical solutions of Maxwell’s equations.
Acknowledgments
This research is supported by the Taiwan National Science Council Grant NSC 99-2112-M-002-016-MY3, and, Devices for Solar Energy, Light Sources, Display, and Sensing under Excellent Research Projects of National Taiwan University grant NTU-ERP-103R89086. The authors thank Yingmin Huang, Benjamin Judkewitz and Changhuei Yang for insightful discussions.
References and links
1. Y. M. Wang, B. Judkewitz, C. A. DiMarzio, and C. H. Yang, “Deep-tissue focal fluorescence imaging with digitally time-reversed ultrasound-encoded light,” Nat. Commun. 3, 928 (2012).
2. X. Xu, H. Liu, and L. V. Wang, “Time-reversed ultrasonically encoded optical focusing into scattering media,” Nat. Photonics 5(3), 154–157 (2011). [CrossRef] [PubMed]
3. K. Si, R. Fiolka, and M. Cui, “Fluorescence imaging beyond the ballistic regime by ultrasound pulse guided digital phase conjugation,” Nat. Photonics 6(10), 657–661 (2012). [CrossRef] [PubMed]
4. J. L. Hollmann, R. Horstmeyer, C. Yang, and C. A. DiMarzio, “Analysis and modeling of an ultrasound-modulated guide star to increase the depth of focusing in a turbid medium,” J. Biomed. Opt. 18(2), 025004 (2013). [CrossRef] [PubMed]
5. S. H. Tseng, “PSTD Simulation of optical phase conjugation of light propagating long optical paths,” Opt. Express 17(7), 5490–5495 (2009). [CrossRef] [PubMed]
6. S. H. Tseng, “Investigating the Optical Phase Conjugation Reconstruction Phenomenon of Light Multiply Scattered by a Random Medium,” IEEE Photon. J. 2(4), 636–641 (2010). [CrossRef]
7. Q. H. Liu, “Large-scale simulations of electromagnetic and acoustic measurements using the pseudospectral time-domain (PSTD) algorithm,” IEEE Trans. Geosci. Rem. Sens. 37(2), 917–926 (1999). [CrossRef]
8. S. D. Gedney, “An anisotropic perfectly matched layer-absorbing medium for the truncation of FDTD lattices,” IEEE Trans. Antenn. Propag. 44(12), 1630–1639 (1996). [CrossRef]
9. Y. Huang, C. Tsai, W. Ting, and S. H. Tseng, “PSTD simulation of the continuous-wave optical phase conjugation phenomenon,” In Review.
10. A. Taflove and S. C. Hagness, Computational Electrodynamics: the finite-difference time-domain method (Artech House, 2000).
11. W. Choi, C. Fang-Yen, K. Badizadegan, S. Oh, N. Lue, R. R. Dasari, and M. S. Feld, “Tomographic phase microscopy,” Nat. Methods 4(9), 717–719 (2007). [CrossRef] [PubMed]
