Open Access
Add to BibSonomy
Abstract
Analog optical computing based on metasurfaces has attracted much attention for achieving high-speed calculating without the electronic processing unit. Wavefront coding imaging systems involve the joint design of an encoded image-capturing module and decoding postprocessing algorithms to obtain a required final image. The decoding postprocessing algorithms, as a typical deconvolution computation, require an additional electronic processing unit to yield a high-quality decoded image. We demonstrate an analog optical deconvolution computing kernel based on nanoantennas metasurfaces for the postprocessing calculation of wavefront coding systems. Numerical simulations are presented to prove that the encoded point spread function can be refocused through a designed optical computing metasurface. The proposed approach opens an opportunity for real-time recovering images in wavefront coding optical systems.
© 2021 Optical Society of America under the terms of the OSA Open Access Publishing Agreement
1. Introduction
Wavefront coding imaging technology involves the joint design of an encoded image-capturing module and algorithms for encoded image restoration to obtain a required final image. For example, a cubic phase mask (CPM) is introduced into the pupil plane to encode wavefront so that the feature of low sensitive to defocus can be produced in optical coded section. The decoded images with larger depth of focus [1,2] can be obtained after the decoded algorithms are employed to magnify the prominent feature. Furthermore, a 3D image can be acquired when the Complementary Kernel Matching (CKM) algorithm is applied in the same wavefront coding system [3]. Encoded image restoration, a considerable influential and fundamental technology for the wavefront coding system, is carried out by deconvolution calculation based on the known encoded information. Deconvolution calculation can be performed through either electronic or optical methods. Although digital computers can deal with complex data, there are limitations in real-time applications due to the speed of computation. Analog optical computing (AOP) [4–12], which uses optical signals to process information directly, provides an alternative way to manipulate large amounts of data and real-time calculation.
Traditionally, analog optical deconvolution computing has been accomplished by using the Fourier approach relying on lenses and filter systems [13]. However, the use of numerous conventional lenses leads to a heavy and complex optical system, not conducive to the development of integrated devices. One option for significantly simplifying the size of an optical system is to utilize metasurfaces [14–21]. It is now well established from a variety of studies that analog optical computing can be performed by metasurfaces, which can be realized by the Pancharatnam–Berry phase [10], the guide mode resonance [21], the surface plasmon [6], etc. Although the existing literature on analog optical computing based on metasurfaces is extensive, more recent attention has focused on the provision of spatial differentiation [20–26]. The calculation of spatial differentiation is usually used for image edge enhancement. However, there are relatively few studies in analog optical deconvolution computing based on metasurfaces.
In this paper, we demonstrate a real-time deconvolution calculation based on metasurfaces. The whole computing is finished after light transmits through two graded index (GRIN) lenses, and one metasurface. The unit cell of the metasurface is made of metal antennas [27,28], which can tune the desired amplitude and phase at the same time. In the simulation, we choose the point spread function (PSF) of a wavefront coding system with CPM as an input function. The recovered airy spot after deconvolution calculation is well-matched with the ideal one. To sum up, the use of the nanoantennas metasurfaces to achieve deconvolution calculation allows for a significant size reduction compared to traditional systems. In addition, this technology opens new doors for real-time encoded image restoration in applications involving wavefront coding.
2. Design structure of the deconvolution computing
2.1 Principle of analog optical deconvolution computing
In liner space invariant system, the input function $f({x,y} )$ and corresponding output function $g({x,y} )$ have a convolution relationship [29] described by
where $h({x,y} )$ is a two-dimensional (2D) impulse response for the purposed analog optical computing system. Equation (1) can be written in the Fourier domain as where $G({{k_x},{k_y}} )$, $H({{k_x},{k_y}} )$, and $F({{k_x},{k_y}} )$ are the Fourier transformation (FT) of their counterparts in Eq. (1), respectively. And $({{k_x},{k_y}} )$ means the 2D frequency variables in the spatial Fourier domain. According to Eq. (1) and Eq. (2), the general AOP system is shown in Fig. 1.We realize optical deconvolution through the same method of Fourier transformation. Figure 2 presents an overview of the whole structure. In our designed system, $f({x,y} )$ and $g({x,y} )$ are cross-polarized monochromatic light, i.e., ${E_x}({x,y} )$ and ${E_y}({x,y} )$, respectively. The transfer function $H({{k_x},{k_y}} )$, also the transmission coefficient $T({{k_x},{k_y}} )$, is associated with the priori information, i.e., the PSF of the wavefront coding system.
The method for both one dimensional (1D) and two dimensional (2D) optical deconvolution computing is similar. Thus, we discuss the 1D system as an example to describe how to construct analog optical deconvolution computing in the following subsections. 2D optical deconvolution computing can be done by following the same steps in 1D deconvolution. Certainly, both the result of 1D and 2D deconvolution computing are discussed in the section of simulation consequence.
2.2 Fourier Transform Block and Inverse Fourier Transform Block
The normal FT block for the 1D analog optical computing is a GRIN lens [28,30]. Its permittivity varies with the coordinate:
Besides the GRIN lens, metalens can be the choice of FT blocks, and is a more effective method to minimize the footprint [20]. However, the GRIN lens has advantages in other ways, like being cheaper and easier to experiment with. After comprehensive consideration, we use GRIN lens instead of metalens in the simulation.
2.3 Control of amplitude and phase with nanoantennas
The transfer function block is carried out by metasurfaces, whose unit cell is a gold V-shaped antenna on the glass substrate. The V-shaped antenna is a simple design that can be easily adjusted to control the intensity and phase of scattering light. Individual one is characterized by the thickness of antenna h, opening angle of the V-shaped antenna $\psi $, the total length of the antenna L, the width of the antenna d, the thickness of glass substrate t, the height of output plane z, and the orientation angle $\theta $ between the x-axis and the bisector of the V-shaped antenna.
In the simulation, we remain the thickness of antenna $h = 100\textrm{nm}$, the width of the antenna $d = 20\textrm{nm}$, the thickness of glass substrate $t = 500\textrm{nm}$, the height of output plane $z = 500\textrm{nm}$ and the orientation angle $\theta \textrm{ = 4}{\textrm{5}^ \circ }$ unchanged. Then we sweep the opening angle of the V-shaped antenna $\psi $ and the total length of the antenna L for desired amplitude and phase control. The characteristics defined in the preliminary scanning process are summarized in Fig. 3(a). In addition, the metasurface pattern is devised for liner incident polarization along the x direction, and y polarized transmitted light. The single unit cell and the whole metasurface are calculated by the finite element method. At a wavelength of 633 nm, the results obtained from the scanning procedure are presented in Fig. 3(c) and Fig. 3(d).
Fig. 3. The sweep characteristic (a) of the unit cell and phase difference (b) for orientation angle. Simulated (c) amplitude and (d) phase for the proposed unit cell with period $p = 300\textrm{nm}$, orientation angle $\theta = {45^ \circ }$.
It is apparent from Fig. 3(d) that all phases $2\pi $ are hard to be found at one orientation angle $\theta $. Here we exploit the Pancharatnam-Berry (PB) phase concept. When the antenna orients at $\theta = {45^ \circ }$ and $\theta ={-} {45^ \circ }$, the $\pi $ phase difference between $E_y^t$ and $- E_y^t$ are highlighted in Fig. 3(b). In addition, the change in amplitude at two angles that differs by 90 degrees is small, within the tolerance of the desired amplitude control. Thus, we can keep the amplitude result in Fig. 3(c) valid at an orientation angle $\theta ={-} {45^ \circ }$ while getting more phase modulation. After collecting all scanning results, we matched the unit cell for both amplitude and phase of transfer function $H({{k_x},{k_y}} )$.
For deconvolution computing, transfer function $H({{k_x},{k_y}} )$ or transmission coefficient $T({{k_x},{k_y}} )$ can be written as
Equation (4) illustrates the relationship between the transfer function and the input function. Transfer functions always depend on the FT of PSFs in different systems. In addition, coordinate in k space here refer to the Cartesian coordinate of metasurfaces because FT happens in the spatial domain.
When the input function is the PSF of the target optical system, the output function is an airy spot after deconvolution computing. Equally, when the target wavefront coding system is used for imaging, the input function is the encoded image. Therefore, the output function can be the decoded image after deconvolution computing.
3. Simulation results
1D deconvolution computing has been simulated to assess the feasibility of the above design. Before undertaking the simulation, the wavefront coding system with CPM is picked out as a prepositive system before analog optical deconvolution computing. Thus, the input function for analog optical deconvolution computing is the PSF of the wavefront coding system with CPM.
According to Fourier optics, the PSF is the Fourier transform of the pupil function in the coherent system. For the incoherent system, the PSF is the absolute square of the PSF in the coherent system. Consequently, for the wavefront coding system under incoherent conditions, the PSF can be written as:
where $P(x )= \exp ({ik\alpha {x^3}} )$ is the pupil function of wavefront coding with CPM. Equation (5) can be further described asFig. 4. The simulation process and result for 1D analog optical deconvolution computing. (a) represents for the first FT block, (c) for metasurfaces, and (e) for the last FT block. (b) is the input function of PSF for CPM. (d) is the transfer function of metasurfaces. (f) is the final result for analog optical deconvolution computing.
The reason for normalizing the input function is to calculate the transmission efficiency easily. From Fig. 4(b), the PSF is zeros out of the radius $r = 3\mathrm{\mu }\textrm{m}$. Consequently, the simulation area can be further decreased, ranging from $- 3\mathrm{\mu }\textrm{m}$ to $3\mathrm{\mu }\textrm{m}$. This width is also the width of the GRIN lens. After the input function is determined, the simulation is checked by exploiting the structure in Fig. 2. Figure 4(a) provides the simulation date of the transmission procedure in the first FT block. The most obvious finding to emerge from the analysis is that light inclines during the whole FT process, making the phase of FT result $FT[{{E_x}(x )} ]$ asymmetric. This feature identified in the former judgment can be seen in the blue solid line in Fig. 4(d), which displays an ideal transfer function for our metasurfaces. This ideal transfer function is calculated and normalized by Eq. (4), providing the guide for matching the unit cell.
The matching procedure is easy to undergo. For example, we first locate the coordinate of each meta-atom. We use the centre of meta-atoms in x-coordinate as match points. Because the period of meta-atoms is 300 nm, we start the matching procedure in $x = 150\textrm{nm}$. The ideal amplitude and phase value are shown in Fig. 5(a). The desired phase is used to match the phase in Fig. 4(d), finding the approximate value with the error of 0.05. Through this step, we can minimize the database to easily find the structure of the purposed meta-atom. Then the desired amplitude is used to match the transmission in Fig. 4(c) within the minimized date set. At last, we can find the structure of the purposed meta-atom with the antenna length $L = 275\textrm{nm}$, opening angle $\psi = {150^ \circ }$ and orientation angle $\theta = {45^ \circ }$ in the location $x = 150\textrm{nm}$. The transmission and phase of this meta-atom, labelled in Fig. 5(b) and Fig. 5(c) by the red cycle, is the best match for the ideal value in Fig. 5(a).
Fig. 5. The example for matching the purposed meta-atom in $x = 150\textrm{nm}$. (a) is the ideal amplitude and phase, which is labelled in the red cycle. The red cycle in (b) is the matched amplitude and in (c) is the matched phase. The red cycle in (b) and (c) also decides the structure of the purposed meta-atom with $L = 275\textrm{nm}$ and $\psi = {150^ \circ }$.
Other meta-atoms’ matching procedure is like the above example. Once the metasurface is permuted, it is necessary to test whether the transmission result agrees with the ideal transfer function. When the input function $f(x )= 1[{\textrm{V/m}} ]$ is imported into the metasurface, the normalized simulated transfer function is shown in Fig. 6. The outcome, standing out in Fig. 6(a), is that the simulated intensity of the transfer function has a similar tendency with the ideal one. The different values that emerge between the two transfer functions may influence the performance of the recovered result. The result, as shown in Fig. 6(b), indicates that the phase of the simulated transfer function is highly consistent in the fitness. An explanation for the better simulation result in phase distribution is that we prioritize to match phase and then look for the amplitude closest to the target. In addition, the database of the unit cell is not large enough to hold all target values, which also causes the imperfect fitness.
Following the confirmation of the metasurface, optical deconvolution computing continues to do computing in the transfer function block and IFT block. Figure 4(c) provides the transmission procedure in the metasurface. The most interesting aspect of this figure is that the output function $FT[{E_y^t(x )} ]$ is nearly unchanged. From Fig. 4(e) we can see the transmission procedure in the IFT block. Light gradually becomes an airy spot, refocusing in the output plane. This is a rather remarkable outcome. We need to emphasize that all the fields of transmission procedure in Fig. 4(a), Fig. 4(c) and Fig. 4(e) are calculated separately due to the low computing power (CPU: i7-8700 K, RAM: 16GB).
The simulated result is shown in Fig. 4(f). The blue dotted line is the ideal airy spot, which is calculated in the FT block for the input function $f(x )= 1[{\textrm{V/m}} ]$. The red solid line is the simulation result for output light $FT[{E_y^t(x )} ]$. Both two lines are normalized, and the maximum intensity of result ${|{E_y^t(x )} |^2}$ is near ${10^{ - 4}}$.
Here, we can define the transmission efficiency in one block as the ratio of the output function’s maximum intensity to the input function’s maximum intensity. Therefore, the transmission efficiency in FT/IFT block is $\textrm{6}\textrm{.25\%}$. In addition, the transmission efficiency in the metasurface is $\textrm{2}\textrm{.89\%}$, proportional to the minimum value of the normalized transfer function. Then, the whole transmission efficiency is $\textrm{0}\textrm{.01\%}$ for our analog optical deconvolution computing. This poor transmission efficiency is related to be related to the low transmission efficiency in the transfer function block and FT/IFT block. This is an important issue for future experimental research. In other words, the energy of output light is so low that input light requires very high energy for an experiment in the future.
A closer inspection of Fig. 4(f) shows the simulation result is in good agreement with the desired airy spot. Excepted that, we can also find the two main differences between the ideal and the simulation result for deconvolution. One is a little shift of the spot’s peak, and another is the noise in the sidelobe. These drawbacks may bring aberration in recovered encoded images.
The final stage of the study comprises 2D deconvolution computing. Before evaluating the performance of 2D analog optical deconvolution computing, we also need to compare the simulated transfer function of the 2D metasurface with the ideal one. The simulated amplitude or phase distributions of the transfer function and the ideal results show a clear contrast in Fig. 7. Both the simulated amplitude distribution in Fig. 7(c) and phase distribution in Fig. 7(d) show a similar tendency with the ideal result in Fig. 7(a) and Fig. 7(b). However, the phase of the simulated transfer function in the 2D metasurface has a decreased performance, when comparing to the transfer function in the 1D metasurface. The unfitness of meta-atoms with target values will increase the influence in 2D computing, due to more meta-atoms in 2D metasurfaces than 1D metasurfaces. And this decreased performance in the 2D metasurface will enlarge the noise in the last recovered result.
Fig. 7. The contrast between ideal transfer function and simulated one. (a) is the 2D amplitude distribution of ideal transfer function. (b) is 2D phase distribution of ideal transfer function. and phase (d) distribution for the 2D metasurface.
Then, 2D analog optical deconvolution computing is used to further validate the usefulness of our method. The simulation process is consistent with the 1D simulation. And the results obtained from the simulation of 2D deconvolution computing are presented in Fig. 8. No significant difference between the 1D and 2D deconvolution computing is evident. The result in Fig. 8(c) and Fig. 8(d) is like a 1D recovered airy spot, fitted well with the desired one. The only inconsiderable difference is that the sidelobe noise in 2D computing is larger than 1D computing.
Fig. 8. The simulation result for 2D analog optical deconvolution computing. (a) is the 2D input function of PSF for CPM. (b) is the final result for analog optical deconvolution computing. The final airy spot is analyzed in the x axis (c) and y axis (d).
The differences between ideal airy disk and simulated result either in the 1D computing or in the 2D computing are caused by a little bit of unfitness between disperse unit cell and continuous transfer function. The most effective method to solve these problems is to enlarge the database by sweeping more characteristics of the unit cell. With the successive increases in the database, the simulated transfer function is further unanimous with the ideal transfer function. Therefore, the output function of metasurfaces will be better and the last recovered airy spot will have fewer drawbacks.
There are always aberrations existing in the imaging system, so we need to evaluate the aberration influence for analog optical deconvolution computing. Wavefront coding system with CPM is particularly sensitive to defocus. When the defocus aberration is present in this system, the PSF will have a shift in the image plane. As a result, the peak position of the PSF deviates from the optical axis [31]. Therefore, we use the offset PSF as an input function, which is shown in Fig. 9(a) and Fig. 9(c). In the same deconvolution metasurface, the 1D and 2D recovered airy spot is present in Fig. 9(b) and 9(d) respectively. The peak position of the recovered airy spot has a large shift from the optical axis. If the defocused coding images are recovered through analog optical deconvolution computing, the recovered images will have shift stripes, just like the recovered result through the computer. The problem of aberration influence may be solved by using optical and electronic computing together.
Fig. 9. Defocus aberration influence for analog optical deconvolution computing. (a) is 1D and (c) is 2D offset PSF in wavefront coding system with CPM. (b) is the 1D and (d) is 2D simulated recovered airy spot.
To further explain the broad feasibility for analog optical deconvolution computing, we choose another input function to compare with the former result in one dimension. The input function can be defined as:
where $a = 2.1/\lambda = 3.3\mathrm{\mu }{\textrm{m}^{ - 1}}$ and $b = {\lambda ^2}/0.9 = 0.445\mathrm{\mu }{\textrm{m}^2}$. This function is like a sinusoidal, which is shown in Fig. 10(a). The result after 1D analog optical deconvolution computing is present in Fig. 10(b). This simulated recovered airy spot offset more distance from the center spot than the recovered result in Fig. 4(f). And the noise for the recovered result is larger. But this result shows the usefulness for analog optical deconvolution computing in other systems, not just the wavefront coding system.Fig. 10. 1D analog optical computing for other input functions. (a) is the input function and (b) is the simulated and ideal recovered airy spot.
4. Conclusion
The present research aims to seek the application of analog optical computing in encoded image restoration, the postprocessing calculation of wavefront coding systems. In this work, our numerical simulations verify the feasibility of optical deconvolution computing, fulfilled by metasurfaces and GRIN lenses, for a wavefront coding system with CPM. The encoded point spread function can be recovered to an ideal airy disk through the designed optical computing metasurfaces. This research of optically encoded image restoration adds to the rapidly expanding field of real-time wavefront coding imaging. However, a limitation of this study is that the transmission efficiency of the system is $\textrm{0}\textrm{.01\%}$. Some methods can effectively deal with the problem of low efficiency. For example, leveraging low-loss, CMOS-compatible dielectric metasurfaces associated with agile amplitude or phase responses achieved from co-polarized light can improve the overall efficiency [20]. The issue of poor transmission efficiency is an intriguing one that could be usefully explored in further research. The performance of analog optical deconvolution computing can also be improved in processing more complicated data sets in the future. Through the design of metasurfaces, like achromatic metalens [32] and multifunction metasurfaces [33,34], RGB images or more imaging processing [35,36] can be dealt with by optical signals.
Funding
National Natural Science Foundation of China (61675002).
Acknowledgments
This work is supported by the National Natural Science Foundation of China (NSFC) (No.61675002).
Disclosures
The authors declare no conflicts of interest.
Data availability
Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.
References
1. E. R. Dowski and W. T. Cathey, “Extended depth of field through wave-front coding,” Appl. Opt. 34(11), 1859–1866 (1995). [CrossRef]
2. H. Zhao, H. Yi, and Y. Li, “Extending the depth-of-field of incoherent imaging system through wave-front coding by using composite phase masks,” in MIPPR 2011: Pattern Recognition and Computer Vision (International Society for Optics and Photonics, 2011), p. 80040F.
3. P. Zammit, G. Carles, and A. R. Harvey, “3D imaging and ranging in a snapshot,” in Optical Systems Design 2015: Computational Optics (International Society for Optics and Photonics, 2015), p. 963004.
4. L. L. Doskolovich, D. A. Bykov, E. A. Bezus, and V. A. Soifer, “Spatial differentiation of optical beams using phase-shifted Bragg grating,” Opt. Lett. 39(5), 1278–1281 (2014). [CrossRef]
5. S. AbdollahRamezani, K. Arik, A. Khavasi, and Z. Kavehvash, “Analog computing using graphene-based metalines,” Opt. Lett. 40(22), 5239–5242 (2015). [CrossRef]
6. A. Pors, M. G. Nielsen, and S. I. Bozhevolnyi, “Analog computing using reflective plasmonic metasurfaces,” Nano Lett. 15(1), 791–797 (2015). [CrossRef]
7. A. Chizari, S. Abdollahramezani, M. V. Jamali, and J. A. Salehi, “Analog optical computing based on a dielectric meta-reflect array,” Opt. Lett. 41(15), 3451–3454 (2016). [CrossRef]
8. A. Youssefi, F. Zangeneh-Nejad, S. Abdollahramezani, and A. Khavasi, “Analog computing by Brewster effect,” Opt. Lett. 41(15), 3467–3470 (2016). [CrossRef]
9. T. Zhu, Y. Zhou, Y. Lou, H. Ye, M. Qiu, Z. Ruan, and S. Fan, “Plasmonic computing of spatial differentiation,” Nat. Commun. 8(1), 1–6 (2017). [CrossRef]
10. T. Zhu, Y. Lou, Y. Zhou, J. Zhang, J. Huang, Y. Li, H. Luo, S. Wen, S. Zhu, and Q. Gong, “Generalized spatial differentiation from the spin hall effect of light and its application in image processing of edge detection,” Phys. Rev. Appl. 11(3), 034043 (2019). [CrossRef]
11. F. Zangeneh-Nejad, D. L. Sounas, A. Alù, and R. Fleury, “Analogue computing with metamaterials,” Nat. Rev. Mater. 6(3), 207–225 (2021). [CrossRef]
12. S. Abdollahramezani, O. Hemmatyar, and A. Adibi, “Meta-optics for spatial optical analog computing,” Nanophotonics 9(13), 4075–4095 (2020). [CrossRef]
13. S. H. Lee, “Optical information processing. Fundamentals,” (1981).
14. M. Khorasaninejad, A. Y. Zhu, C. Roques-Carmes, W. T. Chen, J. Oh, I. Mishra, R. C. Devlin, and F. Capasso, “Polarization-insensitive metalenses at visible wavelengths,” Nano Lett. 16(11), 7229–7234 (2016). [CrossRef]
15. W. T. Chen, A. Y. Zhu, J. Sisler, Y.-W. Huang, K. M. Yousef, E. Lee, C.-W. Qiu, and F. Capasso, “Broadband achromatic metasurface-refractive optics,” Nano Lett. 18(12), 7801–7808 (2018). [CrossRef]
16. Q. Jiang, L. Cao, L. Huang, Z. He, and G. Jin, “A complex-amplitude hologram using an ultra-thin dielectric metasurface,” Nanoscale 12(47), 24162–24168 (2020). [CrossRef]
17. H. Ren, X. Fang, J. Jang, J. Bürger, J. Rho, and S. A. Maier, “Complex-amplitude metasurface-based orbital angular momentum holography in momentum space,” Nat. Nanotechnol. 15(11), 948–955 (2020). [CrossRef]
18. Z.-L. Deng and G. Li, “Metasurface optical holography,” Mater. Today Phys. 3, 16–32 (2017). [CrossRef]
19. Y. Gao, Y. Fan, Y. Wang, W. Yang, Q. Song, and S. Xiao, “Nonlinear holographic all-dielectric metasurfaces,” Nano Lett. 18(12), 8054–8061 (2018). [CrossRef]
20. S. Abdollahramezani, A. Chizari, A. E. Dorche, M. V. Jamali, and J. A. Salehi, “Dielectric metasurfaces solve differential and integro-differential equations,” Opt. Lett. 42(7), 1197–1200 (2017). [CrossRef]
21. A. Saba, M. R. Tavakol, P. Karimi-Khoozani, and A. Khavasi, “Two-dimensional edge detection by guided mode resonant metasurface,” IEEE Photonics Technol. Lett. 30(9), 853–856 (2018). [CrossRef]
22. C. Guo and S. Fan, “Optical image processing using photonic crystal slab,” in Semiconductors and Semimetals (Elsevier, 2019), pp. 93–114.
23. J. Zhou, H. Qian, C.-F. Chen, J. Zhao, G. Li, Q. Wu, H. Luo, S. Wen, and Z. Liu, “Optical edge detection based on high-efficiency dielectric metasurface,” Proc. Natl. Acad. Sci. 116(23), 11137–11140 (2019). [CrossRef]
24. L. Wan, D. Pan, S. Yang, W. Zhang, A. A. Potapov, X. Wu, W. Liu, T. Feng, and Z. Li, “Optical analog computing of spatial differentiation and edge detection with dielectric metasurfaces,” Opt. Lett. 45(7), 2070–2073 (2020). [CrossRef]
25. Y. Zhou, H. Zheng, I. I. Kravchenko, and J. Valentine, “Flat optics for image differentiation,” Nat. Photonics 14(5), 316–323 (2020). [CrossRef]
26. M. Farmahini-Farahani, J. Cheng, and H. Mosallaei, “Metasurfaces nanoantennas for light processing,” J. Opt. Soc. Am. B 30(9), 2365–2370 (2013). [CrossRef]
27. N. Yu, P. Genevet, M. A. Kats, F. Aieta, J.-P. Tetienne, F. Capasso, and Z. Gaburro, “Light propagation with phase discontinuities: generalized laws of reflection and refraction,” Science 334(6054), 333–337 (2011). [CrossRef]
28. G. Yuan, K. S. Rogers, E. T. Rogers, and N. I. Zheludev, “Far-field superoscillatory metamaterial superlens,” Phys. Rev. Appl. 11(6), 064016 (2019). [CrossRef]
29. J. W. Goodman, Introduction to Fourier optics (Roberts and Company Publishers, 2005).
30. A. Silva, F. Monticone, G. Castaldi, V. Galdi, A. Alù, and N. Engheta, “Performing mathematical operations with metamaterials,” Science 343(6167), 160–163 (2014). [CrossRef]
31. S. Chen, Z. Fan, H. Chang, and Z. Xu, “Nonaxial Strehl ratio of wavefront coding systems with a cubic phase mask,” Appl. Opt. 50(19), 3337–3345 (2011). [CrossRef]
32. W. T. Chen, A. Y. Zhu, J. Sisler, Z. Bharwani, and F. Capasso, “A broadband achromatic polarization-insensitive metalens consisting of anisotropic nanostructures,” Nat. Commun. 10(1), 355 (2019). [CrossRef]
33. D. Zhu, Z. Liu, L. Raju, A. S. Kim, and W. Cai, “Building Multifunctional Metasystems via Algorithmic Construction,” ACS Nano 15(2), 2318–2326 (2021). [CrossRef]
34. K. T. P. Lim, H. Liu, Y. Liu, and J. K. W. Yang, “Holographic colour prints for enhanced optical security by combined phase and amplitude control,” Nat. Commun. 10(1), 25 (2019). [CrossRef]
35. J. Chang, V. Sitzmann, X. Dun, W. Heidrich, and G. Wetzstein, “Hybrid optical-electronic convolutional neural networks with optimized diffractive optics for image classification,” Scientific Reports 8, 12324 (2018). [CrossRef]
36. G. Weitzstein, A. Ozcan, S. Gigan, S. Fan, D. Englund, M. Soljacic, C. Denz, D. A. B. Miller, and D. Psaltis, “Inference in artificial intelligence with deep optics and photonics,” Nature 588, 39 (2020). [CrossRef]
