Space-time coupling of the carrier-envelope phase in ultrafast optical pulses

The carrier-envelope phase (CEP) plays an increasingly important role in precise frequency comb spectroscopy, all-optical atomic clocks, quantum science and technology, astronomy, space-borne-metrology, and strong-field science. Hitherto, it has been essentially assumed that CEP is strictly a temporal phenomenon. Here we introduce an approach for space-time calculation of the CEP in the spatially defined region of interest. We find a significant variation of CEP in the focal volume of refracting focusing elements and accurately calculate its value. We discuss the implications and importance of this finding. Our method is particularly suitable for application to complex, real-world, optical systems thereby making it especially useful to applications in research labs as well as in the engineering of innovative designs that rely on the CEP.


INTRODUCTION
Carrier-envelope phase (CEP), φ ce , which is defined as a phase offset between the rapidly oscillating electric field within ultrashort pulse and the peak of the pulse envelope [1], became ubiquitous in many branches of contemporary physics. For instance, in the frequency comb metrology [2], where the frequency comb offset (f 0 ) is linearly proportional to the pulseto-pulse CEP shift [1], this CEP shift, along with pulse repetition rate (f r ) would define a precise frequency ruler -the frequency comb f n = f 0 + nf r [2]. Such frequency combs, capable of spanning very broad spectral range all the way to the extreme ultraviolet [3], already have had transformative impact on precise spectroscopy [4,5], all-optical atomic clocks [2], quantum science and technology [6,7], and even space-borne-metrology [8] and astronomy [9]. CEP also plays a central role in the strong field science [10,11]. Hence, CEP metrology and control are of significant importance and a lot of effort has been dedicated to these goals [12][13][14][15].
Hitherto, in most of these works, albeit with few exceptions [16][17][18], it has been tacitly or explicitly assumed that CEP, φ ce , is strictly a temporal phenomenon, so that E(t) = E en (t) cos(φ ce + 2πf c t), where E en (t) is the pulse envelope and f c -denotes a carrier optical frequency. In the early days of ultrafast science, the treatment of the interaction of femtosecond pulses with the matter was also considered as only a temporal phenomenon. However, it became evident that space-time coupling can't be neglected in the femtosecond [19,20] or attosecond science [21,22]. Within focal volume the temporal profile of the pulse may vary significantly [22][23][24]. Taking into account space-time coupling is critical for analysing experimental data, design, control and avoiding artifacts in the ultrafast science experiments [23,25].
Recently, an analytical approach has been introduced to estimate the CEP distribution in the focal volume [26,27]. While providing valuable insight, these analytical tools are limited only to very simple optical systems and field distributions.
In this work, we demonstrate a new numerical method for space-time calculation of the CEP, φ ce . Our approach is not limited to simple apparatus, rather is well suited for applications to complex, real world, optical systems. We consider several basic key examples of focusing optics, calculate explicitly CEP within the region of interest and discuss implications of our findings.

METHOD
Evolution of ultrashort pulse electric field in space and time, E( r, t), can be decomposed into its spectral components, Ẽ ( r, f ), by a Fourier transform: (1) For the sake of simplicity, we will assume scalar approximation (i.e. assuming that E( r, t) = E( r, t)x) , even though our approach doesn't necessitate this and can be used for analysis of fully vectorial fields. We assume ultrashort pulse propagation along theẑ axis.
Ultrashort laser pulse, E(x, y, z = 0, t), entering an optical system at z = 0 is decomposed into monochromatic fields,Ẽ(x, y, 0, f ), using Eq.1 at the entrance of the system to be analysed. Each of these monochromatic field components of the pulse is propagated using optical code for the propagation of the continuous wave (CW) light fields into spatial region of interest r ∈ R 3 interest . For example, for many scientific and technological applications, such a region of interest is a focal volume of an optical system. Then, for each point r in the region of interested, the complex frequency field components are recomposed together using inverse Fourier transform, E( r, t) = ∞ −∞Ẽ ( r, f )e i2πf t df ≡ F −1 {Ẽ( r, f )}, to obtain a temporal laser field evolution of the optical pulse.
In our work, we have chosen to use commercially available optical codes, intended primarily for lens design, for the complex light field propagation. Using such widely used commercial optical codes has several distinct advantages. First, these codes are usually suitable for analysis of real world optical systems of significant complexity and are not limited to simple systems. For example, field propagation in optical systems consisting of many lenses, mirrors, prisms, gratings, etc. can be easily performed. Second, there are several optical codes available, which have been well tested, debugged and optimized during many years of real-world scientific and industrial applications. Thus we leverage these advantages in our approach. However, these optical codes, capable of physical optics simulations (i.e propagation of complex optical field and not only ray optics) are usually intended for use for the propagation of the continuous wave (CW) rather than ultrashort optical pulses. Decomposing ultrashort pulse into its monochromatic spectral components allows to easily adopt the CW optical codes for the analysis of ultrashort pulses [28]. Specifically, we used OSLO [29] code in this work, but any other code for CW optical field propagation could be used instead, for example Zemax [30] or CodeV [31].
Once the temporal evolution of the pulse is reconstructed in the spatial region of interest, we find the CEP, φ ce ( r), by extracting the phase between the peak of the pulse envelope, at time t max , and the closest peak of the carrier wave, by a conventional formula φ ce ( r) = Arg(E( r, t max )) = arctan(Re(E( r, t max )), Im(E( r, t max ))). The result is unwrapped -whenever the reconstructed CEP jumps between adjacent points in space more than or equal to π radians, the multiples of 2π is added until the jump is less than π.
To facilitate numerical calculations, instead of the continuous Fourier transform (Eq. 1), we use a Discrete Fourier Transform (DFT). Working with DFT has several significant consequences that should be taken into account in the context of our analysis to achieve efficient calculations and to avoid distortion and aliasing artifacts. The spectral amplitudẽ E( r, f ) of the pulse is sampled within the pulse bandwidth at N points with frequencies The sampling bandwidth has to be chosen to represent the spectral bandwidth of the actual pulse accurately. Choosing the bandwidth too narrow will result in pulse distortion in the time domain, and choosing the bandwidth that is much wider than the actual pulse bandwidth would result in inefficient use of the computer resources.
Special care should be taken in choosing the sampling spacing ∆f to minimize the aliasing [32] effect. Using a discrete sampling instead of continuous spectral amplitude, , results in periodic pulses in time domain with a period T = 1/∆f . As the pulse propagates along optical system, it often disperses to some characteristic time scale τ max . The τ max depends both on material dispersion and on the system geometry. In the case of Fourier limited pulse, the shorter the pulse at the entrance of the system, the wider bandwidth it has and more dispersion it will experience that will result in the larger τ max . To avoid the aliasing, it is important to choose T > τ max , implying that the sampling spacing should be chosen so that ∆f < 1/τ max . Strictly speaking, finite pulse bandwidth necessarily implies infinite pulse duration in the time domain.
However, the pulse will decay in time to negligible amplitude. So intelligent judgment has to be made to set the frequency sampling spacing dense enough so that aliasing will become negligibly small and at the same time not to over-sample to avoid wasting the computer resources.
Note that frequency comb spacing (∆f ) used for calculations should not necessarily be identical to the frequency spacing (f r ) of the real laser system frequency comb. For example, the repetition rate of commonly used CEP stabilized femtosecond amplifiers is just of the order of 1 kHz, implying f r = 1 kHz. Such a dense sampling of typical femtosecond pulse spectra would be impractical due to the computer's memory limitations and is not needed to achieve accurate space-time pulse evolution computations and the CEP retrieval.
The frequency spectrum of a typical femtosecond pulse usually spans near-infrared, visible and sometimes ultraviolet (UV) regions. Hence, there is no spectral content in the pulse from zero to about hundreds of Terahertz. To make computations more efficient, it is useful to eliminate this zero-padded region of spectra by shifting the pulse spectral envelope, After space propagation of thẽ E env ( r, f ), applying the convolution theorem we reconstruct the pulse in the time domain, The DFT is implemented in the software using a Fast-Fourier Transform (FFT) algorithm [32]. The FFT algorithms work much more efficiently when the length of the vector to be transformed is a power of 2 of a positive integer number (N = 2 m , m ∈ N). We achieve this condition by slightly reducing spacing in the frequency domain or applying zero-padding in the time domain.

RESULTS
We start by analyzing several simple well-known scenarios in order to confirm the validity of our approach and of the computer simulations.
Good example of such a scenario is focusing of the femtosecond pulse by a lens as shown in the Fig. 1. In this example, the lens is an aspheric Geltech lens made of C0550 Corning glass with effective focal lens f EFL = 3.1 mm. The aspheric surface is chosen to correct spherical aberration, so that chromatic aberration is the main remaining aberration [28].  Focusing is defined as an instance in time with maximum pulse intensity. An emblematic "horseshoe" shape and a forerunner pulse are observed as expected [33]. Figure 1 (b) shows the same pulse 235 fs after the focusing. Also, a typical "v" shape of the femtosecond pulse spatial field distribution after the focus [33] was clearly observed as well.  phase and groups velocities, in collimated beam this CEP shift will be largely the same across the beam, hence no meaningful spatial variation of the CEP is expected. The small observed CEP variation (< 10 −3 rad) results from subtle diffraction of the collimated beam propagating through the slab. In the case of a spherical focusing mirror, the CEP variation is larger than in previous cases but still negligible (∼ 10 −3 rad) for any practical purpose.
However, the situation changes entirely for the refracting focusing element. A Geltech aspheric lens (solid line in the Fig. 2)  As the refractive aspheric lens has exhibited significant CEP variation across the focal plane, we add to our analysis another refractive focusing element often used in the femtosecond optics -an achromatic doublet (for example, in the refractive femtosecond pulse shapers [20]). In this example, we have chosen an uncoated AC254-030 lens by Thorlabs with a focal length f = 30 mm. Further details about this achromatic doublet lens are given in the supplemental document. While for an achromatic doublet the chromatic aberration is notably compensated, a femtosecond pulse propagation through such a lens experiences significant material dispersion. We explore how CEP spatial variations within the focal volume of a lens depend on the focusing conditions. Figure 3 shows CEP variation for the aspheric lens (panel (a)) and for the achromatic doublet (panel (b)) for different NA of the focusing system. For very "relaxed" focusing conditions, such as N A 0.03, which corresponds to the waist size of w=1 mm for the achromatic doublet, there is almost no variation of the CEP phase within the focal volume. However, when the focusing becomes tighter, the variation of CEP within the focal volume increases significantly.
Even at still very gentle focusing conditions of NA=0.13, the variation of the CEP within the focal volume exceeds 2π radians for the achromatic doublet. These variations of CEP are very significant; to put it into perspective, the difference between the cosine-like pulse and sine-like pulses is only π/2 radians. Note that even at small NA values ( example of NA 0.1 is shown in the Fig. 3 (b)) for the short focal length aspheric, meaning the impinging beam width of only few hundred microns), there is still meaningful residual variation of the CEP, particularly near the periphery of the focal volume. Further reduction of the impinging beam width doesn't translate into vanishing CEP phase as occurs for the long focal length lens ( Fig. 3 (a)). We attribute this to the diffraction effect that becomes significant as the beam width approaches the wavelength scale of the incoming laser pulse.
Our results indicate that it is really an interplay between geometrical (focusing) effect, diffraction and chromatic dispersion that leads to CEP variation across the beam-any of these contributions alone would not cause spatial CEP variations. Spatial distributions of the electrical field at the focus for the aspheric lens and achromatic doublet are shown in the Figure 4 for their highest respective NAs. Note that a significant field curvature in both cases accompanies substantial variations in the CEP within a focal volume. The characteristic horseshoe shape [33] in the focus is clearly visible for both lenses.
Similar to the increase of CEP variation with an increase in NA, we observe an increase in the field curvature as well as in the horseshoe shape appearance as we increase the focusing NA for both lenses.

DISCUSSIONS AND CONCLUSIONS
Our results show that there is a significant variation of CEP in the focal volume of the refracting focusing element. Spatial CEP variations strongly increase with a numerical aperture. We have observed a strong correlation between the CEP spatial variations and curvature of the pulse wavefront within the focal region. Use of refractive optics is quite common and ubiquitous in ultrafast science and applications, particularly with lower energy pulses. As we do not observe spatial CEP variations with a focusing mirror, where just a reflection plays a role, our results indicate that it is a combination of geometrical focusing and dispersion that causes such a phenomenon. Our findings imply that if one observes or drives physical phenomena dependent on the CEP with such an optical system, the result may be obscured or washed out by averaging within a focal volume. This means that it is essential to estimate the spatial behaviour of CEP within a region of interest.
To address this challenge, we have introduced an approach for a robust and accurate quantitative evaluation of the spatial behaviour of CEP and the pulse field evolution in space and time. Our method is not limited to simple setups. It can particularly useful as a powerful tool for designing and engineering of the ultrafast optical systems, where CEP matters.
In future work, it will be very interesting and useful to synergistically use both analytical methods [26,27] and our numerical approach. Measuring the spatial distribution of CEP in the lab is a non-trivial [16,17], time and resource-consuming task. Instead, our method can be used, providing the results much faster at a fraction of cost and effort, when compared to the actual lab experiments. Thus, our approach will drive major advances not only in the CEP related experimental work, but also in the development of the theory and the development of new techniques for suppression of the CEP variations in the focal volume.