_{2}) is discussed. The derivations for the nonbandlimited and bandlimited OA signals from many red blood cells (RBCs) are presented. The OA field generated by many RBCs was obtained by summing the OA field emitted by each RBC approximated as a fluid sphere. A Monte Carlo technique was employed generating the spatial organizations of RBCs in two-dimensional. The RBCs were assumed to have the same SO

_{2}level in a simulated configuration. The fractional number of oxyhemoglobin molecules, confined in a cell, determined the cellular SO

_{2}and also defined the blood SO

_{2}. For the nonbandlimited case, the OA signal amplitude decreased and increased linearly with blood SO

_{2}when illuminated by 700 and 1000 nm radiations, respectively. The power spectra exhibited similar trends over the entire frequency range (MHz to GHz). For the bandlimited case, three acoustic receivers with 2, 10, and 50 MHz as the center frequencies were considered. The linear variations of the OA amplitude with blood SO

_{2}were also observed for each receiver at those laser sources. The good agreement between simulated and published experimental results validates the model qualitatively.

## 1.

## Introduction

Blood oxygen saturation (${\mathrm{SO}}_{2}$) monitoring is a useful clinical practice for managing a number of clinical situations. The blood oxygenation level (derived from hemoglobin concentration) is crucially important for understanding brain hemodynamics in response to sensory stimulation.^{1} Moreover, it can provide useful clinical information for numerous medical applications such as evaluating the effects of chemotherapy and radiotherapy on tumors,^{2} monitoring the healing of wounds,^{3} and studying gene expression.^{4} Of all the techniques that have been employed to measure blood oxygenation *in vivo*, pulse oxymetry remains the most commonly used clinical technique. Nevertheless, the poor spatial resolution^{5} and the inability to distinguish between arterial and venous blood due to single point measurements are limitations of this technique. Other imaging modalities such as magnetic resonance imaging,^{6} positron emission tomography,^{7} and single photon emission computed tomography^{8} have also been used to quantify blood oxygenation levels. Each of the aforementioned techniques is resource intensive making the quick assessment of oxygenation levels a complex task.^{9}

Optoacoustics (OAs) is a noninvasive, hybrid imaging modality, which probes the optical and thermoelastic properties of tissue by detecting the pressure waves produced by laser irradiation. The emitted acoustic waves result from the short-pulsed laser-induced transient thermoelastic expansion that occurs when the incident optical energy is absorbed and transforms into heat. Advantages of OA imaging lie in the fact that the generated sound waves scatter two to three orders of magnitude less than optical waves and thus can travel a larger distance due to lower attenuation.^{10} The contrast in OA imaging is directly related to optical absorption, hence optical spectral information of tissue can be obtained. In the recent years, OA has found applications in visualizing carcinomas using gold nanoparticles,^{11} detecting circulating melanoma cells,^{12} single-cell imaging,^{13} and vascular imaging.^{14} In addition, Wang and colleagues have used multi-wavelength OA systems for *in vivo* structural and functional imaging of oxygen saturation in rat brain^{15} as well as imaging hemoglobin concentration and oxygenation.^{16} Esenaliev et al.^{17} studied the feasibility of the OA technique for noninvasive monitoring of blood oxygenation. Further, triple-wavelength OA systems have been used to measure oxygenation from the external jugular vein of *in vivo* sheep.^{18} Recently, studies have also been carried out to generate quantitative parameters by performing frequency domain analysis of the ultrasound radio frequency (RF) signals for tissue characterization.^{19}20.^{–}^{21}

Motivated from these studies, our group recently developed a frequency domain theoretical framework to study the OA signals from a collection of spherical absorbers.^{22}^{,}^{23} The OA field, generated by many absorbers, was obtained by using the linear superposition principle for the OA fields emitted by the individual absorbers. This framework was utilized to study OA signal properties of blood samples with nonaggregating and aggregating red blood cells (RBCs).^{22} It was shown that the OA signal amplitude increased monotonically with hematocrit for the nonaggregating blood and that was similar to published experimental data.^{24} The OA signal amplitude also exhibited a monotonic rise as the RBC cluster size increased for the aggregating blood. It was estimated that the spectral power at 15.6 MHz increased about 11 dB at an aggregation level as compared to that of the nonaggregating blood at the same hematocrit. This prediction was later confirmed in an experimental study.^{21} The effects of RBC oxygenation on OA signals were also examined using this theoretical model.^{23} The OA signals were simulated from mixtures of oxygenated and deoxygenated RBCs. A 6 fold decrease and 5 fold increase of the OA signal amplitude were estimated when blood ${\mathrm{SO}}_{2}$ varied from 0% to 100% for the 700 and 1000 nm illuminating radiations, respectively. Although, slightly nonlinear variations were observed, the overall trends of the simulated results were in accordance with those of published experimental results.^{17}^{,}^{18}

The purpose of the paper is to present a modeling approach to examine the dependence of the OA signals on blood oxygenation. This approach is suitable to model the situations when the individual cells have the same oxygenation level in a blood sample and thus, differs from the previous study reported in Ref. 23. The effect of the finite nature of receiver bandwidth on the OA signals was also investigated in this simulation work. In this study, each RBC was approximated as a fluid sphere. The fractional number of oxyhemoglobin molecules, packaged in a cell, determined the ${\mathrm{SO}}_{2}$ level of that cell. It also defined the blood ${\mathrm{SO}}_{2}$ assuming that ${\mathrm{SO}}_{2}$ states were identical for all constituting cells. The optical absorption coefficient of a RBC is linearly proportional to its ${\mathrm{SO}}_{2}$ level and dictates the OA emission amplitude. A Monte Carlo technique was employed to generate spatially random distributions of RBCs in two-dimensional (2-D). The OA signals were computed from such configurations with different ${\mathrm{SO}}_{2}$ levels. Both the time and spectral domain features were extracted from those signals as a function of blood ${\mathrm{SO}}_{2}$. Our simulated results exhibited good agreement with those of published experimental results.

This paper is organized in the following way: in Sec. 2, we derive the expressions for the nonbandlimited and bandlimited OA signals from a collection of RBCs. The numerical method implemented to compute the OA signals is illustrated in Sec. 3. The simulated signals along with their time and spectral domain features are displayed in Sec. 4. In Sec. 5, the results are discussed and compared to published experimental results. Finally, conclusions of the study have been drawn in Sec. 6.

## 2.

## Theoretical Model

## 2.1.

### Derivation of the Nonbandlimited Signal

The detailed derivation of the theoretical model was presented in Refs. 22 and 23. However, for the sake of completeness, the theoretical framework is discussed here in brief. The analytical expression for the OA field emitted by a spherical absorber at a large distance $r$ from its center due to the uniform illumination by optical radiation is given by,^{25}

## (1)

$${p}_{f}^{\text{single}}(r,{k}_{f})\phantom{\rule{0ex}{0ex}}=\frac{i\mu \beta \text{\hspace{0.17em}}{I}_{0}{\upsilon}_{s}{a}^{2}}{{C}_{P}r}\frac{[\mathrm{sin}\text{\hspace{0.17em}}\widehat{q}-\widehat{q}\text{\hspace{0.17em}}\mathrm{cos}\text{\hspace{0.17em}}\widehat{q}]{e}^{i{k}_{f}(r-a)}}{{\widehat{q}}^{2}[(1-\widehat{\rho})(\mathrm{sin}\text{\hspace{0.17em}}\widehat{q}/\widehat{q})-\mathrm{cos}\text{\hspace{0.17em}}\widehat{q}+i\widehat{\rho}\widehat{\upsilon}\text{\hspace{0.17em}}\mathrm{sin}\text{\hspace{0.17em}}\widehat{q}]}.$$^{25}During the derivation of the wave equation the condition of thermal confinement was imposed (the pressure pulse was launched before the heat conduction took place). Moreover, the solution was obtained by solving the time independent wave equation inside and outside the absorber in spherical polar coordinates and by implementing the appropriate boundary conditions (i.e., the continuity of pressure field and particle velocity at the spherical boundary).

^{22}

^{,}

^{23}

^{,}

^{25}

The pressure field at a large distance from the center of the uniformly illuminated region composed of an ensemble of absorbers with identical physical properties and equal radii becomes,^{22}^{,}^{23}

## (2)

$${p}_{f}^{\text{ensemble}}(r,{k}_{f})\approx \frac{i\mu \beta {I}_{0}{\upsilon}_{s}{a}^{2}}{{C}_{P}r}\phantom{\rule{0ex}{0ex}}\frac{[\mathrm{sin}\text{\hspace{0.17em}}\widehat{q}-\widehat{q}\text{\hspace{0.17em}}\mathrm{cos}\text{\hspace{0.17em}}\widehat{q}]{e}^{i{k}_{f}(r-a)}}{{\widehat{q}}^{2}[(1-\widehat{\rho})(\mathrm{sin}\text{\hspace{0.17em}}\widehat{q}/\widehat{q})-\mathrm{cos}\text{\hspace{0.17em}}\widehat{q}+i\widehat{\rho}\widehat{\upsilon}\text{\hspace{0.17em}}\mathrm{sin}\text{\hspace{0.17em}}\widehat{q}]}\sum _{n=1}^{N}{e}^{-i{\mathbf{k}}_{f}\xb7{\mathbf{r}}_{n}}.$$^{26}

^{,}

^{27}Equation (2) clearly shows that the OA signal amplitude is linearly proportional to the cellular absorption coefficient. The time dependent pressure field for the uniform illumination of the absorbers by a delta function heating pulse can be obtained by employing the Fourier transformation of Eq. (2) as,

## (3)

$${p}_{f}^{\text{ensemble}}(r,t)\approx \frac{i\mu \beta F{\upsilon}_{s}{a}^{2}}{2\pi {C}_{P}r}\phantom{\rule{0ex}{0ex}}{\int}_{-\infty}^{\infty}d\omega \frac{[\mathrm{sin}\text{\hspace{0.17em}}\widehat{q}-\widehat{q}\text{\hspace{0.17em}}\mathrm{cos}\text{\hspace{0.17em}}\widehat{q}]{e}^{i{k}_{f}(r-a-{v}_{f}t)}}{{\widehat{q}}^{2}[(1-\widehat{\rho})(\mathrm{sin}\text{\hspace{0.17em}}\widehat{q}/\widehat{q})-\mathrm{cos}\text{\hspace{0.17em}}\widehat{q}+i\widehat{q}\widehat{\upsilon}\text{\hspace{0.17em}}\mathrm{sin}\text{\hspace{0.17em}}\widehat{q}]}\phantom{\rule{0ex}{0ex}}\sum _{n=1}^{N}{e}^{-i{\mathbf{k}}_{f}\xb7{\mathbf{r}}_{n}},$$^{28}In case of an analytic signal, the imaginary part of the signal is the Hilbert transform of the real part. Therefore, at any time point the real part of the complex valued function provides the instantaneous RF signal and the corresponding signal amplitude can be obtained from its magnitude.

## 2.2.

### Derivation of the Bandlimited Signal

In practice, OA signals are detected using a ultrasonic transducer with finite receiving bandwidth that consequently filters out many frequency components of the signal. The finite nature of the receiver bandwidth can easily be incorporated in the model by introducing a Gaussian function, which can be used to model the frequency response profile of a transducer,^{29} in the frequency domain. Therefore, the time dependent pressure field can be written as,

## (4)

$${p}_{f}^{\text{ensemble}}(r,t)\approx \frac{i\mu \beta F{\upsilon}_{s}{a}^{2}}{2\pi {C}_{P}r}{\int}_{-\infty}^{\infty}d\omega \text{\hspace{0.17em}}\mathrm{exp}[-\frac{{(\omega -{\omega}_{0})}^{2}}{2{\sigma}^{2}}]\phantom{\rule{0ex}{0ex}}\frac{[\mathrm{sin}\text{\hspace{0.17em}}\widehat{q}-\widehat{q}\text{\hspace{0.17em}}\mathrm{cos}\text{\hspace{0.17em}}\widehat{q}]{e}^{i{k}_{f}(r-a-{v}_{f}t)}}{{\widehat{q}}^{2}[(1-\widehat{\rho})(\mathrm{sin}\text{\hspace{0.17em}}\widehat{q}/\widehat{q})-\mathrm{cos}\text{\hspace{0.17em}}\widehat{q}+i\widehat{\rho}\widehat{\upsilon}\text{\hspace{0.17em}}\mathrm{sin}\text{\hspace{0.17em}}\widehat{q}]}\sum _{n=1}^{N}{e}^{-i{\mathbf{k}}_{f}\xb7{\mathbf{r}}_{n}}.$$## 2.3.

### OA Estimation of Blood ${\mathrm{SO}}_{2}$

Blood ${\mathrm{SO}}_{2}$ can be estimated if the OA signals at two laser wavelengths are measured and it can be obtained by evaluating,^{16}

## (5)

$${\mathrm{SO}}_{2}=\frac{{P}_{\mathrm{OA}}^{{\lambda}_{2}}{\u03f5}_{\mathrm{Hb}}^{{\lambda}_{1}}-{P}_{\mathrm{OA}}^{{\lambda}_{1}}{\u03f5}_{\mathrm{Hb}}^{{\lambda}_{2}}}{{P}_{\mathrm{OA}}^{{\lambda}_{1}}\mathrm{\Delta}{\u03f5}_{\mathrm{Hb}}^{{\lambda}_{2}}-{P}_{\mathrm{OA}}^{{\lambda}_{2}}\mathrm{\Delta}{\u03f5}_{\mathrm{Hb}}^{{\lambda}_{1}}}.$$## 3.

## Numerical Method

## 3.1.

### Physical Properties of a Red Blood Cell

RBCs have a biconcave shape under normal physiological conditions. Therefore, it would not be possible to obtain a closed form expression for the OA field generated by an erythrocyte. It is also expected that the OA field amplitude would depend upon the direction of measurement. However, to the best of our knowledge the effect of shape on the angular distribution of OA field amplitude has not been evaluated so far. Therefore, an erythrocyte was approximated as a fluid sphere in this work with a volume of $87\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mu {\mathrm{m}}^{3}$ and radius, $a=2.75\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mu \mathrm{m}$. The density and speed of sound in a cell were taken as ${\rho}_{s}=1092\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{kg}/{\mathrm{m}}^{3}$ and ${\upsilon}_{s}=1639\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{m}/\mathrm{s}$, respectively.^{30} The acoustical properties for the surrounding medium (saline water) were taken as ${\rho}_{f}=1005\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{kg}/{\mathrm{m}}^{3}$ and ${v}_{f}=1498\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{m}/\mathrm{s}$, respectively.^{30} The numerical values of the thermal properties for a RBC were chosen as $\beta =1\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{K}}^{-1}$ and ${C}_{P}=1\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{J}\text{\hspace{0.17em}}\mathrm{Kg}\text{\hspace{0.17em}}{\mathrm{K}}^{-1}$. Further, fluence of the illuminating source was considered to be $F=1\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{J}\text{\hspace{0.17em}}{\mathrm{m}}^{-2}$. Note that thermal parameters and fluence of the optical radiation only control the amplitude of the signals and do not affect their spectral features.

The strength of the OA field emitted by an erythrocyte is governed by the amount of optical energy absorbed by it and thus, is proportional to the absorption coefficient of the cell. The absorption coefficient of a cell depends upon the concentration and oxygenation states of the hemoglobin molecules. Approximately 280 million hemoglobin molecules are packaged in an erythrocyte^{31} and accordingly, the concentration of hemoglobin molecules could be computed to be $5.34\times {10}^{-3}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{moles}/\mathrm{l}$. The optical absorption coefficient of an erythrocyte could be estimated by employing the relationship, $\mu ={c}_{{\mathrm{HbO}}_{2}}{\u03f5}_{{\mathrm{HbO}}_{2}}+\phantom{\rule{0ex}{0ex}}{c}_{\mathrm{Hb}}{\u03f5}_{\mathrm{Hb}}$. Here, ${c}_{{\mathrm{HbO}}_{2}}$ and ${c}_{\mathrm{Hb}}$ are the concentrations of oxyhemoglobin and deoxyhemoglobin molecules in a cell, respectively. A cell with 90% ${\mathrm{SO}}_{2}$ can be interpreted as 90% of the hemoglobin molecules are fully oxygenated (i.e., four oxygen molecules are attached to each hemoglobin molecule) and 10% of them are deoxygenated (i.e., no oxygen molecule is attached to a hemoglobin molecule). The molar extinction coefficients for oxyhemoglobin and deoxyhemoglobin molecules could be found in Refs. 32 and 33 at the relevant optical wavelengths. Table 1 lists the estimated optical absorption coefficients for a cell at various oxygenation states for three incident laser radiations that are commonly used in practice. Note that ${\mathrm{SO}}_{2}$ state of a RBC also determined the blood ${\mathrm{SO}}_{2}$ level because the constituent RBCs were assumed to belong at the same state of oxygenation in the blood sample.

## Table 1

Absorption coefficients of a cell at various oxygenation states for different incident laser wavelengths.

SO2 (%) | μ(m−1) at 700 nm | μ(m−1) at 1000 nm | μ(m−1) at 1064 nm |
---|---|---|---|

0 | 2206.61 | 254.30 | 49.13 |

10 | 2021.61 | 354.80 | 85.01 |

20 | 1836.62 | 455.31 | 120.90 |

30 | 1651.62 | 555.81 | 156.78 |

40 | 1466.62 | 656.31 | 192.67 |

50 | 1281.63 | 756.81 | 228.55 |

60 | 1096.63 | 857.31 | 264.44 |

70 | 911.63 | 957.81 | 300.32 |

80 | 726.64 | 1058.31 | 336.21 |

90 | 541.64 | 1158.82 | 372.09 |

100 | 356.64 | 1259.32 | 407.98 |

## 3.2.

### Generation of the Spatial Organizations of Cells

RBCs are the most numerous blood corpuscles (98% of particles in the blood as RBCs) and their concentration is about 5 million cells per cubic millimeter of blood.^{27} Furthermore, hemoglobin molecules are the main light absorbing chromophores contained in RBCs for the laser wavelengths that are generally used in OA experiments.^{10} Accordingly, RBCs are thought to be the dominant OA sources in a blood sample. Therefore, the contributions from the white blood cells and platelets were neglected in this study. The RBCs were placed in a 2-D space mimicking a blood sample and this area was referred to as the region of interest (ROI). A 2-D ROI was used in this study because 2-D simulations are relatively computationally inexpensive compared to three-dimensional. Our group has used 2-D simulations successfully to interpret experimental results in the context of ultrasound backscattering from tissue.^{34}^{,}^{35} The number of cells that occupied the 2-D ROI was obtained by dividing the product of the size of the ROI and the hematocrit of the blood sample with the largest cross-sectional area of a RBC. The hematocrit is defined as the fractional volume occupied by the RBCs in a blood sample. Nevertheless, for a 2-D sample, the fractional area covered by the RBCs determines the hematocrit. A 45% hematocrit was considered in this study because this corresponds to the normal level of hematocrit in human blood.^{27} The spatial locations of the RBCs within the ROI were generated by implementing a Monte Carlo algorithm known as the random sequential adsorption (RSA) technique forming a 2-D tissue mimicking configuration.^{36} In this technique, coordinates of a RBC were randomly proposed and accepted if the particle did not overlap with the existing cells under the periodic boundary conditions. The proposed move was canceled if the nonoverlapping condition was not satisfied and a new move was initiated. This procedure was repeated until a valid location was identified. In this way, random positions of cells were generated. Figure 2 displays some 2-D simulated tissue realizations with randomly positioned cells. To improve the visibility of the figures, a smaller ROI is presented in each case. Each gray circle in Fig. 2(a) represents an deoxygenated cell for which ${\mathrm{SO}}_{2}=0\%$. Similarly, partially and fully oxygenated cells are indicated by dark and black circles in Fig. 2(b) and 2(c), respectively.

## 3.3.

### Computation of the Nonbandlimited OA Signals

The time dependent wide bandwidth OA signal for a simulated 2-D tissue configuration was obtained by evaluating Eq. (3) and by considering contributions from a wide range of frequencies (MHz to GHz). The trapeziodal rule was used to carry out the integration in Eq. (3) at each time point to obtain the complex valued signal and that consequently provided the instantaneous OA signal and its amplitude. The signal was computed at a large distance ($r=6000\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mu \mathrm{m}$) from the center of the ROI and in the backward direction (in the opposite direction to the direction of light propagation). For this case, the size of the ROI was fixed to $500\times 500\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mu {\mathrm{m}}^{2}$ and that corresponded to $N=4735$. A computer code was written in $C$ and executed on a IBM server (OS: Red Hat Enterprise Linux 5.3EL, Processor: Intel Xeon 3 GHz, RAM: 8 GB, Quad Core, 32-bit). The simulated data were imported in MATLAB R2009b for analysis. For each ${\mathrm{SO}}_{2}$ level, 200 RF lines were generated from 200 different tissue realizations and that took about 9 min and 30 s to run. The associated power spectrum and amplitude histogram, averaged over 200 lines, were obtained to study spectral and time domain features. The best fit Rayleigh distribution curve to the histogram was constructed by using a MATLAB optimization function.^{37}

## 3.4.

### Computation of the Bandlimited OA Signals

Three Gaussian functions modeling diagnostic and high frequency ultrasound transducers were used (with center frequencies, ${f}_{0}=2$, 10, and 50 MHz, and ${\omega}_{0}=2\pi {f}_{0}$) to generate bandlimited OA signals. A 80% receiving bandwidth ($-6\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{dB}$) was considered for each case and the associated $\sigma $ value was computed using the relationship, $2\sigma \sqrt{2\text{\hspace{0.17em}}\mathrm{log}(2)}=\phantom{\rule{0ex}{0ex}}0.8{\omega}_{0}$. Equation (4) was computed to simulate finite bandwidth signals for each ${\mathrm{SO}}_{2}$ level. The first two transducers were chosen to match the diagnostic ultrasound frequency range. Moreover, a 2 MHz transducer was used extensively for conducting *in vitro* and *in vivo* blood ${\mathrm{SO}}_{2}$ measurements.^{17}^{,}^{18} The third transducer was considered because a single element focused transducer with 50 MHz as the center frequency was employed in many OA microscopy and tomography measurements.^{10}^{,}^{15}^{,}^{16} A ROI of $2500\times 500\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mu \mathrm{m}}^{2}$ was chosen for this case. The light beam propagated along the longer dimension. For such a ROI, both constructive and destructive interfering events for acoustic waves took place even in the low frequency range (e.g., 2 MHz signals).

## 4.

## Numerical Results

## 4.1.

### Results for the Nonbandlimited OA Signals

A simulated OA RF line is plotted in Fig. 3(a) for a sample with ${\mathrm{SO}}_{2}=50\%$, when irradiated by the 700 nm laser source. The temporal variation of the OA signal envelope is also shown in Fig. 3(a). The rapid fluctuation of the OA amplitude can be observed in Fig. 3(a) due to the presence of high frequency components. The associated average signal envelope histogram obtained from 200 RF lines and the best fit Rayleigh distribution curve are presented in Fig. 3(b). The Rayleigh distribution provides an excellent fit to the histogram. The mean power spectrum lines, constructed from 200 RF lines for each case, for a series of ${\mathrm{SO}}_{2}$ levels are displayed in Fig. 3(c) over a wide range of frequencies (MHz to GHz). The ultrasound power spectral slopes at different bandwidths are comparable for all spectral lines. This is because, the effect of interference at any frequency (i.e., $\sum _{n=1}^{N}{e}^{-i{\mathbf{k}}_{f}\xb7{\mathbf{r}}_{n}}$) is similar for all samples. However, ultrasound power spectra differ in amplitudes due to different absorption coefficients associated with the blood samples under examination. It can also be noticed over the entire frequency range that spectral intensity decreases with increasing blood ${\mathrm{SO}}_{2}$. The spectral intensity at ${\mathrm{SO}}_{2}=0\%$ is about 36 times more than that of 100%. It can be linked to the fact that for this laser wavelength the absorption coefficient ($\mu $) for a deoxygenated RBC is approximately 6 fold greater than that of a fully oxygenated RBC (see Table 1). Figure 3(d) shows that the OA peak to peak amplitude decreases linearly with blood ${\mathrm{SO}}_{2}$. This is expected since the OA signal amplitude in this model is linearly proportional to $\mu $ [see Eqs. (2) and (3)]. Furthermore, $\mu $ determines the blood ${\mathrm{SO}}_{2}$ and thus the signal amplitude exhibits linear variation with blood ${\mathrm{SO}}_{2}$.

Figure 4(a) demonstrates a computed OA RF line and the associated signal envelope for a blood sample with ${\mathrm{SO}}_{2}=50\%$ and for 1000 nm incident optical radiation. It may be noted that the signal strength has decreased as compared to that of the previous case [Fig. 3(a)]. This is expected because $\mu $ at the 1000 nm is less than that of 700 nm at this ${\mathrm{SO}}_{2}$ level (see Table 1). The reduction of signal amplitude is also clear from Fig. 4(b), where average envelope histogram computed from 200 RF lines is presented; the width of the histogram is reduced in comparison to that in Fig. 3(b). The best fit Rayleigh distribution curve is also plotted in Fig. 4(b) and that provides an excellent fit to the histogram. The mean power spectrum curves are plotted in Fig. 4(c) for various blood samples. The spectral intensity at any frequency increases as the blood ${\mathrm{SO}}_{2}$ level increases. Figure 4(d) illustrates how the peak to peak amplitude varies with blood ${\mathrm{SO}}_{2}$. The peak to peak amplitude at ${\mathrm{SO}}_{2}=100\%$ is nearly 5 times higher than that of 0%. This enhancement can be attributed to the fact that $\mu $ for a fully oxygenated RBC is approximately 5 times greater than that of a deoxygenated RBC (see Table 1).

Figure 5 plots the estimated blood ${\mathrm{SO}}_{2}$ as a function of actual ${\mathrm{SO}}_{2}$. The blood ${\mathrm{SO}}_{2}$ has been assessed by using the simulated OA signals for the 700 and 1000 nm incident laser sources and by evaluating Eq. (5). This figure shows that the estimated value matches the actual ${\mathrm{SO}}_{2}$. Thus, the accurate estimation of blood ${\mathrm{SO}}_{2}$ is possible if the OA signal amplitude varies linearly with blood ${\mathrm{SO}}_{2}$.

## 4.2.

### Results for the Bandlimited OA Signals

A typical OA signal computed for a blood sample with ${\mathrm{SO}}_{2}=50\%$, illuminated by the 700 nm light beam and detected by a receiver with 2 MHz as the center frequency is shown in Fig. 6(a). The signal envelope is also plotted in Fig. 6(a). The two vertical lines drawn in Fig. 6(a) correspond to the front and back surface of the sample. It can be seen that strong OA signals are detected at the two edges due to coherent addition of signals for the cells close to the boundaries. The signals from the proximal and distal edges relative to the transducer are out of phase by $\pi $ radians with respect to each other. Moreover, the signal strength in the central region is reduced due to the destructive interference of waves emitted by the individual cells. The associated average envelope histogram, obtained from 100 OA RF lines, is shown in Fig. 6(b). The counts are maximal in the low amplitude range corresponding to the OA signal from the central region of the ROI. Otherwise, there are almost uniform counts for other amplitudes because of smooth variation of the signal envelope. As expected, the peak to peak amplitude decreases linearly as the blood ${\mathrm{SO}}_{2}$ increases [see Fig. 6(c)]. The second row demonstrates time domain signal features for a 10 MHz transducer. The signal strength as shown in Fig. 6(d) has increased in comparison to that of Fig. 6(a), since the OA signal intensity at 10 MHz is higher than that of 2 MHz at any ${\mathrm{SO}}_{2}$ level [see Fig. 3(c)]. The envelope histogram is presented in Fig. 6(e). The counts are mostly due to the signal from the central region. Figure 6(f) displays that the peak to peak amplitude decreases monotonically with blood ${\mathrm{SO}}_{2}$. Nevertheless, its magnitude at any ${\mathrm{SO}}_{2}$ level is more than that of Fig. 6(c). Similar characteristics can also be observed from the third row, where the signals have been filtered by a 50 MHz detector. It may be mentioned that the OA signal rise time decreases when high frequency components are detected [see Fig. 6(a), 6(d), and 6(g)]. Also, as the high frequency content increases the signal envelope histogram moves toward the Rayleigh distribution [see Fig. 6(b), 6(e), and 6(h)].

The simulated OA signals for three ultrasound detectors and signal properties are shown in Fig. 7, for the 1000 nm heating optical radiation. The RF line processed by a 2 MHz receiver is shown in Fig. 7(a) for a blood sample with ${\mathrm{SO}}_{2}$. The variation of signal envelope is also outlined in Fig. 7(a). The OA signal amplitude distribution is presented Fig. 7(b). The RF line, associated envelope and the envelope histogram look similar to that of the other light source [see Fig. 6(a) and 6(b)]. However, the signal strengths are different for these two cases. Figure 7(c) demonstrates that the peak to peak amplitude increases linearly with blood ${\mathrm{SO}}_{2}$, as has been observed for the nonbandlimited signals [see Fig. 4(d)]. Similar signal properties have been exhibited by the signals filtered by the 10 and 50 MHz receivers as shown in Fig. 7(d)–7(i). The blood ${\mathrm{SO}}_{2}$ has also been assessed (data not shown) using the output signals for these filtering functions and obtained perfect agreement between the estimated and actual value for each receiver.

## 5.

## Discussion

A single particle approach was used in this paper to generate the OA field from a collection of RBCs suspended in a saline solution. One of the assumptions was that the RBCs were the dominant OA sources. The surrounding fluid medium did not absorb light energy and hence, it did not generate any significant OA signal. However, this assumption might not be satisfied in some cases. For example, the light absorption of water increases significantly for wavelengths $>1000\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{nm}$. The absorption of blood serum in the presence of glucose, found in blood samples obtaining from patients with diabetes mellitus, may also become important in the same optical wavelength range.^{38} For such situations, a background signal may arise from the surrounding medium. Therefore, this theoretical framework needs to be modified in order to incorporate the light absorption property of the surrounding medium and to examine theoretically how it affects OA emission of a sample.

Another assumption was that the ${\mathrm{SO}}_{2}$ level of each cell was identical in a blood sample. This assumption should hold if the local oxygen would homogeneously be consumed by the cells. In the process, hemoglobin molecules would become oxygenated and accordingly, determine the cellular absorption coefficient. The model predicts a linear variation of the OA amplitude with the cellular absorption coefficient [see Eqs. (3) and (4)]. Consequently, the OA amplitude varies linearly with blood ${\mathrm{SO}}_{2}$. In another study, it was found that the OA signal amplitude exhibited a slightly nonlinear variation with blood ${\mathrm{SO}}_{2}$, when binary mixtures of cells were considered (such as fully oxygenated and deoxygenated RBCs).^{23} In that work, the blood ${\mathrm{SO}}_{2}$ was increased gradually by replacing the deoxygenated RBCs with the equal number of fully oxygenated RBCs. Therefore, although the proportion of fully oxygenated RBCs linearly increased, the resultant OA signal amplitude was governed by the complex interaction of waves with different strengths from these two populations. In that case, the OA signal did not follow the linear variation [see Eq. (8) of Ref. 23].

The simulated results are consistent with published experimental results. For example, Esenaliev et al.^{17} measured that slope and amplitude of the OA signal increased with incremental change of blood ${\mathrm{SO}}_{2}$ [see Figs. 3(a) and 8(b) of Ref. 17]. This finding could be interpreted as a consequence of enhancement of bulk tissue light absorption coefficient due to oxygenation. In Fig. 8(a), the temporal variation of the simulated OA signal amplitude at the front edge of the ROI is plotted for four representative ${\mathrm{SO}}_{2}$ levels. The OA signals were simulated from a blood sample at 15% hematocrit illuminated by a 1064 irradiating source and detected with a 2 MHz receiver. The associated cellular absorption coefficients are given in Table 1 at various ${\mathrm{SO}}_{2}$ levels. It appears from Fig. 8(a) that the slope and the amplitude of the OA signal increases with increasing blood ${\mathrm{SO}}_{2}$. The signals generated using the theoretical framework presented here also exhibited similar trends.

The same signals were also used to compare the simulation and experimental results. In Fig. 8(b), the experimental data ($\pm \mathrm{SD}$), fitted line at ${\mathrm{SO}}_{2}=100\%$, are plotted against blood ${\mathrm{SO}}_{2}$ [experimental data taken from Fig. 8(c) of Ref. 17]. The simulated data, normalized by the amplitude at ${\mathrm{SO}}_{2}=100\%$, are also shown in Fig. 8(b). It can be seen that the simulated trend is in excellent agreement with the linear fit to the experimental data.

An interesting point to mention is that the signals from the two edges of a sample are out of phase by $\pi $ radians with respect to each other [see Fig. 6(a), 6(d), 6(g) and Fig. 7(a), 7(d), 7(g)]. It is known that the OA pulse emitted by a sphere is N-shaped originating from the nature of thermoelastic expansion.^{25} Therefore, the surface facing toward the detector for each cell close to the front edge of the sample contributes to a positive pressure at the detector (due to the positive sections of the N-shaped curves from the RBCs arriving at the detector in phase). However, the distant surface from the receiver for each cell adjacent to the back edge of the sample creates a negative pressure. These arguments might explain as to why the OA signals arising from the edges maintain a phase difference of $\pi $ radians.

In addition to that, it is clear from Fig. 6(a), 6(d), 6(g) and Fig. 7(a), 7(d), 7(g) that the simulated OA signal amplitudes were detected before the time that it would have taken the OA disturbance to reach the detector from the sample front surface (note that the temporal window corresponding to the front and back edge of the ROI is marked by the two vertical lines in these figures). Therefore, it appears that the OA signals violated the causality condition. Another related observation is that the slope of the OA signal envelope increased as the detected high frequency content increased [see Fig. 6(a), 6(d), and 6(g)]. These might be explained by examination of the effect of the generation of a bandlimited pulse from a N-shaped pulse containing high frequency components emitted by a single sphere (cell). If the contributions from the low frequency components were only considered with a Gaussian amplitude distribution for the detector [e.g., 2 MHz case as in Eq. (4)], the constructed output pulse would become broader in the time domain than that of the actual N-shaped broadband pulse generated from an individual cell. In this case, the leading edge would appear before the arrival of the center of the pulse and and this is what we observed in these simulations. On the other hand, for the second or third Gaussian function corresponding to the high frequency ultrasound detector, the generated pulse would be narrower and thus the leading edge would be closer to the center of the pulse resulting in faster signal rise compared to that of the lower frequency ultrasound detector.

## 6.

## Conclusions

A modeling approach based on the single particle theory is presented to examine the blood ${\mathrm{SO}}_{2}$ dependent OA signal properties. The OA signals were simulated from 2-D tissue configurations corresponding to different blood ${\mathrm{SO}}_{2}$ levels. The fractional number of oxygenated hemoglobin molecules, packaged in a cell, defined the cellular absorption coefficient and also fixed the state of blood ${\mathrm{SO}}_{2}$. Both nonbandlimited and bandlimited signals were generated. For the nonbandlimited case, the OA signal strength decreased and increased linearly with blood ${\mathrm{SO}}_{2}$ for the 700 and 1000 nm laser source, respectively. The signal amplitude at ${\mathrm{SO}}_{2}=100\%$ was found to be nearly 6 times lower and 5 times greater than that of ${\mathrm{SO}}_{2}=100\%$ at these wavelengths, respectively. The power spectral lines also showed similar trends over the entire frequency range (MHz to GHz). The estimated blood ${\mathrm{SO}}_{2}$ level demonstrated a good agreement with the actual ${\mathrm{SO}}_{2}$. For the bandlimited case, three Gaussian functions (with center frequencies 2, 10, and 50 MHz and 80% bandwidth for each receiver) were considered. The OA signal amplitude processed by each filtering function also varied in identical manners for these optical radiations. These findings are consistent with published experimental results validating the model qualitatively. This study also suggests that the accuracy of the OA technique, estimating blood ${\mathrm{SO}}_{2}$ using two light wavelengths, does not depend on the bandwidth of a receiver. Future work will include introducing modifications to the code that will allow taking into account how the optical and ultrasound fields are perturbed by the surrounding tissues by incorporating effects such as inhomogeneous optical distributions and ultrasound frequency dependent attenuation.

## Acknowledgments

R. K. Saha is grateful to the Saha Institute of Nuclear Physics for providing an opportunity to carry out research work in photoacoustics. M. C. Kolios and E. Hysi would like to thank the following funding agencies: 1. the Canadian Institutes of Health Research (CIHR) grant MOP-97959, 2. CIHR grant from the Terry Fox Foundation, 3. the Canada Research Chairs Program awarded to M. C. Kolios, and 4. financial support received from the Natural Sciences and Engineering Research Council of Canada by E. Hysi.

## References

*in vivo*using photoacoustic microscopy,” Appl. Phys. Lett. 90, 053901, 1–3 (2007).APOPAI0003-6935 Google Scholar

*in vivo*,” IEEE Trans. Med. Imaging 24(4), 436–440 (2005).ITMID40278-0062 http://dx.doi.org/10.1109/TMI.2004.843199 Google Scholar

*in vivo*imaging,” Nat. Biotechnol. 24(7), 848–851 (2006).NABIF91087-0156 http://dx.doi.org/10.1038/nbt1220 Google Scholar

*In vivo*Monitoring of blood oxygenation in large veins with a triple-wavelength optoacoustic system,” Opt. Express 15(24), 16261–16269 (2007).OPEXFF1094-4087 http://dx.doi.org/10.1364/OE.15.016261 Google Scholar

*in vivo*by picosecond time-resolved spectrophotometry,” J. Biochem. 109(3), 456–461 (1991). Google Scholar