Probing electron-phonon interactions away from the Fermi level with resonant inelastic x-ray scattering

Interactions between electrons and lattice vibrations are responsible for a wide range of material properties and applications. Recently, there has been considerable interest in the development of resonant inelastic x-ray scattering (RIXS) as a tool for measuring electron-phonon (e-ph) interactions. Here, we demonstrate the ability of RIXS to probe the interaction between phonons and specific electronic states both near to, and away from, the Fermi level. We performed carbon $K$-edge RIXS measurements on graphite, tuning the incident x-ray energy to separately probe the interactions of the $\pi^*$ and $\sigma^*$ electronic states. Our high-resolution data reveals detailed structure in the multi-phonon RIXS features that directly encodes the momentum dependence of the e-ph interaction strength. We develop a Green's-function method to model this structure, which naturally accounts for the phonon and interaction-strength dispersions, as well as the mixing of phonon momenta in the intermediate state. This model shows that the differences between the spectra can be fully explained by contrasting trends of the e-ph interaction through the Brillouin zone, being concentrated at the $\Gamma$ and $K$ points for the $\pi^*$ states, while being significant at all momenta for the $\sigma^*$ states. Our results advance the interpretation of phonon excitations in RIXS, and extend its applicability as a probe of e-ph interactions to a new range of out-of-equilibrium situations.

Interactions between electrons and lattice vibrations are responsible for a wide range of material properties and applications. Recently, there has been considerable interest in the development of resonant inelastic x-ray scattering (RIXS) as a tool for measuring electron-phonon (e-ph) interactions. Here, we demonstrate the ability of RIXS to probe the interaction between phonons and specific electronic states both near to, and away from, the Fermi level. We performed carbon Kedge RIXS measurements on graphite, tuning the incident x-ray energy to separately probe the interactions of the π * and σ * electronic states. Our high-resolution data reveals detailed structure in the multi-phonon RIXS features that directly encodes the momentum dependence of the e-ph interaction strength. We develop a Green's-function method to model this structure, which naturally accounts for the phonon and interaction-strength dispersions, as well as the mixing of phonon momenta in the intermediate state. This model shows that the differences between the spectra can be fully explained by contrasting trends of the e-ph interaction through the Brillouin zone, being concentrated at the Γ and K points for the π * states, while being significant at all momenta for the σ * states. Our results advance the interpretation of phonon excitations in RIXS, and extend its applicability as a probe of e-ph interactions to a new range of out-of-equilibrium situations.

I. INTRODUCTION
Knowledge of the interactions between electrons and phonons is central to understanding a diverse array of condensed matter phenomena. Basic material properties such as low-field electrical transport [1], as well as emergent collective phases including charge density waves [2] and conventional superconductivity [3], are all crucially dependent on the coupling of phonons to electrons near the Fermi surface. Extensive experimental efforts have therefore been devoted to quantifying this coupling through various techniques. Neutron scattering [4], nonresonant inelastic x-ray scattering (IXS) [5], and Raman spectroscopy [6] all probe e-ph coupling via the lattice vancing theoretical understanding of the phonon contribution to the RIXS cross-section [28][29][30][31][32]. Much of the excitement around these developments lies in the momentum-resolution that RIXS offers. The other great strength of RIXS -the ability to probe a particular electronic orbital by tuning the incident energy to a specific absorption edge -has thus far been under exploited.
Here, we demonstrate the power of RIXS to probe the interaction of phonons with two distinct electronic states in graphite: the π * state near the Fermi level, and the σ * state well above it. At the resonances of these two states, we find qualitatively different multiphonon excitations in our RIXS spectra, showing a stark difference in the phonon momenta to which they couple. This behaviour cannot be captured by the currently available theoretical models, which assume dispersionless phonons [19,28,32] or focus on one-phonon processes [29]. We therefore extend a Green's-function approach, previously applied to small molecules [30], to treat the full momentum-dependence of phonons in a crystalline lattice. This model accurately reproduces our experimental spectra at both resonances, revealing contrasting trends of the e-ph interaction strength through the Brillouin zone.

II. RIXS MEASUREMENTS ON GRAPHITE
Graphite is a favorable material for this study, as it has a simple electronic and phononic structure, as well as two well-defined excitations to the low-energy π * and high-energy σ * states at the carbon K edge. Coupling to the electronic states near the Fermi level is known to be dominated by the highest optical phonon band above 0.15 eV [5,6,[34][35][36], which is well resolved by our 47 meV-resolution measurements. There are also no other low-energy excitations (such as magnons) that could obscure the phonon excitations in our RIXS spectra. This allows us to validate our methodology at the π * resonance, before extending our treatment to the σ * states. Given their similarities to graphite [6,34,35], our results are also applicable to the other technologicallyimportant allotropes of carbon.
Our carbon K-edge RIXS measurements on graphite are summarized in Fig. 1(a), which shows approximately zone-center RIXS maps around the π * and σ * resonances (see Appendix A for experimental details). Above the elastic line at 0 eV energy loss, both maps show a series of features of decreasing intensity that are reminiscent of the harmonic progression of phonon excitations seen in previous studies [15][16][17][18][19][20][21][22][23][24][25][26][27]. These features are more intense in the σ * map, suggesting a stronger overall coupling.
We assign contributions to the phonon features through reference to the known phonon dispersion of graphite, shown in Fig. 1(b) [33]. The lowest-energy feature occurs at ∼0.19 eV in both maps, which corresponds to the degenerate transverse (TO) and longitudinal optical (LO) modes at the Γ point [whose displacement pat- Contrasting RIXS response at the π * and σ * resonances. (a) RIXS maps around both resonances, with XAS above for reference (the dashed vertical lines in the XAS mark the peaks of the resonances at 285.6 eV and 298.1 eV respectively). Both maps show a series of phonon features above the elastic line, with the contrasting resonance behavior most apparent for the two-phonon feature between 0.3-0.4 eV. The intense feature above ∼1 eV in the π * map arises from electronic transitions. (b) Phonon dispersion of graphite from Ref. [33], with the transverse (T), longitudinal (L) and outof-plane (Z) acoustic (A) and optical (O) modes marked. The TO mode with significant coupling is marked in red. (c) The 2D-projected Brillouin zone of graphite, with the momenta accessible in our RIXS measurement indicated by the turquoise circle, and high symmetry positions Γ, K, and M marked.
tern is shown in the final state in Fig. 2(a)]. In a purely local picture, successive peaks would then correspond to higher harmonics of this first peak. The next feature around 0.3-0.4 eV energy loss, however, shows an internal structure that differs markedly between the two maps. At the π * resonance (285.6 eV) the feature is split into two (as can also be seen in the linecuts in Fig. 3), while at the σ * resonance (291.8 eV) a single broader asymmetric peak is visible.
To conserve momentum, a single-phonon excitation in a zone-centre spectrum can include only a zone-centre mode in the final state. The final states of a twophonon feature, however, can conceivably consist of pairs of opposite-momenta phonons at any point in the Brillouin zone. The split feature at the π * resonance has peak energies of 0.32 eV and 0.39 eV, which correspond to double the energies of the TO mode at K and Γ respectively (see the Supplemental Material). Meanwhile, the broad feature at the σ * resonance spans double the full bandwidth of the TO mode. Phonons of different momenta therefore appear to be excited at each of the resonances. ωi and wavevector ki excites an electron from the 1s to the π * or σ * bands, leaving a core hole. In the intermediate state, the altered charge density perturbs the positions of surrounding ions and generates nq phonon modes. The excited electron then relaxes via the emission of a photon of energy ω f , leaving n q phonons with total energy (ωi − ω f ) in the final state (the zone-center TO mode is depicted [33]). (b) Electronic band structure of graphite from Ref. [37], with the low-energy π * bands shaded in orange and the high-energy σ * bands in blue.

III. ELECTRON-PHONON INTERACTIONS IN RIXS
We have identified the modes that contribute to the RIXS spectra, but why do phonons of different momenta appear at the π * and σ * resonances? To answer this, we need to understand how phonons are generated in a RIXS measurement. Figure 2(a) shows a schematic of the RIXS process at the carbon K edge. An incident x-ray photon excites an electron from the 1s to the π * or σ * bands, depending on the energy, where it experiences a partially-screened potential from the resulting core hole. The altered charge density of this intermediate state perturbs the positions of the surrounding ions, generating multiple phonons. Finally, the excited electron relaxes to fill the core hole through the emission of a photon, with the difference in incident and scattered photon energies equal to the energy of the phonons remaining in the final state.
The composition of the intermediate electronic state can therefore have a significant impact on the phonon modes that are excited. As both resonances involve the same core-hole, this cannot explain the difference between the RIXS spectra. We instead focus on the band structure of the excited electron, shown in Fig. 2(b). For electrons excited into the low-energy π * states, the steepness of the bands confines them to a small region around the Dirac points at K. These states can only interact with phonons with momenta close to Γ, which scatter electrons within a single Dirac cone, or K, which scatter electrons between different Dirac cones. By contrast, the high-energy σ * bands are flatter and electrons excited into them can span the Brillouin zone, suggesting interactions with a wide range of phonon momenta. This provides an intuitive understanding of our earlier assignment of the RIXS features, where the split two-phonon feature at the π * resonance is dominated by phonons at Γ and K, while the single broad feature at the σ * resonance has contributions from across the TO mode dispersion.
To put this on a more rigorous footing, we need a theoretical model that can account for interactions of the full excitonic intermediate state with multiple phonons of different momenta. Within the exact-diagonalization approaches proposed thus far [19,28,29], this requires a computationally-intractable expansion of the Hilbert space. We therefore employed a Green's-function approach which implicitly accounts for vibronic effects in the intermediate state (full details can be found in Appendix B). We address e-ph interactions within the framework of a linked-cluster (cumulant) expansion using a second-order Fan-Migdal-type self-energy in the time domain. RIXS cross-sections were obtained by considering the non-diagonal elements of an exciton Green's function, as described elsewhere [30]. As well as fully accounting for the phonon dispersion, this approach allows a proper treatment of mode mixing in the intermediate state, which leads to significant renormalization of the spectral features.
In order to focus on the effect of the interactionstrength dispersion, we treat both π * and σ * intermediate states as localised excitons, which is a good approximation over the intermediate-state lifetime of graphite [38]. The intermediate-state lifetime itself was determined using the core-hole lifetime from x-ray photoemission spectroscopy [39], adjusted to mimic the presence of the extra electron (see the Supplemental Material). Given the weak inter-layer coupling of graphite, the optical phonon modes of interest have a minimal dependence on the out- 3. e-ph interactions for low-energy π * excitations. (a) Modeled momentum dependence of G through the 2D Brillouin zone (orange), compared to g determined by IXS (blue) [5], Raman spectroscopy (yellow) [6] and time-resolved ARPES (purple) [36]. (b) Normalized experimental (black points with error bars) and calculated (orange line, including experimental broadening) RIXS spectra on resonance (285.6 eV). The contribution from each phonon generation process is indicated by the vertical orange lines, grouped above by the number of phonons. (c) Incident-energy dependence of the experimental (black points with error bars) and calculated (red to orange lines) spectra, plotted alongside the XAS with arrows indicating the energies.
of-plane momentum [33] and we can safely restrict our analysis to the 2D-projected Brillouin zone. We used the phonon dispersion calculated with density functional perturbation theory, which is in good agreement with IXS measurements [33]. The only free parameter of the model is then the momentum-dependent e-ph interaction strength, G(q), which can be determined by fitting the experimental RIXS spectra.

IV. ELECTRON-PHONON INTERACTIONS NEAR THE FERMI LEVEL
We turn first to our data around the π * resonance. From our earlier assignment of the RIXS features, it appears that significant contributions come only from the TO mode at the K point, and degenerate TO/LO modes at the Γ point (see the Supplemental Material). We therefore include only the TO mode dispersion in our modelling, with a minimal set of momentum points: the Γ and K points, the experimental in-plane momentum transfer q exp ≈ (0.08, 0.00) r.l.u., and K − q exp . Inclusion of the latter two effectively captures the width over which G(q) falls off away from the high-symmetry points. The intensity of the one-phonon peak is sensitive to the fact that q exp = Γ, but for the multi-phonon features it is sufficient to approximate the total momentum transfer as zero. Figure 3(b) shows a fit of this model to the spectrum on resonance, with the corresponding momentum dependence of G(q) plotted in Fig. 3(a). We see good agreement between the experimental and calculated spectra, including the striking split two-phonon feature, despite our simple parameterization of G(q). This confirms that the splitting arises due to the dominant contributions from pairs of phonons at K and Γ, indicated by the vertical lines at 0.32 eV and 0.39 eV respectively. Intensity above 0.45 eV is attributed to overlapping three-and four-phonon features, which are also consistent with contributions coming from pairs of ±K-momenta phonons with additional Γ modes.
We further tested our model through its ability to reproduce spectra with the incident x-ray energy detuned from the maximum of the XAS. Figure 3(c) shows the energy dependence of the experimental and calculated spectra through the resonance. As well as the persistent splitting of the two-phonon feature below the resonance, our model captures the changing relative intensities of the one-and two-phonon features. The non-monotonic energy dependence of these intensities, including the offset of the maximum of the one-phonon peak intensity from zero-detuning [27], arises from a combination of strong e-ph interaction strength and long core-hole lifetime [31]. We emphasise that the only parameter being varied between the spectra in Fig. 3(c) is the incident energy, providing a strong validation of our model.
As the incident energy is raised through the upper tail of the XAS peak, the composition of the intermediate electronic state will change (picking up contributions from the M point, see Fig. 2 [40]. This in turn allows scattering by phonons of different momenta, changing the momentum-dependence of G and explaining the discrepancies between the experimental and calculated spectra above 285.6 eV. Agreement can be improved by allowing G(q) to vary slightly with incident energy (see the Supplemental Material).
We can now compare our results with those from other techniques. Alongside our best-fit G(q), Fig. 3(a) also shows the e-ph interaction strength determined by IXS [5] and Raman spectroscopy [6], as well as a novel time-resolved ARPES technique [36]. We see reasonable agreement between all the techniques. Note, however, that we use a lower-case g for the coupling probed by these other techniques, to highlight that they do not involve a core hole. As discussed above, RIXS measures a related but inequivalent coupling, G, of phonons to an excitonic intermediate state consisting of both the excited electron and core hole [30]. The concept of exciton-phonon coupling is not new, being relevant to Raman and other optical spectroscopies [41,42] (the only difference here being that the hole is often in the valence band rather than core level). Recent work in these fields has shown that, under certain approximations, the two quantities G and g can be related by [43] This expression describes an exciton in initial state n, with centre-of-mass momentum Q, and wave-function |S n = hc A n hc |hc , where h (c) indicate hole (conduction) states. The exciton is scattered by a phonon in mode ν with momentum q, to a final state m with momentum Q + q and wave-function |S m = h c A m h c |h c . It can be seen that Eq. 1 is separated into two terms, the first accounting for electron-phonon scattering with the hole-state unaffected, and the second for hole-phonon scattering with a passive electron, each weighted by the coefficients of the excitonic wavefunction. In cases where the core hole is strongly screened, RIXS reflects the single-particle coupling of the excited electron. When the core hole is weakly screened, RIXS accesses exciton-phonon coupling, which is the quantity of interest for the many out-of-equilibrium situations discussed in the Introduction. From Fig. 3(a), we can also see how the various techniques probe different components of g ν cc (Q, q) and G ν nm (Q, q). Standard ARPES resolves g by electronic momentum Q for bands c below the Fermi energy, but integrates over phonon modes ν and momenta q [7]. By contrast, IXS [5] and neutron scattering [4] resolve the phonon momentum q throughout the Brillouin zone for all modes ν, but average over the equilibrium electronic state. Raman spectroscopy [6] and the novel timeresolved ARPES [36] technique can access low-energy electronic excitations, but are restricted to particular phonon momenta (q = Γ and K for graphite). By contrast, we have shown that RIXS can resolve G by electronic state n, phonon mode ν and phonon momentum q, integrating only over the electronic momentum Q.
By modelling the contributions to the multi-phonon features, we have therefore determined the interactions between the π * electronic state and TO phonons throughout the Brillouin zone, all from a single zone-center spectrum. This method will be invaluable, as the soft x-ray edges at which phonons are enhanced (due to the long intermediate-state lifetimes) and energy resolutions are maximised, also suffer from restricted momentum transfer [see Fig. 1(c)].

V. ELECTRON-PHONON INTERACTIONS AWAY FROM THE FERMI LEVEL
Having successfully reproduced the π * spectra, we now come to the high-energy σ * states. Based on our assignment of the RIXS features, we again restrict our analysis to the TO mode (see the Supplemental Material). In this case, however, the broad profile of the two-phonon peak suggests G(q) is finite over a wide range of momenta. In order to constrain our parameter space, we therefore assume a linear dependence between the high-symmetry points Γ, K, and M . Given the more uniform variation of G(q) around Γ, we can also approximate the experimental momentum transfer as zero for the whole spectrum. Figure 4(a) shows the best-fit dispersion of G(q), which we find to be significant at all the high-symmetry points. The resulting spectrum on resonance is shown in Fig. 4(b), and again we see that a simple parameterization of reciprocal space is able to describe the experimental data well. Unlike the splitting seen at the π * resonance, here our model reproduces the subtle broadening and asymmetry of the two-phonon feature, confirming that these characteristics arise from the varying contributions of pairs of opposite-momenta phonons which disperse gradually through the Brillouin zone.
In order to determine the sensitivity of our fit to the values at the high-symmetry points, we also show spectra in Fig. 4(b) calculated with G(M ) = 0 and G(K) = 0 (G(Γ) must be finite to produce a significant one-phonon peak). We can see that the calculated spectra are most sensitive to the value of G(K), which strongly influences the intensities of all the multi-phonon peaks. While the value of G(M ) has little impact on the three-, four-and five-phonon peaks, we see that a finite value on the order of G(Γ) is required to accurately capture the intensity and shape of the two-phonon peak. We also show in the Supplemental Material how a previous RIXS measurement of G(q) at the σ * resonance of graphite [27], which assumed zero coupling between the high-symmetry points, is unable to reproduce the multi-phonon features.
Our best-fit G(q) is also able to reproduce the changing intensities of the features on detuning below 291.8 eV, as shown in Fig. 4(c). As for the π * states, we expect agreement between the model and experiment to degrade above the resonance due to changes in the intermediate electronic state. Further, just above the σ * 1 resonance, that we are focused on, lies the σ * 2 resonance [see Fig. 2(b)], which arises from zero-point vibrations [40] that are not included in our modelling.

VI. CONCLUSION
We have demonstrated that, through careful choice of the incident photon energy, RIXS is able to probe the interactions between phonons and electrons both close to, and away from, the Fermi level. Our results for graphite reveal that the momentum dependence of the interaction strength, G(q), can be highly distinct in these two regimes: for low-energy π * states it is concentrated in small regions around the Γ and K points, while for high-energy σ * states it is significant across the Brillouin zone. This establishes RIXS as a uniquely versatile probe of e-ph interactions in crystalline materials, applicable to studies of out-of-equilibrium processes alongside ordered phases in correlated systems. Combined with the recent development of ultrafast time-resolved RIXS at xray free-electron lasers [44][45][46][47], our methodology should prove particularly powerful for determining e-ph interactions in pumped systems.
We find that contributions to the multi-phonon features in a RIXS spectrum can come from phonons away from the experimental momentum transfer, allowing us to constrain G(q) throughout the Brillouin zone from a single zone-center spectrum. This is particularly apparent at the π * resonance, where the rapid variation of G with momentum produces a split two-phonon feature.
These contributions are also important, however, for the σ * states with more uniform G(q). Here, the gradually dispersing phonons cause a broadening and asymmetry of the multi-phonon features, as has been seen in other studies [15,22], alongside shifts in their peak energies that could be mistaken for anharmonic effects [22]. As we show in Appendix B, mode mixing in the intermediate state results in the multi-phonon contributions affecting the intensity and resonance behavior of all the phonon features in the spectra, including the one-phonon peak. Neglecting this effect is likely to lead to qualitative errors in the determination of G(q) [27]. These insights from graphite will be important when analysing the phonon excitations in RIXS measurements of more complex materials, such as the much-studied cuprates [14-16, 19, 21, 23, 25, 26].
Our results highlight the importance of developing a robust theoretical description of the complex RIXS process in order to analyse experimental data. The Green's function approach that we employ offers a practical way to incorporate the phonon and e-ph interaction strength dispersions without the need for summation over all possible intermediate states. There are many avenues that remain to be explored, including the effects of itinerancy [32] and changes to the potential energy surface [30] in the intermediate state, and how the core hole affects the coupling probed by RIXS compared to other techniques [30]. The availability of theoretical methods to treat these effects will only become more pressing as the energy resolution of RIXS improves, and it is used to probe phonons in an ever increasing range of topical materials. Our carbon K-edge RIXS measurements were performed at beamline I21 of the Diamond Light Source. A natural graphite single crystal was mounted such that the scattering plane is (h, 0, l), and cleaved in vacuum. We chose an energy resolution of 47 meV (determined from scattering from amorphous carbon tape mounted next to the sample) to balance throughput and resolution. All RIXS spectra were taken with linear horizontal (π) incident x-ray polarisation, at a temperature of 20 K, and with a counting time of 20 min. X-ray absorption spectra in total electron yield were obtained by measuring the sample drain current. Data for energies around the π * resonance were taken with grazing-incident x-rays (θ = 20°) while data around the σ * resonance were taken with normal-incident x-rays (θ = 90°) in order to maximise the intensity of the phonon features. All data have been corrected for incident flux and self-absorption (see the Supplemental Material).

Appendix B: Theoretical model
To model our rich experimental data we employed a phenomenological Green's-function formulation of the RIXS cross-section. In principle, final states involving both phonons and low-energy electronic excitations are allowed, the latter of which would result in a continuum of excitations. As such a continuum is not seen above the background in our data, however, we neglect these weak higher-order process and restrict the final states to those with only phonon excitations present.
We account for electron-electron interactions by constructing an effective quasi-particle exciton Green's function. Due to the different time scales, it is reasonable to separately consider electron-lattice interactions in the presence of an exciton which propagates through the material [30,42]. Expanding the interacting exciton Green's function in a time series, we find two important types of diagrams: closed-loop diagrams where a phonon gets created and destroyed by the same exciton, and those where a phonon connects two excitons and gives rise to the final state phonon population. The former affect the intensities of the RIXS features non-linearly, and can be interpreted as contributions from intermediate state phonons. Evaluating these contributions has proven to be very computationally demanding in exact diagonalization approaches. Here, we account for them using the cumulant representation of the interacting Green's function, which reduces computational cost. The second type of diagram will be treated explicitly up to the number of phonon satellites visible in the experiment.
To simplify the notation, we omit electron-photon matrix elements which affect only the absolute intensity and not the characteristic shape of the phonon features, and drop the index of the excitonic state. The energy-and momentum-dependent cross-section then can be written as where Λ n (ω, q) is the complex off-diagonal part of the exciton Green's function, representing the scattering of an exciton by n phonons of total momentum q = k i −k f , and the presence of phonons in the final state is reflected by the many-body Green's function D n (ω, q).

In the limit of infinite phonon life-time
Here, we account for phonon-phonon interactions by attributing a small lifetime broadening to the phonon propagator.
Previous work has suggested an approximate form of Λ in the limits of weak [29] and intermediate [30] e-ph interactions. For the purpose of this work, however, we use a closed-form solution for Λ in the limit of local excitons coupled linearly to phonons. Although the cross-section in Eq. (B1) is given in the energy domain, from here on we use the time-dependent form. Focusing on the characteristic behavior of the phonon satellites with multiple momentum contributions, we have (B2) The diagonal part of the exciton Green's function is denoted Λ 0 , and is the same dressed propagator that would appear in the absorption process (I XAS ∼ F {Λ 0 }). The bracketed part of Eq. (B2) represents the phonon contributions to the scattering diagram, with the product over n phonons conserving momentum, q = q 1 +...+q n . G(q) is the momentum-dependant electron-phonon interaction strength, and the greater phonon propagator in the zerotemperature limit is simply Thus, the non-diagonal part of the exciton Green's function differs from the diagonal part by the interaction term with a certain number of final state phonons. Despite the simple form, in the limit of small excitons Eq. (B2) accounts for all types of diagrams including those with vertex corrections.
The exciton propagator dressed by the exciton-phonon interaction can be expressed as Λ 0 (t) = L(t)e C(t) , where L is the bare exciton Green's function. We accounted for the exciton lifetime by including a damping factor e −γt/2 in L. Note that this expression contains half of the inverse lifetime, γ/2. We used a value γ/2 = 0.15 eV, slightly adjusted from the core-hole lifetime as determined by x-ray photoelectron spectroscopy mea-surements (0.1 eV) [39] to account for the presence of the extra electron (see the Supplemental Material).
Finally, the cumulant function to second order in the exciton-phonon interaction can be written (t, q).
Substituting the exciton and phonon Green's functions, and taking the time integrals in Eq. (B3) analytically [48], we get C(t) = q G(q) 2 /(ω(q) 2 N )(e −iω(q)t + iω(q)t − 1). It is important to notice that the summation runs over the whole Brillouin zone without restriction, since the total momentum transfer of the contributing diagrams is zero. From the form of Eqs. (B2) and (B3), it is clear that for any final state configuration there are contributions to the cross-section from phonons of all momenta. Therefore, even the intensity of the one-phonon peak depends on G(q) throughout the Brillouin zone, not just at the experimental q point.
To evaluate Eq. (B3), we used the TO phonon dispersion calculated with density functional perturbation theory. Brillouin zone integrations were performed using a reduced set of weighted momentum points which reflect the symmetries of the graphite structure [49]. Having calculated all the contributions to Eq. (B2) in the time domain, we can then perform a Fourier transform at the given incident photon energy and plug it into Eq. (B1) to obtain the RIXS cross-section.

I. SELF-ABSORPTION CORRECTION
Self-absorption -the re-absorption by the sample of previously scattered photons -depends on the energies and angles of both the incident and outgoing x-rays, and must therefore be corrected to reliably compare spectra. We follow previous resonant inelastic x-ray scattering (RIXS) studies [1][2][3] and define a self-absorption correction factor In this expression, u = sin θ/ sin (2Θ − θ) is a geometrical term that depends on the scattering angle 2Θ = 154°, and sample angle θ = 20°(90°) for the π * (σ * ) spectra. The other depends on α 0 , the absorption coefficient before the resonance (taken at 284 eV for the π * resonance and 289 eV for the σ * resonance), and α (ω), the resonant absorption coefficient for x-rays with energy ω and polarization . The absorption coefficients can be obtained from XAS measured in total electron yield, shown in Fig. S1(a) and S2(a). With our incident x-rays polarized in the (h, 0, l) scattering plane (π polarization), we measure the absorption for c ( ⊥ c) using x-rays incident normal (grazing) to the sample surface. From these, we calculate the absorption coefficients as α π (ω i ) = XAS(ω i , ⊥ c) sin 2 θ + XAS(ω i , c) cos 2 θ, The resulting correction factors are shown in Fig. S1(b) and S2(b). The RIXS spectra are then corrected by dividing through by C , (2Θ, θ, ω i , ω o ).
The anisotropic absorption of graphite results in a significant polarization dependence of the correction factors, especially around the π * resonance. The polarization of the outgoing photons will depend on the symmetry of the phonon modes being excited, and was not determined in our experiment. As our model only predicts the relative intensities of the phonon peaks, however, we can normalize the spectra to the intensity of the strongest multi-phonon peak. As shown in Fig. S1(c) and S2(c), after normalisation the difference between the corrected spectra is minimal over the energy range of interest, and does not significantly affect the fits to our model. We therefore correct all spectra assuming that the outgoing photons are π polarised.

II. ASSIGNMENT OF PHONON MODES
In the main text, we assign the contributions to the RIXS spectra based on the known phonon dispersion of graphite [4]. Here we demonstrate this process for the two-phonon features. Fig. S3 shows spectra near the π * and σ * resonances plotted alongside the phonon

III. INTERMEDIATE-STATE LIFETIME
To model our RIXS spectra we adjusted the intermediate-state lifetime from the value γ/2 = 0.1 eV (half-width at half-maximum of the Lorentzian component of the spectral peak) measured by x-ray photoemission spectroscopy [5] due to the presence of the extra excited electron in RIXS. Figure S4 compares best-fit RIXS spectra at the π * and σ * resonances calculated using the bare core-hole lifetime 0.1 eV, to the adjusted value 0.15 eV. We can see that the adjusted lifetime gives better agreement with both experimental spectra, especially for the high-order peaks at the σ * resonance. This adjusted lifetime was used for all calculations shown in the main manuscript, successfully reproducing the energy dependence of the experimental spectra below both resonances.

IV. INDEPENDENT FITTING OF THE ENERGY-DEPENDENT SPECTRA
The calculated spectra shown in Fig. 3(c) of the main manuscript are the result of a simultaneous fit to the experimental spectra around the π * resonance, such that the only parameter of the model that is varied between the spectra is the incident energy.
We have also fitted the spectra at each energy independently to see how this affects the extracted electron-phonon (e-ph) interaction strengths. Figure S5 shows a comparison of the simultaneous and independent fits, with the resulting values of G(q) listed in Table S1.
Both sets of calculated spectra have an acceptable level of agreement with the experimental fits to the experimental RIXS spectra (black points with error bars) around the π * resonance, plotted alongside the XAS with arrows indicating the energies.
TABLE S1. G(q) extracted from simultaneous and independent fits to the RIXS spectra around the π * resonance, at the momenta Γ, q exp ≈ (0.08, 0.00), K − q exp , and K. data given the noise level, but the independent fits do offer better agreement especially for the spectrum at 286 eV. As we note in the main text, we expect composition of the intermediate electronic state and therefore the momentum dependence of G to vary in the high-energy tail of the XAS peak, providing a natural explanation for this. While Table   S1 shows some variation in the extracted interaction strengths, the general trend of G(q) through the Brillouin zone -with G(K) > G(Γ) >> G(K − q exp ) ≥ G(q exp ) -is valid irrespective of the fitting procedure.

V. COMPARISON WITH PREVIOUS RIXS MEASUREMENTS ON GRAPHITE
Feng et al. recently presented lower-resolution RIXS measurements around the σ * resonance of graphite [3], interpreting their data with reference to Raman spectroscopy results.
Their analysis implicitly assumed, however, a coupling of phonons to electrons near the Fermi level, leading them to constrain G(q) to be non-zero only at the Γ and K points.
While this is reasonably accurate for the π * resonance, we show in the main text that RIXS at the σ * resonance probes the interaction of phonons with elections away from the Fermi level, whose strength is significant throughout the Brillouin zone. As well as overlooking the possibilities afforded by tuning the incident photon energy, this serves as a warning that careful consideration of the states involved in the RIXS process is essential for correct interpretation of experimental spectra, especially if the energy resolution is lacking.
We compare our best-fit G(q) at the σ * resonance (blue), to that obtained by Feng et al.
(orange) in Fig. S6(a). From the resulting normalized spectra in Fig. S6(b), we see that the intensities depend on the ratio of G to lifetime). We therefore see that, even for the relatively flat optical modes studied here, accounting for intermediate-state mixing strongly affects the momentum-dependence of G extracted from the data.

VI. MOMENTUM DEPENDENCE OF THE RIXS SPECTRA
The maximum momentum transfer is limited by the low energy of the carbon K edge.
Additionally, the RIXS cross-section falls off rapidly for sample angles away from θ = 0°( grazing incidence) at the π * resonance, and away from θ = 90°(normal incidence) at the σ * resonance. In the main manuscript, we show how these restrictions can be circumvented by analysing the contributions to the multi-phonon peaks, which can come from phonons throughout the Brillouin zone. Figure S7 shows the momentum-dependence of the RIXS spectra at the σ * resonance along the (h, 0, l) direction. The in-plane momentum was varied by changing the sample angle θ while keeping the spectrometer angle 2Θ( = 2 × θ) fixed, so the total momentum transfer q ∝ 4h 2 /(3a 2 ) + l 2 /c 2 is constant. While the overall intensity falls off rapidly with decreasing θ, it can be seen that the relative intensities of the phonon peaks change only minimally, supporting our conclusion that G(q) varies gradually away from the Γ point.
The low intensity of the phonon features at the π * resonance, along with the strong falloff in intensity with increasing θ, precluded measurement of a momentum-dependence with reasonable statistics.