Single-pixel terahertz imaging with enhanced edge detection using angular momentum of structured light Open Access

The fine structure embedded in propagating plane waves offers both theoretical scientific value and practical benefits across numerous photonics applications. Through precise manipulation of electromagnetic field properties—amplitude, wavelength, phase, and polarization—structured light can be generated, continuing to attract attention due to its versatility and the underlying aesthetic principles of wave physics. Within a wide range of notable applications, one can highlight single-pixel terahertz imaging, where engineering the phase and amplitude of illuminating light significantly enhances imaging capabilities.

The ability to generate and detect structured light—electromagnetic waves with a spatial inhomogeneity of amplitude, phase, and polarization—along with understanding its interactions with matter, continues to expand its practical application areas.^1,2 Structured light is successfully employed in photonic communications,^3,4 laser microfabrication,^5,6 imaging improvements,^7–9 tomography,¹⁰ and light sheet microscopy,¹¹ covering frequencies from visible to terahertz (THz).¹² When combined with the optimal placement of optical elements, structured light illumination can enhance imaging metrics, including resolution, sharpness, and contrast, beyond conventional expectations.^13–15 

One of the primary challenges in photonic applications is diffractive signal broadening. This limitation can be addressed through the intelligent structuring of light’s spatiotemporal properties. A good example can be a conical lens, known as an axicon, which enables the generation of the focal line, the so-called nondiffracting Bessel beam,^16–18 which can be further enhanced using a flat metasurface into an optical needle.^19,20 Nondiffracting elliptic Mathieu,²¹ parabolic Weber,²² and self-accelerating Airy beams^23,24 are created similarly using various complex transmission masks. Most intriguingly, the generalization of this concept leads to nondiffracting and nondispersive sub-cycle optical bullets.²⁵ 

The angular momentum of light represents another fascinating discovery linked to the azimuthal (ϕ) wavefront topologies represented by the functional dependence of exp(imϕ), where m labels the so-called topological charge. This wavefront topology creates a phase singularity, better known as an optical vortex,^26–28 due to the helical shape of the wavefront with m nested screw dislocations.^29,30 Frequently observed in electromagnetics, optical vortices carry orbital angular momentum (OAM).³¹ A single electromagnetic vortex is detected as a dark spot of light with a phase singularity at its core. The three-dimensional trajectory of a screw dislocation is a single entity because the time reversal changes the signs of topological charges of single vortices in the beam. For this reason, let us say that a negatively charged optical vortex can be interpreted as a positively charged vortex propagating backward. This has beautiful implications in diffraction, in which the creation and annihilation of positive and negative vortices happen in pairs and, therefore, can be perceived as a dark knot of light.^32–34 This interpretation of the complex evolution of vortical structures with various symmetries has led to the emergence of terms such as knots, braids, and bundles of dark light.^20,35

Terahertz radiation, which bridges the infrared and microwave regions of the electromagnetic spectrum, exhibits remarkable propagation characteristics with minimal losses through numerous dielectric materials and compounds. This property enables noninvasive imaging and nondestructive inspection capabilities in diverse fields, including materials science, biomedical examination, and security applications.^36–38 The raster scan approach³⁹ was first used in optoelectronic THz imaging, and from this implementation, the field of THz imaging techniques experienced a rapid evolution.^40,41 Various innovative techniques have emerged, including holographic imaging,^42–44 real-time imaging,^45–47 compressed detection,^48,49 and computational imaging.^50,51 Ghost imaging has recently garnered significant attention through single-pixel imaging architectures^52,53 that demonstrate promising results in the THz range.^54–56 

The practical implementation of compact THz imaging systems in real-world environments still faces significant technical challenges. Primary among these are the limited power output of THz emitters, the reliability of sensitive THz detectors, compact approaches, and effective designs of passive optical components. Furthermore, there is a lack of optimal solutions for integrating optical elements directly onto a chip, which would eliminate the need for precise optical alignment.⁴¹ Furthermore, image quality is particularly susceptible to sample positioning when using passive optical elements such as mirrors or lenses, especially given that most of the objects studied are relatively large. These positioning issues often result in degraded image quality, manifesting itself as poor contrast, reduced sharpness, or decreased resolution.

Compact photonic components^57–60 offer an efficient alternative to bulky conventional photonic systems. By leveraging diffraction phenomena, these components achieve reduced thickness and weight compared to traditional optical elements.^57,61–63 When sub-wavelength or near-wavelength thickness elements are arranged in specific spatial patterns, they generate phase shifts in transmitted radiation. Through interference, these phase shifts produce desired light patterns. At its fundamental level, this approach utilizes binary or multilevel diffractive elements that create a phase delay of (n − 1)t for the refractive index n and the local thickness t of the material. More sophisticated implementations employ metasurfaces, where local changes are made in the geometric phase of transmitted radiation by a structure composed of individual meta-atoms.^64,65

Fresnel lenses,^57,66 top-hat converters,⁶⁷ and axicons^9,17 are efficiently implemented with these ideas. Flat optical elements efficiently handle the generation of these nondiffracting beams; for instance, Airy beams are created under various conditions⁶⁸ within a wide range of wavelengths extending from optical^24,69 to the THz range.⁷⁰ 

Laser ablation technology can serve as an effective tool for fabricating diffractive zone silicon plates.⁷¹ Offering the advantage of maintaining the same technological conditions across the different fabrication runs, various optical elements can be produced to shape the THz radiation in different ways: by transforming Gaussian modes into Fibonacci structures,⁷² generating Bessel beams,⁹ and creating structured THz light in the form of Airy beams⁷³ in THz imaging experiments. It was shown that beam engineering—in particular, the structuring of THz illumination—has led to significant improvements in imaging capabilities, including enhanced resolution, increased contrast, and greater sensitivity when the material parameters of thin graphene samples have been measured.⁷³ 

Modern cameras and smartphones capture images in single-shot exposure using detector arrays and lenses with millions of pixels. The increase in pixel count beyond 20 megapixels may seem excessive and lead to data storage challenges. An alternative image retrieval method that can consider a single-pixel detector for image acquisition. This technique, which has been in use for over a century, remains particularly valuable in non-visible spectrum applications where detector arrays for specific wavelengths are either costly or unavailable.

Single-pixel imaging is an innovative method for image reconstruction that employs a single, compact, and highly sensitive single-pixel detector. This technique traditionally measures the overlap between the scene and various masks using a single-pixel detector,⁵² which then integrates these measurements with the mask information.⁵³ The effectiveness of this approach has been thoroughly validated through theoretical⁷⁴ and experimental studies.⁷⁵ 

In a classical sense, single-pixel imaging involves the integration of spatial light modulators (SLMs) on the focal plane of the camera lens, which facilitates the modulation of scene images using various masks before measuring light intensities with a single-pixel detector. In contrast, ghost imaging employs diverse structured light distributions, created by SLM, to illuminate the scene, subsequently capturing the reflected or transmitted light intensities.⁷⁶ Instead of capturing individual pixel samples from the observed scene, these methods evaluate the inner products between the scene and a collection of test functions. In particular, the use of random illumination or light collection functions plays a pivotal role in this approach, rendering each measurement a stochastic sum of pixel values that span the entire image.

These techniques are widely adopted in the single-pixel imaging community because they utilize SLMs to spatially modulate either the illumination or the light-collecting schemes. Image scanning techniques can be categorized into three main types: basis scanning, raster scanning, and compressive sampling.⁷⁴ In our research, we implement the raster scanning approach.

Direct implementation of THz imaging requires enhanced functionality, compact optical configurations, and improved user-friendliness.^41,77 An optimal compact solution integrates passive optical elements with detectors on a single chip, thus enabling increased sensitivity while eliminating optical alignment requirements.⁷⁸ However, this approach presents challenges in optical design. When THz illumination is incident on the object near the optical axis, the focus becomes extremely sharp within the substrate’s thickness range, resulting in a high aperture. As the beam waist becomes comparable to the THz radiation wavelength, conventional paraxial optics become ineffective. Consequently, an accurate description requires a nonparaxial illumination approach.

The use of OAM in single-shot microscopy⁷⁹ offers several advantages, including enhanced contrast, edge detection capabilities, shadowing effects, and Fourier filtering. These features typically involve Fourier filters or alternative approaches implemented in the light-collecting (i.e., imaging) portion of the optical system.^80,81

Special underlining within the aforementioned topics deserves edge detection, as it can provide important information for many subsequent visual tasks, such as object recognition, image retrieval, image segmentation, security, or medical imaging, while these visual tasks require extracting object boundaries or perceiving obvious edges.^82,83 Angular momentum and vortex masks routinely contribute in single-shot imaging to enhancements in edge detection of the sample.^84,85 However, the aforementioned methods are mostly pre-detection edge enhancement through physical wavefront engineering. Edge detection can also be performed using post-acquisition digital image processing methods (such as Sobel, Canny, or Laplacian operators in the image space).⁸⁶ 

Vortex phase masks introduce a helical phase structure, which performs a radial Hilbert transform on the light wavefront. This operation redistributes the intensity symmetrically, enhancing the edges in all directions.^87,88 The phase singularity at the center of the vortex mask creates a unique diffraction pattern. The phase changes rapidly at the vortex core, which translates to extremely steep phase gradients in the immediate vicinity of the singularity. These gradients are essentially the rate of change of phase across the beam’s cross section. The steep phase gradients around the singularity lead to localized areas of constructive and destructive interference. Because the phase is changing so rapidly around the singularity, very small changes in the object’s surface or refractive index cause large changes in the interference patterns. This results in high contrast at the edges of the object. Thus, a vortex pattern selectively enhances the high spatial frequencies that correspond to edges in the image. The topological charge of the vortex mask determines the order of spatial differentiation. Higher charges can enhance specific edge features or provide more pronounced edge detection. Depending on the topological charge (orbital angular momentum) of the vortex beam, it can have directional sensitivity, meaning that it can enhance edges in specific orientations.⁷⁹ In a Fourier optics setup, the vortex mask acts as a spatial filter in the frequency domain. Modulating the phase of the Fourier transformed image emphasizes edge information when the inverse Fourier transform is applied.^89,90

However, as far as our knowledge goes, edge detection in single-pixel THz imaging schemes, despite its promising potential in applications, has never been explored up until now.

In this work, the center of gravity is to reveal how the OAM can be viable in extending single-pixel THz imaging capabilities. We have incorporated vortex phases into single-pixel THz imaging through structured illumination using various flat photonic elements designed for THz radiation. Our investigation examines three configurations of the structured THz illumination that incorporate the orbital angular momentum of THz light: a paraxial Fresnel zone plate (ZP), a nonparaxial hyperbolic zone plate, and an axicon. This methodology contrasts with conventional approaches that utilize vortical phase masks in the collection pathway combined with Fourier lenses.⁷⁹ 

Our results demonstrate that structured vortical illumination significantly enhances edge detection in single-pixel raster scan THz imaging. We have verified that the structured vortical illumination detects the edges of the samples in a pronounced way for single-pixel raster scan THz imaging. We have systematically compared edge-detecting setups with their common counterparts and benchmarked them using quantitative metrics such as resolution and contrast. We have also observed significant edge enhancements for a group of near-wavelength line pairs, thus indicating the potential utility of this approach for the characterization of THz materials.

To validate our numerical predictions, we have manufactured numerous photonic elements using two competing technologies: first, a cost-effective 3D printing technology, and second, the high-precision laser ablation technique. Our comparative analysis revealed that both technologies demonstrate comparable performance in generating THz structured light illumination, with laser-ablated elements showing marginally superior results.

We introduce vortex versions of the phase masks described by Eqs. (2)–(4) by combining these elements with Eq. (5). The subsequent procedure is a substitution of the resulting phases into Eq. (1). The designs of the phase masks used in this study are given in Fig. 1.

The numerical implementation of the Rayleigh–Sommerfeld propagator was described by Voeltz.^94,95 In the numerics, we use a mesh grid containing M = 2048 pixels for the x and y coordinates.

Before we proceed to the experiments, we verify our expectations numerically. For experimental verification, we have selected a particular source with a frequency of 253 GHz, giving the wavelength λ = 1.186 mm. This choice was motivated by the 65-nm CMOS (Complementary Metal-Oxide-Semiconductor)-technology-based compact source available (details can be found in Sec. V dedicated to the experimental verification.)

We start our numerical validation with two-phase elements: the first is a standard zone plate [see Eq. (2)], and the other is the thin lens combined with a vortex phase [see Eq. (5)]. We illuminate the phase masks with a large Gaussian beam covering all apertures and present results in Figs. 3(a) and 3(b) for two cross sections—for the longitudinal (yz) and the transverse one (xy). The vortex mask results in an axial intensity drop—a void due to the presence of a phase singularity in the very center of the beam. In both cases, we detect the brightest spots in the longitudinal plane at a distance of z ≈ 23.5 mm, which is closer to the element than the envisioned focal point of f = 30 mm. To understand this discrepancy, we need to compare element sizes and distances with the wavelength. The actual location of the focal point is located $≈19.81λ$ from the phase mask. The numerical aperture of the system, NA = sin θ, can be estimated as tan θ = D/(2z_f), giving us θ = 0.824 rad, leading to NA = 0.734. Such a large value indicates that both phase masks operate within the nonparaxial regime, indicating, therefore, the choice of the Rayleigh–Sommerfeld integral for modeling.

Moreover, as analytical estimates suggest,^96,97 it is wise to ask ourselves how strongly the nonparaxiality of the focusing influences the vectorial field structure. In these simulations, we have assumed incident linear polarization. An appearance of longitudinal polarization is expected together with the cross-polarization, and we can estimate their influence as^98,99 orders 1/(kw) and 1/(kw)², giving ∼18% and 2% maximal contribution to the profile. Thus, in the experiment, we expect to see slightly larger beam profiles than in numerical simulations.

Due to the nonparaxiality of the focusing system, the standard zone plate design gives distorted beam profiles in the longitudinal plane since the design is based on a paraxial approximation.⁹¹ To correct this, we introduce a nonparaxial zone plate design based on Eq. (4) and verify its proper action numerically. The results obtained are shown in Figs. 3(c) and 3(d).

We see an immediate impact on the beam profiles in the longitudinal plane. First of all, the brightest spot is now located at the expected position z_f = f = 30 mm. Next, the beam profiles in the transverse plane are now symmetric for a nonparaxial zone plate (ZP) design and for the nonparaxial ZP with a vortex phase mask. The further location of the focal spot allows us to re-estimate its location as $≈25.45λ$ from the phase mask. The numerical aperture of the system NA can be estimated knowing that tan θ = 0.847, from where the highest angle of a ray that impinges on the focal point is θ = 0.703 rad, so the nonparaxial ZP is NA = 0.646. As the transverse sizes of both beams were approximately the same, we expect the experimentally measured beam sizes to be a bit larger than those obtained numerically.

Finally, we simulate an action of a nonparaxial axicon phase mask with and without combined vortex phase masks; see Figs. 1(e) and 1(f). The numerically simulated beam profiles for the transverse and longitudinal planes are shown in Figs. 3(e) and 3(f). An axicon is an element that does not create the brightest spot behind it like lenses but produces a bright focal line. The length of the focal line is known as a Bessel zone, and the intensity distribution is known as a nondiffracting Bessel beam. Knowing the element size and the axicon phase angle β = 0.4 rad gives us the length of the Bessel zone z_B ≈ 60 mm. The axial intensity distribution of a Bessel beam is approximately I ∝ z² exp(−z/z_B),¹⁰⁰ so it is not uniform and starts with a rising part until it reaches the maximal value at z ≈ 24 mm for both vortical and standard Bessel beams.

The performance of the optical system is often evaluated using a target with alternating light and dark bars.¹⁰¹ Resolution is determined by the finest discernible pattern, expressed as line pairs per millimeter. In single-shot imaging,⁷⁴ each line of the object is blurred, and the line spread function determines the extent of the blur. This blurring significantly affects the fine patterns.

Geometric optics, which treats light as rays, is a cornerstone of optical engineering. It accurately models light propagation in uniform media, bending at interfaces, and curved paths in varying refractive indices. However, it neglects wave phenomena such as diffraction and interference, which are valid when wavelengths are much smaller than interacting structures. Paraxiality, light propagation near the optical axis,¹⁵ is an essential assumption that leads to Gaussian optics.¹⁰¹ 

In THz frequencies, wavelengths become comparable to structure sizes, thus making geometrical optics less accurate. Therefore, high spatial frequencies and nonparaxial propagation become significant, necessitating more sophisticated models. Another important distinction of THz imaging is that conventional single-shot imaging is not available for these wavelengths, and single-pixel imaging with various strategies is used to obtain an image. A detailed review of the differences between single-pixel imaging and single-shot imaging can be found in Ref. 15.

In addition, we briefly recast some arguments and definitions from single-pixel imaging. Although we acknowledge the nonparaxial nature of the imaging setup, for the sake of simplicity, we neglect the influence of longitudinal and cross-polarized electric field components. To evaluate imaging properties, we conduct a numerical virtual experiment using a sample containing a series of three consecutive transparent and non-transparent stripes with progressively decreasing widths, as shown in Fig. 4(a). This sample allows for a direct comparison of object resolution with experimental accuracy, unlike the much finer sample used in our previous numerical studies.^15,92,102

It is important to note that the concept of contrast applies both to the object and to its image. We define image contrast as the percentage of object contrast reproduced in the image space, assuming uniform (plane wave) illumination. However, non-uniform illumination can introduce contrast variations, potentially degrading the image quality.

It is worth recalling that the perfect reproduction of object contrast and resolution is fundamentally limited by diffraction and other inherent physical constraints of optical systems. Even ideal, perfectly designed, and manufactured lens systems cannot completely overcome these limitations, as outlined in Refs. 101 and 103.

By plotting the transferred modulation (contrast) against the number of line pairs per millimeter in the image, we obtain a modulation transfer function (MTF), a common metric for characterizing the resolution capabilities of optical systems.

The limiting resolution of the system is determined by the intersection of the modulation function with a line representing the minimum detectable modulation level of the system sensor. For our analysis, we set this threshold at 20%, assuming it as a rational number for evaluation. This approach provides a quantitative description of the resolving power of the imaging setup.

It is valuable to remember that limiting resolution alone does not fully characterize a system’s imaging performance. Two systems with identical limiting resolutions can exhibit vastly different image qualities. A system with higher modulation at lower spatial frequencies generally produces sharper and more contrasty images.¹⁰⁴ However, when comparing systems directly, the optimal choice often involves a trade-off between high resolution and high contrast at lower frequencies. The ultimate decision depends on the specific imaging task and the relative importance of these factors.

All numerical processes (multiplying by phase functions, applying element phase transformation, and performing FFT) work on the full M × M grid. Compared to direct integration, in the FFT implementation, there are no nested loops or operations that scale with an additional factor of M² (or any higher power), leading to the complexity of $OM6$ because the transformation is not iterated over an additional dimension. Essentially, one operation per pixel (or a set of FFTs over the entire pixel grid) is performed, not a nested operation for each pixel across another independent dimension.

Now, in single-pixel imaging, the configuration is mostly the same except that the sample is being raster-scanned in both x and y directions. There are N_x positions in x and N_y positions in y. For each raster scan position, the field from the plane wave is incident on the sample at a new lateral offset (or equivalently, the sample is “shifted”). After propagation through the lens and onward, instead of reconstructing a full image, we record one intensity value at a single-pixel detector per scan position.

However, even if we are only taking one scalar measurement per scan, the propagation calculation is still performed over an entire M × M grid (because the FFT must be computed to properly account for diffraction, interference, and the phase induced by the elements). That means that the propagation computation [costing roughly $OM2⁡log⁡M$] must be rewritten for every scan position.

Similarly, a direct integration method would have a cost of $ONxNyM4$⁠, which is comparable to the single-shot imaging using direct integration. In the FFT implementation, even if the detector ultimately records only a single value, the propagation calculation still involves an FFT over the entire grid. The nature of the FFT is such that it processes all M² data points together. While it might initially seem that recording a single pixel would save computation, the numerical methods used (such as FFTs) do not compute only one element of the Fourier transform “for free.” The FFT algorithm is designed to efficiently compute the complete transform of the input array. There is not a straightforward way to “skip” computing the rest of the array without essentially performing the full transformation. For this reason, the computational complexity of the single-pixel raster scanning is N_xN_y times higher than that of a single-shot imaging.

With this in mind, we proceed to the evaluation of the performance of individual single-pixel imaging; see Fig. 4.

The sample, which we use for virtual experiments, contains seven groups of three alternating white and black bars of different widths, as presented in Fig. 4(a).

When imaged with a paraxial zone plate, we observe the resolved object; however, the image is blurred, increasing when the bars in the sample become smaller; see Fig. 4(b).

Now we introduce the OAM into a single-pixel imaging setup; see Fig. 4(c). As mentioned in the Introduction, vortical phases are one of the possible ways to detect edges in single-shot imaging. Traditionally, these phase masks are used on the light-collecting side.⁷⁹ However, here, in single-pixel raster scan imaging, we proceed in another way, implementing it for illumination of the sample (we will explain later why we deviate here from the common setup).

As expected, we observe an edge enhancement for the first three groups of pairs of lines with widths varying from 3.32, 2, and 1.416 mm. Note that only the outer edges are resolved for the last four groups with corresponding widths of 1.123, 0.928, 0.781, and 0.684 mm. Another question is the resolution of the vortical zone plate setup. The resolving power of this system depends on its ability to detect not only the space between objects but also the space between edges. Here, as the lines and edges approach each other, the blurs overlap, and the system’s resolving power is reduced and becomes below the acceptable level when compared to the case without the vortex phase. Although the contrast of the largest group for vortex illumination drops from 97% to 89%, the next three groups are resolved with similar contrasts of 82%, 76%, and 40%, respectively. The contrast of vortical illumination drops rapidly for the remaining groups.

To compare and investigate whether the resolution can be further enhanced, we have imaged the sample using a nonparaxial version of the zone plates; see Figs. 4(d) and 4(e). Visually, the nonparaxial zone plate performs better compared to Figs. 4(b)–4(d). The third group looks brighter, and the fourth group is resolved much better than the previous one.

For a better comparison and to quantify the contrast and resolution of the setup, we investigate the cross section of the sample in Figs. 4(h)–4(j). The contrast for the first group remains basically unchanged; however, the second group is resolved now with a contrast of 94%. The contrast of the third group also increased by 10%. The contrast of the fourth group is now at 60%. It is seen that only the last two groups of line pairs are practically resolved with the same power.

A nonparaxial version of the vortical zone plate also performs better in edge detection; see Fig. 4(e). Edges of the first three pairs of lines are resolved with practically the same contrast as without the vortical phase: about 91% for the first two groups, 1% smaller for the third group, and 5% smaller for the fourth and fifth. However, the last groups are not resolved, with only the outer edges of the group becoming highly pronounced.

The last example, which we consider here, is an axicon phase mask; see Figs. 4(f) and 4(g). From our previous study,¹⁵ we know that the axicon produces a modulation transfer function with relatively high contrast for higher spatial frequencies and lower contrast for larger details. Thus, we expect to resolve a larger number of line pairs in comparison with the two previous cases. In fact, visually, the sample looks clearer, and the resolution appears to be increased; see Fig. 4(l). The contrast of the first group is around 90%, and 74% and 66% for the next two groups. The contrast of the fourth group is the same as for the nonparaxial lens at 60%. The next two groups have a contrast of $≈15$% greater than that of the nonparaxial zone plate. Thus, six out of seven groups of line pairs are acceptably resolved, so it is meaningful to ask ourselves whether these behaviors will be translated into a Bessel vortex mask; the results are shown in Fig. 4(g).

They come as a surprise, as we see some “ghosting” in the edge detection and the contrast of even the third-line pair group drops. A possible explanation lies in the structure of a Bessel vortex; see Fig. 3(f). The central intensity ring is surrounded by two additional rings, which can be seen due to the lower intensity of the first ring compared to the Bessel beam; see Fig. 3(g). We plot a cross section of the MTF in Fig. 4(m) and summarize the numerical findings in Table I, demonstrating a worse performance of the axicon element.

The last important note here is regarding the asymmetry in spikes for the last group of stripes in Figs. 4(i), 4(k), and 4(m). This is caused by the numerical implementation of the fast Fourier transform in Ref. 95, which makes the vortex core slightly shifted by a percent of a pixel with respect to the optical axis. For this reason, the target reacts differently to positive and negative vortex phase masks and is also distance z₂ dependent. This effect is especially pronounced for small targets and complex illumination structures, such as a Bessel beam. This effect reduces when the numerical mesh grid has more pixels, M, resulting in an increase in computational time as M². We do mention it here, as this is an indication of what we have to expect in the experiment, as the small structure is much more sensitive to the exact position of a vortex core.

An insightful reader might now ask the next possible good question: What is the design and assembly of the imaging setup? One can naively ask how we can be sure that the positions we choose for the sample and the detector are correct and optimal. The answer to this question is that we have chosen the distance from the illuminating element z₁ to be the place on the z-axis where the illumination is the brightest. For the selection of the position z₂ of the light-collecting element, several simulations were performed with a finer sample, which we have introduced in Ref. 15. These results in the form of three MTF maps are given for a particular selection of elements and frequencies used in the experiment and are presented in Fig. 5. As one can notice, the results obtained are very different, indicating that the best resolution does not coincide with the highest intensity. This conclusion was already presented for smaller numerical apertures and higher THz frequencies.¹⁵ A pair of nonparaxial zone plates provides the largest degree of freedom in the assembly of the optical setup. We selected the distances z₂ as the best common resolution points. These points were used to obtain Fig. 4.

To this end, we plot a cross section of the MTF function for all cases we consider here (Fig. 6). We see that in general, the contrast of resolved line pairs is slightly lower for the vortical illumination than that of nonvortical elements if we consider the internal structure of a group of line pairs; see the 20% level in the graph.

It allows us to get the dependence depicted in Fig. 6(b). For this purpose, we have numerically reevaluated all situations showcased in Fig. 4 using Eq. (10). In doing so, we restricted ourselves to spatial frequencies ranging from 0.45 to 0.75 lp/mm, as they represent the last four groups in the sample in Fig. 4(a). As we read from the graph, the contrast of the whole group of slits is around 50% for the lower spatial frequency representing the fourth group, around 70% for the next spatial frequency in the sample, with a further increase to around 79% in the next one, and reaching the value of $≈92$% for the highest spatial frequencies.

As a small group of alternating transmissive and light-blocking bars can be perceived as some effective medium,^105,106 this consideration makes a strong case for the introduction of such a parameter. Although scalar optics simplifies the interaction of light with gratings by ignoring polarization, the scalar effective medium can be introduced and explained using the Huygens–Fresnel principle. This sensitivity to the outer edges of a structure with a sub-wavelength inner structure might suggest the sensitivity of the current setup to flat boundaries between two materials, which should be investigated in further studies.

At this point, it is rational to recall and better explain the importance of the research by recreating a common single-shot edge detection vortical scenario,^79,81 but with single-pixel raster scan-based THz imaging; see Fig. 7. Here, we present one particular example of many attempts to use the Fourier-based approach with a vortex phase mask at the light-collecting end. Although we might notice some slight horns in some line pair profiles, they are minuscule and disappear when finer features are observed in a sample. An adjustment of the position of the sample, together with an adjustment of the collecting optics, slightly enhanced edge detection for some particular spatial frequencies but, in general, did not resolve edges compared to the structured vortical illumination approach summarized in Fig. 4. Thus, in the single-pixel nonparaxial THz imaging, OAM of light in the light-collecting part does not produce the same effects as observed in single-shot paraxial imaging.

An explanation for this finding comes from the understanding that in single-shot imaging, the edge-enhanced image is obtained as the convolution of the electric field of the image E with a vortical kernel ∝exp iφ.⁷⁹ In single-pixel imaging, lenses are strongly nonparaxial, with NA reaching 0.7. Under these conditions, lenses do not perform Fourier transforms of all objects or image planes.⁹² The second obstacle is the fact that the single pixel is static and located on the axis, so a vortical phase mask on the light-collecting side will result in a lower intensity detected by a detector.

For experimental verification, we were looking at two candidates: the first is relatively inexpensive and fast but possibly not accurate enough, while the other one is more precise and time-consuming in fabrication; however, it is not in the range of low price.

The first technology used to manufacture the vortex beam-shaping elements is the 3D extrusion printing method. This approach is simple and cost-effective, enabling the conversion of computer-generated models into tangible objects. In particular, 3D printing is particularly attractive for THz applications because filament polymers exhibit properties suitable for this range.¹⁰⁷ In this research, High Impact Polystyrene (HIPS) was chosen because of its high transmission of 0.82 and a refractive index of 1.53, comparable to glass commonly used in visible-range applications.

A Creality K1 extrusion printer was used for the printing process, equipped with a 400 μm nozzle and a vertical step of 100 μm—the smallest step available in the manufacturer’s software. The nozzle operated at a temperature of 255 °C, while the construction plate maintained a constant temperature of 100 °C. The printing used a “lines” infill pattern, alternating the infill lines horizontally and vertically in each layer, with a 100% infill density to ensure maximum homogeneity of the components. A single 2-in. optical component required up to 1.5 h to print and consumed ∼20 g of filament; see Fig. 8. Their profiles were measured using an optical profilometer (Sensofar S Neox). Figures 9 and 10 present a selection of diffractive THz elements made by this method.

The profiles we observe are in relatively good shape, although we observe some irregularities in the cross sections—missing steps, wrong depths, and improper widths of the same phase level. Two additional 3D printed elements are nonparaxial zone plates with alternatively embedded vortex phases [see Fig. 10(a)] and a single vortex phase with phase levels N = 8 [see Fig. 10(b)]. We produced these two elements as we tried a two-element setup when we physically attached a nonparaxial zone plate, Fig. 9(a), and a vortex phase, see Fig. 10(a), and when we added two unwrapped phases together; see Fig. 10(b). This was done to select the most experimentally viable option.

The second technology we use to produce THz optical components is the laser ablation of silicon.⁷¹ Due to a high index of refraction (n = 3.46), this material allows us to achieve a high phase modulation in significantly thinner optical components compared to 3D printed polymers. For the production of optical elements, we use a silicon wafer of high resistivity, 525 μm thick. We perform laser ablation using an optical setup employing a 6 W, 220 fs, 1030 nm Carbide laser by Light Conversion Ltd. and an excelliSCAN 14 scanner together with the f–θ lens (f = 100 mm). The beam diameter at the focus using this setup is around 24 μm. We chose the shortest available pulse duration, 1.5 W of power at a frequency of 50 kHz, and set the scanning speed at 200 mm/s. Using these parameters, a single scan of the pattern removes around 1.9 μm of material. To achieve phase modulation of 2π, we need to ablate λ/(n − 1) = 1.2/(3.46 − 1) = 488 μm of the material. For this particular wavelength, laser ablation was a rather slow process—the production of some elements took up to a few days.

For the sake of brevity, only two profiles of elements manufactured by laser ablation are presented; see Fig. 11. As can be seen, for the vortex, the first element is a quantized phase mask with N = 8. Although the heights are smaller now, we observe relatively flat steps of the element. Next, we have a vortical axicon element; see Fig. 11(b). Once again, the precision provided by the laser ablation is good compared to that of the 3D-printed element. This conclusion was reached by importing the profiles obtained using the profilometer and comparing them in a numerical experiment.

Several experiments are performed involving fabricated elements. The THz light intensities produced by the two competing technologies were found to be very similar. In what follows, we briefly recap these experimental findings.

Manufactured elements were characterized with the aim of evaluating their focusing capabilities and comparing the obtained spatial distributions with the numerical expectations of Sec. III. For this purpose, a special THz beam characterization and shaping setup was assembled [Fig. 12(a)]. As a compact THz radiation source, the 65-nm CMOS technology-based differential Colpitts oscillator (Terahercų technologijos MB, Vilnius) was designed for the frequency of 253 GHz. This oscillator employs capacitive and inductive feedback loops to maintain stable fundamental oscillations near 84 GHz, enhancing third harmonic components to deliver an output power of 78 μW at 253 GHz.¹⁰⁸ 

THz radiation was detected using a microbolometer mounted on 3D motorized linear stages. This detector integrates an air-bridged Ti microbolometer coupled with a THz resonant dipole antenna (noise equivalent power in the range of $20−60pW/Hz$⁠, where the geometric parameters of the antenna determine the resonant frequency.^109,110 In addition, in a single-pixel THz imaging setup [Fig. 12(b)], the detector remains fixed while the sample is mounted on motorized stages. The system consists of three lenses; the first lens (L₁) is the photonic component currently under investigation, and the second lens (L₂) is used for the collection of focused radiation. L₂ was manufactured using the HIPS filament extrusion printing method and the diameter of the component is 50.8 mm with a focal distance of 3 cm and is quantized in eight equivalent phase increments. The final lens (L₃) is a plano–convex component with a curvature radius (R = 27 mm), has a focal distance of ∼3 cm, and was also manufactured using 3D printing with HIPS extrusion. The purpose of this lens is to ensure that the detector collects all the radiation. In this configuration, another detector was used, the broadband log-spiral antenna coupled to a nanometric field effect transistor (FET) displaying equivalent noise power reaching $<100pW/Hz$⁠; fabricated by 90-nm CMOS technology.¹¹¹ However, it was not used for beam-shaping characterization due to its top-mounted lens, which enlarges the signal collection area and complicates the scanning. This limitation does not affect single-pixel imaging, as the collecting lens [denoted as L₂ in Fig. 12(b)] mitigates the impact of the enlarged collection area.

The role of a bandpass filter manufactured from steel foil 30 μm with rectangular holes of optimized dimensions was to reduce the possible influence of other harmonics. The filter inserted in front of the detector transmits a narrow frequency band in a specific polarization.¹¹² The filter cells are spaced 780 μm apart horizontally and vertically, with a rectangular cutout measuring 624 μm in length and 92 μm in width. The filter’s central wavelength is 244 GHz with a bandwidth of 36 GHz. The beam profile raster scan and imaging were performed with a step size of 0.1 mm.

For imaging quality evaluation, a modified USAF 1951 resolution target containing both vertical and horizontal slit openings was used. The width of these slits varies from the largest to the smallest as 3, 2, 1.5, 1, and 0.75 mm. The smaller slits are below the THz radiation wavelength (λ = 1.19 mm); therefore, they are not suitable for the spatial resolution estimate.

Exemplified cases are presented in Figs. 12(c)–12(e). First, we see an action of a nonparaxial zone plate given in panel (c). The formation of the focal spot can be seen at a slightly closer distance of z = 29 mm, which might be caused by a small wavefront curvature of the THz Gaussian beam incident on the photonic element. Regardless of the technology used for the optical element manufacturing, the behavior of the spatial profile detected experimentally is very similar. We note an appearance of small interference fringes due to standing waves in the setup.

Next, we tried a variety (recall that we had three combinations of vortex masks and two competing manufacturing technologies) of approaches to generate a focused nonparaxial THz vortex beam, and the best result from a quality point of view is shown in Fig. 12(d). We expected the focal point to form at z = 30 mm, but due to the implementation of the elements, it appears closer to the element in the experimental coordinate frame. The interference fringes due to the parasitic nature are more pronounced as the intensity in the ring of a vortex is less than that in the center of a Gaussian beam.

Moreover, we show the intensity distribution of a vortex Bessel beam; see Fig. 12(e). We have also tried five implementations of the vortex Bessel beam using two technologies and three designs described above. The results are similar, as we observe two interference mechanisms: the first one is caused by parasitic reflections within the element, and the second one is from parasitic interference with adjacent elements in the scheme.

Finally, we conducted imaging experiments for all the cases discussed here to verify numerical expectations and confirm successful edge detection in a sample while maintaining similar contrast and resolution. The best experimental results are presented in Figs. 12(g)–12(j) panels. We chose a nonparaxial zone plate [Fig. 12(g)], as a resolution enhancement over the paraxial zone plate was expected. When the vortex zone plate is placed in the illumination part, it results in a more pronounced edge discernibility in the sample image, as shown in Fig. 12(h). Although the experiment is performed at rather low intensities in a vortex ring, and therefore, the imaging of the sample was a rather challenging task due to the fine-tuning of the positions of the sample, detector, and optical elements, it experimentally confirms the predicted conceptual scheme.

In stark contrast to the virtual numerical experiment, the experimental implementation had a polarization oriented horizontally, so the THz radiation reacted differently to vertical and horizontal slits. To estimate it quantitatively, we have plotted a cross section of the images at one particular position for both vertical and horizontal slits. The results are presented in Figs. 12(i) and 12(j). As expected, the polarization had an influence on the interaction of THz radiation with a sample, with vortex illumination being the most sensitive.

In general, the trend in contrast and resolution is very similar to that predicted numerically—the vortex zone plate allows us to resolve better edges of the studied samples. This effect is more pronounced for horizontal stripes, as the currents induced in the stripes coincide with the orientation of the incoming THz light. A background slope in the intensity distribution that induces asymmetry in edge enhancement is probably caused by the nonuniform spatial profile of the structured light illuminating the sample under test; see Fig. 12(d). As we mentioned in Sec. III, an asymmetry between edges is caused by the very strong sensitivity of a single-pixel detector to the alignment of an optical axis.

Detailed estimates for both measured geometries are given in Table II. The contrast values for vortical cases are lower than expected for groups of 3 and 2 mm, but for a group of slit size 1.5 mm, it is close to the expectation; compare with Table I. We note that for a horizontal group, the contrast seems to be higher for the last two groups than expected when slits are oriented horizontally but much lower when they are oriented vertically, hinting at the fact that the scalar theory does not fully cover the underlying physics here.

As an alternative method for enhanced edge detection in THz imaging, the use of spatial filtering, such as dark-field and phase-contrast methods, can indicate the use of spatial filtering, which selectively modifies frequency components in the Fourier domain.¹¹³ Dark-field filtering removes low-frequency components to enhance scattered signals, while phase-contrast imaging introduces controlled phase shifts to reveal weakly absorbing structures. These methods can improve the contrast by up to 30 dB.¹¹³ While a direct comparison is hardly possible due to differences in experimental targets, the incorporation of vortex phases into single-pixel THz imaging through structured illumination enhances edge contrast without additional filtering directly at the imaging stage. It allows for a reduction of the complexity of the system, making the implementation route easier and simpler.

Conventional single-pixel imaging often relies on basis scanning or compressive sensing techniques,⁷⁴ employing structured illumination or detection patterns derived from Hadamard matrices, Fourier basis, or random patterns.¹¹⁴ These methods typically involve a stationary sample, while the illumination or detection patterns are modulated. An alternative approach involves raster scanning the sample itself across the field of view. Although conceptually it is simpler to implement and requires less computational overhead for image reconstruction, raster scanning can exhibit certain limitations related to the need for a large number of measurements and uncompressed dark noise compared to more computationally intensive methods.⁷⁴ 

The theoretical framework for raster scanning often relies on approximation, where test functions are simplified to delta functions. Although this simplification provides a basic understanding, we have previously revealed that nonparaxial effects significantly alter the imaging characteristics when employing structured illumination and light collection.¹⁵ In a nonintuitive way, we observed notable changes in contrast, resolution, and focal behavior. Moreover, the optimal configuration of diffractive optical elements for illumination and detection is strongly influenced by their specific properties.

In this research, we have focused on the question of whether OAM is a viable option for single-pixel THz imaging. This is not a trivial question, as we found that for single-pixel raster scan imaging, the positioning of the sample with respect to the illumination and detector does not obey laws from single-shot imaging. Although various attack angles on edge detection in single-shot imaging are known, most of them involving vortex phases are performed using several Fourier transforms for the light-collecting part of the setup.

Edge detection for paraxial single-pixel imaging can be assumed to be performed similarly to single-shot imaging using arrays of pixel detectors. However, in the given study, we deal with NA up to 0.7, raising the question of whether this technique can even be transferred to single-pixel imaging. This research revealed that edge detection is not possible in single-pixel raster scan THz imaging when vortical phase masks are placed in the same way as in single-shot paraxial imaging.⁷⁹ It has been shown that the incorporation of vortex phases into single-pixel THz imaging using structured illumination with a few selected examples of flat photonic elements allows one to successfully circumvent limitations created by the strong nonparaxiality of the optical system and detect the edges of samples. It was reached via a comprehensive analysis of three exemplar cases of structured THz illumination with the OAM of light: the paraxial zone plate, the nonparaxial zone plate, and the axicon. This approach is in sharp contrast to common schemes when a vortical phase mask is used in the light collecting part of the optical system. Edge-detecting THz imaging setups and their common counterparts are compared and benchmarked using metrics such as resolution and contrast.

Moreover, as the sizes of fine structures comprising the sample used for imaging reach near wavelength sizes, it is observed that structured vortical illumination increases the intensity ratio between the edge of the sub-wavelength line pair group and the center of the group, where intensity reaches zero as we further decrease the sizes of the line pairs.

This can be perceived as an interaction of structured illumination with vortical phases with the group of line pairs as a whole entity, i.e., as with an effective medium. This finding encourages further investigation of whether this approach is viable for THz material inspection and raises the question of how sensitive it could be to variations in material properties in a sample. Furthermore, exploring the applicability of this method to layered and composite materials could provide more insight into its potential to identify subtle structural variations and detect hidden defects or imperfections.

Numerical predictions were illustrated and supported via various photonic elements manufactured using two competing technologies: low-cost 3D printing technology and high-precision laser ablation technology. We did find that both approaches allow the fabrication of flat optics elements performing similarly when generating THz structured illumination, with laser-ablated elements having a slight edge. It opens a perspective of technological routes and possible trade-offs in further evolution in structured-light-based THz imaging systems. It can also be attributed to the developments in real-time imaging, as structured THz light illumination with OAM can partially increase the recording speed by enhancing edge contrast at the acquisition stage, hence reducing the need for post-processing and extensive image reconstruction.

This research has received the funding from the Research Council of Lithuania (LMTLT): The experimental part of the research was supported via Agreement No. S-MIP-22-76, and the theoretical investigation via Agreement No. S-MIP-23-71.

The authors have no conflicts to disclose.

Sergej Orlov: Conceptualization (lead); Formal analysis (lead); Funding acquisition (equal); Investigation (equal); Methodology (equal); Project administration (equal); Resources (equal); Supervision (equal); Validation (equal); Writing – original draft (lead); Writing – review & editing (equal). Kasparas Stanaitis: Data curation (equal); Investigation (equal); Software (equal); Validation (equal); Visualization (equal); Writing – review & editing (equal). Paulius Kizevičius: Data curation (equal); Formal analysis (equal); Investigation (equal); Methodology (equal); Software (equal); Validation (equal); Visualization (equal); Writing – review & editing (equal). Paulius Šlevas: Data curation (equal); Formal analysis (equal); Investigation (equal); Methodology (equal); Visualization (equal). Ernestas Nacius: Investigation (equal); Methodology (equal); Software (equal). Linas Minkevičius: Data curation (equal); Formal analysis (equal); Funding acquisition (equal); Investigation (equal); Methodology (equal); Project administration (equal); Resources (equal); Supervision (equal); Validation (equal); Writing – review & editing (equal). Gintaras Valušis: Conceptualization (equal); Formal analysis (equal); Funding acquisition (equal); Project administration (equal); Resources (equal); Supervision (equal); Validation (equal); Writing – original draft (equal); Writing – review & editing (equal).

The data that support the findings of this study are available from the corresponding author upon reasonable request.