Multi-dimensional optical imaging systems exploit different degrees of freedom of scattering photons from an object scene, such as polarization, depth, and spectrum, to reveal different information [1–3]. For instance, spectral characteristics reflect the elemental composition, while the polarimetric characteristics contain the surface’s roughness and conductance [4]. The information can be helpful for object inspection and classification in remote sensing and industry applications. However, most existing approaches require mechanical moving parts or multiple modulation processes (e.g., polarizers of different orientations, diffractive gratings, or spectrum filters), which leads to long acquisition time, large volume size, high system complexity, or low sampling resolution.

In recent years, lensless imaging systems have gradually unfolded their advantages compared with traditional lens-based imaging systems. Unlike the point-to-point imaging manner in the latter systems, lensless imaging systems, acting as new paradigms in imaging, are replacing the conventional lenses with encoding masks and directly recording the coded pattern of an object on the sensor. Then after the post-process reconstruction, the information of the object can be recovered. Benefiting from this typical architecture, lensless imaging systems are usually more flexible, light weight, and with less cost than traditional lens-based imaging systems. It has been demonstrated that lensless imaging systems can be used in super-resolution imaging [5,6], three-dimensional imaging [7,8], multispectral [9] and hyperspectral imaging [10], and so on. Diffuser-based scattering imaging has drawn attention since much optical field information can be retrieved during the scattering process. However, there will be a random speckle in the transmittance of conventional diffusers [11], and their optical properties may be unstable with time [12], which makes the scattering effect of the conventional diffuser unpredictable and unrepeatable. Therefore, the time-consuming and tedious characterizations are unavoidable before utilization. Furthermore, due to the limited optical field engineering ability, it is also a great challenge to realize multiple-dimensional imaging for conventional lensless imaging systems owing to the wavelength and polarization insensitivity.

Digital optics could play an important role to conquer these challenges. Digital optics is a new concept that leverages discrete micro-/nano-optical elements to realize optical field manipulation with higher flexibility and resolution. As a core of digital optics, metasurfaces are planar devices comprising arrays of subwavelength meta-atoms. By locally tailoring the geometries of each meta-atom, the metasurface can manipulate the phase, amplitude, and polarization at will [13–24]. The metasurface-based devices have been successfully demonstrated and utilized in numerous applications, e.g., light-field imaging [16], depth sensing [25], planar synthetic aperture [26], polarization detection and imaging [27,28], multispectral imaging [29], and wide field-of-view imaging and detection [18,30]. The concept of the metasurface diffuser has been well established lately to achieve high-resolution bio-imaging [31] and complex optical field imaging [32]. Compared with the conventional diffuser, the predesigned metasurface diffuser significantly simplifies the characterization procedure and shows stable and reliable optical properties. More importantly, by exploiting sensitivities of wavelength and polarization at subwavelength scales, the metasurface has the potential to realize multi-dimensional imaging [33]. Nevertheless, preparing the metasurfaces comprising great nano-pillars or nano-holes is facing enormous challenges, especially in large-area manufacturing. Alternatively, a common and straightforward approach is to develop the geometric phase elements with controllable liquid crystal (LC) orientations [34–38]. On the one hand, the gradual maturity of LC production lines provides low-cost and large-scale production. On the other hand, compared to conventional diffraction optical elements, LCs also revealed unparalleled superiority in terms of operation efficiency, processing difficulty [39], and wavelength-polarization sensitivity.

In this paper, we propose and demonstrate the concept of multi-dimensional computational imaging (MCI) by combining the principles of both lensless computational imaging and metasurface optics. A flat LC-based diffuser fabricated through a standard photoalignment technology and a digital micro-mirror device (DMD) is used to encode spatial–spectral–polarization five-dimensional (5D) object information. The generated point spread functions (PSFs) exhibit linear translation invariance within the memory effect range, which promises the post-processing algorithm to recover two-dimensional information on the $x\text{-}y$ plane. The variation of depth along the $z$ axis will result in phase delays, leading to different PSFs with low correlation, which makes depth information able to be retrieved. Owing to the spin-reversed wavefront coding and dispersive diffraction character of the LC-based geometric phase diffuser, PSFs among different chirality and wavelength are highly irrelevant, benefiting to reveal polarization and spectrum information. Therefore, MCI can be easily achieved by the conventional deconvolution technique. Our proposal has been demonstrated in the visible band, with the merits of snapshot acquisition, multi-dimensional perception ability, simple optical configuration, compact device size, and high fabrication efficiency.

The designed metasurface diffuser is composed of anisotropic liquid crystal molecules (LCMs), based on the Pancharatnam–Berry (PB) phase [40,41], which is also called the geometric phase. The transmission property of the LCMs with a fast axis and slow axis along the $u$ and $v$ directions can be calculated with the Jones matrix, which can be written as $$J_{uv}=\left[\begin{array}{cc}{t}_u& 0\\ 0& t_v\end{array}\right],$$(1)where $t_u$ and $t_v$ are the transmission coefficients along the two main axis directions. Assume the fast axis of the LCM has an orientation angle of $\theta$ with respect to the $x$ axis, the Jones matrix can be represented as $$\begin{gathered} J_{xy}=M·J_{uv}·M^{-1} \\ =\left[\begin{array}{cc}\cos \theta & -\sin \theta \\ \sin \theta & \cos \theta \end{array}\right]\left[\begin{array}{cc}{t}_u& 0\\ 0& t_v\end{array}\right]\left[\begin{array}{cc}\cos \theta & \sin \theta \\ -\sin \theta & \cos \theta \end{array}\right]. \end{gathered}$$(2)

For circularly polarized (CP) incidence, the output field after passing through the LCMs can be calculated as $$\left[\begin{array}{c}{E}_{x\mathrm{out}}\\ E_{y\mathrm{out}}\end{array}\right]=\frac{J_{xy}}{\sqrt{2}}\left[\begin{array}{c}1\\ i\sigma \end{array}\right]=\frac{1}{2\sqrt{2}}\{(t_u+t_v)[\begin{array}{c}1\\ i\sigma \end{array}]+(t_u-t_v)\exp (2i\sigma \theta )[\begin{array}{c}1\\ -i\sigma \end{array}\left]\right\},$$(3)where $\sigma =\pm 1$, corresponding to left-circularly polarized (LCP) and right-circularly polarized (RCP) light, respectively. The second term of Eq. (3) indicates the CP incidence is converted to orthogonal polarization with an additional phase shift equal to twice the orientation angle as $\mathrm{\Delta }\mathrm{\Phi }=2\sigma \theta$. Theoretically, by designing the orientation angle of LCMs between 0 and $\pi$, a $2\pi$ phase-change coverage could be achieved despite wavelength, which leads to the capability to customize the unique phase distribution.

Based on the PB phase, we designed an LC metasurface diffuser with a random phase distribution due to the random wavefront achieving better-reconstructed results and showing more sensitivity to multi-dimensional information during the post-processing algorithm (see Appendix A for details). The designed phase profile was discretized into a four-order phase (i.e., 0, $\pi /2$, $\pi$, $3\pi /2$), which is enough for the spin-reversed wavefront coding. Figure 1(a) shows the designed phase profile of the LC metasurface diffuser, 10.3 mm long and 5.8 mm wide. Figures 1(b) and 1(c) indicate the local phase distribution with four-order phase and the LCMs with four different orientations correspondingly. The LC metasurface diffuser was then fabricated through a DMD-based photoalignment technology. The polarized optical microscope (POM) image of the metasurface diffuser is shown in Fig. 1(d), and Fig. 1(e) shows the local magnification images of Fig. 1(d) under different polarization states.

The functionalities and performances of an imaging system can be quantified by calculating its PSF. Our crucial intuition is to build a diffuser whose PSFs change with chirality, wavelengths, and position so that the associated information can be recorded into a snapshot speckle, and the polarization, spectral, and spatial information can be recovered through the conventional deconvolution technique, shown in Fig. 2. Suppose an object point is at a distance of $z_1$ from the diffuser. The optical field passes through the diffuser resulting in phase changes. Subsequently, the field propagates to the sensor at the distance $z_2$, resulting in the PSF $P_{\lambda ,\sigma ,z}$, which is the squared magnitude of the complex wave field at the sensor plane [42]: $$P_{\lambda ,\sigma ,z}={|F\{A·\exp [i({\mathrm{\Phi }}_{\mathrm{obj}}+{\mathrm{\Phi }}_{\mathrm{LC}}+{\mathrm{\Phi }}_{\text{sensor}})]\}|}^2,$$(4)where $F$ represents the Fourier transform, $A$ is the amplitude of the point source, and ${\mathrm{\Phi }}_{\{\mathrm{obj}/\mathrm{LC}/\text{sensor}\}}$ are the phase delays induced by propagating the object to the LC metasurface diffuser ${\mathrm{\Phi }}_{\mathrm{obj}}=k(z_1+\frac{x^2+y^2}{2z_1})$, the diffuser itself ${\mathrm{\Phi }}_{\mathrm{LC}}=2\sigma \theta$, and propagating from the LC metasurface diffuser to the sensor ${\mathrm{\Phi }}_{\text{sensor}}=k(z_2+\frac{x^2+y^2}{2z_2})$. $k$ is the wavenumber $2\pi /\lambda$, $\sigma$ is the chirality $\pm 1$, and $\theta$ is the orientation angle of LCM. Equation (4) shows that the PSF generated by the diffusion and diffraction of an LC metasurface diffuser depends on the spatial distribution $(x,y,z)$, wavelength $(\lambda )$, and chirality $(\sigma )$ of an object point (see Appendix B for details). These phase terms are summed up at k-space and converted to the spatial domain through Fourier transform. Therefore, this multiple-dimensional information can be reconstructed from PSF analysis in a single image, shown in Figs. 2(c) and 2(d). Information in the ($x,y$) plane can be recovered by deconvolution based on the shift-invariance of PSFs. Information recovery in ($z,\lambda ,\sigma$) spaces relies on the axial-spectral-polarization dependency of the diffuser.

Our Snapshot 5D Imaging system is mainly based on Fresnel diffraction and Fourier optics [43], and for an incoherent imaging system, it is described as $$I(x,y,z,\lambda ,\sigma )=O(x,y,z,\lambda ,\sigma )*\mathrm{PSF}(x,y,z,\lambda ,\sigma ),$$(5)where $*$ denotes convolution and $\mathrm{PSF}(x,y,z,\lambda ,\sigma )$ is the PSF. $O(x,y,z,\lambda ,\sigma )$ and $I(x,y,z,\lambda ,\sigma )$ are the objects and the images (also called the input and output signals). $(x,y,z)$ is the spatial distribution of an object, and $\lambda$ and $\sigma$ are the wavelength and chirality of the incident light. As the description in Section 2.B, in the memory effect range, the shift-invariance PSFs are axial–spectral–polarization dependent, making PSFs uncorrelated under different spatial–spectral–polarization conditions. It can be given by $$\mathrm{PSF}(z_1,{\lambda }_1,{\sigma }_1)\star \mathrm{PSF}(z_2,{\lambda }_2,{\sigma }_2)\approx \{\begin{array}{cc}0,& \mathrm{if}|z_1-z_2|>T_z\vee |{\lambda }_1-{\lambda }_2|>T_{\lambda }\vee {\sigma }_1\ne {\sigma }_2\\ \delta ,& \mathrm{if}|z_1-z_2|<T_z\wedge |{\lambda }_1-{\lambda }_2|<T_{\lambda }\wedge {\sigma }_1={\sigma }_2\end{array},$$(6)where $\star$ denotes the correlation operator, and $\delta$ denotes the impulse response function. $T$ represents the response threshold value of each dimension. The variation of $z$ and $\lambda$ will lead to additional phase delays ${\mathrm{\Phi }}_{\mathrm{\Delta }z}=k(\mathrm{\Delta }z+\frac{x^2+y^2}{2z}·\frac{\mathrm{\Delta }z}{z+\mathrm{\Delta }z})$ and ${\mathrm{\Phi }}_{\mathrm{\Delta }\lambda }=k·\frac{\mathrm{\Delta }\lambda }{\lambda +\mathrm{\Delta }\lambda }·r$, which changes the intensity profile on the detector, decreasing the correlation between PSFs. The polarization selective properties described as in Eq. (3) explain the different phase delay assigned to LCP and RCP incidences, providing significant difference among PSFs.

The speckle pattern with all the spatial distribution, wavelength, and chirality information of the object is the composite response of the LC metasurface diffuser. It can be expressed as $$I={\sum }_{x,y,z,\lambda ,\sigma }I(x,y,z,\lambda ,\sigma )={\sum }_{x,y,z,\lambda ,\sigma }[O(x,y,z,\lambda ,\sigma )*\mathrm{PSF}(x,y,z,\lambda ,\sigma )].$$(7)

Therefore, the object with interesting spatial–spectral–polarization information can be reconstructed from the image overlapping diverse factors by deconvoluting PSF with respect 5D property as follows: $$O(x_0,y_0,z_0,{\lambda }_0,{\sigma }_0)=\text{deconv}[I,\mathrm{PSF}(x_0,y_0,z_0,{\lambda }_0,{\sigma }_0)].$$(8)

We use Wiener deconvolution for image reconstruction, and it is expressed as follows [44]: $$\begin{gathered} O(x_0,y_0,z_0,{\lambda }_0,{\sigma }_0)=\text{deconv}[I,\mathrm{PSF}(x_0,y_0,z_0,{\lambda }_0,{\sigma }_0)] \\ ={\mathrm{FFT}}^{-1}\left\{\frac{\mathrm{FFT}(I)·\mathrm{FFT}{[\mathrm{PSF}(x_0,y_0,z_0,{\lambda }_0)]}^c}{{|\mathrm{FFT}[\mathrm{PSF}(x_0,y_0,z_0,{\lambda }_0,{\sigma }_0)]|}^2+\mathrm{SNR}(f)}\right\}, \end{gathered}$$(9)where $\mathrm{FFT}(·)$ and ${\mathrm{FFT}}^{-1}(·)$ denote the Fourier transform and its inverse, respectively. ${(·)}^c$ is the complex conjugate and $\mathrm{SNR}(f)$ is the signal-to-noise ratio. Using the spatial–spectral–polarization-dependent behavior of the LC metasurface diffuser described in Eq. (6), we can reconstruct the object of interest with high quality.

As an imaging system, the spatial resolution of MCI is measured first. The setup of MCI is shown in Fig. 3(a). A projector (Acer X118H) is used to generate object patterns. The magnification lens of the projector is removed, and an iris is set to filter out background light from the projector. The LC metasurface diffuser then modulates the light that comes from the projector, and a monochromatic camera (Daheng, MER-500-7UM) is used to capture the intensity profile. The exposure time to capture PSFs was set to 80 ms, and that to capture objects was set to 16 ms. The distance from the object patterns to the metasurface is 15 cm, and the camera is placed 2 cm behind the metasurface. To decide the resolution of this system, resolution chart patterns with different gaps, shown in Fig. 3(c), are loaded on the center area of the projector. The recovered images and their profiles are shown in Fig. 3(d), illustrating that the spatial resolution of this system is around 28 μm, as the peak–valley ratio is about 1.52 at this scale. The deconvolution algorithm takes 0.43 s to reconstruct information from the PSF and the corresponding speckle pattern pair. The axial resolution and spectral resolution are discussed in Appendix A, and more experimental results are exhibited in Appendix C.

The same setup shown in Fig. 3(a) is used to generate object patterns for multispectral imaging of MCI. By loading different patterns with different RGB values, objects’ spatial and spectral information could be changed conveniently. In practice, seven different capital letters with different colors are projected respectively. The speckle patterns captured by the camera are shown in Fig. 4(a) (Pseudocolore is applied for convenient exhibition). Then, a single central pixel with the aforesaid colors is lightened up respectively, and the captured speckles are taken as multispectral PSFs of MCI, as shown in Fig. 4(b). The reconstructed multispectral images are finally obtained by deconvolving different spectral PSFs with different speckle patterns, exhibited in Fig. 4(c). The ground truth patterns loaded on the projector are shown in Fig. 4(d). As expected, each speckle pattern is smoothly reconstructed. Figure 4(e) demonstrates the composite multispectral image superimposing with seven individual spectral images of Fig. 4(c).

For verifying the polarization selective ability of MCI, we leverage the setup shown in Figs. 5(a) and 5(b). Different from the setup in Fig. 3, two pairs of the polarizer and the quarter-wave plate were employed. The first pair of the polarizer and quarter-wave plate translates the incident light into CP light. By rotating the quarter-wave plate, LCP light or RCP light can be generated. Then the CP light transmits through the metasurface diffuser and converts incident CP light to the opposite chirality. The LC metasurface diffuser can obtain a high polarization conversion efficiency ($>95\%$) at the designed wavelength, which is about 532 nm, but the efficiency will decrease away from that wavelength. So the second pair of the polarizer and the quarter-wave plate is applied to filter out the co-polarized light remaining in transmitted light. In the experiment, two letters, “H” and “E,” with the same color, were loaded on the projector. Figure 5(c) shows the raw speckle patterns captured by the camera. The rotation arrows in the upper left indicate the polarization of each pattern (blue for RCP, red for LCP). Figure 5(d) shows the PSFs with different polarizations. Figures 5(e) and 5(f) show the reconstructed images by deconvolving the different polarization speckle patterns with different polarization PSFs. There are two rotation arrows in the upper left. The first indicates the polarization of the speckle pattern, and the second indicates the polarization of the PSF used for deconvolving [e.g., the upper picture in Fig. 5(f) has two arrows, blue on the left and red on the right, indicating this picture is the result of the RCP speckle pattern deconvolving with the LCP PSF]. According to these experimental observations, MCI shows high sensitivity to the polarization of incident objects, and objects can be recovered only when the speckle patterns and PSFs are at the same polarization.

To evaluate the performance of our MCI in spatial–spectral–polarization 5D $(x,y,z,\lambda ,\sigma )$ imaging, two objects with different 5D properties are projected and superimposed on the sensor. The experimental setup is shown in Fig. 6(a). The two patterns with different colors are loaded on the left and right areas of the projector. After transmitting through a diaphragm and a polarizer, the beam from the projector becomes polarized light. The left part of the beam is reflected by a mirror to another path. Beams in two paths are polarized into orthogonal circular polarizations while transmitting through quarter-wave plates with different fast axis orientations. Then two beams are superimposed by a beam combiner and modulated by the metasurface diffuser. Finally, the raw speckle pattern is recorded by a camera, shown in Fig. 6(b). Before recovering the actual pattern, a two-point source with different colors was loaded on the projector to generate corresponding 5D PSFs respectively, shown in the upper part of Figs. 6(c) and 6(d). The reconstructed images obtained by deconvolving the PSFs with speckle patterns in Fig. 6(b) were exhibited in the lower parts of Figs. 6(c) and 6(d). The experimental results shown above illustrate that the two objects with different 5D information can be reconstructed from the overlapped speckle pattern. PSFs with different 5D information can retrieve different objects in the superimposed speckle pattern. Note that the quality of reconstructed images in Figs. 6(c) and 6(d) is not as good as shown in Figs. 4(c) and 5(e), which is likely due to the polarization conversion efficiency of our metasurface diffuser. Although the metasurface converts polarization of incident light into the opposite handedness, some co-polarized light remains in transmitted light, causing minor deterioration for recognizing objects with orthogonal polarization.

We demonstrate a lensless snapshot 5D imaging system by employing a computational technique and the metasurface’s spatial–spectral–polarization sensitivity. The 5D information of the object is encoded by the metasurface as speckle patterns, which allow us to apply a deconvolution algorithm to reconstruct the 5D information of interest by corresponding PSFs. Our demonstrations present a compact and inexpensive technique for snapshot 5D imaging that might be promising for material classification and identification, biomedicine, and industry applications. As a general framework of imaging and detection, our proposal is promising to promote the next generation of engineering optics [45,46].

Acknowledgment. Y. Lei thanks Dr. Dongliang Tang for his help in LC metasurface fabrication and experimental design.