Spectral resolution enhancement of hyperspectral imagery by a multiple-aperture compressive optical imaging system

The Coded Aperture Snapshot Spectral Imaging (CASSI) system captures the three-dimensional (3D) spatio-spectral information of a scene using a set of two-dimensional (2D) random-coded Focal Plane Array (FPA) measurements. A compressive sensing reconstruction algorithm is then used to recover the underlying spatio-spectral 3D data cube. The quality of the reconstructed spectral images depends exclusively on the CASSI sensing matrix, which is determined by the structure of a set of random coded apertures. In this paper, the CASSI system is generalized by developing a multiple-aperture optical imaging system such that spectral resolution enhancement is attainable. In the proposed system, a pair of high-resolution coded apertures is introduced into the CASSI system, allowing it to encode both spatial and spectral characteristics of the hyperspectral image. This approach allows the reconstruction of superresolved hyperspectral data cubes, where the number of spectral bands is significantly increased and the quality in the spatial domain is greatly improved. Extensively simulated experiments show a gain in the peak-signal-to-noise ratio (PSNR), along with a better fit of the reconstructed spectral signatures to the original spectral data.


Introduction 1234
Hyperspectral imaging requires sensing a large amount of spatial information across many wavelengths.Traditional hyperspectral imaging techniques scan adjacent zones of the underlying spectral 1 Hoover Fabián Rueda Chacon.Bachelor of Sciences in Computer Science, Master of Sciences in Computer Science and Informatics, Universidad Industrial de Santander, Colombia.Affiliation: Ph.D. student in Electrical and Computer Engineering at the University of Delaware, USA.E-mail: rueda@udel.edu 2 Alejandro Parada Mayorga.Bachelor of Electronic Engineering, Master of Electronic Engineering, Universidad Industrial de Santander, Colombia.Affiliation: PhD student, University of Delaware, USA.E-mail: alejopm@udel.edu 3Henry Arguello Fuentes.Electrical engineer, Master in Electrical Power, Universidad scene and merge the results to construct a hyperspectral 3-Dimensional (3D) data cube.Push-broom spectral imaging sensors, for instance, capture a spectral data cube by using one FPA measurement per spatial line of the scene (Brady, D. J., 2009).Spectrometers based on optical band-pass filters need to scan the Industrial de Santander, Colombia.PhD in Electrical and Computer Engineering, University of Delaware, USA.Affiliation: Associate professor in full-time dedication of the School of Engineering and Computer Systems of the Universidad Industrial de Santander, Colombia.E-mail: henarfu@uis.edu.coscene by tuning band-pass filters in steps (Eismann, M., 2012).These sensing techniques obey the well-known Nyquist criterion, which imposes a severe limit on the required number of samples.More specifically, these methods require scanning a number of zones linearly in proportion to the desired spatial or spectral resolution.As the desired resolution increases, the required number of samples grows considerably such that the cost of sensing a hyperspectral image is extremely high.Recently, a mathematical technique called Compressive Sensing (CS) has allowed signal sampling at rates below the Nyquist rate (Donoho, D. L., 2006).This new technique involves diverse mathematical areas, such as numerical optimization, signal processing, random matrix analysis, and statistics.The enormous potential of CS has been recently applied in areas such as microscopy, holography, tomography and spectroscopy (Willett, Marcia, and Nichols, 2011;Arguello and Arce, 2013).This paper focuses on the application of CS in spectral imaging; this technique has been termed Compressive Spectral Imaging (CSI).CSI senses 2D coded random projections of the underlying scene such that the number of required projections is far less than those in the linear scanning case.CSI exploits the fact that hyperspectral images can be sparse in some basis representations (Candès and Tao, 2011).Formally, suppose that a hyperspectral signal ∈ Թ ேൈெൈ , or its vector representation ∈ Թ ேெ , is ܵsparse on some basis શ, such that ൌ શી can be approximated by a linear combination of ܵ vectors of શ with ܵ ≪ ‫.ܮܯܰ‬ Here, ܰ ൈ ‫ܯ‬ represents the spatial dimensions, and ‫ܮ‬ is the spectral depth of the image cube.CSI allows to be recovered from ݉ random projections with high probability when ݉ ‫‪ሻ‬ܮܯܰ‪݃ሺ‬݈ܵ‬ ≪ ‫.ܮܯܰ‬ The Coded Aperture Snapshot Spectral Imaging (CASSI) system (Wagadarikar, John, Willett, and Brady, 2008;Arguello and Arce, 2011) is a remarkable imaging architecture that effectively implements CSI.Thus, CASSI senses the 3D spectral information of a scene by using 2D random projections, as depicted in Fig. 1(a).The principal components in CASSI include the coded aperture, the dispersive element and the Focal Plane Array (FPA).The coded aperture patterns are the only varying elements in CASSI, while the other optical elements remain fixed during the operation of the instrument.The input-output relation in CASSI can be expressed as ‫ܡ‬ ൌ ۶, where ‫ܡ‬ represents the random projections, ۶ is the transfer function representing the dispersive element and the coded aperture effects, and is the 3D spectral data cube in vector form (Arguello, Correa and Arce, 2013;Arguello, Rueda and Arce, 2013).Given the compressive measurement ‫,ܡ‬ the objective of CS is to recover an estimate of by using an ℓ ଶ െ ℓ ଵ norm-based optimization algorithm, which exploits the sparsity property of the hyperspectral source.
Despite its potential, CASSI faces a limiting trade-off between spatial and spectral resolution, with the total number of recoverable voxels constrained by the size of the FPA.This constraint limits the utility and cost-effectiveness of compressive hyperspectral imaging for many applications.CSI in infrared (IR) wavelengths is an application where FPAs are particularly critical components, because they become very costly when the resolution increases (Arce, Brady, Carin, Arguello, and Kittle, 2014).As a consequence, spectral super-resolution enhancement is a topic of high interest, because high-resolution reconstructions can be attained from lowresolution/low-cost FPA detectors.This paper presents the spectral resolution enhanced multi-aperture CASSI system (SREM-CASSI), which is a generalization of the CASSI system that includes a new multi-aperture section formed by a dispersive element sandwiched with a pair of high-resolution coded apertures.This configuration leads to multiple-coding flexibility of the spatial and spectral characteristics of the hyperspectral scene, thus permitting the reconstruction of highly resolved scenes from multiple-coded low-resolution FPA 2D projections.In particular, the random projections in SREM-CASSI are given by ‫ܡ‬ ൌ ۲۶, where ۶ is the transfer function accounting for the pair of coded apertures and the dispersive element effects and ۲ is a decimation matrix representing the effect of the low-resolution FPA detector.In the following, we introduce the design of the SREM-CASSI optical architecture, along with its optical and matrix model, as well as simulations to evaluate the attainable improvements.

SREM-CASSI System Model
The proposed SREM-CASSI optical architecture is depicted in Figure 1(a).This is composed by eight optical elements: four lenses, two high-resolution coded apertures, a dispersive element (prism or grating) and a low-resolution detector.The spatio-spectral power source density is denoted as ݂ ሺ‫,ݔ‬ ‫,ݕ‬ ߣሻ, where ‫ݔ‬ and ‫ݕ‬ index the spatial domain and ߣ indexes the wavelengths.The source density is first spatially modulated by the coded aperture ܶ ଵ ሺ‫,ݔ‬ ‫ݕ‬ሻ, resulting in a coded field represented as ݂ ଵ ሺ‫,ݔ‬ ‫,ݕ‬ ߣሻ ൌ ܶ ଵ ሺ‫,ݔ‬ ‫ݕ‬ሻ݂ ሺ‫,ݔ‬ ‫,ݕ‬ ߣሻ.Subsequently, the coded field is sheared by the dispersive element, whose output can be expressed as where ݄ሺ‫ݔ‬ െ ‫ݔ‬ ᇱ െ ܵሺߣሻ, ‫ݕ‬ െ ‫ݕ‬ ᇱ ሻ is the optical impulse response of the system, and ܵሺߣሻ represents the dispersion, which occurs only in the horizontal direction.After dispersion, the source density is then modulated by a second coded aperture ܶ ଶ ሺ‫,ݔ‬ ‫ݕ‬ሻ, resulting in the field ݂ ଷ ሺ‫,ݔ‬ ‫,ݕ‬ ߣሻ ൌ ܶ ଶ ሺ‫,ݔ‬ ‫ݕ‬ሻ݂ ଶ ሺ‫,ݔ‬ ‫,ݕ‬ ߣሻ.Finally, the compressive measurements are realized by the integration of the doubly encoded and dispersed data over the detector's spectral range sensitivity.The spectral density just in front of the detector can be expressed as ݃ሺ‫,ݔ‬ ‫ݕ‬ሻ ൌ ‫‬ ݂ ଷ ሺ‫,ݔ‬ ‫,ݕ‬ ߣሻ ݀ߣ.More specifically, ݃ሺ‫,ݔ‬ ‫ݕ‬ሻ can be written as (2) If the optical impulse response of the system is assumed to be linear, Eq. ( 2) can be succinctly expressed as ݃ሺ‫,ݔ‬ ‫ݕ‬ሻ ൌ න ܶ ଶ ሺ‫,ݔ‬ ‫ݕ‬ሻܶ ଵ ሺ‫ݔ‬ ܵሺߣሻ, ‫ݕ‬ሻ݂ ሺ‫ݔ‬ ܵሺߣሻ, ‫,ݕ‬ ߣሻ ݀ߣ. (3) The coded aperture pixel sizes of ܶ ଵ and ܶ ଶ are denoted as ∆ܿ ଵ and ∆ܿ ଶ , respectively.The transmittance functions of both coded apertures are then given by where ‫ݐ‬ ᇲ ᇲ ሺଵሻ and ‫ݐ‬ ᇲ ᇲ ሺଶሻ are binary values accounting for a translucent (1) or blocking (0) element.The term rectሺሻ represents the rectangular step function.In practice, the coded apertures are implemented through the use of digital micro-mirror devices (DMD) or photomasks.
To choose which coded apertures to use, it is important to take care of the throughput of the system.In SREM-CASSI, the transmittances of the coded apertures define the throughput of the system; therefore, both coded apertures are related.More clearly, the transmittance of the new system is the product of the transmittance of the two coded apertures.Although the distribution of the coded aperture entries can be optimized to achieve better reconstruction results, they can be generated completely at random to show the improvement of the SREM system over the traditional CASSI.Furthermore, the use of random distributions entails high incoherence with the signal representation basis, which assures the correct reconstruction of the signal.Figure 2 shows an example of two typical coded aperture realizations with different transmittance levels, where the white pixels represent translucent elements that allow the light to pass through and the black pixels represent blocking elements.Furthermore, assuming the pixel size of the detector is ∆ ௗ , the integration of the continuous field ݃ሺ‫,ݔ‬ ‫ݕ‬ሻ in a single detector pixel can be expressed as Using Eqs.(3-5) in ( 6), the energy captured in the ሺ݊, ݉ሻ ௧ pixel is expressed as where ߱ , represents the noise of the system.Representing the source density ݂ ሺ‫,ݔ‬ ‫,ݕ‬ ߣሻ in discrete form as ݂ ,, , Eq. ( 7) can be succinctly expressed as for is the ratio between the size of the detector and the coded aperture pixels, and ‫ܮ‬ is the number of spectral bands of the data cube.In this paper, it is assumed that Δ ∈ Ժ ା , that is, the detector and coded aperture pixel sizes satisfy ∆ ௗ ൌ ݇ ଵ ∆ܿ ଵ ൌ ݇ ଶ ∆ܿ ଶ , where ݇ ଵ , ݇ ଶ 1 are integers.
Notice that ܰ ᇱ ൌ ே ௱ and ‫ܯ‬ ᇱ ൌ ቒ ெାିଵ ௱ ቓ, where ܰ ൈ ‫ܯ‬ corresponds to the number of pixels in the first coded aperture and ܰ ൈ ሺ‫ܯ‬ ‫ܮ‬ െ 1ሻ in the second coded aperture.A critical requirement to achieve spectral super-resolution is that the pixel sizes of both coded apertures must be smaller than that of the detector, i.e., ∆ܿ ଵ ൏ ∆ ௗ and ∆ܿ ଶ ൏ ∆ ௗ .

SREM-CASSI Matrix Forward Model
The SREM-CASSI FPA measurements given in Eq. ( 8) can be succinctly expressed in matrix notation as where ‫ܭ‬ is the number of captured snapshots, the matrix ۲ represents the decimation originated by the low resolution detector, and are vector representations of ݃ , and ݂ in Eq. ( 8), respectively, ۶ is the projection matrix accounting for the dispersive element and the ݅ ௧ coded apertures, and the vector ࣓ represents the noise of the system.Notice that the coded apertures ܶ ଵ ሺ‫,ݔ‬ ‫ݕ‬ሻ and ܶ ଶ ሺ‫,ݔ‬ ‫ݕ‬ሻ change for every snapshot.Notice also that, in Eq. ( 9), represents the high-resolution spectral source data cube, whereas the vectors correspond to low-resolution measurements.Figure 1(b) shows a sketch of the sensing process to obtain the low-resolution measurements from the high-resolution spectral scene.The snapshots are taken sequentially, and it is assumed that the underlying spectral scene remains static during the integration time of the K snapshots.The optical transmission function of the system is represented by where ‫۾‬ is a ܰሺ‫ܯ‬ ‫ܮ‬ െ 1ሻ ൈ ‫ܮܯܰ‬ matrix representing the dispersive element operation and ‫܂‬ ଵ and ‫܂‬ ଶ are the matrix representations of the coded apertures used in the ݅ ௧ snapshot.Specifically, ‫܂‬ ଵ is a ‫ܮܯܰ‬ ൈ ‫ܮܯܰ‬ block-diagonal matrix of the form where diagሺ‫ܜ‬ ଵ ሻ represents an ‫ܯܰ‬ ൈ ‫ܯܰ‬ matrix with the elements of ‫ܜ‬ ଵ in the diagonal and ൈ is an NM ൈ NM zerovalued matrix.Notice that the function "diag(x)" is defined as a function that places the elements of the vector parameter x in the diagonal of a matrix.
The second coded aperture ‫܂‬ ଶ operation is modeled in the system as an ܰሺ‫ܯ‬ ‫ܮ‬ െ 1ሻ ൈ ܰሺ‫ܯ‬ ‫ܮ‬ െ 1ሻ matrix, with the values of the second coded aperture in its diagonal.Alternately, the dispersive element operation is represented by the matrix ‫,۾‬ which can be written as where ேெ represents an ܰ‫-ܯ‬long one-valued vector.Finally, ‫܌‬ ൌ ൣሾ ∆ ேି∆ ሿ ⊗ ∆ ൧, where ⊗ is the Kronecker matrix product operation and where Notice that the matrix operation ‫܌‬൫ી ୖ ൯ in Eq. ( 13) shifts the columns of ‫܌‬ by ݇ positions circularly to the right.Consequently, the decimation operation resulting from the low-resolution detector can be modeled as For a multiple-snapshot approach, the general model for SREM-CASSI can be written as (16) Furthermore, Eq. ( 16) can be succinctly expressed as where In particular for reconstruction, the hyperspectral signal ۴ ∈ Թ ேൈெൈ , or its vector representation ∈ Թ ே.ெ., is assumed to be ܵ-sparse on some basis શ, such that ൌ શી.Here, ી are the coefficients of the sparse representation.Hence, can be approximated by a linear combination of ܵ vectors from શ with ܵ ≪ ܰ. ‫.ܯ‬ ‫.ܮ‬ Specifically, an estimation መ of the highresolution data cube from the low-resolution measurements can be achieved by solving the optimization problem where ߬ > 0 is a regularization parameter that balances the conflicting tasks of minimizing the least squares residuals and, at the same time, searches for a sparse solution.

Analysis of the Forward Operators
The singular value spectrums for the sensing matrices based on the random selection of the coded apertures for the SREM-CASSI system and the traditional CASSI system are presented in Fig. 3.
is indicated as a measure of ill-posedness, where ߣ ଵ represents the most significant eigenvalue and ߣ the less significant.As ߢ is smaller, the forward operator ۶ is better posed.It can be observed that, although the spread of the singular values behaves in a similar fashion for both architectures regardless of the transmittance level of the coded apertures, the SREM-CASSI condition number is significantly smaller than that of the traditional CASSI.In consequence, the SREM-CASSI optical design leads to better well-conditioned sensing matrices.

Simulations and Results
A high-resolution spectral data cube exhibiting ‫ܮ‬ ൌ 24 spectral bands and ܰ ൌ ‫ܯ‬ ൌ 256 spatial pixels was experimentally obtained using a wide-band Xenon lamp as light source and a visible monochromator that spans between 451 nm and 642 nm (RGB representation in Fig. 4(a)).The image intensity was captured using a CCD camera with a 656x492 pixel resolution and a pixel size of 9.9 µm.A low-resolution spectral data cube was obtained by clustering the 24 bands into 6 bands.The spectral range is the same for both the high-and low-resolution data cubes.The bandwidth of each spectral band in the high-resolution data cube is 8 nm, whereas the low-resolution data cube exhibits 32 nm per band.
The goal of these experiments is to recover the datacube exhibiting 24 bands from the 6-band data cube.To accomplish this, two high-resolution coded apertures with 256x256 and 256x279 pixel resolutions are employed.The entries of these coded apertures are random realizations of Bernoulli random variables, with different levels of transmittance.To obtain an estimation of the highresolution spectral data cube, the optimization problem in Eq. ( 18) is solved by using the Gradient Projection for Sparse Reconstruction algorithm (GPSR) as it exhibits faster computational speed (Figueiredo, Nowak, and Wright, 2007).In addition, the representation basis શ was set to be the Kronecker product of three bases, શ ଵ ⊗ શ ଶ ⊗ શ ଷ , where the combination શ ଵ ⊗ શ ଶ was the 2D-Wavelet Symlet 8 basis and શ ଷ was the Discrete Cosine basis.Due to the random nature of the coded aperture entries, ten trials were performed for each experiment, and the results were averaged.Three different coded aperture/detector pixel ratios Δ (2, 4, 8) were evaluated, along with six different transmittance levels (10%, 20%, 30%, 50%, 80%, and 100%) of the coded apertures.Figure 4 shows the results for different transmittance levels and the corresponding average PSNR of the reconstruction that was achieved.Note that better results are obtained when the transmittance is lower than 50%, with 10%-20% being the best average transmittance ratio interval.It can be noticed that the results worsen when we approach the CASSI architecture (transmittance = 100%).Using the best transmittance level for each experiment, Fig. 5 shows the reconstruction PSNR vs. the number of captured snapshots for different values of Δ, using the SREM-CASSI and the traditional CASSI architectures.There, it is evident that, as the decimation ratio increases, the reconstruction quality decreases.However, capturing more snapshots can alleviate the loss in quality.Thus, a reconstruction PSNR of 28 dB is achieved by using either Δ ൌ 2 and 4 snapshots, Δ ൌ 4 and 8 snapshots, or Δ ൌ 8 and 64 snapshots.Then, if an eight-times smaller resolution detector is available, roughly eight times more snapshots are required to achieve similar reconstruction results.In Fig. 5, it can also be seen how the results from the SREM-CASSI architecture surpass those achieved by the traditional CASSI.In contrast, Fig. 6 shows the reconstruction results of the righthand side object obtained with the traditional CASSI and the proposed architecture when Δ ൌ 4 is used and 128 shots are captured.It can be easily noticed that the SREM reconstruction quality improves on that obtained with the traditional CASSI.
Finally, Fig. 7 shows the comparison between the reconstructed spectrums of three selected points from the data cube for different number of snapshots and Δ ൌ 2. As the number of captured snapshots increases, the spectral signatures approach the original signature.

Conclusions
A spectral resolution enhancement methodology for coded aperture-based multiple-snapshot spectral imaging systems has been developed.The proposed optical architecture exploits the subpixel information from the original hyperspectral signal by means of two high-resolution coded apertures, attaining richer spectral scenes by using a low-resolution detector but at the cost of capturing multiple FPA measurements.The reconstructions attained up to 32.5 dB of PSNR with half the size of a full-resolved FPA (2 dB decay), 31 dB with a detector four times smaller (3.5 dB decay) and 28.5 dB with an eight-times smaller detector (6 dB decay).
Figure 1.The Spectral Resolution Enhanced Multi-aperture CASSI system (SREM-CASSI).(a) SREM-CASSI optical architecture.(b) Representation of the sensing process.First, the ܰ ൈ ‫ܯ‬ ൈ ‫ܮ‬ data cube is spatially modulated by the first coded aperture ሺܶ ଵ ሻ.Subsequently, the coded data cube is spectrally dispersed by a prism in the horizontal direction and modulated again by the second coded aperture ሺܶ ଶ ሻ.Finally, the spatio-spectral coded information is decimated and integrated by the low-resolution detector.The coded aperture pixel sizes ∆ܿ ଵ and ∆ܿ ଶ are smaller than the resolvable pixel at the detector Δ ௗ .

Figure 2 .
Figure 2. Two coded apertures with different transmittance levels.White pixels allow the light to pass through, while black pixels block the light