## Singular value decomposition analysis of a photoacoustic imaging system and 3D imaging at 0.7 FPS |

Optics Express, Vol. 19, Issue 14, pp. 13405-13417 (2011)

http://dx.doi.org/10.1364/OE.19.013405

Acrobat PDF (1210 KB)

### Abstract

Photoacoustic imaging is a non-ionizing imaging modality that provides contrast consistent with optical imaging techniques while the resolution and penetration depth is similar to ultrasound techniques. In a previous publication [Opt. Express **18**, 11406 (2010)], a technique was introduced to experimentally acquire the imaging operator for a photoacoustic imaging system. While this was an important foundation for future work, we have recently improved the experimental procedure allowing for a more densely populated imaging operator to be acquired. Subsets of the imaging operator were produced by varying the transducer count as well as the measurement space temporal sampling rate. Examination of the matrix rank and the effect of contributing object space singular vectors to image reconstruction were performed. For a PAI system collecting only limited data projections, matrix rank increased linearly with transducer count and measurement space temporal sampling rate. Image reconstruction using a regularized pseudoinverse of the imaging operator was performed on photoacoustic signals from a point source, line source, and an array of point sources derived from the imaging operator. As expected, image quality increased for each object with increasing transducer count and measurement space temporal sampling rate. Using the same approach, but on experimentally sampled photoacoustic signals from a moving point-like source, acquisition, data transfer, reconstruction and image display took 1.4 s using one laser pulse per 3D frame. With relatively simple hardware improvements to data transfer and computation speed, our current imaging results imply that acquisition and display of 3D photoacoustic images at laser repetition rates of 10Hz is easily achieved.

© 2011 OSA

## 1. Introduction

### 1.1 Background

1. L. V. Wang, “Ultrasound-mediated biophotonic imaging: a review of acousto-optical tomography and photo-acoustic tomography,” Dis. Markers **19**(2-3), 123–138 (2003-2004). [PubMed]

2. G. J. Diebold, T. Sun, and M. I. Khan, “Photoacoustic monopole radiation in one, two, and three dimensions,” Phys. Rev. Lett. **67**(24), 3384–3387 (1991). [CrossRef] [PubMed]

3. R. A. Kruger, D. R. Reinecke, and G. A. Kruger, “Thermoacoustic computed tomography—technical considerations,” Med. Phys. **26**(9), 1832–1837 (1999). [CrossRef] [PubMed]

4. D. Frauchiger, K. P. Kostli, G. Paltauf, M. Frenz, and H. P. Weber, “Optoacoustic tomography using a two dimensional optical pressure transducer and two different reconstruction algorithms,” Proc. SPIE **4434**, 74–80 (2001). [CrossRef]

9. M. H. Xu and L. V. Wang, “Universal back-projection algorithm for photoacoustic computed tomography,” Phys. Rev. E Stat. Nonlin. Soft Matter Phys. **71**(1), 016706 (2005). [CrossRef] [PubMed]

### 1.2 Singular value decomposition

**g**is a vector that represents the measured data set,

**H**is the imaging operator, and

**f**is a vector that represents a finite-dimensional approximation of the unknown object(s) that produced the data in

**g**. It is generally found that

**H**is singular and cannot be inverted directly. For singular matrices, it can be shown that an

*M*x

*N*matrix,

**H**, can be decomposed by means of Eq. (2):where

**U**is an

*M*x

*M*matrix,

**V**is an

*N*x

*N*matrix, and both are nonsingular. The

*M*x

*N*matrix

**S**is a diagonal matrix with non-zero diagonal entries representing the singular values of the imaging operator. The rows of

**U**and columns of

**V**

^{T}are the orthonormal singular vectors that form a complete linear vector space (in their respective dimensions) describing the measurement space and discretized object space, respectively. It follows that each singular value of

**S**relates the sensitivity of the imaging operator to the corresponding singular vectors in

**U**and

**V**

^{T}.

### 1.3 Estimate of effective singular values

**V**are linearly independent. However, by examining the associated magnitude of the singular values in matrix

^{T}**S**, it is clear not all vectors contribute equally to the overall system response. In fact, some do not effectively contribute at all to the reconstruction of an object [14

14. D. W. Wilson and H. H. Barrett, “Decomposition of images and objects into measurement and null components,” Opt. Express **2**(6), 254–260 (1998). [CrossRef] [PubMed]

15. K. Konstantinides and K. Yao, “Statistical analysis of effective singular values in matrix rank determination,” IEEE Trans. Acoust. Speech Signal Process. **36**(5), 757–763 (1988). [CrossRef]

17. D. J. Kadrmas, E. C. Frey, and B. M. W. Tsui, “An SVD investigation of modeling scatter in multiple energy windows for improved SPECT images,” IEEE Trans. Nucl. Sci. **43**(4), 2275–2284 (1996). [CrossRef] [PubMed]

15. K. Konstantinides and K. Yao, “Statistical analysis of effective singular values in matrix rank determination,” IEEE Trans. Acoust. Speech Signal Process. **36**(5), 757–763 (1988). [CrossRef]

*t*is the effective rank of the matrix. While this threshold is seemingly subjective, it can be evaluated more precisely by stating the rank of the matrix is estimated to be equal to

*t*if,Via the results of [15

15. K. Konstantinides and K. Yao, “Statistical analysis of effective singular values in matrix rank determination,” IEEE Trans. Acoust. Speech Signal Process. **36**(5), 757–763 (1988). [CrossRef]

17. D. J. Kadrmas, E. C. Frey, and B. M. W. Tsui, “An SVD investigation of modeling scatter in multiple energy windows for improved SPECT images,” IEEE Trans. Nucl. Sci. **43**(4), 2275–2284 (1996). [CrossRef] [PubMed]

**S**become exceedingly small when the remaining columns of V

^{T}do not contribute useful information about object space.

**f**in Eq. (1). In this paper, the SVD of the imaging operator,

**H**, was computed and utilized to compute a regularized pseudoinverse solution to Eq. (1). The regularized solution utilized only the singular values and associated singular vectors that were identified as corresponding to useful information as described above.

### 1.4 Experimental objective

10. M. Roumeliotis, R. Z. Stodilka, M. A. Anastasio, G. Chaudhary, H. Al-Aabed, E. Ng, A. Immucci, and J. J. L. Carson, “Analysis of a photoacoustic imaging system by the crosstalk matrix and singular value decomposition,” Opt. Express **18**(11), 11406–11417 (2010). [CrossRef] [PubMed]

13. M. Roumeliotis, P. Ephrat, J. Patrick, and J. J. L. Carson, “Development and characterization of an omnidirectional photoacoustic point source for calibration of a staring 3D photoacoustic imaging system,” Opt. Express **17**(17), 15228–15238 (2009). [CrossRef] [PubMed]

## 2. Methods

### 2.1 Photoacoustic imaging system

10. M. Roumeliotis, R. Z. Stodilka, M. A. Anastasio, G. Chaudhary, H. Al-Aabed, E. Ng, A. Immucci, and J. J. L. Carson, “Analysis of a photoacoustic imaging system by the crosstalk matrix and singular value decomposition,” Opt. Express **18**(11), 11406–11417 (2010). [CrossRef] [PubMed]

^{TM}v8.5. The imaging system utilized a total of 30 ultrasound transducers (model V304, 1” diameter, 2.25 MHz with fractional bandwidth of 65%,

*Panametrics-NDT*, Waltham, Massachusetts) in a staring hemispherical arrangement. Transducers were mounted on custom frames, each supporting 3 transducers at evenly spaced azimuthal angles (36° between columns). Fifteen (15) of the transducers were mounted on frames with zenith angles of 22.5°, 45°, and 67.5°. The remaining 15 transducers were mounted on frames placed azimuthally between the original frames, at zenith angles of 33.75°, 56.25°, and 78.75°. A schematic of the transducer array is shown in Fig. 1(a) .

### 2.2 System calibration scan

10. M. Roumeliotis, R. Z. Stodilka, M. A. Anastasio, G. Chaudhary, H. Al-Aabed, E. Ng, A. Immucci, and J. J. L. Carson, “Analysis of a photoacoustic imaging system by the crosstalk matrix and singular value decomposition,” Opt. Express **18**(11), 11406–11417 (2010). [CrossRef] [PubMed]

^{3}object space and 1.5 mm step-size for a total of 8,000 voxel locations. It was originally intended to complete scans of varying voxel count and step-size. However, calibration scans exceeding 8,000 voxels represented the practical limit of our experimental methodology. For the experimental image reconstruction, a second scan of dimension 16x16x7 mm

^{3}was produced with 30 transducers at a measurement space temporal sampling rate of 10 MHz.

### 2.3 Singular value decomposition

### 2.4 Regularized pseudoinverse image reconstruction

**f**in Eq. (1). The SVD component matrices of Eq. (2) were used to represent the imaging operator in solving for the vector

**f**. The built-in MATLAB function (pinv) was used to compute the pseudoinverse of the matrix

**S**, referred to as

**S**

^{+}. Due to noise present in the imaging operator, singular values with indices greater than 90 percent of the computed matrix rank were set to zero when computing

**S**. The pseudoinverse was then multiplied by the data set measured experimentally in

^{+}**g**, yielding an estimate of the object

**f**.

**g**, and the right side of Eq. (1), which was then added to a master of the image. While this model is computationally expensive, it was deemed instructive to show the results of the simulated data sets being reconstructed via another approach.

### 2.5 Real-time photoacoustic imaging

## 3. Results

### 3.1 Estimate of matrix rank

### 3.2 Image reconstruction of objects with different transducer count

**g**, to produce unique data. These time-dependent pressure measurements were then multiplied by the pseudoinverse obtained from component matrices (

**U**

^{T},

**S**

^{+},

**V**) for each transducer count. Figure 3 shows the centre z-plane of the image reconstruction of each of the sources as a function of transducer count. The second column from the right shows the results of the iterative reconstruction technique applied to each of the 3 object types. A phantom of the objects is displayed in the right-most column.

### 3.3 Image reconstruction of objects with different measurement space temporal sampling rate

### 3.4 Image reconstruction of a point source in real-time (1.4 seconds per frame)

## 4. Discussion

### 4.1 Estimate of matrix rank

19. A. A. Oraevsky, V. G. Andreev, A. A. Karabutov, and R. O. Esenaliev, “Two-dimensional opto-acoustic tomography transducer array and image reconstruction algorithm,” Proc. SPIE **3601**, 256–267 (1999). [CrossRef]

### 4.2 Image reconstruction of objects with different transducer count

### 4.3 Image reconstruction of objects with different measurement space temporal sampling rate

### 4.4 Real-time imaging of a photoacoustic point-like source

21. S. Ashkenazi, “Photoacoustic lifetime imaging of dissolved oxygen using methylene blue,” J. Biomed. Opt. **15**(4), 040501 (2010). [CrossRef] [PubMed]

22. J. Su, A. Karpiouk, B. Wang, and S. Emelianov, “Photoacoustic imaging of clinical metal needles in tissue,” J. Biomed. Opt. **15**(2), 021309 (2010). [CrossRef] [PubMed]

### 4.5 Implications to previous work

12. P. Ephrat, M. Roumeliotis, F. S. Prato, and J. J. Carson, “Four-dimensional photoacoustic imaging of moving targets,” Opt. Express **16**(26), 21570–21581 (2008). [CrossRef] [PubMed]

### 4.6 Technical considerations for faster 3D frame rates

## 5. Conclusion

^{3}. By decomposing the imaging operator and analyzing the magnitude of the singular values, the effective matrix rank of the imaging operator was estimated. The results clearly show a linear increase in the effective matrix rank as both the transducer count and measurement space temporal sampling rate were increased. Image reconstruction of a variety of objects was performed by computing the pseudoinverse of each decomposed imaging operator and multiplying it with the measured data set of each phantom object. As expected, image fidelity was greater in cases with more transducers and increased temporal sampling rate. Image reconstruction was also performed in a period of 1.4 s per frame on a moving point-like source to demonstrate real-time acquisition, reconstruction and display of 3D photoacoustic images. Extension of the real-time 3D photoacoustic imaging to a frame rate of 10Hz should be easily achieved with straightforward improvements to data transfer time and computation speed.

## Acknowledgments

## References and links

1. | L. V. Wang, “Ultrasound-mediated biophotonic imaging: a review of acousto-optical tomography and photo-acoustic tomography,” Dis. Markers |

2. | G. J. Diebold, T. Sun, and M. I. Khan, “Photoacoustic monopole radiation in one, two, and three dimensions,” Phys. Rev. Lett. |

3. | R. A. Kruger, D. R. Reinecke, and G. A. Kruger, “Thermoacoustic computed tomography—technical considerations,” Med. Phys. |

4. | D. Frauchiger, K. P. Kostli, G. Paltauf, M. Frenz, and H. P. Weber, “Optoacoustic tomography using a two dimensional optical pressure transducer and two different reconstruction algorithms,” Proc. SPIE |

5. | M. Frenz, K. P. Kostli, G. Paltauf, H. Schmidt-Kloiber, and H. P. Weber, “Reconstruction technique for optoacoustic imaging,” in |

6. | C. G. A. Hoelen and F. F. M. de Mul, “Image reconstruction for photoacoustic scanning of tissue structures,” Appl. Opt. |

7. | K. P. Kostli, D. Frauchiger, J. J. Niederhauser, G. Paltauf, H. P. Weber, and M. Frenz, “Optoacoustic imaging using a three-dimensional reconstruction algorithm,” IEEE J. Sel. Top. Quantum Electron. |

8. | G. Paltauf, J. A. Viator, S. A. Prahl, and S. L. Jacques, “Iterative reconstruction algorithm for optoacoustic imaging,” J. Acoust. Soc. Am. |

9. | M. H. Xu and L. V. Wang, “Universal back-projection algorithm for photoacoustic computed tomography,” Phys. Rev. E Stat. Nonlin. Soft Matter Phys. |

10. | M. Roumeliotis, R. Z. Stodilka, M. A. Anastasio, G. Chaudhary, H. Al-Aabed, E. Ng, A. Immucci, and J. J. L. Carson, “Analysis of a photoacoustic imaging system by the crosstalk matrix and singular value decomposition,” Opt. Express |

11. | P. Ephrat, L. Keenliside, A. Seabrook, F. S. Prato, and J. J. L. Carson, “Three-dimensional photoacoustic imaging by sparse-array detection and iterative image reconstruction,” J. Biomed. Opt. |

12. | P. Ephrat, M. Roumeliotis, F. S. Prato, and J. J. Carson, “Four-dimensional photoacoustic imaging of moving targets,” Opt. Express |

13. | M. Roumeliotis, P. Ephrat, J. Patrick, and J. J. L. Carson, “Development and characterization of an omnidirectional photoacoustic point source for calibration of a staring 3D photoacoustic imaging system,” Opt. Express |

14. | D. W. Wilson and H. H. Barrett, “Decomposition of images and objects into measurement and null components,” Opt. Express |

15. | K. Konstantinides and K. Yao, “Statistical analysis of effective singular values in matrix rank determination,” IEEE Trans. Acoust. Speech Signal Process. |

16. | K. Konstantinides, “Threshold bounds in SVD and a new iterative algorithm for order selection in Ar models,” IEEE Trans. Signal Process. |

17. | D. J. Kadrmas, E. C. Frey, and B. M. W. Tsui, “An SVD investigation of modeling scatter in multiple energy windows for improved SPECT images,” IEEE Trans. Nucl. Sci. |

18. | P. Ephrat and J. J. L. Carson, “Measurement of photoacoustic detector sensitivity distribution by robotic source placement,” in |

19. | A. A. Oraevsky, V. G. Andreev, A. A. Karabutov, and R. O. Esenaliev, “Two-dimensional opto-acoustic tomography transducer array and image reconstruction algorithm,” Proc. SPIE |

20. | J. Provost and F. Lesage, “The application of compressed sensing for photo-acoustic tomography,” IEEE Trans. Med. Imaging |

21. | S. Ashkenazi, “Photoacoustic lifetime imaging of dissolved oxygen using methylene blue,” J. Biomed. Opt. |

22. | J. Su, A. Karpiouk, B. Wang, and S. Emelianov, “Photoacoustic imaging of clinical metal needles in tissue,” J. Biomed. Opt. |

**OCIS Codes**

(170.0110) Medical optics and biotechnology : Imaging systems

(170.5120) Medical optics and biotechnology : Photoacoustic imaging

**ToC Category:**

Imaging Systems

**History**

Original Manuscript: April 11, 2011

Revised Manuscript: June 3, 2011

Manuscript Accepted: June 4, 2011

Published: June 27, 2011

**Virtual Issues**

Vol. 6, Iss. 8 *Virtual Journal for Biomedical Optics*

**Citation**

Michael B. Roumeliotis, Robert Z. Stodilka, Mark. A. Anastasio, Eldon Ng, and Jeffrey J. L. Carson, "Singular value decomposition analysis of a photoacoustic imaging system and 3D imaging at 0.7 FPS," Opt. Express **19**, 13405-13417 (2011)

http://www.opticsinfobase.org/vjbo/abstract.cfm?URI=oe-19-14-13405

Sort: Year | Journal | Reset

### References

- L. V. Wang, “Ultrasound-mediated biophotonic imaging: a review of acousto-optical tomography and photo-acoustic tomography,” Dis. Markers 19(2-3), 123–138 (2003-2004). [PubMed]
- G. J. Diebold, T. Sun, and M. I. Khan, “Photoacoustic monopole radiation in one, two, and three dimensions,” Phys. Rev. Lett. 67(24), 3384–3387 (1991). [CrossRef] [PubMed]
- R. A. Kruger, D. R. Reinecke, and G. A. Kruger, “Thermoacoustic computed tomography—technical considerations,” Med. Phys. 26(9), 1832–1837 (1999). [CrossRef] [PubMed]
- D. Frauchiger, K. P. Kostli, G. Paltauf, M. Frenz, and H. P. Weber, “Optoacoustic tomography using a two dimensional optical pressure transducer and two different reconstruction algorithms,” Proc. SPIE 4434, 74–80 (2001). [CrossRef]
- M. Frenz, K. P. Kostli, G. Paltauf, H. Schmidt-Kloiber, and H. P. Weber, “Reconstruction technique for optoacoustic imaging,” in Biomedical Optoacoustics II, January 23–24, 2001 (SPIE), 130–137.
- C. G. A. Hoelen and F. F. M. de Mul, “Image reconstruction for photoacoustic scanning of tissue structures,” Appl. Opt. 39(31), 5872–5883 (2000). [CrossRef] [PubMed]
- K. P. Kostli, D. Frauchiger, J. J. Niederhauser, G. Paltauf, H. P. Weber, and M. Frenz, “Optoacoustic imaging using a three-dimensional reconstruction algorithm,” IEEE J. Sel. Top. Quantum Electron. 7(6), 918–923 (2001). [CrossRef]
- G. Paltauf, J. A. Viator, S. A. Prahl, and S. L. Jacques, “Iterative reconstruction algorithm for optoacoustic imaging,” J. Acoust. Soc. Am. 112(4), 1536–1544 (2002). [CrossRef] [PubMed]
- M. H. Xu and L. V. Wang, “Universal back-projection algorithm for photoacoustic computed tomography,” Phys. Rev. E Stat. Nonlin. Soft Matter Phys. 71(1), 016706 (2005). [CrossRef] [PubMed]
- M. Roumeliotis, R. Z. Stodilka, M. A. Anastasio, G. Chaudhary, H. Al-Aabed, E. Ng, A. Immucci, and J. J. L. Carson, “Analysis of a photoacoustic imaging system by the crosstalk matrix and singular value decomposition,” Opt. Express 18(11), 11406–11417 (2010). [CrossRef] [PubMed]
- P. Ephrat, L. Keenliside, A. Seabrook, F. S. Prato, and J. J. L. Carson, “Three-dimensional photoacoustic imaging by sparse-array detection and iterative image reconstruction,” J. Biomed. Opt. 13(5), 054052 (2008). [CrossRef] [PubMed]
- P. Ephrat, M. Roumeliotis, F. S. Prato, and J. J. Carson, “Four-dimensional photoacoustic imaging of moving targets,” Opt. Express 16(26), 21570–21581 (2008). [CrossRef] [PubMed]
- M. Roumeliotis, P. Ephrat, J. Patrick, and J. J. L. Carson, “Development and characterization of an omnidirectional photoacoustic point source for calibration of a staring 3D photoacoustic imaging system,” Opt. Express 17(17), 15228–15238 (2009). [CrossRef] [PubMed]
- D. W. Wilson and H. H. Barrett, “Decomposition of images and objects into measurement and null components,” Opt. Express 2(6), 254–260 (1998). [CrossRef] [PubMed]
- K. Konstantinides and K. Yao, “Statistical analysis of effective singular values in matrix rank determination,” IEEE Trans. Acoust. Speech Signal Process. 36(5), 757–763 (1988). [CrossRef]
- K. Konstantinides, “Threshold bounds in SVD and a new iterative algorithm for order selection in Ar models,” IEEE Trans. Signal Process. 39(5), 1218–1221 (1991). [CrossRef]
- D. J. Kadrmas, E. C. Frey, and B. M. W. Tsui, “An SVD investigation of modeling scatter in multiple energy windows for improved SPECT images,” IEEE Trans. Nucl. Sci. 43(4), 2275–2284 (1996). [CrossRef] [PubMed]
- P. Ephrat and J. J. L. Carson, “Measurement of photoacoustic detector sensitivity distribution by robotic source placement,” in 9th Conference on Photons Plus Ultrasound: Imaging and Sensing 2008, January 20–23, 2008 (SPIE).
- A. A. Oraevsky, V. G. Andreev, A. A. Karabutov, and R. O. Esenaliev, “Two-dimensional opto-acoustic tomography transducer array and image reconstruction algorithm,” Proc. SPIE 3601, 256–267 (1999). [CrossRef]
- J. Provost and F. Lesage, “The application of compressed sensing for photo-acoustic tomography,” IEEE Trans. Med. Imaging 28(4), 585–594 (2009). [CrossRef] [PubMed]
- S. Ashkenazi, “Photoacoustic lifetime imaging of dissolved oxygen using methylene blue,” J. Biomed. Opt. 15(4), 040501 (2010). [CrossRef] [PubMed]
- J. Su, A. Karpiouk, B. Wang, and S. Emelianov, “Photoacoustic imaging of clinical metal needles in tissue,” J. Biomed. Opt. 15(2), 021309 (2010). [CrossRef] [PubMed]

## Cited By |
Alert me when this paper is cited |

OSA is able to provide readers links to articles that cite this paper by participating in CrossRef's Cited-By Linking service. CrossRef includes content from more than 3000 publishers and societies. In addition to listing OSA journal articles that cite this paper, citing articles from other participating publishers will also be listed.

« Previous Article | Next Article »

OSA is a member of CrossRef.