On the interpretation of time-resolved anisotropic diffraction patterns

In this paper, we review existing systematic treatments for the interpretation of anisotropic diffraction patterns from partially aligned symmetric top molecules. Such patterns arise in the context of time-resolved diffraction experiments. We calculate diffraction patterns for ground-state NaI excited with an ultraviolet laser. The results are interpreted with the help of a qualitative analytic model, and general recommendations on the analysis and interpretation of anisotropic diffraction patterns are given.


Introduction
Almost 100 years ago, Friedrich et al [1] published an article in which they described intensity maxima after sending an x-ray beam through crystals, and interpreted these in terms of diffraction of x-rays at the atomic sites. Ever since, x-ray diffraction has been an indispensable tool for determining the atomic structure of matter.
To a good approximation, x-ray diffraction patterns are simple functions of the atomic species and the relative positions of the nuclei with respect to each other; that is, they directly encode the nuclear geometry. This is a major advantage of diffraction over most other techniques (e.g. optical spectroscopy), which usually measure the properties of the electronic structure and infer the nuclear geometry indirectly from models or extensive electronic structure calculations.
For a chemist, the natural next step after resolving the structure of reactants and products in a chemical reaction would be to 'film' the reaction itself, which leads to time-resolved x-ray diffraction (TRXD). In this field, the commonly used technique is pump-probe difference diffraction. A pump laser initiates the chemical reaction, and after a specific delay time, the x-ray probe pulse produces a diffraction pattern. From this time-dependent diffraction pattern, a reference pattern of the unexcited sample is subtracted to yield a difference pattern; this technique removes all contributions that do not depend on the time-dependent nuclear geometry, such as incoherent scattering or the scattering from the major portion of the molecules that have not been excited. The difference patterns, or rather their encoded structures, can then be ordered by the pump-probe delay time to produce a 'molecular movie' of the chemical reaction. While the use of a pump laser means that only photochemical processes can be studied, this is to some extent unavoidable; for most reactions, only laser technology provides a fast enough clocking mechanism.
Unfortunately, so far there has been a distinct lack of suitable x-ray sources for ultrafast TRXD experiments. Existing synchrotron sources produce pulses of about 50 ps; this is fast enough for studying, for example, short-lived transient states in liquid (see [2,3] and references therein) or crystalline samples (see [4] and references therein). However, many processes require a sub-picosecond time resolution. Although such a time resolution can be obtained with slicing [5,6] or plasma sources [7,8], they only yield a low photon flux, making such setups unattractive for many experiments. This situation is expected to change with the availability of free-electron lasers (FEL) in the hard x-ray regime for pump-probe experiments. At the linear coherent light source (LCLS), the first experiments have already been performed [9]- [11]; two other sources, the SCSS [12] and the European XFEL [13], are in the construction stage, and plans for a Swiss FEL have recently been presented [14]. These facilities will provide x-ray pulses of 100 fs duration or less and with intensities exceeding those of current synchrotron facilities by orders of magnitude. As the experimental facilities start to become available, it is the task of the theoretician to provide support in the design of experiments and in the interpretation of results. This aid can take various forms.
On a fundamental level, we can investigate what information is encoded in time-resolved diffraction patterns. While this may look at first like a purely academic exercise, the need for such studies becomes apparent if we recall that standard textbook derivations are only applicable for time-independent systems; for example, they employ time-independent scattering expressions, or assume negligible electric currents (i.e. stationary systems). Studies of the theoretical foundations of TRXD have been conducted by various authors [15]- [23] and mainly differ in the type of system studied, the description of the interaction and the thoroughness of the derivation.
Most of the results for the difference cross section can be cast in the form with σ T being the classical Thomson scattering cross section from a free electron. q is the scattering vector (difference of the wave vectors of incoming and outgoing photons) and R denotes the manifold of nuclear coordinates. As a convention throughout this paper, we notate an n-dimensional vector x as x (n) and drop the superscript for n = 3. N specifies the number of atoms in the molecule. The time-averaged difference density of nuclear geometries ¯ as a function of the pump-probe delay time τ is the convolution of the difference density with the intensity profile I of the x-ray pulse, For brevity, we use the term 'difference density' to denote both of them. Physically, , t)d 3N R is the difference probability of finding a molecule whose molecular geometry is in an infinitesimal volume d 3N R around R (3N ) 0 . The squared molecular form factor |F mol | 2 specifies the signal for a fixed nuclear geometry. Within the independent atom model, it can be written in terms of the atomic form factors f a and bond vectors R ab as Equation (1) has a close resemblance to ordinary x-ray diffraction; the result from the timeindependent theory can essentially be obtained by replacing the time-dependent difference density by the time-independent density of the stationary system. However, there are distinct differences in the interpretation of the results.
In stationary molecules, the electronic potential energy surface usually exhibits one (or more) deep wells that correspond to stable geometries, and the nuclear density is strongly localized at these minima. It is then permissible to replace the density distribution by a deltafunction (or sum thereof); this removes the integration in (1), and the diffraction pattern is 4 basically the molecular form factor (or a sum of such form factors) for the well-defined molecular structure(s). This approximation is not always possible, however. In bulk liquids, for example, the distance between two atoms is not always well defined, and it is necessary to operate with distribution functions (see e.g. [3]).
In general, we expect the approximation of well-defined molecular geometries to fail in TRXD experiments, even at T = 0 K. One reason is apparent from (2): the finite duration of the x-ray pulse smears out all dynamics. On a more fundamental level, however, the wave function of the nuclear coordinates that has been excited by the pump laser is no longer in the minimum of a deep potential well. As a result, we observe quantum-mechanical wavepacket dispersion (see for example [24]- [26]); the difference density becomes delocalized, and an x-ray diffraction measurement yields a whole distribution of molecular structures.
Wavepacket dynamics also play a dominant role in molecular alignment. Using an intense non-resonant laser pulse, molecules can be strongly aligned with respect to the laser polarization axis [27,28], so that angle-resolved information can be directly extracted by a subsequent probe pulse [11,29]- [31].
Alignment also occurs as a side effect when using resonant laser pulses, since the pump laser preferentially excites molecules that are aligned with respect to the polarization axis. Typical difference densities on ultrashort timescales then consist of contributions from the excited wave function and from the 'hole' due to the depletion of the ground state, both of which are weakly aligned with respect to the laser polarization axis.
In section 2, we review and discuss the processing of anisotropic diffraction patterns resulting from such weakly aligned difference densities. Sections 3 and 4 present the details and the results of our calculation on NaI excited from the ground state by an ultraviolet (UV) laser. These calculations serve to illustrate the procedures and caveats detailed before. Our findings are summarized in section 5, and an outlook for further research is presented.

Anisotropic distributions
Anisotropic diffraction patterns have been analysed in several previous studies. For diatomic molecules, the inversion of the diffraction pattern has been studied by Ben-Nun et al [16]. In another paper, Cao and Wilson [17] studied the inversion from a cylindrically symmetric sample assuming constant atomic form factors. Several formulae in these papers are related to those we present here. However, for polyatomic molecules, it is rather difficult to generalize the procedure in [16], or the approximations are too crude, and it is difficult to assign a simple physical picture to the internuclear distribution function in [17].
Alternatively, we note that the formalism of x-ray and electron diffraction is identical within the independent atom model. For a systematic treatment of the class of symmetric top molecules, we can then draw on experience from the field of ultrafast electron diffraction. Kohl and Shipsey [32] calculated the diffraction patterns from molecules selected to be in specific rotational quantum states. Similar to our own research interest, the Zewail group studied the diffraction patterns arising from aligned ensembles [33,34], which resulted in a systematic description of the anisotropies [35]. In the following, we summarize their derivations and discuss the obtained results in some depth.
Throughout this paper, we use the integral notation (1). This is in contrast to [32]- [35], where the averaging procedure is mostly hidden in brackets. Although this introduces some overhead in the notation, we think that the manipulations in this section are more transparent when formally applied to a density distribution. The z-axis in this frame is given by the transition dipole vector D, which is assumed to be parallel to the symmetry axis of the symmetric top.
We restrict our attention to symmetric tops (which include linear molecules as a special case) that are excited by a linearly polarized laser pulse from a thermal ensemble. This already encompasses experiments on a range of symmetric molecules; see for example [2,4,32,35]. Further, we only consider the case where the transition dipole moment (i.e. the molecule-fixed vector that is aligned) is parallel to the symmetry axis. The difference density can then be expanded in a complete basis as ¯ (R (3N ) , τ ) = ¯ (θ, r (3N −6) , τ ) = n P n (cos θ) ¯ n (r (3N −6) , τ ), where P n are the Legendre polynomials, and r denotes the internal coordinates of the molecule (i.e. excluding the centre of mass (CM) and rotational coordinates) in a molecule-fixed frame. For linear molecules, there are five CM and rotational coordinates; to keep the discussion simple, we only consider the general nonlinear case in the following. The angle θ is defined in figure 1. Note that the difference density depends only on one of the three Euler angles. Using the orthogonality of the Legendre polynomials, the inverse relation reads ¯ n (r (3N −6) , τ ) = 2n + 1 2 π 0 dθ sin θ ¯ (θ, r (3N −6) , τ ) P n (cos θ).
Note that here and in the following, we always drop the integration over the two remaining Euler angles, and the CM coordinates. As the difference density is independent of them, these integrations just give rise to additional normalization factors; however, we always have to carry the sinθ -factor from the Jacobi determinant around. As detailed in [32,35], a difference density of the form (4) leads to a difference diffraction , τ ) P n (cos ϑ ab ) j n (qr ab ) (7) 6 in terms of spherical Bessel functions j n , and with r ab as the distance between atoms a and b, and ϑ ab the angle between the ab-vector and the transition dipole moment as defined in figure 1. Note that terms with odd n are antisymmetric under space inversion and give no contribution when calculating the diffraction pattern (1), so we can drop them. At this point, the link between the integration and the variables r ab and ϑ ab is somewhat opaque; we will detail the integration scheme later on.
In the simplest anisotropic case, we only consider the terms n = 0, 2, and the density is a second-order polynomial of cos θ . This is what we get whenever we excite molecules from thermal equilibrium by one-photon absorption [36]. The diffraction pattern (6) is then a simple sum of S 0 (the 'isotropic' curve) and S 2 (the 'anisotropic' curve). If we consider a fixed delay time τ and a fixed length of the scattering vector q, we can assign each detector pixel on the resulting circle a value of P 2 (cos θ q ) and a value of dσ/d . The decomposition of dσ/d into S 0 and S 2 turns into a simple linear fitting problem, whose slope gives the anisotropic and whose y-intercept gives the isotropic curve [16].
Alternatively, we can obtain the isotropic signal by 'magic angle detection' [35]. If the angle α between the incoming x-ray beam and the laser polarization axis fulfils cos 2 α = 1/3, the isotropic signal can be obtained by directly integrating over the azimuthal angle on the detector plate.
Until now, we have merely rewritten the general equation (1). To grasp the full potential of (6) and (7), we now turn our attention to the density components ¯ n and their relation to some 'fundamental' density distributions.

Density distributions
Whenever a molecule interacts with a photon, its rotational quantum number changes. During the laser-matter interaction, a molecule that starts out with a well-defined rotational quantum number becomes an ensemble (coherent superposition) of states with different rotational quantum numbers. These quantum numbers influence the subsequent time evolution through ro-vibrational coupling (e.g. centrifugal distortion).
However, in many practical applications, this coupling can be neglected. That is, the rotational motion of the molecule and its internal evolution are uncoupled, and we can approximate We have to extend this ansatz in the case of multiple species, where we understand species as well-defined, mutually distinguishable configurations that contribute to the signal. For example, there might be a 'hole' term (depletion of the stable ground state) and an excitedstate wavepacket, both of which exhibit different rotational motion due to different moments of inertia (see also the example of NaI in the subsequent sections). We assume in analogy to (8) that for each well-defined species S, the rotational motion and internal evolution are uncoupled, so that we can write 7 Comparison with (4) shows that each density component ¯ n is built up from fundamental 'vibrational' difference densities ¯ S vib that describe the time-dependent structure of the single species, and (real-valued) coefficients c S n that arise from an expansion of the angular distribution in Legendre polynomials, and thus describe the angular distribution of the species. If we do not know the exact angular distribution, these coefficients appear as additional fitting parameters.
With this in mind, we can construct a general analysis procedure for the class of anisotropic diffraction patterns presented here. We first use (6) to decompose the two-dimensional diffraction pattern dσ/d into a set of one-dimensional diffraction curves S n . In a second step, we guess the vibrational difference densities and coefficients, and build up the single density components ¯ n via (9). This guess is compared to the diffraction curves (7) and subsequently improved through some iterative schemes. For the special case of diatomics, the scattering curves can be directly inverted, see (19) below.
We point out that (6) and (7) contain the special case of isotropic densities, where only the n = 0 term contributes. The diffraction pattern (6) is isotropic, and if we write out the Legendre polynomial and spherical Bessel function in (7), the isotropic scattering curve is given by which is the well-known Debye formula for isotropic ensembles [37] (where again the integration is often suppressed). From the case of purely isotropic signals, we can extract a number of relations that are also useful for anisotropic patterns. Using distributions S rot with norm 1 and the orthogonality of the Legendre polynomials, we obtain This is just the trivial result that for an isotropic distribution, we do not have to care about molecular rotation. Furthermore, we can calculate the difference probability density of finding a specific molecular structure, that is, a certain value of r (3N −6) , by integrating ¯ over the remaining angular degree of freedom where we used (5), (9) (with P 0 (x) = 1) and (11). The vibrational difference densities ¯ S vib are thus exactly the difference probability distributions of the molecular structures for the corresponding species S. Equations (11) and (12) are, in fact, valid for arbitrary difference densities. The orthogonality of the Legendre polynomials ensures that all but the isotropic components (c S 0 , ¯ 0 ) vanish during the θ-integration. This leads to an important conclusion: if we decompose the signal according to (6) and focus only on the isotropic component S 0 , we can ignore the rotational coefficients and interpret the resulting density as a (difference) probability distribution; that is, we can process the curve S 0 and interpret the results as if we had an isotropic ensemble in the first place.
In the theory involving isotropic ensembles, (10) is often recast in a form that involves atom-pair-distribution functions (see e.g. [3]), which removes the cumbersome integration over r (3N −6) . We first calculate the difference probability density ab for the distance between atoms a, b being r 0 . In the spirit of (12), we take the probability density of finding some structure r (3N −6) and integrate it over all internal coordinates except the ab distance. Formally, we choose the internal coordinates to be the spherical coordinates ϑ ab , φ ab , r ab of the ab vector, and some 3N − 9 other internal coordinates x ab , and obtain thereby defining a difference pair distribution g ab of the ab distance. The same integration scheme can be used to evaluate (10), and we can express the signal through the pair-distribution functions as Note that pair-distribution functions can be defined in different ways. To reduce the number of terms, the sum sometimes only runs over all atom types, and the distribution function is defined per pair of atom types. In liquid scattering experiments, it is more convenient to convert the pair distribution into a dimensionless quantity called the radial distribution function [3]. The concept of using the spherical coordinates of the ab-vector for the integration and hiding most integrations in the definition of some distribution function can be extended to the anisotropic components. From (7) and (9), we obtain where the species-dependent angle-resolved difference pair distributions g S ab are Note that we now have one pair distribution per species in (15), while we had a single pair distribution in (13). This is caused by the different species having different angular distributions, which we account for by the rotational coefficients. Only if all species have the same angular distribution, or if all rotational coefficients are identical, can we sum up the species-dependent pair distributions to a global pair distribution. Finally, we want to point out that the method presented here might only be feasible for weak alignment, as induced by few-photon absorption. The stronger the alignment, the higher the order of the Legendre polynomials that appear in the decomposition (4) and (6). In particular in the case of multiple species, this leads to a large number of rotational coefficients that, in general, have to be guessed.

Direct inversion
In the following, we study the behaviour of the single density components ¯ n of NaI. NaI has the internuclear distance r as the single internal coordinate, which greatly reduces the complexity. Also, the relevant transition dipole moment is parallel to the molecular axis, so that P n (cos ϑ NaI ) = 1. Choosing real atomic form factors for convenience, (7) simplifies to For a diatomic, it is possible to invert the scattering curves directly to obtain ¯ n (r, τ ). The spherical Bessel functions have an orthogonality relation [38] ∞ 0 dq q 2 j n (uq) j n (vq) = π so that, up to constant factors, (a related formula is (30) of [17]). Since real-world data are only collected for q q max < ∞, the inverted density using [0, q max ] for the integral boundaries exhibits high-frequency oscillations, which are suppressed here by applying the exponential factor with k > 0. Again, applying this to the isotropic contribution (n = 0) gives the well-known inversion formula for isotropic samples [33]. For isotropic distributions, the inversion formula can also be used for studying polyatomic molecules [33]- [35], [39]. Let us assume that all atomic form factors have the same functional form, f a (q) = Z a f (q), with Z a being the number of electrons of atom a, insert (14) into (19) and use (18). We obtain up to constants the expression a,b Z a Z b g ab (r, τ ). This function has a peak or dip whenever r roughly equals some relevant atom-atom distance that has been established (peak) or removed (dip) during the excitation. While the underlying assumptions are too crude for a quantitative analysis, this inversion is still useful for qualitative estimates of the dynamics.
Direct inversion schemes have also been applied previously to anisotropic diffraction patterns. In [34], a mixture of S 0 and S 2 , which is obtained when the laser polarization axis is parallel to the incoming beam, was inverted using (19) with n = 0. In [35], the same procedure was applied to the diffraction pattern along specific lines on the detector (roughly corresponding to different values of θ q ) [35]. As is shown in the references, these procedures can lead to artifacts in the resulting distribution functions, such as predicting the breakage of non-existent bonds.
Therefore, we should always decompose the diffraction pattern into the single curves S n , and then invert the isotropic signal using j 0 , and the anisotropic signal using j 2 , where care has to be taken when interpreting the latter. Comparison with (15) shows that the inversion procedure extracts not only a sum of the angle-resolved difference pair distributions, but also includes contributions from the ϑ-integration and the rotational coefficients, both of which can change the amplitude and sign of the single peaks. Potential energy surfaces and sketch of the relevant nuclear motion of NaI. By absorbing a UV photon, NaI is vertically excited from the X 1 + ground state to the A 1 + first excited state, which exhibits a shallow minimum. The nuclei then vibrate on the excited-state surface with bond distances ranging from approximately 3 to 12 Å. Through the avoided crossing at 7 Å, parts of the wavepacket can make a transition to the electronic ground state and dissociate.

The NaI molecule
To illustrate the general decomposition outlined in the last section, we have calculated the dynamics of NaI excited by a UV photon. The relevant potential energy surfaces and nuclear dynamics of NaI after photo-excitation are shown in figure 2. We chose NaI because it is an extensively studied molecule [24]- [26,40]- [42] and because the excited wavepacket oscillates with large amplitudes, so we can easily distinguish the two species, which are the excited state and the depleted ground state (the 'hole'). We used the WavePacket program package [43] for the calculation, which comes with potential energy curves and dipole moments for NaI based on [41]. The wave function was represented using a pseudospectral basis with a Gauss-Legendre discrete variable representation (DVR) in θ , and an equally spaced grid for the radial coordinate r [44,45]. The initial wave function was chosen to be symmetric under a rotation around the laser polarization axis, and this symmetry is preserved by the Hamiltonian, so we can drop the azimuth angle.
Here, and in the following, we use a vector notation (r, θ, t) = (1/r )(ψ g (r, θ, t), ψe(r, θ, t)) T together with a contact transformation such that ψ g /r specifies the wave function in the electronic ground state and ψ e /r that in the electronic excited state.
In the electric dipole approximation for the molecule-field interaction, the Hamiltonian can be written aŝ where the field-free HamiltonianĤ 0 conserves the rotational quantum number, and the transition dipole momentμ in the second term is responsible for the initial switching of the electronic state with an amplitude depending on the internuclear distance r and the electric field strength ε.
The molecule was initially placed in the electronic, vibrational and rotational ground state by a propagation in imaginary time [46]. It was subsequently excited by a laser with a wavelength of 328 nm (this value was taken from previous calculations [25,26]), a sin 2 shape with 10 fs FWHM and an intensity of 10 11 W cm −2 . The propagation was carried out using the split operator method detailed in [44]. For simplicity, we carried out the propagation nonperturbatively; however, the laser intensity was chosen such that first-order effects dominate; that is, only terms n with n = 0, 2 are significant. The details of the calculation and the scripts used for data processing can be found in the supplementary material (available from stacks.iop.org/NJP/12/113022/mmedia).
We note in passing that the valence electron in NaI undergoes substantial reorganization whenever the internuclear distance passes the avoided crossing [26]. This, however, does not invalidate the independent atom model used to represent the molecular electron density, since the vast majority of the electrons can be considered as core electrons associated with a specific atom. Figure 3 shows the isotropic and anisotropic difference densities of NaI calculated via equation (20). We can immediately identify several features that we discussed earlier. We have two welldefined species: the oscillating excited-state wavepacket and the ground-state 'hole', which shows minor vibrations. Qualitatively, both difference densities are similar and encode the same 'fundamental' time evolution. We can also see that the hole contribution changes sign in ¯ 2 at τ = 1.5 ps. From the discussion in section 2, especially equation (9), we recall that the density components for each species contain some universal vibrational density and a rotational coefficient that describes the orientation of the species. In our calculation, the sign change of the anisotropic density component comes from a sign change of c hole 2 , that is, a change of the rotational distribution of the hole. This does not happen to the excited-state wavepacket, thus emphasizing the need to assign different rotational coefficients to different species.

Dynamics of laser-excited NaI
We can deduce from (21) that the excitation is proportional to cos θ on the wave function and cos 2 θ on a density level. In a simple 'classical' picture, we can interpret this such that the laser excites the molecules ∼cos 2 θ and leaves a 'hole' in the ground state ∼ −cos 2 θ . In the lowest-order perturbation, only the states l l max = 2 contribute, which gives a smallest rotational period of T = 2πh/Bl max (l max + 1) ≈ 45 ps. The sign change after 1.5 ps occurs because the hole is initially not aligned with the laser field, but starts the rotation with a phase that is determined by the details of the excitation process.

Interpretation of the dynamics.
To get a deeper physical insight into this process, we study the NaI dynamics with an analytic model. In a perturbation theory interpretation, the molecule is always subject to the field-free Hamiltonian and interacts with the electric field at a random time once and at two random times for second-order effects (on the wave function level). With this in mind, we can explicitly write down expressions for the difference density. For simplicity, we neglect the non-adiabatic transition, so thatĤ 0 is diagonal in the electronic basis.
If the laser is turned off, the molecule remains in the electronic, vibrational and rotational ground state, that is, Here, ψ GS denotes the vibrational ground-state function, whose sole time dependence is a complex phase. If the laser is turned on and we assume no significant rotation during the laser pulse, then the wave function at the end of the laser pulse t e can be written up to second order as on (R, t e ) = 1 r a 0 ψ GS (r, t e ) + a 2 ψ hole (r, t e ) cos 2 θ a 1 ψ exc (r, t e ) cos θ .
A major portion of the wave function is not excited at all and stays in the original ground state ψ GS . By absorbing a photon, the molecule can switch to the excited electronic state and move on the corresponding surface, which is described by the wave function ψ exc . Due to the form of the interaction Hamiltonian, it acquires an angular factor of cos θ . After the molecule is excited, it may become de-excited by emitting a photon, thereby falling back to the electronic ground state. The resulting wave function ψ hole is the result of a two-photon process and consequently has an orientation dependence of cos 2 θ. For easier comparison of the magnitude of various terms, we have split each contribution into a complex-valued, normalized wave function ψ i and a real-valued coefficient a i . If the onephoton excitation probability is given by p 1, then we would expect a 2 0 ≈ 1 − p, a 2 1 = p, and a 2 2 ≈ p 2 , since the two-photon process arises from de-excitation from the excited subset of the wave function. However, the interference term between ψ hole and ψ GS scales with a 0 a 2 ≈ p and thus has the same order of magnitude as the wave function on the excited-state surface. In fact, this term describes the detailed dynamics of the created hole in the ground state and has to be included for a correct treatment.
As the field-free Hamiltonian does not couple states with different rotational quantum numbers, we rewrite (23) using Legendre polynomials and obtain an expression that is valid for times t > t e on (R, t) ≈ 1 r where ψ hole,0 (r, t e ) = ψ hole,2 (r, t e ) = ψ hole (r, t e ) from (23), and with separate time dependences for each radial wave function. Inserting (22) and (24) into (20) and collecting only terms up to order p, we obtain for the isotropic difference density (25) and for the anisotropic difference density Note that the weights of the ground-and excited-state contributions to the difference densities are of similar magnitude. The last term in (25) and (26) describes the wavepacket corresponding to the excited-state species. It performs the large-amplitude oscillations in figure 3, and one readily observes that the contribution is identical for both difference density components. The remaining terms describe the hole; the first term in (25) describes the static population depletion of the ground state and gives a negative contribution to the difference density around r = 3 Å, while the remaining terms are the coherence or interference terms between the unexcited ground-state wave function and the second-order perturbed wave function. If we assume for the moment that the latter is a real eigenstate ψ x of the unperturbed Hamiltonian, we can write The product of the two stationary wave functions gives some characteristic density distribution, which then oscillates in time as described by the cosine. The hole wave functions can be expressed as the result of an absorption of a photon from the ground state, propagation on the excited-state surface for a certain time and subsequent emission of the photon. Consequently, the phase ϕ x is an arbitrary number that depends on the details of the system and the excitation process.
In general, the hole wave functions ψ hole,i can be written as sums of eigenstates, and the full time dynamics are a superposition of single interference terms of the form (27). In our calculation, we found the most important contributions to ψ hole,i to be the vibrational ground state and first excited state ψ 0 , ψ 1 (we neglect ro-vibrational coupling that would make these states dependent on the rotational quantum number).
Let us first consider the isotropic coherence terms in (25). Since we started in the vibrational ground state, ψ 0 ψ GS = |ψ GS | 2 and E 0 = E GS . Consequently, the coherence term ∼ ψ 0 ψ GS is time independent, has the same form as the depletion of the ground state and effectively only serves to modify the prefactor of the second term in (25). The product ψ 1 ψ GS , in contrast, is a function with one node and localized around r = 3 Å. The energy difference E 1 − E GS is about 32 meV for NaI, corresponding to an oscillation of 130 fs duration. This interference term causes a small oscillation of the difference density around r = 3 Å in figure 3. However, the overall difference density remains negative at all times due to the ground-state depletion.
Essentially the same oscillation is seen in the anisotropic difference density, which we can also trace back to the interference between ψ GS and the first vibrationally excited state in (26). We also find ψ 0 ψ GS = |ψ GS | 2 , which gives a contribution similar to the depletion of the ground state in the isotropic difference density. However, since ψ hole,2 describes the time evolution of the rotationally excited wave function, the time evolution contains an additional rotational energy term, so that E l=2 0 − E GS is now approximately 0.1 meV, corresponding to a slow oscillation of 45 ps duration.
This means that 2 has a contribution of the form c 2 (t)|ψ GS | 2 , where the prefactor c 2 comes from the cosine in (27) and oscillates slowly. From (9), it is then natural to assign −|ψ GS | 2 to the vibrational difference distribution of the hole species and c 2 (t) to the rotational coefficient and describe the slow oscillation as the rotation of the hole.
We can then wonder what the rotational distribution of the hole at the end of the laser pulse is. Obviously, from (27), this is determined by the value of ϕ l=2 0 , which in turn depends on the molecule and the laser parameters. That is, while the hole species rotates with an overall period of 45 ps, the initial orientation is determined by the details of the excitation process. For our calculation, we found a phase shift of slightly less than −π/2, so the hole starts out weakly aligned with the polarization axis, and the sign of the interference term changes for the first time after only 1.5 ps.

X-ray diffraction from laser-excited NaI
For all results, we assumed an x-ray beam with a sin 2 -pulse shape and a half-width of 100 fs, which should roughly correspond to the situation at the LCLS including timing jitter [9]. Atomic form factors were taken from [47].
A typical diffraction pattern is shown in figure 4. As detailed in section 2, this rather complex pattern is built up from two diffraction curves (here calculated from (17)), which encode the internal dynamics of the two species and their rotational coefficients. From section 2, equation (6), we know that for one-photon excitation, the difference signal is up to constant factors given by dσ/d = S 0 (q) − P 2 (cos θ q )S 2 (q). Consequently, if we keep q constant, and look at different pixels corresponding to different θ q (crosses for q = 1.5 Å −1 ), we can plot the cross section as a function of P 2 (cos θ q ) and obtain S 0 (q) and S 2 (q) from the offset and negative slope, respectively.
Using this procedure for all q gives the two scattering curves. In the case of noisy experimental data, a least squares fit should also yield some measure for the error of the diffraction curves. The isotropic curve can be compared with various approximations, such as an azimuthal integration over the detector plane, which gives surprisingly good results (see [35] for a discussion of this and other approximation schemes). However, we want to point out that the correct decomposition itself is simple enough, so one can easily avoid the loss of accuracy inherent in the approximate scheme.
Formally, we can assign each detector pixel a value of q and θ q , where the latter depends on the angle between the wave vector of the incoming x-ray beam and the laser polarization axis. Consequently, we get different diffraction patterns if we vary this angle. However, this  Left: calculated two-dimensional difference diffraction pattern in arbitrary units for NaI at a pump-probe delay time of 500 fs. q z and q ⊥ give the length of the component of q parallel and perpendicular to the laser polarization axis, which was set perpendicular to the wave vector of the x-ray beam. Middle: plot of dσ/d as a function of P 2 (cos θ q ) (q = 1.5 Å −1 , points correspond to those of the left image). Right: isotropic and anisotropic curves S n (q, τ ) encoded in the two-dimensional pattern. Also shown is the result from an azimuthal integration over the detector surface.
only affects how well we can decompose the pattern into the single curves S n ; the curves themselves are independent of the experimental setup. For the figures presented here, we assumed a perpendicular setup where the anisotropy is most pronounced.
In figure 5, we compared the real difference densities to an inversion from the diffraction data using (19). We used a pump-probe delay time of 1650 fs because the vibrational interference pattern is most pronounced in the anisotropic difference density at this time. We assumed LCLS beam parameters of 8 keV photon energy; with a maximum scattering angle of 60 • , only scattering vectors up to q ≈ 4 Å −1 are collected. We find that this range of scattering vectors is too small to resolve-in particular, the vibrational interference termswhereas perfect agreement with the original density distribution is obtained with a larger range of scattering vectors, as indicated in the figure.
Until now, we have only studied NaI in the rotational ground state. In practice, one usually starts with a Boltzmann distribution of initial rotational states. However, our basic findings still hold in this case. The orientation of the hole is still determined by an interference term between the unexcited molecular wave function and a second-order wave function, and the phase between them (i.e. orientation of the hole) depends on the detailed dynamics of the latter. We also find two additional effects that do not occur in the ground-state calculation.
Firstly, we point out that if the molecule absorbs a photon in the rotational ground state, the resulting excited-state wave function ψ exc consists of a single rotational quantum state with l = 1; the density |ψ exc | 2 is always preferentially oriented towards the laser polarization axis. This is a special case; for a nonzero initial rotational quantum number l 0 , ψ exc is composed of contributions with l 0 + 1 and l 0 − 1. As a consequence, its contribution to the difference density also includes rotational interference terms, and the anisotropic density can show amplitude or sign changes that we can interpret as rotation of the excited molecule.
A second effect that we expect to occur is rotational dephasing [28,34,39]. The rotational interference terms for ensemble members with different initial angular momenta oscillate with different timescales. When performing the ensemble average some time after the excitation, the different anisotropic contributions typically cancel. In effect, the diffraction pattern becomes almost isotropic on the timescale of a few picoseconds with alignment revivals at certain times.

Conclusions
Whenever we excite molecules with a linearly polarized laser, we prepare anisotropic ensembles. When we probe the subsequent sub-picosecond dynamics with TRXD, the resulting difference diffraction patterns also show a pronounced anisotropy. It has been demonstrated that such signals can be decomposed into an isotropic and one or more anisotropic diffraction curves. For alignment through one-photon absorption, the decomposition turns into a linear fitting procedure, and gives an isotropic and an anisotropic curve. These curves can be processed independently of each other, and they encode the same basic information about the nonequilibrium molecular structure.
The isotropic curve is described by the usual Debye formula [37] and gives the same signal as we would obtain from an isotropic ensemble. The anisotropic curve contains additional coefficients that describe the rotational distribution of the molecule, which makes data extraction more difficult. On the other hand, it is related to an angle-resolved pair distribution; that is, it also includes explicit information about the orientation of the single bonds with respect to the transition dipole moment, which might be used to supplement the information extracted from the isotropic curve.
We have calculated the diffraction patterns and difference densities for the diatomic NaI molecule excited from the rotational, vibrational and electronic ground states. We have demonstrated both numerically and analytically that the density distribution can be extracted through an exact inversion procedure.
By interpreting the dynamics with a simple analytic model, we could explicitly demonstrate that the orientation of the hole species is strongly dependent on the details of the excitation process. This should serve as a warning that the orientation of the single species cannot be guessed from simple models that omit quantum-mechanical interference terms (for example, equation (32) of [35] assumes that when molecules are preferentially excited along the laser polarization axis, the hole is also oriented in this direction).
Apart from applying the formalism presented here to larger systems, one possible route for expansion is molecular imaging. It has been proposed to align molecules in their stable ground state in the gas phase by a non-resonant laser to extract their structure [22,30,31,48]. Due to imperfect alignment, however, the structural information that we get in the laboratory frame is always blurred to some extent. To try to overcome this limitation, it is intuitively appealing to use the same basic procedure of separating the density into a rotational and an 'internal' component. However, it remains to be seen how useful this concept is in practice and whether it can be used to create a convergent and stable algorithm.