Imaging using terahertz (THz) radiation (0.1–2 THz) is appealing as many visibly opaque materials are transparent in this frequency range [1]. THz imaging enables, for example, the detection of faults in space shuttle panels [2] and investigations into the substructure of paintings [3]. The THz photon energies are non-ionizing, unlike x rays, hence inspections will not damage sensitive electronics [4]. However, THz imaging has two main outstanding challenges. First, THz detector arrays are usually expensive and very narrowband [5,6]. Second, the relatively long wavelengths (0.15–3 mm) impose that diffraction-limited imaging will fail to see micrometer-scale details.

To date, most subwavelength THz imaging techniques use near-field raster scanning techniques, where a sub-diffraction-limited probe tip is raster scanned over a surface, sampling the evanescent field at each location [7,8]. These methods have achieved 10 nm resolution [8]. While impressive, the weak signals emanating from the probe are sensitive to detector noise, necessitating long measurement times, and complex detection schemes are required to improve acquisition rate [9].

Recently, alternative imaging techniques using single-pixel detectors have emerged [10,11]. Single-element detectors are generally cheaper and more robust than detector arrays. These methods spatially pattern the incident beam of radiation (described later). Orthogonal patterns minimize the mean square error [12]. However, unambiguous reconstruction of an $N$-pixel image normally requires $N$ measurements. To circumvent the trade-off between imaging time and resolution, compressive sampling techniques have been developed that make use of assumptions about the nature of the object to reconstruct an image from an undersampled set of measurements [13–16].

Here, we demonstrate near-field THz imaging using an ultra-thin, silicon photomodulator. The ultra-thin wafer allows us to access the THz evanescent fields to achieve a spatial resolution of 9 μm ($\lambda /45$ at 0.75 THz), demonstrated experimentally by imaging a resolution target. We conclude by investigating two different methods that improve the acquisition time by reconstructing images from undersampled sets of measurements: adaptive sampling and compressed sensing.

Single-pixel imaging schemes have been implemented in the THz regime [17–19]. An approach to this employs dynamic spatial patterning of a THz beam using a photo-conductive modulator [20]. This technique exploits the fact that when an optical pump beam is incident onto a photomodulator, such as a silicon wafer, electron–hole pair photo-excitation increases its THz conductivity [21] and reduces THz transmission. By shaping the optical pump beam into a binary intensity pattern, we are able to spatially control which areas of the modulator transmit THz radiation. In this way, the pattern in the optical pump beam is imprinted onto the incident THz beam. Previous studies have shown that the photomodulator thickness limits the achievable resolution, as this determines the amplitude of evanescent fields interacting with the object [20].

Figure 1 illustrates our imaging setup (a detailed schematic is shown in the supplementary material of Ref. [20]). For generation and detection of THz radiation we use a pair of ZnTe crystals in a standard THz time domain spectrometer [22]. The system is pumped by an 800 nm (90 fs) amplified Ti:sapphire laser running at a repetition rate of 1050 Hz. This laser system also provides a photo-excitation beam that is spatially structured by a digital micromirror device. This patterned beam (fluence of $\sim 100{\mathrm{\mu J}/\mathrm{cm}}^2$ per pulse) is projected onto a silicon wafer (6-μm-thick, $8000\mathrm{\Omega }·\mathrm{cm}$ resist) using a $+5\mathrm{cm}$ lens. Note that a 6-μm-thick modulator is the optimal thickness, determined by the interplay between absorption and diffraction—see Supplement 1, Section 1.A.

 

Fig. 1. Illustration of imaging setup: using a digital micromirror device and a lens, a pump pulse is spatially structured and projected onto a silicon wafer. This spatially modulates a coincident THz pulse. This THz pulse then passes through an object and is measured on a single-element THz detector. Inset is an optical image of a resolution test target (cartwheel) manufactured from gold on a 6-μm-thick silicon wafer.

Download Full Size | PDF

In the regions that are illuminated by the optical photo-excitation beam, the silicon conductivity increases, rendering the material response at these locations opaque to THz radiation [21]. Then a single-cycle THz pulse (see Supplement 1, Fig. S1) arrives $\sim 5\mathrm{ps}$ after photo-excitation. Since carrier diffusion can be neglected on these fast timescales, the pattern from the 800 nm beam is transferred to a THz pulse without spatial blurring. The patterned THz field then interacts with an object placed on the exit interface of the Si wafer before being collected onto a single-element detector. We record the peak amplitude of our THz pulse for our signal (see Supplement 1, Fig. S1).

Dynamic spatial encoding of the THz beam enables an image to be reconstructed from a series of measurements with a single-element detector. We measure the total THz field transmitted through an object as it is illuminated with a sequence of spatially patterned beams. Measurement $y_i$ represents the transmission corresponding to the $i$th pattern. To reconstruct the image of the object from the series of measurements, the spatial structure of the set of illuminating patterns can be stored in a sampling matrix $\mathbf{A}$. Here the $i$th row of $\mathbf{A}$ is an $N$-element vector that is a 1D representation of the $i$th projected pattern. Therefore, for an $N$-pixel image we have measurement $y_i=\sum _{j=1}^N{A}_{ij}{x}_j$ where $x_j$ is the $j$th image pixel. This a matrix representation of

As described above, $\mathbf{A}$ represents the basis in which the image is expanded. When the image is fully sampled with $N$ orthogonal measurements, $\mathbf{A}$ satisfies ${\mathbf{AA}}^T=\mathbf{I}$, where $\mathbf{I}$ is the identity matrix. In this case the solution is given by $\mathbf{x}={\mathbf{A}}^T\mathbf{y}/N^2$. Further, it is known that an orthogonal basis minimizes the detector noise [12]. For this reason, when the number of measurements equals the number of image pixels, we use a Hadamard matrix ($\mathbf{A}=\mathbf{H}$) formed from the Paley type-I construction (see Supplement 1, Section 5). $\mathbf{H}$ consists of elements taking the value of 1 or $-1$ in equal number. As our intensity masks cannot represent negative numbers, we use a lock-in amplifier to record the difference in transmission between a mask consisting of 1s (i.e., transmitting radiation) and 0s (i.e., blocking radiation) and its inverse. This measurement also minimizes low-frequency intensity fluctuations of our THz source (see supplementary material of Ref. [20]).

If an object is undersampled, $\mathbf{A}$ has fewer rows than columns, and the problem is underconstrained, with infinitely many solutions that satisfy Eq. (1). Adaptive sampling and compressive sensing are both strategies that make use of assumptions about the nature of the object (such as object sparsity when represented in a particular basis) to choose one solution that most likely represents the object. Depending upon the strength of these assumptions, and the level of noise in the measurements, these approaches have been proven to enable functioning image reconstruction from highly undersampled measurement sets [13–16]. More details of the compressive approaches we use in this work are given below.

In near-field imaging, subwavelength resolution can be achieved due to the interaction of near fields with the object. However, near fields decay exponentially with distance. In our approach, the thickness of the modulator can therefore be expected to play an important role in determining the ultimate resolution of our THz images. To investigate this, we image a subwavelength-sized metallic resolution target (cartwheel) through a silicon photomodulator of varying thickness $h$, as shown in Fig. 2.

 

Fig. 2. (a)–(c) THz images of the cartwheel shown in the inset of Fig. 1 taken through 400-, 110-, and 6-μm-thick silicon wafers, respectively. See Supplement 1, Fig. S4 for a close-up annotated version of part (c). Note, the cartwheels in (a) and (b) have diameters larger than the field of view. The origin of the vertical lines in (a) is discussed in Supplement 1, Section 4. (d)–(f) Calculated diffracted fields of a cartwheel as propagated through 400-, 110-, and 6-μm-thick silicon, respectively. THz polarization is horizontal in the experiment.

Download Full Size | PDF

For a relatively thick modulator ($h=400\mathrm{\mu m}$, of the order of the THz wavelength), we see that the subwavelength features of the cartwheel are not evident in the image [see Fig. 2(a)]. One can understand the resulting image by considering the diffracted field expected for the object when propagated through a thickness $h$ of silicon (refractive $\text{index}=3.44$), plotted in Fig. 2(d) using scalar diffraction theory [23] (see Supplement 1, Section 3 for details). While agreement is imperfect (see below), we see similar blurring to the cartwheel edges, particularly for the high spatial frequency components toward the center of the cartwheel, completely distorting the final image. By reducing the photomodulator thickness to 110 μm, we begin to recover an image resembling a cartwheel [see Figs. 2(b) and 2(e)], with only the center of the image distorted. Only when we reduce the modulator thickness to 6 μm do we recover a reasonably complete image of the cartwheel [see Figs. 2(c) and 2(d)]. The overall trend here is clear: as the thickness of the modulator is reduced, the images sharpen. Hence, due to the increasing spatial frequency content of the cartwheel toward its center, we can estimate our obtained resolution by evaluating the minimal distance for which the cartwheel arms are distinguishable. This leads to resolution estimates of 154, 100, and 9 $(\pm 4)\mathrm{\mu m}$ for $h=400$, 110, and 6 μm, respectively.

We note that there are two main reasons for the discrepancies between the left- and right-hand panels of Fig. 2. First, the polarization of the THz field is important, while our scalar diffraction calculations neglect this. This leads to the breaking of rotational symmetry in the experimental images. Indeed, the effect of polarization can be observed explicitly when we vary the orientation of certain objects (see Supplement 1, Section 2). Second, we must note that the optical pump light has a finite penetration depth in silicon (11 μm for 800 nm [24]), which will influence the diffracted field. This has two effects: first, the modulation efficiency will decrease due to transmitted pump intensity being wasted. Conversely, this finite penetration depth actually decreases the effective thickness of the modulator, increasing resolution. We discuss these effects and the optimal modulator thickness in the Supplement 1, Section 1.A. Further investigation with direct gap modulators is therefore required in order to push the THz image resolution lower than 9 μm using this method.

As with all near-field imaging techniques, small signals go hand in hand with long measurement times. However, as our imaging approach does not rely on raster scanning, we can reduce acquisition time by reducing the total number of independent measurements. It should be noted that this is the first time undersampling with near-field radiation has been performed. In this section, we investigate two strategies to reconstruct images using undersampled sets of measurements: adaptive sampling and compressed sensing.

Using adaptive sampling, we first measure a low-resolution image and then sample regions of interest with progressively higher resolution. In short, identification of coarse edges from this initial low-resolution image determines where to sample with higher resolution, thus reducing the total number of measurements that are made [14]. Edge identification is achieved via a single-tier 2D Haar wavelet decomposition of the low-resolution image. The Haar wavelet transform is a hierarchical structure that highlights the presence of edges at progressively finer scales: edge features yield large wavelet coefficients, while more uniform areas yield low wavelet coefficients [25]. In our experiment, after each higher resolution resampling phase, edge detection is performed on the new image, and the process is repeated until the required resolution is reached. The algorithm is described in detail in Supplement 1, Section 5.

Compressed sensing is an alternative non-adaptive approach that makes assumptions about the object to reduce the total number of measurements. This can include knowledge of the basis in which the representation of the image is sparse [13], or the basis in which the total variance or curvature of the object is expected to be low. Such assumptions typically hold for a wide variety of images, and compressed sensing theory shows how to then recover an image using fewer measurements than the number of pixels in the image [i.e., $\mathbf{A}$ in Eq. (1) has fewer rows than columns] [13]. In this work we make the assumption that the total curvature of the image will be low. We sample the object using a set of random binary patterns, which are chosen due to their high degree of incoherence with respect to a wide range of sparse basis representations. Non-adaptive compressed sensing has the advantage that masks can be designed (and loaded onto a modulator) ahead of time rather than in response to measurement, as is the case for adaptive sampling.

Note that, in both our sampling strategies, we have applied regularization to combat noise in our measurements. This essentially allows the reconstruction algorithm to permit solutions that deviate from the measurements by an amount based on the estimated noise level, while seeking to minimize the level of curvature in the reconstruction. More details are given in Supplement 1, Section 7, alongside unregularized images.

Figure 3 compares reconstructed images of a transmissive object depicting two of Maxwell’s equations as the number of measurements for both adaptive imaging and compressive sensing are reduced. The main features of the object are clearly recovered even when the number of samples is reduced to 35% of the Nyquist limit. These images also display directional contrast modulation due to polarization effects (see Supplement 1, Section 2). Comparing the two methods, we observe that adaptive sampling marginally outperforms compressed sensing, yielding higher contrast features. Note that our aim is to demonstrate the feasibility of these techniques rather than to optimize the algorithms. We anticipate improvements in image quality should one inject more prior knowledge. We also note that there are a wide variety of algorithms available for both adaptive sampling [14,15] and compressive sensing [13,16] that are optimized for particular scenarios.

 

Fig. 3. Compressive $72\times 72\mathrm{THz}$ images of $\nabla ·\mathbf{E}=0$ and $\nabla ·\mathbf{B}=0$ (see Supplement 1 for optical image) with decreasing number of measurements, where the top, middle, and bottom rows, respectively, use 75%, 50%, and 35% measurements as the number of pixels. Because of the low repetition rate of our laser, a full measurement set takes 5 h. Left column: adaptive sampling. Right column: compressed sensing.

Download Full Size | PDF

In conclusion, we have demonstrated a subwavelength THz imaging technique that is compatible with adaptive and compressive sensing algorithms. We photoexcite a 6-μm-thick silicon wafer with a sequence of optical patterns to spatially modulate an incident THz beam. By placing an object at the exit interface of our wafer, we demonstrate imaging at $9(\pm 4)\mathrm{\mu m}$ resolution and observe strong polarization effects at the interface of a conductor. We conclude by comparing two strategies to reconstruct an image from an undersampled set of measurements: adaptive sampling and compressed sensing. To the best of our knowledge, this is the first experimental demonstration, in any spectral regime, of undersampled image reconstruction at highly subwavelength dimensions. It is not obvious that these multipixel sampling approaches would work well in the near field: an implicit assumption is that the transmission through each pattern depends only on the open area and not on aperture shape. Our results demonstrate that this assumption is surprisingly robust down to a scale of $\lambda /45$. Such approaches therefore hold promise for improving signal to noise in other near-field and low-signal imaging techniques.