Next Article in Journal
Theoretical Study on the Ultrafast Selective Excitation of Surface-Enhanced Coherent Anti-Stokes Raman Scattering Based on Fano Resonance of Disk-Ring Nanostructures by Shaped Femtosecond Laser Pulses
Next Article in Special Issue
Simulation Study of Acoustic-Resolution-Based Photoacoustic Microscopy for Imaging Complex Blood Vessel Networks
Previous Article in Journal
A Stabilized Single-Longitudinal-Mode and Wide Wavelength Tunability Erbium Laser
Previous Article in Special Issue
Multimodal In Vivo Imaging of Retinal and Choroidal Vascular Occlusion
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Enhancing Finite Element-Based Photoacoustic Tomography by Localized Reconstruction Method

1
Department of Radiology, Wake Forest University School of Medicine, Winston Salem, NC 27157, USA
2
Comprehensive Cancer Center, Wake Forest University School of Medicine, Winston Salem, NC 27157, USA
3
Department of Chemical, Biological, and Bioengineering, North Carolina Agricultural and Technical State University, Greensboro, NC 27411, USA
4
Center of Excellence in Product Design and Advanced Manufacturing, North Carolina Agricultural and Technical State University, Greensboro, NC 27411, USA
5
Department of Medical Engineering, College of Engineering and Morsani College of Medicine, University of South Florida, Tampa, FL 33612, USA
*
Author to whom correspondence should be addressed.
Photonics 2022, 9(5), 337; https://doi.org/10.3390/photonics9050337
Submission received: 3 April 2022 / Revised: 8 May 2022 / Accepted: 11 May 2022 / Published: 12 May 2022
(This article belongs to the Special Issue Photoacoustic Imaging for Biomedical Applications)

Abstract

:
Iterative reconstruction algorithm based on finite element (FE) modeling is a powerful approach in photoacoustic tomography (PAT). However, an iterative inverse algorithm using conventional FE meshing of the entire imaging zone is computationally demanding, which hinders this powerful tool in applications where quick image acquisition and/or a large image matrix is needed. To address this challenge, parallel computing techniques are proposed and implemented in the field. Here, we present an alternative approach for 2D PAT, which locoregionally reconstructs the region of interest (ROI) instead of the full imaging zone. Our simulated and phantom experimental results demonstrate that this ROI reconstruction algorithm can produce almost the same image quality as the conventional full zone-based reconstruction algorithm; however, the computation time can be significantly reduced without any additional hardware cost by more than two orders of magnitude (100-fold). This algorithm is further applied and validated in an in vivo study. The major vessel structures in a rat’s brain can be imaged clearly using our ROI-based approach, coupled with a mesh of 11,801 nodes. This novel algorithm can also be parallelized using MPI or GPU acceleration techniques to further enhance the reconstruction performance of FE-based PAT.

1. Introduction

Photoacoustic tomography (PAT), also referred to as optoacoustic tomography (OAT), is an emerging non-ionizing and non-invasive imaging modality in the field of biomedical imaging based on the physical principle of photoacoustic effect, where light energy is converted into acoustic energy due to optical absorption and localized thermal expansion in biological tissues [1]. As a hybridized imaging technique, PAT combines both merits of optical imaging modality and the ultrasound imaging modality. By collecting ultrasound waves generated directly from laser-illuminated target tissues (tumor, blood vessels, organs, etc.), it can achieve a much better signal to noise level for greatly enhanced spatial resolution in deep tissues (up to sub 100 µm at cm imaging depths) as ultrasound scattering is two orders less than optical scattering in biological tissues [2,3,4]. Using laser–tissue interaction, the optical properties (absorbed optical energy density, intrinsic absorption coefficient, etc.) and associated physiological/functional information (oxygen saturation, hemoglobin levels, etc.) can be resolved with single or multi-wavelength PAT. These intrinsic optical properties and physiological/functional information may be important for accurate diagnostic decision-making in the early stage of diseases. As a promising imaging tool in life science and clinic application, PAT has shown its great potential to detect breast cancer, visualize joint abnormality, assess vascular and skin diseases, monitor neural activities, including epilepsy in small animals, sound out fluorescent proteins, and evaluate exogenous contrast agents in molecular imaging [5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20].
Thus far, several algorithms are implemented to reconstruct qualitative or quantitative photoacoustic images, i.e., the distribution of absolute absorbed optical energies or initial pressure waves that are proportional to the absorbed optical energies in the tissue. The available algorithms include analytical model-based algorithms and optimization-based iterative algorithms, such as backprojection, Fourier transform, P-transform, the k-wave method, statistical approaches, finite element (FE) based modeling, the deep learning algorithm, etc. [21,22,23,24,25,26,27,28,29,30,31,32,33,34,35,36,37,38] Analytical approaches are more computationally efficient; however, they usually provide qualitative photoacoustic images, and accurate reconstruction using analytical approaches requires idealized imaging scenarios, such as acoustic homogeneity and a regular boundary. In contrast, the current finite element-based PAT algorithms (FE-PAT) used in a number of studies, including FE-PAT in the frequency domain and FE-PAT in the time domain [19,20,21,37,38,39,40,41,42,43,44,45,46], are optimization-based iterative algorithms and seem to provide unrivaled advantages to accommodate geometric irregularity and complex source representations. These algorithms reduce the reconstruction errors in practical situations such as limited datasets or detecting locations, and quantitatively recover multiple parameters (absorbed optical energy, optical absorption coefficient, total hemoglobin levels, oxygen saturation, etc.) from one or a set of discretized photoacoustic equations and optical propagation models [39,40]. Moreover, acoustic heterogeneity in biological tissues can be recovered simultaneously with optical property by the finite element-based algorithms when acoustic heterogeneity is taken into consideration in the photoacoustic wave equation [1,41,42].
However, the FE-based algorithms are iterative optimization algorithms and are computationally demanding using conventional FE meshing of the entire imaging zone, particularly in inverse modeling where the coupling matrix is very dense compared to the sparse matrix in forwarding modeling. The computational cost is even striking for FE-PAT in the time domain, where the differentiation of signals in the time domain is needed instead of a set of linear equations used by the frequency-domain counterpart [38]. It usually requires hours to days for one set of images to be reconstructed with a typical mesh size for in vivo applications such as finger joint imaging or breast imaging, which is unrealistic for potential clinical applications where timely imaging and/or a large image matrix is needed. As the scale increases to improve the spatial resolution, the cost of computational memory and time grows dramatically. To overcome this challenge, parallel computation techniques are proposed and implemented using a Message Passing Interface (MPI) among multiple CPU processors or a Graphics Processing Unit (GPU) [43,44,45]. However, the improvement of computation performance is linearly dependent on the number of processors or threats and, herein, the hardware cost in parallel computation techniques. It is still highly demanded to develop a powerful yet computationally efficient reconstruction algorithm for PAT.
To address the need in the field, here we propose an alternative approach based on a novel algorithm to locoregionally reconstruct the region of interest (ROI) instead of the full imaging zone used in the conventional FE-PAT algorithm. The idea of ROI reconstruction is studied in the field of X-ray-based computed tomography (CT) [46,47]; however, no ROI-based reconstruction algorithms are proposed or discussed in the field of photoacoustic tomography. Furthermore, truncated datasets are used in X-ray-based CT, leading to localized reconstruction in the ROI. In contrast, a cost function with a penalty term defined in the ROI will be used in our localized reconstruction algorithm, which could be used by other optimization-based reconstruction algorithms where a cost function needs to be iteratively minimized to obtain the optimized solution. A two-dimensional (2D) algorithm is implemented to demonstrate the efficiency and accuracy of this proposed algorithm. Simulation data and phantom experimental data are used in the developed 2D algorithm to validate the accuracy and performance of this algorithm. This algorithm is further applied and validated in an in vivo study of a rat’s brain. This novel algorithm can also be parallelized using multiple CPU or GPU accelerated high-performance computing techniques to further enhance the reconstruction performance of FE-based photoacoustic tomography.

2. Materials and Methods

2.1. Reconstruction Algorithm

Our 2D reconstruction algorithm is based on the following Helmholtz photoacoustic wave equation and the absorbing boundary condition in the frequency domain [37].
2 P ( r , ω ) + k 2 P ( r , ω ) = i k c β Φ ( r ) C p P · n ^ = η P + γ 2 P Φ 2
where η = i k 0 3 / 2 ρ + i 3 / 8 k 0 ρ 2 1 i / k 0 ρ ; γ = i / 2 k 0 ρ 2 1 i / k 0 ρ ; ϕ is the angular coordinate at a radial position ρ ; k = ϖ / c is the acoustic wave number defined by the angular frequency ω and the average speed of the acoustic wave in the media c ; β is the thermal expansion coefficient; C P is the specific heat at constant pressure.
The finite element discretization of the above Equation (1) can be written as the following linear equations by choosing triangle elements and the associated basis function ψ i .
[ A ω ] { p } = { B ω } { φ }
where
A i j ω = V ( ψ j ψ i ) d V V k 2 ψ j ψ i d V Γ ( η ψ j + γ 2 ψ j ϕ 2 ) ψ i d s
B i j ω = i k c β C P V ψ j ψ i d V
p = { p 1 , p 2 , , p N } T φ = { φ 1 , φ 2 , , φ N } T
The above Equation (2) is used as the forward model in the FE-PAT algorithm, which will compute the spatial variation of the acoustic pressure wave from the initial guesses or iterative values of the absorbed optical energy distribution in image reconstruction. To obtain a matrix equation capable of inversely solving the unknown distribution of absorbed optical energy, a modified objective function is introduced for localized reconstruction
F ( φ * ) = p c ( φ * ) p m 2 + λ φ * φ * a 2
where φ * = L φ ; φ * a = L φ a ; [ L ] K × N is an introduced filtering matrix; N is the total number of FE nodes in the entire imaging zone; K is the total number of FE nodes in the region of interest (ROI) containing the photoacoustic targets. φ and φ a are the current (calculated) and the actual distribution (unknown) of absorbed optical energy in the entire imaging zone. φ a is usually set equal to the calculated distribution of absorbed optical energy at the previous iteration. L i j = 1 if the ith node in the preselected ROI is the same as the node j in the entire imaging zone; otherwise, L i j = 0 . p c and p m are the calculated and measured acoustic pressures at each boundary site among the total M locations. The edge of the preselected ROI can be estimated in several different ways. It can be chosen by an educated guess from the approximated location and size of the targets. It can be quickly estimated from the photoacoustic sinogram data. It can also be precisely preselected from localizer/scout images, which can be acquired by coarse mesh-based quick reconstruction, fast beamforming methods such as Delay and Sum (DAS), or guidance images such as ultrasonography if an ultrasound-photoacoustic dual-modality imaging system was used.
The optimized solution in the inverse problem is achieved by minimizing the above modified objective function. By applying Taylor’s expansion and considering F / φ * = 0 at the minimum of the function F ( φ * ) (optimized solution), an inverse matrix can be derived
[ J T J + λ I ] { Δ φ ˜ * } = [ J T ] { p m p c }
where Δ φ ˜ * = φ * φ * a is the update vector for the absorbed optical energy in the ROI; JM×K is the Jacobian matrix defined by p i / φ j * at each frequency in the ROI; M is the number of the measurement sites; K is the number of mesh nodes in the ROI; λ is the regularization parameter determined by combined Marquardt and Tikhonov regularization schemes; I is the identity matrix.
In our algorithm, the image formation task is used to update the distribution of absorbed optical energy via the iterative solution of Equations (2) and (4) so that the objective function F ( φ * ) defined in Equation (3) can be minimized. During the iteration, φ * a in the ROI is chosen as the absorbed optical energy in the previous iteration (or initial guess in the first iteration), while the current absorbed optical energy φ * in the ROI can be updated from the solved Δ φ ˜ * in Equation (4). The current absorbed optical energy outside of ROI is considered the background, which is simply approximated by the mean value of φ * around the edge of the ROI. Histogram analysis can be used in more complicated cases. The regularization parameter λ is usually reduced at each iterative step to balance the speed of convergence and accuracy around the minimum of the modified objective function F ( φ * ) .

2.2. Numerical Simulations

The above algorithm was first validated with simulation studies. To illustrate the performance and accuracy of the above algorithm, a 2D circular mesh with 4489 nodes and 8736 triangle elements was used in all the calculations on a desktop PC with an Intel I7-4790 CPU, 3.6 GHz, and 8 GB RAM. The simulated acoustic pressure on the boundary was generated by the forward calculation with Equation (2). The simulated photoacoustic target was centered at x = 3 mm, y = 0 mm with a diameter of 4 mm and the optical absorption coefficient was set as 3 (arbitrary unit). The simulated photoacoustic target was embedded in the background, which had a radius of 15 mm and a uniform optical absorption coefficient of 1 (arbitrary unit). Uniform light illumination was assumed in the simulation, and the light fluence was assumed as I in an arbitrary unit for simplicity. The speed of sound was uniformly assumed as 1495 m/s in the simulated tissue. A total of 60 detectors were used to uniformly cover the full circle on the boundary, 15 mm away from the center. Photoacoustic signals at 50 frequency elements from 30 Khz–1.5 Mhz were generated from the forward model. An additional 5% white noise was applied to simulate the signal noise in real applications.

2.3. Phantom Experiments

Phantom experiments were further conducted, and the experimental data were fed into the developed algorithm to further validate the accuracy and performance of the above algorithm. The 2D mesh and the PC used in the experimental studies were the same as the mesh and PC used in the simulation studies. The experimental data were obtained from tissue mimic phantom scanned by a system previously described in detail. [48] Laser pulses at a selected wavelength (700 nm) from a Surelite OPO Plus (Continuum Inc., San Jose, CA, USA) were expanded to fully illuminate the phantom. The laser had a repetition rate of 20 Hz and a pulse duration of 10 ns. The energy density was controlled below 20 mJ/cm2. The generated acoustic signals were detected by a piezoelectric transducer with a central frequency of 1.0 MHz (Valpey Fisher, Hopkinton, MA, USA). The transducer was rotated around the center at a distance of around 33.5 mm to cover 60 sites on the circular boundary. In the experiment, we embedded a target with a diameter of 3 mm in a 30 mm diameter solid cylindrical background phantom, which was fabricated using a previously reported procedure. The optical absorption coefficient for the background phantom and target were 0.01/mm and 0.07/mm, respectively. The reduced optical scattering coefficient was 1.0/mm for the target and the background.

2.4. In Vivo Experiments

The developed RO-based localized reconstruction was further applied and validated in an in vivo experimental study. A Harlan Sprague Dawley rat was used in this in vivo experiment. The in vivo experimental protocol and procedure were approved by the IACUC committee at the University of Florida. The experiment system and data acquisition procedure were previously described in detail [21]. Briefly, the rat was anesthetized, and its head was placed into the imaging region through an opening covered with a piece of polyethylene membrane at the bottom of a water tank. In the photoacoustic tomography system used in this study, the light source was a pulsed Nd: YAG laser (Altos, Bozeman, MT, USA) with a 4 ns pulse duration and 10 HZ repetition rate. The diameter of the laser beam was expanded to 30 mm by a lens, and the light power was adjusted to 15 mJ cm−2 at the surface of the rat’s brain. An acoustic transducer (Valpey Fisher, Hopkinton, MA, USA) was immersed in the water tank and driven by a motorized rotator to receive light-induced acoustic signals. The central frequency of the transducer was 1.0 MHz, and the bandwidth was 80% at 6 dB. The acoustic signals at 175 detecting locations uniformly distributed along the boundary covering a total of 350 degrees were collected and used for later image reconstruction.

3. Results

3.1. Simulation Results

The simulated data were fed into both full zone-based FE-PAT and localized FE-PAT algorithms. Three different ROIs containing the target were chosen in the localized FE-PAT algorithm: R = 10.0 mm, R = 7.5 mm, and R = 4.0 mm. An initial value of 0.01 was used as the absorbed optical energy at each node for all the reconstruction processes in our simulation studies, including the reconstruction process based on full zone-based FE-PAT and reconstruction process based on localized FE-PAT in various ROI sizes. The reconstructed photoacoustic images were obtained from the fifth iteration when the reconstruction error was low and stable. Furthermore, the reconstructed photoacoustic image (distribution of absorbed optical energies) was spatially filtered to reduce the reconstruction noise associated with finite element modeling. The photoacoustic image reconstructed by conventional full zone-based FE-PAT is shown in Figure 1a, while the photoacoustic images reconstructed by localized FE-PAT in various ROI sizes are shown in Figure 1b–d. The accuracy of the full zone-based FE-PAT in the frequency domain was previously evaluated and reported [37]. As shown in Figure 1a, the simulated target and background were clearly reconstructed by the full zone-based FE-PAT. From the results in Figure 1b–d, we can see that the target can also be clearly reconstructed using localized FE-PAT in various ROI sizes (R = 10.0 mm, 7.5 mm, and 4.0 mm). The targets reconstructed by localized FE-PAT in the ROI, ranging from 10.0 mm to 4.0 mm in radius (Figure 1b–d), were visually identical to the target reconstructed by conventional full zone-based FE-PA (Figure 1a) in terms of the positions and the sizes of the target. The positions and the sizes of the reconstructed targets were also visually in agreement with the actual location and size of the simulated target (D = 4.0 mm centered at x = 3 mm, y = 0 mm). From Figure 1a–d, we can also see that the signal variation in the background (noise) appeared in the entire image zone for the full zone-based FE-PAT algorithm while it only appeared in the chosen ROI for the localized FE-PAT algorithm, as the average signal intensity was assigned to the background that is outside of the ROI. However, the averaged signal intensity in the background (inside and outside of the ROI) appears to be identical for both the full zone-based FE-PAT algorithm and localized FE-PAT algorithm in various ROI sizes.
Further quantitative analysis of the reconstructed results is shown in Figure 2, where the absorbed optical energy (photoacoustic intensity) along a horizontal transect line through the center of the target (x = 3 mm, y = 0 mm) was plotted. As shown in Figure 2, the distribution of the absorbed optical energies in the target reconstructed by localized FE-PAT in various ROI sizes (R = 4.0 mm to R = 10.0 mm in Figure 2) was in great agreement with that reconstructed by the full zone-based FE-PAT (Full-zone in Figure 2). The reconstruction error between the full zone-based FE-PAT algorithm and the localized FE-PAT algorithm in different ROI sizes was calculated and listed in Table 1 by applying the reported error measurement formula along the transect line. The reconstruction errors were less than 5% (in the range of 2.84–4.47%) for ROI sized from 4.0 mm to 10.0 mm in radius, demonstrating the reconstruction accuracy of localized FE-PAT in comparison with the full zone-based FE-PAT. The diameter of the reconstructed target was further analyzed using the full width at half maximum (FWHM) method. The measured diameter of the reconstructed target, the computational time, and the mesh information for different FE-PAT algorithms used in this simulation study are all listed in Table 1. As shown in Table 1, the recovered target size (the diameter) was between 4.08 mm and 4.14 mm using localized FE-PAT in various ROI sizes and has an error between 3.03% and 4.55%, in comparison to the recovered target size (3.96 mm) using the full zone-based FE-PAT algorithm. Compared to the actual diameter (4.0 mm) of the simulated target, it is obvious that the target was reconstructed with high accuracy in terms of the location/size by all the imaging algorithms. The results of the quantitative analysis demonstrated the high reconstruction accuracy of the localized FE-PAT algorithm in various ROI sizes, in terms of the reconstructed photoacoustic intensities (the absorbed optical energies) and the size of the reconstructed target, in comparison to the full zone-based FE-PAT algorithm.
The efficiency of the localized FE-PAT algorithm was further analyzed based on the computational time listed in Table 1. As shown in Table 1, it took only 49–1144 s to generate one full set of images using the localized FE-PAT algorithm, depending on the sizes of the ROI, while it took 7745 s to generate one full set of images using the full zone-based FE-PAT algorithm. The acceleration rate was increased from 7-fold to 24-fold and 158-fold as the radius of the ROI was reduced from 10.0 mm to 7.5 mm and 4.0 mm. It clearly demonstrated that the localized FE-PAT algorithm offers a significant improvement in computation time compared to the full zone-based FE-PAT algorithm, while the reconstruction error remains very low. The acceleration performance is further plotted and fitted in Figure 3. Using a one-term power series model, the acceleration performance can be estimated by α 4.696 or β 2.364 , where α and β are the ratios of the radius and number of nodes between the ROI and the entire imaging field, respectively; β α 2 is the number of nodes and is almost squarely proportional to the radius of the ROI for a relatively uniform 2D mesh. It is worth noting that the computational frame of 3D FE-PAT is based on the forward equations and inverse equations similar to the equations used for 2D FE-PAT (Equations (2) and (4)), although the formula to calculate the stiffness matrix (A in Equation (2)) is different in 3D PAT. As such, this localized reconstruction approach can be further adapted to 3D PAT, and the acceleration performance can be predicted from the power series model established in 2D-based localized reconstruction. The expected acceleration performance for 3D mesh-based localized FE-PAT can be estimated by α 7 , as the number of nodes in a relatively uniform 3D mesh is almost cubically proportional to the radius of the ROI, i.e., 2.364 × 3 ≈ 7. It is expected that the improved imaging speed might even be striking in 3D mesh-based localized FE-PAT, as around 1/3 ROI ( α 1 / 3 ) can lead to more than 1000-fold (37) acceleration.

3.2. Phantom Experimental Results

Further experimental validation with tissue mimic phantom is shown in Figure 4 and Figure 5. The experimentally measured signals on the boundary were remapped to the corresponding virtual detectors on an R = 15 mm mesh with high accuracy, as the far-field condition remained valid for the generated and propagated photoacoustic/ultrasound waves. The signals were transformed into the frequency domain by Fourier transform, and the first 50 frequency elements (30 Khz–1.5 Mhz) were used in the image reconstruction. The speed of sound was assumed uniformly in the phantom during the image reconstruction, and the value 1510 m/s was used for the best reconstruction results. As in the simulation study, an initial value of 0.01 was used as the absorbed optical energy at each node for both the full zone-based FE-PAT and localized FE-PAT. The reconstructed photoacoustic images were obtained from the fifth iteration when the reconstruction error was low and stable. Furthermore, the reconstructed photoacoustic image (distribution of absorbed optical energies) was spatially filtered to reduce the reconstruction noise associated with finite element modeling. A threshold of zero was used to remove the negative values, as the photoacoustic intensity (absorbed optical energy) should be positive.
The reconstructed images of the phantom are shown in Figure 4a,b using the full zone-based FE-PAT algorithm and localized FE-PAT algorithm in a ROI with a radius of 4.0 mm. From the results, we can see that the target can be clearly reconstructed by the localized FE-PAT algorithm (Figure 4b) and was visually identical to the image obtained from the full zone-based FE-PAT algorithm (Figure 4b). Further quantitative analysis of the reconstruction results is shown in Figure 5, where the absorbed optical energy (photoacoustic intensity) along a horizontal transect line through the center of the target (around x = −0.5 mm, y = −2 mm) is plotted. As shown in Figure 5, the distribution of the absorbed optical energies reconstructed by the localized FE-PAT was in great agreement with that reconstructed by the full zone-based FE-PAT. The reconstruction error between the full zone-based FE-PAT algorithm and the localized FE-PAT algorithm was 3.72%, applying the reported error measurement formula along the transect line with values close to zeros excluded. As in the simulation study, the FWHM method was used to determine the diameter of the reconstructed target in the phantom. The measured diameter of the reconstructed target, the computational time, and the mesh information are all listed in Table 2. As shown in Table 2, the recovered target size (the diameter) was determined as 3.22 mm using the localized FE-PAT, and it had an error of 0.31% in comparison to the recovered target size (3.23 mm), which was generated using the full zone-based FE-PAT algorithm. Compared to the actual diameter (3 mm) of the target in the phantom, it seems that the target was reconstructed accurately in terms of the location/size by both the localized FE-PAT and the full zone-based FE-PAT algorithm. The results of the quantitative analysis demonstrated the high reconstruction accuracy of the localized FE-PAT algorithm, in terms of the reconstructed photoacoustic intensities (the absorbed optical energies) and the size of the reconstructed target in the phantom, in comparison to the full zone-based FE-PAT algorithm.
The efficiency of the localized FE-PAT algorithm was further demonstrated based on the computational time listed in Table 2. As shown in Table 2, it took only 49 s to generate one full set of images using the localized FE-PAT algorithm, while it took around 7730 s to generate one full set of images using the full zone-based FE-PAT algorithm. The reconstruction speed was accelerated 158-fold using the localized FE-PAT algorithm (R = 4.0 mm) compared to the conventional full zone-based FE-PAT algorithm. The results clearly demonstrate that the localized FE-PAT algorithm offers a significant improvement in computation time compared to the full zone-based FE-PAT algorithm, while the reconstruction error remains very low for phantom experimental data.

3.3. In Vivo Experimental Results

The localized reconstruction algorithm was further applied and validated in an in vivo experimental study. The measured acoustic signals at 174 locations were transformed into frequency-domain signals by Fourier transform. The first 50 frequency elements (41.5 Khz–2.075 Mhz) and the first 100 frequency elements (41.5 Khz–4.15 Mhz) were used in two separated image reconstruction cases with the goal of comparing the impact of the frequency range on the quality of the reconstructed image. The frequency-domain signals were fed into the developed algorithm using a 2D circular mesh with a radius of 15 mm. The 2D mesh had 11,801 nodes and 23,240 triangle elements. A ROI with a radius of 6 mm was selected based on an educational guess from the location and size of the rat’s brain and confirmed using the quick beamforming method Delay and Sum (DAS). A desktop PC with an Intel I7-4790 CPU, 3.6 GHz, and 8 GB RAM was used in the reconstruction of the in vivo rat’s brain.
The reconstructed images are shown in Figure 6a,b, where a threshold of zero was used to remove the negative values as the photoacoustic intensity (absorbed optical energy) should be positive. As shown in the images, the blood vessels in the rat brain and the middle cerebral artery were clearly imaged by our localized reconstruction algorithm in both frequency ranges (Figure 6a,b). However, there was no visual difference in the reconstructed results between 50 frequency elements (41.5 Khz–2.075 Mhz) and the first 100 frequency elements (41.5 Khz–4.15 Mhz). Moreover, the reconstruction time for the localized reconstruction algorithm with 50 frequency elements was around 11.5 min (685 s), while the reconstruction time for the localized reconstruction algorithm with 100 frequency elements took around 20.5 min (1215 s). The computation performance of the counterpart full zone-based FE-PAT was not available, as the dense mesh (11,801 nodes) used for our in vivo studies is out of the computational capability of full zone-based FE-PAT on a regular PC. The quantitative profile of the reconstructed photoacoustic images across the middle cerebral artery (MCA) further demonstrated that there is no significant difference in the reconstruction results between 50 frequency elements (41.5 Khz–2.075 Mhz) and the first 100 frequency elements (41.5 Khz–4.15 Mhz). The diameter of the reconstructed MCA was further analyzed using the FWHM method, and it was 0.387 and 0.381 mm for 50 frequency elements (41.5 Khz–2.075 Mhz) and the first 100 frequency elements (41.5 Khz–4.15 Mhz), respectively. The reason that the recovered vessel size was larger than the exact size and the increased frequency range did not improve the spatial resolution might be due to the limited spatial resolution of the mesh we used in our reconstruction. The low nominal frequency (1 MHz) of the transducer used in our experiment may result in lower spatial resolution as well.

4. Discussion

Computational cost is one of the major issues that hinder the translation of the powerful finite element-based algorithms (FE-PAT algorithms in frequency domain and time domain) in practical applications of photoacoustic tomography. It usually requires hours to days for one set of images to be reconstructed with a typical mesh size for in vivo studies such as finger joint imaging or breast imaging [37,38,44]. Due to the high computational cost, it is unrealistic to apply this powerful imaging tool to any clinical applications where timely imaging and/or a large image matrix is needed. The proposed method here is based on localized optimization in the region of interest (ROI), where the cost function was modified with a penalty term restricted in the ROI. It is worth noting that only a fraction of pixels/voxels need to be updated in the inverse process using our localized reconstruction algorithm, while all the pixels/voxels in the imaging zone must be updated indistinguishably by currently available full zone-based FE-PAT algorithms as the cost functions are defined in the entire imaging zone. In most practical PAT applications, detectors are usually placed in the far-field to allow for better spatial resolution in the targets. Furthermore, in most commercial PAT scanners, PAT is integrated with ultrasonography, where ultrasonography provides tissue anatomy, and PAT provides molecular/functional information in the region of interest. As such, the ROI containing targets to be photoacoustically imaged is much smaller than the entire imaging zone, and this localized algorithm may significantly reduce the computational cost due to restricted inverse modeling in ROI only. Other computational efficient optimization methods, such as the finite difference method, can also be used to quantitatively reconstruct photoacoustic images. However, they are not as powerful as FEM algorithms, as the finite difference method is based on regular rectangular grids and can easily run into problems handling curved boundaries.
The localized FE-PAT algorithm was implemented and then validated with both numerical simulations and phantom experiments in our studies. In our studies, we focused on the comparison between our localized reconstruction algorithm and conventional full zone-based reconstruction algorithm. The dataset (frequency range, number of detectors, etc.) and the same mesh were used in the localized reconstruction algorithm and conventional full zone-based reconstruction algorithm. As seen in the simulation study, the reconstruction error of the photoacoustic image (absorbed optical energies) by the localized FE-PAT algorithm in various ROI sizes was less than 5% in comparison with the image reconstructed by the conventional full zone-based FE-PAT. The error of the target size recovered by the localized FE-PAT algorithm in various ROI sizes was also less than 5% in comparison to the actual size or the target reconstructed by the conventional full zone-based FE-PAT algorithm. The very low reconstruction error demonstrates the accuracy of the localized FE-PAT algorithm. The phantom experimental study further confirmed the accuracy of the localized FE-PAT algorithm. However, the reconstruction speed can be accelerated by up to 158-fold, and the reconstruction time of the FE-PAT algorithm can be reduced from around 2–3 h to less than 1 min. As shown in Figure 1 and Figure 2, the size of the chosen ROI does not affect the image quality if the targets are contained in the ROI, but it may affect the imaging speed. As such, our algorithm allows for high flexibility in selecting ROI. The localized reconstruction algorithm was further applied and validated in an in vivo experimental study of a rat’s brain. Our results demonstrate that the major vessels can be imaged clearly using our localized reconstruction algorithm coupled with a mesh of 11,801 nodes and 23,240 triangle elements. The total size of the image depends on the nodes used in the reconstruction. As such, the reconstructed image of the in vivo rat brain using 11,801 nodes was larger than the reconstructed image of simulation/phantom experimental studies using 4489 nodes. The size of the middle cerebral artery (MCA) was found to be around 0.381–0.387 mm. Similar to the conventional full zone-based FE-PAT, the spatial resolution of our localized reconstruction algorithm depends on both the size of each finite element and the frequency range of detected acoustic signals. The frequency components chosen in our reconstruction covered the major bandwidth of the generated acoustic signals from the simulated target and ex vivo phantom subject and demonstrated high image quality in our studies compared to the ground truth of the targets. Due to the fact that a transducer with a low nominal frequency was used in in vivo study, both frequency ranges (41.5 Khz–2.075 Mhz and 41.5 Khz–4.15 Mhz) covered the major bandwidth of the detected acoustic signals from the rat’s brain and resulted in similar spatial resolution (0.381 mm vs. 0.387 mm). Similar to the conventional full zone-based FE-PAT, the spatial resolution of our localized reconstruction algorithm can be improved using a denser mesh and a transducer with a higher nominal frequency.
It is worth noting that the targets in our studies were qualitatively reconstructed with high accuracy in terms of the location/size by the localized FE-PAT algorithm. However, the quantitative values of the absorbed optical energies (photoacoustic intensities) in the targets might be underestimated for both the full zone-based FE-PAT algorithm and localized FE-PAT algorithms. For example, the absorbed optical energies in the simulated target were around 2.8, determined using both full zone-based FE-PAT algorithm and localized FE-PAT algorithms, which is lower than the actual absorbed optical energy (3 in arbitrary unit) in the target. As the FE-based reconstruction algorithm is one of the iterative optimization algorithms, the target can be qualitatively reconstructed with high accuracy in the first several iterations, but the quantitative accuracy of the absorbed optical energies in the target may require more iterations and longer computational time. Furthermore, the localized reconstruction algorithm in the ROI is a variant of the conventional finite element-based photoacoustic tomography. Similar to other finite element-based PAT algorithms, it could provide not only qualitative photoacoustic images (pixel intensities proportional to the absorbed optical energy at each location) but also quantitative values (the absolute values of absorbed optical energy or inherent optical absorption coefficient at each location) from one or a set of discretized photoacoustic equations and optical propagation models. It can also potentially accommodate geometric irregularity and complex source representations, reduce the reconstruction errors in practical situations such as limited datasets or detecting locations, and quantitatively recover multiple parameters (optical and acoustic contrasts, total hemoglobin concentrations, oxygen saturation, etc.). Moreover, the localized reconstruction algorithm is capable of significantly accelerating the computational speed of finite element-based photoacoustic tomography without any additional hardware cost. It can provide quantitative images in a timely manner (within minutes, seconds, or real-time). For a small target area or ROI in molecular/functional imaging, our localized reconstruction algorithm may have a temporal resolution of seconds. However, using an additional computational approach such as integration with multiple CPU or GPU-based parallelized computing techniques, our algorithm could also potentially provide real-time imaging capability and provide fast quantitative images for scenarios where large ROI or real-time/video-rate imaging is needed. As our algorithm is a quantitative imaging algorithm, it can provide quantitative information in molecular/functional imaging (concentration of molecular imaging probe, oxygen saturation, total hemoglobin levels, etc.), similar to conventional FE-PAT algorithms. The reconstructed images from our ROI-based approach can be further quantitatively analyzed like other quantitative imaging algorithms. In contrast, most commonly used algorithms such as filtered backprojection usually provide qualitative images, i.e., the initial peak of the generated acoustic signals proportional to the absorbed optical energy in the local tissue, or quantitative information based on the ratio of the qualitative images. The significantly enhanced computational efficiency of using our localized FE-PAT may facilitate the translation of this powerful model-based imaging tool to real-world preclinical/clinical studies with a large number of subjects, including instant tumor detection using handheld ultrasound-photoacoustic probes, the localization of tumor margins during surgery, and image-guided cancer treatments.

5. Conclusions

We presented a localized FE-PAT algorithm and successfully demonstrated its accuracy and performance in comparison with conventional full zone-based FE-PAT algorithms. The computation time was dramatically reduced, and image quality remained almost the same. The work presented here suggests that this localized FE-PAT algorithm can accelerate image reconstruction a hundredfold without any additional hardware cost, making it promising for higher resolution, higher dimension, or real-time imaging. The localized reconstruction algorithm was further applied and validated in an in vivo experimental study of a rat’s brain. Our results demonstrate that the complex vessels in the rat’s brain can be imaged clearly within around 11.5 min using our localized reconstruction algorithm coupled with a mesh of 11,801 nodes and 23,240 triangle elements. Furthermore, this localized reconstruction algorithm can be parallelized using multiple CPU or GPU schemes to further enhance the reconstruction performance of FE-based photoacoustic tomography, which may facilitate the translation of this promising imaging tool into clinical studies with a large number of subjects.

Author Contributions

Y.S. proposed the algorithm and wrote the manuscript. Y.S. designed and performed all the experiments, collected raw data, analyzed the data, and prepared the figures for the manuscript. Y.S. and H.J. conceived and supervised the project. All authors contributed to the critical reading and writing of the manuscript. All authors have read and agreed to the published version of the manuscript.

Funding

This research received no external funding.

Institutional Review Board Statement

The study was conducted in accordance with the Declaration of Helsinki and approved by the IACUC committee at the University of Florida (protocol code D710 and 18 August 2009).

Informed Consent Statement

Not applicable.

Data Availability Statement

The data that supports the findings of this study are available within the article.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Sun, Y.; Jiang, H.; O’Neill, B.E. Photoacoustic Imaging: An Emerging Optical Modality in Diagnostic and Theranostic Medicine. J. Biosens. Bioelectron. 2011, 2, 1000108-1. [Google Scholar] [CrossRef]
  2. Wang, X.; Pang, Y.; Ku, G.; Xie, X.; Stoica, G.; Wang, L.V. Noninvasive laser-induced photoacoustic tomography for structural and functional in vivo imaging of the brain. Nat. Biotechnol. 2003, 21, 803–806. [Google Scholar] [CrossRef] [PubMed]
  3. Wang, L.V.; Yao, J. A practical guide to photoacoustic tomography in the life sciences. Nat. Methods 2016, 13, 627–638. [Google Scholar] [CrossRef]
  4. Cox, B.T.; Laufer, J.; Arridge, S.R.; Beard, P.C. Quantitative spectroscopic photoacoustic imaging: A review. J. Biomed. Opt. 2012, 17, 0612021–06120222. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  5. Oraevsky, A.A.; Karabutov, A.A.; Solomatin, S.V.; Savateeva, E.V.; Andreev, V.A.; Gatalica, Z.; Singh, H.; Fleming, R.D. Laser optoacoustic imaging of breast cancer in vivo. Proc. SPIE 2001, 4256, 6–16. [Google Scholar] [CrossRef]
  6. Manohar, S.; Kharine, A.; Van Hespen, J.C.G.; Steenbergen, W.; van Leeuwen, T. Photoacoustic mammography laboratory prototype: Imaging of breast tissue phantoms. Phys. Med. Biol. 2005, 50, 2543–2557. [Google Scholar] [CrossRef] [PubMed]
  7. Kolkman, R.G.M.; Klaessens, J.H.G.M.; Hondebrink, E.; Hopman, J.C.W.; De Mul, F.F.M.; Steenbergen, W.; Thijssen, J.M.; van Leeuwen, T. Photoacoustic determination of blood vessel diameter. Phys. Med. Biol. 2004, 49, 4745–4756. [Google Scholar] [CrossRef]
  8. Siphanto, R.I.; Thumma, K.K.; Kolkman, R.G.M.; van Leeuwen, T.; De Mul, F.F.M.; Van Neck, J.W.; Van Adrichem, L.N.A.; Steenbergen, W. Serial noninvasive photoacoustic imaging of neovascularization in tumor angiogenesis. Opt. Express 2005, 13, 89–95. [Google Scholar] [CrossRef]
  9. Yang, S.; Xing, D.; Zhou, Q.; Xiang, L.; Lao, Y. Functional imaging of cerebrovascular activities in small animals using high-resolution photoacoustic tomography. Med. Phys. 2007, 34, 3294–3301. [Google Scholar] [CrossRef]
  10. Zhang, E.Z.; Laufer, J.; Beard, P. Three-dimensional photoacoustic imaging of vascular anatomy in small animals using an optical detection system. SPIE 2007, 6437, 64370S. [Google Scholar] [CrossRef]
  11. Hoelen, C.G.A.; De Mul, F.F.M.; Pongers, R.; Dekker, A. Three-dimensional photoacoustic imaging of blood vessels in tissue. Opt. Lett. 1998, 23, 648. [Google Scholar] [CrossRef] [PubMed]
  12. Emelianov, S.Y.; Li, P.; O’Donnell, M. Photoacoustics for molecular imaging and therapy. Phys. Today 2009, 5, 34–39. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  13. Razansky, D.; Distel, M.; Vinegoni, C.; Ma, R.; Perrimon, N.; Köster, R.W.; Ntziachristos, V. Multispectral opto-acoustic tomography of deep-seated fluorescent proteins in vivo. Nat. Photonics 2009, 3, 379. [Google Scholar] [CrossRef]
  14. Zhang, H.F.; Maslov, K.; Stoica, G.; Wang, L. Functional photoacoustic microscopy for high-resolution and noninvasive in vivo imaging. Nat. Biotechnol. 2006, 24, 848–851. [Google Scholar] [CrossRef]
  15. Zhang, H.F.; Maslov, K.; Wang, L.V. In vivo imaging of subcutaneous structures using functional photoacoustic microscopy. Nat. Protoc. 2007, 2, 797–804. [Google Scholar] [CrossRef]
  16. De La Zerda, A.; Zavaleta, C.; Keren, S.; Vaithilingam, S.; Bodapati, S.; Liu, Z.; Levi, J.; Smith, B.; Ma, T.-J.; Oralkan, O.; et al. Carbon nanotubes as photoacoustic molecular imaging agents in living mice. Nat. Nanotechnol. 2008, 3, 557–562. [Google Scholar] [CrossRef]
  17. Bouchard, L.S.; Anwar, M.S.; Liu, G.L.; Hann, B.; Xie, Z.H.; Gray, J.W.; Wang, X.; Pines, A.; Chen, F.F. Picomolar sensitivity MRI and photoacoustic imaging of cobalt nanoparticles. Proc. Natl. Acad. Sci. USA 2009, 106, 4085–4089. [Google Scholar] [CrossRef] [Green Version]
  18. Cheng, R.; Shao, J.; Gao, X.; Tao, C.; Ge, J.; Liu, X. Noninvasive Assessment of Early Dental Lesion Using a Dual-Contrast Photoacoustic Tomography. Sci. Rep. 2016, 6, 21798. [Google Scholar] [CrossRef] [Green Version]
  19. Sun, Y.; Sobel, E.S.; Jiang, H. Quantitative three-dimensional photoacoustic tomography of the finger joints: An in vivo study. J. Biomed. Opt. 2009, 14, 064002. [Google Scholar] [CrossRef]
  20. Sun, Y.; Sobel, E.S.; Jiang, H. First assessment of three-dimensional quantitative photoacoustic tomography for in vivo detection of osteoarthritis in the finger joints. Med. Phys. 2011, 38, 4009–4017. [Google Scholar] [CrossRef] [Green Version]
  21. Zhang, Q.; Liu, Z.; Carney, P.R.; Yuan, Z.; Chen, H.; Roper, S.N.; Jiang, H. Non-invasive imaging of epileptic seizures in vivo using photoacoustic tomography. Phys. Med. Biol. 2008, 53, 1921–1931. [Google Scholar] [CrossRef] [PubMed]
  22. Kruger, R.A.; Reinecke, D.R.; Kruger, G.A. Thermoacoustic computed tomography-technical considerations. Med Phys. 1999, 26, 1832–1837. [Google Scholar] [CrossRef] [PubMed]
  23. Xu, M.; Wang, L.V. Universal back-projection algorithm for photoacoustic computed tomography. Phys. Rev. E 2005, 71, 016706. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  24. Yuan, J.; Xu, G.; Yu, Y.; Zhou, Y.; Carson, P.L.; Wang, X.; Liu, X. Real-time photoacoustic and ultrasound dual-modality imaging system facilitated with graphics processing unit and code parallel optimization. J. Biomed. Opt. 2013, 18, 086001. [Google Scholar] [CrossRef] [Green Version]
  25. Wang, K.; Huang, C.; Kao, Y.-J.; Chou, C.-Y.; Oraevsky, A.; Anastasio, M.A. Accelerating image reconstruction in three-dimensional optoacoustic tomography on graphics processing units. Med. Phys. 2013, 40, 023301. [Google Scholar] [CrossRef]
  26. Anastasio, M.A.; Zhang, J.; Pan, X.; Zou, Y.; Ku, G.; Wang, L.V. Half-time image reconstruction in thermoacoustic tomography. IEEE Trans. Med Imaging 2005, 24, 199–210. [Google Scholar] [CrossRef]
  27. Kunyansky, A.L. Explicit inversion formulae for the spherical mean Radon transform. Inverse Probl. 2007, 23, 373–383. [Google Scholar] [CrossRef] [Green Version]
  28. Zhang, J.; Anastasio, M.A.; La Rivière, P.J. Comparison of iterative reconstruction approaches for photoacoustic tomography. Proc. SPIE 2007, 6437, 256–261. [Google Scholar] [CrossRef]
  29. Li, C.; Wang, L.V. Photoacoustic tomography and sensing in biomedicine. Phys. Med. Biol. 2009, 54, R59–R97. [Google Scholar] [CrossRef]
  30. Ephrat, P.; Keenliside, L.; Seabrook, A.; Prato, F.S.; Carson, J.J.L. Three-dimensional photoacoustic imaging by sparse-array detection and iterative image reconstruction. J. Biomed. Opt. 2008, 13, 054052. [Google Scholar] [CrossRef] [Green Version]
  31. Paltauf, G.; Viator, J.A.; Prahl, S.A.; Jacques, S.L. Iterative reconstruction algorithm for optoacoustic imaging. J. Acoust. Soc. Am. 2002, 112, 1536–1544. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  32. Haltmeier, M.; Nguyen, L.V. Analysis of Iterative Methods in Photoacoustic Tomography with Variable Sound Speed. SIAM J. Imaging Sci. 2017, 10, 751–781. [Google Scholar] [CrossRef]
  33. Treeby, B.E.; Cox, B.T. k-Wave: MATLAB toolbox for the simulation and reconstruction of photoacoustic wave fields. J. Biomed. Opt. 2010, 15, 021314. [Google Scholar] [CrossRef] [PubMed]
  34. Steinberg, I.; Kim, J.; Schneider, M.K.; Hyun, D.; Zlitni, A.; Hopper, S.M.; Klap, T.; Sonn, G.A.; Dahl, J.J.; Kim, C.; et al. Superiorized Photo-Acoustic Non-NEgative Reconstruction (SPANNER) for Clinical Photoacoustic Imaging. IEEE Trans. Med. Imaging 2021, 40, 1888–1897. [Google Scholar] [CrossRef]
  35. Antholzer, S.; Haltmeier, M.; Schwab, J. Deep learning for photoacoustic tomography from sparse data. Inverse Probl. Sci. Eng. 2019, 27, 987–1005. [Google Scholar] [CrossRef] [Green Version]
  36. Poudel, J.; Lou, Y.; A Anastasio, M. A survey of computational frameworks for solving the acoustic inverse problem in three-dimensional photoacoustic computed tomography. Phys. Med. Biol. 2019, 64, 14TR01. [Google Scholar] [CrossRef]
  37. Jiang, H.; Yuan, Z.; Gu, X. Spatially varying optical and acoustic property reconstruction using finite-element-based photoacoustic tomography. J. Opt. Soc. Am. A 2006, 23, 878–888. [Google Scholar] [CrossRef]
  38. Yao, L.; Jiang, H. Finite-element-based photoacoustic tomography in time domain. J. Opt. A Pure Appl. Opt. 2009, 11, 085301. [Google Scholar] [CrossRef]
  39. Yao, L.; Jiang, H. Photoacoustic image reconstruction from few-detector and limited-angle data. Biomed. Opt. Express 2011, 2, 2649–2654. [Google Scholar] [CrossRef]
  40. Yuan, Z.; Jiang, H. Quantitative photoacoustic tomography: Recovery of optical absorption coefficient map of heterogeneous medium. Appl. Phys. Lett. 2006, 88, 231101. [Google Scholar] [CrossRef]
  41. Yuan, Z.; Zhang, Q.; Jiang, H. Simultaneous reconstruction of acoustic and optical properties of heterogeneous media by quantitative photoacoustic tomography. Opt. Express 2006, 14, 6749–6754. [Google Scholar] [CrossRef] [PubMed]
  42. Yuan, Z.; Jiang, H. Simultaneous recovery of tissue physiological and acoustic properties and the criteria for heterogeneous media by quantitative photoacoustic tomography. Opt. Lett. 2009, 34, 1714–1716. [Google Scholar] [CrossRef] [PubMed]
  43. Peng, K.; He, L.; Zhu, Z.; Tang, J.; Xiao, J. Three-dimensional photoacoustic tomography based on graphics-processing-unit-accelerated finite element method. Appl. Opt. 2013, 52, 8270–8279. [Google Scholar] [CrossRef]
  44. Shan, T.; Qi, J.; Jiang, M.; Jiang, H. GPU-based acceleration and mesh optimization of finite-element-method-based quantitative photoacoustic tomography: A step towards clinical applications. Appl. Opt. 2017, 56, 4426–4432. [Google Scholar] [CrossRef] [PubMed]
  45. Sun, Y.; Yuan, Z.; Jiang, H. Enhancing Mesh-based Photoacoustic Tomography with Parallel Computing on Multiprocessor Scheme. Commun. Comput. Phys. 2018, 24, 764–773. [Google Scholar] [CrossRef]
  46. Yu, L.; Zou, Y.; Sidky, E.Y.; Pelizzari, C.A.; Munro, P.; Pan, X. Region of interest reconstruction from truncated data in circular cone-beam CT. IEEE Trans. Med. Imaging 2006, 25, 869–881. [Google Scholar] [CrossRef]
  47. Chityala, R.N.; Hoffmann, K.R.; Bednarek, D.R.; Rudin, S. Region of interest (ROI) computed tomography. Med. Imaging Phys. Med. Imaging 2004, 5368, 534–541. [Google Scholar]
  48. Sun, Y.; Jiang, H. Quantitative three-dimensional photoacoustic tomography of the finger joints: Phantom studies in a spherical scanning geometry. Phys. Med. Biol. 2009, 54, 5457–5467. [Google Scholar] [CrossRef]
Figure 1. Reconstructed photoacoustic images of the simulated target by using full zone-based FE-PAT algorithm (a) and localized FE-PAT algorithm in a ROI with radii of 10.0 mm (b), 7.5 mm (c) and 4.0 mm (d).
Figure 1. Reconstructed photoacoustic images of the simulated target by using full zone-based FE-PAT algorithm (a) and localized FE-PAT algorithm in a ROI with radii of 10.0 mm (b), 7.5 mm (c) and 4.0 mm (d).
Photonics 09 00337 g001
Figure 2. The comparison of reconstructed photoacoustic intensities (absorbed optical energies) across the center of the simulated target by using a full zone-based FE-PAT algorithm and localized FE-PAT algorithm in various ROI sizes (R = 10.0 mm, R = 7.5 mm and R = 4.0 mm).
Figure 2. The comparison of reconstructed photoacoustic intensities (absorbed optical energies) across the center of the simulated target by using a full zone-based FE-PAT algorithm and localized FE-PAT algorithm in various ROI sizes (R = 10.0 mm, R = 7.5 mm and R = 4.0 mm).
Photonics 09 00337 g002
Figure 3. Acceleration performance and fitted curves of the localized FE-PAT algorithm in terms of the relative radius (a) and the relative number of nodes (b) in the ROI. The ordinates are based on a logarithmic scale.
Figure 3. Acceleration performance and fitted curves of the localized FE-PAT algorithm in terms of the relative radius (a) and the relative number of nodes (b) in the ROI. The ordinates are based on a logarithmic scale.
Photonics 09 00337 g003
Figure 4. Reconstructed photoacoustic images of the target in a phantom using the full zone-based FE-PAT algorithm (a) and localized FE-PAT algorithm in a ROI with a radius of 4.0 mm (b).
Figure 4. Reconstructed photoacoustic images of the target in a phantom using the full zone-based FE-PAT algorithm (a) and localized FE-PAT algorithm in a ROI with a radius of 4.0 mm (b).
Photonics 09 00337 g004
Figure 5. The comparison of reconstructed photoacoustic intensities (absorbed optical energies) across the center of the target in the phantom using the full zone-based FE-PAT algorithm and the localized FE-PAT algorithm in a ROI with a radius of 4.0 mm.
Figure 5. The comparison of reconstructed photoacoustic intensities (absorbed optical energies) across the center of the target in the phantom using the full zone-based FE-PAT algorithm and the localized FE-PAT algorithm in a ROI with a radius of 4.0 mm.
Photonics 09 00337 g005
Figure 6. Reconstructed photoacoustic images of the in vivo rat’s brain (a,b) and the qualitative profile across MCA (c) using the localized FE-PAT algorithm in a ROI with a radius of 6 mm. The frequency range used in the reconstruction was 41.5 Khz–2.075 Mhz (a) and 41.5 Khz–4.15 Mhz (b). MCA: middle cerebral artery; V: blood vessels. The transection line was marked as the blue line in (c).
Figure 6. Reconstructed photoacoustic images of the in vivo rat’s brain (a,b) and the qualitative profile across MCA (c) using the localized FE-PAT algorithm in a ROI with a radius of 6 mm. The frequency range used in the reconstruction was 41.5 Khz–2.075 Mhz (a) and 41.5 Khz–4.15 Mhz (b). MCA: middle cerebral artery; V: blood vessels. The transection line was marked as the blue line in (c).
Photonics 09 00337 g006
Table 1. Mesh parameters and performance comparison between conventional full zone-based FE-PAT algorithm and localized FE-PAT algorithm in various ROI sizes for simulation studies.
Table 1. Mesh parameters and performance comparison between conventional full zone-based FE-PAT algorithm and localized FE-PAT algorithm in various ROI sizes for simulation studies.
Radius of ROI (mm)No. of ROI NodesTime Cost (s)Target Size (mm)Error of Target SizeReconstruction Error
*448977453.96--
10.0200911444.083.03%2.84%
7.511193174.134.29%3.67%
4.0323494.144.55%4.47%
* Full zone-based FE-PAT algorithm. The radius of the imaging zone is 15 mm.
Table 2. Mesh parameters and performance comparison between the conventional full zone-based FE-PAT algorithm and the localized FE-PAT algorithm in various ROI sizes for phantom experimental studies.
Table 2. Mesh parameters and performance comparison between the conventional full zone-based FE-PAT algorithm and the localized FE-PAT algorithm in various ROI sizes for phantom experimental studies.
Algorithm *No. of ROI NodesTime Cost (s)Target Size (mm)Error of Target SizeReconstruction Error
Full-zone448977303.23--
Localized328493.220.31%3.72%
* The radius of the full zone-based FE-PAT algorithm is 15 mm. The radius of the ROI used in the localized FE-PAT algorithm is 4.0 mm.
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Sun, Y.; Jiang, H. Enhancing Finite Element-Based Photoacoustic Tomography by Localized Reconstruction Method. Photonics 2022, 9, 337. https://doi.org/10.3390/photonics9050337

AMA Style

Sun Y, Jiang H. Enhancing Finite Element-Based Photoacoustic Tomography by Localized Reconstruction Method. Photonics. 2022; 9(5):337. https://doi.org/10.3390/photonics9050337

Chicago/Turabian Style

Sun, Yao, and Huabei Jiang. 2022. "Enhancing Finite Element-Based Photoacoustic Tomography by Localized Reconstruction Method" Photonics 9, no. 5: 337. https://doi.org/10.3390/photonics9050337

APA Style

Sun, Y., & Jiang, H. (2022). Enhancing Finite Element-Based Photoacoustic Tomography by Localized Reconstruction Method. Photonics, 9(5), 337. https://doi.org/10.3390/photonics9050337

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop