High-harmonic generation (HHG) is capable of producing attosecond pulses, which find widespread usage in observing ultrafast dynamics in atoms, molecules, and solids.1^–8 Given their short pulse width, high peak intensity, and broad spectral width, few-cycle pulses serve as common drivers for generating attosecond pulses.9 Consequently, the accurate characterization of the few-cycle electric field is essential for comprehending and effectively utilizing attosecond pulses in subsequent research.

Conventional pulse diagnosis techniques such as frequency-resolved optical gating (FROG)10 are generally not sensitive to the carrier-envelope phase (CEP) of a laser pulse, failing to access the temporal waveform, i.e., electric field, of the laser pulse to be measured. Over the past decades, several techniques have emerged for diagnosing the electric field of few-cycle pulses. For instance, the method of attosecond streaking camera11 based on photoelectron spectroscopy has the capability to map the vector potential of the streaking field into the delay-dependent momentum shift of the photoelectron. There are also several all-optical methods, such as optoelectronic sampling,12^,13 photon excitation in solids,14^,15 tunneling ionization in a gaseous medium,16 and transient absorption spectroscopy of atoms,17^,18 that offer potential avenues for sampling the waveform of the electric field. Interferometry constitutes another significant branch of technology, utilizing the fundamental property of light interference to enable various high-precision detection methodologies. Notably, it has promoted the detection of gravitational waves19^,20 and the validation of core principles in quantum mechanics.21^,22 This high-precision approach has been extended to the characterization of few-cycle pulses. The approaches involve introducing a perturbing field to influence the spatial interference pattern of the high harmonics (HHs)23 or to alter the frequency shift of the HHs for pulse diagnostics.24^,25 Other noninterferometric approaches based on the perturbing HHG process have also been proposed for ultrafast electromagnetic field characterization.26^,27 Kim et al.26 utilized a noncollinear weak field to perturb the spatial wavefront of HHs, leading to modulation of the propagation angle of HHs, which can be used to reconstruct the waveform of the perturbing field. Nevertheless, although effective in reconstructing temporal data, the aforementioned methods do not assist in capturing spatial information.

As ultrafast science progresses, pulses continue to shorten temporally and broaden spectrally, rendering the temporal and spatial components inseparable. The presence of transmissive and dispersive optical elements along the optical path introduces spatiotemporal coupling effects,28 resulting in the inability to measure the temporal and spatial information of pulses separately for characterization. Consequently, when employing the above techniques for pulse diagnosis, only spatially averaged results are achievable, thus limiting subsequent applications. Given the multidimensional nature of this issue and the inherent difficulties such as spatiotemporal coupling, diagnosing few-cycle pulses in multidimensions remains challenging. Recently, Blöchl et al.29 proposed a near-field imaging method utilizing a nanotip for short-pulse characterization. This technique exploits the local field enhancement effect to sample the electric field with high temporal and spatial resolution. However, it necessitates photocurrent detection devices with a high signal-to-noise ratio. Meanwhile, a method utilizing the multiphoton ionization process in a silicon-based camera for reconstructing the spatiotemporal information of the electric field30 was proposed. However, this approach is only applicable for mid-infrared laser fields with long wavelengths because of the energy band structures of silicon-based materials.

In this study, we propose an all-optical spatiotemporal interferometer for multidimensional characterization of the electric field. Specifically, a signal pulse to be measured is used to modify the temporal phase of an attosecond pulse train from HHG. This phase modulation effectively encodes the complex spatiotemporal structure of the signal pulse, resulting in observable spatiospectral distortions in the HHG spectrum. Through straightforward least squares nonlinear fitting procedures, the spatiotemporal profile of the signal pulse can be retrieved. We have successfully demonstrated the reconstruction of the waveform of a spatiotemporally coupled few-cycle pulse, thereby realizing an all-optical petahertz spatiotemporal oscilloscope.

The principle of the method is shown in Fig. 1. When the driving field interacts with gas atoms, the radiation of HHs will be launched in a time window of around a few hundred attoseconds each half-cycle of the optical field, producing attosecond pulse trains. These attosecond pulses subsequently interfere after radiation. Introducing an additional perturbing field, which carries spatiotemporal coupling effects during the HHG process enables the dipole to introduce an extra spatiotemporal perturbing phase $\sigma (t_e,\tau ,x_f)$.31 This alteration modifies the waveform of the attosecond pulse, as depicted in Fig. 1(b). Due to the time-dependent nature of this perturbing phase, radiations from different time windows experience varying perturbing phases. This will introduce phase differences $\mathrm{\Delta }\varphi (t_e,\tau ,x_f)$ between the attosecond pulses, leading to frequency shifts in each HH after interference. Similarly, the spatial dependence of this perturbing phase results in distinct frequency shifts at different spatial positions of the HHs. Although we discuss the time and spatial dependencies separately, they are fundamentally intertwined and not independent. Finally, by adjusting the delay between the driving field and the perturbing field, the frequency oscillation of the HHs varies across different spatial positions, as depicted in Fig. 1(c).

The HHG mechanism can be described by a three-step model:32 tunnel ionization of bounded electrons, laser acceleration of the freed electron wave packet, and electron recollision with its parent ion. Perturbing phases are primarily introduced during the second step, which can be expressed as $d(t_i,t_e,\tau ,x_f)=\sum _k{d}_0^k(t_i,t_e,x_f)e^{-i{\sigma }^k(t_i,t_e,\tau ,x_f)}$, where the subscript $k$ denotes the corresponding quantum trajectory, $d_0^k$ is the dipole without perturbing fields, $t_i$ is the ionization time of the electrons, $t_e$ is the emission time of the HHs, $\tau$ is the delay between driving field and perturbing field, $x_f$ is the space coordinates in the near field, and $\sigma$ is the perturbing phase. Considering that the ionization and emission times exist in pairs according to the saddle point model,33 the perturbing phase can be denoted as $${\sigma }^k(t_e,\tau ,x_f)=S_1(t_e,\tau ,x_f)-S_0(t_e,x_f),$$(1)where $S_1$ ($S_0$) is the classic action with (without) the perturbing field. Although the number of time windows for HHG is determined by the width of the driving field, for simplicity, we only consider two time slits within one optical cycle. The validity of this model has been confirmed in previous studies.24 When the HHs radiated in two slits are in a condition of constructive interference, the total phase difference is $$\mathrm{\Delta }\varphi (\omega ,\tau ,x_f)={\sigma }_2^k(t_e,\tau ,x_f)-{\sigma }_1^k(t_e,\tau ,x_f)+\omega T_0/2+\pi +\gamma =2m\pi ,$$(2)where $T_0$ is an optical cycle of the driving field, $\omega$ is the frequency of HH, $\gamma$ is the nonadiabatic effect constant,34 and $m$ is an integer. Here, we are considering the frequency oscillation of a specific HH order, and once the order is determined, $t_e$ is also determined. Eventually, we can get the equation for the frequency shift of a specific HH order as (see Section 1 in the Supplementary Material) $$\delta \omega (\tau ,x_f)=\beta (x_f)(\alpha (x_f)E_p(\tau +\mathrm{\Delta }(x_f),x_f)+E_p(\tau ,x_f)),$$(3)where $E_p$ is the perturbing field, $\beta (x_f)$ is the constant, $\alpha (x_f)$ is the intensity ratio of the driving field between two neighboring slits, and $\mathrm{\Delta }(x_f)$ is the time difference between the two attosecond time slits. Equation (3) reveals that the frequency shift of a particular high-order harmonic is proportional to a coherent superposition of two delayed identical waveforms. Figure 2(a) contains all relevant information regarding the perturbing field and can be utilized for spatial and temporal reconstruction of the perturbing field. In the following, we will provide a detailed procedure for waveform reconstruction. Performing a Fourier transform of Eq. (3) with respect to the delay axis yields $$s(\omega ,x_f)={\beta }_1(x_f)s_p(\omega ,x_f)(1+\alpha (x_f)e^{-i\omega \mathrm{\Delta }(x_f)}),$$(4)where $s_p$ is the spatial spectrum of the perturbing field, and ${\beta }_1$ is the constant. Equation (4) indicates that the complex spectrum of the perturbing field is proportional to the experimentally accessible quantity $s$. Given that $|s_p{|}^2$ and $|s{|}^2$ are both experimentally measurable, Eq. (4) can be transformed to $${\left|\frac{s(\omega ,x_f)}{s_p(\omega ,x_f)}\right|}^2={\beta }_2(x_f)(1+\alpha (x_f{)}^2+2\alpha (x_f)\cos (\omega \mathrm{\Delta }(x_f))),$$(5)where ${\beta }_2$ is the constant. Then, we can employ a least squares nonlinear fitting procedure to fit the value of $|\frac{s(\omega )}{s_p(\omega )}{|}^2$ and subsequently determine ${\beta }_2$, $\alpha$, and $\mathrm{\Delta }$ using Eq. (5). With ${\beta }_2$, $\alpha$, and $\mathrm{\Delta }$ known, the complex spectrum of the perturbing field can be calculated via Eq. (4).

To verify the validity of the aforementioned method, we employ the strong field approximation model35 and the thin-slab model36 for simulation. The laser parameters adopted in the simulation are similar to those of the subsequent experiments: the driving field center wavelength is 780 nm, the pulse width is 8 fs, and the spot size in the interaction region is $100\mu \mathrm{m}$, with a peak intensity of $2\times {10}^{14}{\mathrm{W}/\mathrm{cm}}^2$, whereas the peak intensity of the perturbing field is $1\times {10}^{11}{\mathrm{W}/\mathrm{cm}}^2$, carrying a spatial chirp effect, and other parameters of the perturbing field are the same as the driving field. Argon serves as the interacting gas, and two short trajectories are selected for HHG (see Section 2 in the Supplementary Material) while ignoring the propagation effect. Indeed, in our subsequent experiments, we utilize a pulsed nozzle with an aperture size of $200\mu \mathrm{m}$. It injects a very thin layer of gas, with the pressure at the nozzle exit only around a few torr, enabling the neglect of propagation effects. Subsequently, we can use the simulation results to reconstruct the electric field.

Here, we outline the methodology for the multidimensional reconstruction of the spatiotemporal coupling field. First, we extract the center of mass of the 33rd-order harmonic spectrum in the near field as a function of delay and space, as illustrated in Fig. 2(a). When the perturbing field exhibits the spatial chirp effect, the frequency of the oscillations varies across space, presenting a distinctive spatiotemporal structure. Second, a time Fourier transform of Fig. 2(a) yields the spectrum of the frequency oscillation across space. Then, the parameters ${\beta }_2$, $\alpha$, and $\mathrm{\Delta }$ can be obtained via the fitting procedure. For convenience, we present the reconstruction results of the spectrum $s_p$ for two different spatial positions [the red and black dashed lines depicted in Fig. 2(a)], both of which closely match the original spectrum; see Figs. 2(b) and 2(c). Finally, as the reconstructed spectrum of the perturbing field $s_p$ contains all relevant information regarding the perturbing field, its inverse Fourier transform provides both spatial and temporal details of the perturbing field $E_p$, thus achieving the reconstruction of the waveform of the pulses; see Figs. 2(d) and 2(e). Similarly, we present the reconstruction results of the electric field for two different spatial positions to elucidate the details, and these results agree with the original electric field.

To establish the effectiveness and versatility of our approach beyond the electric fields with spatial chirp effect, we also reconstruct the waveform of the perturbing field carrying the pulse front tilt effect, the second-order spatiotemporal coupling effect, different order dispersions, and spatial and temporal chirp effect. The resulting reconstruction, depicted in Section 3 in the Supplementary Material, highlights its consistent performance. These complementary sets of simulation data validate the method’s reliability and universality for the reconstruction of spatiotemporal electric field information.

We conducted experiments to reconstruct spatiotemporally coupled electric fields using the proposed approach; the experimental layout is shown in Fig. 3(a). A sub-ten femtosecond few-cycle near-infrared (NIR) laser pulse, which is generated using a neon-filled hollow-core fiber combined with a set of chirped mirrors, is sent into a Mach–Zehnder interferometer consisting of two 20% reflection beam splitters. The transmitted beam undergoes polarization gating37 optics to reduce the number of attosecond pulses during HHG. An neutral density filter and an aperture are inserted in the reflected beam path to control the intensity of the perturbing field. The gas jet for HHG is filled with argon gas with a backing pressure of 20 bar and is positioned a few millimeters upstream of the laser focus to phase-match the short trajectory.38 The effective peak intensity of the driving field at the focus is estimated to be $3.5\times {10}^{14}\mathrm{W}/{\mathrm{cm}}^2$, whereas the peak intensity of the perturbing field is $\sim 9\times {10}^{11}\mathrm{W}/{\mathrm{cm}}^2$. The delay step size is 400. To create a spatiotemporal coupling field, the perturbing field passes through a pair of misaligned fused silica wedges (with a misalignment angle of 30 deg and an apex angle of 8 deg), inducing a spatial chirp effect at the laser–gas interaction region. This spatial chirp effect, a typical first-order spatiotemporal coupling phenomenon, has been widely studied in previous research.39^,40 We exploit this well-established effect to validate the reliability of our method. Meanwhile, a liftable reflector positioned in front of the nozzle redirects the perturbing field outside the cavity. Utilizing an infrared spectrometer enables the measurement of the spatial spectrum of the perturbing field, revealing a distinctive spatial chirp effect, as depicted in Fig. 3(b). To ensure that the measured spatial spectrum accurately represents the spectrum of the region where the HHG occurs, we have mirrored the optical path of the infrared spectrometer. An objective lens situated along the optical path is used to magnify the focal spot size. The reflected field propagates through both the objective lens and the slit before being detected by the spectrometer. Precise measurement of the spatial spectrum in the near field is achieved by moving the spectrometer via a motorized stage.

Equation (5) implies that the spatial spectrum $|s{|}^2$ and the spatial spectrum $|s_p{|}^2$ must be directly correlated in the near field before utilizing nonlinear fitting of $|\frac{s(\omega )}{s_p(\omega )}{|}^2$. Due to the divergence angle of the HHs, there exists a difference in the spatial width of the HHs in the far field compared with that in the near field. Therefore, we calculate the spatial width of the HHs from the short trajectory in the near field with the equation ${\theta }_s=\frac{{\lambda }_h}{\pi w_h}\sqrt{1+4{\alpha }_s^2{I}_0^2\frac{w_h^4}{w_f^4}}$,41 where ${\theta }_s$ is the divergence angle of the HHs from the short trajectories, ${\lambda }_h$ is the wavelength of the HHs, $w_h$ is the beam waist of the HHs, ${\alpha }_s$ is the coefficient of the short trajectories and is determined by the intensity-dependent dipole phase, $I_0$ is the intensity of the driving field, and $w_f$ is the beam waist of the driving field. In our calculations, we consider both the driving field and the HHs as Gaussian beams, and the HHs are dominated by short trajectories due to their small divergence angle. In addition, the optical path of the infrared spectrometer allows for the detection of the beam waist of the driving field $w_f$, which is $85\mu \mathrm{m}$ in our experiment. By substituting the parameters ${\alpha }_s=-3.8\times {10}^{-14}{\mathrm{W}}^{-1}{\mathrm{cm}}^2$,41${\theta }_s=1.3\mathrm{mrad}$, and $w_f=85\mu \mathrm{m}$, the spatial width of the 25th-order harmonic in the near field is $\sim 42\mu \mathrm{m}$. This procedure establishes a correspondence between spatial points in the far and near field. It should be noted that although we detect the HHs in the far field in our experiments, our simulation findings show that the delay and spatial dependent frequency shift of a specific harmonic order demonstrate a rather similar pattern for both near and far fields (see Section 4 in the Supplementary Material), and therefore, it is legitimate to utilize the detected HH spectrum for pulse reconstruction.

The center of mass of the 25th-order harmonic spectrum as a function of delay and space is then shown in Fig. 4(a). A comparison between Figs. 4(a) and 2(a) illustrates their similar structures, particularly evident in a notable section of small peaks near zero-delay and the spatially inhomogeneous distribution of oscillation frequencies. To present the results of the reconstruction simply and comprehensively, six spatial points are uniformly chosen for reconstruction. Following the identical processing steps mentioned above, we reconstruct the spectrum of different spatial positions of the perturbing field in the interaction region, as illustrated in Fig. 4(b). The reconstructed spectrum shows overall satisfactory agreement with the measured spectrum. The measured and reconstructed spatial spectra are shown in Fig. 4(c). Figure 4(d) shows the reconstructed electric field of the perturbing pulse. Both Figs. 4(c) and 4(d) distinctly demonstrate a profound redshift from top to bottom. The reconstructed electric field demonstrates an obvious spatial chirp effect, consistent with findings in previous research.39^,40 As the perturbing field passes through an extra pair of wedges, its pulse width becomes longer than the driving field. Although FROG cannot access a spatiotemporal waveform, it does provide valuable information such as group delay in a space-averaged manner. We also performed FROG measurement; the results are shown in Section 5 in the Supplementary Material. By comparing the spectral phase reconstructed using space-averaged data in our method with that obtained directly from FROG, we find them to be in good agreement, thereby validating our measurements to some extent.

Notably, we utilize a laser system with unstable CEP, resulting in reconstructed electric fields with averaged CEP. The averaged CEP waveform contains the entire carrier information of the electric field, including high-order dispersion effects. This is similar to an attosecond streaking experiment using a CEP-unstabilized laser.42 Therefore, although the current method is applicable to measure the waveform of a CEP-fixed laser pulse in principle, when a CEP-unstabilized laser is used, the current method can still precisely measure the carrier of the laser field with the envelope information averaged. Due to the CEP being unlocked, we cannot measure it experimentally. Instead, we performed simulations under various CEP conditions (see Section 6 in the Supplementary Material) to demonstrate the sensitivity of our method to CEP. The simulations indicate that our method can accurately retrieve the spatially dependent spectral phase. In the current experiment, the diagnostic pulse frequency reaches up to 0.5 PW, indicating its capability of responding to petahertz electromagnetic fields. Furthermore, although our discussion has thus focused solely on reconstructing linearly polarized electric fields, it can indeed be extended to measure the electric field of a few-cycle pulse with arbitrary polarization as well. As mentioned in a previous study,23 by decomposing the electric field into two orthogonal polarization components and measuring them separately, the spatiotemporal waveform of the vectorial electric field is accessible. In our experimental configuration, the beam size of the HH is enlarged as it propagates to the charge-coupled device camera, and the magnification factor is $\sim 15$ according to the divergence formula ${\theta }_s=\frac{{\lambda }_h}{\pi w_h}\sqrt{1+4{\alpha }_s^2{I}_0^2\frac{w_h^4}{w_f^4}}$. Given that the pixel size of the camera is $20\mu \mathrm{m}$, the spatial resolution of our measurement is $1.3\mu \mathrm{m}$. The temporal information is retrieved via a Fourier analysis of the delay-dependent frequency shift of HH, and the time resolution of our method is then determined by the delay step size, which is 0.4 fs in our experiment. It should be noted that this method utilizes gratings to spectrally resolve the HHs, which already occupy one dimension of the spatial axis. As a result, currently, we can only diagnose one-dimensional spatial information using this approach. In principle, our method can be extended to 2D spatial measurement if the complete spatial spectrum of the far-field HHs can be measured. This will involve technical modification of the current extreme ultraviolet spectrometer and will be subject to future study.

Based on the principle, we consider there are two main influences that can cause the method to break down: intensity and spectral bandwidth. First, when the intensity of the perturbing field is too strong, it will not be possible to analyze it using HH perturbation theory. Second, our approach is based on a few-slit interference model, which will fail if the bandwidth is too large and the pulse duration is too short to span at least two temporal attosecond slits. We have also carried out simulations to test the limitations of our method. The simulation results (see Section 7 in the Supplementary Material) show that our method is still valid for a perturbing field with an intensity ratio of up to 1% and a spectral bandwidth exceeding 500 nm. Furthermore, this method does not require that the wavelength of the perturbing field coincides with that of the driving pulse, that is, simulations (see Section 8 in the Supplementary Material) show that noninteger harmonic or second-harmonic of the fundamental one can still be reconstructed with good accuracy, which suggests the wide applicability of our method.

In summary, we present an all-optical oscilloscope that allows for capturing the waveform of a few-cycle laser pulse in both spatial and temporal domains. The operational principle of this oscilloscope is rooted in the concept that the HHG process can be treated as a time-domain interferometer with few slits. Introducing a perturbing field imposes additional phase information into the HHs, leading to an enriched structure of the HH spectrum in the spatiospectral domain. Upon scanning the delay, the detailed information of the electric field is imprinted into the HH spectrogram, from which the spatiotemporal structure of the perturbing field can be reconstructed. Numerical simulations based on the strong field approximation verify the reliability of this method. Due to the high temporal resolution and fast frequency response rate of the HHG process, this method is applicable for few-cycle pulse measurements, with a frequency response extending up to petahertz. An experiment is performed to successfully retrieve the waveform of a spatially chirped NIR few-cycle laser pulse. This approach provides simple and reliable metrology for the characterization of the waveform of the few-cycle pulses multidimensionally. Owing to its spatiotemporal resolution, the method also has the potential for important applications in probing ultrafast dynamical processes that carry spatiotemporal information.

Qi Zeng received his BS degree from the School of Physics and Electronics, Central South University, China, in 2019. He is currently a PhD student at Huazhong University of Science and Technology (HUST), Wuhan, China. His research interests focus on spatiotemporal measurement and control of the high-order harmonics.

Xinyue Yang received his BS degree from the School of Science, Northwest A&F University, China, in 2023. He is currently a PhD student at Huazhong University of Science and Technology, Wuhan, China. His research interests focus on the control and measurement based on high-order harmonics.

Yimin Deng received his BS degree from Huazhong University of Science and Technology, Wuhan, China, in 2022. He is currently a PhD student at Huazhong University of Science and Technology, Wuhan, China. His research interests focus on element-specific, time-resolved, and transition-channel-resolved measurements of electron dynamics at water window.

Wei Cao is a professor at the School of Physics, HUST, Wuhan, China. He received his PhD in physics from Kansas State University in 2014. From 2014 to 2017, he worked as a postdoctoral fellow at Lawrence Berkeley National Lab in the United States. His research focuses on the generation of the attosecond pulse and its application for ultrafast process diagnosis.

Peixiang Lu is a professor and the vice director at Wuhan National Laboratory for Optoelectronics, HUST, Wuhan, China. He received his PhD from Shanghai Institute of Optics and Fine Mechanics, Chinese Academy of Sciences. He was selected as a Cheung Kong Scholar Chair Professor and the distinguished young scholar of NSFC, and he was selected as OSA fellow in 2016. His research interest is strong ultrafast optics.