Resonant excitation of plasma waves in a plasma channel

We demonstrate resonant excitation of a plasma wave by a train of short laser pulses guided in a pre-formed plasma channel, for parameters relevant to a plasma-modulated plasma accelerator (P-MoPA). We show experimentally that a train of $N \approx 10$ short pulses, of total energy $\sim 1$ J, can be guided through $110$ mm long plasma channels with on-axis densities in the range $10^{17} - 10^{18}$ cm$^{-3}$. The spectrum of the transmitted train is found to be strongly red-shifted when the plasma period is tuned to the intra-train pulse spacing. Numerical simulations are found to be in excellent agreement with the measurements and indicate that the resonantly excited plasma waves have an amplitude in the range $3$ - $10$ GV m$^{-1}$, corresponding to an accelerator stage energy gain of order $1$ GeV.

In the laser wakefield accelerator (LWFA) [1], a short laser pulse propagating through a plasma excites a trailing Langmuir wave, within which the generated electric fields can be of the order E wb = m e cω p /e, where ω p = (n e e 2 /m e ϵ 0 ) 1/2 is the plasma frequency, and n e is the electron density.For electron densities of interest E wb ∼ 100 GV m −1 , some three orders of magnitude greater than is possible in a conventional accelerator.Considerable progress has been made, including, for example, the acceleration of electrons to energies in the GeV range in centimetre-scale accelerator stages [2][3][4][5][6][7][8][9][10], and the application of LWFAs to driving compact light sources [11,12].Recently, free-electron laser gain was demonstrated using laser-accelerated electrons [13,14].
To drive a large amplitude Langmuir (or 'plasma') wave, the duration τ L of the laser pulse must satisfy τ L ≲ T p /2, where T p = 2π/ω p is the plasma period, corresponding to τ L ≲ 100 fs for plasma densities of interest.As a consequence, recent experimental work has been dominated by the use of high energy (joulescale) chirped-pulse-amplification [15] Ti:sapphire lasers.However, this laser material has a high quantum defect (34%) [16] which limits the pulse repetition rate of highenergy systems to f rep ≪ 1 kHz.
An alternative method for driving the plasma wave is to resonantly excite it with a train of low-energy pulses (or a single long, modulated pulse) in which the pulse spacing (or modulation) is matched to T p .An example of this approach is the plasma beat-wave accelerator (PBWA) [1,[17][18][19], in which two long pulses of frequencies ω 1 and ω 2 = ω 1 + ω p are combined to form a pulse modulated at ω p .Beat-wave acceleration of electrons to energies in the 10 MeV range has been reported; of particular relevance to the present work is that by Tochitsky et al. [19], who exploited ponderomotive self-guiding over 3 cm to accelerate electrons to 38 MeV at a gradient of ∼ 1 GV m −1 .
Interest in resonant wakefield excitation has revived [21] with the development of novel laser technologies, such as thin-disk lasers that can generate joule-scale pulses at f rep in the kilohertz range, with high (≳ 10%) wall-plug efficiency [22].The picosecond-duration pulses provided by these systems are too long to drive a plasma wave directly, and a second laser frequency separated by ω p is not currently available to drive a PBWA.A potential solution is the plasma-modulated plasma accelerator (P-MoPA) [23], which comprises three stages: (i) a modulator, in which a long (∼ 1 ps), high-energy (≳ 1 J) laser pulse is spectrally modulated by the low amplitude plasma wave driven by a short (≲ 100 fs), low-energy (≲ 100 mJ) 'seed' laser pulse as they co-propagate in a plasma channel of on-axis density n e,0 ; (ii) a dispersive optical system that converts the spectral modulation to a train of short pulses spaced by T p,0 = 2π m e ϵ 0 /n e,0 e 2 ; (iii) an accelerator stage, also of on-axis density n e,0 , within which the pulse train resonantly drives a large amplitude plasma wave.Numerical simulations [23] show that a 1.7 J, 1 ps driver, with a 140 mJ, 40 fs seed, could accelerate electrons to energies of 0.65 GeV in a 100 mmlong plasma channel with n e,0 = 2.5 × 10 17 cm −3 .
In this Letter we investigate experimentally the accelerator stage of a P-MoPA.We demonstrate guiding of a train of N ≈ 10 short pulses, with a total energy of the order 1 J through 110 mm long plasma channels, equiva- lent to 14 Rayleigh ranges, with n e,0 in the range 10 17 − 10 18 cm −3 .Resonant excitation of a plasma wave within the channel is evidenced by the observation of strong redshifting of the spectrum of the transmitted pulse train when T p,0 was tuned to the pulse spacing in the train.The results are found to be in excellent agreement with numerical simulations, which show that wake amplitudes in the range 3 GV m −1 to 10 GV m −1 were achieved, corresponding to an accelerator stage energy gain of the order 1 GeV.
Figure 1 shows schematically the arrangement employed for these experiments, undertaken with the Astra-Gemini TA3 Ti:sapphire laser at the Rutherford Appleton Laboratory.This laser provides two synchronized beams, here denoted the 'drive' and 'channel-forming' beams, each of central wavelength λ 0 = 800 nm with transform-limited full-width at half-maximum (FWHM) duration of 31 fs.In order to mimic the pulse train employed in the P-MoPA scheme, single laser pulses were converted to a train of short pulses using a Michelson interferometer, as sketched in Fig. 1(a) and described previously [24,25].The temporal intensity profile of the generated pulse train, shown in Fig. 1(b), was determined from single-shot measurements of the spectrum and autocorrelation of the train (see Supplemental Material [20] for further details).
The gas target used in this work was a cell-jet hybrid [10,26], with hydrogen gas pulsed into the target via a solenoid valve and two transducers measuring the pressure on-shot.The laser pulses were coupled into, and out of, the target via a pair of 3 mm radius coaxial pinholes mounted on: (i) the front of the target; and (ii) a motorized plunger that could be moved to adjust the target length L. A relative RMS pressure variation along the laser propagation axis of 4.1 % was measured [20], as shown in Fig. 1

(e).
A hydrodynamic optical-field-ionized (HOFI) channel [27,28] was formed in the target by focusing the channel-forming pulse, of energy ∼ 100 mJ and FWHM pulse duration 80 fs, with an axicon lens of base angle 3.6 • .The transverse intensity profile of the beam produced by the axicon had a central maximum of FWHM spot size (9.8 ± 0.1) µm, as shown in Fig. 1(c).
The pulse train, of total on-target energy E train = (2.5 ± 0.5) J, was focused by an off-axis f /40 paraboloid to the target entrance.The transverse intensity profile of the focused beam [see Fig. 1(d)] was found to have a 1/e 2 intensity radius of (45.5 ± 3.4) µm, a Rayleigh range of z R = (7.9± 0.7) mm, and to contain (64.9 ± 1.5) % of its energy within its FWHM.The delay between the arrival of the channel-forming and drive beams was set to t d = 3.5 ns.After leaving the plasma channel, the energy of the drive beam was reduced, and the beam re-imaged onto a 16-bit camera and a fibre-coupled spectrometer.An example guided mode is shown in Fig. 1(e).
The excitation of plasma waves by the drive pulse was detected through changes in its spectrum [29].The spectra presented in Figs. 2 and  defined as S(λ) = λS meas (λ)/ ∞ 0 λS meas (λ)dλ, where S meas (λ) is the measured spectrum.Figure 2(a) shows S(λ) for an incident pulse train with E train = (2.5±0.5)J and τ = 170 fs, at on-axis densities approximately equal to, and one third of, the resonant value, n e,res ≈ 4.3 × 10 17 cm −3 .As expected, the input spectrum of the pulse train is modulated by the Michelson interferometer to yield N ≈ 10, uniformly-spaced peaks.For the offresonant density, the spectrum of the transmitted train is similar to that of the incident pulse, with some blueshifting apparent in the region λ ≲ 780 nm, likely caused by ionization of the neutral gas collar [30,31] surrounding the HOFI channel and of the gas plumes that extend beyond the target.In contrast, at the resonant density, considerable red-shifting is observed, extending the bandwidth of the input beam by more than 40 nm on the long wavelength side.The new red-shifted light beyond 820 nm is seen to consist of a series of peaks [23]; these arise from spectral modulation of the laser pulse by the wakefield, which generates copies of the input spectrum shifted by ±mω p for integer m.The peaks on the blue side of the spectrum are not visible in Fig. 2, likely due to the additional blue-shift from ionisation.We note that blue-shifting would have predominantly occurred for the first few pulses in the train, and, since the pulse train was negatively chirped, their initial spectra were on the blue side of the mean wavelength.
The density-dependent red-shift seen in Fig. 2(a) strongly indicates resonant plasma wave excitation in the plasma channel.To confirm this, we also measured the transmitted spectra for a temporally-smooth ∼ 1 ps drive pulse of similar energy at on-axis densities match- ing those in Fig. 2(a).As shown in Fig. 2(b), in this case no red-shift was observed, and the spectra were similar for both densities and were dominated by blue-shift of similar magnitude to that observed in Fig. 2(a).
Figure 3 shows the variation with on-axis plasma density of the transmitted spectra when the drive was wellguided [20] by the plasma channel.To quantify the redshift we define the red-shift metric R = ∞ λmin S(λ), where λ min is the longest wavelength in the input spectrum above the noise level.It is evident from Fig. 3(a,b) that the spectra of the pulse train driver exhibit a pronounced red-shift for densities in the range n e,0 = 4 ×10 17 cm −3 to 5 ×10 17 cm −3 , which agrees with the expected resonance density of n e,res = 4.3 × 10 17 cm −3 .For a train of N identical laser pulses, the full-width of the resonance peak is expected [25] to be δn e,0 /n e,res ≈ 8/(3N ), corresponding to δn e,0 ≈ 1.2 × 10 17 cm −3in good agreement with the measured FWHM in R of δn e,0 ≈ 1.6 × 10 17 cm −3 .In contrast, Figs.3(c, d) show that no resonance is observed for the unmodulated drive pulse.Significant red-shifting of the unmodulated drive pulse is observed for n e,0 ≳ 5.5×10 17 cm −3 , likely caused by self-modulation [32][33][34] of the long pulse.
To provide further insight, we compared these measurements with the results of an in-house 2D cylindrical fluid code, benchmarked against the particle-in-cell (PIC) code WarpX [35] (see [20]).The calculations used the retrieved pulse train parameters and modelled the plasma channel as an ideal fully-ionized parabolic waveguide [20].The code ignores the effects of ionization by the laser pulse, and assumes that the temporal envelope of the drive is unchanged by its interaction with the plasma.
Figure 3(b) shows the calculated R for the τ = 170 ns, E train = 2.5 J pulse train in a plasma channel of length L = 110 mm.It can be seen that the position and width of the calculated resonance peak agree closely with those observed in the measurements.For some shots the measured R values reach the calculated curve, but in most cases they are lower.In order to understand this, the energy transmission of the train was measured as a function of the plasma channel length [20].For each cell length the measured energy transmission was found to vary over a wide range, owing to the large pointing jitter of the input pulse train.Shots for which the input beam was well aligned with the channel axis were found to have an input coupling of T 0 = (64 ± 4) %, which is consistent with |c 0 | 2 = (71 ± 5)%, where c 0 is the calculated [20] coupling coefficient between the transverse amplitude profile of the input beam and that of the lowest-order mode of the channel.In contrast, the coupling coefficient deduced from all guided shots is only T 0 = (32 ± 13) %, which reflects the additional losses arising from misalignment with respect to the channel axis.Figure 3(b) shows that if the drive energy is reduced by this factor, i.e. to E train = 800 mJ, the calculated variation of R with density is in excellent agreement with the averaged measurements.At the resonant density, the amplitude of the wakefield driven by the τ = 170 fs pulse train is calculated from the fluid simulation to be 10 GVm −1 (3 GVm −1 ) for E train = 2.5 J (0.8 J).
Further evidence of resonant wakefield excitation is shown in Fig. 4, which shows the measured and calculated variation of R with on-axis density for pulse trains with τ = 200 fs and 170 fs.In this case, E train = (2.5 ± 0.5) J and L = 70 mm.It can be seen that, for both pulse trains, the position, width, and magnitude of the measured variation of R agree well with the calculation assuming E train = 0.8 J.At higher densities, n e,0 ≳ 7 × 10 17 cm −3 , red-shifting arising from self-modulation is again observed.
It has been previously shown that HOFI [27,28] channels achieve higher energy transmission when the wings of the laser pulse have sufficient intensity to ionize the neutral gas collar to form a conditioned [36,37] HOFI channel.PIC simulations [20] of the present experiment indicate that the leading three pulses in the train con- ditioned the HOFI channel, allowing later pulses in the train to be guided with low losses.We note that conditioning of the channel could also be achieved by employing a separate, short pulse immediately ahead of the pulse train [30]; the required energy of the conditioning pulse is ∼ 7 mJ per cm of channel, i.e.only 3% of the drive energy in the present experiment.
In summary we have demonstrated guiding of a train of N ≈ 10 short pulses, with a total pulse train energy of the order 1 J through 110 mm long plasma channels with onaxis densities in the range 10 17 − 10 18 cm −3 .The spectra of the transmitted pulse trains were found to be strongly red-shifted when the plasma period was matched to the pulse spacing in the train.In contrast, no such resonance in the red-shift was observed for an unmodulated drive pulse of the same total energy and duration.Numerical simulations were found to be in excellent agreement with the measurements, and showed that, at resonance, the wake amplitude was in the range 3 − 10 GV m −1 , corresponding to an accelerator stage energy gain of the order 1 GeV.
These results constitute the first demonstration of resonant excitation of a plasma wave by a train of laser pulses guided in a pre-formed plasma channel.The laser and plasma parameters employed in this work are directly relevant to the accelerator stage of the P-MoPA scheme [23], which offers a route to achieving kilohertzrepetition-rate, GeV-scale plasma accelerators driven by plasma modulation of joule-scale, picosecond-duration laser pulses, such as those provided by thin-disk lasers.This research was funded in whole, or in part, by EP-SRC and STFC, which are Plan S funders.For the purpose of Open Access, the author has applied a CC BY public copyright licence to any Author Accepted Manuscript version arising from this submission.
The HOFI plasma channel was generated by the channel-forming beam, reflected by a holed mirror (hole radius 15 mm) to enable coupling of the channel-forming and drive beams into the gas target co-axially.The beam energy remaining at the interaction point was measured using an energy meter to be ∼100 mJ.The channelforming beam was focused to a longitudinally-extended focus by an UVFS axicon lens with 3.6 • base angle, resulting in 1.6 • approach angle of the rays to the axis.The axicon had a hole in the centre of radius 13 mm, and was placed (560±2) mm from the front plate of the gas target.The axicon focus had a Bessel-function transverse profile with a measured full-width at half-maximum (FWHM) spot size of (9.8 ± 0.1) µm.

Pulse train generation
Pulse trains were generated using the same method used in our earlier work [24,25].A Michelson interferometer was installed after the final amplifier in the laser chain, before the grating compressor.A flat (unwedged) compensator plate was used in the 'reflected' arm of the Michelson to account for the one extra pass of the beamsplitter glass made by the 'transmitted' arm.The mirror on the transmitted arm was placed on a motorised stage for remote control of the Michelson arm optical path difference, δX Mich .It was also possible to bypass the Michelson in the laser chain, resulting in a smooth, unmodulated pulse of the same FWHM pulse duration of 1 ps.
The compressor was set for partial compression, to leave a spectral chirp on the laser beam.This resulted in unwanted higher-order phase terms being present in the partially-compressed pulse.These third and fourth order delay (TOD and FOD) spectral terms affect the uniformity of the pulse spacing in the pulse trains.An acousto-optic programmable dispersive filter (AOPDF or 'Dazzler') was used to reduce the magnitude of the TOD and FOD terms in order to generate a pulse train with uniform spacing.
We can extract a measure of uniformity for these pulse trains, U , defined as U ≡ 1 − (τ max − τ min )/2τ avg where τ max is the maximum spacing, τ min is the minimum spacing and τ avg is the average spacing between the pulses.A TOD value of 1 × 10 4 fs 3 was found to correspond to a uniformity measure of U ∼ 0.9 on each train.As the number of pulses in the train was relatively low, on the order of 10, the width of the resonance was expected to be wide enough not to be significantly affected by variations in uniformity on the scale of 10%, hence the uniformity of the pulse trains achieved was deemed to be sufficient.

Pulse train measurement
The temporal intensity profiles of the pulse trains were deduced from single shot autocorrelator (SSA) measurements, in conjunction with information about the settings of the grating compressor and measurements of the input spectrum, as in our earlier work [24,25].
The measurements of the pulse train were conducted as follows.With the laser in a pulsed alignment mode, the pulse train was intercepted prior to its focus and recollimated by a lens with focal length f = −1000 mm.The beam was then sent out of the vacuum chamber via a thin optical window, and redirected to the SSA diagnostic.This consisted of a 50/50 beamsplitter, where one arm was directed via a delay stage and the other via a mirror such that the two beams crossed at an angle 4 • to each other.At the crossing point, a non-linear crystal (BBO Type 1, 38.8 • cut angle) was placed for second harmonic generation, with the generated 400 nm light leaving the crystal perpendicular to the plane of the crystal.The unconverted 800 nm light was dumped.A lens was used to image the blue light in the plane of the crystal on to a CCD.To increase the efficiency of blue light generation, a pair of cylindrical lenses (f = 30 mm) were placed in the path of the two crossing beams, one before and one after the crystal, to focus the 800 nm beams in the plane perpendicular to their crossing angle.

Pulse train retrieval
The compressor grating separation and Michelson spacing were set to generate a pulse train with the desired pulse spacing, and its spectrum was measured.As described above, a controlled amount of TOD was introduced by an AOPDF in order to ensure that the pulse spacing within the train was uniform.First, a model of the compressor grating was used to estimate the TOD term expected for the grating separation.This gave an estimate of the TOD value that needed to be applied with the AOPDF.The uniformity was then optimised by scanning the TOD value about the estimate from the compressor model and analysing the fringe visibility of each pair of peaks in the autocorrelation signal, V = (I max −I min )/(I max +I min ) where I max and I min were the maximum and minimum intensities of each peak.For a pulse train with equal pulse spacing, the autocorrelation peaks from each pair of pulses would appear at precisely the same position on the x-axis of the crystal.
If the pulse spacing varies between pulses, the peaks in the autocorrelation function would be slightly offset from each other in x, resulting in a smearing of the peaks, and hence maximising V maximizes the uniformity of the pulse trains.
Assuming the TOD and FOD terms were zero after applying this method, the only remaining parameter to be determined was the average spacing of the pulses in the train τ = ϕ (2) 2πc/δX Mich where ϕ (2) is the group delay dispersion (second order phase).The value of δX Mich was first estimated from the spectrum of the pulse train, by fitting the modulated spectrum to the spectrum of the fully-compressed pulse multiplied by a modulation function with a period of Ω = 2πc/δX Mich .The value of ϕ (2)  was estimated from a model of the grating compressor.These estimates for δX Mich , ϕ (2) and TOD were used to give the starting points and bounds on the retrieved pulse train properties (e.g. the retrieved value for δX Mich was bounded to be within 10% of the estimated value, and the value for ϕ (2) was bounded to be within 3000 fs 2 of the estimated value).The retrieval method progressed via a numerical optimisation algorithm to find the values for δX Mich , ϕ (2) and TOD that would generate the measured auto-correlation profile, while using the spectrum of the fully-compressed pulse as a constraint.This method enabled the retrieval of δX Mich and ϕ (2) with an uncertainty < 10%.
The calculated values for the spectral phase terms along with the measured power spectrum amplitude give a full description of the pulse train in the spectral domain, and hence the temporal profile could be calculated via a Fourier transform.Figure 6 shows the raw SSA signals and temporal intensity profiles retrieved from the autocorrelation for the two pulse trains used in the experiment.

Drive beam focus
The spot size of the multi-pulse drive beam at focus was calculated using the D4σ method to be (41.5 ± 2.7) × (49.4 ± 2.0) µm 2 , along the minor and major axes of the ellipse respectively.The vacuum Rayleigh range was found to be z R = (7.9± 0.7) mm from a fit of the measured spot size as a function of longitudinal position to w(z The spatial jitter of the drive beam focus was measured to be 25.1 × 37.2 µm 2 RMS, which is on the order of a spot size.The jitter in the channel-forming beam focus position was much smaller than that of the drive beam, measured to be 5.2 × 2.4 µm 2 RMS, owing to the smaller effective f -number of the axicon focus (f # = 1/(2θ) = 11.5)compared to the drive beam focus (f # = 40).FIG. 7. Raw data of the pulse-train input spectrum measured by a spectrometer in the laser area for a series of 6 successive shots, demonstrating the shot-to-shot variation in the position of the spectral peaks.Note that the envelope is not the true spectral intensity envelope of the pulse train as the measurements here have not been corrected for the non-uniform white-light response of the measurement system.

On-shot input spectrum measurement
The spectrum of the modulated drive beam was measured directly after the the pulse-train-generating Michelson interferometer.It was found that the spectral peaks in the drive beam shifted within the spectral envelope of the beam from shot to shot (see Fig. 7).This spectral jitter was believed to arise from small fluctuations in the optical path lengths in the Michelson arms.
In the analysis of the spectral data described in this study, the shot-to-shot spectral shifts on the input beam were accounted for with the following procedure.First, a reference spectrum for the dataset was selected and the peaks in the spectrum identified.Here, a dataset means a series of data taken within a few hour time period, with fixed settings for the pulse-generating Michelson and compressor.Then, on each shot, the peak positions in the input spectrum of that shot were compared to the positions of the reference spectrum to give a shift value for each peak, δλ.The mean of the peaks shift values was calculated, δλ.The output spectrum, f raw (λ) for that shot was then translated as f (λ) = f raw (λ + δλ).Accounting for the shifts in the input spectrum in this way enabled the features in the transmitted spectra to be more directly compared.

Gas target
The cell-jet hybrid gas target [10,26] used in this experiment consisted of a cylindrical chamber of length 110 mm, closed on each end by a plate with a 6 mm diameter circular hole in the centre, to enable the laser pulses to enter and leave.The back plate was mounted on a motorized plunger that allowed the length of the gas target, L, to be varied.The gas inlet was designed in a cone-shape to disperse the gas uniformly in the target.
Hydrogen gas was pulsed into the target via a solenoid valve, opened before the arrival of the laser pulses to ensure that the laser-plasma interaction occurred in a steady-state gas condition.

Longitudinal density profile
The longitudinal plasma density profile depended predominantly on the longitudinal profile of the gas density inside the gas target at the moment of ionisation.During the experiment, the gas pressure was measured at two locations by pressure transducers connected to the top edge of the gas target, one near the entrance (z = 29 mm from the front pinhole) and one near the centre (z = 59 mm).Measured pressures in the target were in the range 10 mbar to 150 mbar and had an uncertainty of 3 mbar.The ratio of these two readings as a function of pressure is shown in Fig. 8(a).The mean ratio of these data is 0.79, suggesting that the centre pressure was consistently higher than that at the front by ∼ 20%.
To gain a better understanding of the full gas profile inside the gas target, calculations were performed using the computational fluid dynamic code code saturne [38,39].Figure 8(b) shows the calculated gas density distribution for the case of 1 bar backing pressure, in steady-state conditions.The simulations show that the gas density directly above the inlet position was greater than in the surrounding gas.In simulation, the ratio of the central pressure to the pressure at the position of the front trans-ducer measured on the top edge of the target shown in Fig. 8(c) was found to be 0.77, close to that found in the experiment.However, the pressure ratio between the same positions along the laser axis, shown in Fig. 8(d), was found to be 0.98, indicating a uniform longitudinal density profile.This was confirmed via separate plasma fluorescence measurements of the longitudinal gas pressure profile along the axis of laser propagation, performed at the Oxford Plasma Accelerator Laboratory [10].These measurements indicated relative gas pressure fluctuations of 4.1 % RMS along the length of the target, consistent with the simulated pressure profile along its central axis [Fig.8(d)].

Selecting well-guided shots
Due to the significant transverse jitter of the drive beam focus position (∼ 30 µm), there were a number of shots for which the drive beam was not coupled into the plasma channel.To distinguish the signal in the data corresponding to well-guided shots, a selection procedure was employed to filter the data before the analysis.The images from the exit mode diagnostic were used to define the selection criteria.When the drive beam was wellguided, the image showed a well-defined spot or group of spots.The images were analysed to give the number of spots, the spot sizes and the average pixel counts within each spot.When the beam was not guided, the image consisted of a low-amplitude, speckled light distribution.The analysis of the images for the non-guided shots identified many low-signal spot regions due to the non-uniform light distribution in the image.For the shot to be considered well-guided, the image was required to satisfy the following criteria; 1.The image must contain fewer than 4 spots; 2. All spots identified in the image must have a spot size smaller than some threshold value for that dataset; 3. The average pixel count for each identified spot must be above some threshold value.
The threshold values for spot size and average pixel counts were determined for each dataset separately.One dataset refers to a series of shots taken with the same input beam properties over a few-hour time period.This accounted for any differences in the alignment of the imaging system to the forward diagnostics between datasets.Histograms of the spot size and average fluence distributions for a representative dataset are shown in Figure 9.A Gaussian mixture model (GMM) [40] was fit to the data in each histogram to identify peaks in the distribution.The GMM returns a sum of n Gaussian components for each distribution, with the i th peak having a mean x i and standard deviation σ i .The number of peaks in the fitted GMM was allowed to vary between 1 and 6.The optimal value of n was selected using the Akaike information criterion [41], a metric used to evaluate GMMs that includes both a goodness-of-fit measure and a preference for fewer components over larger models to avoid over-fitting.
When considering the spot size, the Gaussian peak corresponding to the set with the largest mean spot size was taken as the set of failed shots.Therefore, the threshold spot size was calculated as x S,i=n−1 +3×σ S,i=n−1 , where n is the number of fitted components.When considering the average counts within the guided spot, the Gaussian peak in the GMM corresponding to the lowest mean pixel counts was taken as the set of failed shots.In this case, the threshold value was set at x C,i=1 + σ C,i=1 .This procedure removed between 25% and 50% of the total shots from each dataset.

On-axis plasma density calibration
Previous work has demonstrated that the on-axis density of a HOFI plasma channel is roughly linearly proportional to the initial background gas density, with a proportionality factor, Ξ, of between 5 and 10 [31,36] depending on the delay between the channel-forming and drive pulses.In this experiment, features in the spectrum of the transmitted pulse trains could be used to extract the value of this proportionality factor and allow conversion between pressure and on-axis plasma density of the plasma channels.
The spectral intensity of the input pulse train is modulated at Ω = 2πc/δX Mich .Therefore, a Fourier transform of the spectral intensity, F T [ Ĩ(ω)], has a dominant peak at 2π/Ω = δX Mich /c.case with pulse spacing 170 fs, which exhibits a peak at δX Mich /c ≈ 450 fs.Note that this peak does not occur at the pulse spacing τ = 170 fs, since the pulse spacing is determined by δX Mich and the GDD of the chirped input pulse.Near resonance, the high-amplitude wakefield spectrally modulates the laser pulse train, generating side-bands in the spectrum at ±mω p .In the linear regime, these side-bands have the same spectral shape as the input laser spectrum but have a reduced amplitude and are offset by multiples of ω p .In the case where only the highest amplitude side-bands are visible (with m = 1), we would expect the spectral intensity of the transmitted pulse to consist of the input pulse train spectrum, and two copies of the input pulse train spectrum, shifted by ω p on either side of the central frequency of the pulse train spectrum ω 0 .Therefore, the Fourier transform of the spectral intensity at resonance is expected to consist of the peak at δX Mich /c, as well as an additional peak at 2π/ω p = T p , as demonstrated in Fig. 10(a).
The value of the proportionality constant Ξ was extracted from the data using a fit to the observed peak positions, T p = (161 ± 17) fs, for 80 shots around the observed resonant pressure [see Fig. 10(b)].To perform this analysis the spectral intensity below λ min ≈ 820 nm was set to zero, to focus on the new red-shifted light only, Ĩred (ω), and the position of the highest amplitude peak in |F T [ Ĩred (ω)]| was selected as T p .This is related to the on-axis density via T p = 2π/ω p and hence implies n e,0 = (4.8± 0.5) × 10 17 cm −3 .The factor to convert between the initial electron density n e and the on-axis density n e,0 was then calculated from the fit to the data as Ξ = n e /n e,0 = (8.1 ± 0.6).This value is within the range expected for HOFI plasma channels [31,36].
When the gas target was set shorter than full length, the gas reading on the pressure transducers was reduced for the same backing pressures.This may have been due to a partial covering of the transducer outlet by the gas target plunger and/or a change in the flow dynamics.The on-axis density calibration factor was re-calculated for the shorter gas target length used, L = 70 mm, to account for any differences in the pressure transducer reading.The calibration factor in this case was Ξ = (6.2± 0.5).

Guiding Joule-scale pulse trains over ten centimetres
In this section, additional measurements pertaining to the demonstration of guiding of pulse trains in long plasma channels are presented.Unless stated otherwise, the pulse trains used for this work had a pulse spacing of τ = 170 fs.

Drive beam input mode analysis
It is desirable to match the input beam profile [Fig.11(a)] with the lowest-order mode of the plasma channel, as the lowest-order mode has the longest attenuation length.To estimate the fraction of the input beam energy that is coupled to the lowest-order mode of the plasma channel in a well-aligned case, an overlap integral can be calculated between the input beam intensity profile and the Hermite-Gaussian (HG) modes with a spot size matched to the channel.The transverse electron density profile of the channels were not measured in this experiment and so the modes of the channel could not be calculated directly.Therefore, to estimate the lowest-order mode of the channel for this calculation, the transverse intensity profile of the best-guided output beams were used to estimate the matched spot size: w 0 = (40 ± 3) µm [see Fig. 11(b)].This method assumes that the higher-order modes radiated away completely by the end of the channel, leaving only the lowest-order channel mode remaining.Fifty images of the input focus were used, and in each case the overlap integral of the image with the corresponding HG mode was calculated.HG modes with all combinations of indices n, m varying between 0 and 20 were used.The average value of the coefficients over the 50 events, |a nm | 2 , are shown in Fig. 11(c), with the dominant mode contributions and their RMS fluctuations plotted in Fig. 11(d).This analysis indicated that (71 ± 5) % of the input mode intensity was in the lowest order mode matched to the channel on-average.This is a simplified analysis of the coupling of the input beam into the plasma channel, but represents an estimate of the theoretical maximum coupling efficiency of the input beam.In reality, a number of effects including density ramps at the channel entrance, the effects of conditioning of the channel by the drive beam, and the significant spatial jitter of the drive beam would all be expected to further impact the coupling of the beam into the channel.

Drive beam transmission as a function of plasma channel length
To understand what fraction of the drive beam energy was coupled into the channel, the transmission of the drive beam was measured as a function of the length of the plasma channel by moving the position of the motorised plunger on the gas target.The energy in the transmitted beam was estimated by summing the pixel counts within the image of the exit modes, and comparing these to equivalent measurements of the input beam.For each shot, a correction factor for the spectral response of the CMOS sensor and transport optics within the imaging system was calculated to account for the change in the spectrum of the guided pulse arising from the plasma interaction.
The energy transmission as a function of gas target length is plotted in Fig. 12, for the case of a pulse train with spacing τ = 200 fs at an on-axis density of n e,0 = (2.9 ± 0.3) × 10 17 cm −3 .It can be seen that within each gas target length setting, there is a significant variation in energy transmission.This is mainly a result of the significant pointing jitter of the drive beam, which causes the transverse position of the input pulse to vary with respect to the axis of the plasma channel, reducing the fraction of energy that can be coupled into the guide.
The event with highest transmission for each channel length is expected to correspond to conditions where the drive beam was most closely aligned to the channel axis, and therefore can be used to estimate the maximum coupling efficiency that was achieved in the experiment.The data for these events were fitted to an exponential decay of the form T (z) = T 0 exp(−z/L att ), with T the energy transmission, T 0 the coupling efficiency at the channel entrance, z the propagation length and L att the power attenuation length.The fit yields T 0 = (64 ± 4) % and L att = (124 ± 20) mm.The fitted coupling efficiency, T 0 = (64 ± 4) %, can be compared to |c 0 | 2 , where c 0 is the calculated coupling coefficient between the transverse amplitude profile of the incident beam and that of the lowest-order mode of the channel.This was found [20] to be |c 0 | 2 = (71 ± 5) %.Since higher-order modes will be attenuated rapidly, the value of T 0 calculated by projecting T (z) to the channel input will be approximately |c 0 | 2 , and the good agreement between |c 0 | 2 and the value of T 0 deduced from the well-aligned shots is consistent with this picture.A fit to the average transmission values in Fig. 12 gives a coupling efficiency of T 0 = (32 ± 13)% and L att = (96 ± 20) mm, which takes into account the additional coupling losses due to spatial jitter of the drive beam with respect to the channel axis.

On-axis plasma density dependence
In Fig. 13, the spot size of the transmitted pulse train at the output of the plasma channel is plotted against the on-axis density in the plasma channel, over a range of densities explored in this experiment.The horizontal error bars represent the uncertainty on the measured pressure, combined with the uncertainty of the on-axis density calibration.The vertical error bars represent the standard error of the spot size calculated over all shots within each bin.

2D Fluid Code
The wakefield excited by the pulse train was solved in 2D-cylindrical coordinates following the procedure outlined in Ref. [44].For quasistatic linear wakefields driven by the normalised pulse intensity envelope |a| 2 in axisymmetric channels n e (r), the fluid equation for the perturbation of the normalised potential δϕ = ϕ − 1 ≈ |a| 2 /4 − v z /c and corresponding electron density pertur- bation δn is given by; where ξ = z − ct is the co-moving longitudinal coordinate and k p (r) = ω p (r)/c denotes the local plasma wavenumber.This PDE can be solved numerically for any arbitrary pulse envelope |a| 2 ≪ 1 and axisymmetric plasma channel n e (r).Previous work suggests that the transverse plasma density profiles of HOFI channels are approximately parabolic [31].Hence for this calculation, we assumed that the plasma channel took the form of the matched parabolic channel; n e (r) = n e,0 + ∆n(r/w 0 ) 2 where n e,0 is the on-axis plasma density and ∆n = (πr e w 2 0 ) −1 is the channel depth parameter with r e being the classical electron radius.This channel guides the fundamental Gaussian mode |a| 2 ∼ exp(−2r 2 /w 2 0 ) of spot size w 0 .
According to the paraxial wave equation, assuming that the ponderomotive envelope |a| 2 is in the fundamental Gaussian mode and remains fixed throughout the full propagation in this matched parabolic channel, and that the density perturbation is small relative to the channel depth parameter |δn| ≪ ∆n, the spectral modulation of the pulse train will be given by [42]; )rdr denotes the intensity-weighted transverse average.Equation (4) was used to calculate the spectrum of the pulse train after propagation within the channel over a range of plasma densities, using the pulse envelope extracted from the pulse train retrieval process.The resulting spectra were analysed to give red-shift values, R, using the same method as for the experimental data.

PIC simulations of conditioning effect
Conditioning of the neutral gas collar formed by the HOFI channel was studied using the PIC code WarpX [35], which uses ADK theory [43] to model the tunnel ionisation of neutral atoms in intense laser fields.For simplicity, a neutral gas of atomic hydrogen was assumed.The HOFI channel profile, which is comprised of the ionised electron n e,OFI and ion n H + ,OFI densities from the OFI process and the density of neutral hydrogen atoms n H , was parameterised in the following form n e,OFI (r) = n H + ,OFI (r) = n e,0 exp −r 2 /d 2 0 , n H (r) = n max e −(r−rmax) 2 /d 2 1 , r ≤ r max n ∞ + (n max − n ∞ ) e −(r−rmax) 2 /d 2 2 , r > r max (6) where n e,0 = 4.3 × 10 17 cm −3 , n max = 6.2 × 10 18 cm −3 , n ∞ = 3.2 × 10 18 cm −3 , d 0 = 45 µm, d 1 = 16 µm, d 2 = 8 µm, and r max = 50 µm.An example of the conditioning effect of the pulse train is shown in Fig. 14, which shows that full conditioning of the neutral gas collar is achieved after the first few low-energy pulses have passed.

FIG. 1 .
FIG. 1. Sketch of the experimental layout.(a) Illustration of the pulse train generation scheme.(b) Example single-shot autocorrelator (SSA) measurement for the τ = (170 ± 2) fs pulse train.Upper: comparison between the measured (pink) and retrieved (grey, dashed) SSA signal.Lower: retrieved pulse train intensity profile.(c) Measured axicon focus.(d) Example input mode of the focused multi-pulse drive beam.(e) Example guided mode at the channel exit.All focal spot images are normalized to their maximum.(f) Comparison between the measured and simulated longitudinal gas pressure profile [20].

FIG. 3 .
FIG. 3. Density dependence of S(λ) for: (a,b) a pulse train with τ = 170 fs, Etrain = (2.5 ± 0.5) J; and (c,d) an unmodulated ∼ 1 ps, (2.7 ± 0.5) J pulse.(a) and (c): S(λ), averaged in electron density bins of width ∆ne,0 = 0.24 × 10 17 cm −3 .(b) and (d): mean values of R [black squares] weighted by energy transmission; red circles indicate bins containing data from only a single shot.The ne,0-error bars are a combination of the uncertainties in the measured pressure and the on-axis density calibration.The R-error bars represent the standard error on the mean.In (b), individual data points are plotted (grey circles) and the orange dotted line represents ne,res(τ = 170 fs).Overlaid are the results of the fluid calculations for Etrain = 2.5 J (green) and Etrain = 0.8 J (blue).

FIG. 4 .
FIG. 4. Variation of R with on-axis density for a plasma channel of length L = 70 mm and for pulse trains of energy (2.7 ± 0.5) J and pulse separation: (a) τ = 200 fs; and (b) τ = 170 fs.The results of the fluid calculations, assuming Etrain = 800 mJ, are shown by the blue dashed lines.For each plot the expected resonant density is indicated by the orange dotted line.

FIG. 5 .
FIG. 5. Schematic of the method used to generate a pulse train for the experiment.The chirped laser pulse (left) passes through a Michelson interferometer, with the two arms separated in optical path length by δX Mich .The interference pattern following the Michelson forms a pulse train, with an example intensity temporal profile shown at the top.

FIG. 6 .
FIG. 6. Example of the raw SSA signal (top), measured and retrieved SSA profiles (middle), and deduced temporal intensity profile (bottom); (a) for the case of 170 fs pulse spacing; (b) for 200 fs pulse spacing.

FIG. 8 .
FIG. 8. Gas pressure and density profiles inside the gas target.(a) shows the experimentally measured ratio of the pressure measured by the front transducer to that measured by the centre transducer.(b) shows the central slice of the 3D gas density distribution inside the target calculated using code saturne.(c) and (d) show the gas pressure and density profile lineouts extracted from the 3D fluid simulations for; (c) the top edge of the target (transducer axis); and (d) the centre line of the target (laser axis).

FIG. 9 .
FIG. 9. Histograms showing; (a) the spot size; and (b) the average pixel count in the spot, extracted from the analysis of the exit mode diagnostic images for the dataset with pulse spacing 170 fs and gas target length 110 mm.When multiple spots were detected in the image, the spot with the largest spot size is plotted.The black curve indicates the fitted Gaussian mixture model.The red dashed line shows the median value of the data in each case.The purple dashed line shows the threshold values calculated by the method described in the text.All data with spot size greater than the spot size threshold value and average counts lower than the average counts threshold value were discarded.
Figure 10(a)  shows an example

FIG. 10 .
FIG. 10.Procedure for the on-axis density calibration using signatures in the spectrum of the transmitted pulse train, for an example dataset with laser pulse spacing 170 fs and gas target length 110 mm.(a) shows the absolute value of the Fourier transform of the spectral intensity, |F T [ Ĩ(ω)]| of the input pulse train (grey, dashed) and a transmitted pulse train at resonance (blue, solid).A strong peak in the input spectrum is observed at δX mich /c ∼ 450 fs and a new peak in the resonant spectrum appears at Tp = 2π/ωp.(b) Extracted value of Tp for 80 shots around the resonant pressure, indicating an average peak position of Tp = (161 ± 17) fs.

FIG. 11 .
FIG. 11.Comparison between examples of the (a) input mode and (b) guided mode.(c) Average value of the coefficients |anm| 2 for overlap integrals of 50 input modes with Hermite-Gaussian (HG) modes of spot size w0 = (40 ± 3) µm, as indicated by analysis of the highest throughput guided modes.(d) The dominant contribution is the lowest order mode: |a00| 2 = (0.71 ± 0.05).

FIG. 12 .
FIG. 12. Measured energy transmission of the pulse train, T , with spacing τ = 200 fs as a function of plasma channel length, L, for an on-axis plasma density of ne,0 = 2.9 × 10 17 cm −3 .The transmission of individual events are represented by the scattered grey data points, with the maximum at each value of L plotted in green.Black, square data points show the average transmission for each channel length, with error bars representing the standard error on the mean.

FIG. 13 .
FIG. 13.Measured variation of spot size of the transmitted pulse train (Etrain = (2.5 ± 0.5) J, τ = 170 fs) as a function of on-axis plasma density inside the channel (L = 110 mm).The dashed blue line and shaded region shows the input spot size and RMS jitter.Red markers indicate density bins with only one shot contributing.

FIG. 14 .
FIG.14.Results of a 2D PIC simulation of the propagation of the pulse train and its conditioning of the HOFI channel, whose profile is modelled with the parameterisation given in Eq. (5).The top panel displays the normalised intensity envelope |a| 2 , whilst the bottom panel shows the normalised electron density ne/ne,0 including both the electrons initially ionised by the OFI process and the electrons from ionising the neutral gas collar.