High-order harmonic generation spectroscopy by recolliding electron caustics

Spectral focusing of the recolliding electron in high-order harmonic generation driven by two-color fields is shown to be a powerful tool for isolating and enhancing hidden spectral features of the target under study. In previous works we used this technique for probing multi-electron effects in xenon and we compared our experimental results with time-dependent configuration-interaction singles calculations. We demonstrate here that this technique can be exploited for reconstructing the enhancement factor of the xenon giant dipole resonance and we discuss the sensitivity of this method to macroscopic effects. We then extend the technique to argon in order to test the applicability of this procedure to other targets.


Introduction
Since its discovery, high-order harmonic generation (HHG) has paved the way towards time-resolving ultrafast electronic processes on the attosecond time scale, in particular thanks to its capacity to generate isolated attosecond pulses in the extreme ultraviolet (EUV) energy range and beyond [1][2][3]. The HHG process can be used, however, not only as a source of coherent EUV light, but also as a powerful spectroscopy tool for studying the structure and dynamics of the generating medium, a technique named high-harmonic spectroscopy (HHS). The process is initiated by a strong laser pulse that interacts with the atoms of the target leading to the tunnelionization of an electron wave-packet each half-cycle of the driving field. The electron is then accelerated by the laser field and can recombine with the parent ion resulting in the emission of an EUV field that encodes the information about the target with angstrom spatial resolution and attosecond temporal resolution [4,5]. In this respect, HHS can be viewed as a self-probing spectroscopy, where the target is probed by one of its own electrons [6][7][8].
One of the main challenges of HHS is performing quantitative analysis on the target under study. A quantummechanical treatment of the process based on the strong field approximation (SFA) allows the factorization of the HHG spectrum E(Ω) as the product of a recombining electron wave-packet w(Ω) and a photo-recombination transition dipole moment d PW (Ω) from plane-wave continuum states to the bound state [9]. The dipole d PW (Ω) can then be substituted with exact transition dipole moments d(Ω) which use scattering waves instead of plane waves [10,11], a procedure named quantitative re-scattering theory (QRS). In this framework, in order to collect reliable and significant informations on the target, it is essential to know and control the spectral and spatial properties of the recombining electron wave-packet. This can be calibrated by using another atom with similar ionization potential as a reference, where d(Ω) is assumed to be known [12].
This procedure, in practice, works only in very specific experimental conditions. As a matter of fact, the macroscopic high-harmonic spectrum comes from the coherent build-up of many single-atom emissions and its final shape can substantially differ from the single-atom prediction [13]. A retrieval of the macroscopic response of the target often requires the full knowledge of the sample and pulse properties, in addition to time-consuming calculations. As a consequence, QRS can only be efficiently applied at low pressures and laser intensities, in order to avoid saturation and plasma defocusing [14]. Moreover, QRS cannot reproduce the influence of the ionic state on the propagation of the electron as well as the evolution of the ionic state during this propagation. In order to capture these physical effects one needs to abandon factorization and rely on ab initio procedures. In this case, the picture is complicated by the fact that the spectrum is given by the sum of two main contributions, namely short and long trajectories, which usually have different phase-matching properties. Whilst single-atom predictions show a predominance of the long trajectory contribution, due to the high ionization ratio associated, real experiments show a higher short trajectory contribution, due to their lower divergence and more likely phase-matching conditions [15,16]. This makes it difficult to directly compare experiments with singleatom ab initio calculations, where contributions of different trajectories cannot be easily extracted.
Therefore it becomes necessary to reduce the unknown variables governing the process and/or to get more information from the sample. This can be done by accurately manipulating the generating field, allowing to control one or both of the ionization and recombination processes. This approach, for example, is fundamental if one wants to perform pump-probe-like experiments for studying electron-hole dynamics triggered by tunneling [17,18]. Recently, we have shown that by adding a second harmonic to the fundamental field, it is possible to isolate and enhance spectral features that are normally not accessible in HHG experiments [19]. By varying the phase j and the amplitude ratio α between the two fields, we were able to control the properties of spectral caustics in the cutoff region of the spectrum. Effects within the atomic system that are associated with a specific branch of the caustic can then be probed over a wide spectral range.
We applied this technique to probe the giant dipole resonance (GDR) of xenon, a broad enhancement of the dipole around 100 eV which is due to multi-electronic dynamics involving inner orbitals during the recombination step [12,20]. In this paper we investigate the conditions in which this approach can be considered reliable. We further extend the analysis and demonstrate a good quantitative agreement between our experiment and the time-dependent configuration-interaction singles (TDCIS) calculations. TDCIS is a powerful theoretical framework, since it allows switching on and off the effect of the resonance. We are able in this way to extract the enhancement factor of the GDR. We then apply this technique to another target, argon, in order to show the versatility of the method that can be applied even to systems that do not show multi-electronic resonances.

Method
High-order harmonics (HH) are generated from the interaction of a strong field with a gas target of atoms or molecules. In the simplest case, the driving field t E 1 ( ) is represented by a monochromatic wave of wavelength λ 1 (frequency ω 1 ) with linear polarization in the u ẑ direction: with E 10 being the amplitude of the wave and I E 1 10 2 µ its intensity. In this case, HH comes from the electron-ion recombination of two electron trajectories, namely the short and long trajectories, generated each half-cycle of the laser field. Figure 1(a) represents a simplified picture of the process. Long and short trajectories driven by a laser field of intensity I 1 =6.58×10 13 W cm -2 and associated to the emission of one particular harmonic (HH57) in xenon are reported with black and green arrows respectively. They have been calculated by using the semi-classical three-step model [4]. It can be observed that each half-cycle of the field gives rise to equivalent pairs of short and long trajectories. In this case, the spectrum extends up to the cutoff energy ÿΩ=1.32I p +K, with I p being the ionization potential of the target and K the maximum kinetic energy of the recombining electron, equal to 3.17U p , with U I p 1 2 1 l µ the ponderomotive energy. Around the cutoff, long and short trajectories coalesce into a single cutoff trajectory.
When adding a second field of wavelength λ 2 =λ 1 /2 (frequency ω 2 =2ω 1 ) with the same polarization of the fundamental field, the laser field reads: In this situation, the spectrum is built by the sum of two main contributions, coming from electron trajectories launched during different half cycles of the total electric field. Since each contribution is given by the sum of long and short trajectories, HH are generally produced by the interference of four different classes of trajectories [21,22]. Figure 1(b) represents a simplified picture of this process. However, spectral focusing of the recombining electron wave-packet in correspondence of the cutoff of each contribution allows to selectively enhance and probe regions of the spectra where a single recombination event occurs. Indeed, in correspondence of the two cutoffs associated to each contribution, caustics appear [19,23]. For a given (α, j) combination, the position of the two cutoffs ÿΩ + and ÿΩ − , namely upper and lower semi-classical cutoff, can be calculated using the semi-classical equation ÿΩ ± =I p +K ± , with K + and K − being the maximum kinetic energies reached by the electron recombining from the two different sides of the ion. They give an indication of the two points of the spectrum where short and long trajectories coalesce. The position of the two caustics c W + and c W -associated to the two semi-classical cutoffs ÿΩ + and ÿΩ − , namely upper and lower caustic cutoffs, can be calculated by exploiting the analytic description of HHG by atoms in a two-color field proposed by Frolov et al [24]. Within this framework, the recombining electron wave-packet c w ( ) at the emitted photon frequency ω is represented by the interference of two terms, each one associated to one of the two cutoffs Ω + and Ω − : where t 0,± and t ± denote the ionization and the recombination time, respectively, associated to the semi-classical cutoff trajectory of energy Ω ± , and χ is a constant with respect to the photon energy ω that depends only on the cutoff trajec- is the Airy function of argument: [24]. As for the one color case, the real cutoff is higher than the semi-classical one, which is understood to be due to tunneling [9]. The pioneering definition of the cutoff given by Lewenstein et al [9] corresponds, in the analytical description of the process proposed by Frolov et al [25], to the position where the argument of the Airy function representing the wave-packet goes to zero. Extending the procedure to the two-color case [24] allows to analytically find two new cutoffs L W + and L W -, which are referred to as Lewenstein upper and lower cutoffs: where the dependence of t ± , t 0,± , ÿΩ ± and L W  with α and j has now been made explicit. The two caustic cutoffs c W  correspond instead to the maxima of these two Airy functions, which directly follows from the original definition of caustic given by G B Airy [26]. The solution leads to two caustic cutoff trajectories that depend on the (α, j) combination. By tuning the phase delay j between the two colors, the position of the two caustics can be tuned over a wide range of energies, as shown in figure 1(c) for different values of the ratio α between the two field amplitudes. Such control can be achieved by a pair of wedges, as discussed in detail in [19,27]. In order to highlight the advantages, the limits and the reliability of this method, we will now address the details of how the selection of unique recombination event can be done. We simulated HHG in krypton with the Lewenstein model considering a two-color field with α=0.2. The yield Y(q, j) of each qth harmonic is reported in figure 2  (a) Semi-classical electron displacement of short (green arrows) and long (black arrows) trajectories associated to HH57 generated by the driving field force F 1 (t)=−eE 1 (t); the force F 1 (t), multiplied by 10 3 , is represented by the blue line. (b) Semi-classical electron displacement of short (green arrows) and long (black arrows) trajectories associated to HH57 generated by the driving field force  [24], caustic cutoffs c W  provide a more accurate estimate of the cutoffs. We add to this theoretical argument that c W  can be experimentally measured, since they correspond to a maximum of the spectrum.
In particular, for jä(0.4, 1)π, a strong caustic is predicted in the region of energies between the semi-classical and the Lewenstein lower cutoffs (lower branch). The sketch depicted in the inset 2(c) shows that in this region the cutoff trajectory is enhanced, while trajectories recombining from the opposite direction (hereafter defined opposite trajectories) are suppressed due to lower ionization probabilities [21,22,24]. The selection of one trajectory for each cycle is clear from the comparable intensity of even and odd harmonics in this region. It must be noticed that in the region A indicated by the arrow, spectral focusing is reduced and destructive interference with opposite trajectories produces the observed minimum of harmonic intensity. Also in the region of energies between the semiclassical and Lewenstein upper cutoffs (upper branch) a single recombination event corresponding to the upper cutoff trajectory is selected. Indeed, within the semi-classical framework, opposite trajectories are not present, since their cutoff is lower. This is depicted in the sketch of inset 2(d). Also in this case, the selection of one trajectory for each cycle is clear from the comparable intensity of even and odd harmonics (see figure 2(b)).
In order to assess the reliability of our probe for performing quantitative analysis on the dipole response, we simulated HH from many atoms distributed in a plane perpendicular to the propagation axis. Due to the cylindrical symmetry of the driving field with respect to the propagation axis, the resulting HHG near-field yield Y q r E q r , , , , )| depends only on the radial distance r. Its value at r=0 corresponds to the single-atom yield Y(q, j)=Y n (q, j, 0) reported in figure 2(a). We then propagate E HHG (q, j, r) within the paraxial approximation. The Fresnel diffraction formula [28] leads to the following on-axis far-field field E q, f HHG j ( ) at a distance d: By comparing the far-field yield with the singleatom yield Y(q, j) we are able to evaluate how much the spectrum is reshaped due to the different divergence of trajectories. The plateau changes dramatically, compared to singleatom plateau harmonics. For jä(0.4π, π) the amplitude of the oscillations observed below the lower cutoff reduces, which is attributed to the fact that long trajectories are filtered out due to their higher divergence. For jä(0, 0.4π) the interference pattern observed below the lower cutoff changes with respect to the single-atom one. Before propagation, the pattern comes from the interference of all the four trajectories recombining from opposite sides. After propagation, only the interference between short trajectories survives, which gives rise to the well- known pattern of interference that has been extensively discussed in [29][30][31][32]. Figure 3(b) reports the ratio between far-field and singleatom HH. One interesting result is that in the region of the caustics, the yield is not affected by propagation as it is in the plateau energy region. In this region, maximum and minimum ratios in logarithmic scale are respectively 0.4 and −0.1 and are confined in particular in the region of the lower cutoff around j≈0.2π (region A -black arrow) and in the region of the upper cutoff around j≈0.8π (region B -black arrow). As discussed previously, the interference observed in region A is due to a comparable contribution coming from opposite trajectories. As for region B, here the upper cutoff yield decreases dramatically and becomes comparable to the photon yield of opposite trajectories generated beyond the lower cutoff. In these two regions, after propagation, the relative intensity of interfering trajectories changes, leading to a different modulation amplitude. Apart from these two regions, far-field HH at the caustic resemble the single-atom predictions, which suggests that direct comparison of the cutoff yields with ab initio single-atom predictions is more reliable than comparing HH at the plateau. This observation is of crucial importance. Indeed, one of the main challenges in the original experiment was keeping the gas pressure, and therefore the signal level and contrast, very low in order to avoid propagation effects. The results shown here suggest that the sensitivity of our measurements can be improved, even if more theoretical and experimental investigations should be performed in order to assess up to which extent this can be done. We want to stress that in these simulations only the divergence of the HH is considered. A full propagation taking into account phase-matching in the target should be performed in order to extend these concepts to more general cases.
However, it is clear that the advantage of probing the target at the caustic is the possibility to select only one trajectory, overcoming the problem of determining which trajectory is actually contributing to the process. Moreover, this trajectory selection occurs around relatively high energies, where most of the interesting physics is occurring.

Xenon case: the giant resonance
We applied this technique to probe the giant resonance of xenon. This multi-electron process in HHG can be described within the three-step semi-classical model. During the ionization step, a hole in the external 5p orbital is created due to the tunnel-ionization process. The ionized electron, after being accelerated by the driving field, recombines with the parent ion following two possible channels: it recombines directly with the 5p hole or it exchanges energy with the electrons in the inner 4d orbital, promoting one of these electrons to the 5p shell and recombining with the 4d hole left behind. This results in a broad enhancement of the spectrum starting around 60 eV and having its peak around 100 eV [12,20]. Figure 4 reports the experimental and theoretical results discussed more in details in [19]. In brief, figure 4(a) shows the experimental HHG yield Y exp (ω, j) generated in xenon by a two-color field with λ 1 =1550 nm and α=0.4 as a function of j and the photon energy ÿω. Figure 4   (TDCIS-all). In this last case, the upper branch corresponding to the upper caustic is enhanced, due to the giant resonance phenomenon.
The value of α, the intensity of the fundamental field I 1 =6.58×10 13 W cm -2 , and the absolute value of the phase j were all calibrated from the experiment by comparing the position of the two cutoffs for each j with the observed ones. In our previous work [19], TDCIS and experimental caustic cutoffs were compared with the semi-classical cutoffs calculated from the three-step model. Here, the upper and lower caustic cutoffs For each value of the phase j, we integrate Y exp (ω, j) between 1.25 )and 1.25 ) , in order to filter out the modulation due to the harmonic comb (see figure 5(a)). The results for the two branches named V and C in figure 5(a) ) and (c) respectively. TDCIS simulations including 5p, 5s and 4d orbitals as well as the experimental results have been normalized to the lower energy caustic yield (yield at j≈0.7π and ÿΩ − ≈47.6 eV). This choice is motivated by the fact that, around this point, the presence of a strong caustic and the absence of strong multi-electronic effects suggest agreement with both the TDCIS models. TDCIS simulations with only 5p have been scaled accordingly. Experimental yields below 0.015 are not shown due to their low signal-to-noise ratio.
It can be observed that along the C branch, the experimental yield Y C exp w ( ) agrees very well with the TDCIS calculations only when all the 5p, 5s and 4d are taken into account. This proves that the giant resonance enhances the contributions of trajectories that are suppressed in many other systems, as concluded in [19]. There is a disagreement around 63eV that was already discussed in [19]. Indeed, this region corresponds to region A underlined in figures 2(a) and 3(b).  As discussed previously, single-atom calculations predict here the presence of interferences between the lower cutoff trajectory recombining from one side and the short and long trajectories recombining from the opposite side. This interference is blurred out in the experiment likely due to spatial and temporal averaging effects [33].
As for the V branch, the experimental yield Y V exp w ( ) shows a better agreement with the TDCIS calculations including all the orbitals Y V TDCIS ,all w ( ) with respect to the ones including only the 5p orbital Y V TDCIS ,5p w ( ). However we do not have a satisfactory quantitative agreement. This is probably due to the fact that the caustic corresponding to the V branch is not very pronounced. As a matter of fact, the experimental caustic appears to be broader than the one predicted by theory. This can be attributed again to averaging effects. If the caustic is not pronounced, we expect to be more sensitive to this broadening, that worsens the agreement. This does not happen in correspondence of the C branch, where the caustic is clearly visible.
Since the C branch shows optimal agreement with theory, we can extract useful information from data. Figure 6 reports the ratio of the experimental yield with respect to TDCIS calculations including only the 5p orbital (black line) compared with the ratio of the two TDCIS calculations with and without the 4d and 5s orbital (blue line). This corresponds to the enhancement factor introduced by the giant resonance. Our experiment successfully predicts an enhancement starting at 60 eV and growing up to a factor of 10 at 75 eV. Gray line in figure 6 is the ratio between the experiment and TDCIS including 5p, 5s and 4d orbitals, which is around one apart from the dotted line region where, as said previously, we expect interference of opposite trajectories. The good agreement confirms that controllable probes corresponding to a single recombination event can be used to quantitatively probe multi-electron dynamics.
Within the QRS framework, the HHG yield is described as the product of the recombining electron wave-function and the partial photo-ionization cross section (PICS) Σ(ω) of the target: The recombining electron wave-function c w ( ) is assumed to be not affected by correlations which only modify the PICS value. Due to this factorization, the enhancement factor is given by the ratio of the partial PICS taking into account correlations between 4d, 5s and 5p orbitals (Σ all ) and the one from the 5p orbital alone (Σ 5p ). The result of this calculation obtained within the TDCIS framework by Pabst et al [20] is reported as a red line in figure 6. It can be observed that QRS predicts a smoother enhancement factor that does not fit our experimental observations. Moreover, the steep increase of the yield beyond 60 eV is not reproduced by the QRS model, which means that it can be attributed mostly to the recombining wave-function rather than to a change of the recombination dipole. This is a clear evidence that the recombining electron wave-packet c w ( ) is affected by the correlations of the ionic state, due to the mutual energy exchange between the ion and the wave-packet. The same steep variation of the wave-packet beyond 60 eV has been theoretically predicted in [20]. Our experiment confirms this prediction.
The analysis can be refined by taking into account empirically the observed averaging effects that seem to be responsible of caustic broadening and interference blurring. This analysis is not intended to demonstrate quantitative agreement, but it will confirm that averaging effects are likely the primary reason for the disagreements along the V branch and in region A of the C branch. In order to simulate averaging effects, we applied a 2D-Gaussian filtering to the TDCIS results. Single-atom TDCIS maps Y TDCIS (ω, j) have been convoluted with a 2D-Gaussian function G(ω, j) for reproducing the observed macroscopic broadening. The resulting filtered yield Y TDCIS reads: where σ ω =4 eV and σ j =0.05π. After filtering we performed the same analysis described previously for the two branches V and C. In this way, the new yields Y V TDCIS 5p, It can be observed a better agreement of the V branch experimental yield with the one predicted by TDCIS including all 5p, 5s and 4d orbitals Y V TDCIS all, w ( ). A better agreement through the C branch with Y C TDCIS all, w ( ) is also observed around 63 eV, where the interference predicted by the simulation has been blurred out by the filtering. This confirms that temporal Red line: Σ all (ω)/Σ 5p (ω) (data retrieved from partial PICS shown in [20]). Dotted lines underline the region where we expect a worse agreement. and spatial averaging effects are likely the main reason for the blurring of interferences and the broadening of the caustics in the experiment. Despite filtering, the steep variation of the enhancement factor is still observable both in the theory and in the values retrieved from the experiment, as can be seen from figure 7(c).

Argon case
By applying this method to xenon we were able to tune the caustic in a wide range of energies, from 45 up to 75 eV. This was possible due to the presence of the giant resonance, that enhances the contribution of the upper branch making it clear in the HH spectrum. One could question the applicability of the method to other targets that do not show resonances. For instance, the upper branch in SFA simulations performed in krypton is barely visible, as can be observed from figure 2. In order to assess the versatility of the method, we performed the experiment in argon. In this case, the intensity of the fundamental beam and the parameter α are respectively I 1 =1.15×10 14 W cm -2 and α=0.2. These parameters have been calibrated by a comparison of the experimental data with the expected position of the cutoffs. In this case, backing pressure of the valve was increased to 4 bar in order to evaluate the accuracy of the method even at higher pressures and in order to increase signal. The HH spectra observed are reported in figure 8(a) as a function of j.
Multi-electronic effects are not expected to be strong in the region of energies spanned by the caustic [34]. TDCIS simulations have been performed, by including only the 3p z external orbital. TDCIS simulations taking into account the correlations of all the electrons in the 3p orbital and the ones in the inner 3s orbital have shown to give a negligible contribution in the region considered (not shown). Results for the TDCIS-3p z case are shown in figure 8(b). We applied the same analysis we performed for xenon to the argon results, taking the harmonic yield along the C and V branches underlined in cyan and violet in figure 8. The C branch results reported in figure 8(d) show an optimal agreement between the experimental yield (solid line) and the theoretical one (dotted line) with the exception of the region of energies between 96 and 103 eV, where interference between opposite trajectories is expected. Also in the argon case, the interference is blurred out when theoretical simulations are convoluted with the 2D-Gaussian filter. As for the V branch, a good agreement is observed, both before and after filtering, on the lower caustic. This is probably due to the fact that for lower α values the HHG response in this region is almost flat and the caustic can be detected along all the lower cutoff. Along the upper cutoff, a better agreement is obtained with the Gaussian-filtered theoretical results.
In contrast to xenon, it is interesting to notice that, even if no resonances are present in argon, photons can be detected along all the upper branch. This is due to two reasons: the higher backing pressure used and the high dipole response of argon. Despite a higher pressure is used, and although the broadening of the resonance appears to be higher than the one observed in xenon, a good agreement is observed, even before filtering, in the region of energies where a single recombination event is selected. This confirms that the method can be extended to common targets, even under experimental conditions that usually do not allow direct comparison with ab initio simulations.

Conclusions
We presented a HHS method for selecting single recombination events during HHG in two-color fields. By wisely selecting regions of the spectra corresponding to caustics of the recombining electron wave-packet, we are able to select single recombination events, which allowed us to disentangle the role of the wave-packet enhancement during a multi-electronic resonance. The results showed that ab initio theories are indeed necessary for correctly reproducing the enhancement factor of the GDR in xenon.
Results in argon show that the method can be applied also to targets that do not show resonances, provided that enough XUV photons are generated and collected along the energy region spanned by the caustics. A good agreement is achieved even under experimental conditions for which direct comparison with ab initio theories is usually not reliable. The method constitutes a step forward for quantitative analysis of dipole matrix elements by HHS.
The research leading to these results has received funding from the European