Optoacoustic imaging has emerged as a medical technology capable of providing tomographic reconstructions of tissue structures with submillimeter resolution at depths up to several centimeters.1 Several groups have demonstrated such capabilities in tissue phantoms 2, 3, 4, 5, 6, 7, 8, 9 as well as in in-vitro 10, 11 and in-vivo systems.12, 13, 14, 15 Optoacoustic methods rely on the generation and detection of thermoelastic stress waves (or thermal expansion) generated within a target, following the absorption of pulsed laser radiation. Thus, optoacoustic imaging exploits key benefits of both optical imaging and ultrasound by combining the high optical penetration and high endogenous/exogenous optical contrast of tissue structures in the red/near-infrared with the relatively nondispersive propagation of stress waves in tissues. This provides for the generation of images at moderate sample depths with spatial resolution. These depth and resolution characteristics for functional imaging lie between those provided by optical coherence tomography and diffuse optical imaging.
Numerous methods are available to measure the time-resolved thermoelastic response of a target that can be subsequently used for image reconstruction. The most common approach is to employ piezoelectric transducers placed in acoustic contact with the tissue surface, either directly, or through the use of a coupling liquid. 3, 13, 16, 17, 18, 19, 20, 21 While this approach has been successful, there are some disadvantages. First, the requirement of acoustic contact between the transducer and sample makes detection of the stress waves within the same area as the laser irradiation cumbersome. Second, both the sensitivity of nonfocused piezoelectric transducers, and the horizontal resolution limit of the reconstructed image, scale with the detector diameter.21, 22 As a result, a tradeoff exists between high sensitivity and high spatial resolution when using nonfocused piezoelectric transducers. Alternative signal acquisition methods for optoacoustic tomography have also been implemented, including optical detection of the laser-induced stress waves 4, 7, 8, 23, 24, 25, 26 and interferometric measurement of surface displacement.9, 27, 28, 29 The former, while still requiring acoustic contact, can provide fast acquisition of the stress field over a relatively large field of view. The latter, which forms the basis for this work, utilizes an interferometric system to provide noncontact measurements of time-resolved surface displacement resulting from the absorption of pulsed laser radiation. This approach possesses several advantages over techniques based on piezoelectric or optical stress transducers. First, no acoustic coupling is required between the target and detection system. Second, the probe beam of the interferometer can be easily focused to a dimension of ; a size roughly an order of magnitude smaller than that of the smallest piezoelectric transducers currently employed for optoacoustic imaging. This smaller detection area potentially enables an improvement in lateral resolution as compared to current stress detection approaches.
In this work, we demonstrate and analyze the performance of our optoacoustic technique termed POISe (pulsed optoacoustic interferometric spectroscopy and imaging) in imaging heterogeneous tissues. POISe represents a continuation of the efforts published by Yablon 27 and Payne 28, 29 At the core of the system is an interferometer capable of providing an absolute measurement of the time-resolved surface displacement following irradiation with a short laser pulse. In a previous study we have demonstrated the ability to acquire measurements at multiple locations on the tissue surface with nanosecond temporal resolution and angstrom-level displacement sensitivity.9 Moreover, we have demonstrated that the information content of these measurements is sufficient to provide high-quality 2-D images of tissue-like phantoms containing simulated blood vessels with lateral resolution.9 Here we expand on our previous work by reconstructing 3-D tomographic images of tissue-like phantoms as well as chicken chorio-allantoic membrane vasculature in vivo. We demonstrate the ability of POISe to provide lateral and axial spatial resolutions of and , respectively, at sample depths exceeding ten transport mean free paths in turbid media. For biological tissues that possess a characteristic reduced scattering coefficient of , this corresponds to imaging at tissue depths of .
Materials and Methods
Figure 1 provides a schematic of the POISe system. A pulsed Nd:YAG laser (Quantel Brilliant B) delivers pulses at 1064 or to the sample, generating the thermoelastic expansion that is measured by a modified Mach-Zehnder interferometric system. The interferometer employs a polarized He–Ne laser , whose output is split by an acousto-optic modulator (AOM, ATM-1101A1, IntraAction, Bellwood, IL) driven at . A portion of the He–Ne beam travels undisturbed through the AOM, while the rest is deflected into several frequency-shifted beams. One of the first-order frequency-shifted beams is picked off by a mirror and directed to the optical detection system, and forms the reference arm of the interferometer. The unshifted beam enters the sample arm of the interferometer, and passes through a plate and a focusing lens that delivers the beam to the sample surface. The light reflected from the surface is redirected by the plate and a polarization sensitive beamsplitter to the beamsplitter, where it is recombined with the reference beam. The beams emerging from the beamsplitter are directed to two bandwidth Si PIN photodiodes (S5972, Hamamatsu, Japan) operating in a balanced detection setup. The combined photodiode signal is amplified using a broadband low-noise amplifier (ZFL-500LN, Mini-Circuits, Brooklyn, NY) and digitized by a bandwidth digital oscilloscope (TDS-7054, Tektronix, Wilsonville, OR or WR6051, LeCroy Corp., Chestnut Ridge, NY). The AOM driving signal is also digitized by the oscilloscope. For automated operation, the personal computer controls the pump laser, and the two translation stages (Newport Corporation, Irvine, CA) are used to position the sample relative to the pump and/or probe laser beams through GPIB connections.
When the sample surface is at rest, the interferometer signal is a sinusoidal oscillation at the AOM modulation frequency. Motion of the sample surface modifies the oscillation frequency through the Doppler effect. The time-resolved surface displacement is directly proportional to the instantaneous phase difference between the interferometer signal and the AOM driving signal , and is given by:is the wavelength of the interferometer beam and denotes the Hilbert transform and denotes the raw digitized voltages of the AOM driving signal and the interferometer signal, respectively. Our interferometric surface displacement measurement system provides a time resolution of and a displacement sensitivity of .
To acquire the spatial variation of the time-resolved displacement, we must make measurements at multiple surface locations. At the same time, we must maintain a perpendicular relationship between the probe beam and the sample. Two scanning configurations are possible. One option is to scan the sample using an translation stage such that the probe and pump beams are moving in tandem relative to the sample surface. Alternatively, the probe beam can be scanned while holding the pump beam stationary with respect to the sample. Since scanning the probe beam in both directions is impractical in our current system configuration, this scanning mode is accomplished by mounting the final mirror and lens of the interferometer sample arm on a translation stage that travels in the direction, while the sample and the dichroic that redirects the pump beam toward the sample are mounted on a translation stage that travels in the direction. In this latter case, the scanning area is limited by the clear aperture of the dichroic mirror.
For the 3-D image reconstructions, we acquired data at 350 to 750 distinct locations on the sample surface, with 3 to 5 surface displacement traces acquired at each location and subsequently averaged. The pump laser delivered in continuous operation, resulting in a total acquisition time of at each location. Movement of the probe beam relative to the sample between measurement locations required approximately . This measurement strategy resulted in a total data acquisition time in the range of .
Phantoms were constructed with optical properties to simulate a pair of blood vessels within a turbid tissue in the near-infrared. A 2% Intralipid solution was used as the background medium with a reduced scattering coefficient of , absorption coefficient , and at .30 Simulated vessels were fashioned using hollow polyimide tubes (Cole-Parmer, Vernon Hills, IL) with inner and outer diameters of 120 and , respectively. A solution of India ink with was circulated through the tubes using a syringe pump (Harvard Apparatus, Holliston, MA) to simulate the absorption of whole blood at . Phantoms were constructed by mounting two polyimide tubes on a support frame. This assembly was then submerged within a container filled with the 2% Intralipid solution. A Q-switched Nd:YAG laser (Brilliant B, Quantel, Les Ulis, France) operating at with a pulse duration and beam diameter (after expansion) was used as a pump laser, providing an incident radiant exposure of . At this wavelength, the polyimide tubes are essentially transparent, and their contribution to the displacement signal is negligible. Figure 2 provides a schematic of one of the phantoms that was constructed. Also shown is the typical measurement geometry where the -diam pump-laser beam was delivered to the sample, and surface displacement measurements were taken along a line perpendicular to the longitudinal tube axes in a range approximately relative to the midpoint between the tubes. In all the tube phantom experiments, the pump beam was kept stationary with respect to the sample, while the interferometer beam was scanned to provide measurements at multiple locations. The specific scan areas used to generate the images in each of the experiments are described later.
Several phantom configurations were built. All phantoms contained two polyimide tubes placed within the 2% Intralipid solution at depths of . Three main configurations were employed with tubes placed either side-by-side or one below the other with wall-to-wall distances of . The third configuration consisted of tubes crossing at a angle.
The chicken chorio-allantoic membrane (CAM) vasculature was chosen as the initial in-vivo animal model, as it provides easy access to well-developed vasculature that, when covered by an Intralipid layer, mimics a bulk vascular network. Fertilized chicken eggs old were opened and their contents transferred into transparent plastic cups that were then covered with a breathable film and incubated at for five additional days. At day 9, the CAM was almost fully developed, presenting a surface covered with a vascular network. For the purpose of our experiments, a volume of 0.5% Intralipid solution ( and at ) was poured on top of the membrane to provide a turbid layer in thickness. A Q-switched Nd:YAG laser (Quantel) operating at with a pulse duration and beam diameter was used to irradiate the sample with an incident radiant exposure of . For all CAM measurements, the pump and probe laser beams were kept stationary and the sample was translated in both and directions.
Image Reconstruction Algorithm
Backprojection methods provide a simple approach for reconstructing absorbing heterogeneities that generate displacement signals measured by POISe. These methods take advantage of the direct proportionality between the time of arrival of an acoustic signal at a measurement location and the distance between this location and the point of origin of the signal. Thus, our first approach to image reconstruction was to apply a simple delay and sum beam-forming algorithm to surface displacement measurements acquired at distinct surface locations.3 Figure 3 provides a schematic depiction of this algorithm. The tissue volume is divided into small volume elements (voxels), each of which is evaluated as a potential source from which a thermoelastic disturbance can emanate and be detected by measurements taken at different surface locations. The delay and sum method determines the acoustic source intensity for each voxel by finding the time window within each surface displacement trace that corresponds to the distance between the detector and the interior of the current voxel, and performs a weighted summation of the averaged displacement signals in that window for all the detectors. This is given by the expression:is the location in the imaged volume, is the acoustic source intensity corresponding to location , denotes a sum over all the detection points, is a detector-specific weighting factor, is the time-resolved displacement signal from the ’th measurement location, is the location of the ’th detector, is the speed of sound, and is the voxel size. The operator : denotes the set of discrete displacements acquired within the time window of interest, while ⟨ ⟩ denotes averaging the signal over the specified time interval that, in this case, corresponds to the transit time across the voxel in question.
As depicted in Fig. 3, if a large number of recordings contain prominent displacement features corresponding to a certain voxel location, the image reconstruction algorithm assigns a high source strength to that voxel. However, even for voxels that contain no acoustic sources, their corresponding time windows in the individual traces will contain a displacement feature if acoustic sources lying elsewhere in the sample are equidistant from the probe beam location. Thus, in all likelihood, even “empty” voxels will be assigned a small nonzero acoustic source intensity. This “ghost” acoustic source is the reason why “arcing” artifacts are seen in images reconstructed using simple backprojection algorithms. To mitigate this issue, a set of weights can be used to improve different characteristics of the formed image.
For example, weights can be applied to the measured signal to compensate for directional variations of the detector sensitivity3is the angle between the line connecting the measurement location with the current voxel location and the inward pointing normal passing through the measurement location . Weighting can also be used to compensate for geometrical attenuation of the acoustic wave as it travels from the source within the tissue to the detector location
We have also examined the coherence weighting approach introduced by Liao, Li, and Li31:is the location of the voxel currently being reconstructed, is the number of measurement locations, is the corrected acoustic source intensity corresponding to location , is the coherence factor for the same location, and is the time-averaged displacement signal corresponding to the voxel centered at . The coherence factor has a high value only if most of the time-resolved displacement traces possess a signal indicating the presence of acoustic sources at location .31 Consequently, multiplication by reduces the acoustic source intensity assigned to location by Eq. 3 if the location does not contain a real acoustic source. The end result is a reduction in arcing artifacts and an improvement in the signal-to-noise ratio (SNR) of the reconstructed image.
An additional issue that must be considered is that solution to the thermoelastic wave equation for a small spherical absorber results in a displacement feature whose temporal width corresponds to double the spatial extent of the heterogeneity.32 Consequently, one approach to improve imaging resolution is to sum only the contributions from the recorded traces that exhibit a rising profile within the time window corresponding to the current voxel. Henceforth, we refer to the use of this approach as a “rising-edge-only” reconstruction. While quite effective, this technique introduces high frequency noise, because as much as half of the measurement set is discarded. Thus in measurement sets where the signal-to-noise ratio is intrinsically low, the use of rising-edge compensation is not advantageous.
Results and Discussion
Resolution Limits as Determined from Two-Dimensional Image Reconstruction
2-D images of the parallel-tube phantoms (side-by-side and stacked configurations) were reconstructed from measurements taken at 65 locations spanning the interval at increments along a line segment perpendicular to the tube axes. Three time-resolved displacement traces were acquired at each location. Noting the symmetry of the phantom, image reconstruction was performed in a plane perpendicular to the direction of the polyimide tubes along the line of measurement locations. A area was divided in pixels, and the acoustic source intensity for each pixel was determined using the approach outlined in Sec. 2.4. A directivity correction was applied, equal to the cosine of the angle between the normal to the phantom surface and the line connecting each measurement location with the current pixel being reconstructed. In addition, rising-edge-only compensation was employed.
Figure 4 provides two 2-D images reconstructed from measurements acquired on phantoms in which two parallel polyimide tubes with a wall-to-wall separation of were submerged at depths of 4 and , respectively. Lateral line scans of the reconstructed acoustic-source intensity taken through the perceived centers of the tubes are also shown. Figure 5 shows two images reconstructed from measurements acquired on two phantoms, each containing two stacked polyimide tubes with a wall-to-wall separation of submerged at depths of 4 and , respectively. Both the images and the line scans are normalized to the maximum acoustic intensity in the image.
While the cross sections of the tubes appear somewhat flattened, the overall dimensions agree well with the internal diameter of the tubes, and their positions match their true depths within the phantom. To estimate the resolution limits of our technique, we compared the line scans shown in Fig. 4 with the Rayleigh resolution criterion for optical systems. The Rayleigh criterion states that individual fringes are just resolvable when the combined intensity of both fringes at the saddle point is times the maximum intensity of the combined fringe.33 In our line scans, the acoustic source intensity at the saddle point between the tubes is less than 0.71 of the maximum intensity. This is well within the requirements of the Rayleigh criterion. Thus, POISe image reconstructions clearly resolve objects apart, indicating that our technique provides an intrinsic lateral resolution of better than . At the same time, the vertical (depth) line scans shown in Fig. 5 are characterized by a full width at half maximum (FWHM) of , which agrees very well with the actual inner diameter of the polyimide tubes, thus allowing us to speculate that our axial resolution is .
The resolution of POISe images is comparable to that provided by optoacoustic approaches based on time-resolved stress detection.8, 10, 13, 34, 35 Importantly, unlike POISe, these previous studies have utilized both acoustic contact and either access to the sample or highly absorbing heterogeneities. By comparison, POISe is capable of noncontact imaging of phantoms with realistic tissue optical properties and provides lateral resolution and axial resolution at depths approaching in turbid media.
Three-Dimensional Image Reconstruction
In the previous section we have shown the ability of POISe to provide high resolution images in a 2-D geometry facilitated with the use of phantoms with an axis of symmetry. The next step is to investigate the capability of POISe to provide full 3-D tomographic reconstructions.
Figure 6 is a photograph of a tissue-like phantom in which two polyimide tubes with internal diameter were suspended on a plastic support and made to cross at approximately . The support was placed in a container the filled with a 2% Intralipid solution (as described in Sec. 2.2) until the tubes were approximately below the Intralipid surface. Surface displacement measurements were taken within a rectangular area centered at the crossing point at locations spaced at 1 and increments, respectively. This provided an array of measurement locations. For reconstruction, a volume centered at the tube intersection point was divided into voxels. The acoustic source intensity for each voxel was determined as described in Sec. 2.4 employing coherence factors, rising-edge, and directionality corrections.
Figure 7 shows the top view of the 3-D reconstructed volume using the maximum intensity projection (MIP) method,36 while Fig. 8 shows several reconstructed depth slices for the crossed tube phantom at the positions indicated by the dashed lines in Fig. 7. The high resolution demonstrated in the previous section is maintained in the 3-D image reconstructions, as shown by the line scans in Fig. 8e: the vertical line scan indicates a tube height of , while the horizontal line scans indicate a tube width of thickness. However, one must keep in mind that the tubes cross each vertical section at a angle, thus the horizontal size of the tube is elongated by a factor of . Therefore, a corrected estimate for the horizontal direction is . This is consistent with the expected FWHM of a object convolved with the FWHM horizontal point spread function characteristic of POISe imaging. The edges of the first MIP top view [Fig. 7a] seem to fade away, because the pump laser beam has a Gaussian rather than flat-top profile, and thus less energy reaches the edges of the reconstructed area. Such fading can be compensated by normalizing the final image with respect to the irradiance profile of the pump beam, with the result shown in Fig. 7b. This compensation was also employed, as necessary, in the CAM image reconstructions that follow. Due to an oversight on our part, the dataset used to reconstruct the image shown in Fig. 7 was acquired using an incident radiant exposure 20% higher than the ANSI standards for the maximum allowable exposure for skin. However, it is important to note that such a high exposure is not necessary to obtain images of high quality. This can be appreciated by defining an image contrast-to-noise ratio (CNR) as:is the acoustic intensity within a region of interest and and represent the mean and standard deviation of the acoustic intensity in the nearby background region. We consider two pairs of regions of interest (ROIs) located near the center and toward the edges of the image shown in Fig. 7a: ROIs 2 and 4 are contained within the reconstructed tubes, while ROIs 1 and 3 are adjacent locations that are representative of the image background. Application of Eq. 8 for ROIs 1 and 2 provides a , while analysis of ROIs 3 and 4 provides a . A 20% reduction in pump beam radiant exposure will result in a 20% reduction in CNR that does not significantly degrade image quality.
Figure 9 is a photograph of a CAM target prior to the addition of Intralipid, with blood vessels clearly visible on the surface. A layer of 0.5% Intralipid, approximately thick, was added on top of the CAM as described in Sec. 2.3. Time-resolved surface displacement measurements were acquired in a rectangular area at locations spaced at , as indicated by a dashed rectangle in Fig. 9. This resulted in a array of measurement locations. The acoustic source intensity was determined in a volume whose surface outline is indicated by a solid rectangle in Fig. 9. The volume was divided into voxels. The image reconstruction algorithm described in Sec. 2.4 was used without any rising-edge-only, coherence factor, or weighting factor corrections.
Figure 10 shows depth slices reconstructed at several locations within the scanned area. The slice locations are marked with dotted lines as indicated in Fig. 9b. There is very good agreement between the locations of the main blood vessels in the photograph and in the POISe image reconstructions. Slice (a) is taken at the main branch point at the top of the scanned area, and thus only one vessel is present. Slice (b) is taken below the branching point, and thus two vessels are visible. Slice (c) is taken another down where the smaller vessel branches off. This slice shows one vessel on the left, and two closely spaced vessels on the right. Finally, slices (d) and (e) are taken below the second branching point, and show the cross sections of three vessels.
Figures 11 and 12 show the MIP top view of 3-D reconstructed volumes, together with the matching, zoomed-in regions of photographs of the same samples prior to the addition of the Intralipid. These MIP renderings were normalized to the pump beam profile, as mentioned previously. The resemblance between the POISe image and the photograph is striking, with all the major blood vessels clearly visible in the POISe reconstruction. The smallest vessels visible in Fig. 11 are only in diameter and demonstrate the ability of POISe to image small features at depth in a strongly scattering medium . We also performed a contrast-to-noise ratio (CNR) analysis of these two figures. In Fig. 7, the for the pair of ROIs in the center of the image (indicated by the boxes near location 3 in the figure), while for the pair of ROIs at the left edge of the image (indicated by the boxes above location 9). For Fig. 12, the for the ROIs in the center of the image (close to location 1), while for the ROIs at the right edge of the image (close to location 2). Thus a reduction in the incident radiant exposure by a factor of 4 to comply with the ANSI standard of would still result in images with a reasonable CNR.
It is important to note that these POISe reconstruction results were obtained using a grid of surface displacement measurements where the spatial sampling is fairly coarse, i.e., . We expect that increases in the spatial sampling will improve both the spatial resolution and CNR of the POISe reconstructions at the expense of increased image acquisition time. Finally, while the targets measured in this study presented a flat reflective surface that facilitated alignment of the interferometer, measurement of real tissue surfaces do not present an insurmountable challenge. In an earlier study, Payne 28, 29 demonstrated the use POISe to acquire surface displacement data with this interferometric system from rough surfaces, including that of a human forearm in vivo. Thus we fully expect to be able to translate this system for human use.
We demonstrate the ability of the POISe system to provide 2- and 3-D images of heterogeneous tissue phantoms reconstructed from surface displacement measurements made at multiple locations on the target surface. We demonstrate a lateral/horizontal resolution better than at depths exceeding , while the axial/vertical resolution has been estimated to be . We demonstrate the ability of the POISe system to form high resolution images of vascular networks provided by the chicken chorio-allantoic membrane in vivo through a layer of scattering medium. In the future, we intend to investigate the effect of varying the wavelength of the pump laser source (e.g., by means of an optical parametric oscillator), to provide spectroscopic images and the capacity to probe the physiological state of tissues.
We acknowledge support from the NIH through the Laser Microbeam and Medical Program (P41-RR-01192).