# Bayesian inferencing and deterministic anisotropy for the retrieval of the
molecular geometry $|\Psi(\mathbf{r})|^{2}$ in gas-phase diffraction
experiments

Kareem Hegazy1,2,∗ Varun Makhija3 Phil Bucksbaum1,2,4 Jeff Corbett5 James
Cryan2 Nick Hartmann6 Markus Ilchen2,7,8 Keith Jobe5 Renkai Li9 Igor Makasyuk5
Xiaozhe Shen5 Xijie Wang5 Stephen Weathersby5 Jie Yang9 Ryan Coffee2,6,∗
1Department of Physics, Stanford University, Stanford, California 94305, USA
2Stanford PULSE Institute, SLAC National Accelerator Laboratory, 2575 Sand
Hill Road, Menlo Park, California 94025, USA  3Department of Chemistry and
Physics, University of Mary Washington, Fredericksburg, Virginia 22401, USA
4Department of Applied Physics, Stanford University, Stanford, California
94305, USA  5SLAC National Accelerator Laboratory, 2575 Sand Hill Road, Menlo
Park, California 94025, USA 6Linac Coherent Light Source, SLAC National
Accelerator Laboratory, Menlo Park, California 94025, USA 7European XFEL,
Holzkoppel 4, 22869 Schenefeld, Germany  8Universität Kassel, Heinrich-Plett-
Str. 40, 34132 Kassel, Germany  9Department of Engineering Physics, Tsinghua
University, Beijing 100084, China

###### Abstract

Currently, our general approach to retrieve the molecular geometry from
ultrafast gas-phase diffraction heavily relies on complex geometric
simulations to make conclusive interpretations. In this manuscript, we develop
a broadly applicable ultrafast gas-phase diffraction method that approximates
the molecular frame geometry $|\Psi(\mathbf{r},t)|^{2}$ distribution using
Bayesian Inferencing. This method does not require complex molecular dynamics
simulation and can identify the unique molecular structure. We demonstrate
this method’s viability by retrieving the ground state geometry distribution
$|\Psi(\mathbf{r})|^{2}$ for both simulated stretched NO2 and measured ground
state N2O. Due to our statistical interpretation, we retrieve a coordinate-
space resolution on the order of 100 fm, depending on signal quality, an
improvement of order 100 compared to commonly used Fourier transform based
methods. By directly measuring the width of $|\Psi(\mathbf{r})|^{2}$, we open
ultrafast gas-phase diffraction capabilities to measurements beyond current
analysis approaches since this width is generally only accessible through
simulation. Our method also leverages deterministic ensemble anisotropy; this
provides an explicit dependence on the molecular frame angles. This method’s
ability to retrieve the unique molecular structure with high resolution, and
without complex simulations, provides the potential to effectively turn gas-
phase ultrafast diffraction into a discovery-oriented technique, one that
probes systems that are prohibitively difficult to simulate.

##  I Introduction

Ultrafast molecular gas-phase diffraction Williamson  _et al._ (1997); Ewbank
_et al._ (1992, 1993); Shen  _et al._ (2019); Stankus  _et al._ (2020);
Minitti  _et al._ (2015) is a vital tool for retrieving time dependent
molecular structure. We aim to expand this tool to perform high precision
reconstructions of molecular geometry $|\Psi(\mathbf{r},t)|^{2}$, where
$\mathbf{r}$ is the molecular frame nuclear coordinates, without relying on
molecular dynamics simulations. Since elastic diffraction is a many-to-one
process, directly inverting diffraction patterns for the molecular geometry is
intractable. Typically, we require complex excited state simulations to
provide the molecular geometries which we then validate through comparisons
with measured diffraction patterns or pair-distribution functions (PDFs – a
weighted histogram of pair-wise distances). Consequently, ultrafast gas-phase
diffraction is generally constrained by the limitations of complex excited
state simulations.

Figure 1: Ensemble anisotropy reveals the pair-wise angles
($\theta_{\mu\nu}^{(mf)}$ and $\phi_{\mu\nu}^{(mf)}$) and is measured with
respect to the lab frame Euler angles $\theta^{(\text{lf})}_{\text{I}}$,
$\phi^{(\text{lf})}_{\text{I}}$, and $\chi^{(\text{lf})}_{\text{I}}$. We
illustrate the correspondence between the lab frame (a) and molecular frame
(b) for a single pair-wise distance (the highlighted NO bond). The origin is
defined by the nitrogen atom. The lab frame is defined by the laser
polarization and propagation direction, while the molecular frame is defined
by the molecule’s ground state principal moments of inertia.  Figure 2: The
Axis Distribution Moments (ADMs) encapsulate the ensemble anisotropy which
provides various constraints on the molecular frame as a function of time.
Panel a shows the square norms of the ADMs, while panels b and c show their
time dependence. Panels d and e show the time dependent ensemble anisotropy
probability distribution for $\theta^{(\text{lf})}_{\text{I}}$ and
$\chi^{(\text{lf})}_{\text{I}}$, respectively. Panels f and g show
illustrative line-outs of these Euler angle distributions for
$\theta^{(\text{lf})}_{\text{I}}$ and $\chi^{(\text{lf})}_{\text{I}}$,
respectively, with isotropy indicated by the dashed lines.

A variety of studies sought to reduce reliance on complex simulations. Fourier
transforming the time dependence exposes dissociative and vibronic signals
Ware  _et al._ (2019a, b); Bucksbaum  _et al._ (2020) but it is insensitive to
classes of isomerizations. Methods employing ensemble anisotropy have garnered
much interest Saldin  _et al._ (2010); Starodub  _et al._ (2010); Ho  _et al._
(2009); Elser (2011); Poon  _et al._ (2013); Reckenthaeler  _et al._ (2009);
Pabst  _et al._ (2010a); Yang  _et al._ (2014); Hensley  _et al._ (2012) yet
they struggle to get sub-Angstrom resolution for generic molecular structures.
Optimization methods, while capable of exposing large-scale motion, are
susceptible to local minima Yang  _et al._ (2014). Pattern matching measured
data against sampled isomers Stankus  _et al._ (2019); Natan (2021); Natan
_et al._ (2021) becomes intractable for moderately large geometries due to the
curse of dimensionality. For example, a molecule with $N_{\text{atoms}}$ atoms
has $3N_{\text{atoms}}-6$ degrees of freedom. To independently sample each
degree of freedom 10 times would require $10^{3N_{\text{atoms}}-6}$
geometries, becoming intractable for molecules with 7 or more atoms.
Simulations reduce the geometry-space of isomers to select, but this trade-off
requires previous knowledge Stankus  _et al._ (2019) that potentially imparts
biases.

We employ insights from molecular ensemble anisotropy methods, applied
statistics, and machine learning principles to approximate the molecular
geometry probability density $(|\Psi(\mathbf{r},t)|^{2})$. We access the
molecular frame by decomposing measured data onto anisotropic components.
Then, we iteratively approximate $|\Psi(\mathbf{r},t)|^{2}$ with a statistical
approach uniquely suited for high repetition-rate diffraction facilities where
the resolution strongly improves with signal to noise much faster than
increasing the $q$ range beyond moderate values. Unlike the PDF approach, it
retrieves the molecular distances and angles required to define a unique
molecular geometry. An important new feature of our method is its ability to
quantitatively measure additional parameters, such as the width of
$|\Psi(\mathbf{r},t)|^{2}$. We report order 100 fm resolution without complex
excited state molecular dynamics simulations and require only the ground state
geometry and its transition dipole moment. This method has the potential to
expand ultrafast gas-phase diffraction into a discovery oriented technique,
one that is free of complex excited state simulation limitations and is
applicable to novel molecular systems.

Time dependent ensemble anisotropy accesses the molecular frame through
angular constraints Makhija  _et al._ (2016); Marceau  _et al._ (2017);
Gregory  _et al._ (2021); Sándor  _et al._ (2018a); Mikosch  _et al._ (2013);
Sándor  _et al._ (2019, 2018b). We define anisotropy with the Euler angles
$\theta^{(\text{lf})}_{\text{I}}$ (polar), $\phi^{(\text{lf})}_{\text{I}}$
(azimuthal), and $\chi^{(\text{lf})}_{\text{I}}$ between the lab frame and the
molecular frame (Fig. 1a). An induced rotational wavepacket creates ensemble
anisotropy given by $|\Psi(\theta^{(\text{lf})},\phi^{(\text{lf})},t)|^{2}$.
Axis Distribution Moments (ADMs) Stolow and Underwood (2008); Underwood and
Reid (2000); Gregory  _et al._ (2021) are the coefficients in the Wigner D
matrix expansion of $|\Psi(\theta^{(\text{lf})},\phi^{(\text{lf})},t)|^{2}$

$\mathcal{A}^{l}_{mk}(t)=\frac{2l+1}{8\pi^{2}}\expectationvalue{D^{l}_{mk}\left(\phi_{\text{I}}^{(\text{lf})},\theta_{\text{I}}^{(\text{lf})},\chi_{\text{I}}^{(\text{lf})}\right)}{\Psi(t)}.$
(1)

These ADMs describe the ensemble of molecular frame orientations with respect
to the lab frame. When calculating the ADMs the $l$, $m$, and $k$ are
difference and sum of quantum numbers between rotational eigenstates,
respectively for the total angular momentum, the projection onto the lab frame
z axis, and the projection onto the molecular frame z axis. The ADMs access
the molecular frame by decomposing the measurement into $C_{lmk}(q)$
coefficients, dependent on pair-wise distances and angles
($\theta_{\mu\nu}^{(mf)}$ and $\phi_{\mu\nu}^{(mf)}$) shown in Fig. 1b. The
PDF is not directly sensitive to these angles. Figure 2 illustrates how
transient anisotropy (panels b and c) provides constraints on the these Euler
angles and consequently the molecular frame (panels d-g). For example, at
39.25 ps the anisotropy provides simultaneous constraints on
$\theta^{(\text{lf})}_{\text{I}}$ and $\chi^{(\text{lf})}_{\text{I}}$. At
39.68 ps, $\chi^{(\text{lf})}_{\text{I}}$ (the molecular frame azimuthal
plane) is highly constrained. At 39.85 ps the ensemble is well localized in
$\theta^{(\text{lf})}_{\text{I}}$, resolving measurements along the molecular
frame $\hat{\mathbf{z}}$. Here, P$\left(\phi^{(\text{lf})}_{\text{I}}\right)$
is uniform due to cylindrical symmetry imparted by a linearly polarized pulse.

To approximate $|\Psi(\mathbf{r},t)|^{2}$, we use Bayesian Inferencing with
markov-chain monte carlo (MCMC) techniques to tackle the curse of
dimensionality. Bayesian Inferencing describes a class of statistical
inferencing techniques using Bayes’s Theorem to update one’s model based on
observed data Box and Tiao (2011), allowing us to ultimately approximate
$|\Psi(\mathbf{r},t)|^{2}$ as the probability distribution
$P\left(\mathbf{r},t\right|\left.\bm{\Theta},C\right)$. The
$P\left(\mathbf{r},t\right|\left.\bm{\Theta},C\right)$ distribution is
parameterized by $\bm{\Theta}$, which includes the molecular geometry degrees
of freedom and is retrieved using MCMC techniques. This Bayesian Inferencing
and MCMC approach efficiently samples pertinent geometries consistent with the
measured $C_{lmk}(q)$ in an unbiased fashion. It naturally avoids regions in
our sampling space that are inconsistent with the $C_{lmk}(q)$. Importantly,
since $P\left(\mathbf{r},t\right|\left.\bm{\Theta},C\right)$ is conditioned on
the measured $C_{lmk}(q)$ coefficients and we analytically relate molecular
frame pair-wise distances and angles to the $C_{lmk}(q)$, we actually retrieve
$P\left(\mathbf{r},t\right|\left.\bm{\Theta},C\right)$ with neither the PDF
nor complex molecular dynamics simulations.

In this manuscript, we validate these principles by retrieving
$|\Psi(\mathbf{r})|^{2}$ for the ground states of both simulated NO2 and
measured N2O. Here NO2, an asymmetric top, serves as a test case to show our
method’s broad capabilities and behavior under various experimental
conditions, while the N2O experiment validates such capabilities on measured
results. We chose these molecules to be amenable to conventional methods to
benchmark. We discuss how $P\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$
improves upon the Fourier limited PDF resolution by a factor of 100, and how
this procedure depends more strongly on signal to noise than the measured $q$
range beyond modest thresholds.

##  II Method

Figure 3: To access the molecular geometry term, in the molecular frame, one
must remove the lab frame anisotropy dependence and fit onto the ADMs. For the
NO2 simulation (left) and N2O data (right), we illustrate the analysis steps.
One first measures $\Delta\langle I(\textbf{q})\rangle(t)$, the difference
diffraction pattern, given by Eq. 3 (top row). Removing the detector angular
dependence, one retrieves $B^{m}_{l}(q,t)$ of Eq. 4 (middle row). Which in
turn yields $M_{lmk}(q)$ coefficients Eq. 5 (bottom row) after removing the
time dependent ensemble anisotropy (ADMs). All as described in the text. We
note that in the N2O data (right) we have limited visibility of data due to
experimental limitations illustrated by the hashes.

Our method can be subdivided into three principal concepts. Firstly, we use
ensemble anisotropy, described by the ADMs, to access the molecular frame by
projecting the data onto anisotropic components. Secondly, we select a model,
$P\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$, to approximate
$|\Psi(\mathbf{r})|^{2}$ and develop our statistical approach to solve for
$\bm{\Theta}$ using Bayesian Inferencing. Lastly, we take our statistical
description and use MCMC techniques to solve for $P\left(\bm{\Theta}|C\right)$
to retrieve the optimal $\bm{\Theta}$ parameters ($\bm{\Theta}^{*}$). The code
used for this analysis can be found in Ref. github (2022), which can be run to
reproduce the following results or adapted for other molecules.

###  II.1 Extracting Molecular Frame Information

We describe our analysis procedure for a system given an induced deterministic
ensemble anisotropy under experimental conditions at the SLAC MeV Ultrafast
Electron Diffraction facility (UED) Shen  _et al._ (2019). Our generic pump-
probe setup is similar to most ultrafast diffraction setups, consisting of an
800 nm Ti:Sapphire pump laser and a 120 fs FWHM electron bunch probe. For the
simulated NO2 results, we consider using a single 10 TW/cm2 800 nm pump pulse
to impulsively induce a coherent rotational wave packet and probing it within
a window of high anisotropy variation: [37.5, 41.5] ps. For the measured N2O
sample, a train of 8 identical 800 nm pulses (40 fs duration and $5\times
10^{12}$ W/cm2 irradiance) separated by full quantum revivals induced such
rotational wavepacket. We measured the first field free full quantum revival
over a window of ~3 ps. We masked $q$ regions $[0,3.5]$ Å-1 and above 7.25 Å-1
due to ellipticity in the imaging of the diffraction pattern and poor signal
to noise, respectively. Linearly polarized pump pulses induce azimuthal
symmetry, which sets $m=0$ in Eq. 1
(P$(\phi_{\text{I}}^{(\text{lf})},t)=1/2\pi$), while the Raman excitation of
the wavepacket requires $l$ being even in Eq. 1.

We define anisotropy in two equivalent ways and quantify it through the ADMs.
Firstly, anisotropy is defined by a non-zero projection of the measured
diffraction pattern onto any $Y_{l}^{m}$ with even $l>0$ for a given $\Delta
q$ range. Secondly, anisotropy exists when there is a non-zero
$\mathcal{A}^{l}_{mk}(t)$ for $l>0$. We note other methodologies in Refs.
Stolow and Underwood (2008); Hockett (2015); Underwood and Reid (2000).
Supplementary Section VIII describes our calculation of the ADMs.

We access the molecular pair-wise distances and angles in the molecular frame.
Using the ADMs and the Independent Atom Approximation, we relate measured lab
frame anisotropy in diffraction patterns, $\langle I(\mathbf{q})\rangle(t)$,
to the molecular geometry

$\displaystyle\langle I(\textbf{q})$
$\displaystyle\rangle(t)={}\mathcal{I}\bigg{(}\sum_{\mu}|f^{*}_{\mu}(q)|^{2}+\sum_{\mu,\nu:\mu\neq\nu}\text{Re}\bigg{\\{}f_{\mu}(q)f^{*}_{\nu}(q)\sum_{l}i^{l}8\pi^{2}\sqrt{4\pi\left(2l+1\right)}$
(2) $\displaystyle\times\sum_{m,k}$
$\displaystyle(-1)^{k-m}\underbrace{Y^{-m}_{l}\left(\theta^{(\text{lf})}_{q},\phi^{(\text{lf})}_{q}\right)}_{\text{\small
Lab
Frame}}\expectationvalue{\underbrace{D^{l}_{mk}\left(\phi_{\text{I}}^{(\text{lf})},\theta_{\text{I}}^{(\text{lf})},\chi_{\text{I}}^{(\text{lf})}\right)}_{\text{\small
Ensemble Anisotropy}}\underbrace{j_{l}(q\Delta
r_{\mu\nu})Y_{l}^{-k}\left(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}\right)}_{\text{\small
Molecular Frame Geometry}}}{\Psi(t)}\bigg{\\}}\bigg{)}.$

In Eq. 2, derived in Supplementary Section IX, $f_{\mu}(q)$ is the scattering
amplitude of the $\mu^{\text{th}}$ atom, $j_{l}(qr)$ are the spherical Bessel
functions of the first kind, and the momentum transfer vector is given by
$\mathbf{q}=[q,\theta_{q}^{(\text{lf})},\phi_{q}^{(\text{lf})}]$. The
difference vector
$\Delta\mathbf{r}_{\mu\nu}=\mathbf{r}_{\mu}-\mathbf{r}_{\nu}=[\Delta
r_{\mu\nu},\theta^{(\text{mf})}_{\mu\nu},\phi^{(\text{mf})}_{\mu\nu}]$ are the
molecular frame pair-wise distances and angles between the $\mu^{\text{th}}$
and $\nu^{\text{th}}$ atoms. Figure 1b illustrates the molecular frame pair-
wise distances and angles. This equation concisely shows how the ensemble
anisotropy connects the lab frame to the molecular frame geometry. We note
that directly accessing the molecular frame pair-wise angles
$(\theta^{(\text{mf})}_{\mu\nu},\phi^{(\text{mf})}_{\mu\nu})$ requires
anisotropy and is otherwise inaccessible through the PDF and isotropic
contributions alone. This is evident by isolating the isotropic component
($l=0$, $m=0$, $k=0$) which sets
$Y_{0}^{0}\left(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}\right)=1/(2\sqrt{\pi})$.

For our method, we describe optimal representations of the lab and molecular
frames used in Eq. 2. The molecular frame is defined by the molecule’s
principal moments of inertia before laser excitation with the
$\hat{\mathbf{z}}^{(\text{mf})}$, $\hat{\mathbf{x}}^{(\text{mf})}$, and
$\hat{\mathbf{y}}^{(\text{mf})}$ corresponding to the principle moments of
inertia in decreasing order: A, B, and C respectively. When looking at the
$\Delta\mathbf{r}_{\mu\nu}$ contribution, we isolate the $\mu^{\text{th}}$ and
$\nu^{\text{th}}$ atoms while ignoring other atoms and translate the atom pair
such that $\mathbf{r}_{\nu}$ defines the origin. This is highlighted in Fig.
1b where the nitrogen is translated to the origin. This translation allows us
to define the pair-wise angles and derive Eq. 2. Since we are concerned with a
difference in locations $\Delta\mathbf{r}_{\mu\nu}$, Eq. 2 is invariant under
such molecular frame translations. In the lab frame, the laser polarization
defines $\hat{\mathbf{z}}^{(\text{lf})}$ and the propagation direction of the
probe pulse defines $\hat{\mathbf{y}}^{(\text{lf})}$. The measured signals in
the lab frame, on a 2D detector, are defined by detector parameters
$q=|\textbf{q}|$ and the azimuthal angle $\theta^{(\text{d})}$ defined by
$\hat{\mathbf{z}}^{(\text{lf})}$. Supplementary Section IX describes how to
rewrite $\mathbf{q}$ in terms of the detector coordinates, in particular for
small angle scattering at UED
$\theta^{(\text{lf})}_{q}\approx\theta^{(\text{d})}$ and
$\phi^{(\text{lf})}_{q}\approx 0$ as will be used below.

The primary difficulty of working with Eq. 2 comes from the expectation value
including both the ensemble anisotropy and molecular frame geometry. We want
to separate the ensemble anisotropy into the ADMs. This isolates the time
dependent molecular geometry term that we would like to retrieve. By doing
this, we only require more tractable molecular rotation simulations with
respect to the known ground state geometry in order to retrieve the time
dependent molecular geometry. Otherwise, as Eq. 2 is written, it requires a
priori knowledge of exactly the unknown time dependent geometries for which we
are solving. In this work, we describe various ways to do this under common
experimental conditions.

Focusing on the ground state of NO2, we can separate the ADMs and molecular
geometry contribution in Eq. 2 by applying a rigid rotor approximation:

$\displaystyle\langle I(\mathbf{q})$
$\displaystyle\rangle_{\text{rigid}}(t)={}\mathcal{I}\bigg{(}\sum_{\mu}|f_{\mu}(q)|^{2}+\sum_{\mu,\nu:\mu\neq\nu}\text{Re}\bigg{\\{}f_{\mu}(q)f^{*}_{\nu}(q)\sum_{l}i^{l}8\pi^{2}\sqrt{\frac{4\pi}{\left(2l+1\right)}}$
(3) $\displaystyle\times\sum_{m,k}$
$\displaystyle(-1)^{k-m}\underbrace{Y^{-m}_{l}\left(\theta^{(\text{lf})}_{q},\phi^{(\text{lf})}_{q}\right)}_{\text{\small
Lab Frame}}\expectationvalue{\underbrace{j_{l}(q\Delta
r_{\mu\nu})Y_{l}^{-k}\left(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}\right)}_{\text{\small
Molecular Frame
Geometry}}}{\Psi(0)}\underbrace{\left.\mathcal{A}^{l}_{mk}(t)\right|_{\text{rigid}}}_{\text{\small
Anisotropy}}\bigg{\\}}\bigg{)}.$

Equation 3 is the general form, which we adapt to our specific case by setting
$m=0$ and replacing $\theta^{(\text{lf})}_{q}\approx\theta^{(\text{d})}$ and
$\phi^{(\text{lf})}_{q}\approx 0$. The resulting lab frame measurements are
shown in Fig. 3 (top row).

Figure 4: For simulated NO2 we defined a $|\Psi(\mathbf{r})|^{2}$
distribution, from which we calculated the $C_{lmk}(q)$ under various
experimental conditions. Panel a shows the simulated NO2 distribution that we
use to calculate the simulated NO2 responses ($C_{lmk}(q)$ and $M_{lmk}(q)$).
Panel b shows $M_{lmk}(q)$ for various SNR ratios for the case of an ensemble
temperature of 100 K and kick fluence of 1 J/cm2. Panel c shows the ADMs’
dependence on pump strength (constant ensemble temperature of 25 K) and
temperature (constant pump fluence of 1 J/cm2).

To retrieve $P\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$, we first
isolate the molecular frame geometry terms from Eq. 3 with a series of fits.
The first projects out the measured lab frame anisotropy
$\left(Y^{-m}_{l}\left(\theta^{(\text{lf})}_{q},\phi^{(\text{lf})}_{q}\right)\right)$
from Eq. 3 by fitting the angular dependence of the measured diffraction.

$\displaystyle B^{m}_{l}(q,t)={}$
$\displaystyle\int_{0}^{\pi}\Big{\langle}I\Big{(}\mathbf{q}(q,\theta^{(\text{d})}),t\Big{)}\Big{\rangle}_{\text{rigid}}Y_{l}^{m}\Big{(}\theta^{(\text{lf})}_{q}\left(q,\theta^{(\text{d})}\right),\phi^{(\text{lf})}_{q}\left(q,\theta^{(\text{d})}\right)\Big{)}\sin\Big(\theta^{(\text{lf})}_{q}\left(q,\theta^{(\text{d})}\right)\Big{missing})d\theta^{(\text{d})}$
$\displaystyle\begin{split}={}&\mathcal{I}\sum_{\mu,\nu:\mu\neq\nu}\text{Re}\bigg{\\{}f_{\mu}(q)f^{*}_{\nu}(q)i^{l}8\pi^{2}\sqrt{\frac{4\pi}{\left(2l+1\right)}}\\\
&\hskip
28.45274pt\times\sum_{k}(-1)^{k-m}\expectationvalue{\underbrace{j_{l}(q\Delta
r_{\mu\nu})Y_{l}^{-k}\left(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}\right)}_{\text{\small
Molecular Frame
Geometry}}}{\Psi(0)}\underbrace{\left.\mathcal{A}^{l}_{mk}(t)\right|_{\text{rigid}}}_{\text{\small
Anisotropy}}\bigg{\\}}\end{split}$ (4)

This yields the time $(t)$ and $q$ dependent $B^{m}_{l}(q,t)$ coefficients
shown in Fig. 3 (middle row). The second fit isolates the molecular frame
information through by fitting out the time dependence of $B^{m}_{l}(q,t)$
with the simulated ADMs, $\mathcal{A}^{l}_{mk}(t)$. The resulting
coefficients, $C_{lmk}(q)$, concisely relate measured data to the molecular
frame pair-wise geometry.

$\displaystyle\begin{split}C_{lmk}(q)={}&\mathcal{I}\sum_{\mu,\nu:\mu\neq\nu}\text{Re}\bigg{\\{}f_{\mu}(q)f^{*}_{\nu}(q)(-1)^{k-m}i^{l}8\pi^{2}\\\
\times\sqrt{\frac{4\pi}{\left(2l+1\right)}}&\expectationvalue{\underbrace{j_{l}(q\Delta
r_{\mu\nu})Y_{l}^{-k}\left(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}\right)}_{\text{\small
Molecular Frame Geometry}}}{\Psi(0)}\bigg{\\}}\end{split}$ (5) $\displaystyle
M_{lmk}(q)={}$ $\displaystyle\frac{C_{lmk}(q)}{\sum_{\mu}|f_{\mu}(q)|^{2}}.$
(6)

Here, $M_{lmk}(q)$ are the modified $C_{lmk}(q)$ coefficients that compensate
for the rapid $q^{-4}$ falloff in the electron scattering amplitudes. For the
N2O data, the poor signal to noise precludes all contributions except
$C_{200}(q)$. Figure 3 (bottom row) shows the retrieved $M_{lmk}(q)$ for both
the simulated and measured data. Depending on the data quality and degree of
orthogonality in the ADMs, one may need to employ regularization to return
physical values. Regularization adds a fitting cost to extraneous
coefficients, thus minimize the impact of non-orthogonal ADMs. Supplementary
Section X provides further discussion on fitting the ADMs and regularization.

###  II.2 Applying Bayesian Inferencing

Having isolated the molecular frame terms ($C_{lmk}(q)$), we approximate
$|\Psi\left(\mathbf{r}\right)|^{2}$ with a chosen model and use Bayesian
Inferencing Box and Tiao (2011); Foreman-Mackey  _et al._ (2013) to adapt said
model to the observed data. We again emphasize that $\mathbf{r}$ represents
the molecular frame coordinates of the nuclei, unless otherwise stated. We
approximate $|\Psi(\mathbf{r})|^{2}$ with the probability distribution
$P\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$, which is parameterized
by $\bm{\Theta}$ and conditioned (optimized) on the observed $C_{lmk}(q)$
coefficients. This necessitates two key steps in the analysis. Firstly, one
must choose a functional form of
$P\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$ dependent on the system’s
state and the desired degree of accuracy. For example, one may choose a
multidimensional delta function for a single molecule response, a normal
distribution to model the ground state, or Hermite polynomials to describe
vibronic states. Secondly, given the presumed functional form of
$P\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$, one will find the
globally optimal $\bm{\Theta}$ parameters ($\bm{\Theta}^{*}$) that best
describes the measured $C_{lmk}(q)$ coefficients. Since we are retrieving the
ground state geometry we ultimately, want the normal distribution. The delta
distribution is analogous to calculating the diffraction pattern from a single
geometry and comparing it to one’s measurement. This is important for
intermediate results and in some cases may be the only tractable solution.
Here we explicitly write out
$P\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$ and the $\bm{\Theta}$
parameters

$\displaystyle P\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$
$\displaystyle\approx\left|\Psi\left(\mathbf{r}\right)\right|^{2}$ (7)
$\displaystyle P^{(\delta)}\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$
$\displaystyle=\delta\left(\bm{\Theta}^{(\delta)}-\mathbf{r}\right)$ (8)
$\displaystyle\bm{\Theta}^{(\delta)}$
$\displaystyle=\left[\expectationvalue{\text{NO}^{(1)}},\expectationvalue{\text{NO}^{(2)}},\expectationvalue{\angle\text{ONO}}\right]$
(9) $\displaystyle
P^{(\mathcal{N})}\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$
$\displaystyle=\frac{1}{\sqrt{2\pi}^{N_{\text{dof}}}\prod^{i<N_{\text{dof}}}_{i=0}\bm{\Theta}^{(\mathcal{N})}_{2i+1}}\exp\\{\frac{-1}{2}\sum_{i=0}^{i<N_{\text{dof}}}\left(\frac{\bm{\Theta}^{(\mathcal{N})}_{2i}-\mathbf{r}_{i}}{\bm{\Theta}^{(\mathcal{N})}_{2i+1}}\right)^{2}\Bigg{\\}}$
(10) $\displaystyle\bm{\Theta}^{(\mathcal{N})}$
$\displaystyle=\left[\expectationvalue{\text{NO}^{(1)}},\sigma\left(\text{NO}^{(1)}\right),\expectationvalue{\text{NO}^{(2)}},\sigma\left(\text{NO}^{(2)}\right),\expectationvalue{\angle\text{ONO}},\sigma\left(\angle\text{ONO}\right)\right]$
(11)

The $\bm{\Theta}$ parameters include the $3N_{\text{atom}}-6$ geometric
degrees of freedom $(N_{\text{dof}})$ needed to define the molecular geometry,
and the width parameters in the case of
$P^{(\mathcal{N})}\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$. It is
important to note that $\bm{\Theta}$ has the minimal number of parameters
needed to define $P\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$. Adding
redundant parameters can significantly alter $P\left(\bm{\Theta}|C\right)$.

With our chosen $P\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$ and
measured $C_{lmk}(q)$ we now retrieve $\bm{\Theta}^{*}$, which is the mode of
the posterior distribution $P\left(\bm{\Theta}|C\right)$. It is required to
build $P\left(\bm{\Theta}|C\right)$ and find $\bm{\Theta}^{*}$ in the full
$\bm{\Theta}$-space without marginalizing over individual parameters. This is
because the mean and mode of $P\left(\bm{\Theta}|C\right)$ will disagree due
to correlations between parameters and $P\left(\bm{\Theta}|C\right)$.
Consequently, we must next tackle the curse of dimensionality.

Figure 5: We successfully retrieve the multidimensional posterior
$P^{(\mathcal{N})}\left(\bm{\Theta}|C\right)$ for NO2 and N2O from which we
find $\bm{\Theta}^{*}$. Panel a shows the 1d and 2d projections of the
retrieved $P^{(\mathcal{N})}\left(\bm{\Theta}|C\right)$ distributions for the
simulated NO2 response. The recovered
$P^{(\mathcal{N})}(\mathbf{r}|\bm{\Theta}^{*},C)$ (panel b) is what we compare
to Fig. 4a. The red dashed lines indicate $\bm{\Theta}^{*}$, while the black
“x” and solid black lines indicate the ground truth, respectively. Panel c
shows the 1d and 2d projections of the retrieved
$P^{(\mathcal{N})}\left(\bm{\Theta}|C\right)$ distributions for N2O data,
though only using the $C_{200}(q)$ contribution. The black “x” and solid black
lines indicate previously measured values for N2O Herzberg (1966). For
comparison, panel d shows the simulated PDF from the same $q$ range.

###  II.3 Solving for the high dimensional model parameters $\bm{\Theta}$

We retrieve $P\left(\bm{\Theta}|C\right)$ with the Metropolis Hastings
algorithm (MHA) from the following system of equations:

$\displaystyle C_{lmk}(q)$ $\displaystyle=\int
H_{lmk}\left(q,\textbf{r}\right)\left|\Psi\left(\mathbf{r}\right)\right|^{2}d\textbf{r}$
(12) $\displaystyle C_{lmk}^{(\text{calc})}(q,\bm{\Theta})$
$\displaystyle=\int
H_{lmk}\left(q,\textbf{r}\right)P\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)d\textbf{r}$
(13)
$\displaystyle\begin{split}H_{lmk}\left(q,\textbf{r}\right)&=\mathcal{I}\text{Re}\bigg{\\{}i^{l}(-1)^{k-m}8\pi^{2}\sqrt{\frac{4\pi}{\left(2l+1\right)}}\\\
\times\sum_{\mu,\nu:\mu\neq\nu}|f_{\mu}(q)||&f_{\nu}(q)|j_{l}(q\Delta
r_{\mu\nu})Y_{l}^{-k}\left(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}\right)\bigg{\\}}.\end{split}$
(14)

We note the high dimensionality and complexity of Eq. 13, which is a system of
order 10 equations, each with order 100 terms, embedded in an order
100-dimensional space of $q$ bins. This must be evaluated on a
$N_{\text{dof}}$-dimensional space of all possible molecular geometries and
width parameters. The MHA is chosen for its ability to retrieve probability
distributions from high dimensional integral equations like Eq. 13 Hastings
(1970); Foreman-Mackey  _et al._ (2013).

The MHA is designed to efficiently and preferentially sample regions of
$\bm{\Theta}$-space proportional to the agreement with data, spending the vast
majority of its time sampling regions of high probability. The MHA builds
$P\left(\bm{\Theta}|C\right)$ by accumulating $\bm{\Theta}$ parameters based
on comparing their likelihood to neighboring $\bm{\Theta}^{{}^{\prime}}$
parameters. For instance, if the likelihood of $\bm{\Theta}$ is 1000 times
larger than $\bm{\Theta}^{\prime}$, the MHA will visit $\bm{\Theta}^{\prime}$
once for every 1000 visits to $\bm{\Theta}$. This likelihood is the
probability of observing the data $C_{lmk}(q)$ with a corresponding standard
error of the mean $\sigma_{lmk}(q)$, given the parameters $\bm{\Theta}$:

$\begin{split}&P\left(C|\bm{\Theta}\right)=\left[\prod_{lmk,q}\frac{1}{\sigma_{lmk}(q)\sqrt{2\pi}}\right]e^{K\left(\bm{\Theta}\right)}\\\
&\times\exp{\frac{-1}{2}\sum_{lmk,q}\left(\frac{C_{lmk}(q)-C^{(\text{calc})}_{lmk}(q,\bm{\Theta})}{\sigma_{lmk}(q)}\right)^{2}}.\end{split}$
(15)

Here, $K(\bm{\Theta})$ is used to constrain parameters to physicality e.g.,
$\bm{\Theta}>0\text{ and }\angle\text{ONO}<\pi$. We note the likelihood
function, and hence the MHA, is theory independent and is analogous to a
random walk guided by the relative agreement of neighboring $\bm{\Theta}$
parameters to the data. Reference Foreman-Mackey  _et al._ (2013) gives a
detailed description of combining Bayesian Inferencing and the MHA, as well as
the MHA python package used in this work. Supplementary Section XI describes
our use of the MHA and Bayesian Inferencing in greater detail and how one can
introduce physical intuition, or a priori knowledge, into the MHA.

The measured $q$ range, the induced rotational wavepacket, and the
$\sigma_{lmk}(q)$ are vital in determining the width, shape, and parameter
correlations of $P\left(\bm{\Theta}|C\right)$. To investigate such
dependencies we first define a $|\Psi(\mathbf{r})|^{2}$ distribution for NO2
to calculate $C_{lmk}(q)$. Figure 4a and Table 1 show and describe this
distribution, respectively. Measuring more diffraction patterns increases the
signal to noise ratio (SNR) by reducing $\sigma_{lmk}(q)$ which scales as
$1/\sqrt{N}$. Here, the SNR is the geometric mean of
$C_{000}(q)/\sigma_{000}(q)$ between $0.5<q<4$ Å-1. Figure 4b illustrates the
$C_{lmk}(q)$ coefficients used in this analysis with the following SNRs based
on previous UED Wolf  _et al._ (2019) and x-ray Ware  _et al._ (2019b)
diffraction experiments: 25, 50, 100, 200, 400. To calculate $\sigma_{lmk}(q)$
for NO2 we add Poisson noise to the diffraction patterns and propagate that
noise through the lab frame anisotropy and ADM fits outlined above (see
Supplementary Section XII). For the N2O data, Supplementary Section XII
describes the data processing and retrieval of $\sigma_{lmk}(q)$. Increasing
the average pump fluence by increasing the pump pulse FWHM induces a larger
rotational coherence while increasing the ensemble temperature increases the
spread of initial rotational states. Figure 4c shows how the ADMs’ amplitudes
increase and decrease by increasing the pump fluence and ensemble temperature,
respectively. We also note the ADMs oscillate faster and with more complexity
as both the pump fluence and temperature increase. Unless otherwise stated,
the standard configuration of experimental parameters for our NO2 results is a
$q$ range of [0.5, 10] Å-1, a SNR of 100, a pump fluence of 1 J/cm2 and a 100
K ensemble temperature.

This method ultimately yields the following three results; a distribution of
$\bm{\Theta}$ parameters (the posterior $P\left(\bm{\Theta}|C\right)$), the
optimal set of model parameters ($\bm{\Theta}^{*}$), and a parameterized
probability of molecular geometries $P(\mathbf{r}|\bm{\Theta}^{*},C)$. For
each individual $\bm{\Theta}$ parameter, where the $i^{\text{th}}$ parameter
is denoted as $\Theta_{i}$, we calculate its resolution as the standard
deviation of the projection onto $\Theta_{i}$. This resolution,
$\sigma^{\Theta}$, is the one dimensional standard deviation after
marginalizing over all other parameters, which removes the correlations
between $\bm{\Theta}$ parameters. That is, if one randomly draws some
parameters $\bm{\Theta}$ from $P\left(\bm{\Theta}|C\right)$, the distribution
of parameter $\Theta_{i}$ will have a width of $\sigma^{\Theta}$. The mean and
mode of said marginalized distribution will likely not correspond to
$\bm{\Theta}^{*}$, since $\bm{\Theta}^{*}$ is the mode of the full
$\bm{\Theta}$-space distribution. We find $\bm{\Theta}^{*}$ via a simple mode
search algorithm based on a weighted averaging search (see Supplementary
Section XIII).

##  III Results

Figure 6: Varying experimental parameters affects the resolution (width) of
$P^{(\mathcal{N})}\left(\bm{\Theta}|C\right)$, but our method is most
sensitive to SNR. Panel a shows how the uncorrelated widths of
$P^{(\mathcal{N})}\left(\bm{\Theta}|C\right)$, denoted by $\sigma^{\Theta}$,
changes by increasing $q$ range. Panel b similarly shows the dependence of
$\sigma^{\Theta}$ and $\bm{\Theta}^{*}$ error versus SNR. Panel c shows the
dependence of $\sigma^{\Theta}$ and $\bm{\Theta}^{*}$ error versus pump
fluence (width of the rotational wavepacket). Panel d shows the dependence of
$\sigma^{\Theta}$ and $\bm{\Theta}^{*}$ error versus the molecular ensemble
temperature. Figure 7: Varying the measured $q$ range affects false
correlations in $P^{(\mathcal{N})}\left(\bm{\Theta}|C\right)$ for NO2; a
larger reciprocal space provides more information and dampens false
correlations. Panel a shows the 1d and 2d projections of
$P^{(\mathcal{N})}\left(\bm{\Theta}|C\right)$ for a limited $q$ range of
$[0.5,5]$ Å-1. The red dashed lines illustrate $\bm{\Theta}^{*}$, while the
black “X” and solid lines indicate the ground truth values. Panel b shows the
corresponding $P^{(\mathcal{N})}(\mathbf{r}|\bm{\Theta}^{*},C)$. Similarly,
panel c shows the 1d and 2d projections of
$P^{(\mathcal{N})}\left(\bm{\Theta}|C\right)$ for the broader $q$ range of
$[0.5,20]$ Å-1. Panel d shows the corresponding
$P^{(\mathcal{N})}(\mathbf{r}|\bm{\Theta}^{*},C)$. Panel e shows the
correlation between all $\bm{\Theta}$ parameters as a function of $q$ range.
We note the decrease in correlations with larger $q$, where panels a and b
illustrate how the width and false correlations in
$P^{(\mathcal{N})}\left(\bm{\Theta}|C\right)$ decrease with higher $q$.

To illustrate our method’s efficacy, we retrieve
$P^{(\mathcal{N})}(\mathbf{r}|\bm{\Theta}^{*},C)$ from both simulated NO2
$C_{lmk}(q)$ coefficients under varying experimental conditions and measured
N2O data. As discussed, we first build the posterior
$P^{(\mathcal{N})}\left(\bm{\Theta}|C\right)$ which is shown in Fig. 5 for
simulated NO2 (a) and measured N2O data (c). Panels b and d show
$P^{(\mathcal{N})}(\mathbf{r}|\bm{\Theta}^{*},C)$ for NO2 and the simulated
PDF for N2O, respectively. Tables 1 and 2 give the extracted $\bm{\Theta}^{*}$
and $\sigma^{\Theta}$ for NO2 and N2O, respectively. For NO2,
$P^{(\mathcal{N})}\left(\bm{\Theta}|C\right)$, of width $\sigma^{\Theta}$ ~30
fm, fully encompasses the ground truth values. Despite the largely flat
$\expectationvalue{\angle\text{ONO}}$ distribution, $\bm{\Theta}^{*}$ still
converges on the ground truth values. For N2O data, the retrieved
$P^{(\mathcal{N})}\left(\bm{\Theta}|C\right)$ encompasses the previously
measured results of Ref. Herzberg (1966). The uncorrelated widths of this
distribution $\sigma^{\Theta}$ are of order picometers even with our limited
$q$ range of $[3.5,7.25]$ Å-1 and a very poor SNR. Most notably, the retrieved
$\expectationvalue{\angle\text{N}\text{N}\text{O}}$ is $\pi$ and we resolve
the ~5 pm difference between the N${}^{\text{T}}$N${}^{\text{C}}$ and
N${}^{\text{C}}$O bond distances (Table 2). The retrieved widths
$\sigma\left(\text{N}^{\text{T}}\text{N}^{\text{C}}\right)$ and
$\sigma\left(\angle\text{N}\text{N}\text{O}\right)$ are unphysical due to the
limited $q$ range, as discussed later. The PDF peaks’ ~25 pm widths (Fig. 5d),
are an order of magnitude larger than our resolution. The missing low and high
$q$ components produce ringing artifacts in this inverse Fourier Transform
because of the incomplete Fourier space. This confuses the PDF results as they
are not positive definite and these artifacts can be substantial at large
distances.

Figure 5a and c both show how $\bm{\Theta}^{*}$ does not correspond to the
mean or mode of most $P^{(\mathcal{N})}\left(\bm{\Theta}|C\right)$
projections. This is due to the non-linearity and correlations of
$P\left(\bm{\Theta}|C\right)$ in $\bm{\Theta}$ space. The plotted
$P\left(\bm{\Theta}|C\right)$ projections trace over variables. Tracing over
variables removes correlations that are important to the resolution of
$P\left(\bm{\Theta}|C\right)$ and accuracy of $\bm{\Theta}^{*}$. We must find
$\bm{\Theta}^{*}$ in this correlated space since the geometry parameters are
indeed correlated.

| Input | $\bm{\Theta}^{*}$ | $\sigma^{\Theta}$
---|---|---|---
$\expectationvalue{\text{NO}^{(1)}}$[Å] | 1.35 | 1.3500 | 0.00029
$\sigma\left(\text{NO}^{(1)}\right)$[Å] | 0.03 | 0.03000 | 0.0019
$\expectationvalue{\text{NO}^{(2)}}$[Å] | 1.05 | 1.0500 | 0.00029
$\sigma\left(\text{NO}^{(2)}\right)$[Å] | 0.02 | 0.02000 | 0.0054
$\expectationvalue{\angle\text{ONO}}$[rad] | 2.34 | 2.340 | 0.00047
$\sigma\left(\angle\text{ONO}\right)$[rad] | 0.01 | 0.01010 | 0.015
Table 1: We provide the input $\bm{\Theta}$ parameters for our NO2 simulation, the retrieved $\bm{\Theta}^{*}$ and the uncorrelated widths of $P\left(\bm{\Theta}|C\right)$, $\sigma^{\Theta}$. | $\bm{\Theta}^{*}_{\text{Literature}}$ | $\bm{\Theta}^{*}$ | $\sigma^{\Theta}$
---|---|---|---
$\expectationvalue{\text{N}^{\text{T}}\text{N}^{\text{C}}}$ [Å] | 1.128 | 1.142 | 0.039
$\sigma\left(\text{N}^{\text{T}}\text{N}^{\text{C}}\right)$ [Å] |  | 0.081 | 0.028
$\expectationvalue{\text{N}^{\text{C}}\text{O}}$ [Å] | 1.184 | 1.175 | 0.036
$\sigma\left(\text{N}^{\text{C}}\text{O}\right)$ [Å] |  | 3.08$\times 10^{-8}$ | 0.027
$\expectationvalue{\angle\text{N}\text{N}\text{O}}$ [rad] | 3.142 | 3.142 | 0.061
$\sigma\left(\angle\text{N}\text{N}\text{O}\right)$ [rad] |  | 5.5$\times 10^{-12}$ | 0.062
Table 2: We provide the literature values of $\bm{\Theta}^{*}$ Herzberg (1966)
with the retrieved $\bm{\Theta}^{*}$, along with the uncorrelated widths,
$\sigma^{\Theta}$, from the N2O data.

Our main results are given in Figs. 5 and Tables 1 and 2. We now show how
these results would behave by altering experimental parameters. Figures 6, 7,
and 8 show these exploratory results.

Expanding the measured $q$ range improves $\sigma^{\Theta}$ only until ~8 Å-1,
but continuous to reduce false correlations between $\bm{\Theta}$ parameters
beyond this. Figure 6a summarizes how $\sigma^{\Theta}$ decreases and then
plateaus with increasing $q$ range. Increasing the measured reciprocal range
$q$ provides more information about the system and reduces correlations in
$P\left(\bm{\Theta}|C\right)$ (Fig. 7). Since these correlations diminish with
added information, we refer to them as false correlations. Figures 7a and c
illustrate how the false correlations diminish from a $q$ range of [0.5, 5]
Å-1 to [0.5, 20] Å-1, respectively, under conditions of SNR of 100 and a
fluence of 1 J/cm2. The plotted correlation in Fig. 7e is between all six
$\bm{\Theta}$ parameters.

Our method is strongly sensitive to SNR due to our statistical interpretation.
Figure 6b shows that $\sigma^{\Theta}$ rapidly decreases with increasing SNR
with $q$ of $[0.5,10]$. Increasing SNR by an order of magnitude decreases
$\sigma^{\Theta}$ by an order of magnitude for pair-wise distances and angles.
The shape of $P^{(\mathcal{N})}\left(\bm{\Theta}|C\right)$ does not change
since we are not adding new information given the fixed $q$ range.

Increasing the induced rotational coherence and lowering the ensemble
temperature improves the measured resolution similarly to increasing the SNR.
This is caused by the increased magnitude and variation in the ADMs (Fig. 4c).
Figure 6c and d show how rapidly $P\left(\bm{\Theta}|C\right)$ resolution
improves by increasing rotational coherence (25 K ensemble temperature) and
decreasing ensemble temperature (1 J/cm2 pump fluence), respectively

Figure 8: The $P^{(\delta)}\left(\bm{\Theta}|C\right)$ distribution suffers
from a $q$ dependent systematic error stemming from the false assumption that
a single geometry describes the results measured from an ensemble. Here we
show the 1d projections of $P^{(\delta)}\left(\bm{\Theta}|C\right)$ and
$P^{(\mathcal{N})}\left(\bm{\Theta}|C\right)$ as a function of measured $q$
range. Each column indicates a different $q$ range starting at 0.5 Å-1 with
the end of said $q$ range indicated by the right most border of that column.
The dashed lines are the ground truth values. The bottom plot is the simulated
$C_{200}(q)$ coefficient used for both posteriors and is intersected by black
lines that indicate the upper $q$ range of each column.

Generally, when varying the $q$ range, SNR levels, pump fluence, and ensemble
temperature we find the pair-wise distances’ $\sigma^{\Theta}$ to be of order
100 fm; for the width parameters, $\sigma^{\Theta}$ is order 1 pm. Our
retrieved $\bm{\Theta}^{*}$ values are generally within a relative error of
~$10^{-7}$ and ~$10^{-3}$ from the ground truth values for geometric
parameters and width parameters, respectively. This resolution is often ~100
times better than PDF based methods because our statistical treatment is
highly sensitive to SNR.

As noted above, we run the intermediate delta distribution on the same NO2
simulation (Fig. 4). With
$P^{(\delta)}\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$, we assume a
single molecule response can describe a signal averaged over an ensemble of
geometries. This systematic error introduces an order of magnitude increase in
$\sigma^{\Theta}$. When increasing $q$,
$P^{(\delta)}\left(\bm{\Theta}|C\right)$ seems to converge in an unstable
fashion on the ground truth (Fig. 8), unlike
$P^{(\mathcal{N})}\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$. We note
that this calculation runs ~100 times faster than
$P^{(\mathcal{N})}\left(\bm{\Theta}|C\right)$ since it doesn’t have to sum
over many geometries in Eq. 13. For this reason
$P^{(\delta)}\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$ primarily
serves as the intermediate to quickly test one’s analysis before switching to
the normal distribution. For large molecules,
$P^{(\delta)}\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$ may be the
only tractable method. Supplementary Section XIV provides plots and further
discussion of these results.

##  IV Discussion

We divide the following discussion into three categories covering the two
major points of our analysis and how one may extend it to excited states
dynamics. Firstly, we provide intuition for the need of induced anisotropy to
access the molecular frame geometric angles $\theta^{(\text{mf})}_{\mu\nu}$
and $\phi^{(\text{mf})}_{\mu\nu}$. Secondly, we will compare conventional PDF
based methods to our Bayesian Inferencing approach. We also discuss how
systematic errors and various experimental parameters affect the MHA results.
Finally, we introduce methods to evaluate excited state dynamics by evaluating
only the isotropic component and by introducing a separation of timescale
approximation to Eq. 2. Both these methods do not require any complex excited
state dynamic simulations.

###  IV.1 The Role of Anisotropy

To provide intuition for the distinct angular terms, we condense and label the
important reference frames from Eq. 2

$\expectationvalue{I(q)_{\text{mol}}}=\cdots\sum_{m,k}(-1)^{k-m}\underbrace{Y^{-m}_{l}\left(\theta^{(\text{lf})}_{q},\phi^{(\text{lf})}_{q}\right)}_{\text{Lab
Frame}}\expectationvalue{\underbrace{D^{l}_{mk}\left(\phi_{\text{I}}^{(\text{lf})},\theta_{\text{I}}^{(\text{lf})},\chi_{\text{I}}^{(\text{lf})}\right)}_{\text{Ensemble
Anisotropy}}\underbrace{j_{l}(q\Delta
r_{\mu\nu})Y_{l}^{-k}\left(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}\right)}_{\text{Molecular
Frame Geometry}}}{\Psi(t)}.$ (16)

The molecular frame geometry component separates into distance and angular
terms. The former is the pair-wise distance term ($j_{l}(q\Delta r_{\mu\nu})$)
that governs the $q$ dependence. The latter
$(Y_{l}^{-k}(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}))$, is
the angular decomposition of the molecular geometry and acts as a scaling
parameter. This molecular frame geometry term is coupled to the measured lab
frame anisotropy
($Y^{-m}_{l}(\theta^{(\text{lf})}_{q},\phi^{(\text{lf})}_{q})$) by the
ensemble anisotropy. The ensemble anisotropy,
$D^{l}_{mk}\left(\phi_{\text{I}}^{(\text{lf})},\theta_{\text{I}}^{(\text{lf})},\chi_{\text{I}}^{(\text{lf})}\right)$,
acts as a key from the measured lab frame diffraction pattern to the molecular
frame geometry. Similar derivations in Refs. Baskin and Zewail (2006) and
Parrish and Martinez (2019) exist but do not stress the dependence on
molecular frame coordinates; Ref. Baskin and Zewail (2006) is not treated
fully quantum mechanically as done here in Supplementary Section IX.
Anisotropy is required for our method to have an explicit dependence on the
pair-wise angles. Without anisotropy, $C_{000}(q)$ has no explicit angular
dependence (Eq. 5), just like the PDF.

Stronger impulsive alignment produces a broader coherent rotational wavepacket
which exhibits higher amplitude signals with more variations (Fig. 4c). Larger
amplitude ADMs improve $C_{lmk}(q)$ SNR by lifting higher order coefficients
up out of the noise, resulting in similar resolution improvements to only
increasing SNR, shown in Fig. 6c. Increasing the number of $C_{lmk}(q)$
coefficients improves the $\theta_{\mu\nu}^{(\text{mf})}$ and
$\phi_{\mu\nu}^{(\text{mf})}$ resolution since each $C_{lmk}(q)$ provides a
new angular constraint via
$Y_{l}^{-k}(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})})$ (Eq.
5).

One can produce fast signal variations with an initially broad hot thermal
ensemble. Writing coherence onto hotter molecular ensembles produces weak but
fast varying ADMs, shown in Fig. 4c. Figure 6d shows how quickly the
resolution worsens at higher temperatures. When fitting the ADMs to
$B^{m}_{l}(q,t)$, one ideally measures particular points that include two
separate regions where the ADMs have high variation and sufficiently before
and after the prominent anisotropy signal where their magnitude dampens. We
note, however, one need not strictly measure the entire transient rotational
signal.

When simulating or inducing molecular tumbling is prohibitively difficult, one
may use the induced anisotropy from the photo-excitation. Producing
sufficiently large rotational coherences in large molecules is challenging.
The option of using dipole alignment from the excitation pulse allows this
technique to be more generally applicable to most molecular systems of
interest.

###  IV.2 Bayesian Inferencing and the MHA

We use Bayesian Inferencing to retrieve $P\left(\bm{\Theta}|C\right)$ and
ultimately the $P(\mathbf{r}|\bm{\Theta}^{*},C)$ distribution.
$P(\mathbf{r}|\bm{\Theta}^{*},C)$ is significantly more information rich than
the simple PDF. Our method decomposes anisotropic diffraction into individual
$C_{lmk}(q)$, each of which independently constrains
$P\left(\bm{\Theta}|C\right)$ and consequently the molecular geometries. The
$P(\mathbf{r}|\bm{\Theta}^{*},C)$ distribution, which approximates
$|\Psi(\mathbf{r})|^{2}$, allows us to find the most probable pair-wise
distances and angles to define a most-likely geometry that is unique.
Traditionally, an ideal PDF derived from an infinite $q$ range, being the
inverse Fourier transform of $q\,M_{000}(q)$, is at best a weighted histogram
of unlabeled pair-wise distances. Since our measurements necessarily do not
measure $q$ all the way to 0, and the signal drop-off limits our ability to
measure high $q$, our $q$ range is always limited. These limitations obfuscate
the PDF interpretations by introducing sinusoidal systematics that result in
negative probabilities, shown in Fig 5d where we do not expect any distance
above 2.3 Å. To alleviate the ambiguity of Fig. 5d, we typically simulate
$|\Psi(\mathbf{r})|^{2}$ with a priori knowledge and validate simulation
against the measured PDF. Instead of this traditional approach, our method
uncovers the globally optimal parameters ($\bm{\Theta}^{*}$) from the data for
a given $P\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$. This requires
only the initial ground state geometry, simulations of the coherent rotational
wavepacket when using $C_{lmk}(q)$ for $l>0$, and for excited state dynamics
one additionally needs relevant transition dipole moments. As made clear by
comparing Figs. 5b and d, the $P(\mathbf{r}|\bm{\Theta}^{*},C)$ distribution
is significantly more information rich than the simple PDF, e.g. it provides
the 3d molecular geometry and width of the $|\Psi(\mathbf{r})|^{2}$. This
method thus has the potential to shift ultrafast diffraction to a discovery
method applicable even to systems that extend beyond the scope of theory.

Our method retrieves the labelled pair-wise distances with ~100 times better
resolution than the PDF. This is due to our statistical treatment using
Bayesian Inferencing and the MHA. In the MHA, each $lmk$ and $q$ contribution
is itself a probability distribution; each is an experiment of its own. When
adding a new $\bm{\Theta}$ to $P\left(\bm{\Theta}|C\right)$ the relative ratio
of the likelihood function $P(C|\bm{\Theta})$ (Eq. 15) decreases since each
additional $C_{lmk}(q)$ increases the magnitude of the negative exponent. That
is, with more measurements Eq. 15 becomes smaller and
$P\left(\bm{\Theta}|C\right)$ becomes sharper. Our method therefore heavily
relies on $\sigma_{lmk}(q)$ and $C_{lmk}(q)$. Statistical noise tends to
increase $\sigma_{lmk}(q)$, making $P\left(\bm{\Theta}|C\right)$ wider
($\sigma^{\Theta}$ larger), while systematic errors in $C_{lmk}(q)$ tend to
shift the centriod of $P\left(\bm{\Theta}|C\right)$, as seen in the N2O data.
Supplementary Section XII describes our method for consistently accounting for
both statistical and systematic errors. The PDF error adds in quadrature in
$\sigma_{lmk}(q)$; its scale is set by the largest error bar and
disproportionately suffers from poorly measured data points. Conversely, since
each term of the negative sum in the argument of the exponent in Eq. 15 is
weighted by $1/\sigma_{lmk}(q)$, the MHA reduces the contribution from poorly
measured data points while increasing the effect of high precision
measurements.

Our Bayesian Inferencing approach expands the utility of gas-phase ultrafast
diffraction to measure previously inaccessible variables, such as the width of
$|\Psi(\mathbf{r})|^{2}$. Given
$P\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$ is a generic function
parameterized by $\bm{\Theta}$, one can introduce variables through
$\bm{\Theta}$ by selecting a
$P\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$ that depends on it. Here,
we expanded the measurable parameters of gas-phase ultrafast diffraction to
include the width of $|\Psi(\mathbf{r})|^{2}$.
$P^{(\mathcal{N})}\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$
quantifiably measures $|\Psi(\mathbf{r})|^{2}$ width and reduces the
systematic effects from assuming a single geometry, illustrated in Fig. 8,
even when $|\Psi(\mathbf{r})|^{2}$ is not well modeled by a Gaussian. One can
choose a more accurate $P(\mathbf{r}|\bm{\Theta}^{*},C)$ for their system
e.g., one may choose Hermite polynomials for vibronically excited states. For
the N2O data, our $q$ range of [3.5, 7.25] Å-1 is insufficient to resolve the
width of $|\Psi^{(\text{N}_{2}\text{O})}(\mathbf{r})|^{2}$, making the width a
nuisance fit parameter used to avoid the above-mentioned systematic errors
induced by the delta distribution. Instead, we can use the retrieved distances
from $\bm{\Theta}^{*}$ (using the normal distribution) as input into the delta
distribution in order to approximate $|\Psi(\mathbf{r})|^{2}$.

The MHA performs an unbiased search through $\bm{\Theta}$ space guided by the
$C_{lmk}(q)$ coefficients and correlates each $\bm{\Theta}$ parameter. Our
method is model independent and does not suffer from model bias as might be a
concern for conventional methods. Simultaneously evaluating each $\bm{\Theta}$
parameter artificially introduces correlations between them. Since
$\bm{\Theta}$ is the minimal set of parameters to define
$P\left(\mathbf{r}\right|\left.\bm{\Theta},C\right)$, we expect the parameters
to be uncorrelated. Figure 7 shows how adding information by extending the $q$
range decreases false correlations, which are caused by integrating over
insufficiently large $q$ regions. For the N2O data, we observe these false
correlations, most notably between
$\expectationvalue{\text{N}^{\text{T}}\text{N}^{\text{C}}}$ vs
$\expectationvalue{\text{N}^{\text{C}}\text{O}}$ (Fig. 5c). Simultaneously
evaluating all $\bm{\Theta}$ parameters leverages well-resolved parameters to
constrain poorly resolved parameters. For example, the long OO bond (or
$\angle\text{ONO}$) in our asymmetric NO2 is the best constrained parameter as
it produces the most $q$ oscillations. The MHA removes geometries where the
two NO distances are inconsistent with the well-resolved OO distance. These
correlations similarly help find $\bm{\Theta}^{*}$, as observed with N2O,
where the $P^{(\mathcal{N})}\left(\bm{\Theta}|C\right)$ uncorrelated widths do
not distinguish the $\text{N}^{\text{T}}\text{N}^{\text{C}}$ and
$\text{N}^{\text{C}}\text{O}$ bonds but $\bm{\Theta}^{*}$ does.

The width of $P\left(\bm{\Theta}|C\right)$ ($\sigma^{\Theta}$) relies heavily
on SNR rather than increased $q$ range, which is ideal since it is generally
prohibitively difficult to change the $q$ range at ultrafast diffraction
facilities and easier to reduce the SNR by taking more measurements.
Similarly, Ref. Ma  _et al._ (2020) illustrated the importance of SNR. Below 8
Å-1, $\sigma^{\Theta}$ for the pair-wise distances and angles quickly
decreases as the $q$ range increases, shown in Figs. 6a and 7. Above 8 Å-1,
$\sigma^{\Theta}$ for the pair-wise distances and angles does not change
significantly. In contrast, $\sigma^{\Theta}$ continuously decreases with
improved SNR, shown in Fig. 6b. This is because smaller $\sigma_{lmk}(q)$
makes it less probable for the MHA to visit $\bm{\Theta}$ parameters with
larger residuals. Which is again due to increasing the magnitude of the
exponential’s argument in Eq. 15. For the PDF, the resolution goes as
$2\pi/\Delta q$, or 1.26, 0.63, and 0.31 Å for $q$ ranges of 5, 10, and 20 Å-1
respectively, which is roughly 1000 times larger than our observed resolution
for NO2. Our method therefore lends itself well to high repetition-rate
machines, such as the upcoming LCLS II. We note that increasing the $q$ range
above 8 Å-1 has a larger effect on the width parameters (dashed lines of Fig.
6a) and that one may prioritize higher $q$ over SNR if they are more
interested in the width of $|\Psi(\mathbf{r})|^{2}$.

###  IV.3 Potential Extension to Excited State Dynamics

Our method is broadly applicable to diffraction experiments with laser
excitation, including dynamics from excited electronic states. Here, the laser
excitation imparts one or more units of angular momentum providing at least
the anisotropy terms $C_{20k}(q)$. From N2O data and from low SNR NO2
simulations we see that the $C_{200}(q)$ alone recovers ~4 pm resolution. The
primary difficulty with extending our method to excited states dynamics lies
in isolating the ADMs in rovibronically coupled systems at sufficiently long
timescales. Since the principle moments of inertia change with the geometry,
one must reorient the altered excited state molecular geometry by adding three
molecular frame Euler angles to the $\bm{\Theta}$ parameters. The generally
much wider width of $|\Psi(\mathbf{r},t)|^{2}$ in the excited state dampens
the $C_{lmk}(q)$ coefficients faster than the ground state signal, reducing
the need for extending $q$. Next we discuss two variants, a time separable
method and an isotropic method. Both variations isolate the ADMs as we did in
the rigid-rotor approximation.

The time separable method introduces a separation of time scales since the
vibronic dynamics are typically much faster than the rotational dynamics. Here
we assume the ADMs are relatively stationary during the vibronic motion. This
approximation is analogous to the Born-Oppenheimer approximation and allows us
to separate the rotational and vibronic timescales. Applying this
approximation to an experiment with a single excitation pulse, the measured
results are given by Eq. 17,

$\displaystyle\langle I(\textbf{q})\rangle_{\text{sep}}^{(1)}$
$\displaystyle(t)={}\mathcal{I}\bigg{(}\sum_{\mu}|f_{\mu}(q)|^{2}+\sum_{\mu,\nu:\mu\neq\nu}\text{Re}\bigg{\\{}f_{\mu}(q)f_{\nu}(q)\sum_{l}i^{l}8\pi^{2}\sqrt{\frac{4\pi}{\left(2l+1\right)}}$
(17) $\displaystyle\times\sum_{m_{1},m_{2}}$
$\displaystyle(-1)^{m_{1}-m_{2}}Y^{-m_{2}}_{l}\left(\theta^{(\text{lf})}_{q},\phi^{(\text{lf})}_{q}\right)\expectationvalue{j_{l}(q\Delta
r_{\mu\nu})Y_{l}^{-m_{1}}\left(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}\right)}{\Psi_{\text{vib}}(t)}\mathcal{A}^{l}_{m_{2}m_{1}}(0)\bigg{\\}}\bigg{)}.$

Here, $\mathcal{A}^{l}_{m_{2}m_{1}}(0)$ are the ADMs calculated with respect
to the ground state geometry and evaluated immediately after laser excitation.
The dipole selection rule induces the above anisotropy, and for a single
photon transition provides at least the $C_{20k}(q,t)$ coefficients. This
variation is independent of one’s ability to induce a rotational wavepacket,
particularly for large molecules. However, it requires knowledge of either the
transition dipole moment or the Frank-Condon factor and the ground state
dipole.

To further constrain $P\left(\bm{\Theta}|C\right)$, one can couple to more
$C_{lmk}(q)$ coefficients by introducing a precursor pulse that excites a
rotational wavepacket. This precursor pulse, assumed to be a rotational Raman
impulse, is chosen to have a negligible effect on the vibronic system thus
maintaining consistency with our separation of timescale approximation. The
Raman impulse first induces the rotational coherence. Following the Raman
impulse, the system evolves for a rotational time $\tau$, at this point the
vibronic excitation pulse arrives. One would measure the vibronic dynamics
over a small window $(t\ll\tau)$ for different orientations by also scanning
the delay $\tau$ over an appreciable portion of the rotational evolution. This
window, measured by $t$, is typically of order picosecond or less such that
the $\mathcal{A}^{l}_{mk}(t)$ do not appreciably change. The measured
diffraction images are given by Eq. 18

$\displaystyle\langle I(\textbf{q})\rangle_{\text{sep}}^{(2)}(t,\tau)$
$\displaystyle={}\mathcal{I}\bigg{(}\sum_{\mu}|f_{\mu}(q)|^{2}+\sum_{\mu,\nu:\mu\neq\nu}\text{Re}\bigg{\\{}f_{\mu}(q)f^{*}_{\nu}(q)\sum_{l}i^{l}8\pi^{2}\sqrt{4\pi\left(2l+1\right)}\sum_{m_{1},m_{2}}(-1)^{m_{1}-m_{2}}$
(18) $\displaystyle\times
Y^{-m_{2}}_{l}\left(\theta^{(\text{lf})}_{q},\phi^{(\text{lf})}_{q}\right)\tilde{\mathcal{A}}^{l}_{m_{1}m_{2}}(n,n^{\prime};\tau)\matrixelement{\psi^{n^{\prime}}_{\text{el-
vib}}(t)}{j_{l}(q\Delta
r_{\mu\nu})Y_{l}^{-m_{1}}\left(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}\right)}{\psi^{n}_{\text{el-
vib}}(t)}\bigg{\\}}\bigg{)}$

where $n$ labels the vibronic states, $\ket{\psi^{n}_{\text{el-vib}}(t)}$ is
the vibronic wavefunction (assumed unknown),
$\tilde{\mathcal{A}}^{l}_{m_{1}m_{2}}(n,n^{\prime};\tau)$ are the modified
ADMs, and $t$ is the arrival time of the probe after the second excitation
pulse. These modified ADMs consider the angular momentum transfer by the
photon for the excitation to the vibronic state. This only requires
calculating the ground state transition dipole moment. One then follows the
same analysis procedure outlined above for the ground state for all
$C_{lmk}(q,t)$ coefficients at each probe time $t$. In such an experiment, one
should measure the ensemble anisotropy without the vibronic excitation pulse.
This would provide the best fit ADMs via fits to pulse intensity and
temperature. Supplementary Section IX further describes our separation of
timescale approximation and provides the derivations for Eqs. 17 and 18, and
describes the required molecular frame orientation angles when $l\neq 0$.

The isotropic method uses only the $C_{000}(q,t)$ term similar to conventional
analyses. The $\mathcal{A}^{0}_{00}(t)$ is a constant which eliminates the
need to simulate either the full ADMs or the full modified ADMs. This places
all the time dependence into the molecular frame geometry term of Eq. 16.
Thus, it can be applied to either single (Eq. 17) or double pulse (Eq. 18)
modalities and relies on the Bayesian Inferencing and MHA mapping rather than
anisotropy. One then performs the MHA for each $C_{000}(q,t)$. The
$C_{000}(q,t)$ term only implicitly depends on the pair-wise angles through
$\Delta r_{\mu\nu}$. This is in contrast to the explicit pair-wise angle
dependence in the higher order $C_{lmk}(q)$ terms. Our statistical treatment
likely provides adequate pair-wise angle resolution because we have more pair-
wise distances than are required to uniquely specify a geometry.

When using a Raman inducing precursor pulse, one will likely use a combination
of the isotropic and time separable methods. For fast dynamics one would use
the time separable method for small windows shortly following the rotation
time $\tau$. Longer lived dynamics can be retrieved by the isotropic method.
When retrieving $P\left(\bm{\Theta}|C\right)$ in either case one initiates the
MHA with the ground state $\bm{\Theta}^{*}$ parameters. For each subsequent
time step one initiates MHA with the $\bm{\Theta}^{*}$ parameters from the
previous time step.

Vibronic wavepackets bifurcate into multiple states, e.g. at conical
intersections, causing $P(\mathbf{r},t|\bm{\Theta}^{*},C)$ to bifurcate as
well. We account for these different states by

$P(\mathbf{r},t|\bm{\Theta}^{*},C)=\sum_{i}^{N_{\text{ex}}}c_{i}P(\mathbf{r},\bm{\Theta}^{*}_{i},t|C)$
(19)

where $N_{\text{ex}}$ is the number of excited state distributions with
appreciable population. Conical intersections will induce bifurcations that
spawn a new distribution that adds to $N_{\text{ex}}$. In this way we consider
this method to be fully data-driven since we can change our theoretical
description ($c_{i}$) based on data alone.

In the excited state, $|\Psi(\mathbf{r},t)|^{2}$ often takes on shapes that
are poorly represented by normal distributions Makhija  _et al._ (2021).
Nevertheless, the normal distribution acts as a second order improvement upon
fitting with a single geometry Stankus  _et al._ (2019); Natan (2021). This
improvement dampens the single geometry systematics shown in Fig. 8 and
retrieves a quantifiable measurement of $|\Psi(\mathbf{r},t)|^{2}$ widths.
This is particularly useful for vibronic wavepackets with appreciable width.
One can better describe amorphous $|\Psi(\mathbf{r},t)|^{2}$ Felker and Zewail
(1984) by using more representative distributions, such as Hermite
polynomials, for $P(\mathbf{r},t|\bm{\Theta}^{*},C)$.

Thus far we have only considered diffraction consistent with the independent
atom approximation and all the equations above have been derived under this
approximation. Recently, however, diffraction beyond the independent atom
approximation has been observed in both electron Yang  _et al._ (2020) and
x-ray diffraction Yong  _et al._ (2020). Under such conditions this method
must be modified by either re-deriving the above equations to consider these
effects or by accounting for this signal in the $C_{lmk}(q)$ coefficients. For
MeV electron diffraction, inelastic scattering is limited to the low $q<1$ Å-1
region and can be easily removed from the $C_{lmk}(q)$ coefficients. For
x-rays, diffraction beyond the independent atom approximation contributes a
constant offset after the initial signal turn-on that spans the entire $q$
range. References Yong  _et al._ (2020); Stankus  _et al._ (2019) observed a
constant signal from the electron vacancy after exciting the molecule into a
Rydberg state. Since this signal does not vary appreciably in time, it can
safely be subtracted out before applying our method. Alternatively, one can
alter $P(\mathbf{r},t|\bm{\Theta}^{*},C)$ to include this offset by adding a
corresponding parameter to $\mathbf{\Theta}$.

##  V Conclusion and outlook

We have shown that our method can approximate $|\Psi(\mathbf{r})|^{2}$ with
$P(\mathbf{r}|\bm{\Theta}^{*},C)$ for the ground states of NO2 and N2O. In
simulation, we retrieve ~100 fm resolution for NO2. From measured N2O UED
data, we retrieve ~4 pm resolution despite a short $q$ range of [3.5, 7.25]
Å-1 and high SNR. This method returns the labeled pair-wise distances and the
angles with ~100 times better resolution than PDF-based methods. Furthermore,
we can retrieve the unique molecular geometry more directly from the data. In
spite of similar bond distances and atomic scattering amplitudes for NO2, our
method distinguishes bond distances from the $C_{200}(q)$ coefficient alone.
Adding a comparatively low SNR and a narrow $q$ range of [3.5, 7.25] Å-1 to
the UED data we still observed picometer resolution from the N2O data and can
even begin to resolve the
$\langle\text{N}^{\text{T}}\text{N}^{\text{C}}\rangle$ and
$\langle\text{N}^{\text{C}}O\rangle$ distances. These results are highly
encouraging and illustrate the viability of our Bayesian Inferencing approach
and inspires further expansion into excited states dynamics. The code
repository github (2022) contains the algorithms used for this work and
instructions on how to reproduce these results. It also contains instructions
on how to run this analysis and templates for applying this method to new
molecules.

This Bayesian Inferencing approach is best suited for gas-phase ultrafast
diffraction instruments that have high SNR such as high repetition-rate free
electron facilities, e.g. LCLS-II-HE. Resolution quickly improves with SNR
considerably faster than if one increases $q$ beyond ~8 Å-1. Illustrated in
our measurement of N2O with a limited $q$ range below 8 Å-1, [3.5, 7.25] Å-1,
we demonstrate ~40 times better resolution than the Fourier limit.
Nevertheless, larger $q$ ranges improve resolution for widths of
$|\Psi(\mathbf{r})|^{2}$ and diminish false correlations between $\bm{\Theta}$
parameters.

Our general method has the potential to become common-place for ultrafast gas-
phase diffraction measurements due to its broad applicability and its
independence from complex excited state simulations. In this work, we
validated its use for standard pump probe setup. One can extend this method to
excited state dynamics either with or without anisotropy. Our isotropic method
is well suited for current pump probe setups that generally focus on the
isotropic component. This method however greatly benefits from deterministic
anisotropy that can either be induced by impulsive Raman or by the dipole
moment selection from the excitation pulse. Beyond ultrafast gas-phase
diffraction, one can apply this general framework to other classes of
experiments, e.g. the previously mentioned photo-electron experiments Makhija
_et al._ (2016); Marceau  _et al._ (2017); Gregory  _et al._ (2021); Sándor
_et al._ (2019, 2018b). This is done by deriving the molecular frame response
(Eq. 2) and applying our Bayesian Inferencing approach.

Given its broad applicability, high resolution, amenability to various
measurements, and independence from complex molecular dynamic simulations, our
method has the potential to effectively turn ultrafast gas-phase molecular
diffraction into a discovery-oriented technique. This method can retrieve the
unique molecular geometry of a general molecule to ~100 fm resolution without
relying on complex molecular simulations. Moreover, because our method is
parameterized by $\bm{\Theta}$, we have the opportunity to expand the scope of
ultrafast gas-phase diffraction into previously inaccessible measurements. For
instance, we demonstrated the use of this parameterization to measure the
width of $|\Psi(\mathbf{r},t)|^{2}$; this width is important in the excited
state where single geometries lose their meaning. This method unlocks our
ability to study larger and more complex systems that are currently too
difficult to simulate.

##  VI Data Availability

The UED N2O data used in this analysis will be provided by the corresponding
authors upon reasonable request. The simulated NO2 data, $C_{lmk}(q)$, can be
calculated by the supplied analysis code in Ref. github (2022).

##  VII Code Availability

The code used in this analysis can be found in Ref. github (2022). Here, one
will find a detailed description of the code and how to run it in order to
reproduce the NO2 results. This repository also includes templates for one to
apply this algorithm to new molecules.

###### Acknowledgements.

Use of the Linac Coherent Light Source (LCLS), SLAC National Accelerator
Laboratory, is supported by the U.S. Department of Energy, Office of Science,
Office of Basic Energy Sciences under Contract No. DE-AC02-76SF00515. The UED
work was performed at SLAC MeVUED, which is supported in part by the DOE BES
SUF Division Accelerator and Detector research and development program, the
LCLS Facility, and SLAC under contract Nos. DE-AC02-05-CH11231 and DE-
AC02-76SF00515. Markus Ilchen acknowledges funding by the Volkswagen
foundation for a Peter-Paul-Ewald Fellowship. We thank Markus Gühr for his
help in setting up the SLAC gas phase UED experiment. We thank Theodore
Vecchione for his help operating the SLAC UED facility. We thank Gregory
Stewart for creating Figs. 1 and S2.

## References

  * Williamson  _et al._ (1997) J. Charles Williamson, Jianming Cao, Hyotcherl Ihee, Hans Frey, and Ahmed H. Zewail, “Clocking transient chemical changes by ultrafast electron diffraction,” Nature 386, 159–162 (1997).
  * Ewbank  _et al._ (1992) John Ewbank, Walter Faust, Jing Luo, Jeffrey English, David Monts, David Paul, Qun Dou, and Lothar Schafer, “Instrumentation for gas electron diffraction employing a pulsed electron beam synchronous with photoexcitation,” Review of Scientific Instruments 63, 3352 – 3358 (1992).
  * Ewbank  _et al._ (1993) John Ewbank, Jing Luo, Jeffrey English, Ruifeng Liu, Walter Faust, and Lothar Schafer, “Time-resolved gas electron diffraction study of the 193-nm photolysis of 1,2-dichloroethenes,” The Journal of Physical Chemistry 97 (1993), 10.1021/j100136a016.
  * Shen  _et al._ (2019) X. Shen, J. P. F. Nunes, J. Yang, R. K. Jobe, R. K. Li, Ming-Fu Lin, B. Moore, M. Niebuhr, S. P. Weathersby, T. J. A. Wolf, C. Yoneda, Markus Guehr, Martin Centurion, and X. J. Wang, “Femtosecond gas-phase mega-electron-volt ultrafast electron diffraction,” Structural Dynamics 6, 054305 (2019), https://doi.org/10.1063/1.5120864 .
  * Stankus  _et al._ (2020) Brian Stankus, Haiwang Yong, Jennifer Ruddock, Lingyu Ma, Andres Moreno Carrascosa, Nathan Goff, Sébastien Boutet, Xuan Xu, Nikola Zotev, Adam Kirrander, Michael P Minitti, and Peter M Weber, “Advances in ultrafast gas-phase x-ray scattering,” Journal of Physics B: Atomic, Molecular and Optical Physics 53, 234004 (2020).
  * Minitti  _et al._ (2015) M. P. Minitti, J. M. Budarz, A. Kirrander, J. S. Robinson, D. Ratner, T. J. Lane, D. Zhu, J. M. Glownia, M. Kozina, H. T. Lemke, M. Sikorski, Y. Feng, S. Nelson, K. Saita, B. Stankus, T. Northey, J. B. Hastings, and P. M. Weber, “Imaging molecular motion: Femtosecond x-ray scattering of an electrocyclic chemical reaction,” Phys. Rev. Lett. 114, 255501 (2015).
  * Ware  _et al._ (2019a) Matthew R. Ware, James M. Glownia, Adi Natan, James P. Cryan, and Philip H. Bucksbaum, “On the limits of observing motion in time-resolved x-ray scattering,” Phil. Trans. R. Soc. A (2019a), 10.1098/rsta.2017.0477.
  * Ware  _et al._ (2019b) Matthew R. Ware, James M. Glownia, Noor Al-Sayyad, Jordan T. O’Neal, and Philip H. Bucksbaum, “Characterizing dissociative motion in time-resolved x-ray scattering from gas-phase diatomic molecules,” Phys. Rev. A 100, 033413 (2019b).
  * Bucksbaum  _et al._ (2020) Philip H. Bucksbaum, Matthew R. Ware, Adi Natan, James P. Cryan, and James M. Glownia, “Characterizing multiphoton excitation using time-resolved x-ray scattering,” Phys. Rev. X 10, 011065 (2020).
  * Saldin  _et al._ (2010) D. K. Saldin, V. L. Shneerson, D. Starodub, and J. C. H. Spence, “Reconstruction from a single diffraction pattern of azimuthally projected electron density of molecules aligned parallel to a single axis,” Acta Crystallographica Section A 66, 32–37 (2010).
  * Starodub  _et al._ (2010) Dmitri Starodub, John C. H. Spence, and Dilano K. Saldin, “Reconstruction of the electron density of molecules with single-axis alignment,” in  _Image Reconstruction from Incomplete Data VI_, Vol. 7800, edited by Philip J. Bones, Michael A. Fiddy, and Rick P. Millane, International Society for Optics and Photonics (SPIE, 2010) pp. 188 – 196.
  * Ho  _et al._ (2009) P. J. Ho, D. Starodub, D. K. Saldin, V. L. Shneerson, A. Ourmazd, and R. Santra, “Molecular structure determination from x-ray scattering patterns of laser-aligned symmetric-top molecules,” The Journal of Chemical Physics 131, 131101 (2009), https://doi.org/10.1063/1.3245404 .
  * Elser (2011) Veit Elser, “Three-dimensional structure from intensity correlations,” New Journal of Physics 13, 123014 (2011).
  * Poon  _et al._ (2013) H.-C. Poon, P. Schwander, M. Uddin, and D. K. Saldin, “Fiber diffraction without fibers,” Phys. Rev. Lett. 110, 265505 (2013).
  * Reckenthaeler  _et al._ (2009) Peter Reckenthaeler, Martin Centurion, Werner Fuß, Sergei A. Trushin, Ferenc Krausz, and Ernst E. Fill, “Time-resolved electron diffraction from selectively aligned molecules,” Phys. Rev. Lett. 102, 213001 (2009).
  * Pabst  _et al._ (2010a) Stefan Pabst, Phay J. Ho, and Robin Santra, “Computational studies of x-ray scattering from three-dimensionally-aligned asymmetric-top molecules,” Phys. Rev. A 81, 043425 (2010a).
  * Yang  _et al._ (2014) Jie Yang, Varun Makhija, Vinod Kumarappan, and Martin Centurion, “Reconstruction of three-dimensional molecular structure from diffraction of laser-aligned molecules,” Structural Dynamics 1, 044101 (2014), https://doi.org/10.1063/1.4889840 .
  * Hensley  _et al._ (2012) Christopher J. Hensley, Jie Yang, and Martin Centurion, “Imaging of isolated molecules with ultrafast electron pulses,” Phys. Rev. Lett. 109, 133202 (2012).
  * Stankus  _et al._ (2019) Brian Stankus, Haiwang Yong, Nikola Zotev, Jennifer M. Ruddock, Darren Bellshaw, Thomas J. Lane, Mengning Liang, Sébastien Boutet, Sergio Carbajo, Joseph S. Robinson, Wenpeng Du, Nathan Goff, Yu Chang, Jason E. Koglin, Michael P. Minitti, Adam Kirrander, and Peter M. Weber, “Ultrafast x-ray scattering reveals vibrational coherence following rydberg excitation,” Nature Chemistry 11, 716–721 (2019).
  * Natan (2021) Adi Natan, “Real-space inversion and super-resolution of ultrafast scattering using natural scattering kernels,”  (2021).
  * Natan  _et al._ (2021) Adi Natan, Aviad Schori, Grace Owolabi, James P. Cryan, James M. Glownia, and Philip H. Bucksbaum, “Resolving multiphoton processes with high-order anisotropy ultrafast x-ray scattering,” Faraday Discuss. 228, 123–138 (2021).
  * Makhija  _et al._ (2016) Varun Makhija, Xiaoming Ren, Drue Gockel, Anh-Thu Le, and Vinod Kumarappan, “Orientation resolution through rotational coherence spectroscopy,” (2016), arXiv:1611.06476 [physics.atom-ph] .
  * Marceau  _et al._ (2017) Claude Marceau, Varun Makhija, Dominique Platzer, A. Yu. Naumov, P. B. Corkum, Albert Stolow, D. M. Villeneuve, and Paul Hockett, “Molecular frame reconstruction using time-domain photoionization interferometry,” Phys. Rev. Lett. 119, 083401 (2017).
  * Gregory  _et al._ (2021) Margaret Gregory, Paul Hockett, Albert Stolow, and Varun Makhija, “Towards molecular frame photoelectron angular distributions in polyatomic molecules from lab frame coherent rotational wavepacket evolution,” Journal of Physics B: Atomic, Molecular and Optical Physics 54, 145601 (2021).
  * Sándor  _et al._ (2018a) Péter Sándor, Adonay Sissay, Fran çois Mauger, Paul M. Abanador, Timothy T. Gorman, Timothy D. Scarborough, Mette B. Gaarde, Kenneth Lopata, Kenneth J. Schafer, and Robert R. Jones, “Angle dependence of strong-field single and double ionization of carbonyl sulfide,” Phys. Rev. A 98, 043425 (2018a).
  * Mikosch  _et al._ (2013) Jochen Mikosch, Andrey E. Boguslavskiy, Iain Wilkinson, Michael Spanner, Serguei Patchkovskii, and Albert Stolow, “Channel- and angle-resolved above threshold ionization in the molecular frame,” Phys. Rev. Lett. 110, 023004 (2013).
  * Sándor  _et al._ (2019) Péter Sándor, Adonay Sissay, François Mauger, Mark W. Gordon, T. T. Gorman, T. D. Scarborough, Mette B. Gaarde, Kenneth Lopata, K. J. Schafer, and R. R. Jones, “Angle-dependent strong-field ionization of halomethanes,” The Journal of Chemical Physics 151, 194308 (2019), https://doi.org/10.1063/1.5121711 .
  * Sándor  _et al._ (2018b) Péter Sándor, Adonay Sissay, Fran çois Mauger, Paul M. Abanador, Timothy T. Gorman, Timothy D. Scarborough, Mette B. Gaarde, Kenneth Lopata, Kenneth J. Schafer, and Robert R. Jones, “Angle dependence of strong-field single and double ionization of carbonyl sulfide,” Phys. Rev. A 98, 043425 (2018b).
  * Stolow and Underwood (2008) Albert Stolow and Jonathan G. Underwood, “Time-resolved photoelectron spectroscopy of nonadiabatic dynamics in polyatomic molecules,” in  _Advances in Chemical Physics_ (John Wiley & Sons, Ltd, 2008) Chap. 6, pp. 497–584, https://onlinelibrary.wiley.com/doi/pdf/ 10.1002/9780470259498.ch6 .
  * Underwood and Reid (2000) Jonathan G. Underwood and Katharine L. Reid, “Time-resolved photoelectron angular distributions as a probe of intramolecular dynamics: Connecting the molecular frame and the laboratory frame,” The Journal of Chemical Physics 113, 1067–1074 (2000), https://doi.org/10.1063/1.481918 .
  * Box and Tiao (2011) G.E.P. Box and G.C. Tiao,  _Bayesian Inference in Statistical Analysis_, Wiley Classics Library (Wiley, 2011).
  * github (2022) github, “Bigr github repository,”  (2022).
  * Hockett (2015) Paul Hockett, “General phenomenology of ionization from aligned molecular ensembles,” New Journal of Physics 17, 023069 (2015).
  * Foreman-Mackey  _et al._ (2013) Daniel Foreman-Mackey, David W. Hogg, Dustin Lang, and Jonathan Goodman, “emcee: The MCMC hammer,” Publications of the Astronomical Society of the Pacific 125, 306–312 (2013).
  * Herzberg (1966) Gerhard Herzberg,  _Electronic Spectra and electronic structure of Polyatomic Molecules_ (Nostrand, 1966).
  * Hastings (1970) W. K. Hastings, “Monte Carlo sampling methods using Markov chains and their applications,” Biometrika 57, 97–109 (1970), https://academic.oup.com/biomet/article-pdf/57/1/97/23940249/57-1-97.pdf .
  * Wolf  _et al._ (2019) T. J. A. Wolf, D. M. Sanchez, J. Yang, R. M. Parrish, J. P. F. Nunes, M. Centurion, R. Coffee, J. P. Cryan, M. Gühr, K. Hegazy, A. Kirrander, R. K. Li, J. Ruddock, X. Shen, T. Vecchione, S. P. Weathersby, P. M. Weber, K. Wilkin, H. Yong, Q. Zheng, X. J. Wang, M. P. Minitti, and T. J. Martínez, “The photochemical ring-opening of 1,3-cyclohexadiene imaged by ultrafast electron diffraction,” Nature Chemistry 11, 504–509 (2019).
  * Baskin and Zewail (2006) J. Spencer Baskin and Ahmed H. Zewail, “Oriented ensembles in ultrafast electron diffraction,” ChemPhysChem 7, 1562–1574 (2006), https://chemistry-europe.onlinelibrary.wiley.com/doi/pdf/ 10.1002/cphc.200600133 .
  * Parrish and Martinez (2019) Robert Parrish and Todd Martinez, “Ab initio computation of rotationally-averaged pump-probe x-ray and electron diffraction signals,” Journal of Chemical Theory and Computation 15 (2019), 10.1021/acs.jctc.8b01051.
  * Ma  _et al._ (2020) Lingyu Ma, Haiwang Yong, Joseph D. Geiser, Andrés Moreno Carrascosa, Nathan Goff, and Peter M. Weber, “Ultrafast x-ray and electron scattering of free molecules: A comparative evaluation,” Structural Dynamics 7, 034102 (2020), https://doi.org/10.1063/4.0000010 .
  * Makhija  _et al._ (2021) Varun Makhija, Andrey E. Boguslavskiy, Ruaridh Forbes, Kevin Veyrinas, Iain Wilkinson, Rune Lausten, Michael S. Schuurman, Edward R. Grant, and Albert Stolow, “A quantum molecular movie: polyad predissociation dynamics in the vuv excited 3p$\sigma$2$\sigma$u state of no2,” Faraday Discuss. 228, 191–225 (2021).
  * Felker and Zewail (1984) Peter M. Felker and Ahmed H. Zewail, “Direct observation of nonchaotic multilevel vibrational energy flow in isolated polyatomic molecules,” Phys. Rev. Lett. 53, 501–504 (1984).
  * Yang  _et al._ (2020) Jie Yang, Xiaolei Zhu, J. Pedro F. Nunes, Jimmy K. Yu, Robert M. Parrish, Thomas J. A. Wolf, Martin Centurion, Markus Gühr, Renkai Li, Yusong Liu, Bryan Moore, Mario Niebuhr, Suji Park, Xiaozhe Shen, Stephen Weathersby, Thomas Weinacht, Todd J. Martinez, and Xijie Wang, “Simultaneous observation of nuclear and electronic dynamics by ultrafast electron diffraction,” Science 368, 885–889 (2020), https://www.science.org/doi/pdf/10.1126/science.abb2235 .
  * Yong  _et al._ (2020) Haiwang Yong, Nikola Zotev, Jennifer M. Ruddock, Brian Stankus, Mats Simmermacher, Andrés Moreno Carrascosa, Wenpeng Du, Nathan Goff, Yu Chang, Darren Bellshaw, Mengning Liang, Sergio Carbajo, Jason E. Koglin, Joseph S. Robinson, Sébastien Boutet, Michael P. Minitti, Adam Kirrander, and Peter M. Weber, “Observation of the molecular response to light upon photoexcitation,” Nature Communications 11, 2157 (2020).
  * Zare (1991) Richard N. Zare,  _Angular Momentum: Understanding Spatial Aspects in Chemistry and Physics_ (Wiley-Interscience, 1991).
  * Poulsen  _et al._ (2004) Mikael D Poulsen, Emmanuel Peronne, Henrik Stapelfeldt, Christer Z Bisgaard, Simon S Viftrup, Edward Hamilton, and Tamar Seideman, “Nonadiabatic alignment of asymmetric top molecules: Rotational revivals,” The Journal of chemical physics 121, 783–791 (2004).
  * Underwood  _et al._ (2005) Jonathan G Underwood, Benjamin J Sussman, and Albert Stolow, “Field-free three dimensional molecular axis alignment,” Physical review letters 94, 143002 (2005).
  * Rouzée  _et al._ (2006) A Rouzée, S Guérin, Vincent Boudon, B Lavorel, and O Faucher, “Field-free one-dimensional alignment of ethylene molecule,” Physical Review A 73, 033418 (2006).
  * Takemoto and Yamanouchi (2008) Norio Takemoto and Kaoru Yamanouchi, “Fixing chiral molecules in space by intense two-color phase-locked laser fields,” Chemical Physics Letters 451, 1–7 (2008).
  * Ohshima and Hasegawa (2010) Yasuhiro Ohshima and Hirokazu Hasegawa, “Coherent rotational excitation by intense nonresonant laser fields,” International Reviews in Physical Chemistry 29, 619–663 (2010).
  * Pabst  _et al._ (2010b) Stefan Pabst, Phay J Ho, and Robin Santra, “Computational studies of x-ray scattering from three-dimensionally-aligned asymmetric-top molecules,” Physical Review A 81, 043425 (2010b).
  * Makhija  _et al._ (2012) Varun Makhija, Xiaoming Ren, and Vinod Kumarappan, “Metric for three-dimensional alignment of molecules,” Physical Review A 85, 033425 (2012).
  * Koch  _et al._ (2019) Christiane P Koch, Mikhail Lemeshko, and Dominique Sugny, “Quantum control of molecular rotation,” Reviews of Modern Physics 91, 035005 (2019).
  * Lin  _et al._ (2020) Kang Lin, Ilia Tutunnikov, Junyang Ma, Junjie Qiang, Lianrong Zhou, Olivier Faucher, Yehiam Prior, Ilya Sh Averbukh, and Jian Wu, “Spatiotemporal rotational dynamics of laser-driven molecules,” Advanced Photonics 2, 024002 (2020).
  * Rouzée  _et al._ (2008) A. Rouzée, S. Guérin, O. Faucher, and B. Lavorel, “Field-free molecular alignment of asymmetric top molecules using elliptically polarized laser pulses,” Phys. Rev. A 77, 043412 (2008).
  * Rouzée  _et al._ (2006) A. Rouzée, S. Guérin, V. Boudon, B. Lavorel, and O. Faucher, “Field-free one-dimensional alignment of ethylene molecule,” Phys. Rev. A 73, 033418 (2006).
  * Hamilton  _et al._ (2005) Edward Hamilton, Tamar Seideman, Tine Ejdrup, Mikael D. Poulsen, Christer Z. Bisgaard, Simon S. Viftrup, and Henrik Stapelfeldt, “Alignment of symmetric top molecules by short laser pulses,” Phys. Rev. A 72, 043402 (2005).
  * Holmegaard  _et al._ (2007) Lotte Holmegaard, Simon S. Viftrup, Vinod Kumarappan, Christer Z. Bisgaard, Henrik Stapelfeldt, Edward Hamilton, and Tamar Seideman, “Control of rotational wave-packet dynamics in asymmetric top molecules,” Phys. Rev. A 75, 051403 (2007).
  * Ren  _et al._ (2012) Xiaoming Ren, Varun Makhija, and Vinod Kumarappan, “Measurement of field-free alignment of jet-cooled molecules by nonresonant femtosecond degenerate four-wave mixing,” Phys. Rev. A 85, 033405 (2012).
  * Seideman (2001) Tamar Seideman, “On the dynamics of rotationally broad, spatially aligned wave packets,” The Journal of Chemical Physics 115, 5965–5973 (2001), https://doi.org/10.1063/1.1400131 .
  * Stapelfeldt and Seideman (2003) Henrik Stapelfeldt and Tamar Seideman, “Colloquium: Aligning molecules with strong laser pulses,” Rev. Mod. Phys. 75, 543–557 (2003).
  * Liu  _et al._ (2007) Han Liu, John Lafferty, and Larry Wasserman, “Sparse nonparametric density estimation in high dimensions using the rodeo,” in  _Proceedings of the Eleventh International Conference on Artificial Intelligence and Statistics_, Proceedings of Machine Learning Research, Vol. 2, edited by Marina Meila and Xiaotong Shen (PMLR, San Juan, Puerto Rico, 2007) pp. 283–290.
  * Cumpson and Seah (1992) P. J. Cumpson and M. P. Seah, “Random uncertainties in aes and xps: I: Uncertainties in peak energies, intensities and areas derived from peak synthesis,” Surface and Interface Analysis 18, 345–360 (1992).

Supplemental Information

##  VIII Calculating the axis distribution moments

We first address how to calculate the axis distribution moments (ADMs) for a
simple linear symmetric top, like N2O. A linear rotor has two unique principal
components of inertia, where the single unique moment is much much smaller
than the other two equal components. In the rigid rotor approximation,

$\displaystyle\bra{jm}\ket{\theta,\phi}$
$\displaystyle=Y_{j}^{m}(\theta,\phi)$ (S1) $\displaystyle E_{jm}$
$\displaystyle=Bj(j+1)$ (S2) $\displaystyle B$
$\displaystyle=\frac{\hbar^{2}}{2I}$ (S3)

where $B$ is the rotational constant for N2O and $I$ is the moment of inertia.
Before the alignment pulse ($t<0$) each molecule is in a single rotational
eigenstate $\ket{j,m}$. The alignment pulse launches a rotational wavepacket
by introducing a rotational coherence between eigenstates

$\ket{\psi^{(i)}(t)}=\sum_{j,m}c_{j_{i}m_{i}jm}(t)\ket{j,m}$ (S4)

where $j_{i}$ and $m_{i}$ label the initial ($t<0$) rotational eigenstate for
a single molecule. Our ensemble, however, is made of many molecules in
different rotational eigenstates that are populated according to the Boltzmann
Distribution. This thermal ensemble is represented by the density matrix

$\displaystyle\rho(t)$
$\displaystyle=\sum_{i}p_{i}\ket{\psi^{(i)}(t)}\bra{\psi^{(i)}(t)}$ (S5)
$\displaystyle p_{i}$ $\displaystyle=\frac{\exp{-E_{j_{i}m_{i}}/(k_{b}T)}}{Z}$
(S6)

where we sum over the initial $\ket{j_{i},m_{i}}$ states, $Z$ is the partition
function, $k_{b}$ is the Boltzmann, and $T$ is the temperature.

The ADMs is a sparse basis that encodes the angular relationship between the
lab frame and the distribution of molecular frame. More specifically, the ADMs
$(\mathcal{A}^{l}_{m}(t))$ are the projection of the ensemble anisotropy
(distribution of molecular frames with respect to the lab frame) onto the
Spherical Harmonics basis.

$P(\theta,\phi,t)=\sum_{l,m}\mathcal{A}^{l}_{m}(t)Y_{l}^{m}(\theta,\phi)=\left|\bra{\psi}\ket{\theta,\phi}\right|^{2}$
(S7)

Here $\theta$ and $\phi$ are the polar and azimuthal angles, respectively,
measured from the lab frame to the molecular frame and $P(\theta,\phi,t)$ is
the probability of finding a molecule oriented at these angles. The molecule’s
smallest moment of inertia defines the molecular frame $\hat{\mathbf{z}}$
axis, along the rotor. The laser polarization defines the lab frame z axis
($\hat{\textbf{z}}$). We calculate the ADMs by

$\begin{split}\mathcal{A}^{l}_{m}(t)=&\frac{2l+1}{4\pi}\Tr(\rho Y^{m}_{l})\\\
=&\frac{2l+1}{4\pi}\sum_{i}p_{i}\sum_{j_{1},m_{1}}\sum_{j_{2},m_{2}}c^{*}_{j_{i}m_{i}j_{1}m_{1}}(t)c_{j_{i}m_{i}j_{2}m_{2}}(t)\\\
&\hskip 11.38109pt\times\int
Y_{j_{1}}^{*m_{1}}Y_{l}^{m}Y_{j_{2}}^{m_{2}}\sin\theta d\theta d\phi\\\
=&\frac{2l+1}{4\pi}\sum_{i}p_{i}\sum_{j_{1},m_{1}}\sum_{j_{2},m_{2}}c^{*}_{j_{i}m_{i}j_{1}m_{1}}(t)c_{j_{i}m_{i}j_{2}m_{2}}(t)\\\
&\hskip 11.38109pt\times\sqrt{\frac{(2j_{1}+1)(2l+1)(2j_{2}+1)}{4\pi}}\\\
&\hskip 11.38109pt\times\begin{pmatrix}j_{1}&l&j_{2}\\\
0&0&0\end{pmatrix}\begin{pmatrix}j_{1}&l&j_{2}\\\
-m_{1}&m&m_{2}\end{pmatrix}.\end{split}$ (S8)

Asymmetric rigid rotors have three unique principal axes with $A\neq B\neq C$,
such that $I_{\text{A}}<I_{\text{B}}<I_{\text{C}}$. As a result, they have a
fundamentally different energy level structure. In general, the energy
eigenvalues may be determined analytically for each $J$ using the $D_{2}$
symmetry group of the rigid rotor Hamiltonian. This renders the Hamiltonian
matrix in the $\ket{jmk}$ symmetric top basis block diagonal Zare (1991). Here
$k$ is the angular momentum quantum number corresponding to the projection of
the angular momentum on the molecular frame $\hat{\textbf{z}}$. Writing the
eigenstates in this basis yields,

$\left|jm\tau\right>=\sum_{k}C_{jmk}\left|jmk\right>.$ (S9)

The asymmetric top eigenstates $\ket{jm\tau}$ each correspond to an energy
eigenvalue $E_{jm\tau}$, and the spacing between eigenstates determines the
field-free evolution of the rotational wavepacket excited by the alignment
pulse from an initial state rotational state $(i)$,

$\ket{\psi^{(i)}(t)}=\sum_{jm\tau}c_{jm\tau}\exp\\{\frac{-iE_{jm\tau}t}{\hbar}\bigg{\\}}\ket{jm\tau}.$
(S10)

The coefficients $c_{jm\tau}$ are determined by solving the TDSE for the
asymmetric rigid rotor in a non-resonant time dependent electric field. The
field-matter interaction is typically mediated by the molecular
polarizability, resulting in a series of Raman Transitions. Such a calculation
has been detailed by several authors Poulsen  _et al._ (2004); Underwood  _et
al._ (2005); Rouzée  _et al._ (2006); Takemoto and Yamanouchi (2008); Ohshima
and Hasegawa (2010); Pabst  _et al._ (2010b); Makhija  _et al._ (2012); Koch
_et al._ (2019); Lin  _et al._ (2020), so we do not discuss it here. The
density matrix $\rho_{jm\tau}^{j^{\prime}m^{\prime}\tau^{\prime}}(t)$ can then
be determined as discussed above for the linear molecule. Finally, the ADMs
can be calculated from the density matrix transformed into the $\ket{jmk}$
basis as follows,

$\begin{split}\mathcal{A}^{K}_{QS}(t)=&\frac{2K+1}{8\pi^{2}}\Tr(\rho
D^{K}_{QS})\\\
=&\frac{2K+1}{8\pi^{2}}\sum_{j,m,k}\sum_{j^{\prime},m^{\prime},k^{\prime}}\rho_{jmk}^{j^{\prime}m^{\prime}k^{\prime}}(t)\\\
&\hskip 11.38109pt\times\sqrt{(2j+1)(2j^{\prime}+1)}(-1)^{m-k}\\\ &\hskip
11.38109pt\times\begin{pmatrix}j&j^{\prime}&K\\\
-m&m^{\prime}&Q\end{pmatrix}\begin{pmatrix}j&j^{\prime}&K\\\
-k&k^{\prime}&S\end{pmatrix}.\end{split}$ (S11)

Figure S1: The axis distribution moments vary as a function of pump fluence
and ensemble temperature. The left column varies the pump fluence with a
constant temperature of 25 K. The right column varies the temperature with a
constant fluence of 1 J/cm2.

##  IX Anisotropy Derivation

Figure S2: To relate the ensemble anisotropy to the lab frame observable, we
must rotate between the pair-wise, molecular, and lab frames. Using NO2, we
illustrate these three reference frames. Here we highlight a single NO bond
and moved the oxygen to the origin, as is done in the derivation. On the left
is the lab frame (lf) where the laser polarization and the direction of the
probe beam define $\hat{\mathbf{z}}^{(\text{lf})}$ and
$\hat{\mathbf{y}}^{(\text{lf})}$ respectively. In the molecular frame (mf),
shown in the center, is defined by the molecule’s principal moments of inertia
where the molecular A, B, and C axis define $\hat{\mathbf{z}}^{(\text{mf})}$,
$\hat{\mathbf{y}}^{(\text{mf})}$, and $\hat{\mathbf{x}}^{(\text{mf})}$. Each
pair of two atoms has its own pair-wise frame (pf) where
$\hat{\mathbf{z}}^{(\text{pf})}$ is defined by $\mathbf{r}_{\mu\nu}$.

We now show how deterministic anisotropy allows one to access molecular frame
geometric information by coupling the lab and molecular frame. Using the
Independent Atom model (IAM), the x-ray, or electron, diffraction intensity
from a single molecule is given by

$\displaystyle\begin{split}I(\textbf{q})&=\mathcal{I}\Big{(}\sum_{\mu}|f_{\mu}(q)|^{2}\\\
+&\sum_{\mu,\nu:\mu\neq\nu}\text{Re}\\{f_{\mu}(q)f^{*}_{\nu}(q)\exp\left(i\textbf{q}\cdot\left(\textbf{r}_{\mu}-\textbf{r}_{\nu}\right)\right)\\}\Big{)}\end{split}$
(S12) $\displaystyle\mathcal{I}$
$\displaystyle=\left\\{\begin{array}[]{cl}I_{0}&\text{: x-ray diffraction}\\\
\frac{I_{0}}{R^{2}}&\text{: electron diffraction}\end{array}\right.$ (S15)

where $\mathcal{I}$ is a scaling coefficient, $I_{0}$ is the initial intensity
of the probe, $R$ is the distance between the sample and where the electron
was detected, $\textbf{r}_{\mu}$ is the position of the $\mu^{\text{th}}$
atom, and $f_{\mu}(q)$ is either the electron scattering amplitude or x-ray
form factor of the $\mu^{\text{th}}$ atom. Here $\mathbf{q}$ is the momentum
transfer imparted on either the electron or x-ray after scattering from the
molecule. In the case of x-ray scattering, we assume one has already removed
the anisotropic effects from Thomson Scattering. The difference in the
scattered x-ray and electron wave functions accounts for the $R^{-2}$ factor
in $\mathcal{I}$. The first term $\left(\sum_{\mu}|f_{\mu}(q)|^{2}\right)$ is
independent of the molecule’s geometry and is referred to as the atomic
scattering contribution. The second term depends on the pair-wise distances of
atoms and is known as the molecular diffraction.

Our objective is to represent the lab frame diffraction pattern, parameterized
by the momentum transfer $q=\left|\textbf{q}\right|$ and the detector’s
azimuthal angle $\theta^{(\text{d})}$, in terms of the molecular frame pair-
wise distances and angles
$\Delta\textbf{r}_{\mu\nu}=\textbf{r}_{\mu}-\textbf{r}_{\nu}=\left[\Delta
r_{\mu\nu},\theta^{(\text{mf})}_{\mu\nu},\phi^{(\text{mf})}_{\mu\nu}\right]$.
This derivation focuses on a single $(\mu,\nu)$ pair from the molecular
diffraction sum in Eq. S12, where the $\nu^{\text{th}}$ atom defines the
origin as we rotate between the lab and various body reference frames. Figure
S2 illustrates these various frames serving as an intuitive guide, with the
$\nu^{\text{th}}$ atom translated to the origin. Such translations are allowed
since they cancel in the $\Delta\textbf{r}_{\mu\nu}$ term. For our rotations,
we use the conventions in Ref. Zare (1991). Unless otherwise stated, $\theta$
and $\phi$ represent the polar and azimuthal angles, respectively, in a
spherical coordinate system.

We define the pair-wise frame (pf) such that
$\hat{\mathbf{z}}^{(\text{pf})}=\Delta\hat{\textbf{r}}_{\mu\nu}$, again
emphasizing we translate the molecule such that the $\nu^{\text{th}}$ atom is
at the origin. The pair-wise frame is shown in Fig. S2c. The exponential term
in Eq. S12 is rewritten using the partial wave expansion

$\begin{split}\exp\left(i\textbf{q}\cdot\Delta\textbf{r}_{\mu\nu}\right)=\sum_{l}i^{l}&\left(2l+1\right)j_{l}(q\Delta
r_{\mu\nu})\\\ &\times
Y^{0}_{l}\left(\theta^{(\text{pf})}_{q},\phi^{(\text{pf})}_{q}\right)\end{split}$
(S16)

Here, $j_{l}(q\Delta r_{\mu\nu})$ are the spherical Bessel functions of the
first kind, $Y_{l}^{m}(\theta^{(\text{pf})}_{q},\phi^{(\text{pf})}_{q})$ are
spherical harmonics, and
$\left(\theta^{(\text{pf})}_{q},\phi^{(\text{pf})}_{q}\right)$ are the polar
and azimuthal angles that define $\mathbf{q}$ in the pair-wise frame. In the
above equation, we extract the labelled pair-wise distance $\Delta
r_{\mu\nu}$: one of our parameters of interest.

The molecular frame (mf) is defined by the molecule’s principal moments of
inertia, here the $\hat{\mathbf{z}}^{(\text{mf})}$,
$\hat{\mathbf{y}}^{(\text{mf})}$, and $\hat{\mathbf{x}}^{(\text{mf})}$ axes
correspond to the moments with increasing rotational inertia. Figure S2b shows
the molecular frame for NO2 with the nitrogen translated to the origin. We
passively rotate from the pair-wise frame into the molecular frame, shown in
Fig. S2 as green and orange, respectively.

$\displaystyle\exp\left(i\textbf{q}\cdot\Delta\textbf{r}_{\mu\nu}\right)$
$\displaystyle={}\sum_{l}i^{l}\left(2l+1\right)j_{l}(q\Delta
r_{\mu\nu})\sum_{m_{1}}\left[D^{l}_{m_{1}0}\left(\phi_{\mu\nu}^{(\text{mf})},\theta_{\mu\nu}^{(\text{mf})},0\right)\right]^{\dagger}Y^{m_{1}}_{l}\left(\theta^{(\text{mf})}_{q},\phi^{(\text{mf})}_{q}\right)$
(S17) $\displaystyle={}\sum_{l}i^{l}\sqrt{4\pi\left(2l+1\right)}j_{l}(q\Delta
r_{\mu\nu})\sum_{m_{1}}Y_{l}^{m_{1}}\left(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}\right)Y^{m_{1}}_{l}\left(\theta^{(\text{mf})}_{q},\phi^{(\text{mf})}_{q}\right)$
(S18)

The molecular frame angles $\phi_{\mu\nu}^{(\text{mf})}$ and
$\theta_{\mu\nu}^{(\text{mf})}$ define the orientation of
$\Delta\hat{\textbf{r}}_{\mu\nu}$, where $\chi_{\mu\nu}^{(\text{mf})}=0$ since
$\Delta\textbf{r}_{\mu\nu}$ is a vector. We stress the importance of these
molecular frame geometry angles
$(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})})$ as they are
needed, along with $\Delta r_{\mu\nu}$ to define a unique molecular geometry.
With PDF methods alone, one only has access to unlabeled $\Delta r_{\mu\nu}$
and generally cannot define a unique molecular geometry. These molecular frame
angles are the last two geometric parameters of interest.

To connect our molecular frame calculation to our measurement, we passively
rotate into the lab frame (lf). The lab frame $\hat{\textbf{z}}^{(\text{lf})}$
is defined as the polarization of the alignment laser
$(\hat{\bm{\varepsilon}})$, and $\hat{\textbf{y}}^{(\text{lf})}$ is along the
probe path and normal to the detector.

$\displaystyle\begin{split}\exp\left(i\textbf{q}\cdot\Delta\textbf{r}_{\mu\nu}\right)&={}\sum_{l}i^{l}\sqrt{4\pi\left(2l+1\right)}j_{l}(q\Delta
r_{\mu\nu})\\\ &\hskip
28.45274pt\times\sum_{m_{1}m_{2}}\left[D^{l}_{m_{2}m_{1}}\left(\phi_{\text{I}}^{(\text{lf})},\theta_{\text{I}}^{(\text{lf})},\chi_{\text{I}}^{(\text{lf})}\right)\right]^{\dagger}Y_{l}^{m_{1}}\left(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}\right)Y^{m_{2}}_{l}\left(\theta^{(\text{lf})}_{q},\phi^{(\text{lf})}_{q}\right)\end{split}$
(S19)
$\displaystyle\begin{split}&={}\sum_{l}i^{l}\sqrt{4\pi\left(2l+1\right)}j_{l}(q\Delta
r_{\mu\nu})\\\ &\hskip
28.45274pt\times\sum_{m_{1}m_{2}}(-1)^{m_{2}-m_{1}}D^{l}_{-m_{2}-m_{1}}\left(\phi_{\text{I}}^{(\text{lf})},\theta_{\text{I}}^{(\text{lf})},\chi_{\text{I}}^{(\text{lf})}\right)Y_{l}^{m_{1}}\left(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}\right)Y^{m_{2}}_{l}\left(\theta^{(\text{lf})}_{q},\phi^{(\text{lf})}_{q}\right)\end{split}$
(S20)

Here $\phi_{\text{I}}^{(\text{lf})}$, $\theta_{\text{I}}^{(\text{lf})}$, and
$\chi_{\text{I}}^{(\text{lf})}$ are the conventional Euler angles in the lab
frame that describe the orientation of the molecule’s principal moments of
inertia with respect to the lab frame. Finally, without loss of generality, we
take $m_{1}\rightarrow-m_{1}$ and $m_{2}\rightarrow-m_{2}$, and combine Eqs.
S20 and S12 to get the following general expression for the diffraction
intensity of a single molecule

$\displaystyle\begin{split}I(\textbf{q}&)={}\mathcal{I}\bigg{(}\sum_{\mu}|f_{\mu}(q)|^{2}+\sum_{\mu,\nu:\mu\neq\nu}\text{Re}\bigg{\\{}f_{\mu}(q)f^{*}_{\nu}(q)\sum_{l}i^{l}8\pi^{2}\sqrt{4\pi\left(2l+1\right)}j_{l}(q\Delta
r_{\mu\nu})\\\
&\times\sum_{m_{1},m_{2}}(-1)^{m_{1}-m_{2}}D^{l}_{m_{2},m_{1}}\left(\phi_{\text{I}}^{(\text{lf})},\theta_{\text{I}}^{(\text{lf})},\chi_{\text{I}}^{(\text{lf})}\right)Y_{l}^{-m_{1}}\left(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}\right)Y^{-m_{2}}_{l}\left(\theta^{(\text{lf})}_{q},\phi^{(\text{lf})}_{q}\right)\bigg{\\}}\bigg{)}.\end{split}$
(S21)

We have now expressed the measurable diffraction (Eq. S12) in terms of the
pair-wise molecular frame distances and angles, as well as the lab frame
angles $\left(\theta^{(\text{lf})}_{q},\phi^{(\text{lf})}_{q}\right)$ that
define $\mathbf{q}$.

In gas-phase diffraction experiments one measures an ensemble of molecules at
different orientations, alignments, and possibly differing geometries
depending on the populated rovibronic states. One samples that ensemble at a
variety of times relative to the evolving ensemble anisotropy, revealing the
following observable,

$\displaystyle\langle I(\textbf{q})\rangle(t)$
$\displaystyle={}\mathcal{I}\bigg{(}\sum_{\mu}|f_{\mu}(q)|^{2}+\sum_{\mu,\nu:\mu\neq\nu}\text{Re}\bigg{\\{}f_{\mu}(q)f^{*}_{\nu}(q)\sum_{l}i^{l}8\pi^{2}\sqrt{4\pi\left(2l+1\right)}$
(S22) $\displaystyle\times\sum_{m_{1},m_{2}}$
$\displaystyle(-1)^{m_{1}-m_{2}}Y^{-m_{2}}_{l}\left(\theta^{(\text{lf})}_{q},\phi^{(\text{lf})}_{q}\right)\expectationvalue{D^{l}_{m_{2}m_{1}}\left(\phi_{\text{I}}^{(\text{lf})},\theta_{\text{I}}^{(\text{lf})},\chi_{\text{I}}^{(\text{lf})}\right)j_{l}(q\Delta
r_{\mu\nu})Y_{l}^{-m_{1}}\left(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}\right)}{\Psi(t)}\bigg{\\}}\bigg{)}$

where $\Psi(t)$ is the molecular ensemble wavefunction that describes both the
rotational and vibronic dynamics of the system. This is the general expression
for the scattered intensity from the entire molecular ensemble in lab frame
coordinates, with the explicit dependence on the molecular frame atomic pair-
wise coordinates.

When considering typical experimental setups, we must employ common symmetries
and relate the calculated lab frame diffraction signal to the 2d detector,
parameterized by $q=|\textbf{q}|$ and $\theta^{(\text{d})}$. The detector lies
in the x-z plane of the lab frame where $\theta^{(\text{d})}$ is with respect
to $\hat{\mathbf{z}}^{(\text{lf})}$.

$\displaystyle\alpha$
$\displaystyle=2\sin^{-1}\left(\frac{q\lambda}{4\pi}\right)+\frac{\pi}{2}$
(S23) $\displaystyle\theta^{(\text{lf})}_{q}$
$\displaystyle=\cos^{-1}\left(\sin\left(\alpha\right)\cos\left(\theta^{(\text{d})}\right)\right)$
(S24) $\displaystyle\phi^{(\text{lf})}_{q}$
$\displaystyle=\tan^{-1}\left(\frac{\cos(\alpha)}{\sin(\alpha)\sin(\theta^{(\text{d})})}\right)$
(S25)

Where $\lambda$ is either the deBroglie wavelength of the electron probe, or
the x-ray wavelength, and $\alpha$ is the scattering angle rotated by $\pi/2$.
For the 3.7 MeV electron probe at the SLAC Ultrafast Electron Diffraction
facility Shen  _et al._ (2019) $\lambda=3.0\times 10^{-3}\text{\AA}$ and the
above relations simplify to

$\displaystyle\alpha$ $\displaystyle\approx\frac{\pi}{2}$
$\displaystyle\theta^{(\text{lf})}_{q}$
$\displaystyle\approx\theta^{(\text{d})}$
$\displaystyle\phi^{(\text{lf})}_{q}$ $\displaystyle\approx 0$

For x-ray diffraction at ~10 keV this expression does not simplify due to
larger x-ray scattering angles. Often, one uses a linearly polarized alignment
pump pulse which induces cylindrical symmetry in the ensemble rotation wave
packet, which results in $m_{2}=0$. Equation S21 is derived for an asymmetric
top, for a symmetric top there is symmetry about the molecular frame z axis,
which sets $m_{1}=0$.

It is difficult to extract $\Delta\textbf{r}_{\mu\nu}$ from Eq. S22 in its
current form since rovibronic coupling may affect the time dependent
anisotropy. With rovibronic coupling, to calculate the ensemble anisotropy we
may be required to simulate the excited state with the complex excited state
simulations we do not want to rely on. This coupling, therefore, may render
the anisotropy calculation too difficult. Instead, we consider two methods to
separate the ensemble anisotropy and the molecular frame pair-wise terms,
where we assume the molecular geometry is rigid over the measurement period.
We decompose the ensemble anisotropy into the Axis Distribution Moments (ADMs)
by projecting the ensemble of molecular frame orientations, with respect to
the lab frame (Fig. S2a), onto the Wigner D matrices,

$\displaystyle\mathcal{A}^{l}_{mk}(t)$
$\displaystyle=\frac{2l+1}{8\pi^{2}}\expectationvalue{D^{l}_{mk}\left(\phi_{\text{I}}^{(\text{lf})},\theta_{\text{I}}^{(\text{lf})},\chi_{\text{I}}^{(\text{lf})}\right)}{\Psi(t)}$
(S26) $\displaystyle\left.\mathcal{A}^{l}_{mk}(t)\right|_{\text{rigid}}$
$\displaystyle=\frac{2l+1}{8\pi^{2}}\expectationvalue{D^{l}_{mk}\left(\phi_{\text{I}}^{(\text{lf})},\theta_{\text{I}}^{(\text{lf})},\chi_{\text{I}}^{(\text{lf})}\right)}{\Psi_{\text{rigid}}(t)}.$
(S27)

Simulations of the rotational wavefunction for rigid symmetric and rigid
asymmetric tops Rouzée  _et al._ (2008, 2006); Hamilton  _et al._ (2005);
Holmegaard  _et al._ (2007); Ren  _et al._ (2012); Seideman (2001);
Stapelfeldt and Seideman (2003) produce good agreement with measured alignment
signatures. To extract $\Delta\textbf{r}_{\mu\nu}$ from Eq. S22 we consider
two approximations: the typical rigid rotor approximation and a separation of
time scales.

###  IX.1 Rigid Rotor Approximation

We first consider the rigid rotor approximation, which assumes the molecular
geometry is constant throughout the rotational dynamics. This allows us to
take the expectation value of the molecular geometry, the molecular frame
terms, with respect to the ground state structure at $t=0$. We may also
calculate the ADMs with respect to the ground state geometry, which allows us
to separate the ADMs from the molecular frame terms

$\displaystyle\langle I(\mathbf{q})\rangle_{\text{rigid}}(t)={}$
$\displaystyle\mathcal{I}\bigg{(}\sum_{\mu}|f_{\mu}(q)|^{2}+\sum_{\mu,\nu:\mu\neq\nu}\text{Re}\bigg{\\{}f_{\mu}(q)f^{*}_{\nu}(q)\sum_{l}i^{l}8\pi^{2}\sqrt{\frac{4\pi}{\left(2l+1\right)}}$
(S28) $\displaystyle\times\sum_{m_{1},m_{2}}$
$\displaystyle(-1)^{m_{1}-m_{2}}Y^{-m_{2}}_{l}\left(\theta^{(\text{lf})}_{q},\phi^{(\text{lf})}_{q}\right)\expectationvalue{j_{l}(q\Delta
r_{\mu\nu})Y_{l}^{-m_{1}}\left(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}\right)}{\Psi(0)}\left.\mathcal{A}^{l}_{m_{2}m_{1}}(t)\right|_{\text{rigid}}\bigg{\\}}\bigg{)}.$

This approximation is useful when investigating the ground state structure of
a molecule or when the change in the molecule’s geometry has a negligible
impact on the moments of inertia.

###  IX.2 A Separation of Timescales Approximation for Excited State Dynamics

The second approximation is a separation of time scales between the rotational
and vibronic dynamics. The anisotropy signature, $\mathcal{A}^{l}_{mk}(t)$,
lasts of order one to tens of picoseconds for molecules with a few to tens of
atoms, respectively. When the vibration or dissociation occurs on a much
faster timescale than the change in anisotropy, we can calculate the
rotational dynamics with respect to the known ground state geometry rather
than with the unknown excited state geometry. The disparity in timescales is
very common, and this approximation is analogous to the Born-Oppenheimer
approximation.

We first consider the more general case of a double pump pulse experiment that
first induces a rotational wavepacket and then launches a vibronic wavepacket.
The first pulse increases the ensemble anisotropy and consequently the number
of $C_{lmk}(q,t)$ coefficients. The second pulse further mixes the rotational
states while exciting vibronic modes. Let $\tau$ denote the arrival time of
the second vibration inducing pulse after the first rotation inducing pulse,
and $t$ is the elapsed time after the second pump pulse.

Our initial state is considered to be a thermal ensemble in the vibronic
ground state. We describe a single initial state ($i$) in the Born-Oppenheimer
basis as

$\ket{\psi^{(i)}(0)}=\ket{J_{i}M_{i}K_{i}}\ket{0}$ (S29)

prior to any pulses. After the alignment pulse, and before the vibration
inducing pulse, our coherent rotational state evolves as

$\ket{\psi^{(i)}(\tau)}=\sum_{JMK}c^{(i)}_{JMK}(\tau)\ket{JMK}\ket{0}.$ (S30)

The vibration pump pulse induces the excited state dynamics, while the
photon’s angular momentum mixes the rotational states. We project the
vibronically excited state onto the Born-Oppenheimer basis,

$\ket{\Psi^{(i)}(t,\tau)}=\sum_{n}\sum_{J_{n}M_{n}K_{n}}\ket{J_{n}M_{n}K_{n}}\ket{n}\bra{J_{n}M_{n}K_{n}n}\ket{\Psi^{(i)}(t,\tau)}$
(S31)

where the vibronic and rotational states are mixed by the ground state dipole
moment and its orientation, respectively. To calculate the coefficients we
apply time dependent perturbation theory and assumed an impulsive excitation

$\displaystyle\bra{J_{n}M_{n}K_{n}n}\ket{\Psi^{(i)}(t,\tau)}=$
$\displaystyle\sum_{\gamma}\matrixelement{J_{n}M_{n}K_{n}n}{D^{1*}_{0\gamma}\mu^{1}_{\gamma}}{\psi^{(i)}(t,\tau)}\frac{-i}{\hbar}\int_{0}^{\infty}E_{0}(t^{\prime})e^{-i\Delta
Et^{\prime}/\hbar}dt^{\prime}$ (S32) $\displaystyle=$
$\displaystyle\tilde{E}\sum_{JMK}c^{(i)}_{JMK}(\tau+t)\sum_{\gamma}A(J_{n},J;K_{n},K;M\gamma)\bra{n}\mu^{1}_{\gamma}\ket{0}$
(S33) $\displaystyle A(J_{n},J;K_{n},K;M\gamma)\equiv$
$\displaystyle\bra{J_{n}MK_{n}}D^{1*}_{0\gamma}\ket{JMK}$ (S34)
$\displaystyle=$
$\displaystyle\sqrt{(2J+1)(2J_{n}+1)}(-1)^{\gamma+K-M}\begin{pmatrix}J&1&J_{n}\\\
-M&0&M\end{pmatrix}\begin{pmatrix}J&1&J_{n}\\\ -K&-\gamma&K_{n}\end{pmatrix}.$
(S35)

where $\mu^{1}_{\gamma}$ is the spherical tensor of the transition dipole
moment operator, $E_{0}(t^{\prime})$ is the electric field of the vibration
inducing pulse, $\Delta E$ is the energy difference between the initial
rotational state and the excited state, and

$\tilde{E}={}\int_{0}^{\infty}E_{0}(t^{\prime})e^{-i\Delta
Et^{\prime}/\hbar}dt^{\prime}.$ (S36)

At this point we have only considered a single initial rotational state in our
thermal ensemble. By summing over initial thermal states, weighting them by
the Boltzmann distribution, the complete rovibronic state immediately after
the vibration excitation pulse is given by

$\displaystyle\ket{\Psi(0,\tau)}$
$\displaystyle=\sum_{nJ_{n}MK_{n}}\ket{nJ_{n}MK_{n}}X^{n}_{J_{n}K_{n}}(M;0,\tau)$
(S37) $\displaystyle X^{n}_{J_{n}K_{n}}(M;t,\tau)$
$\displaystyle=\sum_{JK}\left(\frac{1}{Z}\sum_{i}c^{(i)}_{JMK}(\tau+t)p^{(i)}\right)\sum_{\gamma}\bra{n}\mu^{1}_{\gamma}\ket{0}A(J_{n},J;K_{n},K;M\gamma).$
(S38)

Here $Z$ is the partition function of the initial rotational states and
$p^{(i)}$ is Boltzmann weighting factor. Since we are interested in the time
dynamics of the vibronic state, we apply the time translation operator to Eq.
S37.

$\displaystyle\ket{\Psi(t,\tau)}$
$\displaystyle=\sum_{J_{n}MK_{n}}\sum_{n}X^{n}_{J_{n}K_{n}}(M;t,\tau)\ket{\psi^{n}_{\text{el-
vib}}(t)}\ket{J_{n}MK_{n}}$ (S39) $\displaystyle\ket{\psi^{n}_{\text{el-
vib}}(t)}$ $\displaystyle=\hat{U}(t)\ket{n}$ (S40)

Now that we’ve calculated the time dependent rovibronic state in the Born-
Oppenheimer basis, we must apply it to our measurement. We do this by taking
the expectation value of our diffraction observable (Eq. S22) with respect to
our new rovibronic system, Eq. S39. Here, we only look at the molecular
scattering term since it is the only term affected by $\ket{\Psi(t,\tau)}$ and
pull out the scattering amplitudes $f_{\mu}(q)$ due to the independent atom
approximation.

$\displaystyle\begin{split}\langle
I(\textbf{q})\rangle_{\text{sep}}^{(2)}(t,\tau)|_{\text{mol}}&={}\mathcal{I}\sum_{\mu,\nu:\mu\neq\nu}\text{Re}\bigg{\\{}f_{\mu}(q)f^{*}_{\nu}(q)\sum_{l}i^{l}8\pi^{2}\sqrt{4\pi\left(2l+1\right)}\sum_{m_{1},m_{2}}(-1)^{m_{1}-m_{2}}Y^{-m_{2}}_{l}\left(\theta^{(\text{lf})}_{q},\phi^{(\text{lf})}_{q}\right)\\\
\times&\expectationvalue{D^{l}_{m_{2}m_{1}}\left(\phi_{\text{I}}^{(\text{lf})},\theta_{\text{I}}^{(\text{lf})},\chi_{\text{I}}^{(\text{lf})}\right)j_{l}(q\Delta
r_{\mu\nu})Y_{l}^{-m_{1}}\left(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}\right)}{\Psi(t,\tau)}\bigg{\\}}\end{split}$
(S41)
$\displaystyle\begin{split}&={}\mathcal{I}\sum_{\mu,\nu:\mu\neq\nu}\text{Re}\bigg{\\{}f_{\mu}(q)f^{*}_{\nu}(q)\sum_{l}i^{l}8\pi^{2}\sqrt{4\pi\left(2l+1\right)}\sum_{m_{1},m_{2}}(-1)^{m_{1}-m_{2}}Y^{-m_{2}}_{l}\left(\theta^{(\text{lf})}_{q},\phi^{(\text{lf})}_{q}\right)\\\
\times&\sum_{nn^{\prime}}\tilde{\mathcal{A}}^{l}_{m_{1}m_{2}}(n,n^{\prime};t,\tau)\matrixelement{\psi^{n^{\prime}}_{\text{el-
vib}}(t)}{j_{l}(q\Delta
r_{\mu\nu})Y_{l}^{-m_{1}}\left(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}\right)}{\psi^{n}_{\text{el-
vib}}(t)}\bigg{\\}}\end{split}$ (S42)
$\displaystyle\tilde{\mathcal{A}}^{l}_{m_{1}m_{2}}(n,n^{\prime};t,\tau)$
$\displaystyle\equiv\sum_{J_{n}MK_{n}}\sum_{J_{n^{\prime}}^{\prime}M^{\prime}K_{n^{\prime}}^{\prime}}X^{n}_{J_{n}K_{n}}(M;t,\tau)X^{n^{\prime}*}_{J_{n^{\prime}}^{\prime}K_{n^{\prime}}^{\prime}}(M^{\prime};t,\tau)\bra{J_{n^{\prime}}^{\prime}M^{\prime}K_{n^{\prime}}^{\prime}}D^{l}_{m_{2}m_{1}}\ket{J_{n}MK_{n}}$
(S43)

We have again separated the ensemble anisotropy
$(\tilde{\mathcal{A}}^{l}_{m_{1}m_{2}}(n,n^{\prime};t,\tau))$ from the
molecular frame geometry term, which includes all the vibronic dynamics. The
modified ADMs, $\tilde{\mathcal{A}}^{l}_{m_{1}m_{2}}(n,n^{\prime};t,\tau)$,
are analogous to the original ADMs, but now include the coherent rotational
mixing from the vibronic inducing pulse. That is, each vibronic state will
have its own rotational coherence that must be accounted for when calculating
the ensemble anisotropy. Finally, plugging the molecular diffraction term (Eq.
S42) into the full diffraction expression we get

$\displaystyle\langle I(\textbf{q})\rangle_{\text{sep}}^{(2)}(t,\tau)$
$\displaystyle={}\mathcal{I}\bigg{(}\sum_{\mu}|f_{\mu}(q)|^{2}+\sum_{\mu,\nu:\mu\neq\nu}\text{Re}\bigg{\\{}f_{\mu}(q)f^{*}_{\nu}(q)\sum_{l}i^{l}8\pi^{2}\sqrt{4\pi\left(2l+1\right)}\sum_{m_{1},m_{2}}(-1)^{m_{1}-m_{2}}$
(S44) $\displaystyle\times
Y^{-m_{2}}_{l}\left(\theta^{(\text{lf})}_{q},\phi^{(\text{lf})}_{q}\right)\tilde{\mathcal{A}}^{l}_{m_{1}m_{2}}(n,n^{\prime};t,\tau)\matrixelement{\psi^{n^{\prime}}_{\text{el-
vib}}(t)}{j_{l}(q\Delta
r_{\mu\nu})Y_{l}^{-m_{1}}\left(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}\right)}{\psi^{n}_{\text{el-
vib}}(t)}\bigg{\\}}\bigg{)}.$

Due to the difference in timescales between the rotational and vibrational
dynamics, we further simplify Eq. S44. In its current form, Eq. S44 relies on
updating the ensemble anisotropy calculation as the geometry changes with
vibration. This requires us to know what the geometry will be at time $t$,
which is what we are ultimately trying to solve for. Instead, when the change
in ensemble anisotropy is negligible with respect to the timescale of the
vibration we can hold the anisotropy constant

$\displaystyle\tilde{\mathcal{A}}^{l}_{m_{1}m_{2}}(n,n^{\prime};\tau)$
$\displaystyle=\tilde{\mathcal{A}}^{l}_{m_{1}m_{2}}(n,n^{\prime};0,\tau)\approx\tilde{\mathcal{A}}^{l}_{m_{1}m_{2}}(n,n^{\prime};t,\tau)$
(S45) $\displaystyle\begin{split}\langle
I(\textbf{q})\rangle_{\text{sep}}^{(2)}(t,\tau)&\approx{}\mathcal{I}\bigg{(}\sum_{\mu}|f_{\mu}(q)|^{2}+\sum_{\mu,\nu:\mu\neq\nu}\text{Re}\bigg{\\{}f_{\mu}(q)f^{*}_{\nu}(q)\sum_{l}i^{l}8\pi^{2}\sqrt{4\pi\left(2l+1\right)}\sum_{m_{1},m_{2}}(-1)^{m_{1}-m_{2}}\\\
&\times
Y^{-m_{2}}_{l}\left(\theta^{(\text{lf})}_{q},\phi^{(\text{lf})}_{q}\right)\tilde{\mathcal{A}}^{l}_{m_{1}m_{2}}(n,n^{\prime};\tau)\matrixelement{\psi^{n^{\prime}}_{\text{el-
vib}}(t)}{j_{l}(q\Delta
r_{\mu\nu})Y_{l}^{-m_{1}}\left(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}\right)}{\psi^{n}_{\text{el-
vib}}(t)}\bigg{\\}}\bigg{)}.\end{split}$ (S46)

and the ensemble anisotropy and vibronic geometric dependence are completely
separable. We therefore continue to calculate the ensemble anisotropy with
respect to the ground state geometry while the molecular geometry has not
appreciably changed on the timescale of the rotational dynamics.

In some cases, a single pump pulse experiment is preferred over a two pump
pulse experiment when the setup is too difficult or when the anisotropy is
difficult to induce or measure. In such a case, we simplify Eq. S46 by setting
$\tau=0$. Here, the ensemble anisotropy is imprinted immediately after the
pulse by the interaction between the polarized laser and the excitation
dipole.

$\displaystyle\tilde{\mathcal{A}}^{l}_{m_{1}m_{2}}(n,n^{\prime})$
$\displaystyle=\tilde{\mathcal{A}}^{l}_{m_{1}m_{2}}(n,n^{\prime};0,0)\approx\tilde{\mathcal{A}}^{l}_{m_{1}m_{2}}(n,n^{\prime};t,0)$
(S47) $\displaystyle\begin{split}\langle
I(\textbf{q})\rangle_{\text{sep}}^{(2)}(t,\tau)&\approx{}\mathcal{I}\bigg{(}\sum_{\mu}|f_{\mu}(q)|^{2}+\sum_{\mu,\nu:\mu\neq\nu}\text{Re}\bigg{\\{}f_{\mu}(q)f^{*}_{\nu}(q)\sum_{l}i^{l}8\pi^{2}\sqrt{4\pi\left(2l+1\right)}\sum_{m_{1},m_{2}}(-1)^{m_{1}-m_{2}}\\\
&\times
Y^{-m_{2}}_{l}\left(\theta^{(\text{lf})}_{q},\phi^{(\text{lf})}_{q}\right)\tilde{\mathcal{A}}^{l}_{m_{1}m_{2}}(n,n^{\prime};t)\matrixelement{\psi^{n^{\prime}}_{\text{el-
vib}}(t)}{j_{l}(q\Delta
r_{\mu\nu})Y_{l}^{-m_{1}}\left(\theta_{\mu\nu}^{(\text{mf})},\phi_{\mu\nu}^{(\text{mf})}\right)}{\psi^{n}_{\text{el-
vib}}(t)}\bigg{\\}}\bigg{)}.\end{split}$ (S48)

Depending on the system, one may further improve this approximation by
calculating the ensemble dynamics with respect to a reference geometry for
$t>0$. In some cases, the vibronic transience may be on the timescale of the
rotational transience. Once Eqs S45 or S47 no longer hold at some time $t$
there are two options. Firstly, one can use only $C_{000}(q,t)$ which does not
rely on anisotropy and Eq. S44 will be exact. Secondly, one can continue
calculating $\tilde{\mathcal{A}}^{l}_{m_{1}m_{2}}(n,n^{\prime};t,\tau)$ with
respect to a reference geometry. For example, if one knows an excited state
geometry is similar to the ground state one can continue to use
$\tilde{\mathcal{A}}^{l}_{m_{1}m_{2}}(n,n^{\prime};t,\tau)$. One must prove
this through a priori knowledge or through the retrieved geometries at earlier
times. In the case that the dynamics do not deviate from some other known
geometry one may calculate the
$\tilde{\mathcal{A}}^{l}_{m_{1}m_{2}}(n,n^{\prime};t,\tau)$ with respect to
this geometry.

##  X Fitting for $B^{m}_{l}(q,t)$ and $C_{lmk}(q)$, and common mistakes

Our method relies heavily on two fitting procedures that will likely be the
most important steps of the analysis as they define the $C_{lmk}(q)$
coefficients and $\sigma_{lmk}(q)$. Below, we describe how one performs these
fits analytically by minimizing the $\chi^{2}$. These analytical methods,
however, will struggle to fit the measured time dependence with ADMs if there
is not enough anisotropy and/or there is poor SNR. We therefore highly
encourage one to validate the below, or their own, fitting procedure on their
expected simulated data. One can employ L1 regularization techniques to
improve these fits, but will likely have to use gradient descent since the
derivative of $|x|$ is undefined at $x=0$ and we do not know the sign of the
$B^{m}_{l}(q,t)$ and $C_{lmk}(q)$ coefficients. Gradient descent will be much
slower for numerous fits and should be used if the analytical approach is
insufficient.

To minimize the $\chi^{2}$ we cast it as a weighted linear regression problem
with the loss function

$\begin{split}L&=\frac{1}{2}\sum_{\mu}\frac{\left(\sum_{\nu}X_{\mu\nu}F_{\nu}-Y_{\mu}\right)^{2}}{\epsilon_{\mu}}\\\
&=\frac{1}{2}\left(\mathbf{XF}-\mathbf{Y}\right)^{T}\mathbf{W}\left(\mathbf{XF}-\mathbf{Y}\right).\end{split}$
(S49)

Here, $\mathbf{Y}$ is the data vector we wish to fit, the matrix $\mathbf{X}$
are the fit bases (features) that span the columns, $\mu$ sums over the
detector pixels, and $\nu$ sums over the fit bases. The bases are scaled by
the fit coefficients $\mathbf{F}$ and each data point’s contribution to the
fit is weighted by $\mathbf{W}$, where

$\displaystyle\epsilon_{\mu}$ $\displaystyle=\text{Var}\left(Y_{\mu}\right)$
(S50) $\displaystyle\mathbf{W}$
$\displaystyle=\begin{bmatrix}\frac{1}{\epsilon_{0}}&0&\dots&0\\\
0&\frac{1}{\epsilon_{1}}&&\\\ \vdots&&\ddots&\\\
0&&&\frac{1}{\epsilon_{N}}.\end{bmatrix}$ (S51)

Solving Eq. S49 for the optimal fit coefficients yields

$\mathbf{F}=\left(\mathbf{X}^{T}\mathbf{WX}\right)^{-1}\mathbf{X}^{T}\mathbf{WY}.$
(S52)

The first fitting procedure in this method is to fit the measured data,
$I\langle(\mathbf{q}(\theta^{(\text{d})}),t\rangle$, with the spherical
harmonics, $Y_{l}^{m}(\theta^{(\text{lf})}_{q},\phi^{(\text{lf})}_{q})$, to
retrieve the $B^{m}_{l}(q,t)$ coefficients. Where Eqs. S24 and S25 relate
$\theta^{(\text{lf})}_{q}$ and $\phi^{(\text{lf})}_{q}$ in terms of
$\theta^{(\text{d})}$. Although the spherical harmonics are orthonormal, this
orthonormality is broken by our finite sampling due to detector binning. To
account for this now nonzero overlap between different bases, we fit the
spherical harmonics to the data instead of projecting onto them. This is very
important at low $q$ where one often has the best SNR and the fewest bins to
resolve $\theta^{(\text{d})}$. We note that this can still be important even
for the isotropic component due to the Jacobian. In the limit that the number
of bins reaches infinity, the fit (Eq. S52) turns into a weighted projection.
We use the trapezoidal rule to increase the orthonormality of our binned
spherical harmonics

$\displaystyle X_{\mu\nu}$
$\displaystyle=\frac{1}{2}\left(Y_{l}^{\nu}(\theta^{(\text{lf})}_{q\mu},\phi^{(\text{lf})}_{q\mu})+Y_{l}^{\nu}(\theta^{(\text{lf})}_{q\mu+1},\phi^{(\text{lf})}_{q\mu+1})\right)$
(S53) $\displaystyle Y_{\mu}$ $\displaystyle=\frac{1}{2}\left(\langle
I(q,\theta^{(\text{d})}_{\mu},t)\rangle+\langle
I(q,\theta^{(\text{d})}_{\mu+1},t)\rangle\right)$ (S54)
$\displaystyle\epsilon_{\mu}$
$\displaystyle=\frac{1}{2}\left(\text{Var}\left(\langle
I(q,\theta^{(\text{d})}_{\mu},t)\rangle\right)+\text{Var}\left(\langle
I(q,\theta^{(\text{d})}_{\mu+1},t)\rangle\right)\right)$ (S55) $\displaystyle
F_{\nu}$ $\displaystyle=B_{l}^{\nu}(q,t).$ (S56)

In Eq. S52 we are summing over $\theta^{(\text{d})}$ ($\mu)$, but since
$Y_{l}^{m}(\theta^{(\text{lf})}_{q},\phi^{(\text{lf})}_{q})\propto
P_{l}^{m}(\cos(\theta^{(\text{lf})}_{q}))$ we must consider the
$\cos\left(\theta^{(\text{lf})}_{q}\right)$ Jacobian. This can be resolved in
two ways, by rebinning $\theta^{(\text{d})}$ in equally sized
$\cos(\theta^{(\text{lf})}_{q})$ bins, or by introducing the Jacobean into
$\mathbf{W}$. Since rebinning reduces our resolution, we alter the weight
matrix

$\displaystyle\widetilde{\theta^{(\text{lf})}_{q\mu}}$
$\displaystyle=\frac{1}{2}\left(\theta^{(\text{lf})}_{q\mu}+\theta^{(\text{lf})}_{q\mu+1}\right)$
(S57) $\displaystyle\Delta\theta^{(\text{lf})}_{q\mu}$
$\displaystyle=\theta^{(\text{lf})}_{q\mu+1}-\theta^{(\text{lf})}_{q\mu}$
(S58) $\displaystyle\mathbf{W}$
$\displaystyle=\mathbb{1}\begin{bmatrix}\frac{\sin\left(\widetilde{\theta^{(\text{lf})}_{q0}}\right)\Delta\theta^{(\text{lf})}_{q0}}{\epsilon_{0}}\vspace{1mm}\\\
\frac{\sin\left(\widetilde{\theta^{(\text{lf})}_{q1}}\right)\Delta\theta^{(\text{lf})}_{q1}}{\epsilon_{1}}\\\
\vdots\\\
\frac{\sin\left(\widetilde{\theta^{(\text{lf})}_{qN-1}}\right)\Delta\theta^{(\text{lf})}_{qN-1}}{\epsilon_{N_{bins}-1}}\end{bmatrix}$
(S59)

where $\mathbb{1}$ is the identity matrix.

Now we focus on retrieving the $C_{lmk}(q)$ coefficients by fitting the
$\mathcal{A}^{l}_{mk}(t)$ to the $B^{m}_{l}(q,t)$ coefficients. The
$\mathcal{A}^{l}_{mk}(t)$ are likely not orthogonal and may vary strongly in
their magnitude (L2 norm). Consequently, the fit results from ADMs bases with
larger magnitude can easily skew the results of other bases with lower
magnitude. These skews can completely ruin the fit for the lower magnitude
bases, while not being noticeable in the larger magnitude basis fit. This
issue is also mitigated, or exacerbated, by increasing or decreasing the SNR,
respectively. Another way to mitigate this issue is to add regularization
terms to Eq. S49, which will alter Eq. S52. One would ideally like to use L1
regularization for sparsity, but for the reasons mentioned above one would
need to use the gradient descent, which is much slower. To use Eqs. S49 and
S52, one must make the following alterations:

$\displaystyle X_{\mu\nu}$
$\displaystyle=\frac{1}{2}\left(\mathcal{A}^{l}_{m_{2}\nu}(t_{\mu})+\mathcal{A}^{l}_{m_{2}\nu}(t_{\mu+1})\right)$
(S60) $\displaystyle Y_{\mu}$
$\displaystyle=\frac{1}{2}\left(B_{l}^{m_{2}}(q,t_{\mu})+B_{l}^{m_{2}}(q,t_{\mu+1})\right)$
(S61) $\displaystyle\epsilon_{\mu}$
$\displaystyle=\frac{1}{2}\left(\text{Var}\left(B_{l}^{m_{2}}(q,t_{\mu})\right)+\text{Var}\left(B_{l}^{m_{2}}(q,t_{\mu+1})\right)\right)$
(S62) $\displaystyle F_{\nu}$ $\displaystyle=C_{lm_{2}\nu}(q)$ (S63)
$\displaystyle\mathbf{W}$