Multifrequency Digital Terahertz Holography Within 1.39–4.25 THz Range

Terahertz (THz) multifrequency digital holography within the range from 1.39–4.25 THz is demonstrated. Holograms are recorded using an optically pumped molecular THz laser operating at emission lines of 1.39-, 2.52-, 3.11-, and 4.25 THz frequencies, and nanometric field effect transistor with integrated patch antennae as a THz detector. It is revealed that phase-shifting methods allow for qualitative reconstruction of multifrequency THz holograms combined into one “colored” image. It provides more information about the low-absorbing objects with additionally improved quality achieved by removing unwanted information related to the so-called dc term and conjugated beam forming a virtual image. It is shown that the THz holography can be applied for the investigation of low-absorbing objects, and it is illustrated via inspection of stacked graphene layers placed on a high-resistivity silicon substrate.

or monitor the quality of materials via induced phase shifts.A particular role in the family of coherent imaging techniques can be attributed to THz holography due to its fascinating capability to reconstruct 3-D images and due to submicrometer longitudinal resolution extending hence the borders of practical uses into a possible observation of subwavelength-size defects or nonuniformities in inspected objects [9].
Although digital holography [10] is well established at optical wavelengths [11], [12], [13], [14], its implementation in the THz range imposes a challenging task due to the remaining issue of low power of THz radiation sources and a lack of reliable sensitive detectors operating at room temperature.Furthermore, another important concern is the relation between the wavelength, the size of the optical elements, and the scanned area, which for the THz range is unfavorable.Optical setups constructed for the THz range of radiation suffer from large diffraction effects and undesirable influence of limited apertures of optical elements, which are tremendously more significant if compared to a visible range.The lack of reasonably sized detector arrays with relevant pixel pitch and large enough element diameters resulted in point-to-point imaging: focusing the light on the sample and then on the detector with subsequent shifting of the sample.Such configuration enables using only the radiation propagating around the optical axis which diminishes unwanted diffraction effects.It should be underlined that plane-to-plane imaging requires much larger apertures of all used optical elements.In the case of registering THz holograms, both approaches are possible and additional methods [like phase-shifting (PS) techniques [15], filtering in Fourier domain [16], [17], or many others] of removing unwanted components from reconstruction are required.
Up to now, a large variety of interesting efforts and different approaches dedicated to the development of THz holography [18], such as in-line holography [19], [20], off-axis holography [21], [22], [23], near-field holography [24], and time-reversal holographic imaging of hidden objects [25] was demonstrated.Innovations in THz holography via step PS methods allow for improved quality of the recorded holographic images [15], [26], [27].Very recently, lensless Fourier-transform THz digital holography for full-field reflective imaging was Fig. 1.Optical setup for the recording of THz digital holograms at 1.39-, 2.52-, 3.11-, and 4.25-THz.Letter "E" denotes THz emitter: optically pumped molecular THz laser; PM 1 , PM 2 , and PM 3 are parabolic mirrors; C is a chopper; M 1 , M 2 , and M 3 label gold-coated flat mirrors; BS 1 , BS 2 , and BS 3 mark high-resistivity silicon beam splitters; S is the imaged sample; thin napkins, referred to in the rest of this article as paper sheets served as a simple phase shifter; D 1 and D 2 are detectors.The enlarged part from the D 2 illustrates the Si CMOS chip with TeraFET detectors, which consisted of a transistors and integrated patch antennae.A different detector is used for different THz source frequencies, which are marked with yellow squares.A signal intensity dependence on the number of paper sheets, added to a PD is depicted on the left side for four frequencies.A transparent yellow line marks the minimum of the normalized signal.The example of the registered holograms at four different THz frequencies is given in the left coupled panel.A table in the top right corner presents the corresponding phase shift induced by one paper sheet for each frequency.The inset at the bottom-right corner depicts a sample made of aluminum foil with a cut-out letter "E" letter shape.Parts of the sample are covered by different numbers of paper sheets to produce different phase shifts indicated here by the different colors.demonstrated using an optically pumped molecular THz laser operating at 2.52 THz [28].It is worth noting that holographic images, as a rule, are recorded using one wavelength and forming therefore monochromatic images.
In the given work, we enrich the family of THz holographic experiments by demonstrating multifrequency digital holography ranging from 1.39-4.25 THz and unveil its possible application for the low-absorbing objects via inspection of stacked graphene layers.Relying on the setup based on an optically pumped molecular THz laser operating at discrete emission lines of 1.39-, 2.52-, 3.11-, and 4.25-THz frequencies as well as high THz sensitivity of a detector based on nanometric field effect transistors [29] we demonstrate that PS methods [15] allow to qualitatively reconstruct multifrequency THz holograms combined into one "colored" image.The THz holograms display an additionally improved quality, especially, due to the phase retrieval achieved by removing unwanted information related to the so-called dc term and conjugated beam forming a virtual image.Phase values corresponding to the phase shifts introduced by the object were retrieved at different wavelengths providing thus a way for a deeper and more qualitative inspection of the investigated materials.

II. HOLOGRAMS RECORDING SETUPS AND RECONSTRUCTION METHODS
Holograms Recording: Holograms are recorded by employing the experimental optical setup relying on the Mach-Zehnder interferometer depicted in Fig. 1.The coherent THz emission was generated by an optically pumped continuous wave molecular THz laser (FIRL100, Edinburgh Instruments, marked as letter "E" in the setup scheme).The radiation at frequencies of 1.39-, 2.52-, 3.11-, and 4.25 THz within the power range of 2-8 mW was focused using a parabolic mirror with f = 5 cm (PM 1 ) and, after a chopper modulation, was collimated by PM 2 .The collimated THz beam reflected by a flat gold-coated mirror (M 1 ) was then divided into two parts by the 525-μm-thick highresistivity silicon beam splitter (BS 1 ).One part was focused on the pyroelectric detector (D 1 ) by the third parabolic mirror (PM 3 ) and served as the laser power reference (marked in violet), which was later used to compensate for laser mode instabilities.The beam splitter (BS 2 , equivalent to BS 1 ) divided the second part of the beam into two interfering beams marked by yellow and blue colors.The beam in the first arm was reflected from the second flat mirror (M 2 ), passed through the sample (S), and, subsequently, after passing through the third beam splitter (BS 3 ), reached the second detector.The beam in the second arm also was reflected from the third flat mirror (M 3 ), and its phase was shifted by the introduced phase delay element (PD).Consequently, after reflecting from BS 3 , it interfered with the beam from the first arm and also reached the D 2 .The transmitted beam power was registered using a nanometric field-effect transistor-based Thz detector (TeraFET) with integrated patch antennae.Superstrate dielectric and lens are not used as the THz emission is incident directly on the antennae of the detector.A different detector is used for a different THz frequency.All Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.
devices are integrated on a Si CMOS chip (as shown in the inset of Fig. 1) fabricated using a standard 65-nm CMOS technology.The detectors differ in transistor parameters and antennae sizes, which determine their resonant frequencies: 1.39 THz (used with 1.39-and 4.25-THz sources), 0.8 THz (for 2.52 THz source) and 0.98 THz (for measurements at 3.11 THz).The minimum noise equivalent power (NEP) for individual detectors spanned from 32 pW/ √ Hz at 1.39 THz [30] to 75 pW/ √ Hz at 0.98 THz.The operation mechanism and peculiarities of the design of such detectors are described elsewhere in detail [29].Both the direct detection and the sensing at the third-order resonance of the antennae were employed in the signal registration via the lock-in technique.
Holograms Reconstruction: The intensity pattern was recorded as the interference of two beams: 1) the object beam U obj = A obj exp(iφ obj ) and 2) the reference beam Both beams are formed by light field distributions U reaching the hologram plane and both contain information about their amplitude A and phase φ distributions.Adding these two beams forms an interference pattern U obj (x, y) + U ref (x, y) that is later recorded by the detector as intensity distribution I h (x, y).It forms a Fresnel hologram which can be numerically reconstructed.The intensity pattern recorded by the detector hologram is given as which results in forming four terms.Two of them (terms 1 and 2) do not contain any phase information and form the so-called dc term, one of them (term 3) corresponds to forming a virtual image, and the other one (term 4), to a real image.It is worth noting that only the last one is interesting for us in the case of digital reconstruction, because it forms a real image in the numerical plane corresponding to the object plane in a real setup.Thus, all other terms are forming unwanted noise, dc term propagating back and divergent light field distribution coming from the virtual image.The light field distribution forming the real image can be extracted from all registered information using PS methods [31].Such methods assume recording different numbers of successive holograms as interference of constant object beam and the reference beam undergoing particular phase shifts between the following recordings.Depending on the applied method a total or only partial removal of the unwanted terms can be obtained.In this study, a four-step PS method (introducing 0, 0.5π, π, and 1.5π phase shifts in the reference beam) is used as described in our previous investigation [15].For comparison, a single hologram and a two-step PS reconstruction, introducing 0 and π PS in the reference beam, were also conducted.The four-step reconstructing algorithm utilizes for back-propagation the intensity pattern given by the equation where I h1 (x, y), I h2 (x, y), I h3 (x, y), and I h4 (x, y) are the subsequently recorded holograms corresponding to intensity distributions being a result of two interfering beams with introduced PS of 0, 0.5π, π, and 1.5π in the reference beam, respectively, which corresponds to multiplying each hologram by 1, i, −1, and −i values, respectively.Such smartly manipulated intensity pattern I PS4 (x, y) is back-propagated to the image plane (the same where the object was located).The light field distribution in the image plane will consist here only of the term related to the object -U obj = A obj exp(−iφ obj ).
In the case of a two-step PS algorithm, the algorithm simplifies to Here, the reconstructed light field distribution in the image plane will consist of the term related to the object U obj and the term related to a virtual image.It should be taken into account that here the light field distribution forming the virtual image is introducing additional and unwanted noise in the reconstructed image plane.In the case of single hologram reconstruction, only the back-propagation of the recorded intensity pattern forming a hologram is conducted and results in the reconstruction in the image plane of all four terms described by the (1).
As we describe digital holograms it is obvious that the recording is realized in a real experimental setup, while reconstruction is performed numerically.However, one can also perform simulations of the experimental setup, achieving a simulated hologram, to compare the simulations' results with later obtained experimental data.Both holograms (simulated and recorded) form an interference pattern that is later reconstructed numerically.The numerical simulations were conducted for three different methods: 1) single hologram reconstruction, 2) PS2, and 3) PS4 for each of four frequencies used later in the experiment: 1.39-, 2.52-, 3.11-, and 4.25-THz.
The object consisted of the letter "E" with uniform amplitude and different phase changes corresponding to the experimentally recorded one.It should be noticed that the real object had different numbers of papers in three areas of the letter "E" introducing particular phase shifts for each frequency (shown in Fig. 1 and precisely described in the next section).Thus, in the simulation, a set of four objects was created, each having the same amplitude distribution but different phase levels corresponding to the particular phase shifts in the experiment.For each frequency, an appropriate object was illuminated with a plane wave and propagated to the hologram plane (at a distance of 95 mm corresponding to the experimental value) and formed an object beam U obj (x, y).Then, in the hologram plane, a reference beam U ref (x, y) was added.It was assumed that two beams interfere in the hologram plane at an angle of 0. Then, this pattern was converted into intensity, which corresponds to the registration with the detector (registering signal proportional to the intensity of radiation).Such a hologram was then propagated to the image Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.
plane to perform reconstruction.The simulations for a similar setup for the single frequency were already described in [15].The methodology of single frequency is not different from previously described, although all simulations had to be conducted for 214.14-, 118.08-, 96.71-, and 69.72-μm wavelengths and numerical objects designed separately for each frequency.

III. RESULTS AND DISCUSSION
The investigated sample was raster scanned with a spatial resolution of 0.5 mm in x and y directions.A special sample dedicated to illustrating the holographic images was fabricated from metallic foil of 55 × 85 mm 2 size with the "E" letter-shaped aperture (shown at the bottom of the right corner of Fig. 1).The vertical part (marked in white color) was completely hollow.Meanwhile, the horizontal parts were covered with 1 (marked as yellow), 2 (marked as blue), and 3 (marked as purple) paper sheets (thin napkins) that shift the phase 0.2π, 0.36π, 0.41π, and 0.55π for the frequency 1.39-, 2.52-, 3.11-, and 4.25-THz, respectively.Particular thicknesses of one paper sheet were measured with a caliper and thus the refractive index values: 1.177, 1.175, 1.165, and 1.160, respectively, were determined via time-domain spectroscopy measurements.Due to the fact, that in the THz time-domain spectroscopy signal-to-noise ratio decreases for higher frequencies, additional measurements were conducted, and the signal intensity dependence on the number of paper sheets was calculated.A transparent yellow line in Fig. 1 marks the minima of the normalized signal: the corresponding values are: 0.14π, 0.33π, 0.4π, and 0.5π, respectively, for each frequency.
As one can see, as the frequency increases, the PS of a single paper layer also increases, indicating that a few paper sheets are needed to achieve a complete 2π PS.For example, at 1.39 THz 2π PS can be obtained using 14 paper sheets, at 2.52 THz the number of paper sheets decreases to six, at 3.11 THz five paper sheets are needed and, finally, at 4.25 THz only four paper sheets are required.These estimates are given in the table in the top-right corner of Fig. 1.The left panel demonstrates a shift of the local signal minima on the variation of several paper sheets at different frequencies.
A simple PS element fabricated of a different number of the paper sheets (thin napkins) that were pressed together without the use of any binding material, was placed in the second arm of the interferometer (reference beam) and was used to induce a PS, as one paper sheet corresponds to 0.2π, 0.36π, 0.41π, and 0.55π (or 0.14π, 0.33π, 0.4π, and 0.5π, according to the experimental evaluation of the intensity dependence on the number of paper sheets), for the frequency 1.39-, 2.52-, 3.11-, and 4.25-THz, respectively.For simulations we assumed phase-shift values calculated on the basis of TDS measurements.In the experiment for each frequency eight holograms were recorded with different numbers of paper layers in the reference beam.During reconstruction, different sets of four holograms with different phase shifts were taken to verify which of them match the proper phase shifts the most.
Corresponding intensity distributions of holograms are presented in the left coupled panel in Fig. 1 and for all PS in the reference beam in Fig. 2. It is worth noting that during these measurements, no additional PS was placed in the object beam, and the PS was introduced only by the sample parts that were covered with one, two, and three paper sheets.As it is seen, at 1.39 THz frequency, all four different areas of the imaged sample are visible.When the frequency increases to 2.52-and 3.11-THz, only three areas of the sample are detectable in the hologram, while at 4.25 THz the intensity distribution in the hologram becomes less pronounced due to increasing attenuation of the THz radiation intensity with the increase of the number of PS paper sheets.
The experiment was repeated for each frequency, with an increasing PS achieved by adding an extra sheet of paper to the reference beam.The 2-D hologram intensity distributions are shown in Fig. 2(b).These results are in good agreement with the already discussed dependence of the signal intensity on the PS.It is also visible that for the consecutive holograms the attenuation introduced by the paper PS is increasing with a number of the paper sheets and with recording frequency as well.
As it is seen, the "E" letter-shaped sample is well resolved in all recorded holograms except 4.25 THz.Varying the PS, its different parts become more pronounced.When the PS reaches 2π, hologram distributions look similar to that without a PS, but with a strongly reduced intensity due to the introduction of paper sheets which except PS introduce also the attenuation for higher frequencies.This effect correlates with the intensity and its dependence on the PS given in Fig. 2(c)-(f) panels.It should be also underlined that the decrease in the intensity of the reference beam results in a larger difference in the intensities of the interfering beams and, consequently, the lower contrast of the recorded hologram.
Moreover, the contrast values were estimated for each separate area of the imaged sample in the hologram and results are represented in Fig. 2(c)-(f).The contrast was determined by C=(I a -I n )/I n , where I a is the average signal intensity in the specified area (marked with 0, 1, 2, or 3) and I n is the noise or average intensity of the metal covered area.
For 1.39 THz, as indicated previously, the smallest of all four radiation frequency cases PS is observed with an additional 7 paper sheets.For this reason, contrast graphs of one, two, and three areas lie quite close to the zero area graph.For each frequency, all the sample areas exhibit the same tendency.
After recording the holograms in the detector plane, each of them was reconstructed numerically.The four-step PS algorithm was used to reconstruct the information only about the recorded object [27], [32].For comparison, the reconstructions were also conducted for the two-step PS algorithm [16] as well as the single hologram reconstructions.The modified convolution approach [33] was used to back-propagate the light field distribution to the object plane.Each hologram was scanned with 0.5-mm resolution having a total size of 54 × 54 mm.For the reconstruction, the whole scanned hologram area was placed inside a 1024 × 1024 pixel matrix which was rescaled four times to increase the reconstructed image resolution.The digital reconstruction was performed on a 4096 × 4096 pixels calculation matrix with a sampling of 0.125 mm.The registered intensity pattern was recalculated according to the desired algorithm (four-step and two-step PS, or single hologram reconstructions) and propagated back at a distance of 95 mm corresponding to the distance between the object and the hologram planes.The propagation was conducted using an off-axis approach which does not contain any approximations.From each intensity pattern of the hologram, amplitude, and phase distributions in the image plane obtained for all frequencies and reconstruction algorithms are shown in Fig. 3.
Before proceeding further, it is worth recalling that the sample was a transparent letter "E" of which parts were covered with Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.Fig. 3. Amplitude (top part) and phase (bottom part) distributions of the reconstructed holograms at four different frequencies: 1.39-, 2.52-, 3.11-, and 4.25-THz.In each case, the reconstructions in the first row are obtained using a single experimentally measured hologram.The reconstructions in the second row are obtained using the two-step PS method, while the reconstructions in the third row are obtained using the four-step PS method.Colored scales are normalized to the maximum values and are different for each frequency.The amplitude distributions of obtained images are similar in all cases, however, the phase distributions experience strong differences.Due to the removal of unwanted terms in the reconstructed light field distribution: the PS introduced by different areas of the letter "E" starts to be visible and distinguishable.Note that only for the four-step PS method all areas have different phase levels in the reconstructed image, especially for 1.39 and 2.52 THz.In the case of 3.11 and 4.25 THz, reconstructed phase distributions in different areas are less distinguishable.The rippled phase image at 3.11 THz can be associated with the destructive interference effect of the second harmonics of the detector which displays a resonance frequency of 1.4 THz.The phase image at 4.25 THz is blurred because of the THz light scattering comparable with the wavelength roughness of the paper sheets used as the PD in the reference arm.Note that the two-step PS and the four-step PS methods allow to resolve different phase levels.Especially, the latter shows the most significant difference between phase level values for different areas of the sample.a consecutive number of paper sheets.They introduce particular PS, and the paper tissue itself introduced some additional attenuation, for higher frequencies in particular.Thus, in the reconstructed images of our object, one should expect a relatively uniform amplitude pattern: dependent only on the influence of the attenuation of the paper sheets, and the phase values corresponding to the different PS introduced by the sample.For comparison, the contrast of amplitude and phase reconstructions for holograms are recorded in all studied frequencies.The reconstructions were obtained from a single (one-step) hologram and using two-step and four-step PS techniques.Results are presented in Fig. 4. As one can see, the amplitude contrast increases by applying two-step PS and four-step PS for 1.39 THz.With the increase in frequency, the contrast diminishes for all techniques used.Application of two-step PS and four-step PS allows removal of additional interferences and makes the amplitude more uniform, but does not remove interferences from some parts of the setup.Regarding the phase contrast, the four-step PS technique gives the best results.Nonmonotonical phase contrast behavior at 2.52 THz can be associated probably with unwanted interference effects in the optical setup.
To showcase the benefits of multifrequency holography, we have merged images obtained at different wavelengths into a single "colored" image.Initially, theoretical simulations were conducted by generating a virtual sample, simulating its hologram distribution at the selected frequencies, and subsequently performing reconstructions.The calculation matrix was 1024 × 1024 pixels with the sampling equal to 125 μm which corresponded to 128-mm square area.Only the central parts of the reconstructed areas (with the letter "E") are shown in the combined image illustrated in Fig. 5.The results were derived from three frequencies joined into the single image, employing a methodology analogous to the RGB color model.In this model, the signal distribution of the reconstructed hologram at 1.39 THz corresponds to Red, 2.52 THz to Green, and 3.11 THz to Blue channels.The image obtained at 4.25 THz was not taken into account as this wavelength is close to the paper roughness in the sample.As a consequence, the THz radiation becomes diffused here, leading to a blurred image.
The simulations and subsequent RGB combinations were replicated for both amplitude and phase distributions across all three cases: without phase shift (one-step PS), and using two-step and four-step PS techniques.The results are exposed on the left side of Fig. 5.The RGB combinations were also produced from experimentally obtained images and are provided on the right side of the same figure, respectively.As one can see, the simulations and the experiments display very similar results.
In the case of amplitude distributions, when the PS technique was not introduced (one-step PS), areas with varying phase levels exhibit clearly pronounced differences.In simulations, sample areas with zero, one, two, and three paper sheets are sharply defined and highly distinctive in color.The same behavior is observed for reconstruction from experimental data.Moreover, the reconstruction by the two-step PS method also manifests in nonuniform amplitude distribution inside the letter "E".
However, with the introduction of the four-step PS, both in simulations and experiments, the resulting combined image appears quite uniform.This can be explained by considering Figs. 3  and 4. In the case of one-step PS, at all three frequencies areas with different phase levels are notably distinct from each other further manifesting itself in the combined image.Meanwhile, when using the four-step PS, at all frequencies the amplitude distribution in all areas is more uniform, making hence the combined image homogeneous.This result can be expected as the sample consisted of uniform amplitude and areas with different phase levels.Uniform amplitude distribution is visible only in the case of the four-step PS algorithm reconstruction.Different colors in the amplitude for two other methods result from the fact that we observe interference (constructive and destructive) of the beam forming the image of our object with an additional one.In these two cases, reconstructing the hologram creates more than one beam: in the case of the two-step PS additional beam occurs from term 3 in (1); in the case of single hologram reconstruction additionally, three beams are formed corresponding to terms 1, 2, and 3 in (1).Thus, the reconstructed amplitude distributions for one-step and two-step PS are the result of the object beam interfering with additional unwanted beams resulting in different color areas, instead of uniform (white).The proof of such interpretation is visible in Fig. 3 while observing amplitude one-step PS reconstruction for 2.52 THz: the vertical bar of the letter "E" is empty (no paper layer) so theoretically should be the brightest because it introduces no attenuation, but is the darkest from all areas.This supports the interpretation that in one-step and two-step PS methods we see the result of two-and four-beam interference patterns, respectively.
Conversely, the opposite outcome is observed in the case of phase.When employing one-step PS, a blurred and highly uniform image is produced.This could also be explained by the outcome depicted in the graph in Fig. 4, where the phase contrast remains rather consistent across different frequencies with the one-step PS.Simultaneously, when employing the four-step PS, the areas introducing different phase shifts are very distinct, especially in the case of the experiment.
The quality of reconstructions is also affected by the fact that the phase shifts introduced by different numbers of paper layers in the experiment in the reference beam do not correspond to the exact phase values used in the algorithm.Registering holograms with slightly different phase shifts as desired may result in additional disturbances in the reconstructed image.Nevertheless, the reconstruction algorithm is performed numerically and proper phase shifts are used.Moreover, the simulations carried out as comparison (Fig. 5) were assuming exact phase shifts both in recording and reconstruction algorithms.Taking into account the large accordance between simulation and experimental results, it can be assumed that the influence of imperfect phase shifts during recording holograms is rather small and more prominent and disturbing is visible noise in the experimental evaluation.
Combined RGB images, particularly, those representing the four-step PS, offer an advantage over the images acquired using a single frequency because in this case, areas with different phase levels have uniform amplitude, which is expected for very lowabsorbing samples.Simultaneously, the phase distribution at the four-step PS allows more pronounced distinction between the parts of the sample at different frequencies.For instance, in the experimental case, the area with 1 paper sheet is more observable at 2.52 THz (green color in the merged image), the area with 2 paper sheets is more distinguished at 3.11 THz (blue color), and the area with 3 paper sheets is more noticeable at 1.39 THz (pink color).Thus, the phase distribution consists of a pink vertical bar of the letter "E", and each horizontal bar is represented by a distinct color (green, green with blue, and blue with pink) meaning different values of the phase which is also consistent with the used sample.
Overall, the impact of phase information introduced by a low-absorbing sample on the reconstruction process manifests itself in the differences between the three different methods.It appears that the PS4 algorithm is indeed the best choice for reconstructing low-absorbing samples, especially in in-line configuration.
It is also important to acknowledge the presence of visible diffraction effects, despite which the four-step PS technique manages to isolate the phase information.This highlights the algorithm's effectiveness in dealing with low-absorbing samples.
To illustrate the suitability of THz holograms for lowabsorbing materials inspection, we performed a holographic study using samples of stacked graphene layers placed on the Fig. 7. Amplitude (top part) and phase (bottom part) distributions of the reconstructed holograms at two different frequencies: 2.52-and 3.11-THz.The imaged sample consists of a different number of graphene layers (from 1 to 5) placed on the silicon substrate.In each case, the reconstructions in the first row are obtained using a single experimentally measured hologram image.The reconstructions in the second row are obtained using the two-step PS method, while the reconstructions in the third row are obtained using the four-step PS method.The colored scales are normalized to the maximum signal values.It needs to be underlined that using the four-step PS method a phase distribution corresponding to the phase delay map introduced by the sample is visible and recognizable in the phase distribution of the image.At the same time, the amplitude distribution remains at a constant level which stays under the assumptions and the experimental results.
high-resistivity silicon substrate, as shown in Fig. 6(a).As it is known, the thickness of the graphene monolayer is only 0.3 nm, hence, its inspection in the THz range is challenging due to the low absorption in such thin material.It is worth noting that THz imaging together with simultaneous carrier optical modulation can be a rational way to test the material [34]; however, it requires additional optical excitation.Recently, it was demonstrated that THz structured light in compact THz imaging systems can serve as a convenient tool for contactless inspection of stacked graphene layers [35].
Here, we further extend the research boundaries for techniques suitable for nondestructive inspection of 2-D materials.In multispectral digital THz holography experiment, we used five samples with different numbers-1, 2, 3, 4, and 5 layers of graphene-transferred on a high resistivity silicon substrate.The silicon substrate itself (without graphene) served as a reference sample.Results are depicted in Fig. 6(b).As one can see, the optimal frequency for distinguishing numbers of graphene layers can be observed in Fig. 6(b), which corresponds to 3.11 THz.Images of the investigated samples when the PS using two, four, or six paper sheets was introduced to the reference beam are depicted in Fig. 6(c) for 2.52 THz and Fig. 6(d) for 3.11 THz frequencies.It is seen that the latter exhibits better-resolved images in comparison to all others.
For the graphene samples, the four-step and two-step PS algorithms were used for reconstruction, and, for comparison, the single hologram reconstructions were also performed.Each hologram was scanned with 0.5-mm resolution to cover the total size of 18 × 80 mm.For the reconstruction, the whole scanned hologram area was placed inside a 1024 × 1024 pixel matrix which was rescaled four times to increase the reconstructed image resolution.The digital reconstruction was performed on a 4096 × 4096 pixels calculation matrix with a sampling of 0.125 mm.The registered intensity pattern was recalculated according to the desired algorithm (four-step, two-step PS, or single hologram reconstructions) and propagated back at a distance of 95 mm corresponding to the distance between the object and the hologram planes.The amplitude and phase distributions in the image plane for hologram reconstructions for all frequencies and algorithms are shown in Fig. 7.
As one can see, for a frequency of 3.11 THz the reconstructed image the amplitude distribution is sensitive to the PS methods used.It is seen that the amplitude distribution for the four-step PS algorithm is uniform, while its phase distribution for the equivalent case delivers the most distinguishable areas which is straightforwardly related to the introduced PS of the transparent sample with varying thicknesses.In Fig. 7 it can be seen that the amplitude distribution reconstructed from the single hologram for the frequency of 3.11 THz is characterized by variation in amplitude areas corresponding to the different number of the graphene layers.Such a result indicates that the graphene layers form an interferometric pattern that enables them to be distinguished in the amplitude distribution.However, as previously for the letter "E" sample, this is the result of the interference coming out from reconstructing the object beam (having information about the sample) with additional beams treated as noise [terms 1 -3 in (1)].All four beams are reconstructed from the same hologram distribution for the one-step method, and two beams [term 3 and 4 in (1)] are reconstructed using the two-step PS method.It should be noted that the phase distribution corresponding to the phase delay map introduced by the sample is visible and recognizable in the phase distribution of the image reconstructed with the four-step PS algorithm.We expect the sample to be transparent with the introduced particular PS and, therefore, the only four-step PS method gives the expected results without additional interferences.
In the case of the 2.52 THz frequency, the reconstructed amplitude distributions are much more uniform, and the PS algorithms allow to reveal even small phase differences introduced by the sample.The more advanced is the four-step PS, which removes all the unwanted components related to the additional noise and interferences in the reconstructed images.

IV. CONCLUSION
To summarize, multispectral digital THz holography based on an optically pumped molecular THz laser operating at emission lines of 1.39-, 2.52-, 3.11-, and 4.25-THz frequencies and nanometric field effect transistor-based sensors fabricated by 65-nm CMOS technology is demonstrated.It is revealed that the two-step and four-step PS methods allow for qualitative reconstruction of multifrequency THz holograms combined into one "colored" image and furnish more information with additionally improved quality on low-absorbing objects achieved by removing unwanted information related to the so-called dc term and conjugated beam forming a virtual image.It is shown that holographic reconstructions can be applied for the investigation of low-absorbing and thin objects including samples containing graphene layers placed on the high-resistivity silicon substrate.Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.

Fig. 2 .
Fig. 2. Experimentally recorded intensity distribution of "E" letter-shaped sample at different frequencies with varying PS.(a) Scheme of the imaged sample consisting of different areas covered with a varying number of paper sheets.The grey area is a metal foil, the white area marks the hollow part, the yellow area is covered with 1 paper sheet, the green area with 2, and the purple area with 3 paper sheets.(b) Experimentally obtained 2-D intensity distributions at four frequencies: 1.39-, 2.52-, 3.11-, and 4.25-THz.The image was recorded seven times placing an additional paper sheet to the reference beam to introduce a corresponding PS.The colored scale is normalized to the maximum signal.(c)-(f) Contrast dependencies on the number of paper sheets placed in the reference beam (top scale) or corresponding PS (bottom scale) at four different frequencies: (c) 1.39 THz, (d) 2.52 THz, (e) 3.11 THz, and (f) 4.25 THz.In each case, grey, yellow, green, and purple colored lines with symbols represent the contrast values at corresponding areas of the sample marked with the same colors.Shifting the phase in the reference beam by π, the contrast decreases by 96% for the 1.39-THz frequency, by 57% for 2.52 THz, by 59% for 3.11 THz (here, the phase shift is equal to 4π/5) and by 98% for 4.25 THz.

Fig. 4 .
Fig. 4. Contrast of reconstructed amplitude and phase distributions for the sample "E" at different frequencies: 1.39-, 2.52-, 3.11-, and 4.25-THz.Reconstructions were obtained from the single (1-step) hologram and using two-step and four-step PS techniques.The contrast was obtained from the area "0" without paper sheets.The results of amplitude contrast are represented in full circles and straight lines, while the phase contrast results are marked as open circles connected with dashed lines.As one can see, the amplitude contrast deviations are observed in 1 step PS, indicating insufficient accuracy.The unusual behavior in the amplitude contrast 3.11 THz can be associated with the destructive interference effect of the second harmonics of the detector which displays a resonance frequency of 1.4 THz.

Fig. 5 .
Fig. 5. Amplitude (top panel) and phase (bottom panel) distributions of the combined holographic images obtained at three different frequencies-1.39-,2.52-, and 3.11-THz-across all three cases: without phase shift (one-step PS) and using two steps (two-step PS) and four steps (four-step PS) PS techniques.Signal distribution in the reconstructed hologram at 1.39 THz corresponds to Red color, 2.52 THz to the Green one, and 3.11 THz to Blue channels.The distributions depicted on the left side illustrate the outcomes derived from the simulations, whereas the images on the right side portray the experimental case.

Fig. 6 .
Fig. 6.Holographic images of stacked graphene layers placed on a silicon substrate.(a) Photograph of the graphene samples where the different number of graphene layers (from 1 to 5) were transferred on the high-resistivity silicon substrate.The inset on the right depicts part of the experimental setup, where S is the imaged sample, M 2 and M 3 are flat mirrors and BS 2 , BS 3 are beam splitters.(b) Experimentally obtained hologram intensity distributions imaged at four different frequencies without the introduction of the PS in the reference beam.The colored scales are normalized to the maximum values.The contrast value for 1.39 THz is 0.51, for 2.52 THz 0.35, for 3.11 THz 0.83, and for 4.25 THz 0.66.(c) and (d) Hologram intensity distributions at (c) 2.52 THz and (d) 3.11 THz frequencies when the reference beam phase is shifted by placing 0, 2, 4, and 6 paper sheets.The resulting PS is indicated inside each image.Note that the sample details are well resolved at 3.11 THz.

Rusnė
Ivaškeviči ūtė-Povilauskienė received the B.S. degree in physics and M.S. degree in material science and semiconductor physics from the Physics Faculty, Vilnius University, Vilnius, Lithuania, in 2015 and 2018, respectively, and the Ph.D. degree in physical sciences acquired through a joint program between the Center for Physical Sciences and Technology (FTMC), Vilnius, Lithuania, and Vilnius University, in 2023.She is currently a Researcher with THz Photonics lab, FTMC.Her research interests include THz imaging, development and research of compact optical elements, and exploration of 2-D materials.