Three-dimensional momentum imaging of dissociation in flight of metastable molecules

Many processes in molecules happen on fast timescales. For example, rotations and vibrations typically occur on picosecond and femtosecond timescales, respectively. Molecular bond rearrangement and fragmentation often proceed on similarly short timescales, as shown in [1–5] and many others. Hence, femtosecond laser pulses, possessing temporal durations shorter than these timescales, can be used to influence and shed light on molecular dynamics [6–9].

Not all processes in molecules, however, proceed so swiftly. Multiply charged molecular ions can exist in metastable states that lead to fragmentation happening on long timescales that range from picoseconds to even seconds [10–12]. The lifetimes of these transient systems are governed by the relevant potential energy landscape and the mechanisms responsible for decay, which can include tunneling, predissociation, and radiative decay to repulsive states. Investigating the formation, properties, and decay of these metastable molecular ions experimentally and theoretically has been a prominent field of research (see review papers [10–12] and [13–24], for example).

We study the decay dynamics of metastable molecules by employing coincidence three-dimensional (3D) momentum imaging, which provides the complete 3D momenta of the fragments and therefore their kinetic energy release (KER) and angular distributions. This information can in turn facilitate understanding of the dissociation mechanism(s), demonstrated for example in [17, 25, 26]. Hence, the 3D momentum imaging technique has been a powerful tool in studies of molecular fragmentation following ionization by ultrashort intense laser pulses, x-ray (or extreme ultraviolet) photons, or fast ion impact [27, 28], as long as the breakup is prompt. Prompt breakup happens on a sub-picosecond timescale, much shorter than the flight times of the fragments to the detector. The ionization processes in such experiments can readily form multiply charged metastable molecular ions, seen for example in [17, 22, 23, 29–34].

Importantly, in coincidence measurements, a metastable molecular ion may survive beyond the interaction region, traveling through the spectrometer for a non-negligible time before undergoing dissociation in flight. In this unimolecular fragmentation process, which is a subset of delayed dissociation, the survival time of a fragmenting metastable molecule, t_d, is a significant fraction of the time of flight (TOF) of intact metastable ions, t_m. In the present experiments, t_d is on the order of hundreds of nanoseconds to a few microseconds. In general, the observable range of t_d may differ depending on the specific metastable system, as well as the conditions of the experiment.

One of the striking signatures of dissociation in flight in the coincidence time-of-flight (CTOF) spectrum is a long, curved stripe that extends from the prompt breakup region and terminates at the intact metastable ion time of flight. This signature has been noted in CTOF spectra from a myriad of studies [35–48].

Field and Eland developed a method to extract the lifetimes of metastable states decaying in flight by fitting Monte-Carlo simulated time-of-flight-difference distributions N(${t}_{2}-{t}_{1}$) to the corresponding measured time-difference spectrum. Here, t₁ and t₂ are the times of flight of the first and second fragments, respectively. They have demonstrated their technique for a vast array of molecules [36]. Subsequent studies have implemented this method of lifetime retrieval for other metastable molecules [39, 41–44]. Recently, making some simplifying assumptions, Larimian et al [47] calculated t_d kinematically and retrieved the lifetime for deprotonation of the metastable ethylene dication. They also retrieved the momentum distribution of the fragments using Abel inversion [49] of the position image and discussed possible decay routes.

Our aim in this work is to present a more direct approach for extracting information about dissociation in flight of metastable molecules from coincidence momentum imaging measurements, employing basic principles. This method takes advantage of the known symmetries regarding the fragmentation process and allows one to extract not only the lifetime(s) of the metastable molecule but also the momenta, KER, and angular distributions. This information can enable pinpointing of the likely metastable state(s) dissociating in flight, as well as their dissociation mechanisms. Furthermore, this technique is general and can be applied to many different systems that undergo dissociation in flight. While our method is versatile and can provide a wealth of information to deepen understanding of dissociation processes, the focus of this manuscript is the analysis method of retrieving this information from the measurement.

To demonstrate our approach, we examine the same dissociation-in-flight channel as Larimian et al [47]. Namely, we look at deprotonation of metastable ethylene dications, C₂H${}_{4}^{2+}\to$ C₂H₃⁺ +H⁺ (as well as the deuterated equivalent), using the cold target recoil ion momentum spectroscopy (COLTRIMS) technique [27, 28]. Laser pulses with central wavelength of 790 nm, 23 fs duration (full width at half maximum (FWHM) in intensity), and peak intensity of about 3 × 10¹⁴ W cm⁻² are used to doubly ionize the ethylene molecules, introduced as a supersonic jet. The pulse duration was measured using second harmonic frequency-resolved optical gating (SHG FROG) [50], and the intensity was evaluated using the kink in the photoelectron spectrum of neon associated with 2U_p (where U_p is the average quiver energy of the free electron), which represents the transition from direct to rescattered electrons. To determine the 2U_p point, we measured the momenta of Ne⁺ recoil ions at low extraction field, following the method detailed in [51, 52].

As shown in figure 1(a), we observe the well-known signature of dissociation in flight, a long stripe in the CTOF spectrum. In figure 1(b), we note another signature in a density plot of the ion yield as a function of the measured TOF and position, N(TOF, x), where two 'halos' extend from the light and heavy fragments to the small C₂D${}_{4}^{2+}$ spot. Notably, in both the CTOF and N(TOF, x) spectra, the distributions of the two fragments converge to that of the intact metastable dication. Furthermore, as highlighted in figure 1(a), the predictions of Newton's equations of motion for dissociation in flight, which are marked by the open triangles, agree well with the measured data.

Zoom In Zoom Out Reset image size

To accomplish our goal of directly extracting information about dissociation in flight of a metastable molecule from our measurements, we start from the kinematic equations. The coordinate system utilized is depicted in figure 2(a). We employ principles similar to those that have been applied for years to image collision- and laser-induced prompt molecular dissociation in CTOF measurements [27, 28]. The equations of motion become slightly different from those for prompt dissociation to account for the survival time of the dication, t_d. For two-body dissociation in flight of a generic metastable dication, AB${}^{2+}\to$ ${{\rm{A}}}^{+}+$ B⁺, we have the following equations in the x direction, which in our case is along the laser beam propagation and transverse to the spectrometer axis:

$$\begin{eqnarray*}&&{x}_{j}-{x}_{0}={v}_{0x}{t}_{j}+{v}_{{jx}}^{{\prime} }({t}_{j}-{t}_{d})\qquad [j=1,2]\end{eqnarray*}$$

$$\begin{eqnarray}&&\displaystyle \sum _{j=1}^{2}{m}_{j}{v}_{{jx}}^{{\prime} }=0.\end{eqnarray} \tag{ 1 }$$

Here, x_j is the measured position of fragment j, and x₀ is the initial position of the metastable dication. v_0x is the average initial x-component velocity of the dication in the laboratory frame (${v}_{0x}\simeq 0$ for a cold jet), ${v}_{{jx}}^{{\prime} }$ is the x-component dissociation velocity in the AB²⁺ center-of-mass frame, t_j is fragment j's time of flight, and m_j is its mass. Clearly, t_d is needed to properly calculate the transverse momenta of the fragments. Note that the y-component equations (transverse to the spectrometer axis and along the jet flow) are similar to the x equations, except that v_0y, the supersonic jet velocity, is not negligible.

Zoom In Zoom Out Reset image size

At this point, it is interesting to contrast momentum imaging of dissociation in flight to that of prompt breakup, for which ${t}_{d}\simeq 0$. In the case of prompt fragmentation, one can readily see that the transverse equations of motion are decoupled from motion along the z direction, which is parallel to the spectrometer axis. For the problem of dissociation in flight, however, this is not the case, as t_d leaves us with more unknowns than equations in the transverse directions. Thus, we need to determine t_d from the z-component kinematic equations first.

3.1. t_d and lifetime determination

How exactly can we retrieve t_d from the z-direction motion? First, we write the z-component equations of motion:

$$\begin{eqnarray}&&{\ell }-\displaystyle \frac{1}{2}{a}_{m}{t}_{d}^{2}={({a}_{m}{t}_{d}+{v}_{{jz}}^{{\prime} })({t}_{j}-{t}_{d})+\displaystyle \frac{1}{2}{a}_{j}{({t}_{j}-{t}_{d})}^{2}}^{}\qquad [j=1,2].\end{eqnarray} \tag{ 2 }$$

Here, ℓ is the ion flight distance, a_m is the AB²⁺ acceleration, ${v}_{{jz}}^{{\prime} }$ is the z-component dissociation velocity of the jth fragment in the AB²⁺ center-of-mass frame, and a_j is its acceleration. We proceed to write the equations of motion in a more convenient dimensionless form. To that end, we multiply equation (2) by $2/{a}_{m}$:

$$\begin{eqnarray}&&\displaystyle \frac{2{\ell }}{{a}_{m}}-{t}_{d}^{2}=2\left({t}_{d}+\displaystyle \frac{{v}_{{jz}}^{{\prime} }}{{a}_{m}}\right){({t}_{j}-{t}_{d})+\displaystyle \frac{{a}_{j}}{{a}_{m}}{({t}_{j}-{t}_{d})}^{2}}^{}\qquad [j=1,2].\end{eqnarray} \tag{ 3 }$$

Noting that the AB²⁺ TOF is ${t}_{m}=\sqrt{2{\ell }/{a}_{m}}$, we replace the first term on the left-hand side of equation (3) with t_m². Further, dividing both sides by t_m² leads us to dimensionless z-component equations of motion. We also write the equation for momentum conservation in the AB²⁺ center-of-mass frame:

$$\begin{eqnarray}&&1-{t}_{{dm}}^{2}=2\left({t}_{{dm}}+\displaystyle \frac{{v}_{{jz}}^{{\prime} }}{{v}_{m}}\right)({t}_{{jm}}-{t}_{{dm}})+{\eta }_{j}{({t}_{{jm}}-{t}_{{dm}})}^{2}\qquad [j=1,2]\\ &&\displaystyle \sum _{j=1}^{2}{m}_{j}{v}_{{jz}}^{{\prime} }=0.\end{eqnarray} \tag{ 4 }$$

Here, we have defined ${t}_{{dm}}\equiv {t}_{d}/{t}_{m}$, ${t}_{{jm}}\equiv {t}_{j}/{t}_{m}$, and ${\eta }_{j}\equiv {a}_{j}/{a}_{m}$. Also, we have substituted in ${v}_{m}={a}_{m}{t}_{m}$, the velocity of the dication.

The equations above suggest that we can solve for t_dm (and hence t_d), as we have three equations and three unknowns, t_dm, ${v}_{1z}^{{\prime} }$, and ${v}_{2z}^{{\prime} }$. Combining the above equations of motion and the equation for momentum conservation, we eliminate ${v}_{{jz}}^{{\prime} }$ by substitution, resulting in an equation that can be solved for t_dm. Several subsequent algebraic steps lead us to the following quadratic equation:

$$\begin{eqnarray}&&{{at}}_{{dm}}^{2}+{{bt}}_{{dm}}+c=0\\ &&a=[2\beta ({\eta }_{1}-1)+({\eta }_{2}-1)]{t}_{1m}+[\beta ({\eta }_{1}-1)+2({\eta }_{2}-1)]{t}_{2m}\\ &&b=2[\beta (1-{\eta }_{1})+(1-{\eta }_{2})]{t}_{1m}{t}_{2m}+(1+\beta )-(\beta {\eta }_{1}{t}_{1m}^{2}+{\eta }_{2}{t}_{2m}^{2})\\ &&c=(\beta {\eta }_{1}{t}_{1m}+{\eta }_{2}{t}_{2m}){t}_{1m}{t}_{2m}-({t}_{1m}+\beta {t}_{2m}).\end{eqnarray} \tag{ 5 }$$

Here, $\beta ={m}_{1}/{m}_{2}$ is the ratio of the mass of the light fragment to that of the heavy fragment.

While the simplicity of our derivation and resulting equations makes the calculation of t_d and subsequent momentum imaging look quite straightforward, the problem is more convoluted than it initially seems. The quadratic equation for t_d has two solutions, which we denote as ${t}_{d}^{(+)}$ and ${t}_{d}^{(-)}$, where the superscripts correspond to the sign that is chosen in the quadratic formula.

For the vast majority of events, it is not clear which solution is correct, as both are physical based on obvious criteria: t_d must be real and $0\lt {t}_{d}\lt {t}_{m}$. The root of this complication, illustrated in figures 2(c) and (d), is the link between the sign of ${v}_{1z}^{{\prime} }$ and the correct quadratic formula sign choice. Specifically, ${t}_{d}^{(+)}$ always corresponds to ${v}_{1z}^{{\prime} }\lt 0$, meaning the H⁺ (or D⁺) is ejected away from the detector, and ${t}_{d}^{(-)}$ corresponds to ${v}_{1z}^{{\prime} }\gt 0$. Of course, in setting out to solve for t_d, one does not know a priori the sign of ${v}_{1z}^{{\prime} }$. Moreover, both situations—positive ${v}_{1z}^{{\prime} }$ with delay time ${t}_{d}^{(-)}$ and negative ${v}_{1z}^{{\prime} }$ with ${t}_{d}^{(+)}$—lead to the same measured set of time-of-flight values $({t}_{1},{t}_{2})$. This dilemma also thwarts retrieval of ${v}_{1z}^{{\prime} }$ and the z-component momenta, for which one needs to properly evaluate t_d.

As a noteworthy aside, this dual-solution t_d retrieval problem belongs to an extensive family of inverse problems, in which one is trying to retrieve initial conditions from observable parameters. This is an important problem faced in a wide array of fields in science and mathematics, such as medical imaging, x-ray crystallography, optics, geology, acoustics, and many others [53–55]. Frequently, the solution to an inverse problem is not unique, as is the case for dissociation in flight.

To address the inverse problem at hand, we utilize symmetry concepts. For any given t_d, reflection symmetry of ${v}_{1z}^{{\prime} }$ about 0 is expected, as the light fragment is equally likely to be ejected in either the forward or backward directions. Thus, one could select ${t}_{d}^{(+)}$ or ${t}_{d}^{(-)}$ randomly with equal probability and fit an exponential decay function to the resultant N(t_d) distribution to get a lifetime τ.

In the case of dissociation in flight, however, the symmetry of the overall ${v}_{1z}^{{\prime} }$ distribution should be broken to some degree due to the lifetime of the metastable molecule. That is, if τ is the lifetime of the dication, the ratio of the number of molecules that survive for times ${t}_{d}^{(+)}$ and ${t}_{d}^{(-)}$, ${N}^{(+)}$ and ${N}^{(-)}$, respectively, is described by

$$\begin{eqnarray}&&{N}^{(+)}/{N}^{(-)}={{\rm{e}}}^{-[{t}_{d}^{(+)}-{t}_{d}^{(-)}]/\tau }.\end{eqnarray} \tag{ 6 }$$

Thus, one could correct for the symmetry breaking using this factor. The extent of this correction is determined by the magnitude of $| {t}_{d}^{(+)}-{t}_{d}^{(-)}| /\tau$. Since ${t}_{d}^{(+)}$ and ${t}_{d}^{(-)}$ are typically not dramatically different compared to the lifetime τ, this correction is small.

Thus, we use a self-consistent approach in which we start as suggested above, by choosing ${t}_{d}^{(+)}$ or ${t}_{d}^{(-)}$ randomly with equal likelihood for each event. Recall that this choice is exact in the limit $| {t}_{d}^{(+)}-{t}_{d}^{(-)}| /\tau \ll 1$. An exponential decay function N(t_d) $=\,{N}_{0}{{\rm{e}}}^{-{t}_{d}/\tau }$ is then fitted to the resulting N(t_d) distribution to retrieve the lifetime. This lifetime allows computation of the aforementioned factor ${N}^{(+)}/{N}^{(-)}$, given in equation (6), which is then used to weight the choice of ${t}_{d}^{(+)}$ or ${t}_{d}^{(-)}$ in the next iteration. The obtained N(t_d) distribution is again fit with an exponential decay function to retrieve a more accurate lifetime, again allowing calculation of a new weighting factor for the choice of ${t}_{d}^{(+)}$ or ${t}_{d}^{(-)}$. This process is repeated until the lifetime τ converges. Note that for a given iteration, the choice of ${t}_{d}^{(+)}$ or ${t}_{d}^{(-)}$ and the fitting procedure is repeated for multiple trials to account for the finite sample size of our data. Also, the lifetime τ used to compute the weighting factor for the subsequent iteration is the mean value of those obtained in the multiple trials. For more details about our iterative approach, visit appendix B.

When this self-consistent method is applied to our data, we find good convergence within just a few iterations. We performed the analysis on the C₂H₄ and C₂D₄ data to explore the possibility of isotopic effects. Typical trials from the final iteration are shown in figures 3(a) and (b). Two-term exponential decay fits agree well with both measurements, suggesting that at least two metastable states are contributing to the observed dissociation in flight of ethylene dications. The final converged lifetimes for C₂H${}_{4}^{2+}$ are ${\tau }_{1}=202\pm 10$ ns and ${\tau }_{2}=916\pm 40$ ns. For C₂D${}_{4}^{2+}$, we obtain ${\tau }_{1}=269\pm 29$ ns and ${\tau }_{2}=956\pm 83$ ns. The errors here represent the standard deviation of the trials in the last iteration.

Zoom In Zoom Out Reset image size

Note that the weighting factor ${N}^{(+)}/{N}^{(-)}$ has been modified to account for the two lifetimes, ${\tau }_{1}$ and ${\tau }_{2}$, and thus becomes

$$\begin{eqnarray}&&{N}^{(+)}/{N}^{(-)}=\displaystyle \frac{{N}_{0}^{(1)}{{\rm{e}}}^{-{t}_{d}^{(+)}/{\tau }_{1}}+{N}_{0}^{(2)}{{\rm{e}}}^{-{t}_{d}^{(+)}/{\tau }_{2}}}{{N}_{0}^{(1)}{{\rm{e}}}^{-{t}_{d}^{(-)}/{\tau }_{1}}+{N}_{0}^{(2)}{{\rm{e}}}^{-{t}_{d}^{(-)}/{\tau }_{2}}}.\end{eqnarray} \tag{ 7 }$$

As mentioned, the two-term exponential fit indicates that at least two metastable states are responsible for the observed dissociation in flight. The single lifetime for C₂H${}_{4}^{2+}$ reported by Larimian et al, 498 ± 12 ns [47], lies between the two lifetimes that we have measured. This discrepancy could be due to a number of reasons, such as differences in the laser pulse parameters, the method used to compute t_d, or the t_d range chosen for the exponential decay fit. Note that when we perform a single exponential decay fit, shown in figure 3(a), we obtain a lifetime of 491 ± 19 ns, consistent with the previous measurement.

We also note that our measurements suggest a possible small isotopic effect in the shorter lifetime, ${\tau }_{1}$. The difference between the two shorter lifetimes is on the level of 2.2σ, while the longer lifetimes are the same within the measurement uncertainty. As dissociation in flight of ethylene dications is a low-rate channel, we expect that higher statistics data would make the presence or absence of an isotopic effect in the lifetimes more clear cut. Since we currently lack the good quality electronic structure information on these molecular ions needed to understand this isotopic effect and also to keep the focus on the method, we limit this discussion to highlighting the rich information afforded by our technique.

3.2. Momentum imaging

Having retrieved t_d and the lifetimes, we can proceed with performing momentum imaging, another aim of this work. While our solution of choosing ${t}_{d}^{(+)}$ or ${t}_{d}^{(-)}$ with some weighting works on a sample level and is thus a robust method for retrieving the lifetime(s), we note that for a large fraction of the individual events, ${t}_{d}^{(+)}$ or ${t}_{d}^{(-)}$ will be assigned incorrectly. Therefore, this method of computing t_d cannot be used for momentum imaging, as it is done on an event-by-event basis. As such, we need a single value of t_d for each event, even if it is approximate.

To approximate t_d, we neglect ${v}_{1z}^{{\prime} }$ in equations (4), as the term containing this quantity is typically on the order of a few percent compared to the other terms, as further detailed in appendix C. Moreover, as shown in the same appendix, the error that this approximation introduces in the retrieved t_d and the transverse momenta is also estimated to be at most a few percent. Having neglected ${v}_{1z}^{{\prime} }$, the equation for t_d becomes linear and thus has a single solution for each event, given explicitly by

$$\begin{eqnarray}&&{t}_{{dm}}=\displaystyle \frac{({\eta }_{1}{t}_{1m}^{2}-1)(1-{\eta }_{2})-({\eta }_{2}{t}_{2m}^{2}-1)(1-{\eta }_{1})}{2({t}_{2m}-{t}_{1m})(1-{\eta }_{1})(1-{\eta }_{2})}.\end{eqnarray} \tag{ 8 }$$

Employing this approximation of t_d, we have computed the transverse momentum of the D⁺ fragments from C₂D${}_{4}^{2+}$, shown in figure 4. The distribution of the radial momentum, ${p}_{1r}=\sqrt{{p}_{1x}^{2}+{p}_{1y}^{2}}$, shown in figure 4(b), agrees well with the functional form for the projection of an isotropic distribution onto a plane [56]. To verify that the momentum distribution is isotropic, the complete 3D momentum distribution is needed. Therefore, the p_1z component should also be measured.

Zoom In Zoom Out Reset image size

Is there any way to retrieve the z-direction momentum? Recall our initial fundamental problem in evaluating ${v}_{1z}^{{\prime} }$ that is associated with the ${t}_{d}^{(+)}$ and ${t}_{d}^{(-)}$ solutions. Now, we have gone even further and neglected ${v}_{1z}^{{\prime} }$ entirely, eliminating the possibility of recovering p_1z directly from the measurement. It is important to note that the polarization is typically aligned along the spectrometer axis (z direction) in COLTRIMS measurements in order to reduce losses of fast fragments, which are usually ejected along the laser field. This choice leads to equivalent x and y momentum components due to the axial symmetry about the laser field and prevents the direct determination of p_1z.

To retrieve the missing information, i.e., the p_1z momentum component along the laser field, we take advantage of this axial symmetry and align the laser polarization along the y axis. Under these conditions, the measured p_1y distribution is along the laser polarization, while p_1x is transverse. Moreover, the 'lost' p_1z distribution can be recovered from the measured p_1x distribution by taking advantage of the axial symmetry about the laser polarization. Under ideal conditions, this measurement is sufficient to retrieve the complete 3D momentum distributions of the fragments. In many cases, however, imperfections like spatial non-uniformities in the detector response may bias the results.

To circumvent this issue and verify that the momentum distribution is isotropic, we performed two measurements with the polarization along the z and y directions, as illustrated in figures 5(a) and (b), respectively. Note that while the angular distributions drawn in this figure do not resemble the isotropic distribution we measure in the present experiment, they help to better convey the difference between the two measurement schemes.

Zoom In Zoom Out Reset image size

In the first measurement with the polarization along z, the momentum distribution along the polarization cannot be retrieved because ${v}_{1z}^{{\prime} }$ has been neglected. In the second measurement with the polarization along y, the momentum distribution parallel to the laser polarization can be retrieved directly, while the complete transverse momentum distribution can be recovered by using the axial symmetry about the laser field as discussed above.

Let us define φ as a rotation angle in the xy-plane in both measurements. This angle is sketched on figure 4(a). We denote N(${\varphi }_{1}$) as the distribution obtained in the first measurement (polarization along z) and N(${\varphi }_{2}$) as the distribution found in the second measurement (polarization along y). Computing the ratio of these distributions yields the result shown in figure 5(c). Note that the position-dependent detection efficiency cancels out in this ratio, thus eliminating the impact of detector imperfections.

The ratio shown in figure 5(c) is rather flat, directly demonstrating that dissociation in flight yields an isotropic momentum distribution. This distribution is likely the result of t_d being much longer than the rotational timescale of the molecule. Thus, information about the initial alignment of the molecule with respect to the laser polarization is lost, and the resulting distribution is isotropic.

3.3. Kinetic energy release

Finally, just as accessing the z-component momenta is problematic, so too is retrieving the KER on an event-by-event basis. To obtain a KER distribution, we utilize a method based on the onion peeling technique, which has been widely used to analyze photofragment images [57–60]. The transverse momentum distribution, which was obtained using the method described in the previous section, serves as our projected 'onion,' which we slice along the p_x direction. Some sample slices for dissociation in flight of C₂D${}_{4}^{2+}$ are shown in figure 6(a). As we carry out iterative onion peeling subtraction on the slices, the counts that are 'peeled' away are allocated into the appropriate KER bins to accumulate a distribution. For a given iteration being performed on a particular slice, the KER is found using the edges of the slice, as only events with the maximum KER reach this region. We employ this technique for measurements in which the laser polarization is along the y axis, even though the isotropic nature of the distribution would allow one to use data from either of the previously discussed measurement schemes. Further details about this KER retrieval method will be described in a forthcoming publication about the dissociation in flight of metastable carbon dioxide dications [61].

Zoom In Zoom Out Reset image size

The KER distribution obtained using our 'sliced' onion peeling technique for dissociation in flight of C₂D${}_{4}^{2+}$ is presented in figure 6(b). We estimate the uncertainty in the obtained KER to be about 0.3 eV² . The centroid of the KER, at about 4.2 eV, is in good agreement with that obtained by Larimian et al [47]. Finally, as an alternative method of retrieving the KER, we compare the measured $N({p}_{1r})$ with several simulated $N({p}_{1r})$ distributions corresponding to Gaussian KER distributions with different centroids and widths, as illustrated in figure 6(c). As can be seen, the simulated $N({p}_{1r})$ distribution with a KER centroid of 4.15 eV and a FWHM of 0.5 eV agrees well with the measured $N({p}_{1r})$. Furthermore, this result supports that obtained using the 'sliced' onion peeling approach. As mentioned, the KER can supply a great deal of insight into dissociation pathways, e.g., [17, 25, 26]. This pursuit, however, is beyond the scope of this paper, which focuses on the method.

In summary, we have developed a method to study dissociation in flight of metastable molecular ions using coincidence momentum imaging measurements. Our approach, which supplies valuable information about the relevant metastable states, including the lifetime(s) and momentum distributions of the dissociating fragments, has been realized through the application and symmetries of the relevant kinematic equations.

Encountered hurdles such as the inverse problem of choosing ${t}_{d}^{(+)}$ or ${t}_{d}^{(-)}$ and the related problem of retrieving ${v}_{1z}^{{\prime} }$ have been addressed by exploiting symmetries of the fragmentation. The readily expressed forward-backward symmetry breaking in the ${v}_{1z}^{{\prime} }$ distribution was used in a self-consistent manner to obtain the N(t_d) distribution of the sample and hence the lifetimes of the metastable states dissociating in flight. This analysis allowed us to find two lifetimes in the deprotonation of metastable ethylene dications and a possible isotopic effect in the shorter lifetime.

The necessity of a single t_d value for each event to obtain the momenta was fulfilled by neglecting ${v}_{1z}^{{\prime} }$, an assumption which we have shown to be on solid ground, as it introduces minimal error in the calculation of t_d and the momenta. Furthermore, the azimuthal symmetry about the laser polarization was exploited to obtain all the components of the momentum rendered unretrievable by the t_d inverse problem.

Finally, while we have demonstrated this method for the specific case of deprotonation of metastable ethylene dications formed by intense femtosecond laser pulses, this technique is applicable to coincidence measurements on a variety of metastable molecular systems dissociating in flight, which could also be formed via other means, such as x-ray photoabsorption or fast charged particle impact.

We gratefully acknowledge Kanaka Raju P and C W Fehrenbach for assistance with the laser beam. This work was supported by the Chemical Sciences, Geosciences, and Biosciences Division, Office of Basic Energy Sciences, Office of Science, US Department of Energy.

The long, curved dissociation-in-flight stripe in the CTOF spectrum fairly closely follows a third-order polynomial dependence as a function of t₁. To more effectively select true events while suppressing the contribution of random pairs (i.e., improve the signal to 'noise' ratio), we use the coefficients of a third-order polynomial fit to the curved coincidence stripe to straighten it and then apply a simple rectangular gate. More specifically, t_⊥, the perpendicular distance from each (t₁, t₂) data point to the polynomial fit is calculated. A Gaussian function is fit to the $N({t}_{\perp })$ distribution, and a ±3σ gate is applied to the data, as shown in figure A1(a). Figure A1(b) shows the straightened stripe and gate around it. Note that for C₂H₄, the left and right gating bounds for t₁ were chosen to be 825 ns and 2400 ns, respectively, and for C₂D₄, they were 1200 ns and 2400 ns, respectively. These data selection schemes allow one to avoid the prompt breakup region, as well as the area of the CTOF spectrum near the end of the dissociation-in-flight stripe, where 'noise' and other channels dominate.

Zoom In Zoom Out Reset image size

Recall that due to the inverse problem related to the sign of ${v}_{1z}^{{\prime} }$ (and associated sign choice in the t_d equation), we employ an iterative procedure that starts with a choice of ${t}_{d}^{(+)}$ or ${t}_{d}^{(-)}$ that is equally likely. Then, we use the values of τ in subsequent iterations to weight the choice of ${t}_{d}^{(+)}$ or ${t}_{d}^{(-)}$, using equation (7). As this method utilizes a random number generator to choose ${t}_{d}^{(+)}$ or ${t}_{d}^{(-)}$, it is imperative that the sample size be large enough to ensure true randomness.

As we are applying this analysis to experimental data of a relatively weak channel for which statistics are limited, one may worry about the robustness of the proposed method that relies on random number generation. To address this issue, as mentioned in section 3.1, each iteration of the lifetime determination procedure consists of multiple trials. That is, for each iteration, the analysis is simply repeated multiple times (each time with a randomly selected seed). Each trial uses the same solution choice weighting scheme. An exponential decay function is fit to the resulting $N({t}_{d})$ distributions to retrieve lifetimes for each trial. At the end of an iteration, the amplitudes ${N}_{0}^{(1)}$ and ${N}_{0}^{(2)}$ and lifetimes ${\tau }_{1}$ and ${\tau }_{2}$ used to calculate the weighting factor in equation (7) are the average values from all the trials in that iteration.

As mentioned, when applied to our ethylene data, the lifetimes and amplitudes converge within just a few iterations. Plots illustrating the convergence of these quantities are shown in figures B1(a)–(d). Note that in the first iteration, ${t}_{d}^{(+)}$ and ${t}_{d}^{(-)}$ are selected with equal probability, and each iteration consists of 1000 trials.

Zoom In Zoom Out Reset image size

To obtain a single value of t_d needed for evaluating the momentum of each event, an approximation is necessary. We start with the first expression in equations (4) and divide both sides by ${t}_{{jm}}-{t}_{{dm}}$ to obtain

$$\begin{eqnarray}&&\displaystyle \frac{1-{t}_{{dm}}^{2}}{{t}_{{jm}}-{t}_{{dm}}}=2{t}_{{dm}}+2\displaystyle \frac{{v}_{{jz}}^{{\prime} }}{{v}_{m}}+{\eta }_{j}({t}_{{jm}}-{t}_{{dm}})\qquad [j=1,2].\end{eqnarray} \tag{ 9 }$$

This equation is solved for t_dm after neglecting the second term on the right-hand side that contains ${v}_{{jz}}^{{\prime} }$. How valid is neglecting this term? To explore this question, we perform simulations using a few typical input random functional distributions, including a single-term exponential decay $N({t}_{d})$ distribution with $\tau =900$ ns, a Gaussian KER distribution centered at 4 eV with a 1 eV FWHM, and an isotropic angular distribution. Let us denote the terms on the right-hand side of equation (9) as '1st term,' '2nd term,' and '3rd term,' in left-to-right order. We evaluate the quantity R ≡∣2nd term∣/(1st term + 3rd term), shown in figure C1(a). For 99% of the events, $R\lt 10$%. Furthermore, the error that neglecting the second term introduces into the recovered t_d and p_1x, shown in figures C1(b) and (c), is also reasonably small. For t_d, about 77% of the simulated events lie below the 5% error level, and for p_1x, about 98% of the events lie below the same error level.

Zoom In Zoom Out Reset image size

Since the validity of neglecting ${v}_{1z}^{{\prime} }$ depends on the magnitude of ${v}_{1z}^{{\prime} }$, and for our simulations we have assumed values of this quantity approximately matching the measured ones, it is reasonable to explore how large ${v}_{1z}^{{\prime} }$ can be before the approximation breaks down. Thus, we performed simulations with a few larger KER values (and hence larger maximum values of ${v}_{1z}^{{\prime} }$). Even for a high KER of 20 eV (Gaussian distribution with 1 eV FWHM and the same lifetime and angular distribution as before), 99% of the simulated events have $R\lt 20$%, as shown in figure C2(a). Moreover, as shown in figures C2(b) and (c), 75% of the events have <20% error in the retrieved t_d, and about 98% of the events are below the same error level in p_1x. Given this extreme example, we are thus assured that neglecting ${v}_{1z}^{{\prime} }$ for momentum computation is a reasonable approximation for our case.

Zoom In Zoom Out Reset image size

It is noteworthy that in general it is not merely the KER (and hence the maximum ${v}_{1z}^{{\prime} }$) that is important for consideration but the ratio ${v}_{{jz}}^{{\prime} }/{v}_{m}$, as can be readily seen in equation (9). Recall that v_m is the velocity of the dication. Therefore, in certain cases, it may also be desirable to increase the spectrometer voltage in the experiment to increase v_m and thereby improve the validity of this approximation.

Due to the nature of our method for evaluating t_d, it is not straightforward to propagate errors in order to determine the uncertainty in the recovered lifetime(s). Here we demonstrate how simulations aid us in pinpointing important sources of error.

As measured TOF values are used to calculate t_d, a rather pertinent question is what effect the uncertainty in these quantities has on the retrieved lifetimes. Finite time resolution, number truncation by the time-to-digital converter (TDC) unit used (which has a 25 ps resolution), and uncertainty in the absolute time (i.e., the error associated with the exact time of the laser-molecule interaction) all influence the accuracy of the retrieved lifetime. While we have performed simulations to investigate each of these effects, for the sake of brevity, we use those exploring truncation to illustrate the impact that this type of uncertainty can have.

To examine the effects of truncation of the measured TOF values, we simulate events with the same parameters as in the previous appendix (4 eV KER centroid). Thus, each simulated event has an associated t_d and ${v}_{1z}^{{\prime} }$. The $[j=1,2]$ equations in (4) are rearranged into quadratic expressions that can be solved for t₁ and t₂:

$$\begin{eqnarray}&&0={\eta }_{1}{t}_{1}^{2}+2\left((1-{\eta }_{1}){t}_{d}+\displaystyle \frac{{v}_{1z}^{{\prime} }{t}_{m}}{{v}_{m}}\right){t}_{1}+({\eta }_{1}-1){t}_{d}^{2}-\displaystyle \frac{2{v}_{1z}^{{\prime} }{t}_{d}{t}_{m}}{{v}_{m}}-{t}_{m}^{2},\end{eqnarray} \tag{ 10 }$$

$$\begin{eqnarray*}&&0={\eta }_{2}{t}_{2}^{2}+2\left((1-{\eta }_{2}){t}_{d}-\displaystyle \frac{\beta {v}_{1z}^{{\prime} }{t}_{m}}{{v}_{m}}\right){t}_{2}+({\eta }_{2}-1){t}_{d}^{2}+\displaystyle \frac{2\beta {v}_{1z}^{{\prime} }{t}_{d}{t}_{m}}{{v}_{m}}-{t}_{m}^{2}.\end{eqnarray*}$$

In solving for t₁ and t₂, we choose the positive root in the quadratic formula because choosing the negative root yields negative t₁ and t₂ values or ${t}_{1}\lt {t}_{d}$, which is unacceptable, in contrast to the t_d equation for which both roots can make physical sense. The values of t₁ and t₂ are truncated to varying levels of precision in ns, simulating the possible digitizer accuracy. Then a random fraction is added to the truncated number, as is done for the measured (digitized) data, to match the original number of digits. For example, if t_i is a time-of-flight value truncated to n decimal place(s), the new TOF after adding the random fraction is ${t}_{i}^{\prime} ={t}_{i}+r\times {10}^{-n}$, where r is a random number between 0 and 1. The new values of t₁ and t₂ are then used to compute t_d and obtain a lifetime, which is then compared to the original input lifetime. A few examples of the effect of truncation on the recovered lifetime are shown in figures D1(a)–(c). As can be seen in the figure, as more of the true digits are initially truncated, the retrieved lifetime veers away from the original lifetime of 900 ns, and the range of $N({t}_{d})$ that can be used for fitting starts to deteriorate. As mentioned, similar tests were performed to examine the effects of finite resolution and absolute time uncertainty. These simulations yielded results comparable to the number truncation tests.

Zoom In Zoom Out Reset image size

While our simulations have proven instructive in identifying influential error sources such as those related to the measured time-of-flight values, note that the errors presented in the body of this manuscript reflect those evaluated by statistical means.