## Abstract

We report on experimental demonstration of photoacoustic tomography for simultaneous recovery of the optical and acoustic properties of heterogeneous media. Photoacoustic images are reconstructed from a series of tissue-like phantom experiments using a finite element-based reconstruction algorithm coupled with a scanning photoacoustic imaging system. The results show that the images obtained are quantitative in terms of the shape, size, location and optical and acoustic properties of the heterogeneities examined.

©2006 Optical Society of America

## 1. Introduction

Non-invasive laser-induced photoacoustic tomography (PAT), as a hybrid imaging method, overcomes the resolution disadvantage of pure optical imaging and contrast/speckle disadvantages of pure ultrasound imaging [1–7]. It has been shown that PAT is a powerful imaging method for visualizing both the structural and functional information of biological tissues with excellent spatial resolution and satisfactory imaging depth [8, 9]. PAT is concerned with an inverse problem where a single short-pulsed light beam illuminates an object and the photoacoustic waves excited by thermo elastic expansion are measured using wideband ultrasound transducers in multiple locations around the object. The geometry of the object and spatial distribution of the absorbed optical energy density and acoustic property are recovered from the measured scattered fields using a reconstruction algorithm. To date various PAT reconstruction algorithms have been developed and applied to the detection of breast cancer, skin cancer, vascular diseases and brain tumors in small animals [10–14]. These algorithms, however, almost all rely on analytical solutions to the photoacoustic wave equation in a regularly shaped measurement configuration without appropriate boundary conditions applied. The major assumption made in these linear algorithms is that biological tissues are acoustically homogeneous which is not true in reality. To overcome these limitations, we have developed an iterative nonlinear reconstruction algorithm that is based on the finite element solution to the full Helmoholtz-like photoacoustic wave equation without the homogeneous acoustic property approximations [15]. Together with an iterative Newton method, our algorithm is able to precisely solve the Helmholtz wave equation and fulfill reliable inverse computation for an arbitrary measurement configuration [16–18]. Further, using simulated data [15] we show that our algorithm has an additional advantage compared to the existing ones: the images of both acoustic speed and absorbed optical energy density can be recovered simultaneously. It is known that there exist significant differences in acoustic properties between normal and tumor tissues and that tissue acoustic property images have the potential to differentiate benign from malignant breast lesions [19, 20]. Thus our PAT method offers a unique ability of reconstructing both tissue optical and acoustical properties in a single modality for a more complete diagnosis of diseases. In this paper, we demonstrate experimental evidence that the optical and acoustic images of heterogeneous media can indeed be obtained simultaneously using our reconstruction method.

## 2. Methods and materials

Since our reconstruction algorithm has been described in detail elsewhere [15–18], we only give a brief outline here for context. The frequency-domain wave equation for PAT imaging in an acoustically heterogeneous medium is written as [15],

where *p* is the acoustic pressure; *k*
_{0} = *ω*/*c*
_{0} is the wave number described by the angular frequency, *ω* and the speed of acoustic wave in a reference or coupling medium, *c*
_{0}; *β* is the thermal expansion coefficient; *C*_{p}
is the specific heat; * Φ* is the absorbed optical energy density that is the product of optical absorption coefficient and incident optical fluence;

*O*is a coefficient that depends on both acoustic speed and attenuation as follows [15]:

in which *c* is the speed of acoustic wave in the medium/tissue; *α* is the acoustic attenuation coefficient. The finite element discretization of Eq. (1) is stated as,

where the elements of the matrix * A* and

*are expressed as*

**b**$${B}_{i}=-{\mathit{ik}}_{0}\frac{{c}_{0}\beta}{{C}_{p}}\u3008\sum _{k}{\Phi}_{R,k}{\psi}_{k}{\psi}_{i}\u3009+{k}_{0}\frac{{c}_{0}\beta}{{C}_{p}}\u3008\sum _{l}{\Phi}_{I,l}{\psi}_{l}{\psi}_{i}\u3009$$

in which * Φ* ,

*O*and

*are expressed by their real and imaginary components,*

**p***ψ*is the basis function, <> indicates integration over the problem domain and the second-order absorption boundary conditions are employed here [16],

where $\eta =\frac{-\mathit{ik}-\frac{3}{2\rho}+\frac{i3}{{8\mathit{k\rho}}^{2}}}{\frac{1-i}{\mathit{k\rho}}},\gamma =\frac{\frac{-i}{{2\mathit{k\rho}}^{2}}}{\frac{1-i}{\mathit{k\rho}}}$ and *θ* is the angular coordinate at a radial position *ρ*. To obtain a matrix equation capable of inverse solution, a combined Marquardt and Tikhonov regularization scheme is used for the inverse calculation [15],

where **p***o* = (${p}_{1}^{o}$,${p}_{2}^{o}$,…,${p}_{M}^{o}$
)^{T}, * p^{c}* = (${p}_{1}^{c}$,${p}_{2}^{c}$,…,${p}_{M}^{c}$
)

^{T}, and ${p}_{i}^{o}$ and ${p}_{i}^{c}$ are observed and computed complex acoustic field data for i=1, 2…,

*M*boundary measurement locations,

*Δχ*is the update vector for the optical and acoustic properties , ℑ is the Jacobian matrix formed at the boundary measurement sites,

*λ*is a scalar and

*is the identity matrix. Thus here the image formation task is to update optical and acoustic property distributions via iterative solution of Eqs. (3) and (6) so that a weighted sum of the squared difference between computed and measured acoustic pressure can be minimized.*

**I**The experimental system used was based on a mechanical scanning mechanism [16–18]. In this system, pulsed light with incident fluence 10 mJ/cm^{2} well below the safety standard from a Nd: YAG laser (wavelength: 532nm, pulse duration: 3–6ns) were coupled into the phantom via an optical subsystem and generated acoustic pressure wave. A wide-bandwidth, 1 MHz transducer was used to receive the acoustic signals. The transducer and the phantom were immersed in a water tank. A rotary stage rotated the receiver relative to the center of the tank. One set of data was taken at 120 positions when the receiver was scanned circularly over 360°. The complex wavefield signal was first amplified by a preamplifier (Gain: 17dB, 0.03KHz~3MHz, Onda, Corporation, Sunnyvale, CA), and then amplified further by a Pulser/Receiver (GE Panametrics, Waltham, MA). A data acquisition board converted it into digital one that was fed to a PC. The radius of the receiver motion path was 25.4mm in the phantom experiments conducted.

## 3. Results and discussion

In the first two experiments, we embedded one glycerol powder-containing target (10mm in diameter) in a 30 mm-diameter solid cylindrical phantom. We then immersed the object-bearing solid phantom in a 50.8mm-diameter water background, as displayed in Fig. 1(a). The acoustic speed in the region of the target containing 10% (phantom 1) and 28% (phantom 2) glycerol was 1540 and 1620m/s, respectively, while the acoustic speed for the background was 1480m/s. For the third experiment, a 5mm- and a 10mm-diameter target (each contained 28% glycerol) were embedded into a 30 mm-diameter solid cylindrical phantom. Again the entire solid phantom was immersed in a 50.8mm-diameter water background, as shown in Fig. 1(b). In the final experiment, we embedded two nanoparticle-containing objects (3mm in diameter) in a 50 mm-diameter solid cylindrical phantom, as plotted in Fig. 1(c). Each target had 0.5nM nanoparticles. We then immersed the object-bearing solid phantom in a 110.6 mm-diameter water background. For all the experiments, the phantom materials used consisted of Intralipid as scatterer and India ink as absorber with Agar powder for solidifying the Intralipid and India ink solution. The optical absorption and scattering coefficients for the phantom were 0.004mm^{-1} and 1.0mm^{-1}, respectively. The optical absorption contrast between the target and phantom is about 10 for experiments 1–3, and 20 for experiment 4.

Figures 2(a)–2(b) show the reconstructed acoustic and optical images for experiment 1, while Figs. 2(c)–2(d) present the recovered acoustic and optical images for experiment 2. We see that the object in each case is clearly detected. By estimating the full width at half maximum (FWHM) of the acoustic property, as displayed in Fig. 3, we found that the recovered size of these objects ranged from 9.1 to 9.3 mm, in good agreement with the actual object size of 10mm. In our experiments, the 10 and 28% glycerol containing target used gave an acoustic speed of 1540 and 1620m/s, respectively. As can be seen from Fig. 3, we are able to obtain quantitatively resolved images in terms of the acoustic property as well as the position, shape and size of the objects.

The optical images obtained (i.e., absorbed optical energy density) are dependent on the incident optical fluence, thus giving just a relative measure of the optical absorption property of the media. It is noted the existing PAT reconstruction methods can image only the distribution of absorbed optical energy density that is the product of both the local optical absorption coefficient and fluence distribution within the irradiated medium. It is well known that it is the tissue absorption coefficient that directly correlates with tissue structural and functional information such as blood oxygenation. We have developed a reconstruction algorithm that is based on the finite element solution to the photoacoustic wave equation and photon diffusion equation for recovering the optical properties of heterogeneous media. For the developed algorithm, the absorption coefficient can be iteratively reconstructed from the recovered absorbed energy density and the optical fluence that is computed by the diffuse optical tomography model using the radiant exposure as the source [18].

Figures 4(a) and 4(b) display the reconstructed acoustic speed and optical images for the third experiment having two targets. By estimating the FWHM of the acoustic profiles plotted in Fig. 5, we found that the recovered size of the two objects was 4.6 and 9.8mm, respectively. Again, both size values were in good agreement with the exact values of 5 and 10mm. As shown in Fig. 4, the reconstructed acoustic properties for these cases are also accurate relative to the exact value, which is further confirmed by the quantitative analysis shown in Fig. 5.

Figures 6(a)–6(b) present the reconstructed acoustic and optical properties of two objects having 0.5nM nanoparticles within each of the objects. By estimating the full width half maximum (FWHM) of the acoustic speed plotted in Fig. 7, we found that the recovered diameter size of these objects ranges from 2.8mm to 3.3mm, in good agreement with the actual object size of 3mm. In addition, it is noted the imaging resolution in Fig. 4 is not as good as that in Fig. 6 due to the low optical contrast between the targets and phantoms for experiments 1–3. Further, from Figs. 2, 4 and 6, we can also see that the image quality of both optical and acoustic images appears to be similar when different acoustic contrasts between the target and the background were used.

In summary, we have demonstrated experimental evidence that both optical and acoustic property images can be obtained using our photoacoustic imaging method. This capability of reconstructing two entirely distinct tissue parameters in a single modality may bring in significant clinical value for decision-making.

## References and links

**1. **G. Paltauf, J. Viator, S. Prahl, and S. Jacques, “Iterative reconstruction algorithm for optoacoustic imaging,” J. Acoust. Soc. Am. **112**, 1536–1544(2002). [CrossRef] [PubMed]

**2. **S.J. Norton and T. Vo-Dinh, “Optoacoustic diffraction tomography: analysis of algorihms,” J. Opt. Soc. Am. A **20**, 1859–1866(2003). [CrossRef]

**3. **A.A. Karabutov, E. Savateeva, and A. Oraevsky, “Imaging of layered structures in biological tissues with opto-acoustic front surface transducer,” Proc. SPIE **3601**, 284–295(1999). [CrossRef]

**4. **R.A. Kruger, D. Reinecke, and G. Kruger, “Thermoacoustic computed tomography-technical considerations,” Med. Phys. **26**, 1832–1837(1999). [CrossRef] [PubMed]

**5. **Roy. G.M. Kollkman, John H.G.M. Klaessens, Erwin Hondebrink, Jeroen CW Hopman, Frits F M de Mul, Wiendelt Steenbergen, Johan M Thijssen, and Ton G van Leeuwen, “Photoacoustic determination of blood vessel diameter,” Phys. Med. Biol. **49**, 4745–4756(2004). [CrossRef]

**6. **R. I. Siphanto, K. K. Thumma, R.G.M. Kolkman, T.G. van Leeuwen, F.F.M. de mul, J.W. van Neck, L.N.A. van Adrichem, and W. Steenbergen, “Serial noninvasive photoacoustic imaging of neovascularization in tumor angiogenesis,” Opt. Express **13**, 89–95 (2005). [CrossRef] [PubMed]

**7. **B. Yin, D. Xing, Y. Wang, Y. Zeng, Y. Tan, and Q. Chen, “Fast photoacoustic imaging system based on 320-element linear transducer array,” Phys. Med. Biol. **49**, 1339–1346 (2004). [CrossRef] [PubMed]

**8. **J.A. Viator, G. Au, G. Paltauf, S. Jacques, S. Prahl, H. Ren, Z. Chen, and J. Nelson, “Clinical testing of a photoacoustic probe for port wine stain depth determination,” Lasers Surg. Med. **30**, 141–148(2002). [CrossRef] [PubMed]

**9. **X. Wang, Y. Xu, M. Xu, S. Yokoo, E. Fry, and L. Wang, “Photoacoustic tomography of biological tissues with high cross-section resolution: Reconstruction and experiment,” Med. Phys. **29**, 2799–2805 (2002). [CrossRef]

**10. **P. Liu, “The P-transform and photoacoustic image reconstruction,” Phys. Med. Biol. **43**, 667–674 (1998). [CrossRef] [PubMed]

**11. **K.P. Kostli, D. Frauchiger, J. Niederhauser, G. Paltauf, H. Weber, and M. Frenz, “Optoacoustic imaging using a three-dimensional reconstruction algorithm,” IEEE J. Sel. Top. Quantum Electron **7**, 918–923(2001). [CrossRef]

**12. **X. Wang, Y. Pang, G. Ku, X. Xie, G. Stoica, and L.V. Wang, “Non-invasive laser-induced photoacoustic tomography for structural and functional imaging of the brain in vivo,” Nat. Biotech. **21**, 803–806(2003). [CrossRef]

**13. **Joel J. Niederhauser, Michael Jaeger, Robot Lemor, Peter Weber, and Martin Frenz, “Combined ultrasound and optoacoustic system for real-time high-contrast vascular imaging in vivo,” IEEE Trans. Biomed. Eng. **24**, 436–440 (2005).

**14. **C.G.A. Hoelen, F.F. de Mul, R. Pongers, and A. Dekker, “Three-dimensional photoacoustic imaging of blood vessls in tissue,” Opt. Lett. **23**, 648–650 (1998). [CrossRef]

**15. **H. Jiang, Z. Yuan, and X. Gu, “Spatially varying optical and acoustic property reconstruction using finite element-based photoacoustic tomography,” J. Opt. Soc. Am. A **23**, 878–888 (2006). [CrossRef]

**16. **Z. Yuan, C. Wu, H. Zhao, and H. Jiang, “Imaging of small nanoparticle-containing objects by finite element-based photoacoustic tomography,” Opt. Lett. **30**, 3054–3056 (2005). [CrossRef] [PubMed]

**17. **Z. Yuan, H. Zhao, C. Wu, Q. Zhang, and H. Jiang, “Finite element-based photoacoustic tomography: phantom and chicken bone experiments,” Appl. Opt. **45**, 3177–3183 (2006). [CrossRef] [PubMed]

**18. **Z. Yuan and H. Jiang, “Quantitative photoacoustic tomography: recovery of optical absorption coefficient maps of heterogeneous media,” Appl. Phys. Lett. **88**, 231101 (2006). [CrossRef]

**19. **J.F. Greenleaf and R.C. Bahn, “Clinical imaging with transmissive ultrasound computerized tomography,” IEEE Trans. Biomed. Eng. **28**, 177–185 (1981). [CrossRef] [PubMed]

**20. **S. Schreiman, J.J. Gisvold, J.F. Greenleaf, and R.C. Bahn, “Ultrasound transmission computed tomography of the breast,” Radiol. **150**, 523–530 (1984).