Translator Disclaimer
2 January 2018 Quantitative photoacoustic imaging of two-photon absorption
Author Affiliations +
Abstract
Photoacoustic tomography (PAT) is a hybrid imaging modality where we intend to reconstruct optical properties of heterogeneous media from measured ultrasound signals generated by the photoacoustic effect. In recent years, there have been considerable interests in using PAT to image two-photon absorption, in addition to the usual single-photon absorption, inside diffusive media. We present a mathematical model for quantitative image reconstruction in two-photon photoacoustic tomography (TP-PAT). We propose a computational strategy for the reconstruction of the optical absorption coefficients and provide some numerical evidences based on synthetic photoacoustic acoustic data to demonstrate the feasibility of quantitative reconstructions in TP-PAT.

1.

Introduction

Photoacoustic imaging (PAT) is a recent biomedical imaging modality that can provide high-resolution images of optical contrast of heterogeneous media such as biological tissues.113 In a typical PAT experiment, a short pulse of near-infrared (NIR) photons is sent into the medium to be probed. Photons then propagate inside the medium and a portion of them gets absorbed during the propagation process. The energy absorbed by the medium leads to the heating of the medium, and the heating then forces the medium to expand. The medium cools down when the remaining photons exit, and the cooling process leads to the contraction of the medium. The expansion and contraction of the medium initializes a pressure change inside the medium, which then propagates in the form of ultrasound waves. The ultrasound signals on the surface of the medium are then measured. These measurements are used to infer information on the optical properties of the medium.1431

We study here the problem of using PAT to image two-photon absorption properties of tissue-like heterogeneous media.3242 Here by “two-photon absorption,” we mean the absorption event where an electron transfers to an excited state after simultaneously absorbing two photons whose total energy exceed the electronic energy band gap.4345 Even though it occurs less frequently in normal biological tissues than its single-photon counterpart (i.e., the absorption event where an electron transfers to an excited state after absorbing the energy of a single photon), two-photon absorption is extremely useful in practice. In recent years, various types of materials with strong two-photon absorptions have been proposed and engineered as exogenous contrast agents for different optical imaging modalities.4649 Many such materials can be tuned to be associated with specific molecular signatures. Therefore, they can be used to visualize particular cellular functions and molecular processes inside biological tissues.

There have been extensive experimental investigations on measuring two-photon absorption properties of various materials using PAT.33,34,3942,50,51 Even though many of these studies use special experimental techniques, they do show collectively the feasibility of experimentally detecting contributions from two-photon absorptions to measured photoacoustic signals. That is, these studies showed that it is indeed possible to have strong enough photoacoustic effect from two-photon absorption that can be experimentally detected. However, it has not been satisfactorily demonstrated so far how to, if it is possible at all, separate the photoacoustic effect due to single-photon absorption from that due to two-photon absorption. In the rest of this paper, we demonstrate, computationally, through a model-based reconstruction algorithm, that it is possible to get quantitative reconstructions of both single-photon and twophoton absorptions and therefore separate them, if indeed the two-photon absorption contribution to the photoacoustic effect can be detected in the measured acoustic data.

2.

Mathematical Models

The main physical processes involved in a two-photon photoacoustic tomography (TP-PAT) experiment are the propagation of NIR photons and the propagation of ultrasound signals in the underlying medium. In optically heterogeneous media such as the biological tissues, it is well established now that the propagation of NIR photons can be modeled with a diffusion equation for the local density of photons.17,23,30,52 The main difference between TP-PAT and the regular PAT is that two-photon absorption, in addition to single-photon absorption, needs to be considered in the model for light propagation. Let us denote by ΩRd (d2) the medium to be probed and denote by u(x) the density of photons at position xΩ. We then have that u(x) solves the following semilinear diffusion equation:53

Eq. (1)

·γ(x)u(x)+σ(x)u(x)+μ(x)|u|u(x)=0,in  Ωu+κγuν=g(x),on  Ω,
where is the usual gradient operator with respect to the spatial variable x, and the function g(x) models the incoming NIR illumination source on the boundary Ω. The function γ(x) is the diffusion coefficient of the medium, and the functions σ(x) and μ(x) are the single-photon absorption and the (intrinsic) two-photon absorption coefficients, respectively. The unit outer normal vector at point x on the boundary Ω is denoted by ν(x), and the notation uν=ν·u is used in the Robin boundary condition. The coupling parameter κ in the boundary condition is the rescaled extrapolation length. Its value depends on many parameters and can be explicitly calculated in specific settings.54

The main difference between the nonlinear diffusion equation (1) and the classical linear diffusion equation for modeling light propagation in PAT is the extraterm μ(x)|u|u(x) that models the two-photon absorption mechanism. It is not the objective of this work to carefully derive the current model from first principles for light propagation in diffusive media. However, we do want to offer the following insights. First, in probabilistic representation, the solution to the diffusion equation, the photon density u(x) (after being normalized by the source strength), is the probability of finding a photon at the location x. Therefore uu is the probability of finding two photons simultaneously at x. This probability, multiplied by the absorption rate μ, is the total two-photon absorption at x. The reason we use |u|u, not uu, is to ensure that the solution to the resulting nonlinear diffusion equation, that is u, is nonnegative, as is required by physics; see more discussions in Ref. 53. Second, one can indeed derive our nonlinear diffusion model from the nonlinear radiative transfer model used in Ref. 32 if the underlying medium is highly scattering media, following classical results in kinetic theory.54 Therefore, Eq. (1) is not a completely artificial model.

The nonlinear term μ|u|u makes the model Eq. (1) harder to solve computationally than the classical linear diffusion model. It is important to emphasize that this nonlinear diffusion model is indeed a well-posed mathematical model that admits a unique solution for a given illumination source g under classical assumptions on regularities of the coefficients and the domain. Classical numerical discretization schemes, such as finite element and finite difference methods, can be used to discretize the equation. Iterative schemes such as the Newton’s method can be used to solve the resulting nonlinear algebraic system; see more detailed discussions in Ref. 53.

The initial pressure field generated by the photoacoustic effect in TP-PAT is the product of the Grüneisen coefficient of the medium, Γ, and the total energy absorbed locally by the medium, σu+μ|u|u. Note that here the total absorbed energy consists of two components, the contribution from single-photon absorption, σu, and the contribution from two-photon absorption, μ|u|u. Therefore, we write the initial pressure field as53

Eq. (2)

H(x)=Γ(x)[σ(x)u(x)+μ(x)|u|u(x)],xΩ,
where the Grüneisen coefficient is nondimensionalized, and it describes the efficiency of the photoacoustic effect of the underlying medium.

The change of pressure field generates ultrasound waves that propagate following the standard acoustic wave equation, the same model equation for ultrasound propagation in the regular PAT11

Eq. (3)

1c2(x)2pt2Δp=0,in  (0,+)×Rdp(t,x)=χΩH(x),pt(t,x)=0,in  {t=0}×Rd,
where p(t,x) is the pressure field, and c(x) is the speed of the ultrasound waves. In most biological applications of PAT, the ultrasound speed c is assumed known and is often taken as the speed of ultrasound in water since most biological tissues behave like water to ultrasound waves. The function χΩ is the characteristic function of the domain Ω. It should be understood as the extension operator that extends the initial pressure field inside the medium Ω to the whole space Rd, that is,
χΩH(x)={H(x),xΩ0,xRdΩ.

The acoustic datum measured in TP-PAT is the ultrasound signal on the surface of the medium, p|(0,T]×Ω, for time T sufficiently long, and very often, we need to measure data that are generated from multiple illumination sources. From the measured data, we are interested in reconstructing the physical coefficients (Γ,γ,σ,μ) of the underlying medium. Note that among all the coefficients, the two-photon absorption coefficient μ is the only new coefficient that appears in TP-PAT. The coefficients (Γ,γ,σ) are also quantities to be reconstructed in the regular PAT.3,5,13,17,30,52,5561 Mathematical analysis in Ref. 53 shows that one can not simultaneously reconstruct all four coefficients (Γ,γ,σ,μ) even from data collected from multiple illuminations, if all illumination sources have the same optical wavelength. We will therefore only focus on the two absorption coefficients (σ,μ) in the rest of the paper. The reconstruction of all four coefficients using multispectral data following the ideas in Refs. 3, 21, 55, 6263.64.65.66 will be the subject of a future work.

3.

A Two-Step Reconstruction Method

We now present a numerical method for the reconstruction of the absorption coefficients (σ,μ). We follow the standard two-step procedure in quantitative PAT image reconstructions. In the first step, the qualitative step, we reconstruct the initial pressure field, H in Eq. (2), from measured acoustic data using the wave equation (3). In the second step, the quantitative step, we reconstruct the absorption coefficients from the initial pressure field H using the nonlinear diffusion equation (1). We emphasize that the reason for choosing this two-step reconstruction strategy, instead of the more recent one-step algorithms such as those in Refs. 6768.69, is that the two-step method allows us to avoid solving the nonlinear diffusion equation (1) in the quantitative reconstruction step in this specific setup; see more discussions in Sec. 3.2.

3.1.

Qualitative Step: Reconstructing Initial Pressure Fields

In the qualitative step of TP-PAT, we aim at reconstructing the initial pressure field H from measured datum p|(0,T]×Ω. This step is the same as that in the regular PAT, which has been extensively studied in the past. Many efficient algorithms have been proposed; see for instance Refs. 4, 27, 31, 7071.72.73.74.75.76.77.78.79.80.81.82 for an incomplete list of works in this direction. We implement here a simple least-square-based algorithm for the reconstruction.

To simplify the presentation, let us denote by A the linear operator that takes the initial pressure field H(x) to the acoustic field on the boundary Ω, i.e.,

Eq. (4)

p(t,x)|(0,T]×Ω=AH.
Our objective is to invert the operator A to find H for a given measurement p*(t,x)|(0,T]×Ω. We solve this problem in the least-square sense, that is, we search for H as the minimizer of the misfit functional

Eq. (5)

ϕ(H)0TΩ(pp*)2dxdtAHp*L2((0,T]×Ω)2.
Standard least-square theory then implies that the minimizer H solves the normal equation

Eq. (6)

AAH=Ap*,
where A denotes the L2-adjoint of the operator A. We therefore need to invert the self-adjoint operator AA to find H.

We solve the normal equation (6) using the conjugate gradient method.83 In a nutshell, the method seeks a solution of the normal equation by iteratively choosing “conjugate” or “AA-orthogonal” directions, and minimizing the magnitude of the residual, AHp*L2((0,T]×Ω)2, in each of these conjugate directions.

More precisely, let Hk be the value of H at iteration k and let {hi}i=1k be the set of “AA-orthogonal” directions constructed in the first k iterations. The directions {hi}i=1k satisfy the AA orthogonality relation,   2ik

hi,AAhjL2(Ω)=0,  1ji1.

We now search for an update of Hk in the direction hk such that the residual is minimized after the update. That is, we minimize the residual Ψ(α) over α with

Ψ(α)=A(Hk+αhk)p*L2((0,T]×Ω)2=A(Hk+αhkH),A(Hk+αhkH)L2((0,T]×Ω)=Hk+αhkH,AA(Hk+αhkH)L2(Ω),
where it was recalled p*=AH. The optimality condition immediately gives that the step length at iteration k is
αk=hk,AA(HkH)L2(Ω)hk,AAhkL2(Ω)=hk,skL2(Ω)hk,AAhkL2(Ω),with  sk=A[A(HkH)].
Note that if we define rk=A(HkH)=AHkp* as the residual of the original problem at step k, then sk=Ark is simply the so-called normal residual corresponding to rk. The updated value Hk+1 is then obtained as
Hk+1=Hk+αkhk,
while the normal equation residual sk is updated as
sk+1=skαkAAhk.

The conjugate gradient method updates the search direction following:

hk+1=sk+1+sk+1L2(Ω)2skL2(Ω)2hk.
We summarize the conjugate gradient method in Algorithm 1 following the routine in Ref. 83, with an accuracy tolerance parameter ϵ>0 and the maximal number of iteration K.

Algorithm 1

CG algorithm for qualitative reconstruction.

1: Set parameters ϵ and K; set k=0
2: Set initial guess H=0
3: Evaluate the residual r=p*AH and the normal residual s=Ar
4: Set initial search directions h=s
5: Evaluate the size of normal residual γ=sL2(Ω)2
6: whilekK and γ/Ap*L2(Ω)2>ϵdo
7: g=Ah
8: α=γ/gL2((0,T]×Ω)2
9: H=H+αh
10: r=rαg, s=Ar
11: β=sL2(Ω)2/γ
12: γ=sL2(Ω)2
13: h=s+βh
14: k=k+1
15: end while

It is often the case that a regularization term is added to the misfit functional ϕ(H) in Eq. (5). In our implementation, we did not include a regularization term in ϕ(H). When it is needed, the two algorithmic parameters ϵ and K can both serve as mechanisms to regularize the reconstruction. We did not pursue in this direction in this study, but we understand that tuning regularization can refine some of the reconstruction results that we show in the next section.

Let us mention that even though the operator A and its adjoint operator A are called in each iteration of the algorithm, these operators are never explicitly formed in the numerical implementation. We only need to know the actions of these operators on given vectors. For instance, to evaluate Ah for a given h(x), we solve the acoustic wave equation (3) with initial condition p(0,x)=h(x) and record the solution on the boundary of Ω: p(t,x)|(0,T]×Ω. To evaluate Ar for a given r(t,x), we first solve the following adjoint wave equation:

Eq. (7)

1c2(x)2vt2Δv=0,in  (0,T)×Ωv(t,x)=0,vν(t,x)=r,in  (0,T)×Ωv(t,x)=0,vt(t,x)=0,in  {t=T}×Ω.
We then take Ar=vt(0,x). The derivation of Eq. (7) is straightforward and has been documented previously,67,84 so we omit the details here.

3.2.

Quantitative Step: Reconstructing Absorption Coefficients

The second step, the quantitative step, is to reconstruct the optical coefficients from the initial pressure field H recovered in the first step. In recent years, this step was the subject of many computational studies in the case of the regular PAT; see, for instance, Refs. 3, 14, 16, 17, 23, 30, 52, 5556.57.58.59.60.61, 8485.86 for a partial list of references.

Our main objective here is to develop an algorithm to reconstruct the absorption coefficients (σ,μ) in TP-PAT to show that we can separate two-photon absorption from single-photon absorption. We assume that both the Grüneisen coefficient Γ and the diffusion coefficient γ are known already, for instance from a regular PAT reconstruction.

We assume that we have reconstructed initial pressure fields generated from J2 illuminations sources. We denote by {Hj=Γ[σuj+μ|uj|uj]}j=1J those initial pressure fields, where uj denotes the solution of the nonlinear diffusion equation with illumination sources gj (1jJ).

We first reconstruct from the initial pressure fields {Hj}j=1J, using the fact that Γ and γ are known, the quantities

Eq. (8)

σuj+μ|uj|uj=HjΓ,1jJ.
This allows us to replace the term σuj+μ|uj|uj in the nonlinear diffusion equation (1) for source j to obtain the following linear diffusion equation for uj (1jJ)

Eq. (9)

·(γuj)=HjΓ  in  Ω,uj+κujν=gj  on  Ω.
We can solve this linear elliptic equation to reconstruct uj, again since Γ and γ are known. Therefore, we can reconstruct the quantities

Eq. (10)

σ+μ|uj|=HjΓuj,1jJ.
Therefore, at each point xΩ, we have the following system to determine σ and μ
(1|u1|1|uJ|)(σμ)=(H1Γu1HJΓuJ).
We then reconstruct (σ,μ) by solving this small linear system, in least-square sense, at each point xΩ, to get

Eq. (11)

(σμ)=[(11|u1||uJ|)(1|u1|1|uJ|)]1(11|u1||uJ|)(H1Γu1HJΓuJ)=(Jj=1J|uj|j=1J|uj|j=1J|uj|2)1(j=1JHjΓujj=1JHj|uj|Γuj).
We have assumed here that the small 2×2 matrix
(Jj=1J|uj|j=1J|uj|j=1J|uj|2)
in Eq. (11) is invertible at each point xΩ. Theoretical analysis in Ref. 53 shows that one can indeed invert this matrix if the illuminations are selected carefully, that is, the illuminations are sufficiently different from each other. In our numerical experiments, we observe that this matrix is invertible for almost all illuminations that we have tried.

We now summarize the quantitative reconstruction step in Algorithm 2.

Algorithm 2

Noniterative algorithm for quantitative reconstruction.

1: forj1,Jdo
2: Reconstruct the quantities σuj+μ|uj|uj following Eq. (8)
3: end for
4: forj1,Jdo
5: Solve the diffusion equation (9) to reconstruct uj
6: end for
7: forj1,Jdo
8: Reconstruct the quantities σ+μ|uj| following Eq. (10)
9: end for
10: for each point xΩdo
11: Evaluate ω=j=1J|uj(x)|, θ=j=1J|uj(x)|2, ξ=j=1JHjΓuj and ζ=j=1JHj|uj|Γuj
12: Evaluate (σ,μ) using the formula: [σ(x)μ(x)]=(Jωωθ)1(ξζ)
13: end for

Let us emphasize two important features of the quantitative reconstruction algorithm. First, even though the reconstruction of the coefficients (σ,μ) from initial pressure field H is a nonlinear inverse problem, our reconstruction method is noniterative. Therefore, there is no convergence issues at all. The method is guaranteed to give the correct reconstruction result. Second, the reconstruction algorithm is computationally cheap. The major computational cost of the reconstruction algorithm is the solution of the J linear diffusion equations in Eq. (9). The cost in dealing with the algebraic calculations in the rest of the algorithm, in Eqs. (8), (10), and (11), is almost negligible.

4.

Numerical Implementations

We now provide some details on the numerical implementation of the two-step algorithm for quantitative TP-PAT image reconstructions. We limit ourselves to two-dimensional simulations. Nothing changes in three-dimensional case besides the increasing of the computational cost. We use the notational convention x=(x,y) for the spatial variable.

Since the units of the physical quantities in the nonlinear diffusion model (1) and the acoustic wave Eq. (3) are very different, we first normalize the problems by taking the following convention. We take the spatial domain Ω to be the unit square Ω=[0,1]2 and set the ultrasound speed c=1. We set the time interval where the measurements are taken as (0,T] with T=3. This convention means that if we take the size of Ω in units of cm2, the ultrasound speed at 1.5×105  cm/s, then T=3 equals 20  μs. We observed in our numerical simulations (see discussions in the next section) that these choices of Ω, c, and T for the wave equation are sufficient to capture all of the physical wave signals generated by the initial conditions H(x) we have tested, at least up to the numerical discretization errors. For the diffusion problem, we set U0=1011 as the characteristic photon density involved in the system, and normalize the solution u and the boundary illumination against U0.

The wave equation (3) is posed on R2, not inside Ω. We therefore have to make a truncation to have a finite domain for the wave simulation. We do this by using the technique of perfectly matched layers (PML).87,88 We surround our physical domain Ω with a PML region of thickness 0.2 to have the computational domain Ω˜=[0.2,1.2]2. We use a split-field PML scheme (see, e.g., Refs. 87 and 88). This scheme reduces to the undamped wave equation in the physical domain Ω, coupled with a damped wave split-field scheme in the PML region Ω˜Ω. Ultimately, we end up solving the system of equations, assuming again that c=1

{pxt+τxpx=vxx,pyt+τypy=vyy,vxt+τxvx=px,vyt+τyvy=py,
where p=px+py, τx(x) and τy(x) are absorptive terms supported only in the PML region Ω˜Ω. In our simulations, we use τx(x)=τx(x,y)=χx>1(x)α(x1)2+χx<0(x)αx2 with α a given constant. Similarly, τy(x,y)=τx(y,x). Initial and boundary conditions can be transformed into this first-order formulation in a straightforward way.

We discretize these equations using standard second-order finite differences in space, and first-order finite differences in time on uniform spatial-temporal grids. The spatial grid covering Ω˜ consists of 141×141 spatial points

{(xi,yj):xi=i/100,yj=j/100,20i,j120},
and the temporal grid covering [0,T]=[0,3] consists of 3001 grid points
{tk:tk=k/1000,0k3000}.
The velocity fields vx and vy are solved for at staggered half-time steps, i.e., at times {tk+1/2}, using values of the split pressure fields px and py at the usual time steps {tk}. The pressure field p=px+py is then updated using the velocity fields at these staggered times. Hence, the scheme is known as a leapfrog scheme and reduces to standard second-order finite difference time stepping in the physical domain Ω where τx=τy=0. With our spatial step size h=1/100 and our temporal step size Δt=1/1000, we clearly satisfy the CFL stability condition
(cΔt)2h2=104106=102<1,
which is necessary for the explicit finite difference scheme we implemented to be stable.

To solve the adjoint wave equation (7), we first perform the change of variable t=Tt to transform the equation into an initial value (instead of final value) problem. We then apply the same type of spatial-temporal discretizations to the new equation.

The nonlinear diffusion equation (1) and the linear diffusion equations involved in the reconstruction process, mainly in Eq. (9), are all discretized using a standard first-order finite element method with about 12,000 elements on a triangular mesh of Ω. The nonlinear algebra system resulting from the discretization of Eq. (1) is solved with the Newton’s method. In the forward simulation, the initial pressure field H that is needed in the acoustic wave equation (3) is linearly interpolated from the quadrature points of the triangular elements. In the reconstruction process, the initial pressure field H reconstructed in the qualitative step, i.e., the first step, is interpolated back to the quadrature points of the triangular elements as datum for the quantitative reconstruction step. These interpolation processes induce additional noise in the reconstruction process besides the artificial white noise we add to the acoustic data that we discuss below.

To generate synthetic data, we solve the nonlinear diffusion equation (1) with the true physical coefficients to generate H and then solve the acoustic wave equation (3) to produce p(t,x)|(0,T]×Ω. To add noise to the synthetic data, we use the following strategy. At each point xΩ where the ultrasound signal is measured, we generate an independent Gaussian “white-noise” process wx(t), t(0,T], that satisfies

E[wx(t)]=0,and,E[wx(t)wx(s)]=δ(ts),
where E denotes the operation of taking expectation. We then scale the white noise wx(t) according to the power of the signal p(t,x), to generate noisy data p˜(t,x) with a specified noise-to-signal ratio (NSR) η
p˜(t,x)=p(t,x)+η[0Tp2(t,x)dt0Twx2(t)dt]1/2wx(t),t(0,T].

In our numerical simulations in the upcoming section, we set the tolerance level ϵ=106 in Algorithm 1 and run this algorithm for a maximum number of K=1000 iterations. The quantity AHkp*L2((0,T]×Ω)2 is usually guaranteed to decrease monotonically [see, e.g., the discussion in (Ref. 83, Sec. 7.4)]. However, it is clear from our description above that our implementation of the operators A and A constitute only approximate adjoints of one another due to errors associated with the finite difference approximations and PML region in the wave solver. Therefore, we also force Algorithm 1 to exit if nonmonotonic behavior of AHkp*L2((0,T]×Ω)2 is encountered.

5.

Numerical Simulations

We now present some numerical simulations to demonstrate the performance of our quantitative reconstruction strategy. Our main objective is to show that the photoacoustic data measured in TP-PAT allow quantitative separation of the single-photon and the two-photon absorption coefficients. In all the simulation results below, we set the Grüneisen coefficient Γ=1 and the diffusion coefficient γ(x)=0.02+0.01sin(2πy). Moreover, we use data collected from four different illumination sources. The first source function takes a constant value the top and right sides of the boundary and is zero everywhere else. That is,

g1(x)={1,x(0,1)×{1}{1}×(0,1)0,otherwise.
The second to fourth sources, g2, g3, and g4, are obtained by rotating g1 by π/2, π, and 3π/2, respectively, along the boundary. Note again that g1 is normalized against U0 already.

To measure the quality of the reconstructions, we use relative L2 errors. Let f be a quantity to be reconstructed, ft its true value, and fr the reconstructed value. Then, the relative L2 error of the reconstruction, denoted by EL2(f), is the ratio between the size of the error in the reconstruction and the size of the true quantity. That is,

EL2(f)=ftfrL2(Ω)ftL2(Ω).

5.1.

Experiment I

In the first group of numerical simulations, we attempt to reconstruct the absorption coefficients (σ,μ) shown in Fig. 1. We first perform reconstructions, using Algorithm 1, on the initial pressure field H generated by the four illumination sources {gi}i=14 from acoustic data of different noise levels. The quality of the reconstructions, in terms of the relative L2 errors, is summarized in Table 1, third column. We observed, as has been confirmed by many works in the PAT community, the qualitative reconstruction is of high quality. We show in Figs. 2 and 3 some reconstructions with illuminations g1 and g2, respectively. Shown are the true initial pressure field H, the reconstructed H using clean data (NSR η=0.0) and noisy data with NSR η=0.1. The relative L2 errors of the reconstruction in Fig. 2 are 0.04% for the case of η=0.0 and 2.7% for the case of η=0.1, whereas those for the reconstructions in Fig. 3 are 0.06% for the case of η=0.0 and 2.6% for the case of η=0.1.

Fig. 1

The true absorption coefficients, (a) σ (cm1) and (b) μ×U0 (cm1), used to generate synthetic data in experiment I.

JBO_23_1_016002_f001.png

Table 1

Quality of reconstructions in experiment I. Shown are relative L2 errors in the reconstructions of various initial pressure fields (third column) and the corresponding absorption coefficients in Fig. 1 (fourth column) from ultrasound data with different noise levels (controlled with the NSR η).

NSRIlluminationEL2(H)[EL2(σ),EL2(μ)]
η=0.00g14.05×104(0.46,3.33)×102
g26.66×104
g36.22×104
g44.96×104
η=0.01g17.13×103(1.71,4.08)×102
g27.30×103
g38.25×103
g48.00×103
η=0.05g11.78×102(3.80,6.33)×102
g21.80×102
g31.77×102
g41.93×102
η=0.10g12.68×102(5.15,8.01)×102
g22.63×102
g32.48×102
g42.55×102

Fig. 2

(a) The initial pressure field H(x) generated from illumination g1 using the true absorption coefficients in Fig. 1, (b) as well as the reconstructions of H using ultrasound data with NSR η=0.0 (clean data) and (c) NSR η=0.1.

JBO_23_1_016002_f002.png

Fig. 3

The same as Fig. 2 except that the illumination source used is g2.

JBO_23_1_016002_f003.png

Let us mention that even though we can see clearly the two-photon absorbing inclusions in the true initial pressure field H and the reconstructed H, the true absorption coefficients in Fig. 1 are very different from the H in Figs. 2 and 3. In other words, knowing H does not provide us enough information about the true absorption coefficients unless we perform the next step, the quantitative step, of the reconstruction.

In Fig. 4, we show the reconstructions of the coefficient pair (σ,μ) in Fig. 1 from the initial pressure fields we obtained in the qualitative step, using Algorithm 2. Shown, from left to right, are reconstructions using noisy data with NSR η=0.00, η=0.05, and η=0.10, respectively. The quality of the reconstructions is very high with relative L2 errors (1.71,4.08)×102, (3.80,6.33)×102, and (5.15,8.01)×102, respectively; see the fourth column of Table 1 for the reconstruction result using data with η=0.01, which we did not show here since it is too similar to the case with η=0.00.

Fig. 4

The absorption coefficients (a) σ and (b) μ×U0 reconstructed using noisy data with NSR η=0.00, η=0.05, and η=0.10 (from left to right).

JBO_23_1_016002_f004.png

The reconstructions in Fig. 4 show that by performing quantitative reconstructions, we can separate the two-photon absorption coefficient from the single-photon absorption coefficient from the initial pressure field. This is clearly important for practical applications of TP-PAT where two-photon absorption is the main quantity of interest.

5.2.

Experiment II

In the second group of numerical simulations, we study the reconstruction of the absorption coefficients (σ,μ) shown in Fig. 5. In Fig. 6, we present the true initial pressure filed H computed with illumination source g1, and the reconstructions of this H with clean ultrasound data (NSR η=0.00) and noisy data (NSR η=0.10). By visual inspection, we can see the presence of both the single-photon absorption and the two-photon absorption inclusions in H. The reconstructions are impressively good, with relative L2 errors 9.30×104 and 2.51×102 for η=0.00 and η=0.10, respectively, even when H is this complicated. We have also performed similar reconstructions for H generated from the other three illumination sources g2, g3, and g4. The relative errors in the reconstructions are summarized in Table 2, third column. Despite its slight degeneration as noise level increases, the quality of the reconstructions of H remains high at moderate noise levels.

Fig. 5

The true absorption coefficients, (a) σ (cm1) and (b) μ×U0 (cm1), used to generate synthetic data in experiment II.

JBO_23_1_016002_f005.png

Fig. 6

(a) The initial pressure field H(x) generated from illumination g1 using the true absorption coefficients in Fig. 5, (b) as well as the reconstructions of H using ultrasound data with NSR η=0.00 (clean data) and (c) NSR η=0.10.

JBO_23_1_016002_f006.png

Table 2

Quality of reconstructions in experiment II. Shown are relative L2 errors in the reconstructions of various initial pressure fields (third column) and the corresponding absorption coefficients in Fig. 5 (fourth column) from ultrasound data with different noise levels.

NSRIlluminationEL2(H)[EL2(σ),EL2(μ)]
η=0.00g19.30×104(2.19,5.26)×102
g21.16×103
g31.17×103
g49.38×104
η=0.01g18.09×103(2.77,5.78)×102
g28.25×103
g37.79×103
g47.68×103
η=0.05g11.72×102(4.31,7.73)×102
g21.82×102
g31.82×102
g41.78×102
η=0.10g12.51×102(5.60,9.34)×102
g22.48×102
g32.40×102
g42.79×102

To separate μ from σ in the initial pressure fields, we perform quantitative reconstructions using Algorithm 2. In Figs. 7(a)7(c), we show the reconstructions from data with NSR η=0.00, η=0.05, and η=0.10, respectively. The relative L2 errors in the reconstructions are (2.19,5.26)×102, (4.31,7.73)×102, and (5.60,9.34)×102, respectively. Once again, we see good separation of the two different absorption coefficients, which were mixed together in the initial pressure fields H in Fig. 6. The quantitative reconstruction results are summarized in Table 2, fourth column.

Fig. 7

Reconstruction of the absorption coefficient pair (a) σ and (b) μ×U0 in Fig. 5 using data at different noise levels (NSR η=0.00, η=0.05, and η=0.10 from left to right).

JBO_23_1_016002_f007.png

6.

Concluding Remarks

We studied in this paper quantitative image reconstructions in TP-PAT, aiming at reconstructing the single-photon absorption and the two-photon absorption coefficients of biological tissues from measured ultrasound signals generated by the photoacoustic effect of light absorption. We introduced a nonlinear diffusion equation as the model for light propagation in TP-PAT and presented a two-step image reconstruction strategy, including a noniterative quantitative reconstruction step, based on this model. We showed, with computational simulations, that while single-photon absorption and two-photon absorption are mixed in the images of the initial pressure fields, they can be stably separated from each other through the quantitative reconstruction step, using Algorithm 2, even when the ultrasound data contain relatively high level of random noises. Our numerical simulations confirm the results of mathematical analysis of the problem in a previous publication.53

Compared to the case in the regular PAT, quantitative image reconstruction in TP-PAT is far less investigated, theoretically or computationally, to date. Our numerical simulations show great promise in the quantitative imaging of the two-photon absorption. However, there are still a lot of issues that need to be addressed. For instance, it would be very interesting to test the two-step reconstruction method we proposed against experimentally measured data to see what types of resolution and contrast we can get for the two-photon absorption coefficient. It would also be important to develop efficient algorithms to reconstruct the diffusion coefficient γ in addition to the absorption coefficients. Last but not the least, reconstructing the Grüneisen coefficient with multispectral data, following for instance the ideas in Refs. 3, 21, 55, 6263.64.65.66, could also be extremely useful as well.

Disclosures

The authors have no relevant financial interests in this article and no potential conflicts of interest to disclose.

Acknowledgments

We would like to thank Prof. John C. Schotland (University of Michigan) for useful discussions on the two-photon absorption photoacoustic imaging. This work was partially supported by the National Science Foundation through grants DMS-1321018 and DMS-1620473.

References

1. 

P. Beard, “Biomedical photoacoustic imaging,” Interface Focus, 1 602 –631 (2011). http://dx.doi.org/10.1098/rsfs.2011.0028 Google Scholar

2. 

P. Burgholzer, H. Grun and A. Sonnleitner, “Photoacoustic tomography: sounding out fluorescent proteins,” Nat. Photonics, 3 378 –379 (2009). http://dx.doi.org/10.1038/nphoton.2009.109 NPAHBY 1749-4885 Google Scholar

3. 

B. T. Cox, J. G. Laufer and P. C. Beard, “The challenges for quantitative photoacoustic imaging,” Proc. SPIE, 7177 717713 (2009). http://dx.doi.org/10.1117/12.806788 PSISDG 0277-786X Google Scholar

4. 

P. Kuchment and L. Kunyansky, “Mathematics of thermoacoustic tomography,” Eur. J. Appl. Math., 19 191 –224 (2008). http://dx.doi.org/10.1017/S0956792508007353 0956-7925 Google Scholar

5. 

C. Li and L. Wang, “Photoacoustic tomography and sensing in biomedicine,” Phys. Med. Biol., 54 R59 –R97 (2009). http://dx.doi.org/10.1088/0031-9155/54/19/R01 PHMBA7 0031-9155 Google Scholar

6. 

A. A. Oraevsky, S. L. Jacques and F. K. Tittel, “Measurement of tissue optical properties by time-resolved detection of laser-induced transient stress,” Appl. Opt., 36 402 –415 (1997). http://dx.doi.org/10.1364/AO.36.000402 APOPAI 0003-6935 Google Scholar

7. 

A. A. Oraevsky, A. A. Karabutov, “Optoacoustic tomography,” Handbook of Optical Biomedical Diagonstics, SPIE Press, Bellingham, Washington (2002). Google Scholar

8. 

S. K. Patch and O. Scherzer, “Photo- and thermo-acoustic imaging,” Inverse Prob., 23 S1 –S10 (2007). http://dx.doi.org/10.1088/0266-5611/23/6/S01 INPEEY 0266-5611 Google Scholar

9. 

O. Scherzer, Handbook of Mathematical Methods in Imaging, Springer-Verlag, Berlin (2010). Google Scholar

10. 

E. M. Strohm, M. J. Moore and M. C. Kolios, “Single cell photoacoustic microscopy: a review,” IEEE J. Sel. Top. Quantum Electron., 23 137 –151 (2016). http://dx.doi.org/10.1109/JSTQE.2015.2497323 Google Scholar

11. 

Photoacoustic Imaging and Spectroscopy, Taylor & Francis, Oxford (2009). Google Scholar

12. 

L. V. Wang, “Photoacoustic tomography,” Scholarpedia, 9 10278 (2014). http://dx.doi.org/10.4249/scholarpedia.10278 1941-6016 Google Scholar

13. 

J. Xia, J. Yao and L.V. Wang, “Photoacoustic tomography: principles and advances,” Electromagn. Waves, 147 1 –22 (2014). http://dx.doi.org/10.2528/PIER14032303 PELREX 1043-626X Google Scholar

14. 

S. Arridge et al., “Accelerated high-resolution photoacoustic tomography via compressed sensing,” Phys. Med. Biol., 61 8908 –8940 (2016). http://dx.doi.org/10.1088/1361-6560/61/24/8908 PHMBA7 0031-9155 Google Scholar

15. 

J. R. Cook, W. Frey and S. Emelianov, “Quantitative photoacoustic imaging of nanoparticles in cells and tissues,” ACS Nano, 7 1272 –1280 (2013). http://dx.doi.org/10.1021/nn304739s ANCAC3 1936-0851 Google Scholar

16. 

B. Cox et al., “Quantitative spectroscopic photoacoustic imaging: a review,” J. Biomed. Opt., 17 061202 (2012). http://dx.doi.org/10.1117/1.JBO.17.6.061202 JBOPFO 1083-3668 Google Scholar

17. 

B. T. Cox et al., “Two-dimensional quantitative photoacoustic image reconstruction of absorption distributions in scattering media by use of a simple iterative method,” Appl. Opt., 45 1866 –1875 (2006). http://dx.doi.org/10.1364/AO.45.001866 APOPAI 0003-6935 Google Scholar

18. 

X. Gao et al., “Quantitative imaging of microvasculature in deep tissue with a spectrum-based photo-acoustic microscopy,” Opt. Lett., 40 970 –973 (2015). http://dx.doi.org/10.1364/OL.40.000970 OPLEDP 0146-9592 Google Scholar

19. 

W. Poon et al., “Determination of biodistribution of ultrasmall, near-infrared emitting gold nanoparticles by photoacoustic and fluorescence imaging,” J. Biomed. Opt., 20 066007 (2015). http://dx.doi.org/10.1117/1.JBO.20.6.066007 JBOPFO 1083-3668 Google Scholar

20. 

D. Razansky and V. Ntziachristos, “Hybrid photoacoustic fluorescence molecular tomography using finite-element-based inversion,” Med. Phys., 34 4293 –4301 (2007). http://dx.doi.org/10.1118/1.2786866 MPHYA6 0094-2405 Google Scholar

21. 

D. Razansky, C. Vinegoni and V. Ntziachristos, “Multispectral photoacoustic imaging of fluorochromes in small animals,” Opt. Lett., 32 2891 –2893 (2007). http://dx.doi.org/10.1364/OL.32.002891 OPLEDP 0146-9592 Google Scholar

22. 

J. Ripoll and V. Ntziachristos, “Quantitative point source photoacoustic inversion formulas for scattering and absorbing media,” Phys. Rev. E, 71 031912 (2005). http://dx.doi.org/10.1103/PhysRevE.71.031912 Google Scholar

23. 

P. Shao, T. Harrison and R. J. Zemp, “Iterative algorithm for multiple illumination photoacoustic tomography (MIPAT) using ultrasound channel data,” Biomed. Opt. Express, 3 3240 –3249 (2012). http://dx.doi.org/10.1364/BOE.3.003240 BOEICL 2156-7085 Google Scholar

24. 

T. Tarvainen et al., “Reconstructing absorption and scattering distributions in quantitative photoacoustic tomography,” Inverse Prob., 28 084009 (2012). http://dx.doi.org/10.1088/0266-5611/28/8/084009 Google Scholar

25. 

C. Tian et al., “Imaging and sensing based on dual-pulse nonlinear photoacoustic contrast: a preliminary study on fatty liver,” Opt. Lett., 40 2253 –2256 (2015). http://dx.doi.org/10.1364/OL.40.002253 OPLEDP 0146-9592 Google Scholar

26. 

P. K. Upputuri, M. Krisnan and M. Pramanik, “Microsphere enabled subdiffraction-limited optical-resolution photoacoustic microscopy: a simulation study,” J. Biomed. Opt., 23 045001 (2017). http://dx.doi.org/10.1117/1.JBO.22.4.045001 JBOPFO 1083-3668 Google Scholar

27. 

K. Wang and M. A. Anastasio, “A simple Fourier transform-based reconstruction formula for photoacoustic computed tomography with a circular or spherical measurement geometry,” Phys. Med. Biol., 57 N493 (2012). http://dx.doi.org/10.1088/0031-9155/57/23/N493 PHMBA7 0031-9155 Google Scholar

28. 

K. Wang et al., “Investigation of iterative image reconstruction in three-dimensional optoacoustic tomography,” Phys. Med. Biol., 57 5399 –5423 (2012). http://dx.doi.org/10.1088/0031-9155/57/17/5399 PHMBA7 0031-9155 Google Scholar

29. 

L. Yao, Y. Sun and H. Jiang, “Transport-based quantitative photoacoustic tomography: simulations and experiments,” Phys. Med. Biol., 55 1917 –1934 (2010). http://dx.doi.org/10.1088/0031-9155/55/7/009 PHMBA7 0031-9155 Google Scholar

30. 

R. J. Zemp, “Quantitative photoacoustic tomography with multiple optical sources,” Appl. Opt., 49 3566 –3572 (2010). http://dx.doi.org/10.1364/AO.49.003566 APOPAI 0003-6935 Google Scholar

31. 

B. Zhu and E. M. Sevick-Muraca, “Reconstruction of sectional images in frequency-domain based photoacoustic imaging,” Opt. Express, 19 23286 –23297 (2011). http://dx.doi.org/10.1364/OE.19.023286 OPEXFF 1094-4087 Google Scholar

32. 

P. Chantharasupawong, R. Philip and J. Thomas, “Simultaneous optical and photoacoustic measurement of nonlinear absorption,” Appl. Phys. Lett., 102 041116 (2013). http://dx.doi.org/10.1063/1.4789870 APPLAB 0003-6951 Google Scholar

33. 

Y.-H. Lai et al., “Nonlinear photoacoustic microscopy via a loss modulation technique: from detection to imaging,” Opt. Express, 23 525 –536 (2014). http://dx.doi.org/10.1364/OE.22.000525 OPEXFF 1094-4087 Google Scholar

34. 

G. Langer et al., “Two-photon absorption-induced photoacoustic imaging of rhodamine B dyed polyethylene spheres using a femtosecond laser,” Opt. Express, 21 22410 –22422 (2013). http://dx.doi.org/10.1364/OE.21.022410 OPEXFF 1094-4087 Google Scholar

35. 

G. Langer et al., “Two-photon absorption-induced photoacoustic and luminescence imaging employing a femtosecond laser,” Proc. SPIE, 8943 89431V (2014). http://dx.doi.org/10.1117/12.2037313 PSISDG 0277-786X Google Scholar

36. 

G. J. Tserevelakis et al., “Hybrid multiphoton and optoacoustic microscope,” Opt. Lett., 39 1819 –1822 (2014). http://dx.doi.org/10.1364/OL.39.001819 OPLEDP 0146-9592 Google Scholar

37. 

B. E. Urban et al., “Investigating femtosecond-laser-induced two-photon photoacoustic generation,” J. Biomed. Opt., 19 085001 (2014). http://dx.doi.org/10.1117/1.JBO.19.8.085001 JBOPFO 1083-3668 Google Scholar

38. 

P. W. Winter et al., “Two-photon instant structured illumination microscopy improves the depth penetration of super-resolution imaging in thick scattering samples,” Optica, 1 181 –191 (2014). http://dx.doi.org/10.1364/OPTICA.1.000181 Google Scholar

39. 

Y. Yamaoka, M. Nambu and T. Takamatsu, “Frequency-selective multiphoton-excitation-induced photoacoustic microscopy (MEPAM) to visualize the cross sections of dense objects,” Proc. SPIE, 7524 75642O (2010). http://dx.doi.org/10.1117/12.841258 PSISDG 0277-786X Google Scholar

40. 

Y. Yamaoka, M. Nambu and T. Takamatsu, “Fine depth resolution of two-photon absorption-induced photoacoustic microscopy using low-frequency bandpass filtering,” Opt. Express, 19 13365 –13377 (2011). http://dx.doi.org/10.1364/OE.19.013365 OPEXFF 1094-4087 Google Scholar

41. 

Y. Yamaoka and T. Takamatsu, “Enhancement of multiphoton excitation-induced photoacoustic signals by using gold nanoparticles surrounded by fluorescent dyes,” Proc. SPIE, 7177 71772A (2009). http://dx.doi.org/10.1117/12.808367 PSISDG 0277-786X Google Scholar

42. 

C. S. Yelleswarapu and S. R. Kothapalli, “Nonlinear photoacoustics for measuring the nonlinear optical absorption coefficient,” Opt. Express, 18 9020 –9025 (2010). http://dx.doi.org/10.1364/OE.18.009020 OPEXFF 1094-4087 Google Scholar

43. 

H. Mahr, “Two-photon absorption spectroscopy,” Quantum Electronics, Academic Press, Cambridge, Massachusetts (2012). Google Scholar

44. 

M. Rumi and J. W. Perry, “Two-photon absorption: an overview of measurements and principles,” Adv. Opt. Photonics, 2 451 –518 (2010). http://dx.doi.org/10.1364/AOP.2.000451 AOPAC7 1943-8206 Google Scholar

45. 

A. C. Tam and C. K. Patel, “Two-photon absorption spectra and cross-section measurements in liquids,” Nature, 280 304 –306 (1979). http://dx.doi.org/10.1038/280304a0 Google Scholar

46. 

W. Denk, J. Strickler and W. Webb, “Two photon laser scanning fluorescence microscopy,” Science, 248 73 –76 (1990). http://dx.doi.org/10.1126/science.2321027 SCIEAS 0036-8075 Google Scholar

47. 

P. T. C. So, “Two-photon fluorescence light microscopy,” Encyclopedia of Life Sciences, Nature Publishing Group, London (2002). Google Scholar

48. 

J. Ying, F. Liu and R. R. Alfano, “Spatial distribution of two-photon-excited fluorescence in scattering media,” Appl. Opt., 38 224 –229 (1999). http://dx.doi.org/10.1364/AO.38.000224 APOPAI 0003-6935 Google Scholar

49. 

W. R. Zipfel, R. M. Williams and W. Webb, “Nonlinear magic: multiphoton microscopy in the biosciences,” Nat. Biotechnol., 21 1369 –1377 (2003). http://dx.doi.org/10.1038/nbt899 NABIF9 1087-0156 Google Scholar

50. 

Y. Bae, J. J. Song and Y. B. Kim, “Photoacoustic study of two-photon absorption in hexagonal ZnS,” J. Appl. Phys., 53 615 –619 (1982). http://dx.doi.org/10.1063/1.329967 JAPIAU 0021-8979 Google Scholar

51. 

Y. Yamaoka et al., “Photoacoustic microscopy using ultrashort pulses with two different pulse durations,” Opt. Express, 23 17063 –17072 (2014). http://dx.doi.org/10.1364/OE.22.017063 OPEXFF 1094-4087 Google Scholar

52. 

G. Bal and K. Ren, “Multi-source quantitative PAT in diffusive regime,” Inverse Prob., 27 075003 (2011). http://dx.doi.org/10.1088/0266-5611/27/7/075003 INPEEY 0266-5611 Google Scholar

53. 

K. Ren and R. Zhang, “Nonlinear quantitative photoacoustic tomography with two-photon absorption,” SIAM J. Appl. Math., 78 (2018). SMJMAP 0036-1399 Google Scholar

54. 

R. Dautray and J.-L. Lions, Mathematical Analysis and Numerical Methods for Science and Technology, 6 Springer-Verlag, Berlin (1993). Google Scholar

55. 

J. Laufer et al., “Quantitative determination of chromophore concentrations from 2d photoacoustic images using a nonlinear model-based inversion scheme,” Appl. Opt., 49 1219 –1233 (2010). http://dx.doi.org/10.1364/AO.49.001219 APOPAI 0003-6935 Google Scholar

56. 

A. Pulkkinen et al., “A Bayesian approach to spectral quantitative photoacoustic tomography,” Inverse Prob., 30 065012 (2014). http://dx.doi.org/10.1088/0266-5611/30/6/065012 INPEEY 0266-5611 Google Scholar

57. 

K. Ren, H. Gao and H. Zhao, “A hybrid reconstruction method for quantitative photoacoustic imaging,” SIAM J. Imaging Sci., 6 32 –55 (2013). http://dx.doi.org/10.1137/120866130 Google Scholar

58. 

K. Ren, R. Zhang and Y. Zhong, “Inverse transport problems in quantitative PAT for molecular imaging,” Inverse Prob., 31 125012 (2015). http://dx.doi.org/10.1088/0266-5611/31/12/125012 INPEEY 0266-5611 Google Scholar

59. 

T. Saratoon et al., “A gradient-based method for quantitative photoacoustic tomography using the radiative transfer equation,” Inverse Prob., 29 075006 (2013). http://dx.doi.org/10.1088/0266-5611/29/7/075006 INPEEY 0266-5611 Google Scholar

60. 

J. Zhang and M. A. Anastasio, “Reconstruction of speed-of-sound and electromagnetic absorption distributions in photoacoustic tomography,” Proc. SPIE, 6086 608619 (2006). http://dx.doi.org/10.1117/12.647665 PSISDG 0277-786X Google Scholar

61. 

X. Zhang et al., “Forward-backward splitting method for quantitative photoacoustic tomography,” Inverse Prob., 30 125012 (2014). http://dx.doi.org/10.1088/0266-5611/30/12/125012 INPEEY 0266-5611 Google Scholar

62. 

G. Bal and K. Ren, “On multi-spectral quantitative photoacoustic tomography in diffusive regime,” Inverse Prob., 28 025010 (2012). http://dx.doi.org/10.1088/0266-5611/28/2/025010 INPEEY 0266-5611 Google Scholar

63. 

B. T. Cox, S. R. Arridge and P. C. Beard, “Estimating chromophore distributions from multiwavelength photoacoustic images,” J. Opt. Soc. Am. A, 26 443 –455 (2009). http://dx.doi.org/10.1364/JOSAA.26.000443 JOAOD6 0740-3232 Google Scholar

64. 

D. Razansky et al., “Multispectral opto-acoustic tomography of deep-seated fluorescent proteins in vivo,” Nat. Photonics, 3 412 –417 (2009). http://dx.doi.org/10.1038/nphoton.2009.98 NPAHBY 1749-4885 Google Scholar

65. 

P. Shao, B. Cox and R. J. Zemp, “Estimating optical absorption, scattering, and Grueneisen distributions with multiple-illumination photoacoustic tomography,” Appl. Opt., 50 3145 –3154 (2011). http://dx.doi.org/10.1364/AO.50.003145 APOPAI 0003-6935 Google Scholar

66. 

Z. Yuan and H. Jiang, “Simultaneous recovery of tissue physiological and acoustic properties and the criteria for wavelength selection in multispectral photoacoustic tomography,” Opt. Lett., 34 1714 –1716 (2009). http://dx.doi.org/10.1364/OL.34.001714 OPLEDP 0146-9592 Google Scholar

67. 

T. Ding, K. Ren and S. Vallelian, “A one-step reconstruction algorithm for quantitative photoacoustic imaging,” Inverse Prob., 31 095005 (2015). http://dx.doi.org/10.1088/0266-5611/31/9/095005 INPEEY 0266-5611 Google Scholar

68. 

A. Pulkkinen et al., “Direct estimation of optical parameters from photoacoustic time series in quantitative photoacoustic tomography,” IEEE Trans. Med. Imag., 35 2497 –2508 (2016). http://dx.doi.org/10.1109/TMI.2016.2581211 ITMID4 0278-0062 Google Scholar

69. 

K. Ren and F. Triki, “A global stability estimate for the photoacoustic inverse problem in layered media,” (2017). Google Scholar

70. 

P. Burgholzer et al., “Exact and approximative imaging methods for photoacoustic tomography using an arbitrary detection surface,” Phys. Rev. E, 75 046706 (2007). http://dx.doi.org/10.1103/PhysRevE.75.046706 Google Scholar

71. 

B. T. Cox, S. R. Arridge and P. C. Beard, “Photoacoustic tomography with a limited-aperture planar sensor and a reverberant cavity,” Inverse Prob., 23 S95 –S112 (2007). http://dx.doi.org/10.1088/0266-5611/23/6/S08 INPEEY 0266-5611 Google Scholar

72. 

D. Finch, M. Haltmeier and Rakesh, “Inversion of spherical means and the wave equation in even dimensions,” SIAM J. Appl. Math., 68 392 –412 (2007). http://dx.doi.org/10.1137/070682137 SMJMAP 0036-1399 Google Scholar

73. 

M. Haltmeier, T. Schuster and O. Scherzer, “Filtered backprojection for thermoacoustic computed tomography in spherical geometry,” Math. Methods Appl. Sci., 28 1919 –1937 (2005). http://dx.doi.org/10.1002/(ISSN)1099-1476 Google Scholar

74. 

M. Haltmeier and G. Zangerl, “Spatial resolution in photoacoustic tomography: effects of detector size and detector bandwidth,” Inverse Prob., 26 125002 (2010). http://dx.doi.org/10.1088/0266-5611/26/12/125002 INPEEY 0266-5611 Google Scholar

75. 

Y. Hristova, “Time reversal in thermoacoustic tomography-an error estimate,” Inverse Prob., 25 055008 (2009). http://dx.doi.org/10.1088/0266-5611/25/5/055008 INPEEY 0266-5611 Google Scholar

76. 

B. Huang et al., “Improving limited-view photoacoustic tomography with an acoustic reflector,” J. Biomed. Opt., 18 110505 (2013). http://dx.doi.org/10.1117/1.JBO.18.11.110505 JBOPFO 1083-3668 Google Scholar

77. 

L. Kunyansky, “Thermoacoustic tomography with detectors on an open curve: an efficient reconstruction algorithm,” Inverse Prob., 24 055021 (2008). http://dx.doi.org/10.1088/0266-5611/24/5/055021 INPEEY 0266-5611 Google Scholar

78. 

L. V. Nguyen, “A family of inversion formulas in thermoacoustic tomography,” Inverse Probl. Imaging, 3 649 –675 (2009). http://dx.doi.org/10.3934/ipi Google Scholar

79. 

G. Paltauf et al., “Photoacoustic tomography with integrating area and line detectors,” Photoacoustic Imaging and Spectroscopy, 251 –263 CRC Press, Boca Raton, Florida (2009). Google Scholar

80. 

J. Qian et al., “An efficient Neumann-series based algorithm for thermoacoustic and photoacoustic tomography with variable sound speed,” SIAM J. Imaging Sci., 4 850 –883 (2011). http://dx.doi.org/10.1137/100817280 Google Scholar

81. 

R. W. Schoonover and M. A. Anastasio, “Compensation of shear waves in photoacoustic tomography with layered acoustic media,” J. Opt. Soc. Am. A, 28 2091 –2099 (2011). http://dx.doi.org/10.1364/JOSAA.28.002091 JOAOD6 0740-3232 Google Scholar

82. 

B. E. Treeby, “Acoustic attenuation compensation in photoacoustic tomography using time-variant filtering,” J. Biomed. Opt., 18 036008 (2013). http://dx.doi.org/10.1117/1.JBO.18.3.036008 JBOPFO 1083-3668 Google Scholar

83. 

A. Björck, Numerical Methods for Least Squares Problems, SIAM, Philadelphia, Pennsylvania (1996). Google Scholar

84. 

H. Ammari et al., “Mathematical modelling in photo-acoustic imaging of small absorbers,” SIAM Rev., 52 677 –695 (2010). http://dx.doi.org/10.1137/090748494 SIREAD 0036-1445 Google Scholar

85. 

G. Bal, K. Ren, “Non-uniqueness result for a hybrid inverse problem,” Tomography and Inverse Transport Theory, 559 29 –38 American Mathematical Society, Providence, Rhode Island (2011). Google Scholar

86. 

A. V. Mamonov and K. Ren, “Quantitative photoacoustic imaging in radiative transport regime,” Comm. Math. Sci., 12 201 –234 (2014). http://dx.doi.org/10.4310/CMS.2014.v12.n2.a1 Google Scholar

87. 

J.-P. Bérenger, “A perfectly matched layer for the absorption of electromagnetic waves,” J Comput. Phys., 114 185 –200 (1994). http://dx.doi.org/10.1006/jcph.1994.1159 Google Scholar

88. 

F. Collino and C. Tsogka, “Application of the perfectly matched absorbing layer model to the linear elastodynamic problem in anisotropic heterogeneous media,” Geophysics, 66 294 –307 (2001). http://dx.doi.org/10.1190/1.1444908 GPYSA7 0016-8033 Google Scholar

Biographies for the authors are not available.

© 2018 Society of Photo-Optical Instrumentation Engineers (SPIE) 1083-3668/2018/$25.00 © 2018 SPIE
Patrick Bardsley, Kui Ren, and Rongting Zhang "Quantitative photoacoustic imaging of two-photon absorption," Journal of Biomedical Optics 23(1), 016002 (2 January 2018). https://doi.org/10.1117/1.JBO.23.1.016002
Received: 5 July 2017; Accepted: 5 December 2017; Published: 2 January 2018
JOURNAL ARTICLE
11 PAGES


SHARE
Advertisement
Advertisement
Back to Top