Reprogrammable Graphene-based Metasurface Mirror with Adaptive Focal Point for THz Imaging

Recent emergence of metasurfaces has enabled the development of ultra-thin flat optical components through different wavefront shaping techniques at various wavelengths. However, due to the non-adaptive nature of conventional metasurfaces, the focal point of the resulting optics needs to be fixed at the design stage, thus severely limiting its reconfigurability and applicability. In this paper, we aim to overcome such constraint by presenting a flat reflective component that can be reprogrammed to focus terahertz waves at a desired point in the near-field region. To this end, we first propose a graphene-based unit cell with phase reconfigurability, and then employ the coding metasurface approach to draw the phase profile required to set the focus on the target point. Our results show that the proposed component can operate close to the diffraction limit with high focusing range and low focusing error. We also demonstrate that, through appropriate automation, the reprogrammability of the metamirror could be leveraged to develop compact terahertz scanning and imaging systems, as well as novel reconfigurable components for terahertz wireless communications.

Graphene shows great promise for the implementation of THz metasurfaces. This is due to its widely demonstrated plasmonic properties that provide subwavelength behavior required in the unit cell building block 25,26 . The reasonably low loss of graphene plasmons at THz frequency band has been the main motivation of several wave steering devices 14,27 . The plasmon tunability of graphene is achievable through chemical doping or electrostatic biasing. The tunablity of graphene opens the door to globally reconfigurable metasurfaces with fine-grained control of the amplitude-phase response 28 , or it can provide switching between diverse electromagnetic functionalities in particular setups 29 . Plasmon tunability can also be leveraged locally, in which case graphene becomes a natural candidate for the implementation of (re)programmable THz metasurfaces. In spite of their potential, graphene-based programmable designs have remained relatively unexplored, and current proposals are mainly limited to dynamic beam manipulation for diffusion 30 , encryption 31 , or vorticity control 17 .
In this work, a (re)programmable metamirror is proposed as shown in Figure 1 The device is conceived as a 2-bit coding metasurface that leverages the tunability of its graphene-based unit cells to control both the position and depth of focus. The shift of focus is achieved through a modification of the metasurface phase profile, that is achieved by simply changing the bias applied to individual unit cells. The generalized Snell's law of reflection allows to calculate the exact phase profile required for each target focal point, thus opening the door to fast precise focusing. If the change in the phase profile is automatized, it can be used to develop novel THz scanning and imaging devices. Figure 1 shows the proposed metamirror schematically. A coding metasurface composed of graphene-loaded unit cells is programmed through an FPGA to focus the reflected wave in an arbitrary position. Next, we show the details concerning the design of the unit cell, the coding metasurface, and an imaging device based on the systematic reprogramming of the coding metasurface.

Results
Unit cell. A reflective metasurface operating as a mirror calls for the use of unit cells with a high reflection amplitude and a wide phase range. Fabry-Pérot resonant structures composed of metallic patches over a square grounded substrate can provide such functionality by adjusting the size or position of the metallic patch to achieve the required phase 15 . This approach, however, offers no reconfigurability. In this work, instead, we achieve reconfigurability through changes in the electrostatic bias applied to uniformly sized graphene patches. This modifies the complex conductivity of the graphene sheet and, by extension, the phase of the reflected wave. Figure 2(a) shows the proposed unit cell. It maintains the Fabry-Pérot cavity structure with a graphene sheet over a grounded silicon substrate, yet adding a second graphene sheet over a high-density polyethylene (HDPE) layer. HDPE is chosen due to its particularly low loss behavior in THz band. The unit cells and the graphene patches are dimensioned to provide proper performance at the operation frequency, i.e. 2 THz in this work. The approximated circuit model shown in Figure 2(b) is developed in order to provide a theoretical validation of the unit cell behavior. As observed in Figure 2(c,d), the theoretical and simulation results show a close agreement. More details on the characteristics of the unit cell and its validation are provided in the Methods section.
Compared to the the single-layer graphene structures, our proposed unit cell based on dual-layer graphene enjoys more flexibility to achieve a wide range of phase response and a high reflection amplitude, simultaneously. The dual layer structure, modeled as two parallel RLC cells, provides an additional degree of freedom as each graphene patch can be biased independently. Figure 2(e,f) shows the reflection amplitude and phase response of www.nature.com/scientificreports www.nature.com/scientificreports/ the unit cell at 2 THz, respectively, as functions of chemical potential of the top and bottom graphene layers (μ c,1 and μ c,2 , respectively). A quick exploration reveals a complete 2π phase range as well as large design space areas with high reflectivity.
The coding set is built by picking points with as high amplitude as possible, and phases equal to multiples of π 2 /2 N bit where N bit is the number of bits. We will see that two bits are enough for our purpose. We choose μ c,1 = {0.6, 1.3, 0.1, 0.4} eV and μ c,2 = {0, 0.6, 0.1, 0.1} eV corresponding to the bit combinations B = {00, 01, 10, 11}, respectively. As observed in Figure 2(g,h), a phase shift of π/2 is maintained over the range of 1.9-2.1 THz with an amplitude close to 0.7. This high reflection amplitude for all states guarantees a high focusing efficiency for the designed metamirror.
Coding metasurface reflector in the THz range. Consider a reflective metasurface under illumination of an incident plane wave at elevation angle θ i and azimuth angle ϕ i . It is used to focus the energy in a given point defined by The coordinate origin is at the geometric center of the metasurface. The generalized Snell's law of reflection is used to formulate the phase profile over the metamirror. The required phase profile Φ(x′, y′) at an arbitrary point on the metamirror can be written as where λ 0 is the wavelength in free space, k 0 is the incident wavevector, and r′ = (x′, y′) is the vector position of an arbitrary point on the metasurface.
When the metamirror is illuminated by a normally incident wave ( The actual phase profile of the coding metasurface is obtained by spatially discretizing Φ(x′, y′), and rounding off the phase values to those offered by the closest coding states. Figure 3(a) shows a coding metasurface programmed to focus the beam at x = 0, y = 0, and z = 750 μm with four states. The size of the metasurface is 68 × 68 unit cells or 9λ × 9λ at 2 THz. As observed in Figure 3(b), the metasurface achieves effective focusing of the electrical field very close to the desired spot. The actual point of maximum field is (0, 0, 705), which implies an error of around 6%. By moving away from the focal point, the field spreads out as expected. As demonstrated in Figure 3(c-e), side lobes appear around the target focal spot due to the finite and discrete nature of the metasurface. Both the focal point accuracy and the side lobe amplitude www.nature.com/scientificreports www.nature.com/scientificreports/ can be improved with larger metasurfaces and higher resolution in spatial phase profile. Note that Figure 3(c-e) demonstrate close agreement between the numerical results and the theoretical approach, which is detailed in the Methods section. Figure 3(f-h) provide an assessment of the beam waist as a function of the focal length, and the aperture size of the metasurface D. To give an instance, the mentioned metasurface with focal length of 750 μm and aperture size of 1360 μm, the waist diameter of the half-maximum energy density is 2w = 104.9 μm, whereas the diffraction-limited full beam waist is given by . It is worth noting that the proposed metamirror operates close to the diffraction limit even for long focal length, which is generally difficult to achieve. Moreover, larger metasurfaces are apparently required to achieve a larger focal point range, because of the necessity to operate below the Fraunhofer distance Figure 4(a) shows how the 68 × 68 coding metasurface would be programmed to laterally shift the focus to x = 136 μm, y = 272 μm, z = 750 μm. Figure 4(b-d) confirm that energy is focused around the desired target with an excellent agreement with theory. The actual position of maximum field is less than 30 μm away from the target focal point. Maintaining the focal depth at z = 750 μm and with a lateral size of D = 1.36 mm, we perform a beam waist analysis analogous to that of Figure 3 but in the lateral dimensions. Charts (e-f) from Figure 4 show the beam waist in the X direction and its relative difference to the diffraction limit as a function of the positions of the focal point, respectively. It is observed that as the focal point is moved away from the metasurface center, especially in the Y direction, the beam waist in the X direction tends to increase since. A complementary behavior is observed for the beam waist in Y direction, as shown in charts (g-h) of Figure 4. The reason for such behavior is that the achieved asymmetric phase profile of the metasurface with specified dimensions results in an elliptical beam instead of a circular beam. As a matter of fact, a metasurface with finite dimensions is not capable to produce the necessary complete phase profile to achieve a circular beam pointed away from the metasurface center. As expected, one can operate close to the diffraction limit once being close to the metasurface center.
To note, the electric field around focal point decreases when the Z distance of focal point increase and also when the focal point position moves away from the Z axis which results in a decrease in the focusing efficiency. www.nature.com/scientificreports www.nature.com/scientificreports/ Programmable metamirror for 3D confocal terahertz imaging. The results obtained above demonstrate that the proposed metamirror can set the focus to any point in a distance from the metamirror. Due to the unequivocal relation between the phase profile and focal point given in Equation 2, it is relatively easy to develop algorithms that scan the focal point of the metamirror to perform three-dimensional confocal imaging. Figure 5 exemplifies the confocal imaging concept by the proposed metamirror. Let us consider three spherical objects, as point scatterers, placed within the area of the metasurface influence. Let us also assume that the diameter of the spheres D s is commensurate to the beam waist. By means of the principle of reversibility of the optical path, a strong reflected field should be observed once a plane wave is focused onto a sphere by the metasurface. On the contrary, weak reflected field should be observed when the sphere is absent on the focal point. Thus, one can make a 3D map (image) of the reflected field response when the focal point is scanned in space. Differences in the reflected field response due to the scan of focal point creates the confocal image.
Charts (a-l) of Figure 5 show the reflected field distributions at the far-field of the metasurface with assuming that sphere I is at the focal point (first column), sphere II is close to the focal point (second column), and sphere III is far from the focal point (third column). Each row assumes different locations for all the spheres. From different plots, it is clear that strong and weak reflected field distribution can be observed when the spheres are present and absent at the focal point, respectively. It confirms that the confocal imaging scheme works well in a variety of cases: (A) focus in the center of the metasurface, (B) focus not centered, (C) focus not centered and sphere II closer, and (D) focus not centered and sphere II in a different depth.

Discussion
In this work, we have presented for the first time the design of a fully reprogrammable metamirror that can be used in terahertz confocal imaging. The metamirror acts as a focusing reflectarray that leverages the tunability of graphene to achieve the required phase range and programmability at the unit cell level. At the metasurface level, we adopted a coding metasurface approach that allows to dynamically shift the focal point by adjusting the center and slope of a concentric phase gradient. This is in stark contrast to the design by Huang et al. 32 , wherein a single graphene layer is laid over the whole metasurface area and tuned using a single bias source. Such approach only allows to control the focal depth with limited range. Focusing with programmable coding metasurfaces has been only demonstrated in the GHz range with devices that use PIN diodes for reconfiguration. This limits the amount of achievable unit cell states, reducing the precision and focusing range.
We demonstrated analytically and numerically that the proposed metamirror can shift the focal point both laterally and in depth. With 2 bits per unit cell and assuming a total size of 9λ × 9λ, energy is focused around the desired focal point with less than 6% deviation in the evaluated points. Moreover, the metamirror can operate close to the diffraction limit as long as the focal point is located around the focal depth and centered with respect to the metasurface. We also showed that the programmability of the focal point can be applied to the development of compact devices for 3D terahertz confocal imaging.
As a final note, it is worth highlighting that this work brings the concept close to realization by making reasonable assumptions. First, we consider a graphene quality leading to a relaxation time of τ = 0.6 ps, which is well achievable with Chemical Vapor Deposition (CVD) techniques and encapsulated in hexagonal boron nitride www.nature.com/scientificreports www.nature.com/scientificreports/ (h-BN) 33,34 . Those techniques in fact produce significantly higher relaxation times, which would lead to lower loss and better overall performance. Moreover, considering the electrically thin boron nitride layer has a negligible impact in the simulation results. To mention, when the graphene has a lower relaxation time, the plasmonic losses within the structure increase. Therefore, it results in a device with a lower efficiency. Regarding the biasing scheme, different gate-controllable designs have been simulated 35 and experimentally validated 36 , which could be adapted to local addressing of unit cells. Before actually addressing fabrication, a challenge to be addressed in future work would consist on minimizing the interference imposed by the metallic biasing structure on the response of the metasurface. It is worthy to note that bilayer graphene, a material consisting of two layers of graphene, exists in the various forms such as Bernal-stacked form 37 . Depending on the relative position of layers, different quantum effects can happen in the spacing. But what we have here is different. In our proposed dual-layer structure, the insulating layer (here HDPE) between two graphene patches is sufficiently large for quantum effects to be negligible 38 . Therefore, crystallographic alignment is not required, and micrometric precision of current lithographic techniques is enough to achieve the results reported in the paper. Therefore, the process of fabrication of dual-layer metasurface structure is generally similar to one-layer structure; and the existing transfer printing methods could be employed 36,39 .

Methods
Unit cell modeling. The unit cell structure shown in in Figure 2(a) is implemented in CST Microwave Studio 40 . The length of the unit cell is d u = 20 μm and graphene patches are squared with side d g = 12 μm. The operation frequency is set to 2 THz. We consider bulk silicon as the substrate with relative permittivity ε r = 11.9, loss tangent tan(δ) = 2.5 × 10 −4 , and thickness h si = 10 μm. The extra layer is composed of a layer of HDPE with relative permittivity ε r = 2.37, loss tangent tan(δ) = 0.011, and thickness h HDPE = 4 μm.
The graphene layers are modeled as infinitesimally thin sheets with surface impedance σ ω = Z 1/ ( ), where σ(ω) is the frequency-dependent complex conductivity of graphene. The conductivity is given by where e,  and k B are constants corresponding to the charge of an electron, the reduced Planck constant and the Boltzmann constant, respectively 41 . Variables T, τ and μ c correspond to the temperature (T = 300 K in this paper), the relaxation time and the chemical potential of the graphene layer. Note that this expression neglects effects at the graphene edges and considers that the Drude-like intraband contribution dominates, which are www.nature.com/scientificreports www.nature.com/scientificreports/ experimentally validated assumptions at the sizes and frequencies considered in this work 42 . In all cases, the relaxation time of graphene τ = 0.6 ps, which is well achievable with current fabrication and encapsulation techniques 33 . The chemical potential of the top and bottom layers are initially left as parameters for exploration and then fixed to discrete values leading to the adequate amplitude and phase.
Simulation methods. We implement the unit cells and the complete metasurface in CST Microwave Studio 40 .
For the simulation of unit cells, a Floquet port is applied to excite an x-polarized normal plane wave onto the unit cell and the back-scattered wave is then measured. This way, we obtain the amplitude and phase of the reflection coefficient of the proposed unit cells while taking into consideration mutual coupling between adjacent unit cells. After obtaining the response of individual unit cells, the complete metamirror is modeled as an array of such unit cells with different phases. The metasurface is illuminated by a normal plane wave with x-polarization to obtain the scattered fields above the metasurface and, thus, the field distribution around the focal point. Finally, we simulate the THz imaging device by estimating the scattering response of the spheres in different positions. To this end, we again consider x-polarized normal illumination. To identify the presence of a sphere in the focal point of the metamirror, we compare the far field response of the metamirror with and without the spherical objects. Strong variations indicate the presence of a sphere in the focal point of the metamirror. To note, for the results shown in Figures 3 and 4, the reflected fields i.e. the difference between the total wave and incident wave are extracted to survey the achieved focal points while for the results shown in Figure 5, the difference between the reflected fields with and without the sphere are obtained to identify the presence of the spherical object in the focal point of the metamirror.
Theoretical validation. To gain physical insight into the scattering behavior of the proposed unit cell, an approximate circuit model has been developed. As sketched in Figure 2(b), the graphene patches are modeled as RLC series circuits (R 1 , L 1 , C 1 ; R 2 , L 2 , C 2 ) in parallel with equivalent transmission lines representing the HDPE and grounded silicon layers. The values of resistance, inductance, and capacitance of the graphene patches have been extracted from the reflection parameter at the boundary between air and unit cell, which only depend on the chemical potentials of the top and bottom graphene layers for a fixed patch size. For the results shown in Figure 2(c,d), (R 1 , L 1 , C 1 ; R 2 , L 2 , C 2 ) are (20.32 Ω, 0.013 pH, 0.107 fF; 35.97 Ω, 26.26 pH, 3.38 fF).
To validate the response of the metamirror around the focal point, we need to use analytical models compatible with near field conditions. Fresnel diffraction theory cannot be used because of the small focal length to lens diameter ratio of our device. Alternatively, we apply the Huygens' principle to model the currents excited at the metasurface as electrically small sources of radiation. Consider an incident wave polarized along the x-direction that excites an x-directed current density J x on all unit cells. Therefore, we can write the inhomogeneous potential wave equation as