Hyperspectral imaging devices are common in remote sensing reconnaissance as well as biomedical applications.¹ Typically, these instruments are point or slit imagers and must scan the scene temporally to build up a two-dimensional image, or use bandpass filters to scan the scene spectrally. More specific types, including Hadamard encoding slit spectrometers,² liquid crystal and acousto-optic tunable filters,³ Fourier transform spectrometers,⁴ spectral-temporal scanning,⁵ and more recently volume holographic methods⁶; all scan to build up the spatial-spectral cube of information. In this article the computed tomographic imaging spectrometer (CTIS) is used. It captures both spatial and spectral information in a single frame. Its acquisition time is only constrained by the exposure time of the camera (typically on the order of milliseconds) and required signal-to-noise ratio. There are no moving parts or narrow-band filters, and nearly all collected light (about 70%) is passed to the detector at all times.

Although CTIS has been used in the past for stellar classification,⁷ fluorescence microscopy,⁸ and geology measurements,⁹ this report is on the groundbreaking use of CTIS for functional mapping of the human retina. This application is ideal for exploring the true potential for CTIS technology since the eye is constantly moving, requiring snapshot camera operation. Hyperspectral imaging offers functional rather than structural information, such as optical coherence tomography, and provides complementary information for the clinician. Our results show how hemoglobin spectral signatures provide both qualitative and quantitative oxygen saturation maps. This data offers the possibility of monitoring retinal ischemia from either systemic diseases such as diabetes or from localized retinal arterial and vascular occlusions—the leading causes of untreatable blindness.

Retinal hyperspectral imaging offers the ability to capture in vivo metabolic and physiologic information using chromophore spectra to classify tissues and quantify cellular metabolites.^{10, 11, 12} All of the software tools developed for remote sensing as well as signal and image processing can be applied to the hyperspectral data obtained. Clinical trials are underway at the Doheny Eye Institute, Keck School of Medicine at the University of Southern California to establish a significant patient database.

Current ophthalmic hyperspectral imaging systems^{13, 14} typically require long exposure times $\gtrsim 5\mathrm{s}$ , and therefore are plagued with motion artifacts from rapid eye movements (saccades). Both major and minor saccades occur frequently enough to significantly affect both band-sequential and line-scanning instruments. These data then require significant and laborious postprocessing and registration. Image registration compromises either spatial or spectral resolution, depending on the type of system used. One does not know which wavelengths are affected, cannot be registered, or simply miss the relevant spatial features used for registration in a clinical setting until after the patient has left. For retinal imaging especially, postacquisition registration has artifacts from both changes in the lens-eye geometry and tissue movement from arterial pulses between successive images. Furthermore, patients cannot fixate long enough or accept the high illumination levels required to keep exposure times short enough. The true snapshot capability of the CTIS imager avoids all these problems and provides an easy, patient-friendly way to obtain retinal spectral data. Guidelines of the Association for Research in Vision and Ophthalmology for human investigation were followed, and our Human Investigation Committee approved the study protocol.

CTIS captures the spatial and spectral information of the retina by imaging the scene through a two-dimensional grating disperser as illustrated in Fig. 1 . This produces multiple, spectrally dispersed images of the retina that are recorded by a focal plane array (FPA). From the captured intensity pattern, computed-tomography algorithms are used to reconstruct the scene into a “cube” of spatial ( $x$ and $y$ ) and spectral (wavelength) information. The image cube in wavelength space is then reconstructed from a single image. Note that the basic CTIS design uses just two lenses and a focal plane detector. The CTIS instrument concept originated in Japan^{15, 16} and Russia¹⁷ and has been significantly advanced to maturity by a group at the Jet Propulsion Laboratory (JPL) in Pasadena^{18, 19, 20} and one at the University of Arizona. ^{21, 22, 23, 24} The physics and the mathematics describing CTIS imaging and data cube reconstruction are similar to that of positron-emission tomography, hence the name “computed tomography” imaging spectrometer.

The specific type of diffractive optical element used is called a computer-generated hologram (CGH). It operates at multiple orders, creating an array of images on the sensor. Development of techniques for fabricating the grating with electron-beam lithography has been the main driver in the development of this instrument.²⁵ The CGHs are composed of cells of square pixels (typically $8\times 8$ pixels/cell) that are arrayed to form a 2-D grating. The pixel depths are designed using a CTIS-specific iterative optimization algorithm such that the dispersion order intensities are nearly equal after distribution across the focal plane array. This maximizes the overall system signal-to-noise ratio by avoiding saturation in a single diffraction order. An analog direct-write electron-beam lithography technique was developed to accurately fabricate a CGH. This technique requires only a single exposure followed by iterative development. It allows the production of fine-featured analog surface-relief profiles in thin films of E-beam resist with a depth error less than 5%.

Figure 2 shows how the various orders are dispersed across the focal plane as well as spectra derived for an image scene. It is important to note that each image is not simply composed of single wavelengths; spatial and spectral information from each object pixel is multiplexed over the entire array. Figure 2 shows how the spectrum of a single pixel is distributed by the diffractive disperser. Hence, a single image contains all the information required to reconstruct the spectral image cube. Note that there is a zero-order image for focusing—a difficult task for many spectral imagers.

A calibration matrix was necessary to perform the retinal reconstruction; it was obtained by measuring the location on the image plane of pixels from the object plane at different wavelength with a movable fiber optic coupled to a monochromator. The efficiencies, spot centroids, and point spread functions of all the diffraction orders are measured at all wavelengths in the operating band. This data is used to build a “system transfer matrix” of field-stop-to-detector connection coefficients. Once the system transfer matrix was known, the images of the retina were reconstructed to obtain their spectra. Following Descour ²¹ the iterative expectation-maximization (EM) reconstruction algorithm was used. This algorithm was developed for reconstructing positron-emission tomography medical images,^{26, 27} but other emission tomography algorithms are applicable as well. This estimation approximates the spatial and spectral information using voxel basis functions, which define the solution space for accurate reconstructions. Although there are various forms of iterative algorithms, EM is used due to its robust solutions and proven reliability. The reconstruction starts with an initial condition for the scene which results in a predicted detector image. Corrections are then made to the scene based on back projection of the ratio of the measured detector image and predicted detector images. This process is repeated iteratively until the predicted detector image matches the measured detector image to a desired error tolerance. Typically 10 to 20 iterations were required for good retinal reconstructions. National Institute of Standards and Technology (NIST) traceable reflectance color targets were used to verify the accuracy of the CTIS reconstruction process. Comparison of the CTIS spectral reflectance shows agreement to the NIST standards within a few percent over 450 to $750\mathrm{nm}$ .

The latest version of CTIS reported here replaces the usual monochrome camera with a digital color CCD (Qimaging $\copyright$ Retiga 4000R) as the detector. Imposing spatial or spectral modulation on the image scene significantly improves the spectral reconstructions.²⁰ The Bayer color filters on the camera pixels that provide each RGB plane act to provide spectral modulation. It was found that using the color image as the initial condition for the reconstructions significantly reduces reconstruction time by decreasing the number of iterations required. This is due to the fact that the color image already contains some spectral information as the initial “guess.” Figure 3 shows spectral calibration data for the CTIS taken with the color camera system. Results are shown in Fig. 3a for NIST traceable color targets in which the CTIS results match the standard to a few percent across the spectral range. A second set of standards (rare earth), shown in Fig. 3b, allows verification with high-frequency spectral components.

A Zeiss $\mathrm{ff}450+$ fundus camera was used as the objective assembly to image the retina of the human test subjects. An example of a raw image from the imager is shown in Fig. 4 . A filter matched to the camera flash bulb spectra optimized the dynamic range of the silicon detector from 500 to $600\mathrm{nm}$ . The fundus objective had a large 24 $\times 36\mathrm{mm}$ image plane. Demagnification $(1∕7.2\times )$ was needed in order to compensate for this large size and was accomplished using a $2\times$ tele-extender along with a 180-mm commercial-off-the-shelf (COTS) lens assembly for collimation. The field stop was set at $11\times 11\mathrm{mm}$ and combined with a 50-mm COTS re-imaging lens reduced the field of view to $18\mathrm{deg}$ with a Fundus camera setting of $50\mathrm{deg}$ , with minimal vignetting. This field of view converts to an image plane at the detector across $208\times 208$ pixels. The Fundus camera has three different settings for the field of view and each maps a higher resolution image across the pixel space of the zeroth order CTIS. The Fundus camera settings of 50, 30, and $20\mathrm{deg}$ correspond to 18, 10.8, and $7.2\mathrm{deg}$ , respectively, with our optics and field stop. Small vessels ( $<75\mathrm{\mu }\mathrm{m}$ in diameter) are best visualized at the 7.2-deg setting while larger vessels ( $>150\mathrm{\mu }\mathrm{m}$ in diameter) can be easily imaged at the 18-deg setting. Although it has a smaller field of view, the $18\text{-}\mathrm{deg}$ setting still images both the optic disk and macula together, allowing the clinician to easily compare with typical retinal images. Due to the large numerical aperture of the fundus objective, a field lens is implemented to counterbalance the spatial variance due to the vignetted high-angle rays. The field lens was placed just outside the field stop so that any dust on the lens would be out of focus, making a negligible change to the demagnification.

This article reports on the co-registered image and associated oximetry unmixed reflectance spectra for different parts of the human retina. Incident illumination is spectrally modified due to scattering and absorption contributions within the eye as approximately described by Lambert-Beer’s law. Previously measured extinction coefficients are used for oxy-hemoglobin and deoxy-hemoglobin²⁸ to calculate absorption and scattering coefficients due to the turbid media of the layered retinal structure. For this calculation both scattering and absorption due to the aqueous and vitreous humor are neglected. Our model assumes all scattering and absorption are only due to the vascular surface and underlying RPE layers of the retina. To a first approximation, the vascular structure’s absorption is given by the linear model, which simply multiplies the absorption and scattering, respectively, per unit length and per unit concentration with the absorption and scattering concentration through some retinal tissue depth. The goal is to find the percentage of oxygen saturation based on a weighted linear model. This model is constructed based on previous literature using typical molar concentration (150 g/l) and molecular weight (64,500 g/mole) for healthy individuals. In the generalized linear oximetry model the natural logarithm is used with the ratio of the signal versus background response for capillary vessels. This is purely based on Lambert-Beer’s law where each response is added in the exponent of the exponential function.

Oxygen saturation algorithms were used based on previous multiwavelength oximetry models from slit scanning spectrometers. Specifically, a three-wavelength oximetry model was used to localize the saturation values spatially for each scene. A totally unsupervised approach was taken in order to show how saturation values, localized at positions known to be a vessel, can be found for each dataset. The model is described elsewhere.^{29, 30, 31} The output of this function allows a spatial-saturation map. The local oxygen saturation is related to the health of the particular area of an individual’s retina.

Reflectance from wavelengths shorter than $500\mathrm{nm}$ is adversely affected by light scattering in ocular media and erythrocytes.³² For wavelengths longer than $585\mathrm{nm}$ , a drastic decrease in hemoglobin extinction coefficients causes light scattering and absorption to dominate in deeper retinal tissues and intervening ocular media. In order to properly perform three-wavelength oximetry, the following optimization procedure was employed. Three wavelengths were randomly chosen from the imaged wavelength interval of 454 to $646\mathrm{nm}$ . Subsequently these three wavelengths were modified according to a multivariate optimization algorithm with a spectral resolution of $1\mathrm{nm}$ (the imaged hyperspectral data cube is in increments of $4\mathrm{nm}$ between hyperspectral planes). Linear interpolation of the hyperspectral planes from the three wavelengths provided hyperspectral data at in-between wavelengths for each iteration step of the optimization process. Each hyperspectral plane used for the interpolation was normalized and white-corrected (divided by the measured spectra of the flashbulb). The respective normalization factors were subject to the same multivariate optimization process used for the wavelengths, resulting in a $3+6$ parameter optimization. The fitness function for this optimization process was the necessary condition that in order to properly perform three-wavelength oximetry calculations, the oximetry value for each mapped pixel in the image had to be within the physically meaningful interval [0,1]. Upon convergence of the optimization procedure, the optimized set of the three wavelengths together with the six normalization factors for the hyperspectral planes necessary for the interpolation were determined and the oximetry values for the entire image were calculated.

To suppress background scatter a Monte Carlo optimization routine was used in combination with the three-wavelength model already described. Starting again with the oxy-hemoglobin and deoxy-hemoglobin extinction coefficients and three wavelengths, a stochastic search of the local molar concentration space is given by

where $p$ is the probability of an event occurring at the correct value and $y$ is the perturbing search vector spanning the region of molar concentration of interest. The output yields the oxygen saturation map given by $\phi$ , since Eq. 1 is performed upon the entire cube simultaneously. Furthermore, the angular dependency of the bulb illumination is taken out using a Gaussian illumination approximation. The Gaussian function is adjusted based on the tendency of the bulb illumination to vary spatially for each image around the edges. This is due to operator pointing motion and the complex dynamic function of the human lens-eye system.

Two volunteers were both considered healthy but with an age difference of 30 years. Figure 5 shows a 12-deg image of the optic disk for volunteers 1 and 2, respectively. The results show a clear distinction between veins, arteries, and background. Some mixture of saturation value distribution is expected as has been found in past line-scanning instruments.^{32, 33} Regions within vessel capillaries agree well with the 30 to 35% oxygen saturation difference expected for healthy veins and arteries. ^{32, 33, 34, 35} The saturation for most of the background spatial locations in between the capillary regions shows a tendency to be within the 90 to 100% regime. This is consistent with the subjects being healthy.

The method for acquisition of retinal spectral images demonstrated in this paper has wide applications in any biomedical imaging where patient motion, arterial pulses, or patient comfort must be considered, such as bedside care. One implementation of a snapshot CTIS is inexpensive, is hand held, and does not require a computer to take data. This opens the door for the routine use of spectral imaging in areas outside of biomedicine, for example, forensics, art conservation, color, and field use in archaeology and agriculture.