Benjamin T. Cox

k-Wave: MATLAB toolbox for the simulation and reconstruction of photoacoustic wave fields

A new, freely available third party MATLAB toolbox for the simulation and reconstruction of photoacoustic wave fields is described. The toolbox, named k-Wave, is designed to make realistic photoacoustic modeling simple and fast. The forward simulations are based on a k-space pseudo-spectral time domain solution to coupled first-order acoustic equations for homogeneous or heterogeneous media in one, two, and three dimensions. The simulation functions can additionally be used as a flexible time reversal image reconstruction algorithm for an arbitrarily shaped measurement surface. A one-step image reconstruction algorithm for a planar detector geometry based on the fast Fourier transform (FFT) is also included. The architecture and use of the toolbox are described, and several novel modeling examples are given. First, the use of data interpolation is shown to considerably improve time reversal reconstructions when the measurement surface has only a sparse array of detector points. Second, by comparison with one-step, FFT-based reconstruction, time reversal is shown to be sufficiently general that it can also be used for finite-sized planar measurement surfaces. Last, the optimization of computational speed is demonstrated through parallel execution using a graphics processing unit.

Two-dimensional quantitative photoacoustic image reconstruction of absorption distributions in scattering media by use of a simple iterative method

Photoacoustic imaging is a noninvasive biomedical imaging modality for visualizing the internal structure and function of soft tissues. Conventionally, an image proportional to the absorbed optical energy is reconstructed from measurements of light-induced acoustic emissions. We describe a simple iterative algorithm to recover the distribution of optical absorption coefficients from the image of the absorbed optical energy. The algorithm, which incorporates a diffusion-based finite-element model of light transport, converges quickly onto an accurate estimate of the distribution of absolute absorption coefficients. Two-dimensional examples with physiologically realistic optical properties are shown. The ability to recover optical properties (which directly reflect tissue physiology) could enhance photoacoustic imaging techniques, particularly methods based on spectroscopic analysis of chromophores.

Modeling nonlinear ultrasound propagation in heterogeneous media with power law absorption using a k-space pseudospectral method

The simulation of nonlinear ultrasound propagation through tissue realistic media has a wide range of practical applications. However, this is a computationally difficult problem due to the large size of the computational domain compared to the acoustic wavelength. Here, the k-space pseudospectral method is used to reduce the number of grid points required per wavelength for accurate simulations. The model is based on coupled first-order acoustic equations valid for nonlinear wave propagation in heterogeneous media with power law absorption. These are derived from the equations of fluid mechanics and include a pressure-density relation that incorporates the effects of nonlinearity, power law absorption, and medium heterogeneities. The additional terms accounting for convective nonlinearity and power law absorption are expressed as spatial gradients making them efficient to numerically encode. The governing equations are then discretized using a k-space pseudospectral technique in which the spatial gradients are computed using the Fourier-collocation method. This increases the accuracy of the gradient calculation and thus relaxes the requirement for dense computational grids compared to conventional finite difference methods. The accuracy and utility of the developed model is demonstrated via several numerical experiments, including the 3D simulation of the beam pattern from a clinical ultrasound probe.

Photoacoustic tomography with a limited-aperture planar sensor and a reverberant cavity

Biomedical photoacoustic tomography (PAT) is a soft-tissue imaging modality which combines the high spatial resolution of ultrasound (US) with the contrast and spectroscopic opportunities afforded by imaging optical absorption. Planar US arrays composed of piezoelectric or optical detector elements with small element sizes and fast acquisition times are readily available, making them an attractive option for imaging applications. An exact and efficient, FFT-based PAT reconstruction algorithm, that converts acoustic measurements recorded over a plane to a PAT image, is known. However, to capture sufficient data for an exact PAT reconstruction with a planar geometry requires an infinitely wide array. In practice it will be finite, resulting in a loss of resolution and introducing artefacts into the image. To overcome this limitation it is proposed that acoustic image sources, provided by enclosing the target in a reverberant cavity, are used to generate a periodically repeating sound field. Measurements of this periodic sound field can be used to reconstruct a PAT image exactly from measurements of reverberation made over a finite aperture. The existing FFT-based PAT reconstruction algorithm with only minor additional modifications can be used to generate the image in this case.

A gradient-based method for quantitative photoacoustic tomography using the radiative transfer equation

Quantitative photoacoustic tomography (QPAT) offers the possibility of high-resolution molecular imaging by quantifying molecular concentrations in biological tissue. QPAT comprises two inverse problems: (1) the construction of a photoacoustic image from surface measurements of photoacoustic wave pulses over time, and (2) the determination of the optical properties of the imaged region. The first is a well-studied area for which a number of solution methods are available, while the second is, in general, a nonlinear, ill-posed inverse problem. Model-based inversion techniques to solve (2) are usually based on the diffusion approximation to the radiative transfer equation (RTE) and typically assume the acoustic inversion step has been solved exactly. Here, neither simplification is made: the full RTE is used to model the light propagation, and the acoustic propagation and image reconstruction are included in the simulations of measured data. Since Hessian- and Jacobian-based minimizations are computationally expensive for the large data sets typically encountered in QPAT, gradient-based minimization schemes provide a practical alternative. The acoustic pressure time series were simulated using a k-space, pseudo-spectral time domain model, and a time-reversal reconstruction algorithm was used to form a set of photoacoustic images corresponding to four illumination positions. A regularized, adjoint-assisted gradient inversion using a finite element model of the RTE was then used to determine the optical absorption and scattering coefficients.

Artifact Trapping During Time Reversal Photoacoustic Imaging for Acoustically Heterogeneous Media

Several different reconstruction algorithms have been proposed for photoacoustic tomography, most of which presuppose that the acoustic properties of the medium are constant and homogeneous. In practice, there are often unknown spatial variations in the acoustic properties, and these algorithms give, at best, only approximate estimates of the true image. The question as to which approach is the most robust in these circumstances is therefore one of practical importance. Image reconstruction by ¿time reversal¿-using a numerical propagation model with a time-varying boundary condition corresponding to the measured data in reversed temporal order-has been shown to be less restrictive in its assumptions than most, and therefore a good candidate for a general and practically useful algorithm. Here, it is shown that such reconstruction algorithms can ¿trap¿ time reversed scattered waves, leading to artifacts within the image region. Two ways to mitigate this effect are proposed.

Reconstructing absorption and scattering distributions in quantitative photoacoustic tomography

Quantitative photoacoustic tomography is a novel hybrid imaging technique aiming at estimating optical parameters inside tissues. The method combines (functional) optical information and accurate anatomical information obtained using ultrasound techniques. The optical inverse problem of quantitative photoacoustic tomography is to estimate the optical parameters within tissue when absorbed optical energy density is given. In this paper we consider reconstruction of absorption and scattering distributions in quantitative photoacoustic tomography. The radiative transport equation and diffusion approximation are used as light transport models and solutions in different size domains are investigated. The simulations show that scaling of the data, for example by using logarithmic data, can be expected to significantly improve the convergence of the minimization algorithm. Furthermore, both the radiative transport equation and diffusion approximation can give good estimates for absorption. However, depending on the optical properties and the size of the domain, the diffusion approximation may not produce as good estimates for scattering as the radiative transport equation.

The frequency-dependent directivity of a planar fabry-perot polymer film ultrasound sensor

A model of the frequency-dependent directivity of a planar, optically-addressed, Fabry-Perot (FP), polymer film ultrasound sensor is described and validated against experimental directivity measurements made over a frequency range of 1 to 15 MHz and angles from normal incidence to 80deg. The model may be used, for example, as a predictive tool to improve sensor design, or to provide a noise-free response function that could be deconvolved from sound-field measurements in order to improve accuracy in high-frequency metrology and imaging applications. The specific question of whether effective element sizes as small as the optical-diffraction limit can be achieved was investigated. For a polymer film sensor with a FP cavity of thickness d, the minimum effective element radius was found to be about 0.9d, and that an illumination spot radius of less than d/4 is required to achieve it

Quantitative photoacoustic imaging: fitting a model of light transport to the initial pressure distribution

Photoacoustic imaging, which generates a map of the initial acoustic pressure distribution generated by a short laser pulse, has been demonstrated by several authors. Quantitative photoacoustic imaging takes this one stage further to produce a map of the distribution of an optical property of the tissue, in this case absorption, which can then be related to a physiological parameter. In this technique, the initial pressure distribution is assumed to be proportional to the absorbed laser energy density. A model of light transport in scattering media is then used to estimate the distribution of optical properties that would result in such a pattern of absorbed energy. The light model used a finite element implementation of the diffusion equation (with the delta-E(3) approximation included to improve the accuracy at short distances inside the scattering medium). An algorithm which applies this model iteratively and converges on a quantitative estimate of the optical absorption distribution is described. 2D examples using simulated data (initial pressure maps) with and without noise are shown to converge quickly and accurately.

Measurement of Broadband Temperature-Dependent Ultrasonic Attenuation and Dispersion Using Photoacoustics

The broadband ultrasonic characterization of biological fluids and tissues is important for the continued development and application of high-resolution ultrasound imaging modalities. Here, a photoacoustic technique for the transmission measurement of temperature-dependent ultrasonic attenuation and dispersion is described. The system uses a photoacoustic plane wave source constructed from a polymethylmethacrylate substrate with a thin optically absorbent layer. Broadband ultrasonic waves are generated by illuminating the absorbent layer with nanosecond pulses of laser light. The transmitted ultrasound waves are detected by a planar 7-mum high-finesse Fabry-Perot interferometer. Temperature- induced thickness changes in the Fabry-Perot interferometer are tracked to monitor the sample temperature and maintain the sensor sensitivity. The measured -6-dB bandwidth for the combined source and sensor is 1 to 35 MHz, with an attenuation corrected signal level at 100 MHz of -10 dB. The system is demonstrated through temperature-dependent ultrasound measurements in castor oil and olive oil. Power law attenuation parameters are extracted by fitting the experimental attenuation data to a frequency power law while simultaneously fitting the dispersion data to the corresponding Kramers-Kronig relation. The extracted parameters are compared with other calibration measurements previously reported in the literature.