Quantitative Photoacoustics Using the Transport Equation Simon - PowerPoint PPT Presentation

Quantitative Photoacoustics Using the Transport Equation Simon Arridge 1 Ben Cox 3 Teedah Saratoon 3 Joint work with: Tanja Tarvainen 1 , 2 1 Department of Computer Science, University College London, UK 2 Department of Physics and Mathematics, University of Eastern Finland, Finland 3 Department of Medical Physics, University College London, UK Colóquio Brasileiro de Matemática, July 29-August 2, 2013 S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 1 / 54

Outline Introduction 1 PhotoAcoustics 2 3 PhotoAcoustic Forward Model Quantitative PhotoAcoustic Tomography 4 5 Summary Acknowledgements 6 S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 2 / 54

Outline Introduction 1 PhotoAcoustics 2 3 PhotoAcoustic Forward Model Quantitative PhotoAcoustic Tomography 4 5 Summary Acknowledgements 6 S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 3 / 54

Introduction Outline Photoacoustic Imaging outline of photoacoustic imaging Photoacoustic image reconstruction Spectroscopic photoacoustic imaging Artefacts in photoacoustic imaging Quantitative Photoacoustic Imaging Models of light transport Multispectral reconstructions Unknown scattering: diffusion-based inversions Unknown scattering: using radiative transfer equation S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 4 / 54

Outline Introduction 1 PhotoAcoustics 2 3 PhotoAcoustic Forward Model Quantitative PhotoAcoustic Tomography 4 5 Summary Acknowledgements 6 S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 5 / 54

PhotoAcoustic Tomography PhotoAcoustic Signal Generation Spatially varying chromophore concentrations (naturally occuring or contrast agents) give rise to optical absorption in the medium. The absorption and scattering coefficients µ a and µ ′ s determine the fluence distribution Φ , and thence the absorbed energy distribution H . This energy generates a pressure distribution p 0 via thermalisation, which because of the elastic nature of tissue, then propagates as an acoustic pulse. The pulse is detected by a sensor resulting in the measured PA time series p ( t ) . S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 6 / 54

PhotoAcoustic Tomography PhotoAcoustic Imaging : 3 Modes S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 7 / 54

PhotoAcoustic Tomography PhotoAcoustic Signal Generation S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 8 / 54

PhotoAcoustic Tomography PhotoAcoustic Spherical BackProjection S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 9 / 54

PhotoAcoustic Tomography Optical part of the direct problem Optical part of the direct problem H ( r ) = µ a ( r )Φ( r ) absorbed absorption light energy density coefficient fluence = heat per unit volume S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 10 / 54

PhotoAcoustic Tomography Acoustic part of the direct problem Acoustic part of the direct problem p ( r ) | t = 0 = Γ( r ) H ( r ) = Γ( r ) µ a ( r )Φ( r ) Grüneisen parameter � � c 2 ∇ 2 − ∂ 2 p = 0 ∂ t 2 S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 11 / 54

PhotoAcoustic Tomography Fabrey-Perot Detector S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 12 / 54

PhotoAcoustic Tomography 3D PhotoAcoustic Scanner S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 13 / 54

PhotoAcoustic Tomography PAT Acoustic Inversion (Image Reconstruction) Initial value Problem � � c 2 ∇ 2 − ∂ 2 p = 0 ∂ t 2 p | t = 0 = Γ µ a Φ � ∂ p � = 0 � ∂ t � t = 0 S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 14 / 54

PhotoAcoustic Tomography PAT Acoustic Inversion (Image Reconstruction) Initial value Problem � � c 2 ∇ 2 − ∂ 2 p = 0 ∂ t 2 p | t = 0 = Γ µ a Φ � ∂ p � = 0 � ∂ t � t = 0 Boundary value Problem (t running backwards from T to 0) � � c 2 ∇ 2 − ∂ 2 p = 0 ∂ t 2 p ( r , t ) | t = T = 0 p obs ( r s , t ) p ( r , t ) | ∂ Ω = S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 14 / 54

PhotoAcoustic Tomography PAT Acoustic Inversion (Image Reconstruction) Initial value Problem � � c 2 ∇ 2 − ∂ 2 p = 0 ∂ t 2 p | t = 0 = Γ µ a Φ � ∂ p � = 0 � ∂ t � t = 0 Boundary value Problem (t running backwards from T to 0) � � c 2 ∇ 2 − ∂ 2 p = 0 ∂ t 2 p ( r , t ) | t = T = 0 p obs ( r s , t ) p ( r , t ) | ∂ Ω = S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 14 / 54

PhotoAcoustic Tomography Heterogeneous Sound Speed S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 15 / 54

PhotoAcoustic Tomography Spectroscopic PAT S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 16 / 54

PhotoAcoustic Tomography Spectroscopic PAT absorption at different wavelengths gives spectral images S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 17 / 54

PhotoAcoustic Tomography Spectroscopic PAT absorption at different wavelengths gives spectral images but fluence is also different at different wavelengths S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 17 / 54

PhotoAcoustic Tomography Spectroscopic PAT absorption at different wavelengths gives spectral images but fluence is also different at different wavelengths tumour type LS174T S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 17 / 54

PhotoAcoustic Tomography Spectral Distortion Spectral Distortion Spectrum corrupted by wavelength dependence of fluence S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 18 / 54

PhotoAcoustic Tomography Structural Distortion Structural Distortion Structural distortion due to non-uniform internal light fluence Structural distortion at each wavelength = spectral distortion at each point S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 19 / 54

Outline Introduction 1 PhotoAcoustics 2 3 PhotoAcoustic Forward Model Quantitative PhotoAcoustic Tomography 4 5 Summary Acknowledgements 6 S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 20 / 54

PhotoAcoustic Forward Model Second order wave equation for homogeneous media � � ∂ 2 ∇ 2 − 1 ∂ p ( r , t ) = 0 , p ( r , 0 ) = p 0 ( r ) , ∂ t p ( r , 0 ) = 0 (1) c 2 ∂ t 2 0 � � � p ( r , t ) = 1 g ( r , t | r 0 , t 0 ) ∂ p 0 ( r ) − p 0 ( r ) ∂ g ( r , t | r 0 , t 0 ) d r 0 (2) c 2 ∂ t 0 ∂ t 0 0 � sin ( c 0 kt ) c 2 0 e i k · ( r − r 0 ) d k g ( r , t | r 0 , t 0 ) = k = | k | (3) ( 2 π ) 3 c o k � � 1 p 0 ( r ) cos ( c 0 kt ) e i k · ( r − r 0 ) d k d r 0 p ( r , t ) = (4) ( 2 π ) 3 = F − 1 {F { p 0 ( r ) } cos ( c 0 kt ) } (5) Simple numerical algorithm using wave propagator cos ( c 0 kt ) (Cox and Beard 2005) S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 21 / 54

PhotoAcoustic Forward Model Linear Lossless Acoustic Equations p ( r , t ) = F − 1 {F { p 0 ( r ) } cos ( c 0 kt ) } (6) Cannot input a time-varying pressure, so no use for time-reversal imaging FFT ⇒ periodic boundary conditions (wave wrapping) Instead: solve equivalent first-order system ∂ u − 1 = ∇ p linearised momentum conservation (7) ∂ t ρ 0 ∂ρ = − ρ 0 ∇ · u linearised mass conservation (8) ∂ t c 2 p = P ρ linearised equation of state (9) p acoustic pressure, u particle velocity, ρ acoustic density S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 22 / 54

PhotoAcoustic Forward Model k-Space Acoustic Propagation Model S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 23 / 54

PhotoAcoustic Forward Model Modelling Acoustic Absorption Photoacoustic waves may contain frequencies much higher than conventional ultrasound imaging (tens of MHz) Acoustic absorption in soft tissue over ranges of interest (1-50 MHz or so) typically takes the form α = α 0 ω y , 1 ≤ y ≤ 1 . 5 What wave equations account for absorption like this? Can they be time-reversed to correct for absorption during image reconstruction? S.Arridge (University College London) QPAT using RTE IMPA, Rio de Janiero 27-07-13 24 / 54