1.

Introduction

Hyperspectral imaging (HSI) is a powerful tool widely used in a diverse range of applications, such as remote sensing,^1,2 food testing,³ pharmaceutical testing,⁴ forensic testing,⁵ and disease diagnosis,⁶ as it can acquire spatially resolved spectral information. In HSI, hyperspectral images are acquired by obtaining the light intensities in the spatiospectral three-dimensional (3D) $(x,y,\lambda )$ data cube via various methods, such as whisk broom sensing (point-scanning to record spectral information of each pixel with a diffraction grating and a linear sensor array),⁷ push broom sensing [line-scanning to record a one-dimensional (1D) spectral image with a two-dimensional (2D) camera],⁸ staring sensing (spectral filtering of incoming light to record the images of each spectral band with a 2D camera),⁹ snapshot spectral imaging (a nonscanning method to record the remapped 3D data cube with a 2D camera),¹⁰ and compressed sensing (CS)-based coded aperture snapshot spectral imaging (CASSI)^11,12 as shown in Figs. 1(a)–1(c). Among them, the staring method is most widely used in modern HSI applications because it is compatible with on-the-fly spectral imaging, as its image acquisition is easily performed with a 2D camera after spectral filtering by a tunable filter.

Fig. 1

Data acquisition of a 3D data cube in HSI. (a) Various scanning methods to obtain a 3D data cube in HSI. The pixels measured during a detector integration period are depicted for each scanning method. (b) Snapshot spectral imaging. (c) CASSI. (d) Data acquisition of CS-powered HSI. The data points in the 3D data cube are partially sampled and then processed to reconstruct the complete data set.

However, the staring method often suffers from low light transmittance and therefore limited signal-to-noise ratio (SNR) especially in the infrared (IR) region due to the lack of sensitive image sensors. Also the wide-field illumination of the staring method is an impediment to nonlinear optical measurements, which require high optical intensity at the sample to produce a detectable signal. Similarly, CASSI cannot be applied to nonlinear HSI due to its wide-field illumination. To overcome these issues, it has been demonstrated that scanning in the spatiotemporal space ($x,y,\tau$), or time-domain HSI, provides a higher SNR by virtue of its multiplex advantage (also known as Fellgett’s advantage) compared to direct measurements in the spectral domain. Unfortunately, the major drawback of time-domain HSI is its low image acquisition rate, primarily due to slow optical scanners for multidimensional spectral image acquisition. Even with high-speed optical scanning using a resonant scanner, only up to 2.4 frames per second (fps) with a pixel resolution of $100\mathrm{pixels}\times 100\mathrm{pixels}$ has been realized.¹³ To boost the image acquisition speed, one possible approach is to combine time-domain HSI with CS^14,15 [Fig. 1(b)], which provides a technique to recover a complete data set from fewer measurements than what is required by the Nyquist–Shannon sampling theorem based on two assumptions: (i) the signal from the target object is sparse in a suitable domain and (ii) the measurement process is incoherent; the sampling intervals along the $\tau$ axis are randomized in the case of Fourier-transformation-based measurements like time-domain HSI. Although CS has been demonstrated for time-domain HSI, the methods implemented in the previous demonstrations are not applicable to high-speed hyperspectral image acquisition because the requisite sampling randomization was implemented during postmeasurement processing¹⁶ or was paired with a lower camera frame rate.¹⁷

In this paper, we report a CS-powered method to achieve high-speed time-domain HSI. Specifically, we used broadband (200 to $1600{\mathrm{cm}}^{-1}$) Fourier-transform coherent anti-Stokes Raman scattering (FT-CARS), a type of nonlinear Raman process,^13,18,19 as a platform to evaluate this method. We sampled time-domain interferograms sparsely by utilizing the 3D Lissajous scanning method²⁰ to boost the hyperspectral image acquisition rate and then reconstructed hyperspectral images by solving an inverse problem regularized with spatiospectral total variation (SSTV).²¹ We demonstrated hyperspectral image acquisition at a frame rate of 25 fps with a pixel resolution of $100\mathrm{pixels}\times 100\mathrm{pixels}$ under experimentally feasible conditions, which is 10 times higher than that of the previous record.¹³ We compared the original sinusoidal Lissajous scanning and the triangular Lissajous scanning for CS-based high-speed time-domain HSI in terms of noise robustness and achievable compression ratios, showing that Raman hyperspectral images were reconstructed with reasonable quality in both methods. Because the present scheme is applicable to many realistic Raman or HSI methods and these simulations are based on existing experimental demonstrations that used commercially available equipment, our method offers a powerful platform for the video-rate visualization of rapid chemical processes in various systems, such as inhomogeneous chemical reactions and living cells.

2.

Methods

Figure 2(a) schematically shows the optical setup of our simulated CS-based FT-CARS imaging system. The light source is a femtosecond laser with a repetition rate of 80 MHz. In FT-CARS measurements, molecular vibrations are coherently excited by a pump pulse and interrogated by a time-delayed probe pulse.²² Pump–probe pulse pairs are generated by a Michaelson interferometer, where the optical path length of one of the arms is modulated by a resonant scanner. The incident beams are then spatially scanned by the two galvanometric scanners for interrogating different spatial points of the sample. A long-pass and short-pass filter pair, placed before and after the sample, respectively, isolates the blueshifted CARS signals. Just as in previously demonstrated FT-CARS measurements,²² the CARS signal is acquired as a function of the pump–probe delay, which is Fourier-transformed to obtain the Raman spectrum as a function of Raman shift or wavenumber.

Fig. 2

CS-based FT-CARS imaging. (a) Schematic of the simulated experimental setup. AL, achromatic lens; APD, avalanche photodiode; CoM, concave mirror; L, lens; LPF, long-pass filter; P, polarizer; PBS, polarizing beam splitter; RS, resonant scanner; S, sample; Sc, scanner; SPF, short-pass filter; and $\lambda /4$, quarter-wave plate. (b), (c) Conceptual illustration of sinusoidal (b) and triangular (c) Lissajous scanning trajectories in 3D space spanned by positions ($x,y$) and optical delay $\tau$. Typical sampling patterns at one spatial point are shown in the caption bubbles. (d), (e) Distributions of the number of sampled points with sinusoidal scanning (d) and triangular scanning (e) for$xy$ scanning. The upper bound is set at 300 for both images.

In the present scheme, the ($x,y,\tau$) data cube was scanned in a 3D Lissajous trajectory²⁰ such that the sampling interval along the $\tau$ axis became pseudorandom at each spatial point [Fig. 2(b)]. To this end, there are several options for scanning functions to generate the 3D Lissajous scanning trajectory, such as triangular and sinusoidal functions, especially for the $xy$ scanning. Each of them has advantages and disadvantages in the context of CS recovery. As discussed in the previous research,²⁰ triangular functions are more advantageous than sinusoidal functions in terms of the uniformity of the number of sampled points in the $xy$ plane, as shown in Fig. 2(c). On the other hand, sinusoidal functions realize shorter pixel dwell time than triangular functions because sinusoidally oscillating scanners, such as a resonant scanner and a microelectrochemical system scanner, can operate at higher frequencies. Because shorter pixel dwell time leads to a more randomized sampling pattern along the $\tau$ axis at each pixel, as shown in the caption of Fig. 1(b), more reliable CS recovery is expected.¹⁴ We compared the performances of triangular and sinusoidal scanning functions in the following numerical simulations.

The procedure to reconstruct the hyperspectral image from the sparsely sampled data is introduced as follows. First, we define $N_x,N_y$, and $N_b$ as the dimensions of the hyperspectral image, $N=N_x{N}_y{N}_b$ as the full size of data points in the spectral domain, and $M$ as the total number of sampled points in the time domain. Let $\mathbf{u}=[{\mathbf{u}}_{\mathrm{1,1}}^{\mathrm{T}},\cdots ,{\mathbf{u}}_{N_x,N_y}^{\mathrm{T}}{]}^{\mathrm{T}}\in {\mathbb{C}}^N$ be the hyperspectral image made by stacking the spectrum of each spatial position on top of one another and $\mathbf{f}=[{\mathbf{f}}_{\mathrm{1,1}}^{\mathrm{T}},\cdots ,{\mathbf{f}}_{N_x,N_y}^{\mathrm{T}}{]}^{\mathrm{T}}\in {\mathbb{R}}^M$ be a stack of coarsely measured time-domain interferograms at each spatial position. A hyperspectral image $\mathbf{u}$ was reconstructed from the coarse time-domain signals $\mathbf{f}$ obtained with Lissajous scanning by solving the following optimization problem with a regularization term called the anisotropic SSTV:²¹

Fig. 3

Selection of the hyperparameter values. RMSE values as a function of $\lambda$ are shown. The measurement duration is 40 ms. Dashed lines represent the optimal values of $\lambda$.

3.

Results

As a proof-of-concept demonstration, we simulated hyperspectral image reconstruction with a simple artificial hyperspectral image, shown in Fig. 4(a). We first determined the sinusoidal Lissajous scanning trajectory by three sinusoidal functions with scanning frequencies of $f_x=\mathrm{12,705}\mathrm{Hz}$, $f_y=\mathrm{12,505}\mathrm{Hz}$, and $f_{\tau }=\mathrm{12,520}\mathrm{Hz}$, and the triangular Lissajous trajectory by two triangular functions with $f_x=1005\mathrm{Hz}$ and $f_y=1060\mathrm{Hz}$ and a sinusoidal function with $f_{\tau }=\mathrm{12,520}\mathrm{Hz}$, with a measurement duration of 40 ms so that the trajectory covered a large number of spatial positions and never repeated itself and the restricted isometry constant of the measurement matrix was small. This measurement duration was 10 times shorter than the corresponding raster scan measurement. We assumed that one kind of chemical with five Raman vibrational modes was spread in a droplet-like concentration distribution shown in Fig. 4(a). Then we sampled corresponding time-domain interferograms in the generated Lissajous trajectory. Figures 4(b)–4(d) show the results of the hyperspectral image reconstruction with and without CS. Simple pixel-wise NDFT of undersampled data failed to recover the spectra due to aliasing, as shown in Fig. 4(b), whereas five input peaks were recovered in the spectral image reconstructed via CS [Figs. 4(c) and 4(d)].

Fig. 4

Numerical simulation of recovering undersampled data. (a) Assumed spectrum of a chemical and its concentration map. (b) Typical spectrum and intensity map at $1000{\mathrm{cm}}^{-1}$ of sparsely sampled data obtained by pixel-wise NDFT without CS. (c), (d) Typical spectra and intensity maps at $1000{\mathrm{cm}}^{-1}$ of the hyperspectral images reconstructed via CS from time-domain signals coarsely sampled with (c) sinusoidal Lissajous scanning (d) and triangular Lissajous scanning. The arrows in (b)–(d) indicate the points whose spectra are displayed.

To quantitatively evaluate the reconstruction performance, we reconstructed spectral images from the sparsely sampled time-domain interferograms with different noise levels, different compression ratios, and different scanning schemes (faster sinusoidal and slower triangular). The reconstruction performances were evaluated by the RMSE between the ground-truth time-domain interferograms and the pixel-wise inverse Fourier transform of the estimated spectral image. First, we investigated the measurement-duration dependence of the reconstruction performance in the sinusoidal and triangular scanning schemes [Fig. 5(a)]. To simulate the reconstruction with different measurement durations, the measurement duration was simply changed with the same scanning frequencies in the sinusoidal case, whereas the scanning frequencies were modified so that the number of samples at each pixel became more uniform in the triangular case. These simulation conditions are consistent with the realistic experimental situation, where the frequency of the resonant scanner cannot be significantly changed, while that of the galvanometric scanner can flexibly be changed within a certain frequency range. As seen in Fig. 5(a), the sinusoidal scanning provides lower RMSE values with different measurement durations than that of the triangular scanning, presumably because the sinusoidal Lissajous trajectory has more randomized sampling intervals along the $\tau$ axis at each pixel. Figures 5(b) and 5(c) show typical intensity maps at $1000{\mathrm{cm}}^{-1}$ and typical spectra for measurement durations of 30 ms with sinusoidal Lissajous scanning and triangular Lissajous scanning, respectively. Next, we studied the noise dependence of the reconstruction performance by adding Gaussian noise to the input time-domain interferograms [Fig. 5(d)]. The noise amplitude was varied from zero to one with respect to the amplitude of the $1000{\mathrm{cm}}^{-1}$ mode of the sample. Although there are slight differences in the performances of these two scanning functions, both methods work well, as shown in Figs. 5(b), 5(c), 5(e), and 5(f).

Fig. 5

Performance of CS-powered time-domain HSI with different scanning functions. (a) Performance of CS-powered time-domain HSI with different compression ratios. The top axis represents the measurement time normalized by $T_{\mathrm{raster}}$, the measurement time in the raster-scan method. (b), (c) Intensity maps at $1000{\mathrm{cm}}^{-1}$ and typical spectra with measurement durations of 30 ms with sinusoidal Lissajous scanning and triangular Lissajous scanning, respectively. (d) Performance of CS-powered time-domain HSI with different values of noise amplitudes and a measurement duration of 40 ms. (e), (f) Typical intensity maps at $1000{\mathrm{cm}}^{-1}$ and typical spectra with a noise amplitude of 0.5 with sinusoidal Lissajous scanning and triangular Lissajous scanning, respectively. Arrows in (b), (c), (e), and (f) indicate the points whose spectra are displayed.

To further evaluate the performance of our method in practical applications, we also demonstrated the reconstruction of a biological hyperspectral image in a semiexperimental manner using a Raman spectral image of a Euglena gracilis cell taken by a spontaneous Raman microscope (Renishaw inVia). A $31\mathrm{pixel}\times 33\mathrm{pixel}$ Raman hyperspectral image of a wild-type E. gracilis cell was obtained using 74 mW of 532 nm excitation, a $50\times$ objective, and a per-pixel measurement duration of 300 ms. The image was denoised by truncated singular value decomposition using the top-10 components, and then spectral background was removed using the rolling-ball filter method.²⁵ We generated the time-domain interferograms by inverse Fourier transforming the spectrum at each pixel after resizing it to $100\mathrm{pixels}\times 100\mathrm{pixels}$ and applying a window function. The sinusoidal Lissajous sampling proposed above was used with a measurement duration of 60 ms. Figures 6(a)–6(j) show typical intensity maps and spectra of the reconstructed spectral image and the Fourier transform of the generated full-size time-domain interferograms. Although the spatial resolution is degraded due to the use of total variation regularization, the same molecular distributions between the reconstructed [Figs. 6(a), 6(c), and 6(e)] and ground-truth [Figs. 6(b), 6(d), and 6(f)] spectral images and peak distributions among their corresponding spectra [Figs. 6(g)–6(i)] are evident. Figure 6(k) shows the performance of CS-based time-domain HSI for different measurement durations.

Fig. 6

Recovery of a spectral image of an E. gracilis cell by CS. (a), (c), (e) Intensity maps at 1000, 1300, and $1440{\mathrm{cm}}^{-1}$ of hyperspectral images reconstructed via CS, respectively. (b), (d), (f) Intensity maps at 1000, 1300, and $1440{\mathrm{cm}}^{-1}$ of hyperspectral images obtained by Fourier transform of generated time-domain interferograms, respectively. (g), (i) CS-reconstructed spectra at the points indicated by arrows in (a) and (c), respectively. (h), (j) Fourier transform of the generated time-domain interferograms at the points indicated by arrows in (b) and (d), respectively. (k) Performance of CS-based time-domain HSI for different measurement durations.

4.

Discussion

It is worthwhile to compare our scheme with other CS-based HSI methods that operate in the frequency domain.^20,26 To date, the frequency-domain approach has used random masks²⁶ or Lissajous scanning²⁰ to randomly sample frequency or spatial components; these reconstructions are only reliable when the spectra are nonsparse because CS works well under the condition that signal sampling occurs in a nonsparse domain.²⁷ On the other hand, in our time-domain approach, spectra composed of sharp peaks are preferred because the energy of such sharp peaks is delocalized in the time domain. Thus FT-CARS spectroscopy is a suitable platform for CS because its spectra have a sparse nature originating from the narrow molecular lines spread in a broad spectral range and the absence of background signals such as fluorescence.

In summary, we have numerically demonstrated CS-based time-domain HSI with a spectral image acquisition rate of $10\times$ higher than that of conventional raster scanning. By combining the Lissajous-scanning–based sampling method and spectral image reconstruction by SSTV, reliable Raman spectral image reconstruction was realized for a simple synthetic bead image and a cell image taken by a spontaneous Raman microscope. We compared different scanning methods, faster sinusoidal scanning, and slower triangular scanning and showed that both schemes are feasible for acquiring Raman hyperspectral images at a video rate (with a measurement duration of 30 ms). The present method requires only slight modifications to the existing FT-CARS imaging setup and is extendable to other time-domain HSI methods, such as FT-IR and FT-two-photon excitation.²⁸