2.1.

System Description

Figure 2(a) shows the design of the time-resolved imaging system used in this study and is described in detail in Ref. 15. Briefly, the platform incorporates a tunable Ti-Sapphire laser (Mai Tai HP, Spectra-Physics, USA) used in conjunction with a closed-loop power control system (providing up to 1.8 W power at 800 nm). The time-resolved measurements of the temporal point spread function (TPSF) are made using a gated intensified charge coupled device (ICCD) camera (Picostar HR, LaVision GmbH, Germany). The image intensifier is operated at 300 ps gatewidth at a temporal resolution of 40 ps over a 4 ns time window with images recorded using a 12-bit CCD camera. It should be noted that the noise characteristics of the ICCD are highly dependent on the voltage applied across the image intensifier.¹⁵ In the studies described herein, the intensifier was operated at 570 V. Thus the gated measurements having less than 200 counts were excluded from reconstruction to avoid noise artifacts.¹⁶ The fluorescence measurements are acquired using discrete filter sets (Semrock, USA). A noteworthy feature of the system employed in this study is the incorporation of a Pico projector module (Optoma PK101) as the pattern generation mechanism. The in-built light source in the projector is replaced by the output from the laser. The Pico projector allows the easy implementation of grayscale patterns (via a USB connection using PowerPoint presentation software) with a compact light engine as opposed to the Discovery 1100 kit implemented in the previous incarnation of the system.

Fig. 2

(a) Time-resolved wide-field optical tomography system; (b) grayscale gradient pattern; (c) pattern projected on the imaging stage. The detected pattern was normalized to the full-field uniform pattern to account for nonuniform intensity characteristics; (d) intensity profiles comparing projected pattern with the ground-truth; (e) through (f) temporal characteristics (IRF FWHM and $t_0$) of grayscale patterns are invariant to the gray level generated.

The spatial and temporal characteristics of the system when using the Pico projector were established using a grayscale gradient pattern shown in Fig. 2(b). The objective of this investigation was to determine the range of gray levels that can be accurately projected and detected by the system. Moreover, the temporal characteristics of the instrument response function (IRF) of the system when using the grayscale pattern were evaluated (as the gray levels are generated on the digital micro-mirror device using modulation frequencies in the kHz range) in order to establish the range of gray levels that can be accurately implemented in the optimization procedure without introducing temporal errors during time-gated acquisition. The temporal response of the pattern was therefore acquired for the grayscale pattern and a uniform full-field pattern using a diffusing white paper. The maximum intensities recorded for each pixel for the gradient pattern were normalized by the corresponding intensity values in the uniform pattern. This was done to mitigate error in grayscale estimates due to nonuniform intensity distribution in the projected pattern owing to the projector optics. The normalized grayscale gradient pattern recorded by the camera is shown in Fig. 2(c), and the profile across the center of the pattern shown in Fig. 2(d) shows an accurate ($<10\%$ error) generation of grayscale intensities for levels greater than 0.25. Figure 2(e) to 2(f) compare the photon arrival times ($t_0$) and the IRF width (at half maximum) across the grayscale pattern and the uniform pattern. It is worth noting that the generation of gray levels has a negligible effect on the system IRF. The spatial variation in the IRF characteristics are, however, nonnegligible for both patterns and must be accounted for in the forward model while setting up the inverse problem. Herein, the IRF recorded during pre-experiment calibration is convolved with the time-resolved forward model to accurately model the light propagation in this system. Note that even though the Pico projector is a consumer-grade product with low-quality optics; it allows one to implement a wide-field optical tomography at very low cost ($<\$200$) and can replace advantageously galvo-scanning based systems even in the case of time-resolved studies.

2.2.

Reconstruction Scheme

As stated previously, in this study we employ time-gates as the datatype for reconstruction. Owing to the complex geometry and wide range of optical properties encountered in small animal models, we implemented a time-resolved Monte Carlo (MC) based method to generate an accurate photon propagation forward model. The use of an MC-based method further ensures the accuracy of the model when using early time-gates as datatypes for solving the inverse problem.¹⁷ The MC model used in this study has been described in detail in Ref. 16 and is discussed briefly below.

Consider a volume $\mathrm{\Omega }\in R^3$, with photons incident over an area $A\in R^2$, then the time-resolved signal detected at the point detector at $r_d$ is given by Eq. (1).

Here $r_s$ is a point in $A$, $\eta (r)$ is the effective quantum yield for a fluorophore with quantum yield ($\varphi$), $J^{\mathrm{fluo}}(r_s,r_d,r,t)$ is the spatial and temporal sensitivity function with respect to $\eta (r)$ for the time gate $t$ and total absorption coefficient ${\mu }_a^x$ at the excitation wavelength with a contribution ${\mu }_{\mathrm{af}}^x$ from the fluorophore. The absorption coefficient of the fluorophore is linearly related to the fluorophore concentration via the extinction coefficient. However, the uncertainty in the values of quantum yield and extinction coefficient of a fluorophore in vivo, which may change due to variations in the microenvironment biochemical properties, introduces uncertainties in the reconstruction of the absolute fluorophore concentration biodistribution. In this work, we therefore report relative effective quantum yield as the quantitative reconstructed parameter. For a fluorophore in the medium that has a mono-exponential decay with lifetime $\tau$, under the assumption of equal absorption and scattering coefficients at the excitation and emission wavelengths, the fluorescence signals can be simulated by convolving the temporal signals generated at the excitation wavelength by the exponential time decay of the fluorophore given by $e^{-t/\tau }$. The Jacobians for FMT, $J^{\mathrm{fluo}}(t)$ can therefore be calculated using Eq. (2).

$J^{\mathrm{exc}}(t)$ is the time-gated weight function at the excitation wavelength and $\mathrm{IRF}(t)$ is the experimentally recorded IRF. It should be noted that the experimentally recorded patterns are used to calculate $J^{\mathrm{exc}}(t)$ ensuring that the spatial characteristics of the system are accurately modeled by the MC model. The forward model $J^{\mathrm{exc}}$ is computed in a massively parallel computing environment (Blue Gene, Center for Nanotechnology Innovations, RPI) using 1024 nodes allowing the computation of the time-gated Jacobian in less than 7 min for each pattern and all detectors.

To cast the inverse problem, we employed a Born formulation where the experimental time-domain emission measurements are normalized by the continuous wave (CW) excitation flux at the same position. This normalization efficiently mitigates the dependence of the detected fluorescent signal on the optical properties of the examined tissue.¹⁸ Note also that this formulation employs the CW excitation flux at the same position to alleviate the unavoidable time errors associated with drift and jitter in time-resolved studies. The time-gated Jacobians are computed using the perturbative MC technique, which is the most computationally efficient when large number of detectors are considered.¹⁹ The inverse problem is solved using a least-squares solver, lsqr (MATLAB, Natick, Massachusetts). The solver was stopped after 200 iterations for all studies herein.

2.3.

Optimization Algorithm

The iterative correction of excitation pattern in AWFT modifies the spatial distribution of incident power by modulating the gray levels in the pattern. The optimization procedure therefore employs the measurements from the detectors spanning the area of the pattern on the surface of the model. As a part of pre-optimization calibration, the excitation pattern is measured using a white diffusing paper, and its spatial characteristics are determined. Next a transformation matrix mapping the image used to define the full field pattern to the measured full field is computed. This is done to account for spatial deformations in the projected pattern due to the projector optics.

Consider the animal model shown in Fig. 3(a) where the white border represents the excitation pattern. The optimization of the patterns described here is based on time-resolved data, and the process begins with the acquisition of temporal measurements for each of the $N$ patterns in the user-defined set with the laser operating at power ${\mathrm{Pow}}^i$ ($i=1,2,3\dots$) at the $i$th iteration. As the dynamic range limits the maximum photons recorded in each time-resolved curve, for each pattern $P_n$ ($n=1\dots N$), an excitation map ($E_n$) is constructed from the maximum measured counts (peak counts of the TPSF) at each detector. The underlying principle of the optimization algorithm considers each pixel in the pattern independently and corrects the grayscale level based on the excitation signal in transmittance. Therefore, at each pixel at $(x,y)$ in the $n$th pattern at the $i$th iteration the multiplicative update value $\alpha$ is computed and the corresponding pixel is updated using the formula given in Eq. (3). It should be noted that in Eq. (3), the calculation of $\alpha$ includes the weighted photon counts within a neighborhood of radius $r$, in order to incorporate the effect of diffusion on transmitted photon counts.

Fig. 3

(a) Optimization scheme; (b) algorithm flowchart for single pattern optimization.

In Eq. (3), $\chi$ is the maximum desired photon counts, $(x_k,y_k)$ is the position of the $k$th neighbor within a radius $r$, $d_k$ is Euclidean distance of the $k$th neighbor from $(x,y)$, and $K$ is the total number of detectors included in the neighborhood. The selection of detectors within distance $r$ incorporates the effect of photon diffusion when computing the update value for the pattern. Furthermore, we define $\beta$ as the maximum power amplification allowed when acquiring the TPSF. This ensures that the algorithm does not exceed the maximum permissible exposure limit of laser power during optimization. The iteration is completed after each pixel in the pattern has been updated. The update values obtained after this step for pattern $P_n$ (referred to as $A_n$) represent the relative change in laser power required at each position on the pattern to achieve desired excitation signal level. The experimental implementation of this updated pattern is achieved by normalizing the pattern to 1 and increasing the laser power in the subsequent step by a factor of the maximum update value, $A_n^{\max }$.

This process is repeated until $A_n^{\max }$ is less than or equal to 1 (indicating convergence with all detectors reaching $\chi$ counts). Furthermore, we include a second stopping criterion to ensure the termination of the optimization scheme when $A_n^{\max }$ converges to a constant value greater than 1, wherein the iterations are terminated when the maximum change in the new pattern (${\parallel P_n^{i+1}-P_n^i\parallel }_{\infty }$) is less than 10%. The above process is repeated for all patterns in the base pattern set, and the final set of patterns is referred to as the transmittance-optimized patterns. This optimization procedure is outlined in Fig. 3(b).

2.4.

Grayscale Limits During Pattern Optimization

As noted in the system description, the experimental implementation of the adaptive optimization of the patterns was done using the Pico-projector system. The nonlinearity of the projected gray levels as shown in Fig. 2(d) results in inaccurate projection of gray levels below 0.25. In order to mitigate the effect of this hardware characteristic, the optimization of patterns was limited to grayscale values above 0.25 with all pattern positions having values smaller than 0.25 remaining unchanged during optimization. Furthermore, all patterns obtained during the optimization process were acquired post-hoc and used in the Monte Carlo simulations in order to ensure the accuracy of the model during reconstruction. It should also be noted that the in silico validation of the algorithm applied no limits on the pattern.

2.5.

Neighborhood Radius

As stated above, the radius of the neighborhood incorporates the effects of diffusion when updating the gray level for a pixel on the pattern. The ideal radius selection for each pixel would typically include the height of the model at the pixel and the local optical properties, leading to subject dependent radius and thus necessitating a priori knowledge of the anatomy. In this work, we select the radius based on the distance of the edge of the model from the pattern to ensure that all photons exiting the model are accounted for in the optimization process. It should be noted that a radius of 2.5 mm was selected for the in silico studies and 4 mm for the experimental validation. Note that simulations with radii of 3 and 8 mm were also carried out leading to identical convergence rate as reported here when using 2.5 mm (results not shown).

2.6.

Batch Optimization of Patterns

The optimization procedure shown in Fig. 3(b) will result in each pattern converging to the optimality condition at a different rate based upon the region of excitation and will require different power amplification. For instance, photons from a pattern probing the thoracic cavity will be highly attenuated due to highly absorbing and scattering organs (e.g., heart and lungs), while the signal from another pattern in the same set probing the abdominal area will encounter a thicker volume with low absorbing tissue (e.g., stomach, intestines etc.). This leads to variation in the signal attenuation profiles recorded on the surface for each pattern and necessitates the individual optimization of each pattern in the set.

In this implementation of the above algorithm, all patterns were batch optimized to ensure experimental efficiency (acquisition time). In other words, the temporal measurements for all patterns in the set were acquired consecutively and followed by the optimization of each pattern in the set. Specifically, the optimization algorithm was implemented in two steps. First, the pattern update value was unbounded ($\beta =\infty$) during optimization and the optimal patterns were computed. As stated previously, each pattern will have been assigned different power amplification for the next iteration. In the second step, the minimum power amplification among all patterns in the set is put as the upper bound on the value of $\alpha$ [$\beta =\min (A_n^{\max })$] where, $n=1,\dots ,N$. This ensured that the laser power employed during the acquisition step will not saturate the ICCD for any pattern in the set. Moreover, the consistent source power across all patterns allowed a rapid acquisition of the complete pattern set with high experimental efficiency (typically less than 5 min for the optimization of a 36 patterns set). It should be noted that this algorithm can also be applied to CW and FD datatypes without loss of generality.

2.7.

In Silico Validation

The efficacy of the optimization algorithm for whole-body imaging of small animals was validated using a synthetic mouse model (Digimouse).²⁰ Figure 4 shows the model employed, where seven major organs in the animal torso were modeled to simulate the propagation of photons in small animals. The optical properties at 730 nm computed using values given in Refs. 21 and 22 were assigned to each organ and have been compiled in Table 1.

Fig. 4

(a) Synthetic mouse model used for in silico validation. Three $3\times 3\times 3{\mathrm{mm}}^3$ fluorescent inclusions are shown in red. The profile marked in black and the dotted line represents the region of excitation and the discrete detectors, respectively; (b) uniform bar illumination pattern; (c) Battle-Lemarie wavelet pattern; (d) sinusoidal illumination pattern.

Table 1

Optical properties of major organs in a small animal model at 730 nm.

Furthermore, two $3\times 3\times 3{\mathrm{mm}}^3$ fluorophore inclusions were positioned as shown in Fig. 4(a) to simulate the effect of optical properties and/or geometry on the fluorescence signal intensity. An area of illumination spanning the whole body (as shown in Fig. 4) was selected as the base pattern, and the excitation field was recorded at discrete detectors at 1 mm separation in transmittance mode. Three different illumination patterns were considered for evaluation. First, a uniform bar pattern,⁶ second a wavelet pattern (Battle-Lemarie, Scale 1),⁴ and last a sinusoidal pattern (spatial frequency of $0.125{\mathrm{mm}}^{-1}$) as shown in Fig. 4(b) to 4(d), respectively.^9,10 It is worth noting that the Monte Carlo simulations of time-resolved excitation and fluorescence signals also exhibit the Poisson noise characteristics of time-resolved imaging instrumentation. Thus no noise was explicitly added to the simulated measurements.

The optimization scheme was employed to achieve $\chi =4000$ counts at each detector, which is less than the maximum number of photon counts that can be detected by the ICCD camera used in the study (12-bits). After convergence, the fluorescence forward model was computed for the uniform and optimized pattern to determine the corresponding changes in fluorescence signal intensity.¹⁶

2.8.

In Vitro Validation

The improvement in performance in tomographic reconstruction was investigated using an agarose phantom mimicking a small animal model. The baseline pattern set employed in this tomographic study are binary bar-shaped patterns described in Ref. 6. The performance of this pattern set in a nonplanar geometry is tested using a truncated hemi-cylindrical phantom shown in Fig. 5 with an absorption coefficient of $0.1{\mathrm{cm}}^{-1}$ and reduced scattering coefficient of $13{\mathrm{cm}}^{-1}$. Furthermore, a cylindrical absorber with six times absorption ($0.6{\mathrm{cm}}^{-1}$) simulating absorbing organs in the small animal torso (e.g., liver) is placed along the central axis of the phantom. The geometry of the phantom coupled with the absorber results in a high dynamic range across the surface in transmittance.

Fig. 5

Tissue-mimicking phantom used for in vitro validation. The uniform full-field pattern (recorded during the experiment) projected on the phantom is shown on the xy-plane at $z=0$. The green tubes represent the fluorophores contained in glass capillaries located in the center of the phantom. The red tube shows the position of an absorbing perturbation occluding the central fluorescent inclusion. The black dotted lines represent the projection of the tubes on the excitation field.

Two 10 mm long capillary tubes (1.5 mm inner diameter) containing 14 pmol of Cardiogreen (Sigma Aldrich, MO) in 10 μL phosphate-buffered saline were placed at depths of 11 and 9 mm. Moreover, the deeper inclusion was placed below the absorber simulating the occlusion of fluorescent markers by absorbing organs in vivo.

The tomographic imaging session employed 36 bar patterns ($34\times 24\mathrm{mm}$) described previously. The temporal measurements at excitation wavelength (${\lambda }_{\mathrm{ex}}=780\mathrm{nm}$) were acquired at 40 ps interval spanning 2 ns (50 gates) using 300 ps gates ($\mathrm{MCP}\text{Voltage}=570\mathrm{V}$; $\text{integration time}=25\mathrm{ms}$). The measurements at excitation wavelength for 36 patterns were acquired in 3 min. The fluorescence measurements (at ${\lambda }_{\mathrm{em}}=832\mathrm{nm}$) were acquired for the optimized pattern set ($\mathrm{MCP}\text{Voltage}=570\mathrm{V}$; $\text{integration time}=500\mathrm{ms}$). The acquisition of measurements spanning 3 ns (75 time-gates) was completed in 18 min. The fluorescence measurements using the uniform pattern set was also acquired for comparison of reconstruction performance.

The tomographic reconstruction of the effective quantum yield of the two inclusions in three-dimensional (3-D) was performed using measurements at 231 point detectors uniformly sampling the surface of the phantom at 2 mm intervals. The time-gate datatype derived at each detector was defined by the number of photons measured relative to the peak photon count and was used to construct the time-gated measurement vector. The effective quantum yield was reconstructed by solving the inverse problem with measurements derived from multiple-gates spanning the rising portion of the TPSF. Furthermore, an x-ray CT image of the phantom was acquired following the optical imaging protocol using small animal CT imaging system (Scanco VivaCT40). The CT data was acquired at an isotropic resolution of 0.038 mm. The images were subsequently resampled to 1 mm isotropic resolution and segmented into three regions—agarose, absorber, and fluorophore—using Amira (VSG US, Burlington, Massachusetts). It should, however, be noted that beside the geometrical boundary information, no a priori information on the location of the fluorescence inclusion was enforced in the inverse problem. The 3-D image of the volume (obtained by combining the three regions specified above) was used to generate a voxel-based model geometry with 1 mm voxels. The 3-D model was registered with the optical measurements using fiducial markers positioned outside the phantom. The fully segmented model was subsequently used for validation of the tomographic performance of the pattern optimization scheme.

The dimensions of the 50% isovolume of the reconstructed quantum yield and the position of the centroid of the reconstructed inclusions are used as the resolution metrics to compare the adaptive pattern optimization scheme with the uniform pattern set.

3.1.

Transmittance Optimized Patterns

The optimization of photons transmitted through the synthetic animal model using the full-field pattern shown in Fig. 4 was completed in six iterations for the bar pattern, seven iterations for the wavelet pattern, and eight iterations for the sinusoidal pattern. Figure 6(a) and 6(b) shows the change in the two stopping criterion parameters over the optimization process and for all patterns the process was terminated due to less than 10% change in pattern characteristics. Figure 7(a), 7(c), and 7(e) shows the iterative optimization in the pattern spatial intensity distribution and the corresponding change in transmitted signal intensity on the surface for each of the three patterns. It should be noted that the gray levels on the pattern are adjusted to match the geometry and organ distribution in the animal model with the maximum intensity retained in the thickest portion of the body comprising of the highly absorbing organs like the liver and spleen.

Fig. 6

Stopping criteria for in silico validation of optimization scheme: (a) value of $A_n^{\max }$ (pattern update parameter) for each iteration; (b) change in pattern for each step during optimization. The dotted lines in (a) and (b) indicate the target value of respective parameter for stopping condition.

Fig. 7

Pattern optimization: (a) optimization of uniform bar pattern over five iterations; (b) excitation signal (log scale) for each step of optimization shown in (a); (c) optimization of wavelet pattern over six iterations; (d) excitation signal (log scale) for each step of optimization shown in (c); (e) optimization of sinusoidal pattern over seven iterations; (f) excitation signal (log scale) for each step of optimization shown in (e). All excitation signal maps have been normalized to 4000 counts to account for increase in laser power during practical implementation.

Figure 7(b), 7(d), and 7(f) shows the transmitted excitation field measured on the surface of the animal model. The excitation measurements obtained using the baseline patterns indicate a dynamic range of $\sim 4$ orders of magnitude across the area of illumination due to high fluence near the edge of the model and negligible transmission through the central portion. Following the characteristics of the optimized patterns described above, the optimization process leads to $\sim 2$ orders of magnitude improvement in transmitted light for the central portion of the torso reducing the effective dynamic range in the excitation field to less than 2 orders of magnitude.

3.2.

Fluorescence Signal Optimization

In order to test the hypothesis that the optimization of transmitted excitation field improves the fluorescence measurements, we considered the maximum photon fluence at each inclusion shown in Fig. 4. Figure 8 shows the relative change in the number of photons generated by each of the two inclusions during the optimization process. It is worth noting that the adaptive illumination intensity resulted in more than 600% increase in signal generated at inclusion 1 while retaining the signal level from inclusion 2 at the original level. A comparison of the fluorescence signal measured on the surface at an early time gate (10% of peak) and the maximum gate in Fig. 9(a) to 9(f) demonstrates the effect of nonoptimized measurements in fluorescence tomography with significantly lower signal from deeply embedded fluorescent sources like inclusion 1. Additionally, the optimization process resulted in an improvement in signal from inclusion 1 at both gates suggesting the advantages of this approach for fluorescence tomography. Specifically, the early gate measurements were improved with more than 1.5 order of magnitude increase in signal from inclusion 1 compared with $\sim 1$ order of magnitude increase at the maximum gate. Moreover, the signal from inclusion 2 is not modified for both gates considered and retains the high signal level obtained using the original patterns. The effect of signal improvement for tomographic imaging was quantified by measuring the number of detectors which had a useful signal level (above threshold of 200 counts) at both gates [c.f. Fig. 9(g) to 9(h)]. For all patterns considered, the optimization process resulted in $\sim 60\%$ increase in number of detectors at the early gate and $\sim 75\%$ increase for the maximum gate.

Fig. 8

Relative change in the number of photons generated by inclusion 1 (solid line) and inclusion 2 (dotted line) for the three patterns.

Fig. 9

Improvement in fluorescence signal at the early gate (10% of peak) and maximum gate measured on the surface after excitation optimization for (a) to (b) uniform bar pattern; (c) to (d) wavelet pattern; and (e) to (f) sinusoidal pattern. (g) to (h) Number of detectors on the surface having signal higher than threshold of 200 photons for the early gate and maximum gate, respectively.

An important characteristic of the measurements shown in Fig. 7 is the near equivalence of the signal obtained using all three types of patterns considered. Specifically, the improvement in signal resulting from pattern optimization (and the corresponding increase in useful detector readings) is the same for all patterns indicating the general applicability of this approach for fluorescence imaging. It suggests also that, due to the highly diffusing and heterogeneous nature of the preclinical models, different patterns will have similar performances when employed in transmittance.

3.3.

Tomographic Reconstruction Using Optimized Patterns

The experimental validation of AWFT using the small animal phantom was performed by batch optimization of the 36 patterns in the baseline set. The iterative scheme was terminated in 4 iterations with the power amplification factor for pattern 31 reaching a value of 0.65 as shown in Fig. 10(a) while the change in pattern characteristics [Fig. 10(b)] was significantly higher than stopping condition, indicating a successful convergence of the optimization scheme.

Fig. 10

Stopping criterion parameters for each of the 36 patterns over four iterations: (a) pattern update parameter $A_n^{\max }$; (b) change in patterns over iterations. The black arrows denote the pattern with lowest value of $A_n^{\max }$ for the respective iteration. This value is used as $\beta$ in the second step of the batch optimization. The optimization process was stopped when value of $A_n^{\max }$ for pattern 31 reached the target condition of $A_n^{\max }<1$.

All patterns were simultaneously optimized in $\sim 3\min$ (including $\sim 1\min$ of pattern setup time per iteration). The optimization procedure of the entire pattern set (4 iterations), including the acquisition of temporal measurements was completed in 24 min (note that the full time-gated data was acquired for each pattern at every iteration). Figure 11(a) and 11(b) compare the structure of a baseline half-field pattern and the corresponding transmittance-optimized grayscale patterns. It can be inferred that the optimization procedure removed the sections of pattern beneath the region with smaller thickness to correct for the dynamic range while retaining the high intensity across the central axis to compensate for the absorber. Figure 11(c) and 11(d) compare the excitation field acquired for the original and optimized patterns for the maximum gate. As expected, the original pattern results in a higher signal towards the edges of the phantom [c.f. Fig. 11(c)] demonstrating the role of model geometry as a source of high dynamic range. Conversely, the transmittance optimized field provides a high signal across the area of the pattern [c.f. Fig. 11(d)] with maximum of 4000 photons across the area of the pattern representing the successful optimization of the transmitted signal.

Fig. 11

In vitro demonstration of uniform bar pattern optimization: (a) baseline bar pattern illuminating the phantom shown in Fig. 5; (b) excitation optimized pattern after four iterations. The blue profile represents the detector area shown in the measurements in (c) to (j), the orange outline identifies the full-field of illumination, and the green lines are projections of the fluorescent inclusions; (c) to (d) comparison of excitation signal at maximum gate using the baseline and optimized pattern; (e) to (f) comparison of fluorescence signal at an early gate (10% of maximum) using the baseline and optimized pattern; (g) to (h) comparison of fluorescence signal for the early gate at 50% of maximum using the baseline and optimized pattern; (i) to (j) comparison of fluorescence signal at maximum gate using the baseline and optimized pattern. All measurements shown here represent raw photon counts with a threshold at 200 counts. The measurements for the remaining 35 patterns are shown in Video 1 (MOV, 8.18 MB) [URL: http://dx.doi.org/10.1117/1.JBO.18.3.036006.1].

A similar comparison of the fluorescence field shown in Fig. 11(e) to 11(j) demonstrates that the optimization of the excitation field improves the fluorescence signal from inclusion 1. This demonstrates an improved contrast and higher signal over a larger number of detectors when using AWFT. The results of optimization of all the patterns used in the study are shown in Video 1.

3.4.

Improved Information Content

The impact of the improvement in the overall signal level demonstrated in the measurements on the quality of the tomographic information recorded is quantified using a two-step process. First the percentage increase in number of source-detector pairs at different gates along the temporal measurement which are above the noise threshold level is used to quantify an increase in the number of usable measurements The noise threshold for this system was empirically determined to be 200 photon counts. Figure 12(a) and 12(b) shows $\sim 200\%$ increase in the number of s-d pairs for excitation signal and $\sim 80\%$ in the fluorescence signal at the early gates, with a smaller yet significant improvement for the maximum gate ($\sim 25\%$ for fluorescence field, 15% for excitation field). The improvement in the time-gated signal at the early gates and the late gates establishes the advantages of this imaging technique for high-resolution optical reconstruction^23–25 and lifetime multiplexing applications, respectively.¹⁴ It is also worth noting that the transmittance optimized patterns result in a higher increase for the late gates when compared with the early gates as the optimization procedure corrects for the effects of the highly absorbing local perturbation in the experimental phantom which reduces signal at the late gates.

Fig. 12

(a) to (b) Increase in number of source-detector pairs with signal above 200 counts after pattern optimization at excitation and emission wavelengths, respectively; (c) increase in nonredundant information with increase in s-d pairs shown in (a) and (b).

Next, the analysis of tomographic information content in the measurements is extended to the quantification of increase in nonredundant information using a model based analysis. Following the work by Freiberger et al., we define the orthogonality of the different rows in the fluorescence weight matrices (obtained by discretizing the fluorescence sensitivity function $J^{\mathrm{fluo}}$) as a measure of the information content.²⁶

For a set of $N$ patterns, let $M$ be the total number of measurements (s-d pairs) above the set threshold. Consider a set of measurements $S_n$ for the $n$th pattern. Then the nonredundant information $h_n$ for the $n$th pattern is defined as the average orthogonality of each measurement in $S_n$ with respect to the remaining measurements ($M\setminus S_n$). $h_n$ is computed using Eq. (4):

Here, $j_m$ and $j_n$ are the $m$th and $n$th rows from $J^{\mathrm{fluo}}$. The average nonredundant information content for the entire pattern set is therefore given by

The average information content value of 1 represents maximal nonredundancy in the measurement suggesting maximum information content in the measurements. Figure 12(c) shows the comparison of average nonredundant information content for the uniform and the transmittance-optimized pattern set at individual time-gates along the TPSF.

It is observed that the optimization of transmittance excitation field leads to $\sim 6\%$ increase in the information content at the early gates and $\sim 16\%$ increase for gates along the decaying section of the TPSF. This closely follows the trend in increase in number of s-d pairs observed in Fig. 12(a).

3.5.

Reconstruction Results

In order to determine the impact of increased information content on the quality of reconstructions in the image space, the effective quantum yield was reconstructed using three time gates spanning the rising portion of the TPSF (early gates at 20%, 50%, and the maximum gate). These gates were selected on the rising part of the TPSF as it correlates to resolution improvement in the reconstruction, whereas late gates are critical in lifetime multiplexed studies, which we do not consider herein.^16,27 Figure 13(a) and 13(b) shows the 50% iso-volumes of the effective quantum yield obtained by solving the inverse problem by least-squares minimization. The lower number of source-detector pairs available when using the uniform pattern set (more than 36% reduction) leads to reduced tomographic information in the measurements as discussed previously. This results in poor discrimination of the two objects in the reconstruction and their incorrect localization [c.f. Fig. 13(a)]. Specifically, both objects are reconstructed 5 mm above the expected position along the $z$-axis. Moreover, the lower resolution of the reconstructed volume is evidenced by the connected reconstructed objects. On the contrary, the reconstructed volume obtained using the transmittance optimized patterns allows us to discriminate and accurately localize the two inclusions as shown in Fig. 13(b). The resolution metrics comparing the reconstructed 50% isovolumes with the ground-truth values obtained from the x-ray CT are compiled in Table 2. It is worth noting that the two inclusions are accurately localized with $\sim 1\mathrm{mm}$ error.

Fig. 13

Reconstructed effective quantum yield using three gates spanning the rising portion of the TPSF at 10%, 50%, and 100% of maximum value using (a) uniform patterns and (b) transmittance optimized patterns. 50% isovolume of reconstructed effective quantum yield when using uniform patterns and transmittance optimized patterns, respectively. The 50% iso-contours across $x=28\mathrm{mm}$ and $z=12\mathrm{mm}$ are projected along the corresponding planes.

Table 2

Comparison of resolution metrics of reconstructed quantum yield with ground truth values derived from x-ray CT.

Also, the separation between the two objects was found to be 9 mm (expected—8.3 mm) resulting in less than 1 mm localization error. Furthermore, the reconstructed dimension of the inclusions had $\sim 1$ and $\sim 2\mathrm{mm}$ error along the $x$-axis and $y$-axis, respectively. The resolution of inclusion 2 along the $z$-axis also had less than 2 mm error; however, inclusion 1 had significantly lower z-resolution due to the limited projection angles in the planar imaging configuration. Note that the imaging volume was discretized with $1{\mathrm{mm}}^3$ voxels, and hence the reported errors are on the same scale as the smallest element of volume (voxel).