Application of Matched-Filter Concepts to Unbiased Selection of Data in Pump-Probe Experiments with Free Electron Lasers

Abstract

Pump-probe experiments are commonly used at Free Electron Lasers (FEL) to elucidate the femtosecond dynamics of atoms, molecules, clusters, liquids and solids. Maximizing the signal-to-noise ratio of the measurements is often a primary need of the experiment, and the aggregation of repeated, rapid, scans of the pump-probe delay is preferable to a single long-lasting scan. The limited availability of beamtime makes it impractical to repeat measurements indiscriminately, and the large, rapid flow of single-shot data that need to be processed and aggregated into a dataset, makes it difficult to assess the quality of a measurement in real time. In post-analysis it is then necessary to devise unbiased criteria to select or reject datasets, and to assign the weight with which they enter the analysis. One such case was the measurement of the lifetime of Intermolecular Coulombic Decay in the weakly-bound neon dimer. We report on the method we used to accomplish this goal for the pump-probe delay scans that constitute the core of the measurement; namely we report on the use of simple auto- and cross-correlation techniques based on the general concept of “matched filter”. We are able to unambiguously assess the signal-to-noise ratio (SNR) of each scan, which then becomes the weight with which a scan enters the average of multiple scans. We also observe a clear gap in the values of SNR, and we discard all the scans below a SNR of 0.45. We are able to generate an average delay scan profile, suitable for further analysis: in our previous work we used it for comparison with theory. Here we argue that the method is sufficiently simple and devoid of human action to be applicable not only in post-analysis, but also for the real-time assessment of the quality of a dataset.

1. Introduction

Recently we described the results of a pump-probe experiment, in which the lifetimes of doubly excited states of neon dimers were measured [1]. The dimers were excited by absorption of two EUV photons from the Free Electron Laser (FEL) FERMI-1 [2] and probed via ionization by a UV laser pulse.

The excited dimers decayed by Interatomic Coulombic Decay to stable dimer cations Ne${}_2^{+}$, which were detected by a time-of-flight (TOF) mass spectrometer. Ionization by the UV pulse led to a repulsive state of the dimer, which dissociated, so that the yield of dimer ions was reduced, and the yield of Ne${}^{+}$ increased. The dimer sample used in that work was very dilute with a large atomic background, as it was produced in a supersonic expansion of neon gas, with a yield of about 1%. The data were therefore noisy, and a substantial effort was needed to analyze and filter them. In this paper, we describe the methods used for that analysis.

The data consisted of a number of temporal scans, each lasting approximately 30 min, in which the ion-TOF signal was measured as a function of the delay between pump and probe pulses, and these scans had poor signal-to-noise ratios. Simply averaging them did not give good results, probably because some of them suffered from poor FEL conditions or poor alignment of the pump and probe pulses. One should reject those scans that do not contribute significant signal, while avoiding human bias. A possible approach is to exploit the difference between signal and noise with respect to auto- and cross-correlation. This is an instance of the more general concept of “matched filter”, i.e., the linear filter maximizing the signal-to-noise ratio (SNR) of a measured noisy sample [3]. Under the assumption of white noise, the matched filter is the complex-conjugate time-reversal of the signal sought, i.e., application of the matched filter returns the autocorrelation of the signal [4]. The concept is widely applied in signal processing, where one is primarily interested in extracting a burst of periodic signal from a noisy sample. As far as peak detection is concerned, a popular field of application is chromatography, but the general results obtained there apply equally well to our case; Ref. [5] explicitly discusses the use of a matched filter to determine amplitude and time shift of the peak being sought. We consider the use of auto- and cross-correlation for three purposes:

• validate the scans to be included or excluded from averaging
• determine the weight with which each scan should enter the average
• determine, if desired, by how much to temporally shift a scan prior to averaging

2. Results

For a set of delays $\left\{d_i\right\}$, $i=0..n$, let us consider two delay scans $\mathbf{R}=\{R_0,\cdots ,R_n\}$ and $\mathbf{S}=\{S_0,\cdots ,S_n\}$ each consisting of a sequence of Ne${}_2^{+}$ ion-TOF signals, reduced as explained in Section 4, and padded with zeros outside of the delay range scanned. We will consider $\mathbf{R}$ our reference scan, which for convenience is assumed to be noise-free.

The definition of cross-correlation without normalization, ${\mathrm{\Gamma }}_{\mathbf{RS},q}$ is:

$${\mathrm{\Gamma }}_{\mathbf{RS},q}=\sum _i{R}_i{S}_{i+q}\equiv {(\mathbf{R}★\mathbf{S})}_q$$

(1)

Note that $\mathrm{\Gamma }$ is linear in the two sequences (scaling either one by a factor $\alpha$ scales $\mathrm{\Gamma }$ by the same factor), and that ${\mathrm{\Gamma }}_{\mathbf{SS},q}$ is the autocorrelation of $\mathbf{S}$, which has a maximum for $q=0$:

$${(\mathbf{S}★\mathbf{S})}_0=\sum _i{\left(S_i\right)}^2$$

(2)

We now assume that any scan $\mathbf{N}=\{N_0,\cdots ,N_n\}$ consisting of pure noise has zero (negligible) cross-correlation with any other scan (including itself except, obviously, at zero-shift), thus:

$$\begin{array}{ccc}\hfill {(\mathbf{N}★\mathbf{N})}_q& =& {\delta }_{q0}\sum _i{\left(N_i\right)}^2\hfill \end{array}$$

(3)

$$\begin{array}{ccc}\hfill {(\mathbf{R}★\mathbf{N})}_q& =& 0\hfill \end{array}$$

(4)

with ${\delta }_{q0}$ the Kronecker delta; we use Equation (3) to define a variance ${\sigma }^2={(\mathbf{N}★\mathbf{N})}_0={\sum }_i{\left(N_i\right)}^2$. Equations (3) and (4) strictly hold for white noise and infinite n; in a real-life situation we can expect that noise correlation just decays much faster than signal correlation: thus all scans whose autocorrelation is sharply peaked are probably pure noise, and will also have poor cross-correlation with the reference scan.

Any noisy scan $\mathbf{S}$ can be written in terms of the reference scan $\mathbf{R}$ and a pure-noise sequence $\mathbf{N}$ as: $\mathbf{S}=\alpha \mathbf{R}+\mathbf{N}$, with $\alpha$ a real number. In a worse-case scenario, $S_i$ may be shifted (for simplicity by an integer index r), that is: $S_i=\alpha R_{i+r}+N_i$. Then:

$${(\mathbf{R}★\mathbf{S})}_q=\alpha {(\mathbf{R}★\mathbf{R})}_{q-r}+{(\mathbf{R}★\mathbf{N})}_q=\alpha {(\mathbf{R}★\mathbf{R})}_{q-r}$$

(5)

note that $\mathbf{R}★\mathbf{S}$ is proportional to $\mathbf{R}★\mathbf{R}$ but shifted by r as expected; this is the argument invoked in Ref. [5] to associate the position of the maximum of the cross-correlation function to the shift of the peak being sought. In our work [1] it was not necessary to include a shift; while Equation (5) could be used to estimate the signal amplitude $\alpha$, it is preferable to use the autocorrelation instead:

$${(\mathbf{S}★\mathbf{S})}_q={\alpha }^2{(\mathbf{R}★\mathbf{R})}_q+{(\mathbf{N}★\mathbf{N})}_q={\alpha }^2{(\mathbf{R}★\mathbf{R})}_q+{\sigma }^2{\delta }_{q0}$$

(6)

(note that the shift by r cancels, here). Equation (6) tells us in particular that all autocorrelation sequences should be the same except for a scale factor ${\alpha }^2$ and a sharp noise peak at $q=0$. Visual inspection shows that this is true for a number of scans which we consider good; other scans exhibit a narrower, or structured, sometimes negative, autocorrelation, probably indicating correlated noise, i.e., a drift during the measurement; the rest only exhibits the noise peak, indicating no signal at all (Figure 1). From now on we will consider $(\mathbf{R}★\mathbf{R})$ to be noiseless and, as empirically found, satisfactorily approximated by a gaussian; we derive the width of the gaussian from our best scan (Figure 1).

Our choice of a gaussian is purely empirical, and our method does not critically depend on it: to the extent that all good scans have the same shape of the autocorrelation function, Equation (6), the best fit parameters will be the same for all good scans, except for a scale factor for the height. What is important is that the fitting function provides a reasonable approximation of $\alpha$ in Equation (6), i.e., of the peak value of the signal component of the autocorrelation.

It is nevertheless instructive to discuss some limiting cases for the shape of the autocorrelation function: we begin by noting that the expected shape of a pump-probe signal such as that of our experiment is an exponential decay (we ignore the possible complication brought on by the presence of more than one decay constant) convoluted with the instrumental resolution (in our case a gaussian, coming from the finite duration of the FEL and UV pulses). The respective autocorrelation functions are an exponential and a gaussian. For the case of non-negligible instrumental resolution the autocorrelation will resemble a gaussian near $q=0$, but will have broader wings decaying as $exp(-|q\left|\right)$ rather than $exp(-q^2)$. An early truncation of the scan will clip the wings of the autocorrelation function. In our case an early truncation does somewhat contribute to determining the shape of the autocorrelation function, but does so equally for all scans (the pump-probe delay was scanned in reverse, and point index $i=0$ in Figure 1 consistently corresponds to the maximum value of the delay).

Let us finally consider two possible sources of artefacts, namely a constant, or a linearly drifting, baseline (remembering that the scans are padded with zeros outside of the delay range scanned). This would primarily contribute a slowly decaying component to the autocorrelation signal (of width comparable to that of the scan itself). Let us however note that because the scans are acquired by rapid double-background subtraction (Equation (9) and Figure 4), we expect a complete baseline cancellation. Finally, a drifting signal would alter the shape of the peak and consequently of its autocorrelation function: assuming an even drift that causes a loss of the optimal experimental conditions (temporal or spatial overlap; quality of the focus; resonance wavelength) one can speculate that its main effect would be a reduction of the measured width.

3. Discussion

Given a set of scans $\left\{{\mathbf{S}}^{\left(k\right)}\right\}$, their weighted average is

$$\langle \mathbf{S}\rangle =\frac{{\sum }_k{w}_k{\mathbf{S}}^{\left(k\right)}}{{\sum }_k{w}_k}$$

(7)

$\langle \mathbf{S}\rangle$ is the quantity reported in our work ([1] note for the sake of exactness that in Figure 2a therein, the unity baseline was not subtracted); we want to determine which scans to include in the average, and their weights $w_k$.

From Equation (6) and Figure 1b we can estimate the signal-to-noise ratio as:

$${\mathit{SNR}}_{\mathit{k}}=\sqrt{\frac{{({\mathbf{S}}^{\left(k\right)}★{\mathbf{S}}^{\left(k\right)})}_0}{{({\mathbf{N}}^{\left(k\right)}★{\mathbf{N}}^{\left(k\right)})}_0}-1}=\sqrt{\frac{{\alpha }_k^2{(\mathbf{R}★\mathbf{R})}_0}{{({\mathbf{N}}^{\left(k\right)}★{\mathbf{N}}^{\left(k\right)})}_0}}=\frac{{\alpha }_k}{{\sigma }_k}$$

(8)

A simple analysis of its trend over the course of the experiment reveals some regularities that we exploit to qualitatively classify scans, and to define a quantitative criterion that we adopt to accept or reject them. When the scans are ordered in the sequence they were acquired (Figure 2a) no obvious trend is visible for the $SNR$, but one does note a large number of scans with $SNR=0$: they are either those for which a gaussian fit of the autocorrelation was not successful, or those which have been excluded a priori (e.g., because the scan was aborted). The ordering of the scans by increasing $\mathit{SNR}$ (Figure 2b) reveals a gap between scans with $\mathit{SNR}<0.25$ and $\mathit{SNR}>0.45$. We cannot find an obvious reason for this gap; because this was one of the first resonant two-photon experiments performed at FERMI we can speculate a threshold behavior of some of the less-controlled parameters of the FEL (peak intensity; second harmonic content). In any case, we decided to use the condition $\mathit{SNR}>0.45$ as a discriminant to include a scan in the average.

Let us now come to the weight with which each accepted scan enters the average. We show in Appendix A that the weight which maximizes the SNR of $\langle \mathbf{S}\rangle$ is $w_k={\alpha }_k/{\sigma }_k^2={\mathit{SNR}}_k/{\sigma }_k$; in Ref. [1] we used $w_k={\mathit{SNR}}_k$, which gives a slightly worse result, although the difference is not significant. We presume that the latter fact depends on the limited number of samples, the predominance of few of them (see Figure 2a), and the low dispersion of the ${\sigma }_k$. In both cases we observe an improvement of the signal-to-noise ratio of $\langle \mathbf{S}\rangle$ by a factor ≈2.25 relative to the best single scan.

Unfortunately at the time of the experiment we were not anticipating the need for this test, so we are not able to perform further checks on which parameters and events mostly affected our experiment. Likewise our primary goal was to identify a simple unbiased method of data selection, and we do not attempt a systematic analysis of its merits and limitations, for which we refer the reader to the vast existing literature ([6,7] and references therein); we simply note that the method is most easily applicable to the case of white noise and single- or well-separated peaks, although generalizations to other noise [6] or multiple peaks [8] have been discussed. Let us finally note that because the autocorrelation is the Fourier Transform of the power spectrum (Wiener-Khinchin Theorem), the same determination of the signal-to-noise ratio could have been accomplished by Fourier Transform of each scan (in the frequency domain, the white-noise component appears as a constant).

Despite its simplicity and limitations, we believe that our simple test can be further characterized and profitably applied in future experiments. Let us note that the original and most immediate application of the method addresses a situation common to the beginning of many experiments, namely the need of making a weak signal visible when the expected signal shape is not precisely known in advance.

4. Materials and Methods

4.1. Experimental

The experiment was performed at the Low Density Matter (LDM) beamline [9] at FERMI FEL-1 [2]. The EUV pulses had an average energy of 16 $\mathsf{\mu }$J, duration estimated between 60 and 80 fs FWHM, circular polarization, focal size 30 $\mathsf{\mu }$m FWHM, and were tuned to the wavelength resonant with the two-photon doubly excited target state, 75.65 nm. The UV pulses had an energy of 35 $\mathsf{\mu }$J, estimated duration 200 fs, focal size 80 $\mathsf{\mu }$m FWHM, fixed wavelength (261 nm). The Ne dimers were produced by adiabatic expansion of Ne gas at a pressure of 0.8 MPa and temperature of 190 K through a 100 $\mathsf{\mu }$m nozzle. The target ions were detected by the time-of-flight mass spectrometer of the LDM endstation [10].

4.2. Data Acquisition

The Free Electron Laser FERMI currently operates at a repetition rate of 50 Hz (10 Hz at the time of the experiment). This is determined by the operation rate of the LINAC that accelerates the electron bunches used to generate the FEL pulses. Single-shot ion TOF traces (Figure 3) are acquired, tagged with a progressive integer (“bunchnumber”) and stored along with a wide selection of experimental parameters and machine parameters, for post-processing. The relevant single-shot quantities for us are the FEL energy per pulse ($I0$), the UV energy per pulse ($I_{\mathrm{UV}}$), and the integrated area a of the TOF peak of interest (specifically the ${}^{20}$Ne${}_2^{+}$ peak).

A preset number of single-shot values $a_j$ at nominally identical conditions are acquired and aggregated into a datapoint $A_i$; one in three shots is a blank shot (no gas sample), which is subtracted from the average to eliminate spurious signals (Equation (9)). Although in reality the bunchnumber is never reset, i.e., each bunchnumber is unique since the inception of FERMI, for the purpose of this work we will think of it as spanning a range $j=\{1..m\}$ for each datapoint. A set of datapoints at a sequence of delay values constitutes a delay scan. Delay scans are averaged into a global delay scan, where each scan enters with the weight determined in Equation (8).

4.3. Data Reduction

From the single-shot values we build a datapoint as follows:

$$A_i=\frac{{\displaystyle \sum _{j=\mathrm{gas}}}{a}_j}{{\displaystyle \sum _{j=\mathrm{gas}}}{\left(I{0}_j\right)}^p}-\frac{{\displaystyle \sum _{j=\mathrm{blank}}}{a}_j}{{\displaystyle \sum _{j=\mathrm{blank}}}{\left(I{0}_j\right)}^p}$$

(9)

with the new index i labelling sets of nominally identical measurements at a series of pump-probe delays $d_i$. ${\tilde{A}}_i$ is the analogue of $A_i$ with UV off. We consider equally spaced delays, so $d_i=\Delta t\times i(i=0..n)$. Let us note that single-shot data can be excluded from the sum by applying in post-analysis suitable filter conditions (e.g., a threshold on $I0$, or on the quality of the FEL spectrum); for this reason the number of single-shots aggregated into a datapoint is not strictly constant. For each i we define the quantity $S_i=A_i/{\tilde{A}}_i-1$: this is the base datapoint of a delay scan (Figure 4). The integer p in Equation (9) is the normalization order; for the data presented in Ref. [1] we chose $p=0$, i.e., no normalization. Let us also note that the intensity of the UV pulse was kept constant and is very stable, so it does not appear in the data processing. Let us finally note that at negative delays, and large positive delays we expect the UV to have no effect, i.e., $S_i=0$.

Figure 4. Schematic representation of the generation of a delay scan $\left\{S_i\right\}$, as a function of pump-probe delay $\left\{d_i\right\}$ from single-shot data $a_j$. The green (red) dots symbolize shots that pass (fail) the filter conditions.

Figure 4. Schematic representation of the generation of a delay scan $\left\{S_i\right\}$, as a function of pump-probe delay $\left\{d_i\right\}$ from single-shot data $a_j$. The green (red) dots symbolize shots that pass (fail) the filter conditions.

5. Conclusions

Although rather simple, the method proposed has many advantages and applications. It provides a quick and reliable method to evaluate the quality of a measurement in real time, even when the shape of the expected signal is unknown. In fact, in the approximation of white noise, and of an unstructured pump-probe signal, it provides a first coarse estimate of the signal shape, and thus of the delay range that needs to be scanned. In post-analysis it provides an unbiased criterion to include each scan in an averaged measurement.