Expand this Topic clickable element to expand a topic
Skip to content
Optica Publishing Group

Fast single photon avalanche photodiode-based time-resolved diffuse optical tomography scanner

Open Access Open Access

Abstract

Resolution in diffuse optical tomography (DOT) is a persistent problem and is primarily limited by high degree of light scatter in biological tissue. We showed previously that the reduction in photon scatter between a source and detector pair at early time points following a laser pulse in time-resolved DOT is highly dependent on the temporal response of the instrument. To this end, we developed a new single-photon avalanche photodiode (SPAD) based time-resolved DOT scanner. This instrument uses an array of fast SPADs, a femto-second Titanium Sapphire laser and single photon counting electronics. In combination, the overall instrument temporal impulse response function width was 59 ps. In this paper, we report the design of this instrument and validate its operation in symmetrical and irregularly shaped optical phantoms of approximately small animal size. We were able to accurately reconstruct the size and position of up to 4 absorbing inclusions, with increasing image quality at earlier time windows. We attribute these results primarily to the rapid response time of our instrument. These data illustrate the potential utility of fast SPAD detectors in time-resolved DOT.

© 2015 Optical Society of America

1. Introduction

Diffuse optical tomography (DOT) and fluorescence mediated tomography (FMT) are emerging as important tools in biomedical research [1–10]. DOT utilizes multiple measurements between source and detector pairs through bulk biological tissue, and by solving the resulting inverse problem retrieves the internal distribution of intrinsic optical properties. It is well-recognized that the high-degree of light scatter in biological tissue is the primary reason for the relatively poor imaging resolution achievable for both DOT and FMT [11]. As such, improving this resolution – through either new instrument designs or image reconstruction strategies - is a persistent area of interest in the field.

Among the reported measurement modes for DOT, time-resolved (TR) measurement of transmitted light from a pulsed laser is known to be the most information rich, but also the most technically complex. Several strategies for utilizing TR data have been reported in the literature [12–18]. Of particular interest for this paper, a number of groups [19–24] (including ours) have studied the improvement in imaging resolution obtained by measuring “early-arriving photons”. These have traversed preferentially less-scattered paths between the source and detector pair, and therefore yield narrower “photon density sensitivity functions” (PDSFs; alternatively known as the Jacobian or “banana” functions). This in turn results in a better-conditioned inverse problem and consequently improved imaging resolution.

Our group has studied the early photon (EP) effect in detail in recent years, and in particular the influence of the measurement time, optical properties, and instrument configuration on the width of the PDSF. We demonstrated that (under otherwise equivalent conditions) the reduction in PDSF width relative to time-integrated, quasi-continuous wave (CW) measurements were dependant on the instrument temporal impulse response function (TIRF) width [25]. In particular, instruments with shorter TIRFs allowed measurements closer to the theoretical PDSF width predicted by time-resolved Monte Carlo photon propagation models [26]. Our recent survey of the time-resolved DOT and FMT literature [25] indicated that reported instruments have intrinsic TIRFs on the order of 80-800 ps, the fastest of which are achievable with multi-channel plate (MCP) detectors or streak cameras. These are generally expensive, and in the case of MCPs, extremely sensitive to photo-damage. Most recently, we studied the performance of a new fast single photon avalanche photodiode (SPAD; ID-100-50; ID Quantique, Geneva, Switzerland) detector. This detector type resulted in an overall instrument TIRF full width at half maximum (FWHM) of 59 ps, which is significantly faster than that of most TR-DOT and TR-FMT instruments. As we and other authors have noted [27–29], it also has a number of other apparent advantages, including, i) high quantum efficiency (of 35%), ii) resistance to photo-damage compared to, e.g. photomultiplier tubes (PMTs) and MCPs, iii) relatively low cost, and iv) relatively compact size (8 x 6.1 x 3.9 cm). We demonstrated that the experimentally measured PDSF was reduced by 65% relative to the quasi-CW case using this SPAD. This was the narrowest PDSF we have able to measure with any instrument configuration that we have tested previously. As such, this type of detector has significant promise for TR-DOT and potentially in FMT imaging. As we also discuss, SPADs have a number of disadvantages, including the small detector area (50 μm x 50 μm) which results in relatively poor geometric collection efficiency for diffuse (unfocused) light, and the long after-pulsing tail which results from diffusion processes inside the semi-conductor.

In this paper, we have extended this proof-of-concept work and constructed a complete SPAD-based TR small animal DOT system optimized for early-photon measurements. This instrument consisted of three fast SPADs, a femtosecond Titanium:Sapphire laser and multiple rotation stages to provide complete 360 degree angular sampling of the sample. An orthogonally mounted charge coupled device (CCD) camera was used to image the profile of the sample and extract its outline. We first describe the design of this instrument and then demonstrate its performance in imaging tissue-like optical phantoms with multiple embedded attenuating targets. As we show this instrument is able to clearly resolve up to 4 small attenuating targets embedded deeply in the phantom. We attribute this performance primarily to the rapid TIRF of the instrument. This work validates the potential of fast SPAD detectors in TR-DOT instrument design. We plan to extend this work to fluorescence measurements in the future to test performance in a TR-FMT configuration.

2. Materials and methods

2.1 Instrument design

A schematic diagram and photograph of the SPAD-based time-resolved tomographic imager is shown in Figs. 1(a) and 1(b), respectively. We used a femtosecond Titanium: Sapphire laser (MaiTai XF-1, Newport corporation, Irvine, CA) operating at 750 nm with 80 MHz repetition rate. To control the power at the sample, we used a variable attenuator (VA) configuration, which consisted of a polarizing beam-splitter and rotatable half wave plate; the laser power transmitted or reflected at 90° could be controlled by changing the angle of the half-wave plate. This configuration also allowed us to use the reflected portion of the beam to generate a synchronizing trigger signal with an optical constant fraction discriminator (OCFD) unit (Becker & Hickl, Berlin, Germany)

 figure: Fig. 1

Fig. 1 (a) schematic diagram and (b) photograph of the SPAD-based DOT imager. See text for details.

Download Full Size | PDF

Samples to be imaged were mounted on an inner rotation stage controlled by stepper motor (PK245-01AA, Velmex Inc). The transmitted light from the sample was collimated with a pair of planoconvex lens (f = 25 mm; Edmund Optics) and measured simultaneously with three SPAD detectors which were mounted on a larger motorized rotation stage (PK266-03A-P1, Velmex). The detectors were separated by fixed angles of 35° (this was the minimum spacing we could physically achieve). However, by rotating both rotation stages, an arbitrary number of source-and-detector measurement combinations could be generated. The electronic outputs of the SPADs were connected to a router (HRT-82; Becker & Hickl) which was coupled into a time correlated single photon counting (TCSPC) module (SPC-130; Becker & Hickl) to perform parallel acquisition of time-resolved signals on the 3 channels. This router is capable of parallel acquisition of up to 8 channels, so that additional detectors could in principle be added in the future.

In order to recover the boundaries of the sample, we also mounted a digital monochromatic camera (Chameleon USB 2.0; Point Grey Research, Inc) at the side of the stage, oriented at 90° to the laser axis. This was focused on the sample surface. As the sample was rotated the camera acquired images for each lateral view. As discussed below, these images were analyzed in post processing to recover the outline of the sample.

2.2 Optical phantoms

We first tested the performance of the instrument using a custom made optical phantom with symmetrical, cylindrical geometry as shown in Figs. 2(a)-2(c). We fabricated a 25 mm diameter × 100 mm long resin phantom with TiO2 (Sigma-Aldrich) and India ink added to approximately match the optical properties of biological tissue at near infrared wavelengths. Specifically, the reduced scattering coefficient (µ's) was ~15 cm−1 and the absorption coefficient (µa) was ~0.1 cm−1 [30]. The inside of the cylinder was removed (machined) so that the phantom had a hollow core with 22 mm inner diameter. A set of 12, 2 mm diameter holes were drilled in the solid base, distributed at different positions as shown in Fig. 2(a). Absorbing rods (approximately 30 times the absorption coefficient of the background) could be inserted into specific combinations of holes. The core of the phantom was then filled with a liquid solution containing 1% Intralipid (Liposyn II; Hospira), 50 ppm India ink, so as to generate an approximately homogenous cylinder with small absorbing inclusions. Here, we used either black absorbing rods (with approximately 30 times background attenuation), or 2 mm diameter straws filled with intralipid solution with additional India Ink added to create either 2, 4, 6, or 8 times the background µa. This design allowed us to flexibly select the position and number of absorbing inclusions for testing and characterization of our instrument and image reconstruction algorithm. In future work, this will also allow us to test the effects of background attenuation, autofluorescence, and wavelength dependent optical properties.

 figure: Fig. 2

Fig. 2 (a) diagram showing the distribution and number of holes for the cylindrical phantom, as well as (b) top and (c) side view of the phantom. (d) diagram, (e) top view, and (f) side view of the irregularly-shaped optical phantom.

Download Full Size | PDF

Second, to test the ability of the system to image irregularly shaped samples, a similar approach was used to fabricate a phantom with the geometry shown in Figs. 2(d)-2(f). In this case, four holes were drilled in the base for mounting absorbing rods. The phantom was filled with matching intralipid solution as above.

2.3 Data collection

Measurement of Instrument Temporal Response. The overall TIRF of the instrument was measured by placing a thin diffusive plate at the location of the sample. The laser power was set to 0.2 mW, which produced a count rate of approximately 106 counts / s at the detector. The TCSPC hardware was configured to acquire data for 2 s with a 6.25 ns time range and 6 ps electronic resolution with a temporal resolution of 6 ps.

Phantom Measurements. Optical phantoms were placed on the sample stage and rotated through 360° with 1° increments. After each full sample rotation, the outer stage was rotated by ± 10° (for a total of 3 positions) and the acquisition was repeated. The position of the illumination laser spot was not moved during the experiment, but could be scanned in the future. Therefore, the detectors were positioned at a total of 9 positions corresponding ± 45°, ± 35°, ± 25°, ± 10° and 0° angles from the optical axis, and a total of 9 x 360 rotations = 3240 measurements generated. The total scan time was about 36 minutes. However, the scanning procedure could be reduced in the future, for example, by rotating the inner stage with 5° increments, thereby resulting in a scan time of less than 10 minutes. For these measurements, the laser power was 20 mW, resulting in an approximate photon count rate of 106 counts / s. As above, the acquisition time for each measurement was 2 s, resulting in approximately 104 photon counts per detector channel at maximum of the transmitted time-resolved curve.

2.4 Tomographic image reconstruction

2.4.1 Sample boundary reconstruction

We first analyzed the monochrome CCD camera images to extract the sample (phantom) physical boundary. For each rotation image (corresponding to the 360 rotation angles), we analyzed pixel values along the scanning axial plane and located the two edges of the phantom. We set all pixels between the edges (i.e. inside the phantom) to 1, and all pixels outside the phantom to zero. We then applied an inverse radon transform to reconstruct the boundary of the target phantom. The reconstructed image was then down-sampled onto the same grid used in the diffuse optical tomography reconstruction. Specifically, we used a 51 x 51 grid, for a 2.5 x 2.5 cm2 image, corresponding to a 0.5 x 0.5 mm2 pixel size.

2.4.2 Diffuse optical tomography reconstruction

Forward Modeling: Forward light modeling in diffuse optical tomography is normally performed using some approximation to the Boltzmann Transport Equations (BTE), or by numerical (Monte Carlo) implementation of the BTE. In the case of imaging with “quasi-CW” (i.e. time-integrated) photons, the time-independent diffusion approximation is frequently used, and is generally considered to yield acceptable accuracy in macroscopic imaging geometries (sample length > 10 / (3(µ's + µa))), when optical scattering greatly exceeds absorption (i.e. µ's >> µa), and away from object boundaries. In this case, the “weight matrix” W, between a source and detector pair for the linearized image reconstruction problem is given by:

W(rs,rd)=U0(r,rs)G(rd,r)dr3
where, rs and rd are the positions of the source and detector, respectively and r is a coordinate within the diffusive medium [31]. U0 is the photon density in a homogeneous medium from a source, and G is the Green’s function describing light propagation to the detector [32]. U0 and G can be represented by the solution of diffusion equation for a point source in an infinite medium:
Φ(r)=14πDrexp(μaDr)
where Φ is the photon fluence rate and D = [3(μ's + μa)]−1.

Modeling of time-dependent photon propagation at early times from a pulsed laser source is a more complex problem. Time-resolved Monte Carlo (TR-MC) simulations [33] have been used previously because they are extremely accurate for arbitrary geometries, but are very computationally expensive, particularly at early time points where very small numbers of photons are detected. Simulation times are therefore very long even on modern computers with hardware acceleration, when Jacobians for many source-and-detector combinations must be computed. For this reason, we and others have alternatively used the time-dependent diffusion approximation to the BTE. We have noted previously that this results in slight overestimation of the full-width at half maximum of the Jacobian at early times [34], however the error is expected to be small (less than 10%). However, use of the diffusion approximation is by no means inherent to the system presented here and time-resolved Monte Carlo could be performed in the future. Here W for a given time gate t for a source and detector position (rs and rd, respectively) is given by:

W(rs,rd,t)=0tU0(r,rs,τ)G(rd,r,tτ)dτdr3
where τ is an integration variable for time, and the time-dependent solution of BTE has the form expressed is:
Φ(r,t)=c(4πDct)3/2exp(r24Dctμact)
Image Reconstruction: For each set of time-resolved measurements through the phantom, we analyzed both the early photon data and the quasi-CW (time integrated) data as in Fig. 3(a). For the quasi-CW case, the data was summed over time bins up to 10%-of-peak on the falling edge of the transmitted time-resolved curve, so as to remove the contribution of the semiconductor diffusion tail from the data (see below). For early-photon data types, we summed data in three time windows corresponding to, i) 1% to 15% of the peak on the rising edge of the time-resolved transmission curve (which we refer to as the EP1 data type in this paper), ii) 1%-50% on the rising edge (EP2 data type), and iii) 1%-100% on the rising edge (EP3 data type). The rationale for this windows is as follows: EP1 represents “very early” photons measurable with the SPAD detectors, producing narrower PDSFs than we were able to achieve previously, and EP2 represents a broader time gate that – according to our previous work approximates the effect of the longer TIRF of a PMT-based instrument, with correspondingly wider PDSF [25, 26]. EP3 represents a later gate which trades the benefits of early-photons and the signal-to-noise properties of quasi-CW data.

 figure: Fig. 3

Fig. 3 (a) Example measured time-resolved normalized transmitted Photon density curve and time points used in calculating the EP1, EP2 and EP3 data types. (b) The overall measured system TIRF showing the rapid instrument response (c) Example measured time-resolved normalized transmitted photon density (log scale). Carrier processes in the semiconductor detector result in a “diffusion tail” (red arrow) which is independent from the light diffusion through scattering media. (d) The measured SNR for the different data types used in these analyses.

Download Full Size | PDF

We then used a linear image reconstruction strategy where the forward problem is posed as a system equations: [y] = [W][δμa]. Here, y is the set of calibrated source-detector pair measurements. We used the calibration method described by Li et. al. [35], which required measurement of a ‘homogeneous’ phantom, i.e. one where no absorbing targets were included. This was scanned with the identical geometry as the target phantoms. W is the weight matrix and δμa is the absorption due to the inclusions. For each data type, corresponding weight matrices were computed using Eq. (1) and (3) above. The image was then computed by inverting this system of equations used randomized algebraic reconstruction technique (r-ART) [36]. The number of iterations and regularization parameter λ were held identical in all cases and were 300 and 0.8, respectively.

2.4.3 Reconstruction quality–dice similarity coefficient

To quantitatively evaluate the quality of the reconstructed images, we introduced the Dice similarity coefficient (DSC) which is an index of the spatial overlap between two images [37]. In our case, we define

DSC(A,B)=2X(AB)/(A+B)
where A is the binary image of the true phantom configuration (known target locations), and B is a binary image generated from the reconstructed axial image. The operator ∩ represents the overlap between the two images. In order to convert the reconstructed images to binary images, we used a threshold γ=0.25 of maximum, so that all pixels exceeding γ were set to 1, and all pixels less than γ were set to 0. Use of this threshold is somewhat arbitrary but results in inclusion of more of the reconstructed image (as well as associated image noise and artifacts) in the analysis than, for example, a threshold of 0.5 of maximum. The value of DSC ranges from 0 (indicating no spatial overlap) to 1 (indicating complete spatial overlap).

3. Results and discussion

3.1 Basic instrument characterization

An example measured instrument TIRF for our system is shown in Fig. 3(b). The TIRF full width at half maximum (FWHM) was 59ps, which is in agreement with our previous work with the same detector but with a different instrument configuration [26]. As noted previously, short instrument TIRF is important for improved imaging performance using early photons, since it is inversely correlated with the width of the photon density sensitivity function (PDSF) between source and detector pairs in DOT [25].

An example measured transmitted time-resolved curve through our 2.5 cm thick, cylindrical tissue-mimicking phantom (Fig. 2(a)) is shown in Fig. 3(c). These data are plotted on a logarithmic scale for clarity. As shown, when measured with the SPAD, a relatively long tail at late times is observed, which is a measurement artifact due to carrier diffusion effects in the semiconductor [38]. This well-understood effect distorts the decay portion of the transmitted time-resolved curve, which may be of concern for some applications, e.g. time-resolved measurement of optical properties fluorescence lifetime measurements or imaging with late-photons [28]. However, in the analysis presented here, we are primarily interested the rise portion of the curve, where distortion due to this effect is not observed. As noted, for the quasi-CW we summed only to the 10%-of-peak point on the falling edge of the TR curve (point V, Fig. 3(a)) to remove this artifact. This point corresponds approximately to the change in slope from this effect, both in the TIRF and transmitted TR data here and in our previous work [26].

We next analyzed the signal to noise ratio (SNR) for the 4 different data types that we consider in this paper: EP1, EP2, EP3 and time-integrated quasi-CW data. We took 100, 2 s measurements, and the corresponding SNRs (defined as 20 log(I/σ), where I is the mean signal, and σ is the standard deviation) are shown in Fig. 3(d). As anticipated, the early photon data types have between 28 dB (for EP1) and 36 dB (for EP3) lower SNRs compared to quasi-CW data. These data are generally consistent with our earlier results and re-iterate the noise-resolution trade-off inherent to imaging with early photons. Nevertheless, even at the earliest time window we were able to obtain an SNR of 28 dB, which, as we show, resulted in good quality image reconstructions.

3.2 Phantom experiments

Transmission Measurements at Different Time Points: We tested several combinations of absorbers (position and number) in our cylindrical tissue-mimicking phantom. Three example configurations are shown in Fig. 4. These include combinations of two (Fig. 4(a)), three (Fig. 4(b)) and four (Fig. 4(c)) absorbing rods placed at relatively deep positions in the cylinder. For the two-rod phantom configuration, the rods were separated by 7 mm (center-to-center). Rod separations for the example 3- and 4- rod phantom configurations ranged between 6.9 and 13.7 mm center-to-center, depending on the pair.

 figure: Fig. 4

Fig. 4 (a-c) Example 2, 3 and 4-absorber configurations for the cylindrical phantom. Normalized measured data are shown as a function of the rotation angle for EP1, EP2, EP3 and Quasi-CW data types is shown for (d-f) the detector at location D in the diagrams (0 degrees to the optical axis, and for (g-i) the detector at location D’ in the diagrams (45 degrees from the optical axis).

Download Full Size | PDF

Phantoms were scanned with the SPAD imager for a single axial slice in each case. Example normalized projection data are shown for the 3 phantom examples. Traces corresponding to quasi-CW data (blue line), EP1 data (red line), EP2 data (purple line) and EP3 data (green line) as a function of the rotation angle are shown in Figs. 4(d)-4(i).

Figures 4(d)-4(f) shows example data for the SPAD detector position aligned with the optical axis (i.e. at 0°; labeled as D in Figs. 4(a)-4(c)) and data from a detector at a position offset 45° from the optical axis (labeled as D’ in Figs. 4(a)-4(c)) are shown in Figs. 4(g)-4(i). For each configuration, the earlier time data shows more structure from the absorber attenuation than the quasi-CW case. For example, small “dips” are evident for EP data at approximately 250° rotation angle in Fig. 4(h) and approximate 130° angle in Fig. 4(i) which are not evident for the CW case. These are likewise more pronounced for earlier EP data types (i.e. EP1 vs. EP2 vs. EP3). These additional features are evident in the raw data, and directly contribute to improved image reconstructions in the next section. As anticipated, the data corresponding to earlier windows (e.g. EP1) had more evident noise, as would be expected from our SNR data (Fig. 3(d)). We have previously presented strategies for trading of noise and resolution properties in time-resolved DOT data sets [39] as have others [40].

Cylindrical Phantom Reconstructions: Fig. 5 shows the image reconstructions obtained for the three phantom configurations shown in Fig. 4. Figures 5(a)-5(c) show the reconstructions obtained from quasi-CW data. Figures 5(d)-5(f), 5(g)-5(I), and Figs. 5(j)-5(l) show reconstructions obtained using the EP3, EP2 and EP1 data types, respectively. The true location and size of the absorbing heterogeneities are indicated with green circles for each reconstruction for comparison. As shown, the EP data types yield significantly better performance in image domain compared to quasi-CW data. For phantoms with either two and three absorbing heterogeneities (Figs. 5(a), 5(d), 5(g), 5(j) and Figs. 5(b), 5(e), 5(h), 5(k)), the number and approximate position of the targets could be retrieved with quasi-CW data. However, significantly worse resolution, localization and separation are evident compared to the EP data. Likewise, image reconstructions were more accurate with EP1 data than EP3 data. The advantages of EP data are most evident when viewing the phantom with four absorbing inclusions (Figs. 5(c), 5(f), 5(i), 5(l)). The two deeper absorbing inclusions merge into an extended blurred object for quasi-CW data (Fig. 5(c)), whereas they were more clearly separated for EP data, and in particular with the EP1 data type. We also performed Dice similarity coefficient (DSC [37];) analysis as shown in Fig. 5(m). As indicated a strong increase in DSC was observed for all phantoms when using EP data types, and in general the earlier the data, the higher the DSC value. In particular, DSC increased by between 48% and 165% when considering images obtained using EP1 data versus quasi-CW data. In one specific case (4-target) the EP1 data type produced a marginally worse DSC value than EP2. We note that DSC is an imperfect metric in evaluating DOT images because small offsets in the center of reconstructed targets often result in significantly reduced DSC value. Visual inspection of Figs. 5(i) and 5(l) indicate that the EP1 image is at least comparable quality.

 figure: Fig. 5

Fig. 5 Reconstructed images of the example phantoms with 2, 3 and 4 embedded absorbing targets for: (a-c) quasi-CW data, (d-f) EP3 data, (g-i) EP2 data, and (j-l) EP1 data. (m) Dice similarity coefficient data for the 3 phantom configurations and 4 data types is also shown.

Download Full Size | PDF

Sensitivity at different time gates: To test the sensitivity of the instrument at different time gates, we used a series of phantoms with two, 2 mm diameter cylindrical inclusions (straws) filled with between 2 and 8 times the background absorption. The results are shown in Fig. 6. These data indicate that, for this particular geometry, the minimum contrast that could be resolved with CW and EP3 data was 2:1, whereas this was around 4:1 for EP2 and EP1. At lower contrast (2:1), the heterogeneities could be detected, but their position and concentration could not be correctly recovered. This is a well-known effect, since earlier time gates involve rejection of large numbers of transmitted photons, and subsequent reduction of signal to noise (as in Fig. 3(d), and in [39]). The exact sensitivity for any given sample is highly dependent on the size, position and number of the inhomogeneities. However, in general the earliest time gates are expected to yield less sensitivity.

 figure: Fig. 6

Fig. 6 Reconstructed images for 2-target phantom with either 8:1 (a-d), 4:1 (e-h), or 2:1 (j-l) contrast for CW (a,e,i), EP3 (b,f,j), EP2 (c,g,k) and EP1 (d,h,l) data types. The intensity versus contrast is shown (m-p) normalized to the highest contrast case for, By inspection, the position or contrast of the targets could not be recovered for 2:1 contrast for EP1 and EP2 data types (dashed rectangle).

Download Full Size | PDF

Irregularly Shaped Phantom Reconstructions: To verify that the system is capable of imaging asymmetrical targets (such as small animals) we scanned our irregularly-shaped phantom, and the results are shown in Fig. 7. The outline of the phantom was extracted using the boundary detection and inverse-radon method described in section 2.41 above, as shown in Fig. 7(a). The reconstructed shape was subsequently down-sampled to the 0.5 x 0.5 mm2 DOT reconstruction grid and filled in to generate the axial slice, as shown in Fig. 7(b). Reconstructions were then performed with quasi-CW data (Fig. 7(c)), EP3 data types (Fig. 7(d)), EP2 data type (Fig. 7(e)) and EP1 data type (Fig. 7(f)). The true location and size of the absorbing targets (green circles) are also indicated in each panel. As above, the EP data type resulted in a significantly improved image reconstruction compared to that produced with the quasi-CW data. DSC data for each data type is shown in Fig. 7(g). Here, the DSC metric increased by 185% for reconstructions performed using the EP1 data type compared to the quasi-CW data type. We note that the quasi-CW reconstruction in this case was relatively poor, and in general could potentially be improved, e.g. through alternate scanning geometries, forward model calculations or inversion routines. However, the reconstructions shown here were performed using the same methodology and indicate the improvement obtained in this case.

 figure: Fig. 7

Fig. 7 (a) The boundary reconstruction of the irregularly-shaped phantom shown in Figs. 3d-f. (b) The down-sampled mask of axial slice used in the reconstructed images. Reconstructions obtained with (c) quasi-CW data, (d) EP3 data, (e) EP2, and (f) EP1 data are shown. (g) DSC data for the 3 data types.

Download Full Size | PDF

4. Discussion and conclusions

In summary, we designed, constructed and validated a novel SPAD-based time-resolved DOT system. We were primarily interested in the use of fast SPAD detectors with respect to measurement of early-arriving photons in DOT and FMT [26]. As we have demonstrated previously, the fast temporal response time results in correspondingly narrower instrument PDSF widths and an improved inverse image reconstruction problem. In the context of fluorescence imaging, convolution of the excitation light with the fluorescence exponential decay results in temporal broadening of the transmitted time-resolved curve. However, in principle, measurement with a faster instrument TIRF would result in less distortion on the rising edge of the curve during measurement [25] independent of this process.

Our previous data [26] indicated that at very early times (less than 15% of maximum on the rising edge of the transmitted time-resolved curve), SPAD detectors resulted in reduction of PDSF width by about 10% compared to PMTs. The relationship between this narrowness and imaging resolution is complex (due to the ill-posedness of the reconstruction problem), however, our earlier work [39] implies that these detectors would yield a further ~10% improvement in resolution versus PMTs under identical imaging conditions. We re-iterate that in addition to their rapid TIRFs, SPADs have a number of practical advantages including their low cost and insensitivity to photo-damage compared to, for example, PMTs and multi-channel plate (MCP) detectors. Potential limitations for specific DOT and FMT applications include the small detector active area (50 µm x 50 µm) and the carrier “diffusion tail” which can distort the decay portion of the curve.

We tested our instrument in imaging custom fabricated phantoms with both cylindrical and irregular geometries. These phantoms allowed testing of multiple combinations of target absorber configurations. The reconstructions shown in Figs. 5–7 demonstrate the high-accuracy DOT images that can be obtained with this instrument, particularly with early photon data types. We re-iterate that accurate localization of four distinct heterogeneities is extremely difficult in diffusive optical imaging. A literature search failed to reveal another example where this has been reported for purely optical systems, either at the clinical or small-animal scale. Evaluation of our reconstructions using the Dice similarity coefficient (DSC) showed that reconstructed images with the earliest photons increased the DSC between 48% and 185% compared to quasi-CW photons over all the phantoms we imaged. Moreover, the DSC generally increased when earlier time windows were used (EP3 vs EP2 vs EP1). As we have noted, DSC is an imperfect metric in evaluating DOT images, however it does provide a useful metric in quantifying the improvement that is evident by inspection of the images. We attribute this imaging performance primarily to the rapid response time of the SPAD detectors. Likewise, use of early photons represents an inherent tradeoff since it entails rejection of large numbers of transmitted photons. In practice, this results in corresponding reduction in SNR and ultimately imaging sensitivity.

In combination, these data underscore the potential of fast SPAD detectors for time-resolved DOT. We next plan to incorporate fluorescence measurements into the instrument, and apply it to small animal FMT imaging.

Acknowledgments

This research was supported by the National Institutes of Health (R01EB012117; NIBIB).

Reference and links

1. J. P. Culver, R. Choe, M. J. Holboke, L. Zubkov, T. Durduran, A. Slemp, V. Ntziachristos, B. Chance, and A. G. Yodh, “Three-dimensional diffuse optical tomography in the parallel plane transmission geometry: evaluation of a hybrid frequency domain/continuous wave clinical system for breast imaging,” Med. Phys. 30(2), 235–247 (2003). [CrossRef]   [PubMed]  

2. N. Deliolanis, T. Lasser, D. Hyde, A. Soubret, J. Ripoll, and V. Ntziachristos, “Free-space fluorescence molecular tomography utilizing 360° geometry projections,” Opt. Lett. 32(4), 382–384 (2007). [CrossRef]   [PubMed]  

3. X. Montet, J. L. Figueiredo, H. Alencar, V. Ntziachristos, U. Mahmood, and R. Weissleder, “Tomographic fluorescence imaging of tumor vascular volume in mice,” Radiology 242(3), 751–758 (2007). [CrossRef]   [PubMed]  

4. S. Patwardhan, S. Bloch, S. Achilefu, and J. Culver, “Time-dependent whole-body fluorescence tomography of probe bio-distributions in mice,” Opt. Express 13(7), 2564–2577 (2005). [CrossRef]   [PubMed]  

5. R. B. Schulz, J. Ripoll, and V. Ntziachristos, “Experimental fluorescence tomography of tissues with noncontact measurements,” IEEE Trans. Med. Imaging 23(4), 492–500 (2004). [CrossRef]   [PubMed]  

6. B. J. Tromberg, B. W. Pogue, K. D. Paulsen, A. G. Yodh, D. A. Boas, and A. E. Cerussi, “Assessing the future of diffuse optical imaging technologies for breast cancer management,” Med. Phys. 35(6), 2443–2451 (2008). [CrossRef]   [PubMed]  

7. G. Zacharakis, J. Ripoll, R. Weissleder, and V. Ntziachristos, “Fluorescent protein tomography scanner for small animal imaging,” IEEE Trans. Med. Imaging 24(7), 878–885 (2005). [CrossRef]   [PubMed]  

8. K. M. Tichauer, R. W. Holt, K. S. Samkoe, F. El-Ghussein, J. R. Gunn, M. Jermyn, H. Dehghani, F. Leblond, and B. W. Pogue, “Computed tomography-guided time-domain diffuse fluorescence tomography in small animals for localization of cancer biomarkers,” J. Vis. Exp. 65, e4050 (2012). [PubMed]  

9. F. Leblond, S. C. Davis, P. A. Valdés, and B. W. Pogue, “Pre-clinical Whole-body Fluorescence Imaging: Review of Instruments, Methods and Applications,” J. Photochem. Photobiol. B 98(1), 77–94 (2010). [CrossRef]   [PubMed]  

10. F. Leblond, K. M. Tichauer, R. W. Holt, F. El-Ghussein, and B. W. Pogue, “Toward whole-body optical imaging of rats using single-photon counting fluorescence tomography,” Opt. Lett. 36(19), 3723–3725 (2011). [CrossRef]   [PubMed]  

11. S. R. Arridge, “Optical tomography in medical imaging,” Inverse Probl. 15(2), R41–R93 (1999). [CrossRef]  

12. B. Das, F. Liu, and R. Alfano, “Time-resolved fluorescence and photon migration studies in biomedical and model random media,” Rep. Prog. Phys. 60(2), 227–292 (1997). [CrossRef]  

13. D. Hall, G. Ma, F. Lesage, and Y. Wang, “Simple time-domain optical method for estimating the depth and concentration of a fluorescent inclusion in a turbid medium,” Opt. Lett. 29(19), 2258–2260 (2004). [CrossRef]   [PubMed]  

14. S. Keren, O. Gheysens, C. S. Levin, and S. S. Gambhir, “A comparison between a time domain and continuous wave small animal optical imaging system,” IEEE Trans. Med. Imaging 27(1), 58–63 (2008). [CrossRef]   [PubMed]  

15. A. T. Kumar, S. B. Raymond, G. Boverman, D. A. Boas, and B. J. Bacskai, “Time resolved fluorescence tomography of turbid media based on lifetime contrast,” Opt. Express 14(25), 12255–12270 (2006). [CrossRef]   [PubMed]  

16. A. T. Kumar, S. B. Raymond, A. K. Dunn, B. J. Bacskai, and D. A. Boas, “A time domain fluorescence tomography system for small animal imaging,” IEEE Trans. Med. Imaging 27(8), 1152–1163 (2008). [CrossRef]   [PubMed]  

17. M. S. Patterson, B. Chance, and B. C. Wilson, “Time resolved reflectance and transmittance for the non-invasive measurement of tissue optical properties,” Appl. Opt. 28(12), 2331–2336 (1989). [CrossRef]   [PubMed]  

18. J. Wu, L. Perelman, R. R. Dasari, and M. S. Feld, “Fluorescence tomographic imaging in turbid media using early-arriving photons and Laplace transforms,” Proc. Natl. Acad. Sci. USA 94(16), 8783–8788 (1997). [CrossRef]   [PubMed]  

19. K. Chen, L. T. Perelman, Q. Zhang, R. R. Dasari, and M. S. Feld, “Optical computed tomography in a turbid medium using early arriving photons,” J. Biomed. Opt. 5(2), 144–154 (2000). [CrossRef]   [PubMed]  

20. F. Leblond, H. Dehghani, D. Kepshire, and B. W. Pogue, “Early-photon fluorescence tomography: spatial resolution improvements and noise stability considerations,” J. Opt. Soc. Am. A 26(6), 1444–1457 (2009). [CrossRef]   [PubMed]  

21. M. Niedre and V. Ntziachristos, “Comparison of fluorescence tomographic imaging in mice with early-arriving and quasi-continuous-wave photons,” Opt. Lett. 35(3), 369–371 (2010). [CrossRef]   [PubMed]  

22. M. J. Niedre, R. H. de Kleine, E. Aikawa, D. G. Kirsch, R. Weissleder, and V. Ntziachristos, “Early photon tomography allows fluorescence detection of lung carcinomas and disease progression in mice in vivo,” Proc. Natl. Acad. Sci. U.S.A. 105(49), 19126–19131 (2008). [CrossRef]   [PubMed]  

23. G. M. Turner, A. Soubret, and V. Ntziachristos, “Inversion with early photons,” Med. Phys. 34(4), 1405–1411 (2007). [CrossRef]   [PubMed]  

24. L. Zhao, H. Yang, W. Cong, G. Wang, and X. Intes, “L p regularization for early gate fluorescence molecular tomography,” Opt. Lett. 39(14), 4156–4159 (2014). [CrossRef]   [PubMed]  

25. N. Valim, J. Brock, M. Leeser, and M. Niedre, “The effect of temporal impulse response on experimental reduction of photon scatter in time-resolved diffuse optical tomography,” Phys. Med. Biol. 58(2), 335–349 (2013). [CrossRef]   [PubMed]  

26. Y. Mu, N. Valim, and M. Niedre, “Evaluation of a fast single-photon avalanche photodiode for measurement of early transmitted photons through diffusive media,” Opt. Lett. 38(12), 2098–2100 (2013). [CrossRef]   [PubMed]  

27. A. Puszka, L. Di Sieno, A. D. Mora, A. Pifferi, D. Contini, G. Boso, A. Tosi, L. Hervé, A. Planat-Chrétien, A. Koenig, and J.-M. Dinten, “Time-resolved diffuse optical tomography using fast-gated single-photon avalanche diodes,” Biomed. Opt. Express 4(8), 1351–1365 (2013). [CrossRef]   [PubMed]  

28. M. Mazurenka, L. Di Sieno, G. Boso, D. Contini, A. Pifferi, A. D. Mora, A. Tosi, H. Wabnitz, and R. Macdonald, “Non-contact in vivo diffuse optical imaging using a time-gated scanning system,” Biomed. Opt. Express 4(10), 2257–2268 (2013). [CrossRef]   [PubMed]  

29. Y. Mu and M. Niedre, “A fast SPAD-based small animal imager for early-photon diffuse optical tomography,” in Engineering in Medicine and Biology Society (EMBC),201436th Annual International Conference of the IEEE(IEEE, 2014), pp. 2833–2836.

30. M. J. Niedre, G. M. Turner, and V. Ntziachristos, “Time-resolved imaging of optical coefficients through murine chest cavities,” J. Biomed. Opt. 11, 064017 (2006).

31. V. Ntziachristos, C.-H. Tung, C. Bremer, and R. Weissleder, “Fluorescence molecular tomography resolves protease activity in vivo,” Nat. Med. 8(7), 757–761 (2002). [CrossRef]   [PubMed]  

32. L. V. Wang and H.-i. Wu, Biomedical optics: principles and imaging (John Wiley & Sons, 2012).

33. J. Chen, V. Venugopal, and X. Intes, “Monte Carlo based method for fluorescence tomographic imaging with lifetime multiplexing using time gates,” Biomed. Opt. Express 2(4), 871–886 (2011). [CrossRef]   [PubMed]  

34. N. Valim, J. Brock, and M. Niedre, “Experimental measurement of time-dependent photon scatter for diffuse optical tomography,” J. Biomed. Opt. 15(6), 065006 (2010). [CrossRef]   [PubMed]  

35. C. Li and H. Jiang, “A calibration method in diffuse optical tomography,” J. Opt. A, Pure Appl. Opt. 6(9), 844–852 (2004). [CrossRef]  

36. R. J. Gaudette, D. H. Brooks, C. A. DiMarzio, M. E. Kilmer, E. L. Miller, T. Gaudette, and D. A. Boas, “A comparison study of linear reconstruction techniques for diffuse optical tomographic imaging of absorption coefficient,” Phys. Med. Biol. 45(4), 1051–1070 (2000). [CrossRef]   [PubMed]  

37. K. H. Zou, S. K. Warfield, A. Bharatha, C. M. C. Tempany, M. R. Kaus, S. J. Haker, W. M. Wells 3rd, F. A. Jolesz, and R. Kikinis, “Statistical validation of image segmentation quality based on a spatial overlap index,” Acad. Radiol. 11(2), 178–189 (2004). [CrossRef]   [PubMed]  

38. A. Tosi, A. Dalla Mora, F. Zappa, A. Gulinatti, D. Contini, A. Pifferi, L. Spinelli, A. Torricelli, and R. Cubeddu, “Fast-gated single-photon counting technique widens dynamic range and speeds up acquisition time in time-resolved measurements,” Opt. Express 19(11), 10735–10746 (2011). [CrossRef]   [PubMed]  

39. Z. Li and M. Niedre, “Hybrid use of early and quasi-continuous wave photons in time-domain tomographic imaging for improved resolution and quantitative accuracy,” Biomed. Opt. Express 2(3), 665–679 (2011). [CrossRef]   [PubMed]  

40. S. S. Hou, W. L. Rice, B. J. Bacskai, and A. T. Kumar, “Tomographic lifetime imaging using combined early- and late-arriving photons,” Opt. Lett. 39(5), 1165–1168 (2014). [CrossRef]   [PubMed]  

Cited By

Optica participates in Crossref's Cited-By Linking service. Citing articles from Optica Publishing Group journals and other participating publishers are listed here.

Alert me when this article is cited.


Figures (7)

Fig. 1
Fig. 1 (a) schematic diagram and (b) photograph of the SPAD-based DOT imager. See text for details.
Fig. 2
Fig. 2 (a) diagram showing the distribution and number of holes for the cylindrical phantom, as well as (b) top and (c) side view of the phantom. (d) diagram, (e) top view, and (f) side view of the irregularly-shaped optical phantom.
Fig. 3
Fig. 3 (a) Example measured time-resolved normalized transmitted Photon density curve and time points used in calculating the EP1, EP2 and EP3 data types. (b) The overall measured system TIRF showing the rapid instrument response (c) Example measured time-resolved normalized transmitted photon density (log scale). Carrier processes in the semiconductor detector result in a “diffusion tail” (red arrow) which is independent from the light diffusion through scattering media. (d) The measured SNR for the different data types used in these analyses.
Fig. 4
Fig. 4 (a-c) Example 2, 3 and 4-absorber configurations for the cylindrical phantom. Normalized measured data are shown as a function of the rotation angle for EP1, EP2, EP3 and Quasi-CW data types is shown for (d-f) the detector at location D in the diagrams (0 degrees to the optical axis, and for (g-i) the detector at location D’ in the diagrams (45 degrees from the optical axis).
Fig. 5
Fig. 5 Reconstructed images of the example phantoms with 2, 3 and 4 embedded absorbing targets for: (a-c) quasi-CW data, (d-f) EP3 data, (g-i) EP2 data, and (j-l) EP1 data. (m) Dice similarity coefficient data for the 3 phantom configurations and 4 data types is also shown.
Fig. 6
Fig. 6 Reconstructed images for 2-target phantom with either 8:1 (a-d), 4:1 (e-h), or 2:1 (j-l) contrast for CW (a,e,i), EP3 (b,f,j), EP2 (c,g,k) and EP1 (d,h,l) data types. The intensity versus contrast is shown (m-p) normalized to the highest contrast case for, By inspection, the position or contrast of the targets could not be recovered for 2:1 contrast for EP1 and EP2 data types (dashed rectangle).
Fig. 7
Fig. 7 (a) The boundary reconstruction of the irregularly-shaped phantom shown in Figs. 3d-f. (b) The down-sampled mask of axial slice used in the reconstructed images. Reconstructions obtained with (c) quasi-CW data, (d) EP3 data, (e) EP2, and (f) EP1 data are shown. (g) DSC data for the 3 data types.

Equations (5)

Equations on this page are rendered with MathJax. Learn more.

W( r s , r d )= U 0 (r, r s )G( r d ,r)d r 3
Φ(r)= 1 4πDr exp( μ a D r )
W( r s , r d ,t)= 0 t U 0 ( r , r s ,τ)G( r d ,r,tτ)dτd r 3
Φ(r,t)= c (4πDct) 3/2 exp( r 2 4Dct μ a ct )
DSC( A,B ) = 2 X ( AB )/ ( A+B )
Select as filters


Select Topics Cancel
© Copyright 2024 | Optica Publishing Group. All rights reserved, including rights for text and data mining and training of artificial technologies or similar technologies.