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

Multispectral imaging via nanostructured random broadband filtering

Open Access Open Access

Abstract

It is a challenge to acquire a snapshot image of very high resolutions in both spectral and spatial domain via a single short exposure. In this setting one cannot trade time for spectral resolution, such as via spectral bands scanning. Cameras of color filter arrays (CFA) (e.g., the Bayer mosaic) cannot obtain high spectral resolution. To overcome these difficulties, we propose a new multispectral imaging system that makes random linear broadband measurements of the spectrum via a nanostructured multispectral filter array (MSFA). These MSFA random measurements can be used by sparsity-based recovery algorithms to achieve much higher spectral resolution than conventional CFA cameras, without sacrificing spatial resolution. The key innovation is to jointly exploit both spatial and spectral sparsity properties that are inherent to spectral irradiance of natural objects. Experimental results establish the superior performance of the proposed multispectral imaging system over existing ones.

© 2020 Optical Society of America under the terms of the OSA Open Access Publishing Agreement

1. Introduction

Although digital color photography of very high spatial resolution is taken for granted nowadays, it still remains a challenge to acquire an image of simultaneously high spectral and spatial resolutions in a single, short exposure. The task is difficult because of the fundamental imcompatibility between the 2D organization of pixel sensors in modern CMOS and CCD imaging devices and the 3D nature of the multispectral images. A multispectral image of $K$ bands can be viewed as a three-dimensional $N_x \times N_y \times K$ data cube, with a two-dimensional image plane (sampled on an $N_x \times N_y$ lattice) plus the third dimension of spectrum. As this 3D volume of pixels cannot be captured directly by the 2D pixel sensor array, some indirect apparatus becomes necessary to circumvent the incompatibility.

A straightforward solution is to use an optical beam splitter, such as prism, to separte light into $K$ spectral bands, each of which is captured by a 2D pixel sensor array [1]. But this camera architecture of $K$ sensor arrays is clumsy and expensive; a more serious handicap is that splitting incoming lights into $K$ parts rapidly reduces light influx for each spectral band down to an unacceptable level.

A common approach of acquiring multispectral images of $K$ bands with a 2D sensor array is temporal multiplexing: taking $K$ exposures of the imaged scene one at a time, one exposure per spectral band. Each of the $K$ spectral bands is acquired through switchable filtering elements installed in the light path (tunable filters [2,3], tunable illumination [4,5], etc). But the temporal multiplexing strategy must employ complex, delicate moving parts of optics, and hence it has high cost and poor portability. Another serious drawback of temporal multiplexing is its susceptibility to motion artifacts [6].

A more practical approach of multispectral image acquisition is spatial multiplexing via the use of an $N_x \times N_y$ multispectral filter array (MSFA) [613]. $K$ spectral bands are subsampled in space and the resulting $K$ types of spectral samples are interleaved in a mosaic pattern in the 2D sensor array (illustrated in Fig. 1(a)). From these MSFA samples the underlying multispectral image needs to be reconstructed by solving an ill-posed inverse problem. The multispectral imaging system of spatial multiplexing is schematically described in Fig. 2. The most common spatial multiplexing scheme for multispectral photography (in visible spectrum) is that of Bayer color filter array (CFA) [14], which can be found in all types of consumer electronics and industrial applications. The spatial multiplexing-based multispectral cameras, compared with their temporal multiplexing counterparts, enjoy the advantages of simpler hardware architecture, lower cost and higher portability [6,15,16]. But they suffer from reduced spatial resolution of each spectral band due to the strategy of trading off between the spatial and spectral resolutions. This problem becomes more acute if the target number of spectral bands $K$ increases, as the sampling rate for the underlying $N_x \times N_y \times K$ 3D data cube of the multispectral image decreases at rate $1/K$. The main thrust of this research is to improve the reconstruction of multispectral images despite the low sampling rate.

 figure: Fig. 1.

Fig. 1. A multispectral filter array and the constituent narrowband filters with the peak wavelength evenly shifted.

Download Full Size | PDF

 figure: Fig. 2.

Fig. 2. A typical multispectral imaging system of spatial multiplexing type.

Download Full Size | PDF

In existing CFA designs, the $K$ types of filter elements have a Gaussian-shaped, narrowband spectral response with the peak wavelength evenly shifted. This is essentially a scheme of uniformly sampling the spectrum. Uniform sampling is less efficient than random sampling according to the theory of compressive sensing (CS). CS is an extensively-studied and proven technique for reducing the sampling rate necessary to recover a signal [1719], as long as the signal is sparse in some space. Indeed, one of the earliest applications of CS is the coded acquisition of single-band (monochrome) images at sampling rates much lower than conventional. The CS theory requires the samples to be random linear measurements of the image. In the case of CS-based monochrome imaging, the DMD (digital micromirror device) technology is used to make random linear measurements [20]. However, in the main stream designs of MSFA, the spectrum domain is sampled approximately uniformly; in other words, $K$ narrowband spectral filters of evenly spaced peak wavelength are used, as illustrated in Fig. 1(b). For this matter, existing temporal multiplexing schemes for multispectral imaging also perform approximately uniform sampling in the spectral domain. The CS theory presents a strong critique against the existing practice of approximately uniformly sampling the spectrum, suggesting that significant room still exists to improve the performance of existing multispectral cameras by random sampling.

A pioneer work on CS-like multispectral imaging is the work published by Wagadarikar and Gehm et al. [21,22], called coded aperture snapshot spectral imaging (CASSI). The CASSI system utilizes a coded aperture and one or more dispersive elements to modulate a spectral scene. A sensor array captures a multiplexed 2D projection of the 3D data cube representing the scene. The specific sampling kernel of the modulation depends on the relative positioning of the coded aperture and the dispersive element(s) in the system. One of problems with CASSI is that the proposed optical mechanism of dispersion cannot produce strictly linear measurements of the spectral signal; hence, the CASSI system is ill suited for a CS-based recovery of multispectral images.

In this paper, we propose a novel multispectral imaging system of spatial multiplexing. The system uses a new type of $N_x \times N_y$ MSFA to take $N_x N_y$ random linear measurements in the spectrum, one per pixel location. This special MSFA can be built by the technology of nanostructured random broadband filtering based on surface plasmonic polariton (SPP). From these random linear measurements of MSFA, the multispectral image is recovered by an optimal reconstruction algorithm that jointly exploits both spectral and spatial sparsity via the minimization of a so-called $\ell _{1,2}$-norm of the Laplacian matrix, which is assembled by the spatial second derivatives of all pixels. The proposed method combines random measurement of the spectrum and joint spatial-spectral sparse approximation together, and can recover highly local, distinctive features of the signal in the spectral domain, such as spikes, convexity, inflection on the spectral response curve. These distinctive local features in spectrum are vital for many multispectral imaging applications.

The rest of the paper is organized as follows. In Section 2, we design the random broadband multispectral filters based on SPP and model the random sampling process of the spectrum. An multispectral image formation model is reviewed in Section 3. Based on this model, we propose a new recovery algorithm based on joint spatial-spectral sparsities to push for higher restoration performance. Experiment results are reported in Section 4; they corroborate our analysis and demonstrate the superiority and effectiveness of the proposed method. Finally, Section 5 concludes the paper.

2. Random sampling via broadband MSFA

2.1 Multispectral sampling process

To generate spectral images, a MSFA camera is equipped with sensors of different types to record the spectral irradiance of a real world scene. In existing MSFA cameras of $K$ spectral bands, MSFA is of a pre-determined pattern of $K$ different narrowband filters that are spatially arranged in a periodic manner. Up to now the research on MSFA design focuses either on improving spatial arrangement of the filters [7,8,11], or on broadening the spectral range (LWIR [12,23], NIR [13] and UV [24], etc.). To the best of our knowledge, no work has been reported on designing the type of MSFA that can generate random linear measurements of the spectrum with broadband filters. This void is quite a surprise given the common knowledge that random sampling is more efficient for signal reconstruction than regular, uniform sampling. The existing method to increase spectral resolution is expedient: partitioning the target spectrum range into narrower bands and covering them with corresponding narrowband filters. Such a design is to merely trade spatial resolution for spectral resolution. Alternatively, if different broadband filters of mutually overlapping, multimodal spectral responses can be realized in MSFA, which will be demonstrated shortly in the next subsection, then the pixel values become random linear measurements of the spectrum, carrying far more information on the spectral signatures of the imaged scene. With these random measurements, one can take the recovery approach of compressive sensing to achieve higher spectral resolution.

In the following, we first model the process of random sampling via a broadband MSFA. Let $x_{i}(\lambda )$ ($i=1, 2, \ldots , N$) be the spectral reflectance of the surface point pertaining to pixel $i$, with $\lambda$ being the wavelength and $N$ being the number of pixels. Pixel $i$ in the MSFA camera has a broadband filter with a spectral response function $\gamma _i(\lambda )$. Under the illumination of a light source with strength spectrum $l(\lambda )$, the spectral distribution of the reflected light, recorded by pixel $i$, is given by $l(\lambda ) x_{i}(\lambda )$. Therefore, the sample value at pixel $i$ is given by

$$y_i = \int_\Lambda \gamma_i(\lambda) l(\lambda) x_{i}(\lambda) d\lambda$$
where $\Lambda =[\lambda _{min},\lambda _{max}]$ is the spectral range of the MSFA camera.

Assuming that the underlying multispectral image is to be reconstructed to a resolution of $K$ bands, we partition the imaged spectral range $\Lambda$ equally into $K$ subranges, i.e., $\Lambda =\bigcup _{k=1}^{K}\Lambda _k$; consequently, Eq. (1) can be rewritten as

$$\begin{aligned} y_i & =\sum_{k=1}^{K} \int_{\Lambda_k} \gamma_i(\lambda) l(\lambda) x_{i}(\lambda) d\lambda\\ & =\sum_{k=1}^{K} x_{k,i}\int_{\Lambda_k} \gamma_i(\lambda) l(\lambda) d\lambda,\\ \end{aligned}$$
where $x_{k,i}$ is the $\gamma _i(\lambda ) l(\lambda )$-weighted average of $x_{i}(\lambda )$ in subrange $\Lambda _k$. Letting
$$a_{k,i} = \int_{\Lambda_k} \gamma_i(\lambda) l(\lambda) d\lambda ,$$
Eq. (2) becomes
$$y_i= \sum_{k=1}^{K} a_{k,i} x_{k,i}.$$
As the weight $a_{k,i}$ is decided by the spectral response $\gamma _i(\lambda )$ of the filter at pixel $i$, one can make $y_i$ a random measurement of the spectrum of pixel $i$ by randomizing $\gamma _i(\lambda )$. In the next subsection, we introduce a new nanostructured MSFA design based on surface plasmonic polariton that can realize different random spectral responses $\gamma _i(\lambda )$ at different pixels.

2.2 Random broadband filter based on SPP

Recently, plasmonic optical filters have emerged as an attractive alternative to main stream color filters of organic dyes, which can be coupled with image sensors for multispectral imaging [25,26]. The above design goal of embedding different broadband linear filters into MSFA can be achieved by plasmonic filters because they offer high tunability in transmission spectra. The type of filters is based on surface plasmonic polariton. SPP can be excited by an incident electromagnetic wave if their wavelength vectors match. This is usually achieved by nanopatterning a metal film. The resonant frequency of SPP is determined by the metal materials, dielectric materials, profiles and dimensions of the patterns, etc. As a result, the tunability of SPP enables its application as a spectral filter [25]. By changing the material and the structure dimensions of the nanopatterns, both the resonance wavelength and the band width can be tuned. In principle, an unlimited number of filters with different spectral responses can be designed.

For a proof of the idea, we design a MSFA with nine broadband filters in a $3\times 3$ pattern as shown in Fig. 3. The filters consist of square holes in a square lattice array in an aluminium film with a thickness of 100nm, where the aluminium film is sandwiched between semi-unlimited SiO$_2$ lower cladding and 100nm Si upper cladding. The square holes are filled with Si. When light illuminates this structure, extraordinary transmission (EOT) [27] can be observed at different wavelengths due to surface plasmon resonances at Al-SiO$_2$ and Al-Si interface. The resonance wavelength can be estimated as [28]

$$\begin{aligned} &\lambda_{peak}= \\ &\frac{2\pi \sqrt{\frac{\epsilon_m\epsilon_d}{\epsilon_m+\epsilon_d}}}{\sqrt{(k_{||}\cos\phi+i\frac{2\pi}{a_x}+j\frac{2\pi}{a_x})^{2}+(k_{||}\sin\phi-i\frac{2\pi}{a_y}+j\frac{2\pi}{a_y})^{2}}}, \end{aligned}$$
where $k_{||}$ is the in-plane wave vector magnitude and $\phi$ is the azimuthal angle of incident light, $a_x$ and $a_y$ are the lattice dimensions, $i$ and $j$ are the scattering orders of the array. $\epsilon _m$ and $\epsilon _d$ are the permittivities for the metal and dielectric medium respectively. Using finite-difference time-domain (FDTD) simulation, transmission spectra are calculated numerically and the results are shown in Fig. 4. By tuning the grating width $a$ in a range of 200$\sim$400nm and the period $p$ in a range of 450$\sim$650nm (see Fig. 3), the spectral response is randomly distributed in the visible band, which is ideal for the spectral imaging. These filters have been experimentally demonstrated by either electron-beam lithography [29] and nanoimprint [30] techniques. Furthermore, the plasmonic filters have been successfully integrated in CMOS image sensors [31,32], which encourages the practical application of the proposed spectral imaging technique.

 figure: Fig. 3.

Fig. 3. Concept schematic of nanostructured MSFA with SPP-based broadband filters.

Download Full Size | PDF

 figure: Fig. 4.

Fig. 4. Nine SPP-based broadband filters.

Download Full Size | PDF

3. Recovery based on joint spectral-spatial sparsities

In this section, we investigate how to recover multispectral images from random spectral measurements of the SPP-based broadband MSFA. In designing our recovery algorithm a key innovation is to explore and exploit joint spectral-spatial sparsities that are intrinsic of imaging natural scenes.

Using the bidirectional reflectance distribution function (BRDF) that models how light is reflected at an opaque surface, the spectral reflectance $x_{k,i}$ in Eq. (2) is determined by two factors: the material and the geometry of the object surface, and it can be modeled as

$$x_{k,i}= \tau_{k,i} \cos \alpha_i + \kappa_{k,i} \eta \cos^{s} \beta_i,$$
where the term $\cos \alpha _i$ accounts for Lambertian diffuse reflection, with $\alpha _i$ being the angle between the incident light and the surface normal at the surface point corresponding to pixel $i$; the term $\eta \cos ^{s} \beta _i$ approximates specular reflection, where $s>1$ is the degree of surface shininess, and $\beta _i$ is the angle of the reflection ray with the mirror direction of the incident angle. $\tau _{k,i}$ and $\kappa _{k,i}$ are the diffuse reflectance coefficient and the specular reflectance coefficient of pixel $i$ in band $k$, which are both determined by the spectral signatures of the surface material but independent of the geometric attributes.

It should be mentioned that deriving joint spectral-spatial sparsities from the BRDF formulation, as will be conducted below, is an indirect analysis because the BRDF model is an approximation. A direct and potentially more accurate method is to collect and analyze the statistics of spectral reflectance from a large natural spectral dataset, which requires a great deal of efforts.

Now consider a surface area $\Omega$ consisting of a uniform material in which all pixels $i \in \Omega$ have the same spectral reflectance coefficients $\tau _k$ and $\kappa _k$ in band $k$. In surface area $\Omega$, band $k$ of the multispectral spectral image ${{\mathbf {x}}}$ can be written as

$$\begin{aligned} {{\mathbf{x}}}_k & = (x_{k,1}, x_{k,2}, \ldots, x_{k,N})^{T}\\ & =\tau_k (g_{1}, g_{2}, \ldots, g_{N})^{T} +\kappa_k (h_{1}, h_{2}, \ldots, h_{N})^{T} \\ & =\tau_k {{\mathbf{g}}} +\kappa_k {{\mathbf{h}}}, \end{aligned}$$
where ${{\mathbf {g}}}=(\cos \alpha _1, \cos \alpha _2, \ldots , \cos \alpha _N)^{T}$ and ${{\mathbf {h}}}=(\eta \cos ^{s} \beta _1, \eta \cos ^{s} \beta _2, \ldots , \eta \cos ^{s} \beta _N)^{T}$ are the geometric attribute vectors of all $N$ pixels of ${{\mathbf {x}}}$. Note that $\tau _k$ and $\kappa _k$ both remain a constant for all pixels $i \in \Omega$, because it is decided only by the spectral signature of $\Omega$, independent of the geometric attributes of $\Omega$ such as surface curvature and roughness. In contrast, $\cos \alpha _i$ and $\eta \cos ^{s} \beta _i$ vary in surface area $\Omega$ as angles $\alpha _i$ and $\beta _i$ change from pixel to pixel. In other words, all variations in band ${{\mathbf {x}}}_k$ in surface area $\Omega$ are attributed to those of ${{\mathbf {g}}}$ and ${{\mathbf {h}}}$. By the analysis in [33,34], ${{\mathbf {g}}}$ and ${{\mathbf {h}}}$ can be modeled by piecewise linear functions on surface $\Omega$ of small curvature, or even by piecewise constant functions on flat surfaces. Therefore, ${{\mathbf {x}}}_k =\tau _k {{\mathbf {g}}}+ \kappa _k {{\mathbf {h}}}$ can also be approximated by piecewise linear functions, and consequently the Laplacian (the second derivative) $\nabla ^{2} {{\mathbf {x}}}_k$ is zero at most points in $\Omega$ and has nonzero values only at positions of surface discontinuities, i.e., $\nabla ^{2} {{\mathbf {x}}}_k$ being a sparse representation of ${{\mathbf {x}}}_k$. Moreover, as the geometry terms ${{\mathbf {g}}}$ and ${{\mathbf {h}}}$ is the same in all spectral bands, the discontinuities in images ${{\mathbf {x}}}_k$, $k=1, 2, \ldots , K$, are spatially aligned with those of ${{\mathbf {g}}}$ and ${{\mathbf {h}}}$. Thus the Laplacians $\nabla ^{2} {{\mathbf {x}}}_k$ are not only each sparse, but also they take on nonzero values simultaneously in $\Omega$.

Another type of discontinuity in ${{\mathbf {x}}}_k$ is directly caused by changes in $\tau _k$ and $\kappa _k$. This corresponds to the boundaries of two objects of different spectral signatures or to the locations where surface materials switch. Clearly, in these cases, the Laplacians $\nabla ^{2} {{\mathbf {x}}}_k$ will also have spatially aligned nonzero values.

To summarize, the Laplacian images $\nabla ^{2} {{\mathbf {x}}}_k$, $k=1, 2, \ldots , K$, tend to agree with one the other in the positions of zero and nonzero values, regardless whether the nonzero values are due to geometry discontinuities or to change of surface materials. Let

$$D({{\mathbf{x}}})= \left[ \nabla^{2} {{\mathbf{x}}}_{1}, \nabla^{2} {{\mathbf{x}}}_{2}, \ldots, \nabla^{2} {{\mathbf{x}}}_{K} \right]$$
be the $N \times K$ matrix whose columns are the $K$ Laplacians of ${{\mathbf {x}}}={({{\mathbf {x}}}_1, {{\mathbf {x}}}_2, \ldots , {{\mathbf {x}}}_K)}$. As aforementioned, most rows of $D({{\mathbf {x}}})$ are zero vector, and only few rows include non-zero elements. This property can be imposed to a multispectral image ${{\mathbf {x}}}$ by minimizing the $\ell _{0,2}$-norm of the matrix $D({{\mathbf {x}}})$,
$$\lVert D({{\mathbf{x}}}) \rVert_{0,2}=\# \{i | \lVert {{\mathbf{r}}}_i \rVert_2 \neq 0, i=1, 2, \ldots, N\},$$
where ${{\mathbf {r}}}_i$ is pixel $i$ or row $i$ of matrix $D({{\mathbf {x}}})$. Since the $\ell _{0,2}$-norm minimization is a non-convex, NP-hard combinatorial optimization problem, it is often relaxed to the convex $\ell _{1,2}$-norm minimization, that is
$$\lVert D({{\mathbf{x}}}) \rVert_{1,2}= \sum_{i=1}^{N} \lVert {{\mathbf{r}}}_i \rVert_{2}.$$
In (10), the so-called $\ell _{1,2}$-norm of a matrix, denoted by $\lVert \cdot \rVert _{1,2}$, is the sum of the $\ell _2$-norm of all rows of the matrix. Minimizing the $\ell _{1,2}$-norm encourages nonzero elements to occur in few rows, which is called simultaneous sparsity [35,36]. This sparsity allows us to pose multispectral reconstruction as the following problem of constrained minimization:
$$\begin{aligned} \min ~~ &\lVert D({{\mathbf{x}}}) \rVert_{1,2} \\ s.t. ~~& y_i= \sum_{k=1}^{K} a_{k,i} x_{k,i},~~~~i=1, 2, \ldots, N. \end{aligned}$$
Next we explore the low rank property of multispectral images. A well-known property of irradiance spectra is that they are essentially lower-dimensional signals. It has been shown that most spectra may be represented by coefficients corresponding to as few as 6 basis vectors [37,38]. Typically there are high correlations among the spectrum domain, because the surface can be decomposed into few number of local patches with distinct spectral reflectance. This simple observation justifies that multispectral images are low rank. Therefore, we add a low rank constraint on ${{\mathbf {x}}}$ into (11) and finally arrive at
$$\begin{aligned} \min ~~ &\lVert D({{\mathbf{x}}}) \rVert_{1,2} \\ s.t. ~~& y_i= \sum_{k=1}^{K} a_{k,i} x_{k,i},~~~~i=1, 2, \ldots, N, \\ &Rank({{\mathbf{x}}})\leq R, \end{aligned}$$
where $R$ is the maximal rank of ${{\mathbf {x}}}$. Usually, the above minimization problem with low rank constraints can be written as a nuclear norm minimization [39,40]
$$\begin{aligned} \min ~~ &\lVert D({{\mathbf{x}}}) \rVert_{1,2}+\omega\lVert {{\mathbf{x}}} \rVert_* \\ s.t. ~~& y_i= \sum_{k=1}^{K} a_{k,i} x_{k,i},~~~~i=1, 2, \ldots, N. \end{aligned}$$
where the nuclear norm $\lVert \cdot \rVert _*$ is defined as the sum of the absolute values of the singular values of the matrix, and $\omega$ is a positive constant.

Finally, it should be mentioned, because the camera recorded surface irradiance is the product of reflectance and illuminance, if illuminance is constant or sufficiently smooth in spectral domain, then all the above reflectance based sparsity properties clearly hold for irradiance as well. As imaging is mostly carried out under sun lighting or sufficiently white indoor lighting, the spectral profile of the illuminance is indeed flat or very smooth.

The minimization problem (13) is convex and can hence be solved by any of many existing optimization techniques. We choose the iterative shrinkage-thresholding algorithm based on the alternating direction method [39,41]. The algorithm is detailed in the following pseudo-code.

4. Experimental results

The proposed multispectral image recovery algorithm is implemented and tested on simulated data of the SPP-based MSFA of broadband random spectral filters. The experimental results are reported in this section and compared with those of the existing demosaicking method using narrow band filters.

In our experiments, we generate our ground truth images of 9 spectral bands by uniformly partitioning the 31-bands Columbia University image set [42,43] . The color versions of the sixteen sample images in this test set are listed in Fig. 5. In the experiments, we use the MSFA with nine random broadband filters (shown in Fig. 4) arranged in a $3 \times 3$ pattern to perform random sampling of the spectrum. These nine random broadband filters are selected out of 100 ones that are randomly produced by the COMSOL Multiphysics software, under the selection criterion of minimal pairwise correlations. For a fair comparison against the existing method of uniform spectral sampling, nine narrowband Gaussian filters, which have evenly spaced peak wavelengths, are also arranged in the same $3\times 3$ pattern to form a conventional MSFA. The two different MSFAs, the one of random broadband filters and the other of uniformly spaced narrowband filter, are used to make two equally-sized sets of spectral samples, which are then processed by Algorithm 1, respectively. In other words, the minimization problem (13) is solved for each of the two MSFA spectral sampling schemes.

 figure: Fig. 5.

Fig. 5. Color image versions of the test set in Table 1.

Download Full Size | PDF

Tables Icon

Table 1. PSNR results of the test set in Fig. 5.

In order to make our simulation more realistic and the results indicative of a real physical device, we account for the effects of two types of noises: 1) the errors in the spectral response functions of random broadband filters; 2) the sensor noises. The former is the inaccuracy of the filter model used by the reconstruction algorithm to approximate the true spectral responses of a broadband filter. Although filter model errors can be reduced or removed by making precise measurements of true spectral response curves of all broadband filters via careful calibration, we still run our simulations with a 5% random Gaussian error embedded into $\{a_{k,i}\}_{k=1, 2, \ldots , K, i=1, 2, \ldots , N}$ in (13). In addition, we add a Gaussian white noise of strength $\sigma$ into sensor readings $y_i$, $i=1, 2, \ldots , N$, in (13); the simulations are run for $\sigma = 0$ (no sensor noises), $\sigma = 5$ and $\sigma = 10$, respectively.

The comparison results of the two MSFA spectral sampling schemes are presented in Table 1, in which the peak signal to noise ratio (PSNR) of the reconstructed multispectral images is included. In the reported experiments without the filter model errors, the proposed method, on average, enjoys PSNR performance gains over the conventional system by 2.61 dB, 4.04 dB and 5.47 dB at $\sigma = 0$ , $\sigma = 5$ and $\sigma = 10$. Furthermore, when there are the 5% filter model errors in the experiments, the performance gains of the proposed method become 2.45 dB, 3.41 dB and 5.67 dB, respectively. These results show that this new approach has the better performance no matter whether there are sensor noises and model errors or not, and is actually practically feasible.

We also evaluate the spectral accuracy of competing methods. The difference between the target and reconstructed spectral signatures can be measured by mutual correlation of the two spectral signals; the higher the mutual correlation the higher reconstruction quality in spectral domain. The comparison results of the two competing MSFA spectral sampling schemes are presented in Table 2, in which the averaged mutual correlation $\rho$ between the spectral images ${{\mathbf {x}}}$ and ${{\mathbf {y}}}$ is calculated by

$$\rho =\frac{1}{N}\sum_{i=1}^{N} \frac{\sum_{k=1}^{K} x_{k,i} y_{k,i}}{\sqrt{\sum_{k=1}^{K} x_{k,i}^{2}}\sqrt{\sum_{k=1}^{K} y_{k,i}^{2}}}.$$
As shown by the table, the proposed random MSFA scheme achieves higher mutual correlation with the ground truth spectrum than the regular MSFA scheme, on every test 9-bands multispectral image. On average, the increase in mutual correlation by the proposed method is 0.038, which is considered very significant in the field of multispectral imaging.

Tables Icon

Table 2. Averaged mutual correlations of the narrow band MSFA method and the proposed random MSFA method (no filter model errors, and $\sigma = 5$) on the test set in Fig. 5.

In addition, we also compute the MSE of each of the 9 spectral bands, and plot in Fig. 6 the per band results for the image of Fig. 5(a). The comparison results for all 16 images in the test set are tabulated in Table 3. As shown by the table, the proposed random MSFA scheme achieves lower MSE than the regular MSFA scheme for each of the 16 test images.

 figure: Fig. 6.

Fig. 6. MSE of each of the 9 spectral bands in the test multispectral image Fig. 5(a).

Download Full Size | PDF

Tables Icon

Table 3. MSE results for all 16 images in Fig. 5.

Next, let us qualitatively inspect the spectra of different methods. In many applications of multispectral imaging, users scrutinize and compare distinctive signatures in the spectra of imaged objects. It is, therefore, important to compare the abilities of different multispectral imaging systems to recover subtleties in the spectral curves, such as the points of reflection, peak and valley, the order of continuity, etc. In Fig. 7 and Fig. 8, we plot the spectra of selected image patches recovered by our method and the conventional method so that the reader can inspect the behaviors of the two methods in relation to the ground truth spectra. The higher spectral accuracy of the proposed method (blue curve) over the conventional method (red curve) is clear through the following observations, which are also marked by the green boxes in these figures.

  • • For instances, in Fig. 7(a), the ground truth spectrum (black curve) from 550 to 700nm is a monotonically increasing cubic curve, but the spectrum recovered by the conventional method (red curve) erroneously becomes a nonmonotonic polyline and creates a false peak point near 660nm. The proposed method (blue curve) almost coincides with the ground truth, preserving all distinct spectral signatures of the original, such as peak/valley locations and convexity.
  • • In Fig. 7(b), the conventional method produces an erroneous valley point near 660nm; it also changes the shape of spectral curve from 400nm to 550nm from quadratic to triangular.
  • • In Fig. 8(a), from 475nm to 550nm the original spectral curve is cubical like but it is flattened by the conventional method to linear; the original valley point near 510nm is removed. This moves the original valley point from 512nm to 475nm. Again, the proposed method has none of these errors.
  • • Similarly, in Fig. 8(b), the original cubic-shaped spectral curve from 590nm to 700nm with two reflection points is simply straightened by the conventional method, but not by the proposed method.

 figure: Fig. 7.

Fig. 7. Recovered spectra of two selected patches in the test multispectral image Fig. 5(d).

Download Full Size | PDF

 figure: Fig. 8.

Fig. 8. Recovered spectra of two selected patches in the test multispectral image Fig. 5(h).

Download Full Size | PDF

In all the above cases, the proposed method faithfully reproduces the original spectra. The difference between the two methods in reconstruction quality should be expected, as random broadband spectral sampling carries more information on the original spectral signal than narrowband uniform spectral sampling.

 figure: Fig. 9.

Fig. 9. Recovered spectra of two selected patches in the test multispectral image Fig. 5(i).

Download Full Size | PDF

In order to demonstrate the potential of the proposed multispectral imaging system in spectral analysis and recognition tasks, we present two examples of distinguishing between fake and real products (see Fig. 9 and Fig. 10). In each photograph, there are a real produce and a plastic model that are difficult to distinguish by naked eyes. One has to compare the spectral signatures of the real and fake objects to tell them apart. As shown in Fig. 9 and Fig. 10, the proposed method reproduces the spectral curves of both real and fake objects with so high a precision that it can make reliable identification of the imaged object or substance.

 figure: Fig. 10.

Fig. 10. Recovered spectra of two selected patches in the test multispectral image Fig. 5(k).

Download Full Size | PDF

5. Conclusion

We proposed a novel snapshot multispectral imaging system that consists of a nanostructured random broadband MSFA and an accompanying recovery algorithm that exploits joint spectral-spatial sparsities. The proposed system can significantly increase the precision of other multispectral imaging systems using narrowband MSFA.

Funding

National Key Research and Development Program of China (2019YFA0706604, 2019YFB2203402); National Natural Science Foundation of China (61976169, 61871304, 61836008, 61875157, 61632019); Natural Sciences and Engineering Research Council of Canada.

Disclosures

The authors declare no conflicts of interest.

References

1. B. A. Spiering, “Multispectral imaging system,” U.S. Patent, No. 5,900,942 (1999).

2. F. König and W. Praefcke, “Practice of multispectral image acquisition,” Proc. SPIE 3409, 34–41 (1998). [CrossRef]  

3. N. Gat, “Imaging spectroscopy using tunable filters: a review,” Proc. SPIE 4056, 50–64 (2000). [CrossRef]  

4. M. Parmar, S. Lansel, and J. Farrell, “An led-based lighting system for acquiring multispectral scenes,” Proc. SPIE 8299, 82990P (2012). [CrossRef]  

5. M. C. Lin, C. W. Tsai, and C. H. Tien, “Spectral image reconstruction by a tunable led illumination,” Proc. SPIE 8870, 88700C (2013). [CrossRef]  

6. P. J. Lapray, X. Wang, J. B. Thomas, and P. Gouton, “Multispectral filter arrays: Recent advances and practical implementation,” Sensors 14(11), 21626–21659 (2014). [CrossRef]  

7. J. Brauers and T. Aach, “A color filter array based multispectral camera,” in Proceedings of Workshop Farbbildverarbeitung, (Ilmenau, 2006), (2006), pp. 5–6.

8. R. Ramanath, W. E. Snyder, and H. Qi, “Mosaic multispectral focal plane array cameras,” Proc. SPIE 5406, 701–712 (2004). [CrossRef]  

9. L. Miao, H. Qi, and W. E. Snyder, “A generic method for generating multispectral filter arrays,” in Proceedings of IEEE International Conference on Image Processing (IEEE, 2004), pp. 3343–3346.

10. L. Miao and H. Qi, “The design and evaluation of a generic method for generating mosaicked multispectral filter arrays,” IEEE Trans. on Image Process. 15(9), 2780–2791 (2006). [CrossRef]  

11. S. W. Wang, C. Xia, X. Chen, W. Lu, M. Li, H. Wang, W. Zheng, and T. Zhang, “Concept of a high-resolution miniature spectrometer using an integrated filter array,” Opt. Lett. 32(6), 632–634 (2007). [CrossRef]  

12. J. Mercier, T. Townsend, and R. Sundberg, “Utility assessment of a multispectral snapshot lwir imager,” in Proceedings of Workshop on Hyperspectral Image and Signal Processing: Evolution in Remote Sensing (IEEE, 2010), pp. 1–5.

13. Y. M. Lu, C. Fredembach, M. Vetterli, and S. Süsstrunk, “Designing color filter arrays for the joint capture of visible and near-infrared images,” in Proceedings of IEEE International Conference on Image Processing (IEEE, 2009), pp. 3797–3800.

14. B. E. Bayer, “Color imaging array,” U.S. Patent, No. 3,971,065 (1976).

15. S. H. Baek, I. Kim, D. Gutierrez, and M. H. Kim, “Compact single-shot hyperspectral imaging using a prism,” ACM Trans. Graph. 36(6), 1–12 (2017). [CrossRef]  

16. X. Lin, Y. Liu, J. Wu, and Q. Dai, “Spatial-spectral encoded compressive hyperspectral imaging,” ACM Trans. Graph. 33(6), 1–11 (2014). [CrossRef]  

17. E. J. Candès, J. Romberg, and T. Tao, “Robust uncertainty principles: Exact signal reconstruction from highly incomplete frequency information,” IEEE Trans. Inf. Theory 52(2), 489–509 (2006). [CrossRef]  

18. E. J. Candès and J. Romberg, “Quantitative robust uncertainty principles and optimally sparse decompositions,” Found Comput. Math. 6(2), 227–254 (2006). [CrossRef]  

19. D. L. Donoho, “Compressed sensing,” IEEE Trans. Inf. Theory 52(4), 1289–1306 (2006). [CrossRef]  

20. M. F. Duarte, M. A. Davenport, D. Takhar, J. N. Laska, T. Sun, K. E. Kelly, and R. G. Baraniuk, “Single-pixel imaging via compressive sampling,” IEEE Signal Process. Mag. 25(2), 83–91 (2008). [CrossRef]  

21. M. Gehm, R. John, D. Brady, R. Willett, and T. Schulz, “Single-shot compressive spectral imaging with a dual-disperser architecture,” Opt. Express 15(21), 14013–14027 (2007). [CrossRef]  

22. A. Wagadarikar, R. John, R. Willett, and D. Brady, “Single disperser design for coded aperture snapshot spectral imaging,” Appl. Opt. 47(10), B44–B51 (2008). [CrossRef]  

23. J. Hershey and Z. Zhang, “Multispectral digital camera employing both visible light and non-visible light sensing on a single image sensor,” U.S. Patent, No. 7,460,160 (2008).

24. H. K. Aggarwal and A. Majumdar, “Multi-spectral demosaicing technique for single-sensor imaging,” in Proceedings of National Conference on Computer Vision, Pattern Recognition, Image Processing and Graphics (IEEE, 2013), pp. 1–4.

25. Q. Chen, J. He, X. Shi, and Y. Ma, Application of Surface Plasmon Polaritons in CMOS Digital Imaging (INTECH Open Access Publisher, 2012).

26. W. L. Barnes, A. Dereux, and T. W. Ebbesen, “Surface plasmon subwavelength optics,” Nature 424(6950), 824–830 (2003). [CrossRef]  

27. H. Ghaemi, T. Thio, D. Grupp, T. W. Ebbesen, and H. Lezec, “Surface plasmons enhance optical transmission through subwavelength holes,” Phys. Rev. B 58(11), 6779–6782 (1998). [CrossRef]  

28. K. Walls, Q. Chen, S. Collins, D. R. Cumming, and T. D. Drysdale, “Automated design, fabrication, and characterization of color matching plasmonic filters,” IEEE Photon. Technol. Lett. 24(7), 602–604 (2012). [CrossRef]  

29. Q. Chen and D. R. Cumming, “High transmission and low color cross-talk plasmonic color filters using triangular-lattice hole arrays in aluminum films,” Opt. Express 18(13), 14056–14062 (2010). [CrossRef]  

30. Q. Chen, C. Martin, and D. Cumming, “Transfer printing of nanoplasmonic devices onto flexible polymer substrates from a rigid stamp,” Plasmonics 7(4), 755–761 (2012). [CrossRef]  

31. Q. Chen, D. Das, D. Chitnis, K. Walls, T. Drysdale, S. Collins, and D. Cumming, “A cmos image sensor integrated with plasmonic colour filters,” Plasmonics 7(4), 695–699 (2012). [CrossRef]  

32. Q. Chen, D. Chitnis, K. Walls, T. D. Drysdale, S. Collins, and D. R. Cumming, “Cmos photodetectors integrated with plasmonic color filters,” IEEE Photon. Technol. Lett. 24(3), 197–199 (2012). [CrossRef]  

33. X. Wu and G. Zhai, “On sparse representations of color images,” in Proceedings of IEEE International Conference on Image Processing (IEEE, 2011), pp. 1229–1232.

34. D. Gao, X. Wu, G. Shi, and L. Zhang, “Color demosaicking with an image formation model and adaptive pca,” J. Vis. Commun. Image Represent. 23(7), 1019–1030 (2012). [CrossRef]  

35. J. A. Tropp, A. C. Gilbert, and M. J. Strauss, “Algorithms for simultaneous sparse approximation,” Signal Process. 86(3), 572–588 (2006). [CrossRef]  

36. J. Mairal, F. Bach, J. Ponce, G. Sapiro, and A. Zisserman, “Non-local sparse models for image restoration,” in Proceedings of IEEE International Conference on Computer Vision (IEEE, 2009), pp. 2272–2279.

37. D. H. Marimont and B. A. Wandell, “Linear models of surface and illuminant spectra,” J. Opt. Soc. Am. A 9(11), 1905–1913 (1992). [CrossRef]  

38. M. Parmar, S. Lansel, and B. Wandell, “Spatio-spectral reconstruction of the multispectral datacube using sparse recovery,” in Proceedings of IEEE International Conference on Image Processing (IEEE, 2008), pp. 473–476.

39. J. F. Cai, E. J. Candès, and Z. Shen, “A singular value thresholding algorithm for matrix completion,” SIAM J. on Optim. 20(4), 1956–1982 (2010). [CrossRef]  

40. S. Gu, L. Zhang, W. Zuo, and X. Feng, “Weighted nuclear norm minimization with application to image denoising,” in Proceedings of IEEE Conference on Computer Vision and Pattern Recognition (IEEE, 2014), pp. 2862–2869.

41. C. Li, “Compressive sensing for 3d data processing tasks: applications, models and algorithms,” Ph.D. thesis, Rice University (2011).

42. F. Yasuma, T. Mitsunaga, D. Iso, and S. Nayar, “Generalized assorted pixel camera: Post-capture control of resolution, dynamic range and spectrum,” Tech. rep., Department of Computer Science, Columbia University (2008).

43. http://www.cs.columbia.edu/CAVE/databases/multispectral/.

Cited By

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

Alert me when this article is cited.


Figures (10)

Fig. 1.
Fig. 1. A multispectral filter array and the constituent narrowband filters with the peak wavelength evenly shifted.
Fig. 2.
Fig. 2. A typical multispectral imaging system of spatial multiplexing type.
Fig. 3.
Fig. 3. Concept schematic of nanostructured MSFA with SPP-based broadband filters.
Fig. 4.
Fig. 4. Nine SPP-based broadband filters.
Fig. 5.
Fig. 5. Color image versions of the test set in Table 1.
Fig. 6.
Fig. 6. MSE of each of the 9 spectral bands in the test multispectral image Fig. 5(a).
Fig. 7.
Fig. 7. Recovered spectra of two selected patches in the test multispectral image Fig. 5(d).
Fig. 8.
Fig. 8. Recovered spectra of two selected patches in the test multispectral image Fig. 5(h).
Fig. 9.
Fig. 9. Recovered spectra of two selected patches in the test multispectral image Fig. 5(i).
Fig. 10.
Fig. 10. Recovered spectra of two selected patches in the test multispectral image Fig. 5(k).

Tables (3)

Tables Icon

Table 1. PSNR results of the test set in Fig. 5.

Tables Icon

Table 2. Averaged mutual correlations of the narrow band MSFA method and the proposed random MSFA method (no filter model errors, and σ = 5 ) on the test set in Fig. 5.

Tables Icon

Table 3. MSE results for all 16 images in Fig. 5.

Equations (14)

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

y i = Λ γ i ( λ ) l ( λ ) x i ( λ ) d λ
y i = k = 1 K Λ k γ i ( λ ) l ( λ ) x i ( λ ) d λ = k = 1 K x k , i Λ k γ i ( λ ) l ( λ ) d λ ,
a k , i = Λ k γ i ( λ ) l ( λ ) d λ ,
y i = k = 1 K a k , i x k , i .
λ p e a k = 2 π ϵ m ϵ d ϵ m + ϵ d ( k | | cos ϕ + i 2 π a x + j 2 π a x ) 2 + ( k | | sin ϕ i 2 π a y + j 2 π a y ) 2 ,
x k , i = τ k , i cos α i + κ k , i η cos s β i ,
x k = ( x k , 1 , x k , 2 , , x k , N ) T = τ k ( g 1 , g 2 , , g N ) T + κ k ( h 1 , h 2 , , h N ) T = τ k g + κ k h ,
D ( x ) = [ 2 x 1 , 2 x 2 , , 2 x K ]
D ( x ) 0 , 2 = # { i | r i 2 0 , i = 1 , 2 , , N } ,
D ( x ) 1 , 2 = i = 1 N r i 2 .
min     D ( x ) 1 , 2 s . t .     y i = k = 1 K a k , i x k , i ,         i = 1 , 2 , , N .
min     D ( x ) 1 , 2 s . t .     y i = k = 1 K a k , i x k , i ,         i = 1 , 2 , , N , R a n k ( x ) R ,
min     D ( x ) 1 , 2 + ω x s . t .     y i = k = 1 K a k , i x k , i ,         i = 1 , 2 , , N .
ρ = 1 N i = 1 N k = 1 K x k , i y k , i k = 1 K x k , i 2 k = 1 K y k , i 2 .
Select as filters


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