Ashkan Javaherian

Ashkan Javaherian, Sean Holman

Inspired by the recent advances on minimizing nonsmooth or bound-constrained convex functions on models using varying degrees of fidelity, we propose a line search multigrid (MG) method for full-wave iterative image reconstruction in photoacoustic tomography (PAT) in heterogeneous media. To compute the search direction at each iteration, we decide between the gradient at the target level, or alternatively an approximate error correction at a coarser level, relying on some predefined criteria. To incorporate absorption and dispersion, we derive the analytical adjoint directly from the first-order acoustic wave system. The effectiveness of the proposed method is tested on a total-variation penalized Iterative Shrinkage Thresholding algorithm (ISTA) and its accelerated variant (FISTA), which have been used in many studies of image reconstruction in PAT. The results show the great potential of the proposed method in improving speed of iterative image reconstruction.

Ashkan Javaherian, Sean Holman

Quantitative photo-acoustic tomography (QPAT) seeks to reconstruct a distribution of optical attenuation coefficients inside a sample from a set of time series of pressure data that is measured outside the sample. The associated inverse problems involve two steps, namely acoustic and optical, which can be solved separately or as a direct composite problem. We adopt the latter approach for realistic acoustic media that possess heterogeneous and often not accurately known distributions for sound speed and ambient density, as well as an attenuation following a frequency power law that is evident in tissue media. We use a Diffusion Approximation (DA) model for the optical portion of the problem. We solve the corresponding composite inverse problem using three total variation (TV) regularised optimisation approaches. Accordingly, we develop two Krylov-subspace inexact-Newton algorithms that utilise the Jacobian matrix in a matrix-free manner in order to handle the computational cost. Additionally, we use a gradient-based algorithm that computes a search direction using the L-BFGS method, and applies a TV regularisation based on the Alternating Direction Method of Multipliers (ADMM) as a benchmark, because this method is popular for QPAT and direct QPAT. The results indicate the superiority of the developed inexact Newton algorithms over gradient-based Quasi-Newton approaches for a comparable computational complexity.

Ashkan Javaherian, Sean Holman

We present an optimization framework for photo-acoustic tomography of brain based on a system of coupled equations that describe the propagation of sound waves in linear isotropic inhomogeneous and lossy elastic media with the absorption and physical dispersion following a frequency power law using fractional Laplacian operators. The adjoint of the associated continuous forward operator is derived, and a numerical framework for computing this adjoint based on a k- space pseudospectral method is presented. We analytically show that the derived continuous adjoint matches the adjoint of an associated discretised operator. We include this adjoint in a first-order positivity constrained optimization algorithm that is regularized by total variation minimization, and show that the iterates monotonically converge to a minimizer of an objective function, even in the presence of some error in estimating the physical parameters of the medium.

Ashkan Javaherian, Seyed Kamaledin Setarehdan

The acoustic wave equation governs wave propagation induced by either volumetric radiation sources, or by surface sources of monopole or dipole type. For surface sources, boundary value problems yield wavefield representations via the Kirchhoff-Helmholtz or Rayleigh-Sommerfeld integrals. This study begins by establishing an equivalence between the analytic expressions of the associated monopole and dipole integral formulations and their full-waveform approximations. Leveraging this equivalence, we introduce reception operators that map free space pressure wavefields-obtained by solving the wave equation-onto measured fields restricted to the boundary. Building on this trace mapping, we derive the adjoint of the forward operator. We show that, under the common practical assumption of Dirichlet-type boundary data, the adjoint operator coincides-up to a constant factor-with the interior-field time-reversed form of the dipole integral formula, evaluated on the receiver surfaces. This study aims to advance the approximation of forward problems and the solution of inverse problems in acoustics, with a particular focus on applications that require accurate amplitude modeling, including therapeutic ultrasound optimization, attenuation reconstruction, and photoacoustic tomography.

Ashkan Javaherian

We extend the forward-adjoint operator framework derived in our previous study to photoacoustic tomography (PAT). In that earlier work, the acoustic forward operator included a reception operator that maps, at each time step, the pressure wavefield in free space onto the boundary (receiver surface). It was shown that this reception operator serves as a left-inverse of an emission operator that maps the pressure restricted to the boundary (emitter surface) onto free space, perfectly complying with the reciprocity law of physics. In this study, we define the full PAT forward operator as a composite mapping composed of an acoustic forward operator equipped with a scaled variant of the previously proposed reception operator, and an operator describing the photoacoustic source. Singularities arising both in the reception step (due to the boundary restriction) and in the photoacoustic source (due to its instantaneous nature) are regularized using regularized Dirac delta distributions. The resulting PAT forward-adjoint operator pair satisfies an inner-product relation, which we verify through numerical experiments on a discretized domain. The effectiveness of the proposed operator pair is further demonstrated using an iterative minimization framework that yields both qualitatively and quantitatively accurate reconstructions of an initial pressure distribution from the corresponding Dirichlet-type boundary data.