Electron spin dynamics during microwave pulses studied by 94 GHz chirp and phase-modulated EPR experiments

Electron spin dynamics during microwave irradiation are of increasing interest in electron paramagnetic resonance (EPR) spectroscopy, as locking electron spins into a dressed state finds applications in EPR and dynamic nuclear polarization (DNP) experiments. Here, we show that these dynamics can be probed by modern pulsed EPR experiments that use arbitrary waveform generators to produce shaped microwave pulses. We employ phase-modulated pulses to measure Rabi nutations, echoes, and echo decays during spin locking of a BDPA (1,3-bisdiphenylene-2-phenylallyl) radical at 94 GHz EPR frequency. Depending on the initial state of magnetization, different types of echoes are observed. We analyze these distinct coherence transfer pathways and measure the decoherence time T_2ρ, which is a factor of 2–3 times longer than T_m. Furthermore, we use chirped Fourier transform EPR to detect the evolution of magnetization profiles. Our experimental results are well reproduced using a simple density matrix model that accounts for T_2ρ relaxation in the spin lock (tilted) frame. The results provide a starting point for optimizing EPR experiments based on hole burning, such as electron–nuclear double resonance or ELectron–electron DOuble Resonance (ELDOR)-detected NMR.

1 Introduction

Long microwave (MW) pulses with a duration of several microseconds are building blocks for pulse sequences in electron paramagnetic resonance (EPR) spectroscopy. They serve as preparation pulses to drive forbidden EPR transitions in hyperfine spectroscopy experiments like ELectron–electron DOuble Resonance (ELDOR)-detected nuclear magnetic resonance (EDNMR) (Schosseler et al., 1994) and can also lock spins into a dressed state, characterized by a distinct reorganization of their energy states (Cohen-Tannoudji and Reynaud, 1977; Jeschke, 1999). This enables electron–nuclear polarization transfer via cross-polarization (CP) under modified Hartmann–Hahn matching conditions as NOVEL or eNCP (Henstra et al., 1988; Weis et al., 2000; Weis and Griffin, 2006; Pomplun et al., 2008; Rizzato et al., 2013; Bejenke et al., 2020). Spin locking (SL) can also be harnessed to decouple spins from their surroundings, resulting in longer decoherence times of electron spins in EPR and optically detected magnetic resonance due to decoupling of nuclear interactions (Jeschke and Schweiger, 1997; Rizzato et al., 2023). SL experiments were recently proposed for electron–electron distance determination, where the upper limit of detectable distances is determined by the electron decoherence time (Wili et al., 2020). The altered energy landscape of electron spins locked by a MW field and under relaxation results in distinct spin dynamics, which are not yet fully explored.

Figure 1 Overview of the pulse sequences applied throughout this work and connection to the computational approach. While the spin lock experiments (a) were used as an analytical tool to monitor spin dynamics during MW irradiation, chirp echo EPR provided insights into the spin state following long, rectangular MW pulses. The results were analyzed using density matrix simulations. Additional calculations based on the Spinach library (Hogben et al., 2011) and the Bloch equations (Bloch, 1946) were used for validation and visualization (c).

One important aspect is the effect of decoherence during spin locking (i.e., T_2ρ relaxation), which has recently been discussed in the context of EPR spectroscopy by Wili et al. (2020). Some examples of treating this process exist in the literature. Most of these approaches rely on a phenomenological description of decoherence by rate equations but differ in the chosen interaction frame, which in turn determines the spin states that are affected by relaxation. For example, in a theoretical treatment of CP-ENDOR (Bejenke et al., 2020), we approximated the final state after a long SL pulse by propagating the spin density matrix in the dressed state yet without relaxation terms, followed by a heuristic implementation of T_2ρ that eliminates off-diagonal matrix elements. Hovav et al. (2010) described relaxation during MW irradiation in dynamic nuclear polarization experiments using density matrix propagation in the eigenframe of the equilibrium Hamiltonian without MW irradiation. Similarly, Hajduk et al. (1993) and Cox et al. (2017) proposed to model the effect of NMR shaped pulses or EDNMR preparation pulses with the Bloch equations in the rotating frame, including T_m relaxation. In contrast, Desvaux et al. (1994), De Luca et al. (1999), and Michaeli et al. (2004) expressed relaxation rates during irradiation in the eigenframe of the complete Hamiltonian including an irradiation term (i.e., the tilted frame) when analyzing the effect of dipolar interaction, molecular motion, and chemical exchange to the overall coherence decay rate. Very recently, Jeschke et al. (2025) calculated electron T_2ρ decays not using rate equations but from both the analytical pair product approximation and cluster correlation expansion. Their results show that a rigorous, quantitative treatment of decoherence during SL is currently challenging. Therefore, phenomenological treatment of relaxation using rate equations is still highly relevant to describe SL experiments.

This paper presents our efforts to gain insight into the electron spin dynamics of long MW pulses ($T_{\mathrm{2}\mathit{\rho }}\lesssim t_{\mathrm{pulse}}\lesssim T_{\mathrm{1}\mathit{\rho }}$). We focused on rationalizing decoherence effects and, ultimately, the trajectories of spins during MW irradiation to determine correct excitation profiles in EPR hole-burning experiments. Simultaneously, the experiments showcase the applicability of shaped pulses at our commercial W-band spectrometer. In the literature, most shaped-pulse experiments were performed at MW frequencies lower or equal to Q-band frequencies (34 GHz) (Endeward et al., 2023). The use of shaped MW pulses at high frequencies has been so far constrained to inversion pulses (Bahrenberg et al., 2017; Kuzhelev et al., 2018) or highly specialized experimental setups (Kaminker et al., 2017; Subramanya et al., 2022; Quan et al., 2023).

Experimentally, we followed two main routes. On the one hand, we performed pulse experiments during SL by coherently manipulating the locked spins using periodic phase modulations (PM) of the locking field (Grzesiek and Bax, 1995; Jeschke, 1999; Chen and Tycko, 2020; Wili et al., 2020). We report 94  GHz SL Rabi nutations, echoes, and echo decays measured with pulse sequences established by Wili et al. (2020) that were analyzed with electron spin dynamic simulations (Fig. 1a). On the other hand, we applied chirp echo Fourier transform (FT) EPR spectroscopy (CHEESY) (Wili and Jeschke, 2018) to measure the excitation profiles of long, rectangular MW pulses in the regime between Rabi nutation behavior and steady-state conditions (Fig. 1b). While the first class of experiments allows for observation of the dynamics during MW irradiation, the second class enables observation of the spin state immediately after a MW pulse. Therefore, the two approaches provide complementary information that we can use as a basis for density matrix simulations of an isolated electron spin under MW irradiation. Using the T_2ρ values obtained from the PM experiments, we were able to simulate the CHEESY excitation profiles. Comparison with simulations using the well-established Spinach library (Hogben et al., 2011) showed consistent results when using the suited interaction frame. Together with calculations based on the Bloch equations, we demonstrate that, for MW irradiation periods on the order of T_2ρ, calculating in the correct frame of reference is crucial for describing the experimental data.

2 Theoretical description

In this section, we provide the spin Hamiltonians used for the simulations of PM experiments and chirp echo profiles in their respective interaction frames. As will become evident from the results, a theoretical description based on the simple model of an isolated $S=\mathrm{1}/\mathrm{2}$ electron spin under MW irradiation, including a MW frequency offset, is sufficient to reproduce our experimental results. In the following, we refer to a spin system during free evolution as the bare state in contrast to spins subject to continuous MW irradiation, which are in the dressed or SL state. We label laboratory frame operators without a prime, rotating frame operators with a single prime (^′), and tilted frame (eigenframe during microwave irradiation) operators with a double prime (^′′). Operators in the nutating frame (interaction frame for PM pulses during spin lock) are denoted with an asterisk (^*). Furthermore, we indicate how the time evolution of the density matrix in the different frames can be computed under the effect of decoherence. All equations are expressed in angular frequency units (ω), while experimental results are reported in frequency units ($\mathit{\nu }=\mathit{\omega }/\mathrm{2}\mathit{\pi }$) as recorded.

2.1 Spin Hamiltonian

The laboratory frame Hamiltonian $\widehat{\mathcal{H}}$ for an isolated electron spin under irradiation with linearly polarized MW radiation in the x direction contains the electron Zeeman and the microwave irradiation terms:

$$\begin{array}{}\text{(1)}& \widehat{\mathcal{H}}\left(t\right)={\mathit{\omega }}_{\mathrm{0}}{\widehat{S}}_z+\mathrm{2}{\mathit{\omega }}_{\mathrm{1}}\cos ({\mathit{\omega }}_{\mathrm{MW}}t+{\mathit{\varphi }}_{\mathrm{MW}}){\widehat{S}}_x,\end{array}$$

where ω₀ is the Larmor frequency; ω_MW and ϕ_MW are the MW frequency and phase, respectively; and ω₁ is the Rabi frequency proportional to the MW B₁ field. To eliminate the time dependence, $\widehat{\mathcal{H}}\left(t\right)$ is transformed into a rotating frame with ${\widehat{U}}_{\mathrm{1}}\left(t\right)=\exp (-i{\mathit{\omega }}_{\mathrm{MW}}t{\widehat{S}}_z)$ (Abragam, 1961; Eckardt, 2017). Subsequent application of the rotating-wave approximation (RWA) yields the time-independent rotating frame Hamiltonian:

$$\begin{array}{}\text{(2)}& {\displaystyle }{\widehat{\mathcal{H}}}^{\prime }={\widehat{U}}_{\mathrm{1}}^{\mathit{†}}\left(t\right)\widehat{\mathcal{H}}\left(t\right){\widehat{U}}_{\mathrm{1}}\left(t\right)-i{\widehat{U}}_{\mathrm{1}}^{\mathit{†}}\left(t\right){\displaystyle \frac{\mathrm{d}{\widehat{U}}_{\mathrm{1}}\left(t\right)}{\mathrm{d}t}}\text{(3)}& {\displaystyle }\stackrel{\mathrm{RWA}}{\approx }{\mathrm{\Omega }}_S{\widehat{S}}_{z^{\prime }}+{\mathit{\omega }}_{\mathrm{1}}\cos \left({\mathit{\varphi }}_{\mathrm{MW}}\right){\widehat{S}}_{x^{\prime }}+{\mathit{\omega }}_{\mathrm{1}}\sin \left({\mathit{\varphi }}_{\mathrm{MW}}\right){\widehat{S}}_{y^{\prime }},\end{array}$$

where ${\mathrm{\Omega }}_S={\mathit{\omega }}_{\mathrm{0}}-{\mathit{\omega }}_{\mathrm{MW}}$ is the frequency offset of the MW to the electron Larmor frequency. Unless stated otherwise, MW frequency irradiation is assumed to be applied along y^′ (i.e., ϕ_MW=90°), which is consistent with our experimental SL field.

2.1.1 PM pulses and nutating frame

In the rotating frame, the Hamiltonian ${\widehat{\mathcal{H}}}_{\mathrm{PM}}^{\prime }$ during a PM pulse is defined as in Eq. (4), including a periodic variation of the phase ϕ_MW(t) (Eq. 5) that is matched with ω₁ (Wili et al., 2020).

$$\begin{array}{}\text{(4)}& {\displaystyle }{\widehat{\mathcal{H}}}_{\mathrm{PM}}^{\prime }={\mathrm{\Omega }}_S{\widehat{S}}_{z^{\prime }}+{\mathit{\omega }}_{\mathrm{1}}\cos \left({\mathit{\varphi }}_{\mathrm{MW}}\left(t\right)\right){\widehat{S}}_{x^{\prime }}+{\mathit{\omega }}_{\mathrm{1}}\sin \left({\mathit{\varphi }}_{\mathrm{MW}}\left(t\right)\right){\widehat{S}}_{y^{\prime }},\text{(5)}& {\displaystyle }{\mathit{\varphi }}_{\mathrm{MW}}\left(t\right)={\mathit{\varphi }}_{\mathrm{0}}+a_{\mathrm{PM}}\cos \left({\mathit{\omega }}_{\mathrm{PM}}t+{\mathit{\varphi }}_{\mathrm{PM}}\right),\end{array}$$

where ϕ₀ is the spin lock phase, a_PM the modulation depth of the phase modulation, and ω_PM its frequency. ϕ_PM is the phase of what will later be defined as the PM pulse. Assuming a small modulation amplitude a_PM and strong driving conditions (i.e., ω₁≫Ω_S), this expression can be transformed into the nutating frame by a rotation with ${\widehat{U}}_{\mathrm{PM}}\left(t\right)=\exp \left(-i{\mathit{\omega }}_{\mathrm{PM}}t{\widehat{S}}_{y^{\prime }}\right)$, yielding Eq. (6) after application of the rotating wave approximation (Grzesiek and Bax, 1995; Wili et al., 2020).

$$\begin{array}{}\text{(6)}& {\widehat{\mathcal{H}}}_{\mathrm{PM}}^{*}\stackrel{\mathrm{RWA}}{\approx }{\mathrm{\Omega }}_{\mathrm{d}}{\widehat{S}}_{y^{*}}+{\displaystyle \frac{{\mathit{\omega }}_{\mathrm{1}}{a}_{\mathrm{PM}}}{\mathrm{2}}}\left(-\cos \left({\mathit{\varphi }}_{\mathrm{PM}}\right){\widehat{S}}_{x^{*}}+\sin \left({\mathit{\varphi }}_{\mathrm{PM}}\right){\widehat{S}}_{z^{*}}\right)\end{array}$$

Here, ${\mathrm{\Omega }}_{\mathrm{d}}={\mathit{\omega }}_{\mathrm{1}}-{\mathit{\omega }}_{\mathrm{PM}}$ is a nutating frame offset that is unavoidable due to MW inhomogeneity. In contrast to the experiments performed by Wili et al. (2020), our experiments do not fulfill the abovementioned strong driving condition. Therefore, we used Eq. (4) for numerical simulations. Nevertheless, Eq. (6) provides an effective way of understanding PM experiments, demonstrating the analogy to a rotating frame experiment (see Eq. 3) with pulses along x^* or z^* and with a nutation frequency of $\frac{{\mathit{\omega }}_{\mathrm{1}}{a}_{\mathrm{PM}}}{\mathrm{2}}$. The validity of Eq. (6) under our experimental conditions is discussed in the Results section (see Sect. 4.1.1).

2.1.2 Tilted frame

During constant MW irradiation, a spin ensemble in the rotating frame is subject to an effective field ${\mathit{\omega }}_{\mathrm{eff}}={\mathrm{\Omega }}_S{\mathit{z}}^{\prime }+{\mathit{\omega }}_{\mathrm{1}}{\mathit{y}}^{\prime }$ described by the vector sum of the spin offsets and the microwave field. Therefore, ${\widehat{\mathcal{H}}}^{\prime }$ is not diagonal. Diagonalization of ${\widehat{\mathcal{H}}}^{\prime }$ leads to a rearrangement of the energy states due to the altered quantization axis (Rizzato et al., 2013) and is achieved by a rotation of ${\widehat{\mathcal{H}}}^{\prime }$ with ${\widehat{U}}_{\mathrm{2}}$ around x^′ into a so-called tilted or spin-locked frame (Eq. 7). Note that, in this case, we apply the rotations on the operator and not the coordinate frame so that ${\widehat{\mathcal{H}}}^{\prime \prime }\left(t\right)$ is obtained by Eq. (8) (Ernst et al., 1998).

$$\begin{array}{}\text{(7)}& {\displaystyle }{\widehat{U}}_{\mathrm{2}}=\exp \left(-i\mathit{\theta }{\widehat{S}}_{x^{\prime }}\right),\text{(8)}& {\displaystyle }{\widehat{\mathcal{H}}}^{\prime \prime }\left(t\right)={\widehat{U}}_{\mathrm{2}}{\widehat{\mathcal{H}}}^{\prime }\left(t\right){\widehat{U}}_{\mathrm{2}}^{\mathit{†}}.\end{array}$$

In Eq. (7), $\mathit{\theta }={\mathrm{Tan}}^{-\mathrm{1}}({\mathit{\omega }}_{\mathrm{1}}/{\mathrm{\Omega }}_S)$ is the tilt angle of the z^′′ axis in the Hamiltonian eigenframe with respect to z^′, defined such that $\mathit{\theta }\in [\mathrm{0}\mathit{°},\mathrm{180}\mathit{°}]$ for ω₁>0. Tan⁻¹() is the four-quadrant inverse tangent (e.g., implemented in MATLAB as the atan2() command) that returns the correct quadrant of θ by using ω₁ and Ω_S as separate input values (Bejenke et al., 2020). The transformation yields the tilted-frame Hamiltonian ${\widehat{\mathcal{H}}}^{\prime \prime }$ where z^′′ is aligned with ω_eff (Eq. 9).

$$\begin{array}{}\text{(9)}& {\widehat{\mathcal{H}}}^{\prime \prime }={\mathit{\omega }}_{\mathrm{eff}}{\widehat{S}}_{z^{\prime \prime }}\end{array}$$

This transformation is necessary to implement spin relaxation as discussed in the following.

2.2 Density operator and time evolution

The spin system evolution can be calculated using the density operator formalism. Here, the Hamiltonian ($\widehat{\mathcal{H}}$, ${\widehat{\mathcal{H}}}^{\prime }$, ${\widehat{\mathcal{H}}}^{\prime \prime }$, or ${\widehat{\mathcal{H}}}^{*}$) acts on an ensemble of electron spins described by a time-dependent 2×2 density operator $\widehat{\mathit{\rho }}\left(t\right)$ in the same frame of reference (i.e., $\widehat{\mathit{\rho }}\left(t\right)$, ${\widehat{\mathit{\rho }}}^{\prime }\left(t\right)$, ${\widehat{\mathit{\rho }}}^{\prime \prime }\left(t\right)$, or ${\widehat{\mathit{\rho }}}^{*}\left(t\right)$). For the frame transformation of the density matrix, the previously defined transformation matrices are used with the same sense of rotation as for the corresponding Hamiltonian (Bejenke et al., 2020) (e.g., ${\widehat{\mathit{\rho }}}^{\prime \prime }\left(t\right)={\widehat{U}}_{\mathrm{2}}{\widehat{\mathit{\rho }}}^{\prime }\left(t\right){\widehat{U}}_{\mathrm{2}}^{\mathit{†}}$). The time evolution of any $\widehat{\mathit{\rho }}$ is computed from the solution of the Liouville–von Neumann equation in the respective frame and under assumption of time-independent Hamilton operators; for example, in the rotating frame,

$$\begin{array}{}\text{(10)}& {\widehat{\mathit{\rho }}}^{\prime }\left(t\right)={\widehat{U}}^{\prime }\left(t\right){\widehat{\mathit{\rho }}}^{\prime }\left(\mathrm{0}\right){{\widehat{U}}^{\prime }}^{\mathit{†}}\left(t\right),\end{array}$$

where ${\widehat{U}}^{\prime }\left(t\right)=\exp (-i{\widehat{\mathcal{H}}}^{\prime }t)$ is a time-dependent propagator, and ${\widehat{\mathit{\rho }}}^{\prime }\left(\mathrm{0}\right)$ is the initial density operator. To include relaxation into the time evolution of ${\widehat{\mathit{\rho }}}^{\prime }$, calculations are performed in Liouville space (Hovav et al., 2010). For this purpose, the Hamiltonian in Hilbert space is converted to the corresponding Liouvillian superoperator $\widehat{\widehat{L}}$ as shown in Eq. (11) for ${\widehat{\mathcal{H}}}^{\prime }$. The uppercase T denotes the transposed matrix and, 𝟙 is the 2×2 unit matrix (Kuprov, 2023).

$$\begin{array}{}\text{(11)}& {\widehat{\widehat{L}}}^{\prime }=\mathbb{1}\otimes {\widehat{\mathcal{H}}}^{\prime }-{{\widehat{\mathcal{H}}}^{\prime }}^{\mathrm{T}}\otimes \mathbb{1}\end{array}$$

The respective density operator in Hilbert space is converted into a four-element density vector operator $\widehat{\mathit{\rho }}\left(t\right)$. The ordering of the matrix elements into a vector is shown in Eq. (12) where α and β are spin eigenstates (Ernst et al., 1987; Gyamfi, 2020).

$$\begin{array}{}\text{(12)}& {\widehat{\mathit{\rho }}}^{\prime }\left(t\right)=\left(\begin{array}{c}{\widehat{\mathit{\rho }}}_{\mathit{\alpha }\mathit{\alpha }}^{\prime }\left(t\right)\\ {\widehat{\mathit{\rho }}}_{\mathit{\alpha }\mathit{\beta }}^{\prime }\left(t\right)\\ {\widehat{\mathit{\rho }}}_{\mathit{\beta }\mathit{\alpha }}^{\prime }\left(t\right)\\ {\widehat{\mathit{\rho }}}_{\mathit{\beta }\mathit{\beta }}^{\prime }\left(t\right)\end{array}\right)\end{array}$$

The evolution of the density vector operator is then evaluated using a propagation superoperator ${\widehat{\widehat{U}}}^{\prime }\left(t\right)$ and the solution of the Liouville-von Neumann equation in Liouville space:

$$\begin{array}{}\text{(13)}& {\widehat{\mathit{\rho }}}^{\prime }\left(t\right)={\widehat{\widehat{U}}}^{\prime }\left(t\right){\widehat{\mathit{\rho }}}^{\prime }\left(\mathrm{0}\right)=\exp \left(-i{\widehat{\widehat{L}}}^{\prime }t\right){\widehat{\mathit{\rho }}}^{\prime }\left(\mathrm{0}\right).\end{array}$$

Starting from a Boltzmann-populated equilibrium density operator $\widehat{\mathit{\rho }}\left(\mathrm{0}\right)={\widehat{\mathit{\rho }}}^{\prime }\left(\mathrm{0}\right)\propto {\widehat{S}}_{z^{\prime }}$, Eq. (13) allows for numerical calculation of electron spin ensemble evolution under MW irradiation using matrix exponentials with standard mathematics software such as MATLAB. For better readability, we define the initial density matrix as positive (i.e., ${\widehat{\mathit{\rho }}}^{\prime }\left(\mathrm{0}\right)={\widehat{S}}_{z^{\prime }}$) so that the sign of the quantum mechanical expectation values is the same as of corresponding macroscopic magnetization vectors (i.e., $\mathit{M}\left(\mathrm{0}\right)=(\mathrm{0},\mathrm{0},M_{\mathrm{0}}{)}^{\mathrm{T}}$). Note that, due to the negative gyromagnetic ratio of the electron, the equilibrium density matrix in the high-temperature approximation is ${\widehat{\mathit{\rho }}}^{\prime }\left(\mathrm{0}\right)=\mathbb{1}/\mathrm{2}-\left|R\right|{\widehat{S}}_{z^{\prime }}$, where R is a thermodynamic population factor (Bejenke et al., 2020). $\mathbb{1}/\mathrm{2}$ does not evolve and is neglected in the following. In the tilted frame, the equilibrium density matrix transforms to ${\widehat{\mathit{\rho }}}^{\prime \prime }\left(\mathrm{0}\right)={\widehat{U}}_{\mathrm{2}}{\widehat{\mathit{\rho }}}^{\prime }\left(\mathrm{0}\right){\widehat{U}}_{\mathrm{2}}^{\mathit{†}}=\cos \left(\mathit{\theta }\right){\widehat{S}}_{z^{\prime \prime }}-\sin \left(\mathit{\theta }\right){\widehat{S}}_{y^{\prime \prime }}$.

To describe relaxation during MW irradiation, relaxation terms are added in the tilted frame and, thus, acting during propagation of the ${\widehat{\mathit{\rho }}}^{\prime \prime }\left(t\right)$ density matrix. For this, a relaxation superoperator is defined using a phenomenological decoherence rate (Feintuch and Vega, 2017). In our simulations, we consider only transverse relaxation since the effect of longitudinal relaxation was found to be almost negligible on the timescale of our experiments. As relaxation occurs in the tilted frame, we name it T_2ρ in agreement with the literature (Desvaux et al., 1994; Michaeli et al., 2004; Wili et al., 2020), and the corresponding relaxation superoperator is ${\widehat{\widehat{R}}}_{\mathrm{2}\mathit{\rho }}$.

$$\begin{array}{}\text{(14)}& {\widehat{\widehat{R}}}_{\mathrm{2}\mathit{\rho }}=\left(\begin{array}{cccc}\mathrm{0}& & & \\ & -\frac{\mathrm{1}}{T_{\mathrm{2}\mathit{\rho }}}& & \\ & & -\frac{\mathrm{1}}{T_{\mathrm{2}\mathit{\rho }}}& \\ & & & \mathrm{0}\end{array}\right)\end{array}$$

Evolution of the spin system can then be calculated by solving the Liouville-von Neumann equation (Eq. 13) with a modified propagation superoperator ${\widehat{\widehat{U}}}^{\prime \prime }\left(t\right)$ that includes ${\widehat{\widehat{R}}}_{\mathrm{2}\mathit{\rho }}$ (Feintuch and Vega, 2017).

$$\begin{array}{}\text{(15)}& {\widehat{\widehat{U}}}^{\prime \prime }\left(t\right)=\exp \left(\left(-i{\widehat{\widehat{L}}}^{\prime \prime }+{\widehat{\widehat{R}}}_{\mathrm{2}\mathit{\rho }}\right)t\right)\end{array}$$

Note that the effect of ${\widehat{\widehat{R}}}_{\mathrm{2}\mathit{\rho }}$ in Eq. (15) is the exponential dampening of the coherent terms ${\widehat{\mathit{\rho }}}_{\mathit{\alpha }\mathit{\beta }}^{\prime \prime }\left(t\right)$ and ${\widehat{\mathit{\rho }}}_{\mathit{\beta }\mathit{\alpha }}^{\prime \prime }\left(t\right)$ in the density vector operator ${\widehat{\mathit{\rho }}}^{\prime \prime }\left(t\right)$ of Eq. (12). Only if the density operator is in the eigenbasis of the full spin Hamiltonian will these entries correspond to pure coherence. Otherwise, they correspond to mixtures of populations and coherences, which changes the relaxation behavior (vide infra).

3 Materials and methods

3.1 EPR sample

A standard powder sample of 0.1 % ($m/m$) 1,3-bisdiphenylene-2-phenylallyl (BDPA) in a polystyrene matrix was used for EPR experiments throughout this work. At W-band frequencies, the powder spectrum of this radical features a single Gaussian line with a full width at half maximum (FWHM) of approximately 8.8 G or 24 MHz, due to unresolved hyperfine couplings below 10 MHz (see Fig. A7b). Measurements were done at temperatures between 50 and 100  K.

3.2 Spectrometer configuration and characterization

All experiments were performed on a commercial Bruker E680 W-band EPR spectrometer equipped with an ENDOR resonator (model EN600-1021H) and a liquid helium cryostat (Oxford Instruments). Pulse sequences were set up with the pulse programming language PulseSPEL implemented in Bruker's Xepr operating software. Shaped pulses were generated using a Bruker SpinJet arbitrary waveform generator (AWG) with a 1.6 GS s⁻¹ (where GS represents gigasamples) sampling rate and ±400 MHz bandwidth. In our setup, W-band frequency is reached via a three-step, heterodyne mixing scheme. First, a constant frequency ν_x in the range of 9.6±0.4 GHz is generated by an X-band source. Second, this carrier wave is mixed with the two AWG output channels I(t) and Q(t) for the real and imaginary parts, respectively, using an IQ mixer. Third, this shaped X-band pulse is mixed with the output of a phase-locked oscillator at ν_PLO=84.5 GHz to reach the pulse shape at W-band frequencies (see Sect. A1). As the AWG frequency is centered at zero, active compensation of the input amplitude imbalance and filtering are applied to minimize the contribution of ${\mathit{\nu }}_{\mathrm{LO}}={\mathit{\nu }}_x+{\mathit{\nu }}_{\mathrm{PLO}}$ (“LO leakage”) and mirror frequencies. The pulses are amplified with a solid-state power amplifier (SSPA, Quinstar) with 2 W output power and transferred to the resonator through a circulator. Spin echoes are recorded after demodulation with ν_PLO and subsequent amplification via a low-noise amplifier. Quadrature detection is done using ν_x as a continuous reference. To characterize our experimental conditions, we analyzed the phase stability of our setup following Endeward et al. (2023). The experiment (see Sect. A2) showed that the data suffer from phase noise that can, however, be averaged out using multiple scans. We acquired the MW resonator profile in a range of 600 MHz using standard Rabi nutation experiments measured at different local oscillator frequencies ν_LO (see Fig. A3a). At 50 K, we obtained a representative resonator bandwidth of approximately 140 MHz and a maximal Rabi frequency of ν₁≈15 MHz. Both values depend on the tuning conditions. Likewise, the amplification curve of the SSPA was measured with Rabi nutation experiments for the whole range of input power levels. This experiment showed clear nonlinearity of the MW amplification that was compensated for during pulse shape generation (see Fig. A3b).

3.3 Experiments with phase-modulated pulses

If not denoted otherwise, all phase-modulated (PM) experiments during SL were initiated by a rectangular $\mathit{\pi }/\mathrm{2}$ pulse immediately followed by a 90° phase-shifted SL pulse. During this SL, the PM pulses were applied (see Fig. A4). The locked spins were read out by a rectangular π pulse in the rotating frame, after a delay τ, which led to a spin echo at 2τ after the SL pulse. Integration of this echo yielded the detected signal.

For each data point of a parameter sweep (e.g., t_PM or τ₂, vide infra) during a PM experiment, a pulse shape file ($I/Q$ values), including the spin lock as well as all PM pulses, was generated with a home-written MATLAB script, indexed, and uploaded to the AWG. The complete experiment, including rectangular detection pulses and delays, was programmed in PulseSPEL. Analogously, PM pulse phase cycling (PC) was achieved by generating individual pulse shapes for each combination of PM phases. Due to the memory limitations of the SpinJet AWG, the pulse shape files were generated with 10 ns time resolution of individual $I/Q$ values, which was well below the AWG sampling rate but sufficient to encode the PM pulses for the maximal modulation frequency with only minimal distortions of the pulse shape (see Fig. A4). PM pulses were optimized with the experimental procedure reported by Wili et al. (2020). The modulation frequency ω_PM was calibrated by monitoring the effect of a long PM pulse with low modulation depth a_PM=0.04 on the intensity of the final spin echo as a function of ω_PM. The PM frequency with the strongest effect was chosen for the following experiments as this corresponds to the matching condition ω_PM=ω₁ (see Sect. A5). For the SL echoes, a_PM=0.4 was used, and the PM pulse lengths were determined from the SL Rabi nutations.

3.4 Generation of chirp pulses

The chirp echo EPR experiments were set up based on the work by Doll and Jeschke (2014) and Wili and Jeschke (2018). Additional information can be found in the work by Endeward et al. (2023) and Spindler et al. (2016). Chirp pulses with a frequency range of 200 MHz and lengths of 1 µs and 0.5 µs were generated with a home-written MATLAB routine. The pulse shapes featured a nonlinear frequency sweep as a compensation for the previously measured resonator profile. This ensured approximately constant critical adiabaticity Q_crit during the frequency sweep (Baum et al., 1985; Doll and Jeschke, 2014) (see Sect. A6 for details). Additionally, the pulses were treated with a wideband, uniform rate, smooth truncation (WURST) amplitude modulation to reduce Fourier transform (FT) artifacts caused by the pulse edges (Kupce and Freeman, 1995). Compensation of the nonlinearity of the SSPA yielded the final pulse shape. These pulse shape modifications are also implemented in the pulse function of MATLAB's EasySpin toolbox (Stoll and Schweiger, 2006). Experimental inversion profiles of representative chirp pulses demonstrated that these pulses could achieve broadband excitation and inversion (Fig. A6).

We are aware that there are multiple possibilities for further optimizing the chirp pulses (Endeward et al., 2023). Most importantly, we measured the resonator profile only once and at a different temperature than used for most of the experiments shown in this work. Compensation based on this profile will lead to imperfections as the resonator profile is dependent on temperature and tuning. Furthermore, there are other more advanced compensation schemes that take into account the complete impulse response function of the spectrometer and not just its amplitude response (Doll et al., 2013; Endeward et al., 2023). However, Doll and Jeschke (2017) showed that the effect of the latter is usually small. As this work is constrained to narrow-lined EPR samples where possible distortions caused by the chirp pulses are limited, we refrained from any further correction steps. Nevertheless, we expect that improved corrections might become necessary for samples with larger spectral widths.

3.5 CHEESY experiments

Chirp echo EPR spectroscopy (CHEESY) experiments were set up following the Böhlen–Bodenhausen scheme. Namely, 1 and 0.5 µs chirp pulses were used as $\mathit{\pi }/\mathrm{2}$ and π pulses, respectively, to achieve simultaneous refocusing of spins with different resonance frequencies (Böhlen et al., 1989). An eight-step phase cycle of the chirp pulses ensured that only the desired signal was detected (see Sect. B4). BDPA spectra were measured at a frequency offset ${\mathit{\nu }}_{\mathrm{off}}\approx -\mathrm{30}\mathrm{MHz}$ to avoid zero-frequency artifacts from LO leakage and baseline imperfections. The Rabi frequency of the chirp pulses was optimized by chirp nutations. For this purpose, the intensity and shape of the Fourier-transformed chirp echo were monitored as a function of ν₁ (Doll and Jeschke, 2014). The ν₁ value where the integral of the FT spectrum was maximal was chosen as the optimal setting. To measure the spectral intensity across the whole 200 MHz excitation bandwidth of the chirp pulses, a frozen sample of 4-hydroxy-2,2,6,6-tetramethylpiperidinyloxyl (TEMPOL) in H₂O / glycerol was used as it yielded a broad and structured FT spectrum. First, the Rabi frequency of the $\mathit{\pi }/\mathrm{2}$ pulse was optimized using a non-optimized π pulse. The experiment showed an optimum at ν₁=7.4 MHz which, integrated over the resonator profile, corresponded to Q_crit≈0.6. This value was close to the theoretical optimum of $Q_{\mathrm{crit}}=\mathrm{2}\ln \left(\mathrm{2}\right)/\mathit{\pi }\approx \mathrm{0.44}$ predicted by Jeschke et al. (2015) for a $\mathit{\pi }/\mathrm{2}$ chirp pulse (see Fig. A7a). Second, the optimized ν₁ was set for the $\mathit{\pi }/\mathrm{2}$ pulse, and the power of the π pulse was varied in the same way. The Rabi frequency with maximal signal was ν₁≈5.6 MHz. Considering that this pulse was shorter than the $\mathit{\pi }/\mathrm{2}$ pulse, its low optimal MW intensity was surprising. As discussed by Jeschke et al. (2015), this effect can be attributed to a combination of a dynamic Q_crit-dependent phase shift and ν₁ inhomogeneity. This observation suggested that the inversion efficiency of the pulse was low, presumably causing significant reduction of echo signal. Despite these nonideal conditions, we could achieve good agreement between the chirp echo FT spectrum of BDPA and the corresponding echo-detected (ED) EPR experiment (Fig. A7b). This demonstrated that we can detect narrow-lined EPR spectra using FT spectroscopy without creating distortions or artifacts. Additionally, we measured the offset dependence of the FT BDPA spectrum, which showed that the FT intensity decreases with larger frequency offsets, approximately following the shape of the resonator profile (Fig. A7c, d). Similarly, the broad EPR spectrum of TEMPOL was perturbed by the FT approach. Within the excitation bandwidth of the chirp pulses, the approximate convolution of the EPR spectrum and the resonator profile was obtained, with some additional distortions caused by pulse imperfections (see Fig. A8).

CHEESY experiments to determine pulse inversion profiles of rectangular pulses consisted of the analyzed pulse, followed by the optimized chirp echo. To generate a hole profile narrower than the EPR spectrum of BDPA, a low-power MW pulse (ν₁≈0.7 MHz) was used. The CHEESY spectra were obtained by subtracting the background echo without an initial pulse from the detected signal. The spectra were normalized to the background, and unity was added to afford the z^′ magnetization profiles, as demonstrated by Wili and Jeschke (2018). In contrast to the detection of pulse profiles with ELDOR (Hovav et al., 2015; Zhao et al., 2023), this procedure yielded the complete hole shape in a single experiment without any parameter sweep. Similar experiments for the analysis of inversion profiles via FT were reported by Tait and Stoll (2017) yet without the use of chirped pulses.

3.6 R-CHEESY

We performed refocused (R)-CHEESY experiments with a modified pulse sequence consisting of two chirped π pulses of equal length, which refocused coherence generated by the preceding analyzed pulse. In analogy to the Böhlen–Bodenhausen scheme for the standard chirp echo, the second refocusing chirp pulse was required to refocus spins with different offsets at the same time (Böhlen et al., 1989; Jeschke et al., 2015). Both pulses of this sequence were set with an optimized MW power, as described in the previous section. To filter all undesired echoes, a 32-step phase cycle was used (Bodenhausen et al., 1977; Cano et al., 2002) (see Sect. B5). Without a background signal for reference, all spectra of one data set were normalized to the global maximum. As the FT provided a complex signal, the coherence's x^′ and y^′ components could be separated into the real and imaginary signal parts. However, in contrast to the z^′ profiles where phasing could be done with respect to the Gaussian line of the unperturbed BDPA spectrum, the $x^{\prime }{y}^{\prime }$ profiles only contained the coherence generated by the analyzed pulse. This prevented an absolute phase determination, as there was no reference phase for the data processing. Therefore, the phase was adjusted for the individual FT spectra to maximize the absorptive (i.e., axially symmetric to zero) signal in the real channel. In this way, we achieved agreement with simulations (vide infra).

3.7 Fourier transformation of chirp echoes

Chirp echoes obtained from CHEESY spectroscopy were corrected for offset with a zeroth-order baseline correction. A symmetric FT window around the echo signal was chosen, and the experimental signal was apodized with a cosine window function. The result was filled with zeros, with the number of zeros equal to the number of data points, and shifted using the fftshift command in MATLAB. This time-domain signal was Fourier transformed using the fast Fourier transformation algorithm implemented in MATLAB. Subsequently, a constant phase correction with ϕ₀ was applied to afford the phased, complex FT spectrum. If not denoted otherwise, ϕ₀ was selected to maximize the absorptive contribution in the signal's real part.

For the CHEESY z^′ profiles, this process could be automated using a home-written MATLAB routine. In the case of R-CHEESY, the FT window had to be selected manually for each data set, because the automatized selection of the time-domain signal maximum was unstable. If not adjusted manually, this led to linear phase drifts in the resulting FT spectra.

3.8 Spin density matrix simulations

The home-written simulations reported in this work were performed based on the density operator formalism in Liouville space. The initial density operator ${\widehat{\mathit{\rho }}}^{\prime }\left(\mathrm{0}\right)$ for a $S=\mathrm{1}/\mathrm{2}$ electron spin system was defined according to Sect. 2.2. The Hilbert space Hamiltonian for each piecewise constant element of the simulated sequence was set up in the selected frame of reference. If necessary, ${\widehat{\mathit{\rho }}}^{\prime }\left(\mathrm{0}\right)$ was transferred to the same reference frame. Both the initial density matrix and the Hamiltonians were then transformed into their respective Liouville space representations, and the latter was used to calculate the corresponding propagators (Eqs. 11–13). The density matrix was propagated by consecutive application of the propagation superoperators (Eq. 13), and the final output was obtained as the expectation value of the Cartesian operator of interest. Simulation parameters were set in agreement with their experimental counterparts.

SL Rabi nutations (Sect. 4.1.1) and SL echoes (Sect. 4.1.2) were calculated in the rotating frame using ${\widehat{\mathit{\rho }}}^{\prime }\left(\mathrm{0}\right)={\widehat{S}}_{y^{\prime }}$ or ${\widehat{S}}_{z^{\prime }}$ for SL and without-preparation (WOP) SL experiments, respectively, as well as assuming Ω_S=0. While the delays between the PM pulses were simulated in a single step using ${\widehat{\mathcal{H}}}^{\prime }$ as defined in Eq. (3), the PM pulses were incremented with a step size Δt=1 ns as ${\widehat{\mathcal{H}}}_{\mathrm{PM}}^{\prime }$ is time dependent. In this way, a constant Hamiltonian could be assumed within these increments, as demonstrated by Laucht et al. (2016). Relaxation was neglected during these calculations. For each value of t_PM or τ₂, the result was obtained as the expectation value $\langle {\widehat{S}}_{y^{\prime }}\rangle$ (see Fig. 2c). In the case of the SL echo experiments, MW inhomogeneity was included in the SL echo simulations by repeating the calculation for different ν₁ values. The sum of these individual traces weighted with an experimental ν₁ profile yielded the depicted simulation results (details in Sect. A5). To calculate the effect of PC, the simulations were repeated with different phases ϕ_PM2, and the results were combined in the same way as in the experiment.

The CHEESY profiles were simulated starting from ${\widehat{\mathit{\rho }}}^{\prime }\left(\mathrm{0}\right)={\widehat{S}}_{z^{\prime }}$. A linear array of spin offsets ${\mathrm{\Omega }}_S/\mathrm{2}\mathit{\pi }$ was generated corresponding to the experimental ν−ν_off axis. For each value of Ω_S, the initial density matrix ${\widehat{\mathit{\rho }}}^{\prime }\left(\mathrm{0}\right)$ was transformed into the tilted eigenframe of the MW Hamiltonian ${\widehat{\mathcal{H}}}^{\prime }$ (Eqs. 3, 7, and 8). The final state after a pulse with length t_p was calculated using a single propagator ${\widehat{\widehat{U}}}^{\prime \prime }$ calculated from ${\widehat{\mathcal{H}}}^{\prime \prime }$ (Eq. 9) and including T_2ρ relaxation (see Eqs. 14 and 15). Again, MW inhomogeneity was included as the weighted sum of the results for different ν₁ values using the experimental ν₁ profile shown in Sect. A5. After back-transformation of ${\widehat{\mathit{\rho }}}^{\prime }\left(\mathrm{0}\right)$ into the rotating frame, the final result was obtained as the expectation values $\langle {\widehat{S}}_{z^{\prime }}\rangle$ for the CHEESY z^′ profiles and as $\langle {\widehat{S}}_{x^{\prime }}\rangle$ and $\langle {\widehat{S}}_{y^{\prime }}\rangle$ for the $x^{\prime }{y}^{\prime }$ profiles (see Fig. 5e).

3.9 Spinach simulations

Simulations using the Spinach library (Hogben et al., 2011) were done based on an adapted version of the holeburn.m function. The CHEESY spectrum was calculated for an anisotropically broadened powder EPR line (see Sect. A14) in three steps. First, the initial equilibrium spin system was propagated under a soft MW SL pulse including transversal relaxation, in the tilted frame. Second, an ideal, non-selective $\mathit{\pi }/\mathrm{2}$ rotation was applied that transformed all residual z^′ components into coherences. Third, simulation of the powder-averaged free induction decay (FID) and FT yielded the CHEESY spectrum. In order to filter the FID of the SL pulse, a two-step PC of the SL was applied in the simulations (i.e., ${\mathit{\varphi }}_{\mathrm{SL}}=[\mathrm{0}\mathit{°},\mathrm{180}\mathit{°}]$ and ${\mathit{\varphi }}_{\det }=[+,+]$). Due to the increased computational cost of the calculations, ω₁ inhomogeneity was not included. The results were compared to calculations using the unmodified holeburn.m function and rotating frame relaxation implemented as t1_t2.m (see Sects. 4.2.1 and A14).

4 Results

4.1 Spin-locked EPR experiments

We performed SL EPR experiments to probe the spin dynamics during MW irradiation.

4.1.1 SL Rabi nutations

Figure 2a shows the pulse sequence for SL Rabi nutation experiments, which can be used to determine the optimal PM pulse lengths for further SL EPR experiments. The sequence consists of an initial rectangular $\mathit{\pi }/\mathrm{2}$ pulse (x^′) that prepares electron coherence, which is subsequently locked by the SL pulse along y^′. The preparation and SL pulses were performed with ν₁=15 MHz, which was considerably less than the 100 MHz used by Wili et al. (2020) at 34 GHz and not sufficient to lock the entire BDPA spectrum with a width of ≈24 MHz in the transversal plane. The effects of this will be discussed later in this section.

Figure 2 (a) Pulse sequence for the SL Rabi nutation. The echo was integrated as a function of the PM pulse length t_PM. (b) Experimental SL Rabi nutation traces for ν₁=15 MHz with different modulation amplitudes a_PM. (c) Simulation approach for the SL Rabi nutation experiments by stepwise evolution of the rotating frame density operator during a phase-modulated pulse. (d, e) FT spectra of the experimental nutation traces (red) in the low-frequency region (d) and enlarged in the high-frequency region (e). Simulated nutations are shown in black, accurately reproducing the experimental peak positions. In (d), dotted lines mark the positions where a signal is expected for the respective SL Rabi frequency ${\mathit{\nu }}_{\mathrm{1}}\cdot a_{\mathrm{PM}}/\mathrm{2}$. The dotted line in (e) marks ν=2ν₁. For experimental parameters, see Sect. B1.

During SL and after an initial SL delay τ₀, a PM pulse was applied by periodic modulation of the SL phase for a variable length t_PM and with a modulation amplitude a_PM. The overall SL length τ₀+τ₁ was constant. The effect of the PM pulse was monitored by a bare-state (rotating frame) echo produced by a π pulse after a delay τ that refocused coherence in the rotating frame at 2τ after the SL pulse. Figure 2b shows the SL Rabi nutation traces for three different modulation amplitudes a_PM. The SL Rabi traces contained a dominant slow oscillation with a superimposed weaker and faster modulation. FT of the time traces (Fig. 2d) showed that the dominant, slower oscillation corresponded to a SL Rabi frequency of ${\mathit{\nu }}_{\mathrm{1}}{a}_{\mathrm{PM}}/\mathrm{2}$. This agrees with the theoretical predictions from the nutating frame transformation (Eq. 6). The faster oscillation frequency appeared at 2ν₁ plus two sidebands (Fig. 2e). We simulated the experiments using density matrix calculations (see Fig. 2c and Sect. 3.8), neglecting a bare-state offset by setting Ω_S=0. Figure 2d and e show almost quantitative agreement between the Fourier transformed experimental and simulated time traces. Interestingly, the fast oscillation also directly resulted from the simulation. We attribute this effect to a deviation from the RWA, leading to a minor yet visible modulation with the counter-rotating wave that precesses with 2ν₁. The simulations also qualitatively reproduced the sidebands of 2ν₁, which might be caused by frequency mixing. Comparison with the results by Laucht et al. (2016) shows that they reported similar fast frequency oscillations in their simulation for tilted-frame excitation by modulation of the frequency offset instead of phase modulations.

This good agreement without simulating any bare-state offset was unexpected, as the strong driving condition ω₁≫Ω_S is not met in our experiments. Simulations including a bare-state offset showed that severe distortions occurred for a single, sharp offset ΔΩ>ω₁ (see Fig. A10a). However, simulating a Gaussian distribution of Ω_S values that was centered at 0 MHz yielded a nearly unperturbed FT spectrum (Fig. A10b). This suggests that offset effects are averaged during the experiment, so the final spectrum agrees with the predictions from the ideal ω₁≫Ω_S case.

Figure 3 Pulse sequence for standard (a) and WOP (b) SL echoes. Real components of the quadrature detection signal of the SL echo (c) and WOP SL echo (d) spectra measured with different phase cycles. The detection phase was set to maximize y^′ magnetization in the real channel. Experimental traces are in different shades of red, and corresponding simulations are in black. Signals (i)–(iv) were assigned in (e) and belong to evolution pathways that are refocused at the beginning of PM3. In (c), the asterisks denote the residual signals (ii), (iii), and (iv). The imaginary part of the spectra with signals refocused at the end of τ₃ can be found in Fig. A11. For experimental parameters, see Sect. B2.

4.1.2 SL echoes

In the next step, we examined the feasibility of observing spin echoes during SL. Using a PM $\mathit{\pi }/\mathrm{2}$ pulse followed by a delay and a PM π pulse, the analog of a Hahn echo can be created during SL (Wili et al., 2020). The refocused SL coherence is read out by a third PM $\mathit{\pi }/\mathrm{2}$ pulse that realigns the magnetization with the SL field, which can be refocused with a bare-state echo (see Fig. 3a). We monitored the refocusing of SL coherence by variation of the delay τ₂ between the refocusing (PM2) and readout pulse (PM3). The SL pulse length was kept constant by simultaneous variation of τ₃. For ${\mathit{\tau }}_{\mathrm{2}}={\mathit{\tau }}_{\mathrm{1}}=\mathrm{3}\mathrm{\mu }\mathrm{s}$, refocusing was observed as a strong echo (Fig. 3c), which indicates that PM experiments are feasible under W-band experimental conditions, i.e., weaker MW pulses and considerable offset contributions. In addition to this, weak signals (Fig. 3c, asterisks) whose position depended on the length of the delay τ₀ were observed. According to Wili et al. (2020), these are caused by incomplete preparation of the coherent state $S_{y^{\prime }}$, leaving spins with residual $S_{z^{\prime }}$ components unaligned with the spin lock field.

In the same way as the $\mathit{\pi }/\mathrm{2}$ preparation pulse, the PM pulse rotation angles can deviate from the desired $\mathit{\pi }/\mathrm{2}$ and π values, so several signals can be produced by different effective rotation angles of the individual PM pulses. We used two- and four-step phase cycling (PC) of the phase-modulated inversion pulse PM2 to eliminate the small crossing echoes and two-step phase cycling of the phase-modulated readout pulse PM3 to filter non-locked components. A two-step PC of PM2 (${\mathit{\varphi }}_{\mathrm{PM}\mathrm{2}}=[\mathrm{0}\mathit{°},\mathrm{180}\mathit{°}]$) with opposite detection phase (${\mathit{\varphi }}_{\det }=[+,-]$) selects coherence transfer pathways in which PM2 changes the coherence order by one, i.e., only when it acts as a $\mathit{\pi }/\mathrm{2}$ pulse. In contrast, a four-step PC of PM2 (${\mathit{\varphi }}_{\mathrm{PM}\mathrm{2}}=[\mathrm{0}\mathit{°},\mathrm{90}\mathit{°},\mathrm{180}\mathit{°},\mathrm{270}\mathit{°}]$) with alternating detection phase (${\mathit{\varphi }}_{\det }=[+,-,+,-]$) selects a transfer pathway, where PM2 causes a change in coherence order by two, i.e., where it acts as a π pulse. As expected, the two-step PC of PM2 eliminated the main echo signal, whereas it was retained by the four-step PC (see Fig. 3c). Additionally, the small crossing echoes were affected by the different phase cycles, demonstrating that they arise from different coherence transfer pathways.

As incomplete preparation is supposed to be significant under W-band conditions, we analyzed the SL coherence pathways by omitting the $\mathit{\pi }/\mathrm{2}$ preparation step. This type of experiment without preparation (WOP) was already implemented for CP-ENDOR by Bejenke (2020), but the spin dynamics were not reported. Figure 3b depicts the modified pulse sequence, and panels (c) and (d) show the comparisons of the SL echo with and without preparation obtained by variation of τ₂ while keeping τ₀=2 µs and τ₁=3 µs fixed. Both experiments were performed using the same experimental parameters.

The WOP SL echo experiment produced signals at 1 µs (ii), 2 µs (iii), and 5 µs (iv) (Fig. 3d), which differs significantly from those observed in (c). Notably, the primary echo (i) is missing, while the observed echo positions coincide with the artifacts observed in (c). To rationalize this observation, we analyzed the refocusing times τ₂ and the effect of PC. The strongest WOP SL echo (iii) is observed at τ₂=τ₀, suggesting that it is created by a coherence that starts to evolve during τ₀ and refocuses as a stimulated echo. This is corroborated by the PC, as the signal is retained by a two-step PC of PM2 and vanishes upon four-step PC. With the same rationale, signals (ii) and (iv) can be assigned to a refocused echo and Hahn echo, respectively. Coherence diagrams visualizing the idealized evolution pathways (i.e., Ω_S=0) responsible for each signal are shown in Fig. 3e, along with their respective refocusing condition. These results are consistent with the spins having a Boltzmann distribution in the rotating frame corresponding to z^′ magnetization prior to SL. Upon start of the SL pulse, i.e., under MW irradiation, this becomes a coherence in the tilted frame, which immediately evolves with ν₁.

The peak positions (i)–(iv) as well as their line width and the effect of PC could be reproduced by density matrix simulations following the procedure described in Sect. 3.8 and Fig. 2c. The phase and intensity of the WOP signals agreed only partially. Repetition of the experiment under nominally the same conditions showed that the phase of the individual signals was not reproducible. Detection of the out-of-phase component by setting ${\mathit{\varphi }}_{\mathrm{PM}\mathrm{3}}=\mathit{\pi }/\mathrm{2}$ underlined that there is a phase instability for the WOP echoes (see Sect. A9). We suppose that this is the case because the turning angles of the PM pulses were not optimized for these specific coherence pathways.

Altogether, our results are consistent with the assumption that equilibrium z^′ magnetization becomes a coherence during MW irradiation. Hence, these results show that, under our experimental conditions, conventional spin locking with a preparation step and long MW pulses acting on spins initially in thermal equilibrium are closely related. By application of PM pulses, the spin states can be transferred from SL coherence to SL populations and vice versa.

4.1.3 Relaxation during spin lock

To measure relaxation times during SL, the delays τ₁ and τ₂ of the SL echo sequence (see Fig. 3a) were both increased so that refocusing was achieved throughout the experiment, while the evolution time in the transverse plane was increased. Simultaneously, τ₃ was decreased to ensure a constant SL length. The comparison of the decays obtained with ν₁=15 MHz and the bare-state Hahn echo decay can be found in Fig. 4a. The experimental traces were fitted using both a mono-exponential and a stretched-exponential function. We calculated the residuals, and statistical errors were determined by bootstrapping (for details, see Sect. A10). With this, we could compare the two decay models and evaluate the systematic error introduced when using the simplified mono-exponential model. The decay of coherence was approximately 2–3 times slower in the SL state than in the rotating frame with mono-exponential decay times of 3.6±0.2 and 1.46±0.02 µs, respectively (red and blue curves in Fig. 4a). Fitting with a stretched-exponential decay function led to $T_{\mathrm{2}\mathit{\rho }}=\mathrm{5.0}\pm \mathrm{0.3}$ and $T_{\mathrm{m}}=\mathrm{2.19}\pm \mathrm{0.03}\mathrm{\mu }\mathrm{s}$. Analysis of the residuals (bottom part in Fig. 4a) showed that there is a systematic deviation in the mono-exponential fit (red curve). Namely, the intensity is overestimated for t_p<10 µs and underestimated for t_p>10 µs. This systematic error was reduced when using the stretched-exponential fit (dark red curve), which is consistent with the literature concerned with bare-state echo decays, where stretched-exponential functions are often used. (Mims et al., 1961; Jahn et al., 2022).

Figure 4 SL relaxation measurements obtained with the standard SL echo sequence as depicted in Fig. 3 by variation of τ₁+τ₂. (a) SL relaxation trace measured for ${\mathit{\nu }}_{\mathrm{1}}={\mathit{\nu }}_{\mathrm{PM}}=\mathrm{15}\mathrm{MHz}$ (red) compared to the T_m relaxation trace in the bare state (blue). Experimental time traces are in gray, and the envelope function used for fitting the SL state decay is depicted as black dots. Mono-exponential fits for T_2ρ and T_m are in red and blue, respectively. Residual of the mono-exponential fit for T_2ρ in the bottom part (red) compared to the residuals from the stretched-exponential fitting (dark red). (b) T_2ρ data measured for different SL powers and comparison with T_m values obtained from Hahn echo decays with different ν₁ values and pulse lengths. The corresponding traces and fits can be found in Fig. A13. Error bars correspond to the statistical error, while the systematic error might be larger. For experimental parameters, see Sect. B3.

Investigation of the SL power dependence of T_2ρ (Figs. 4b and A13) showed no clear trend, and an overall increase in the transverse relaxation time by a factor of 2 to 4 compared to T_m=1.4 µs was obtained from a mono-exponential fit (red). The stretched-exponential fitting yielded T_2ρ values between 4 and 8 µs. These results show a similar trend as the low-field results by Wili et al. (2020) that were performed at much higher MW fields ν₁≈100 MHz, reporting an increase in T_2ρ relative to T_m by a factor of 4.5 for stretched-exponential fitting. These observations show potential for designing new pulse sequences that mitigate short decoherence times. The stretching exponents ξ for the stretched-exponential fits were between 0.9 and 1.7 (see Sect. A10). Importantly, there is a strong correlation between the T_2ρ value and ξ (Pearson correlation within the bootstrap samples above 0.9 for all traces). Therefore, the choice of the stretching coefficient affected the obtained decay time.

Analysis of Hahn echo decays measured with lower ν₁ (i.e., longer and more selective pulses) revealed increased T_m times (blue in Fig. 4b), indicating a contribution of instantaneous diffusion to the bare-state decoherence. Additionally, we measured the echo decay with the three- and four-pulse Carr–Purcell (CP) sequences (Carr and Purcell, 1954). Both CP3 and CP4 led to an increase in decoherence time by around 50 % relative to the two-pulse measurement with the same ν₁ (see Sect. A11).

Notably, the SL echo decays were perturbed by an oscillation for small τ₁+τ₂ (gray curve), which, to our knowledge, has not been reported in the literature. We attribute this to the edge effects of the PM pulses. A detailed analysis and quantitative explanation of these artifacts can be found in Sect. A12 but is of little consequence for the following results. We circumvented this effect by fitting the envelope function of the time trace (black dots in Fig. 4a). This method is also used to analyze Hahn echo decay traces perturbed by electron spin echo envelope modulations (ESEEM) (Soetbeer et al., 2021). However, the periodic oscillations in the residuals for ${\mathit{\tau }}_{\mathrm{1}}+{\mathit{\tau }}_{\mathrm{2}}<\mathrm{10}\mathrm{\mu }\mathrm{s}$ are indicative of small artifacts caused by this envelope fitting. Overall, these experiments provided us with an estimate for the decoherence time during SL over a broad range of MW powers. In the following sections, this will be used to model the spin dynamics during MW irradiation by calculating trajectories of our $S=\mathrm{1}/\mathrm{2}$ model system. For this, a transverse relaxation time of T_2ρ≈4 µs will be assumed as an approximation of the mono-exponential T_2ρ time. Due to the strong correlation of T_2ρ and ξ as well as the limitation of our simulation method to mono-exponential decays, we refrain from using the stretched-exponential results for the following calculations.

4.2 CHEESY detection of pulse excitation profiles

CHEESY experiments were performed to analyze the spin state immediately after long, rectangular MW pulses.

Figure 5 (a, b) CHEESY pulse sequence and schematic data processing. Inversion (z^′) profiles were obtained by subtracting a chirp echo background spectrum without the first pulse from the normalized hole profiles and vertically shifting the results by one (Wili and Jeschke, 2018). (c) Experimental (red) and simulated (black) inversion profiles for a low-power MW pulse (ν₁=0.77 MHz). For better visibility, the latter are vertically shifted by −0.1. Simulations were performed using the simulation routine depicted in (e) with a variation of t_p and Ω_S (see Sect. 3.8). Results from an analogous simulation omitting the tilted-frame transformation are shown in light blue. (d) Rabi nutation obtained from CHEESY inversion profiles at $\mathit{\nu }-{\mathit{\nu }}_{\mathrm{off}}=\mathrm{0}$ (red) and from the tilted-frame simulation (black). For experimental parameters, see Sect. B4.

4.2.1 z^′ profiles

We probed the residual z^′ magnetization after a soft MW pulse of varied length t_p with a chirp echo (see Fig. 5a). For each t_p, FT of the chirp echo transient and subtraction of a background chirp echo without additional pulse yielded the inversion profile as a hole burned in the initial z^′ magnetization (Fig. 5b). Figure 5c shows the obtained hole shape as a function of t_p (red).

For short irradiation times, the profiles were sinc-shaped, as expected for a rectangular pulse. A distinct Rabi oscillation was detected, i.e., a periodic variation of the hole depth at the center frequency $\mathit{\nu }-{\mathit{\nu }}_{\mathrm{off}}=\mathrm{0}$ as a function of t_p (Fig. 5d, red). FT of the Rabi oscillation trace yielded a nutation frequency ν₁=0.77 MHz. The magnitude of the oscillation decreased with increasing pulse length. Comparison with the experimental T_2ρ relaxation traces and simulations (see Fig. A17) showed that this decay was faster than the SL coherence decay. We attribute this fast decay to MW inhomogeneity, resulting in a broadened distribution of nutation frequencies. Simulations including both the experimentally determined T_2ρ relaxation rate and MW inhomogeneity could quantitatively reproduce the Rabi trace (black in Fig. 5d; see Sect. A13). For long MW pulses with t_p≳10 µs, a steady state was reached where the excitation profile no longer depended on t_p and resembled a Lorentzian-shaped saturation profile with FWHM≈1.6 MHz. This saturation hole shape agrees with observations of the central hole in standard EDNMR experiments with rectangular high-turning angle pulses and a full width at half maximum of 2ν₁ (Nalepa et al., 2014). The simulations also reproduced the transition from a t_p-dependent, sinc-shaped hole to a steady-state Lorentzian as well as the excitation bandwidth. Moreover, the dampening of the Rabi oscillation curve was simulated quantitatively. Crucially, the simulation results changed when the transformation in the tilted frame was omitted. In this case, the simulations for long t_p no longer matched the experiment, and an increased width of the simulated hole was observed (light blue in Fig. 5c). We attribute this to the fact that, if not in the eigenstate of the full MW Hamiltonian, the relaxation term in the Liouville space propagator does not act on pure coherences but on mixed coherence and population density matrix elements. This leads to incorrect dampening of population terms during decoherence.

The most apparent difference between experimental and simulated profiles in the tilted frame is the resolution in the frequency dimension, i.e., the degree to which the narrow sinc lobes for intermediate pulse lengths are resolved. While the oscillations are clearly visible in the simulations up to t_p=3 µs, they are hardly visible in the experimental spectra with t_p>1.6 µs. This is presumably caused by the intrinsic experimental FT resolution limit posed by the finite FT window during data processing.

Figure 6 CHEESY z^′ profiles simulated with the Spinach library (Hogben et al., 2011) using propagation in the tilted frame (black) and comparison with the experimental profiles (red). For better visibility, the simulations are vertically shifted by −0.1. Analogous calculations in the rotating frame using the holeburn.m function are in light blue. For experimental parameters, see Sect. B4.

We compared our results to simulations using the Spinach library (Hogben et al., 2011), where we explicitly calculated the hole burning in a powder-averaged EPR spectrum by Fourier transformation of a time-domain CHEESY signal (see Sects. 3.9 and A14). As in the home-written simulation routine, transformation in the tilted frame during propagation of the SL pulse was necessary to reproduce the experimental results (compare black and light blue lines in Fig. 6). Compared to Fig. 5, the line width agreed more accurately because the FT was explicitly included in the simulation (see Sect. 3.9). However, the dampening of the Rabi oscillation obtained from the intensity at $\mathit{\nu }-{\mathit{\nu }}_{\mathrm{off}}=\mathrm{0}$ was not quantitatively reproduced (compare Figs. A18 and 5d), because increased computational time prevented incorporation of MW inhomogeneity.

Comparison of tilted-frame simulations with different T_2ρ values showed that using a shorter relaxation time T_2ρ=1.4 µs could not reproduce the data as the decay of the side lobes was too fast. Using a longer relaxation time T_2ρ=10 µs did not affect the visual shape of the pulse profiles. In contrast, completely omitting T_2ρ could not reproduce the experiment as no steady state was reached, even for t_p=100 µs (see Fig. A19). So far, our simulations neglected hyperfine coupling. This is a strong simplification, as nuclear coupling, in particular its pseudo-secular component, is known to introduce additional evolution pathways during irradiation (Jeschke, 1996). Using the Spinach implementation, nuclear coupling could be introduced straightforwardly. Simulation with up to four coupled protons with A<10 MHz led to no change in the CHEESY profiles when using an appropriately expanded relaxation superoperator in the tilted frame (see Fig. A20c).

Overall, the Spinach simulations validated the home-written density matrix simulations by demonstrating that, under our conditions, including nuclear couplings has no visible effect on the z^′ profiles. At the same time, they highlighted that omitting MW inhomogeneity leads to an underestimation of the Rabi nutation decay. Including the FT step in the simulations improved the agreement between experimental and simulated z^′ profiles and allowed for a qualitative assessment of the effect of T_2ρ on reaching dampening of the sinc-shaped features of the profiles. However, the Spinach simulation requires longer computation time, preventing fast calculation and inclusion of MW inhomogeneity. This is the advantage of the home-written density matrix simulation, which we used as the main computational tool in the following.

4.2.2 $x^{\prime }{y}^{\prime }$ profiles

The standard CHEESY approach is limited to observing z^′ magnetization and hole-burning profiles. To observe the coherence generated by a MW pulse, we used a modified chirp sequence where the chirp echo was replaced with two chirped π pulses that form, together with the analyzed pulse, a refocused echo sequence (R-CHEESY; see Fig. 7a). This experiment is closely related to the ABSTRUSE sequence developed by Cano et al. (2002). FT of the hereby generated chirp echo yielded a complex signal. Its phase was adjusted so that the absorptive signal was maximized in the real component.

The results of the R-CHEESY experiment are shown in Fig. 7c (red). Oscillating components were observed for short pulse lengths t_p that changed to a steady state for longer irradiation in the same manner and timescale as observed in the CHEESY z^′ profiles. Monitoring the intensity of the real component at $\mathit{\nu }-{\mathit{\nu }}_{\mathrm{off}}=\mathrm{0}$ showed a continuous Rabi oscillation (Fig. 7d). FT of this trace yielded the Rabi frequency ν₁=0.63 MHz (these experiments were performed with slightly different tuning of the microwave cavity and different power settings compared to the z^′ profiles in the previous section). After long MW irradiation (t_p≫T_2ρ), the steady state exhibited negligible signal in the real part, yet a clear dispersive signal was detected in the imaginary component. This corresponded to no detectable magnetization perpendicular to the MW axis but to a distinct magnetization component parallel to ν_eff. To rationalize this, we simulated the $x^{\prime }{y}^{\prime }$ profiles using the procedure described in Fig. 5e and Sect. 3.8. Again, there was a good agreement between experimental and simulated data for both the initial Rabi oscillation period and the steady state (Fig. 7c, d, black). As in the experiment, no contribution remained in the real part of the steady-state R-CHEESY signal, while a dispersive Lorentzian line was retained in the imaginary part. This can be explained by the T_2ρ relaxation that leads to a decay of coherence in the tilted frame, i.e., of magnetization in the $x^{\prime \prime }{y}^{\prime \prime }\approx x^{\prime }{z}^{\prime }$ plane. Therefore, only the component parallel to ν_eff which decays with the longer T_1ρ is visible in the steady state.

Figure 7 (a) R-CHEESY pulse sequence. (b) Pulse sequence for measuring longitudinal SL relaxation T_1ρ. (c) Experimental (red) and simulated (black) R-CHEESY spectra for different t_p values with their real and imaginary components. (d) Experimental (red) and simulated (black) Rabi nutation obtained from the R-CHEESY excitation profiles for varied pulse lengths t_p as the real intensity at $\mathit{\nu }-{\mathit{\nu }}_{\mathrm{off}}=\mathrm{0}$. FT of the experimental data yielded ν₁=0.63 MHz. (e) Intensity of the dispersive steady-state signal from (c) and comparison with the T_1ρ trace measured with ν₁≈0.8 MHz. For experimental parameters, see Sect. B5.

This is further validated by the decay of the experimental steady-state signal along y^′ for t_p≳10 µs. Overlay of the maximum y^′ signal intensity as a function of t_p with a T_1ρ decay measured at a similar spin lock power showed good agreement (see Fig. 7b and e). Thus, this decay could be attributed to T_1ρ relaxation. As no longitudinal relaxation was included in our theoretical framework, this effect was not reproduced in the simulations. Simulating this effect is complicated by the fact that standard longitudinal relaxation ultimately approaches a non-zero equilibrium distribution. How to account for this in the case of T_1ρ is disputed in the literature and goes beyond the scope of this work (e.g., Slichter, 1990, p. 243; Levitt, 2013, p. 316). Development of a simulation procedure and further benchmark experiments should instead be the objective of future work.

Figure 8 (a, b) Experimental steady-state z^′ and y^′ profiles (t_p=10 µs) and analytical solution from Eq. 17. Experimental MW powers (a: ν₁=0.77 MHz, b: ν₁=0.63 MHz) were used for calculating the analytical solution. (c) Evolution of M(t) calculated analytically from the tilted-frame solution of the Bloch equations (Eq. 16, black) and comparison with the numerical solution for the standard Bloch equation in the rotating frame (Eq. A19, light blue). Effective MW field vector ν_eff is shown as a gray vector, and its components ν₁=0.77 and $\mathit{\nu }-{\mathit{\nu }}_{\mathrm{off}}=\mathrm{1.7}\mathrm{MHz}$ are shown as dotted gray lines along y^′ and z^′; the y^′ and z^′ components of M(t_p) are shown as black and light blue dotted lines after t_p=10 µs for the tilted- and rotating-frame solutions, respectively. (d) FWHM of the calculated inversion profiles using both the tilted- (black) and the rotating-frame (light blue) expression as a function of t_p. The horizontal and vertical lines mark the values $\mathrm{FWMH}=\mathrm{2}{\mathit{\nu }}_{\mathrm{1}}=\mathrm{1.54}\mathrm{MHz}$ and $t_{\mathrm{p}}=T_{\mathrm{2}\mathit{\rho }}=\mathrm{4}\mathrm{\mu }\mathrm{s}$, respectively. Oscillations for small t_p are caused by the oscillating side lobes of the initial sinc-shaped profiles. For additional parameters, see Sect. B6.

4.2.3 Comparison with Bloch equations

As we assume an isolated $S=\mathrm{1}/\mathrm{2}$ electron spin, in principle, one can represent the spin system by a magnetization vector in three-dimensional Cartesian space and calculate its evolution from the Bloch equations (Bloch, 1946). Being defined in the rotating frame, this is analogous to the density matrix calculations without transformation into the tilted frame. Thus, it leads to errors when implementing relaxation in the case of long MW pulses (see below). Redfield (1955) proposed a modification of the Bloch equations in the limit of strong irradiation in which the $x^{\prime }{y}^{\prime }$ relaxation was divided into a fast process perpendicular to the irradiation field and a slower one for the parallel component. Indeed, this modification resulted in line narrowing compared to the steady-state solution of the standard Bloch equations under strong driving conditions (Schenzle et al., 1984; DeVoe and Brewer, 1983). Although this model is similar to our tilted-frame approach, it is different in that there are now two axes with slow relaxation and that it does not directly depend on the spin offset. In contrast to this, we analyzed the Bloch propagation in the tilted frame, in analogy to our density matrix calculations. Solution of the Bloch equations transformed in the tilted frame yielded the following rotating-frame magnetization vector (for the derivation, see Sect. A15):

$$\begin{array}{}\text{(16)}& {\mathit{M}}^{\prime }\left(t\right)=\left(\begin{array}{l}\sin \left({\mathit{\omega }}_{\mathrm{eff}}t\right)\sin \left(\mathit{\theta }\right)M_{\mathrm{0}}\cdot e^{-t/T_{\mathrm{2}\mathit{\rho }}}\\ \left[-\cos \left({\mathit{\omega }}_{\mathrm{eff}}t\right)\sin \left(\mathit{\theta }\right)\cos \left(\mathit{\theta }\right)\cdot e^{-t/T_{\mathrm{2}\mathit{\rho }}}\right.\\ \left.+\sin \left(\mathit{\theta }\right)\cos \left(\mathit{\theta }\right)\right]M_{\mathrm{0}}\\ \left[\cos \left({\mathit{\omega }}_{\mathrm{eff}}t\right){\sin }^{\mathrm{2}}\left(\mathit{\theta }\right)\cdot e^{-t/T_{\mathrm{2}\mathit{\rho }}}+{\cos }^{\mathrm{2}}\left(\mathit{\theta }\right)\right]M_{\mathrm{0}}\end{array}\right).\end{array}$$

For the limiting case of infinite microwave irradiation (t→∞), all terms including an exponential decay vanished, and an analytical steady-state solution was obtained:

$$\begin{array}{}\text{(17)}& \underset{t⟶\mathrm{\infty }}{lim}{\mathit{M}}^{\prime }\left(t\right)=\left(\begin{array}{l}\mathrm{0}\\ \sin \left(\mathit{\theta }\right)\cos \left(\mathit{\theta }\right)M_{\mathrm{0}}\\ {\cos }^{\mathrm{2}}\left(\mathit{\theta }\right)M_{\mathrm{0}}\end{array}\right)=\left(\begin{array}{l}\mathrm{0}\\ \frac{\left(\mathit{\nu }-{\mathit{\nu }}_{\mathrm{off}}\right){\mathit{\nu }}_{\mathrm{1}}}{{\mathit{\nu }}_{\mathrm{1}}^{\mathrm{2}}+{\left(\mathit{\nu }-{\mathit{\nu }}_{\mathrm{off}}\right)}^{\mathrm{2}}}\\ \mathrm{1}-\frac{{\mathit{\nu }}_{\mathrm{1}}^{\mathrm{2}}}{{\mathit{\nu }}_{\mathrm{1}}^{\mathrm{2}}+{\left(\mathit{\nu }-{\mathit{\nu }}_{\mathrm{off}}\right)}^{\mathrm{2}}}.\end{array}\right)\end{array}$$

This steady-state expression was already proposed by Abragam (1961) and is in agreement with the recent literature on EDNMR (Nalepa et al., 2014; Cox et al., 2017). Comparison with the experimental steady-state y^′ and z^′ profiles showed quantitative agreement (Fig. 8a and b; compare to Fig. 5c). Additionally, in contrast to the approaches by Abragam (1961) and Nalepa et al. (2014), the complete spin trajectory during a MW pulse including relaxation could be calculated analytically using Eq. (16). The results for a spin offset ${\mathrm{\Omega }}_S/\mathrm{2}\mathit{\pi }=\mathrm{1.7}\mathrm{MHz}$ and ν₁=0.77 MHz are shown as a black curve in Fig. 8c for a pulse length t_p=10 µs. Starting from ${\mathit{M}}^{\prime }=(\mathrm{0},\mathrm{0},M_{\mathrm{0}}{)}^{\mathrm{T}}$, the magnetization vector precesses around ν_eff (gray arrow) and approaches a steady state (y^′ and z^′ components marked by dashed black line) that lies on ν_eff. The whole evolution occurs in a plane perpendicular to ν_eff that cuts the z^′ axis at M₀. Again, we compared our results to propagation in the rotating frame by numerical solution of the Bloch equations in the rotating frame (for details, see Sect. A15). Analysis of this trajectory (Fig. 8c, blue) showed that there is a decay of the vector components parallel to the effective field in addition to the precession around ν_eff. This effect is already clearly visible for t_p=10 µs and culminates in $\mathit{M}(t\to \mathrm{\infty })=\mathbf{0}$ for even longer pulses (not shown). In agreement with the density matrix simulations in the rotating frame, this leads to an increased hole depth for a certain offset Ω_S and, thus, to a broadening of the burned hole. Quantitative analysis of the width of the calculated hole profiles by their full width at half maximum (FWHM) value (Fig. 8d) showed that, while the tilted-frame solution approached a steady state at the theoretical value of 2ν₁ (horizontal black line), the width of the hole calculated in the rotating frame continuously increased. This effect becomes significant for t_p≳T_2ρ (vertical black line). The experimental FWHM (from the data in Fig. 5) shows only a slight increase in the FWHM for long t_p that is not consistent with the rotating-frame results. We attribute the comparably small increase in the FWHM to spectral diffusion as observed, for example, by Hovav et al. (2015).

5 Conclusions

In this work, we demonstrated the feasibility of shaped-pulse experiments using periodic phase modulations and chirp pulses on a commercial Bruker E680 W-band spectrometer. These experimental tools were then applied to the analysis of electron spin dynamics during spin locking. In particular, we reported that spin dynamics during MW irradiation obey principles similar to the bare state, i.e., Larmor, precession around effective fields and transverse relaxation in the eigenframe of the complete spin Hamiltonian. We performed tilted-frame experiments to measure SL relaxation times and to analyze the relation between bare and SL dynamics. Inversion and excitation profiles of low-power MW pulses were used to benchmark density matrix simulations of these SL dynamics in a reduced spin system. Quantitative agreement was reached using experimental parameters without modification. Additional Spinach-based simulations, including electron–nuclear coupling, powder averaging, and explicit detection, suggested that our findings can be generalized to more complex spin systems and simulation routines. Comparison of these results with simulations based on the Bloch equations showed that, in both cases, transformation into the tilted eigenframe of the MW Hamiltonian is needed for pulse lengths exceeding the SL decoherence time T_2ρ to achieve quantitative agreement with the experimental data.

While we performed the experiments with a BDPA radical which has a narrow and unstructured EPR spectrum, we expect some of the results to be transferable to other samples, including the frequently used nitroxide radicals. First chirp EPR experiments of nitroxide radicals are a promising starting point for using FT EPR for broadband experiments on our W-band spectrometer, and the work in this direction is still ongoing.

Multiple further questions still need to be addressed. Knowledge of the dependence of T_2ρ on the chemical environment could increase mechanistic understanding of this relaxation process and, ultimately, enable tuning of the relaxation properties in the SL state to increase the decoherence time in pulsed EPR experiments. Furthermore, the effect of T_1ρ was not analyzed. A similar analysis of the longitudinal SL relaxation would be interesting if an appropriate benchmark experiment could be designed. Additionally, a detailed, quantitative analysis of the hyperfine decoupling in the SL state could be used to estimate the applicability of hyperfine spectroscopic experiments in the tilted frame.

Appendix A: Supplementary data and analysis

A1 Generation of shaped pulses

Figure A1 Schematic description of the pulse generation in a Bruker E680 W-band EPR spectrometer equipped with an AWG. See Sect. 3.2 for more details.

A2 Phase stability of the spectrometer

To assess the phase stability of the pulses produced by the AWG, we performed a spin echo experiment with transient detection as described by Endeward et al. (2023). The results are displayed in Fig. A2. Analysis of the data shows that strong phase noise is observable, likely due to a deteriorating performance of the 84.5 GHz MW source. When looking at the average of a hundred shots (filled circles), the phase agrees reasonably well with the theoretical value (red “x” symbols), suggesting that phase cycling is possible whenever large numbers of scans are used. This is in agreement with our experiments where we observe efficient phase cycling. However, this procedure will lead to a loss in the signal-to-noise ratio due to partial cancellation.

A3 Resonator profile and amplifier characterization

Figure A3a shows the ν₁ profile of the Bruker EN600-1021H ENDOR resonator measured with BDPA at a temperature of 50 K. The profile is clearly asymmetric and shows some oscillations, which are indicative of reflections in the resonator (Doll and Jeschke, 2014; Endeward et al., 2023).

Additionally, Rabi nutation experiments were performed at different AWG output voltages to measure the amplification curve of the SSPA. This was achieved by variation of the amplitude parameter (Amp) in the Xepr software that scales the input IQ voltage values by a factor between 0 % and 100 %. Due to saturation effects, the amplification curve was clearly nonlinear. Therefore, all nominal MW powers were corrected with this curve in order to yield the target ν₁.

A4 Generation of PM pulses

PM pulses were generated with a home-written MATLAB routine following Eq. (5) (Wili et al., 2020). Figure A4a shows two PM pulses with phases of ϕ_PM1=0 and ${\mathit{\varphi }}_{\mathrm{PM}\mathrm{2}}=\mathit{\pi }/\mathrm{2}$, respectively. The I(t) and Q(t) values are calculated from ϕ_MW using Eq. (A1):

$$\begin{array}{}\text{(A1)}& I\left(t\right)=\cos \left({\mathit{\varphi }}_{\mathrm{MW}}\right),Q\left(t\right)=\sin \left({\mathit{\varphi }}_{\mathrm{MW}}\right).\end{array}$$

A representative pulse shape including two PM pulses is shown in Fig. A4b. Note that some distortions caused by the comparably large time increment of 10 ns are visible for ν_PM=15 MHz. However, choosing a smaller time increment was impossible due to the limited number of IQ values for a single pulse shape supported by the Xepr software.

A5 Calibration of ν_PM and measurement of the ν₁ profile

The frequency-stepped SL resonance spectrum was measured to set ν_PM=ν₁ (Wili et al., 2020). In this experiment, a single long PM pulse with low modulation amplitude (t_PM=1 µs and a_PM=0.04) was applied during SL (Fig. A5a). This caused a signal decrease proportional to the degree of matching of ν₁ and ν_PM. The frequency ν_PM with a maximal signal decrease was used in the PM pulse experiments (i.e., ν_set; see Fig. A5b, black).

To obtain a discrete probability distribution function of ν₁ (i.e., P(ν₁)) from this experiment, the spectrum was vertically shifted and normalized. The abscissa was scaled by ν_set to allow for transfer between experiments with different ν_set values. For simulation of MW inhomogeneity, all data points exceeding $\mathrm{0.3}P({\mathit{\nu }}_{\mathrm{1}}/{\mathit{\nu }}_{\mathrm{set}}=\mathrm{1})$ were selected, and the spectrum was again normalized to a sum of one (Fig. A5b, red). Calculations were repeated for each of these ${\mathit{\nu }}_{\mathrm{1}}/{\mathit{\nu }}_{\mathrm{set}}$ values and the results weighted by $P({\mathit{\nu }}_{\mathrm{1}}/{\mathit{\nu }}_{\mathrm{set}})$. For the simulation of SL echoes, linear interpolation was used between the selected data points to increase the resolution of the ${\mathit{\nu }}_{\mathrm{1}}/{\mathit{\nu }}_{\mathrm{set}}$ axis.

Figure A2 (a) Pulse sequence to measure the phase setting and stability using the AWG. The phase of the rectangular $\mathit{\pi }/\mathrm{2}$ pulse was varied, while the π pulse was kept with a constant phase. (b) Real and imaginary contributions of the echo transients obtained with the pulse sequence in (a). The small dots mark a hundred individual scans for each phase ${\mathit{\varphi }}_{\mathit{\pi }/\mathrm{2}}$ belonging to one color. The average of these results is marked by a larger filled circle of the respective color. The red “x” symbols mark the theoretical position. Experimental parameters: t_π=60 ns, and τ=1 µs.

Figure A3 (a) Resonator profile used for the correction of the chirp pulses. Acquisition was done at 50 K by measuring two-pulse Rabi nutations of BDPA at different local oscillator frequencies. (b) Amplification curve of the solid-state amplifier used for correction of the chirp pulses. The data were obtained by room temperature measurements of BDPA Rabi nutations at different attenuation values of the AWG output voltage. The latter are given as relative amplitude (Amp) values that can be between 0 % and 100 %. A second-order polynomial $f\left(x\right)=a+bx+c{x}^{\mathrm{2}}$ was used for fitting. Pulse sequence: $\mathrm{pulse}-T-\mathit{\pi }/\mathrm{2}-\mathit{\tau }-\mathit{\pi }-\mathit{\tau }-\mathrm{echo}$. Experimental parameters: t_π≈40–400 ns, T≈10 µs, and τ≈0.3–1 µs.

A6 Generation of chirp pulses

Shapes of chirp pulses were generated using a home-written MATLAB routine. After defining the pulse length t_p and sweep width ν_SW, an incremented frequency array of ν_AWG was generated that spanned the whole sweep range. For each increment i, the length Δt(i) was adapted using the experimental resonator profile in Fig. A3a so that the critical adiabaticity Q_crit(i) (Eq. A2) was constant throughout the pulse and equal to its mean value $\stackrel{\mathrm{‾}}{Q_{\mathrm{crit}}}$ (Baum et al., 1985; Jeschke et al., 2015).

$$\begin{array}{}\text{(A2)}& {\displaystyle }{Q}_{\mathrm{crit}}\left(i\right)={\displaystyle \frac{\mathrm{2}\mathit{\pi }\left[{\mathit{\nu }}_{\mathrm{1}}\right(i\left){]}^{\mathrm{2}}\mathrm{\Delta }{t}_{\mathrm{p}}\right(i)}{\mathrm{\Delta }\mathit{\nu }}}=\stackrel{\mathrm{‾}}{Q_{\mathrm{crit}}},\text{(A3)}& {\displaystyle }\stackrel{\mathrm{‾}}{Q_{\mathrm{crit}}}={\displaystyle \frac{\mathrm{2}\mathit{\pi }\stackrel{\mathrm{‾}}{{\mathit{\nu }}_{\mathrm{1}}^{\mathrm{2}}}{t}_{\mathrm{p}}}{{\mathit{\nu }}_{\mathrm{SW}}}}.\end{array}$$

Here, ν₁(i) is the Rabi frequency for increment i as obtained from the resonator profile. Δν corresponds to the constant frequency increment ${\mathit{\nu }}_{\mathrm{AWG}}(i+\mathrm{1})-{\mathit{\nu }}_{\mathrm{AWG}}\left(i\right)$. Δt_p is the length of time increment i, and $\stackrel{\mathrm{‾}}{{\mathit{\nu }}_{\mathrm{1}}^{\mathrm{2}}}$ is the mean value of ${\mathit{\nu }}_{\mathrm{1}}^{\mathrm{2}}$ as obtained from the resonator profile.

Figure A4 (a) MW phase ϕ_MW during a representative pulse shape including PM pulses (red) and corresponding carrier wave (after Wili et al., 2020). The following parameters were used: ${\mathit{\nu }}_{\mathrm{1}}={\mathit{\nu }}_{\mathrm{PM}}=\mathrm{15}\mathrm{MHz}$, ${\mathit{\varphi }}_{\mathrm{0}}=\mathit{\pi }/\mathrm{2}$, a_PM=0.4, τ₀=100 ns, ${\mathit{\tau }}_{\mathrm{1}}={\mathit{\tau }}_{\mathrm{2}}=\mathrm{200}\mathrm{ns}$, t_PM1=80 ns, t_PM2=160 ns, ϕ_PM1=0, and ${\mathit{\varphi }}_{\mathrm{PM}\mathrm{2}}=\mathit{\pi }/\mathrm{2}$. A time increment of 1 ns was used for better visualization. Note the shift of PM2 with respect to the carrier wave due to its PM phase. (b) I(t) and Q(t) values of a representative pulse shape including PM pulses as uploaded to the spectrometer. The following parameters were used: ${\mathit{\nu }}_{\mathrm{1}}={\mathit{\nu }}_{\mathrm{PM}}=\mathrm{15}\mathrm{MHz}$, ${\mathit{\varphi }}_{\mathrm{0}}=\mathit{\pi }/\mathrm{2}$, a_PM=0.4, ${\mathit{\tau }}_{\mathrm{0}}={\mathit{\tau }}_{\mathrm{2}}=\mathrm{500}\mathrm{ns}$, τ₁=1000 ns, t_PM1=80 ns, t_PM2=160 ns, ϕ_PM1=0, and ${\mathit{\varphi }}_{\mathrm{PM}\mathrm{2}}=\mathit{\pi }/\mathrm{2}$. A time increment of 10 ns was used in agreement with the experiments.

Figure A5 (a) Pulse sequence for measuring the frequency-stepped SL resonance spectrum and result for ν_set=12 MHz (b, black). The processed data used for simulations are shown in red. Experimental parameters: t_π=48 ns, τ₀=1 µs, τ₁=2 µs, t_PM=1 µs, a_PM=0.04, τ=0.3 µs, SPP=5, and four-step PC of the π detection pulse.

The nonlinear frequency sweep was interpolated to fit the constant AWG sampling rate of Δt_AWG=0.625 ns and converted to a phase modulation by Eq. (A4) (Doll and Jeschke, 2017).

$$\begin{array}{}\text{(A4)}& {\mathit{\varphi }}_{\mathrm{AWG}}\left(t\right)=\underset{\mathrm{0}}{\overset{t}{\int }}\mathrm{2}\mathit{\pi }{\mathit{\nu }}_{\mathrm{AWG}}\left(t^{\prime }\right)\mathrm{d}{t}^{\prime }+{\mathit{\varphi }}_{\mathrm{AWG}}\left(\mathrm{0}\right)\end{array}$$

This phase was inserted into Eq. (A1), yielding the AWG input I(t) and Q(t). A wideband, uniform rate, smooth truncation (WURST) amplitude modulation function was applied to these I(t) and Q(t) values in order to reduce pulse edge artifacts (Kupce and Freeman, 1995):

$$\begin{array}{}\text{(A5)}& \begin{array}{rl}{I}_{\mathrm{WURST}}\left(t\right)=& I\left(t\right)\cdot \left(\mathrm{1}-{\left|\cos \left({\displaystyle \frac{\mathit{\pi }t}{t_{\mathrm{p}}}}\right)\right|}^{n_{\mathrm{W}}}\right)\\ & \mathrm{where}t\in [\mathrm{0},t_{\mathrm{p}}],\end{array}\end{array}$$

where n_W is a parameter that reflects the steepness of the WURST modulation function. For the experiments shown in this work, n_W=50 was used. The final pulse shapes were obtained after correction with the amplifier profile (see Fig. A3b). In order to assess the performance of these pulses, we measured inversion profiles of representative pulse shapes. We used a spin echo experiment as illustrated in Fig. A6a, b and compared the echo intensity with and without the pulse of interest (Doll and Jeschke, 2014). In contrast to the work by Bahrenberg et al. (2017), we varied the frequency of the detection sequence and the B₀ field, while keeping the analyzed pulse at a constant frequency centered at the resonator dip with ν_resonator=ν_LO. In this way, the shaped pulses stayed fixed with respect to the resonator profile, and we could monitor the effect of the resonator bandwidth on the inversion properties of the pulses. In the case where the detection frequency is fixed at the position of the resonator dip, we would expect overestimation of the pulse performance because inversion is always measured at the position with maximal ν₁.

With this procedure, two factors had to be taken into account. First, ν₁ of the detection pulses depends on their position with respect to the resonator profile. In order to keep the pulse length constant, the AWG output voltage (set by the amplitude parameter in the Xepr software) of the detection pulses was scaled with the inverse resonator profile (Doll and Jeschke, 2014). Second, signal detection by the spectrometer is done in a quadrature detection scheme by demodulation with ν_PLO and ν_x. Thus, the resulting echo was oscillating with ν_AWG that was different for each data point. We demodulated the echoes using a complex exponential function oscillating with ν_AWG (see Fig. A6b). ϕ₀ was optimized for the background spectrum without the analyzed pulse and then used for the echo with pulse without further modification. Figure A6c shows the inversion profiles of representative chirp pulses with a pulse length of t_p=1 µs in a frequency range of ±80 MHz. While the π pulse with ν_SW=100 MHz shows good inversion efficiency over the whole sweep range, some artifacts are visible in the case of ν_SW=200 MHz. The pulse with ν₁ optimized for a $\mathit{\pi }/\mathrm{2}$ pulse shows a reasonably flat profile with $I/I_{\mathrm{0}}\approx \mathrm{0}$ as expected.

Figure A6 (a) Experiment used for acquisition of inversion profiles together with the used pulse sequence (b), following the procedure of Doll and Jeschke (2014). (c) Representative inversion profiles of 1 µs WURST chirp pulses corrected for the resonator profile. Inversion (inv.) pulses were performed with maximal ν₁, while the ν₁ value of the $\mathit{\pi }/\mathrm{2}$ pulse was optimized with a chirp nutation (see Fig. A7). The frequency value in the legend gives the sweep range the pulses were optimized for. Experimental parameters: t_π=72 ns, τ=0.6 µs, and SPP=5.

A7 Characterization of chirp echoes

Figure A7a shows the chirp nutations to determine the ideal ν₁ value for both the $\mathit{\pi }/\mathrm{2}$ and the π chirp pulse (Doll and Jeschke, 2014). Despite the low Q_crit of the chirp π pulse (see Sect. 3.5), the FT of a chirp echo of BDPA showed quantitative agreement with the corresponding ED EPR spectrum, as it is shown in Fig. A7b. For the latter, the field axis B was converted into a frequency axis ν using $\mathit{\nu }={\mathit{\nu }}_{\mathrm{off},\mathrm{FS}}-g_{\mathrm{e}}{\mathit{\mu }}_{\mathrm{B}}hB$; ν_off,FS was selected manually to align the ED EPR with the FT spectrum. Figure A7c shows the phased chirp echo FT spectrum of BDPA as a function of the offset relative to the center of the resonator profile (Endeward et al., 2023). The corresponding integral and the phase ϕ₀ of the signal are shown in (d). The intensity profile closely follows the resonator profile. Fig. A8 shows the analogous comparison for a TEMPOL radical. The spectral width exceeds the sweep range of the chirp pulses by more than 200 MHz, which prevented the detection of the full EPR spectrum by FT spectroscopy. Within the sweep range, the FT spectrum is the approximate convolution of the echo detected EPR and the resonator profile. Additionally, there are some further distortions in the line shape, presumably caused by the pulse edges or reflections in the MW resonator. Although this result shows the limitations of our experimental setup to completely detect broad EPR spectra, it suggests that FT EPR might be useful for other applications to partially benefit from the multiplex advantage of the FT approach.

Figure A7 (a) Representative chirp echo nutations for the $\mathit{\pi }/\mathrm{2}$ and π chirp pulses using a 1 µs and 0.5 µs chirp echo with a sweep width of 200 MHz. The $\mathit{\pi }/\mathrm{2}$ nutation was performed with ${\mathit{\nu }}_{\mathrm{1},\mathit{\pi }}=\mathrm{10}\mathrm{MHz}$, while the π nutation was acquired using the optimized value of ${\mathit{\nu }}_{\mathrm{1},\mathit{\pi }/\mathrm{2}}=\mathrm{7.7}\mathrm{MHz}$. (b) Comparison of the BDPA EPR spectrum obtained by FT of a chirp echo (red) and by echo-detected EPR (light blue). The asterisks denotes a zero-frequency artifact in the FT. (c) FT EPR spectra of BDPA as a function of the frequency offset is shown. Acquisition was done by variation of the magnetic field. (d) Integral of each FT spectrum is shown as a function of the frequency offset (red) as well as its corresponding phase (black). The dashed lines mark the sweep range of the chirp pulses. Note that the center of the integral profile in (d) is shifted to negative ν_off values by approximately 10 MHz, which was presumably caused by tuning imperfections. Pulse sequence: $\mathit{\pi }/{\mathrm{2}}_{\mathrm{chirp}}-\mathit{\tau }-{\mathit{\pi }}_{\mathrm{chirp}}-\mathit{\tau }-\mathrm{echo}$. Experimental parameters: $t_{\mathit{\pi }/\mathrm{2}}=\mathrm{1}$, t_π=0.5, $\mathit{\tau }\approx \mathrm{1}-\mathrm{2}\mathrm{\mu }\mathrm{s}$, SPP≈20, and eight-step PC.

Figure A8 Comparison of the TEMPOL EPR spectrum obtained by FT of a chirp echo (red) and by echo-detected EPR (light blue). For the latter, the field axis was converted into the corresponding frequency axis and horizontally shifted to overlap with the FT spectrum. FT pulse sequence: $\mathit{\pi }/{\mathrm{2}}_{\mathrm{chirp}}-\mathit{\tau }-{\mathit{\pi }}_{\mathrm{chirp}}-\mathit{\tau }-\mathrm{echo}$. Experimental parameters: $t_{\mathit{\pi }/\mathrm{2}}=\mathrm{1}$, t_π=0.5, τ=2 µs, SPP=2000, and eight-step PC.

A8 SL Rabi nutations

The signals observed in the FT of Fig. 2, centered at the SL Rabi frequency $\frac{a_{\mathrm{PM}}{\mathit{\nu }}_{\mathrm{1}}}{\mathrm{2}}$ and at 2ν₁, were also retained when varying ν₁=ν_PM (see Fig. A9). Again, simulations agreed well with the experimental data and the theoretical signals.

Figure A9 FT spectra of the experimental SL Rabi nutation traces (red) obtained with a fixed a_PM=0.4 and different ν₁ values. (a) Low-frequency region on the order of the SL Rabi frequency $\frac{a_{\mathrm{PM}}{\mathit{\nu }}_{\mathrm{1}}}{\mathrm{2}}$. (b) Complete frequency region. FTs of simulated SL Rabi nutations with the same MW power are displayed in black. Vertical lines mark the positions where a signal is expected, i.e., at $\mathit{\nu }=\frac{a_{\mathrm{PM}}{\mathit{\nu }}_{\mathrm{1}}}{\mathrm{2}}$ in (a) and at ν=2ν₁ in (b). For experimental parameters, see Sect. B1.

If a single, bare-spin offset Ω_S was introduced to the simulation, a strong distortion of the FT spectra occurred for Ω₁≳ω₁ (Fig. A10a). In contrast, simulation of a Gaussian distribution of Ω_S values with FWHM=24 MHz did not yield these artifacts (Fig. A10b).

Figure A10 FTs of SL Rabi nutations (red, as in Fig. A9) compared to simulations with a single electron spin offset ${\mathrm{\Omega }}_S/\mathrm{2}\mathit{\pi }=\mathrm{10}\mathrm{MHz}$ (a) and with a Gaussian distribution of Ω₁ values (b). For the latter, the weighted sum of 2001 traces calculated for ${\mathrm{\Omega }}_S/\mathrm{2}\mathit{\pi }$ between −24 and 24 MHz was calculated using weighting factors from a Gaussian distribution function centered at 0 MHz and with FWHM=24 MHz.

A9 Phase of SL echo experiments

During a SL echo experiment, we can distinguish between two kinds of phase information. On the one hand, the final integrated echo has a phase depending on the magnetization phase at the end of the SL pulse. On the other hand, the SL echo itself (i.e., the one generated by the PM pulses) has a phase relative to the readout PM3 pulse. As both of these phases carry additional experimental information, these will be discussed in the following section. The imaginary signal component of the SL echo experiments in Fig. 3, which refers to the phase of the integrated echo, is shown in Fig. A11. Comparison of real and imaginary component shows that the signals are clearly separated in the τ₂ dimension. While the real signals can be assigned to spins that are refocused at the end of τ₂, the imaginary signals all stem from spins refocused at the end of τ₃. Namely, (i), (ii), and (iii) can all be explained by stimulated echoes resulting from different combinations of delays. Signal (iv) is a combination of a refocused echo followed by a stimulated echo. Please note that only signals where PM3 acts as a $\mathit{\pi }/\mathrm{2}$ pulse are visible because of the two-step PC of PM3 in all experiments. The signal in (a) that is marked with an asterisk is a contribution from the real part caused by phase drifts between the experiments.

Figure A11 (a, b) Imaginary signal (red) of the SL echo experiments shown in Fig. 3c and d for the SL echo and WOP SL echo, respectively. For comparison, the real parts of the signals (see Fig. 3) are shown in gray. The asterisk marks an artifact arising from phase drifts between the experiment. (c) Coherence diagrams explaining the signals (i) to (iv). All spectra were obtained with an additional two-step PC on PM3. Thus, only signals where PM3 acted as a $\frac{\mathit{\pi }}{\mathrm{2}}$ pulse are visible. For experimental parameters, see Sect. B2.

This phase separation of the SL echo signals differs from the phase of the SL echoes at the position of PM3. In order to measure this phase, quadrature detection was applied on the readout pulse PM3. Namely, its two-step phase cycle that was applied in all SL echo experiments was shifted by 90° so that, instead of ${\mathit{\varphi }}_{\mathrm{PM}\mathrm{3}}=[\mathrm{0}\mathit{°},\mathrm{180}\mathit{°}]$, ${\mathit{\varphi }}_{\mathrm{PM}\mathrm{3}}=[\mathrm{90}\mathit{°},\mathrm{270}\mathit{°}]$ was used. The results for the different phases are displayed in Fig. A12. Although the experiment was performed under nominally the same conditions as the one shown in Fig. 3d, the relative phase of the signals with ${\mathit{\varphi }}_{\mathrm{PM}\mathrm{3}}=[\mathrm{0}\mathit{°},\mathrm{180}\mathit{°}]$ is different in this spectrum. As the first peak is caused by a refocused echo and the second by a stimulated echo, an inverted sign would be expected for the two, which is the case in Fig. A12 but not in Fig. 3d. This already shows that the phase of the echoes is not stable in our experiment. Additionally, the imaginary signal of the quadrature detection experiment (i.e., ${\mathit{\varphi }}_{\mathrm{PM}\mathrm{3}}=[\mathrm{90}\mathit{°},\mathrm{270}\mathit{°}]$) shows significant absorptive contributions for some of the signals. This again demonstrates that the observed phase in this experiment deviates from the theoretically expected one.

Figure A12 WOP SL echo measured with standard (dark red) and 90° phase-shifted (blue) PM3 pulse. The former is obtained under identical conditions as the bottom trace in Fig. 3d. For experimental parameters, see Sect. B2.

A10 T_2ρ measurements

The SL relaxation decays and their mono-exponential and stretched-exponential least-square fits (Eq. A6) corresponding to the T_2ρ values in Fig. 4 are shown in Fig. A13a and b. Equations (A6) and (A7) were used for the fitting:

$$\begin{array}{}\text{(A6)}& {\displaystyle }I=I_{\mathrm{0}}\cdot \exp \left(-{\displaystyle \frac{{\mathit{\tau }}_{\mathrm{1}}+{\mathit{\tau }}_{\mathrm{2}}}{T_{\mathrm{2}\mathit{\rho }}}}\right),\text{(A7)}& {\displaystyle }I=I_{\mathrm{0}}\cdot \exp \left(-{\left({\displaystyle \frac{{\mathit{\tau }}_{\mathrm{1}}+{\mathit{\tau }}_{\mathrm{2}}}{T_{\mathrm{2}\mathit{\rho }}}}\right)}^{\mathit{\xi }}\right).\end{array}$$

Oscillations are visible in the beginning of some of the experimental traces. Following the literature concerned with the fitting of T_m traces that include electron spin echo envelope modulation (ESEEM) (Soetbeer et al., 2021), we fitted the envelope function of the oscillation maxima. Specifically, we generated a time window with a length of the inverse of the main oscillation frequency for each data point. The data point with maximal intensity was selected for each of these windows. After removing the redundancies in this selected data set, only the data points at a maximum of the dominant oscillation were obtained. A detailed analysis of the oscillation frequencies can be found in Sect. A12. Statistical analysis was performed using the bootstrapping method (e.g., Press et al., 1992). For each N point data set, N_boot=1000 synthetic data sets were generated by randomly taking N points from the data set with replacement and thus, allowing points to be selected multiple times. Each of these bootstrap data sets was fitted and the final result and its error were obtained as the mean and standard deviation of the fit parameters from all synthetic data sets. Errors are given as the 95 % confidence intervals (i.e., as 2σ). Note that this only takes into account statistical errors and cannot give an estimate of systematic errors introduced by the fit model. These can instead be analyzed using the residuals of the fitting (Fig. A13c). The analogous data for the Hahn echo decay and CPMG experiments can be found in Fig. A13d–f.

Figure A13 Experimental SL echo decay traces (gray) and the corresponding mono-exponential (a) and stretched-exponential (b) fits. Data points selected for the envelope fitting are displayed as black dots. The corresponding T_2ρ values can be found in Fig. 4. (c) Residuals of both fits, scaled by a factor of 3 for better visualization. (d–f) Analogous data for the Hahn echo and the Carr–Purcell decay measurements, where N=0 for Hahn echo decays and $N=\mathrm{1},\mathrm{2}$ for CP3 and CP4, respectively. The corresponding T_m values can be found in Figs. 4 and A15. For experimental parameters, see Sect. B3.

Note that the oscillations in the residuals for small τ₁+τ₂ indicate that the removal of oscillating components in the decay traces is incomplete (see Sect. A12 for more detail). Analysis of the Pearson correlation coefficient for T_2ρ and ξ in the bootstrap samples was performed using the corrcoef command in MATLAB (see Fig. A14b). For all values of ν₁, the Pearson correlation of the fitting parameters ρ(ξ,T_2ρ) was above 0.9, which shows that the parameters are highly correlated.

Figure A14 (a) Stretching exponents ξ for the T_2ρ fits (dark red) and for T_m (dark blue) in Fig. A13b and e. (b) Pearson correlation coefficient for ξ and T_2ρ of the different fits in Fig. A13b and e calculated from the fitting results of the bootstrap samples. For experimental parameters, see Sect. B3.

A11 Carr–Purcell sequence

We measured decoherence times using the Carr–Purcell sequence with three and four pulses (CP3 and CP4). The results are shown in Fig. A15, and the decay traces can be found in Fig. A13d–f. For both experiments, the decays were well reproduced both by a mono-exponential fit and a stretched-exponential fit, and an increase in decoherence time was observed compared to the two-pulse Hahn echo experiment.

Figure A15 (a) Pulse sequence for the Carr–Purcell (CP) experiment with N repetitions of the central refocusing block. (b) Decoherence times were obtained from the fitting in Fig. A13d–f for the Hahn echo decay (N=0) and the CP3 and CP4 experiments ($N=\mathrm{1},\mathrm{2}$, respectively). For $N=\mathrm{1},\mathrm{2}$, a 16-step PC of the last two π pulses was applied. All experiments were performed with ν₁=15 MHz. For further experimental parameters, see Sect. B3.

A12 Oscillations in T_2ρ traces

Distortions of the SL state relaxation traces, consisting of oscillations for small τ₁+τ₂ values, were observed in most tilted-frame echo decay experiments. As deuteration of the used BDPA radical showed no effect on the oscillation frequency (data not shown), ESEEM-like modulation caused by residual hyperfine couplings could be excluded. Moreover, a phase shift of the oscillation by π was observed when changing ϕ_PM by $\frac{\mathit{\pi }}{\mathrm{2}}$ for all three PM pulses simultaneously (see Fig. A16a), which suggested an effect depending on an oscillation with $\left|\cos \right({\mathit{\varphi }}_{\mathrm{PM}}\left)\right|$. FT of the experimental T_2ρ decays for different ν₁ and τ₁+τ₂ increments (red in Fig. A16c and d) showed seemingly random but well-defined oscillation frequencies.

Considering the principle of PM pulses (Fig. A16b; Wili et al., 2020), it is apparent that, depending on the time t between the beginning of the spin lock and the PM pulse and on the PM phase ϕ_PM, there is a jump in ϕ_MW at the beginning and end of the PM pulses (blue and black “x” symbols in Fig. A16b). Its magnitude depends on $\left|\cos \right({\mathit{\omega }}_{\mathrm{PM}}t+{\mathit{\varphi }}_{\mathrm{PM}}\left)\right|$. As such abrupt changes in pulses are known to produce artifacts, these were a plausible explanation for the observed deviations from an ideal echo decay. Because the 80 ns increment of τ₁+τ₂ in our experiments was larger than the Nyquist threshold of half of the inverse PM frequency $\frac{\mathrm{1}}{\mathrm{2}{\mathit{\nu }}_{\mathrm{PM}}}\approx \mathrm{30}\mathrm{ns}$, not the ν_PM oscillation but a low-frequency alias would be detected in this case (see Fig. A16e). Thus, we simulated the oscillations by selecting values from $\left|\cos \right({\mathit{\omega }}_{\mathrm{1}}t\left)\right|$ with the experimental time increment and subsequently Fourier transformed this undersampled oscillation (black in Fig. A16c and d). The frequencies obtained from this FT agreed quantitatively with some of the signals observed in the experimental FT spectra of the decay traces. Additional signals could be reproduced when repeating the procedure with a $\left|\cos \left(\frac{{\mathit{\omega }}_{\mathrm{1}}}{\mathrm{2}}t\right)\right|$ function (blue in Fig. A16c and d). This is reasonable because, while the position of PM3 shifted with τ₁+τ₂, the position of PM2 varied with $\frac{{\mathit{\tau }}_{\mathrm{1}}+{\mathit{\tau }}_{\mathrm{2}}}{\mathrm{2}}$. Hence, a combination of the pulse edges caused by both PM2 and PM3 could explain most of the observed signals. Importantly, this method also predicts that there is no oscillation in the trace for ν₁=12.5 MHz, which agrees well with the experiment. In the future, we could detect the decay traces with a time increment τ₁+τ₂ adjusted so that no oscillations are present in the traces. This might improve the results of the fitting for various ν₁.

Figure A16 (a) Excerpt of T_2ρ nutation traces obtained with ν₁=15 MHz and ${\mathit{\varphi }}_{\mathrm{PM}}=\mathrm{0},\frac{\mathit{\pi }}{\mathrm{2}}$. (b) Schematic drawing of arbitrary PM pulses for ${\mathit{\varphi }}_{\mathrm{PM}}=\mathrm{0},\frac{\mathit{\pi }}{\mathrm{2}}$ (black, blue). The first points of the PM pulses are marked with an “x ” symbol and a continuous cosine wave as a reference for ϕ_PM=0 is marked in gray. FTs of the experimental T_2ρ traces for different MW powers ν₁ (c) and increments of $t={\mathit{\tau }}_{\mathrm{1}}+{\mathit{\tau }}_{\mathrm{2}}$ (d) show distinct frequencies (red). Assuming undersampling of oscillations $\left|\cos \right({\mathit{\omega }}_{\mathrm{PM}}t\left)\right|$ and $\left|\cos \left(\frac{{\mathit{\omega }}_{\mathrm{1}}}{\mathrm{2}}t\right)\right|$, the FT of the resulting traces yielded alias frequencies in quantitative agreement with the oscillations observed in the experiment (blue and black in c and d). (e) Illustration of alias frequencies for a $\left|\cos \left(\frac{{\mathit{\omega }}_{\mathrm{1}}}{\mathrm{2}}t\right)\right|$ oscillation with ν_PM=15 MHz and different sampling increments as used in (d). For experimental parameters, see Sect. B3.

A13 Simulation of MW inhomogeneity

Figure A17 shows the experimental Rabi nutation from the CHEESY z^′ profiles (see Fig. 5d) in red. Additionally, the simulations with and without MW inhomogeneity (black and blue, respectively) are shown. Comparison with an experimental T_2ρ decay (gray) shows that experiment and simulation including the ν₁ distribution agree well and decay significantly faster than the T_2ρ trace. In contrast, the simulation without MW inhomogeneity decays slower and closely resembles the T_2ρ decay. The ν₁ distribution for calculation of the black curve is taken from the frequency-stepped tilted-frame resonance spectrum without further modification.

Figure A17 Experimental Rabi nutation obtained from the CHEESY z profiles (red) and comparison with the simulation with (black) and without (blue) including MW inhomogeneity (ih.) using the ν₁ distribution function P(ν₁). For both simulations, T_2ρ was set to 4 µs. For comparison, the experimental T_2ρ decay (ν_PM=15 MHz) is shown in gray. For experimental parameters, see Sect. B4.

A14 Spinach simulation of the CHEESY z profiles

We chose an electron spin with a rhombic g tensor (see Fig. A20b) as a model for our Spinach simulations which differs from the almost isotropic g tensor of the BDPA radical. Setting the spin parameters to the literature values known for BDPA could not reproduce the experimental Gaussian line shape, because this shape is caused by the coupling to a large number of nearly equivalent protons in the BDPA molecule. Calculating the line shape accurately would require including all these nuclei into the simulations, which is computationally unfeasible. Hence, an artificial g anisotropy was introduced to generate a broad enough EPR line to monitor the whole inversion profiles. However, this caused a line shape different from a Gaussian EPR line. Namely, without any coupled nuclei, this leads to a distinct rhombic powder pattern with a singularity at $\mathit{\nu }-{\mathit{\nu }}_{\mathrm{off}}=\mathrm{0}$ (black curve in Fig. A20a). This caused small artifacts in the CHEESY simulations, as the CHEESY intensity scales with the EPR intensity, leading to profiles slightly too narrow in the center. Simulation of MW inhomogeneity was prevented by the increased computation time of the Spinach simulations. Although the shape of the thus simulated profiles agreed well with the experiments upon visual comparison, the Rabi nutation at $\mathit{\nu }-{\mathit{\nu }}_{\mathrm{off}}=\mathrm{0}$ decayed more slowly in the simulation (see Fig. A18). As described in Sect. A13 the decay of the experimental Rabi nutation is a combination of decoherence and MW inhomogeneity, which cannot be reproduced when omitting the MW inhomogeneity.

Figure A18 Experimental Rabi nutation obtained from the CHEESY z profiles (red) and comparison with the Spinach simulation without MW inhomogeneity (black). The simulation results closely match the simulation results with the home-written density matrix approach without MW inhomogeneity (light blue curve in Fig. A17). For experimental parameters, see Sect. B4.

By repeating the simulations with different T_2ρ times, the sensitivity of the simulations with regard to the relaxation time was analyzed (Fig. A19a). Using the bare-state transverse relaxation time T_m=1.4 µs (light red) did not reproduce the experimental profiles as the decay of the sinc oscillations was overestimated. This is most clearly visible for t_p=1.1 and 1.4 µs where significant features of the experimental profiles are missing. In contrast, the visual difference between T_2ρ=4 and 10 µs is small. Therefore, deriving T_2ρ solely from the CHEESY profiles is not easily possible. Importantly, omitting T_2ρ completely could not reproduce the experimental data (light blue), because no steady state was reached (compare t_p=10, 50, and 100 µs). The analogous simulations with T_m in the rotating frame are shown in Fig. A19b. While the profiles without relaxation (light blue) agree quantitatively with the tilted-frame results, the profiles with relaxation show broadening for large t_p values, which is not observed in the experiment. This deviation increased with decreasing decoherence time. Again, this demonstrates that relaxation is the reason why calculations need to be performed in the tilted frame.

Figure A19 (a) Experimental CHEESY z^′ profiles (red) compared to results from Spinach simulations in the tilted frame with different values of T_2ρ. The spectra were shifted vertically for the sake of clarity. (b) Analogous simulations performed in the rotating frame using the holeburn.m function of Spinach together with the t1_t2.m model for relaxation.

Including nuclei led to distinct changes in the simulated EPR spectra (Fig. A20a). Namely, coupling to one proton nucleus with the known hyperfine coupling tensor of BDPA (A tensor taken from Weis et al., 2000; see Fig. A20b) led to a distinct splitting (dark red curve), whereas multiple protons resulted in a more homogeneous line shape. Again, this irregular EPR line shape slightly distorted the simulated CHEESY profiles. Apart from these differences, the CHEESY profiles calculated with a different number of coupled protons are identical (Fig. A20c). This suggests that nuclear coupling does not affect the inversion profiles and can be neglected in the following considerations.

For these expanded spin systems, the relaxation superoperator ${\widehat{\widehat{R}}}_{\mathrm{2}\mathit{\rho }}$ had to be adapted. Comparison of the performance of different relaxation matrices (data not shown) suggested that using a standard relaxation matrix for transversal relaxation where the diagonal elements R_ij,ij are $T_{\mathrm{2}\mathit{\rho }}^{-\mathrm{1}}$ if $|i\rangle ⟶|j\rangle$ is an electron spin transition is not sufficient to describe all relaxation pathways. Instead, $R_{ij,ij}=T_{\mathrm{2}\mathit{\rho }}^{-\mathrm{1}}$ was set for all elements i≠j with all other elements of ${\widehat{\widehat{R}}}_{\mathrm{2}\mathit{\rho }}$ being zero. This is equivalent to defining an equal decoherence rate for all types of SL coherence. We assume that this becomes necessary as states in the tilted frame are mixed and the nuclear spin partially decoupled. A more thorough investigation of this effect is beyond the scope of this work and requires more analysis in the future.

Figure A20 (a) Normalized EPR spectra of an electron spin with g anisotropy coupled to up to four proton nuclei simulated using the Spinach library (Hogben et al., 2011). The experimental FT spectrum of BDPA (see Fig. A7b) is displayed as a dashed gray line. (b) Parameters used for the simulations. Elements of the g tensor were set to arbitrary values to achieve sufficient spectral width. Hyperfine couplings are taken from Weis et al. (1999) for ²H-BDPA and corrected by the ratio of the gyromagnetic ratios of ¹H and ²H. Protons one to four were added successively to obtain the spectra with zero to four coupled nuclei. (c) CHEESY profiles for the spin systems with zero to four nuclei (same color code as in a) simulated with the Spinach library. For clarity, the spectra are vertically shifted.

A15 Frame transformation of the Bloch equations

In order to perform calculations with the Bloch equations in the vector analogue of the tilted frame, we start with the Bloch equations in the rotating frame without relaxation terms (Bloch, 1946):

$$\begin{array}{}\text{(A8)}& {\dot{\mathit{M}}}^{\prime }=\left(\begin{array}{lll}\mathrm{0}& -{\mathrm{\Omega }}_S& {\mathit{\omega }}_{\mathrm{1}}\\ {\mathrm{\Omega }}_S& \mathrm{0}& \mathrm{0}\\ -{\mathit{\omega }}_{\mathrm{1}}& \mathrm{0}& \mathrm{0}\end{array}\right){\mathit{M}}^{\prime }\end{array}$$

where ${\mathit{M}}^{\prime }=(M_{x^{\prime }},M_{y^{\prime }},M_{z^{\prime }})$ is the magnetization vector and ${\dot{\mathit{M}}}^{\prime }$ its time derivative. We apply a rotation around x^′, i.e., perpendicular to the MW field along y^′, by $\mathit{\theta }={\mathrm{Tan}}^{-\mathrm{1}}({\mathit{\omega }}_{\mathrm{1}}/{\mathrm{\Omega }}_S)$. In the following, we will use the following shorthand notation:

$$\begin{array}{}\text{(A9)}& {\displaystyle }{R}_x\left(\mathit{\theta }\right)=\left(\begin{array}{lll}\mathrm{1}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& \cos \left(\mathit{\theta }\right)& -\sin \left(\mathit{\theta }\right)\\ \mathrm{0}& \sin \left(\mathit{\theta }\right)& \cos \left(\mathit{\theta }\right)\end{array}\right)\equiv R,\text{(A10)}& {\displaystyle }{R}_x(-\mathit{\theta })=R_x^{-\mathrm{1}}\left(\mathit{\theta }\right)=\left(\begin{array}{lll}\mathrm{1}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& \cos \left(\mathit{\theta }\right)& \sin \left(\mathit{\theta }\right)\\ \mathrm{0}& -\sin \left(\mathit{\theta }\right)& \cos \left(\mathit{\theta }\right)\end{array}\right)\equiv R^{-\mathrm{1}}.\end{array}$$

Multiplication from the left with R and insertion of $R^{-\mathrm{1}}R=\mathbb{1}$ yields the evolution of the magnetization vector in the tilted frame ${\mathit{M}}^{\prime \prime }=R{\mathit{M}}^{\prime }$.

$$\begin{array}{}\text{(A11)}& \begin{array}{rl}{\dot{\mathit{M}}}^{\prime \prime }& =R\left(\begin{array}{lll}\mathrm{0}& -{\mathrm{\Omega }}_S& {\mathit{\omega }}_{\mathrm{1}}\\ {\mathrm{\Omega }}_S& \mathrm{0}& \mathrm{0}\\ -{\mathit{\omega }}_{\mathrm{1}}& \mathrm{0}& \mathrm{0}\end{array}\right){\mathit{M}}^{\prime }\\ & =R\left(\begin{array}{lll}\mathrm{0}& -{\mathrm{\Omega }}_S& {\mathit{\omega }}_{\mathrm{1}}\\ {\mathrm{\Omega }}_S& \mathrm{0}& \mathrm{0}\\ -{\mathit{\omega }}_{\mathrm{1}}& \mathrm{0}& \mathrm{0}\end{array}\right)R^{-\mathrm{1}}R{\mathit{M}}^{\prime }\\ & =R\left(\begin{array}{lll}\mathrm{0}& -{\mathrm{\Omega }}_S& {\mathit{\omega }}_{\mathrm{1}}\\ {\mathrm{\Omega }}_S& \mathrm{0}& \mathrm{0}\\ -{\mathit{\omega }}_{\mathrm{1}}& \mathrm{0}& \mathrm{0}\end{array}\right)R^{-\mathrm{1}}{\mathit{M}}^{\prime \prime }.\end{array}\end{array}$$

To obtain the new evolution matrix, we calculate the matrix product:

$$\begin{array}{}\text{(A12)}& \begin{array}{rl}& R\left(\begin{array}{lll}\mathrm{0}& -{\mathrm{\Omega }}_S& {\mathit{\omega }}_{\mathrm{1}}\\ {\mathrm{\Omega }}_S& \mathrm{0}& \mathrm{0}\\ -{\mathit{\omega }}_{\mathrm{1}}& \mathrm{0}& \mathrm{0}\end{array}\right)R^{-\mathrm{1}}\\ & =\left(\begin{array}{lll}\mathrm{1}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& \cos \left(\mathit{\theta }\right)& -\sin \left(\mathit{\theta }\right)\\ \mathrm{0}& \sin \left(\mathit{\theta }\right)& \cos \left(\mathit{\theta }\right)\end{array}\right)\left(\begin{array}{lll}\mathrm{0}& -{\mathrm{\Omega }}_S& {\mathit{\omega }}_{\mathrm{1}}\\ {\mathrm{\Omega }}_S& \mathrm{0}& \mathrm{0}\\ -{\mathit{\omega }}_{\mathrm{1}}& \mathrm{0}& \mathrm{0}\end{array}\right)\\ & \left(\begin{array}{lll}\mathrm{1}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& \cos \left(\mathit{\theta }\right)& \sin \left(\mathit{\theta }\right)\\ \mathrm{0}& -\sin \left(\mathit{\theta }\right)& \cos \left(\mathit{\theta }\right)\end{array}\right)\\ & =\left(\begin{array}{lll}\mathrm{1}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& \cos \left(\mathit{\theta }\right)& -\sin \left(\mathit{\theta }\right)\\ \mathrm{0}& \sin \left(\mathit{\theta }\right)& \cos \left(\mathit{\theta }\right)\end{array}\right)\\ & \left(\begin{array}{lll}\mathrm{0}& -{\mathrm{\Omega }}_S\cos \left(\mathit{\theta }\right)-{\mathit{\omega }}_{\mathrm{1}}\sin \left(\mathit{\theta }\right)& -{\mathrm{\Omega }}_S\sin \left(\mathit{\theta }\right)+{\mathit{\omega }}_{\mathrm{1}}\cos \left(\mathit{\theta }\right)\\ {\mathrm{\Omega }}_S& \mathrm{0}& \mathrm{0}\\ -{\mathit{\omega }}_{\mathrm{1}}& \mathrm{0}& \mathrm{0}\end{array}\right)\\ & =\left(\begin{array}{lll}\mathrm{0}& -{\mathrm{\Omega }}_S\cos \left(\mathit{\theta }\right)-{\mathit{\omega }}_{\mathrm{1}}\sin \left(\mathit{\theta }\right)& -{\mathrm{\Omega }}_S\sin \left(\mathit{\theta }\right)+{\mathit{\omega }}_{\mathrm{1}}\cos \left(\mathit{\theta }\right)\\ {\mathrm{\Omega }}_S\cos \left(\mathit{\theta }\right)+{\mathit{\omega }}_{\mathrm{1}}\sin \left(\mathit{\theta }\right)& \mathrm{0}& \mathrm{0}\\ {\mathrm{\Omega }}_S\sin \left(\mathit{\theta }\right)-{\mathit{\omega }}_{\mathrm{1}}\cos \left(\mathit{\theta }\right)& \mathrm{0}& \mathrm{0}\end{array}\right)\\ & =\left(\begin{array}{lll}\mathrm{0}& -{\mathit{\omega }}_{\mathrm{eff}}& \mathrm{0}\\ {\mathit{\omega }}_{\mathrm{eff}}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& \mathrm{0}& \mathrm{0}\end{array}\right),\end{array}\end{array}$$

where the last transformation step is performed with Ω_S=ω_effcos (θ) and ω₁=ω_effsin (θ). Insertion into Eq. (A11) yields the differential equations for the tilted-frame magnetization. In this tilted-frame expression, we add the transversal relaxation rates on the matrix diagonals, thus defining the decay of $M_{x^{\prime \prime }}$ and $M_{y^{\prime \prime }}$.

$$\begin{array}{}\text{(A13)}& \begin{array}{rl}{\dot{\mathit{M}}}^{\prime \prime }& =\left(\begin{array}{lll}-\mathrm{1}/T_{\mathrm{2}\mathit{\rho }}& -{\mathit{\omega }}_{\mathrm{eff}}& \mathrm{0}\\ {\mathit{\omega }}_{\mathrm{eff}}& -\mathrm{1}/T_{\mathrm{2}\mathit{\rho }}& \mathrm{0}\\ \mathrm{0}& \mathrm{0}& \mathrm{0}\end{array}\right){\mathit{M}}^{\prime \prime }\\ & =\left(\begin{array}{l}-\mathrm{1}/T_{\mathrm{2}\mathit{\rho }}\cdot M_{x^{\prime \prime }}-{\mathit{\omega }}_{\mathrm{eff}}\cdot M_{y^{\prime \prime }}\\ {\mathit{\omega }}_{\mathrm{eff}}\cdot M_{x^{\prime \prime }}-\mathrm{1}/T_{\mathrm{2}\mathit{\rho }}\cdot M_{y^{\prime \prime }}\\ \mathrm{0}\end{array}\right).\end{array}\end{array}$$

This set of differential equations is solved by the standard approach of defining $M_{{+}^{\prime \prime }}=M_{x^{\prime \prime }}+i{M}_{y^{\prime \prime }}$.

$$\begin{array}{}\text{(A14)}& \begin{array}{rl}{\dot{M}}_{{+}^{\prime \prime }}& ={\dot{M}}_{x^{\prime \prime }}+i{\dot{M}}_{y^{\prime \prime }}={\mathit{\omega }}_{\mathrm{eff}}\left[-M_{y^{\prime \prime }}+i{M}_{x^{\prime \prime }}\right]\\ & -\mathrm{1}/T_{\mathrm{2}\mathit{\rho }}\cdot \left[M_{x^{\prime \prime }}+i{M}_{y^{\prime \prime }}\right]\\ & =\left[i{\mathit{\omega }}_{\mathrm{eff}}-\mathrm{1}/T_{\mathrm{2}\mathit{\rho }}\right]M_{{+}^{\prime \prime }}\end{array}\end{array}$$

The solution of this differential equation is an exponential of the following form:

$$\begin{array}{}\text{(A15)}& \begin{array}{rl}{M}_{{+}^{\prime \prime }}\left(t\right)& =\left[M_{x^{\prime \prime }}\left(\mathrm{0}\right)+i{M}_{y^{\prime \prime }}\left(\mathrm{0}\right)\right]\\ & \cdot \left[\cos \left({\mathit{\omega }}_{\mathrm{eff}}t\right)+i\sin \left({\mathit{\omega }}_{\mathrm{eff}}t\right)\right]\cdot e^{-t/T_{\mathrm{2}\mathit{\rho }}}.\end{array}\end{array}$$

After factorization and separating the expression for $M_{{+}^{\prime \prime }}\left(t\right)$ into its real and imaginary parts, an expression for M^′′(t) is obtained.

$$\begin{array}{}\text{(A16)}& {\mathit{M}}^{\prime \prime }\left(t\right)=\left(\begin{array}{l}\left[M_{x^{\prime \prime }}\left(\mathrm{0}\right)\cos \left({\mathit{\omega }}_{\mathrm{eff}}t\right)\right.\\ \left.-M_{y^{\prime \prime }}\left(\mathrm{0}\right)\sin \left({\mathit{\omega }}_{\mathrm{eff}}t\right)\right]\cdot e^{-t/T_{\mathrm{2}\mathit{\rho }}}\\ \left[M_{y^{\prime \prime }}\left(\mathrm{0}\right)\cos \left({\mathit{\omega }}_{\mathrm{eff}}t\right)\right.\\ \left.+M_{x^{\prime \prime }}\left(\mathrm{0}\right)\sin \left({\mathit{\omega }}_{\mathrm{eff}}t\right)\right]\cdot e^{-t/T_{\mathrm{2}\mathit{\rho }}}\\ M_{z^{\prime \prime }}\left(\mathrm{0}\right)\end{array}\right)\end{array}$$

The initial conditions in the tilted frame are calculated from the rotating-frame magnetization at t=0, which is ${\mathit{M}}^{\prime }\left(\mathrm{0}\right)=(\mathrm{0},\mathrm{0},M_{\mathrm{0}})$, so that, after back-transformation into the rotating frame, we end up with an analytical expression for M^′(t).

$$\begin{array}{}\text{(A17)}& {\displaystyle }\begin{array}{rl}& {\mathit{M}}^{\prime \prime }\left(\mathrm{0}\right)=R\cdot M^{\prime }\left(\mathrm{0}\right)=\left(\begin{array}{lll}\mathrm{1}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& \cos \left(\mathit{\theta }\right)& -\sin \left(\mathit{\theta }\right)\\ \mathrm{0}& \sin \left(\mathit{\theta }\right)& \cos \left(\mathit{\theta }\right)\end{array}\right)\left(\begin{array}{l}\mathrm{0}\\ \mathrm{0}\\ M_{\mathrm{0}}\end{array}\right)\\ & =\left(\begin{array}{l}\mathrm{0}\\ -\sin \left(\mathit{\theta }\right)M_{\mathrm{0}}\\ \cos \left(\mathit{\theta }\right)M_{\mathrm{0}}\end{array}\right),\\ & {\mathit{M}}^{\prime \prime }\left(t\right)=\left(\begin{array}{l}\sin \left({\mathit{\omega }}_{\mathrm{eff}}t\right)\sin \left(\mathit{\theta }\right)M_{\mathrm{0}}\cdot e^{-t/T_{\mathrm{2}\mathit{\rho }}}\\ -\cos \left({\mathit{\omega }}_{\mathrm{eff}}t\right)\sin \left(\mathit{\theta }\right)M_{\mathrm{0}}\cdot e^{-t/T_{\mathrm{2}\mathit{\rho }}}\\ \cos \left(\mathit{\theta }\right)M_{\mathrm{0}}\end{array}\right),\\ & {\mathit{M}}^{\prime }\left(t\right)=R^{-\mathrm{1}}{\mathit{M}}^{\prime \prime }\left(t\right)=\left(\begin{array}{lll}\mathrm{1}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& \cos \left(\mathit{\theta }\right)& \sin \left(\mathit{\theta }\right)\\ \mathrm{0}& -\sin \left(\mathit{\theta }\right)& \cos \left(\mathit{\theta }\right)\end{array}\right)\\ & \left(\begin{array}{l}\sin \left({\mathit{\omega }}_{\mathrm{eff}}t\right)\sin \left(\mathit{\theta }\right)M_{\mathrm{0}}\cdot e^{-t/T_{\mathrm{2}\mathit{\rho }}}\\ -\cos \left({\mathit{\omega }}_{\mathrm{eff}}t\right)\sin \left(\mathit{\theta }\right)M_{\mathrm{0}}\cdot e^{-t/T_{\mathrm{2}\mathit{\rho }}}\\ \cos \left(\mathit{\theta }\right)M_{\mathrm{0}}\end{array}\right)\\ & =\left(\begin{array}{l}\sin \left({\mathit{\omega }}_{\mathrm{eff}}t\right)\sin \left(\mathit{\theta }\right)M_{\mathrm{0}}\cdot e^{-t/T_{\mathrm{2}\mathit{\rho }}}\\ \left[-\cos \left({\mathit{\omega }}_{\mathrm{eff}}t\right)\sin \left(\mathit{\theta }\right)\cos \left(\mathit{\theta }\right)\cdot e^{-t/T_{\mathrm{2}\mathit{\rho }}}+\sin \left(\mathit{\theta }\right)\cos \left(\mathit{\theta }\right)\right]M_{\mathrm{0}}\\ \left[\cos \left({\mathit{\omega }}_{\mathrm{eff}}t\right){\sin }^{\mathrm{2}}\left(\mathit{\theta }\right)\cdot e^{-t/T_{\mathrm{2}\mathit{\rho }}}+{\cos }^{\mathrm{2}}\left(\mathit{\theta }\right)\right]M_{\mathrm{0}}\end{array}\right).\end{array}\end{array}$$

In the limit of long evolution times, all terms with $e^{-t/T_{\mathrm{2}\mathit{\rho }}}$ decay to zero, and a reduced magnetization vector remains:

$$\begin{array}{}\text{(A18)}& {\displaystyle }\begin{array}{rl}\mathit{M}(t⟶\mathrm{\infty })& =\underset{t⟶\mathrm{\infty }}{lim}\left(\begin{array}{l}\sin \left({\mathit{\omega }}_{\mathrm{eff}}t\right)\sin \left(\mathit{\theta }\right)M_{\mathrm{0}}\cdot e^{-t/T_{\mathrm{2}\mathit{\rho }}}\\ \left[-\cos \left({\mathit{\omega }}_{\mathrm{eff}}t\right)\sin \left(\mathit{\theta }\right)\cos \left(\mathit{\theta }\right)\right.\\ \left.\cdot e^{-t/T_{\mathrm{2}\mathit{\rho }}}+\sin \left(\mathit{\theta }\right)\cos \left(\mathit{\theta }\right)\right]M_{\mathrm{0}}\\ \left[+\cos \left({\mathit{\omega }}_{\mathrm{eff}}t\right){\sin }^{\mathrm{2}}\left(\mathit{\theta }\right)\right.\\ \left.\cdot e^{-t/T_{\mathrm{2}\mathit{\rho }}}+{\cos }^{\mathrm{2}}\left(\mathit{\theta }\right)\right]M_{\mathrm{0}}\end{array}\right)\\ & =\left(\begin{array}{l}\mathrm{0}\\ \sin \left(\mathit{\theta }\right)\cos \left(\mathit{\theta }\right)M_{\mathrm{0}}\\ {\cos }^{\mathrm{2}}\left(\mathit{\theta }\right)M_{\mathrm{0}}\end{array}\right)\\ & =\left(\begin{array}{l}\mathrm{0}\\ \frac{\left(\mathit{\nu }-{\mathit{\nu }}_{\mathrm{off}}\right){\mathit{\nu }}_{\mathrm{1}}}{{\mathit{\nu }}_{\mathrm{1}}^{\mathrm{2}}+{\left(\mathit{\nu }-{\mathit{\nu }}_{\mathrm{off}}\right)}^{\mathrm{2}}}\\ \mathrm{1}-\frac{{\mathit{\nu }}_{\mathrm{1}}^{\mathrm{2}}}{{\mathit{\nu }}_{\mathrm{1}}^{\mathrm{2}}+{\left(\mathit{\nu }-{\mathit{\nu }}_{\mathrm{off}}\right)}^{\mathrm{2}}}\end{array}\right),\end{array}\end{array}$$

where the last step of the rearrangement was done using Mathematica.

A16 Solution of Bloch equations in the rotating frame

When adding the transversal relaxation terms in the rotating frame, Eq. (A19) is obtained (Bloch, 1946):

$$\begin{array}{}\text{(A19)}& \begin{array}{rl}{\dot{\mathit{M}}}^{\prime }& =\left(\begin{array}{l}{\dot{M}}_{x^{\prime }}\\ {\dot{M}}_{y^{\prime }}\\ {\dot{M}}_{z^{\prime }}\end{array}\right)\\ & =\left(\begin{array}{lll}-\mathrm{1}/T_{\mathrm{2}}& -{\mathrm{\Omega }}_S& {\mathit{\omega }}_{\mathrm{1}}\\ {\mathrm{\Omega }}_S& -\mathrm{1}/T_{\mathrm{2}}& \mathrm{0}\\ -{\mathit{\omega }}_{\mathrm{1}}& \mathrm{0}& \mathrm{0}\end{array}\right)\left(\begin{array}{l}{M}_{x^{\prime }}\\ M_{y^{\prime }}\\ M_{z^{\prime }}\end{array}\right).\end{array}\end{array}$$

As the analytical solution of Eq. (A19) is rather cumbersome and provides little intuitive insight into the physical processes (Torrey, 1949), we used a numerical differential equation solver implemented in MATLAB (ode45) to obtain the results displayed in Fig. 8.

Appendix B: Additional information on the figures in the main text

In the following, the experimental conditions and data processing steps of the figures in the main text are given. This includes a reference to the experimental data sets and processing scripts that are uploaded in a separate data repository (https://doi.org/10.25625/B11CUC, Lenjer, 2024). The experimental data sets are provided as .DTA and .DSC files and were processed with the eprload command implemented in EasySpin (Stoll and Schweiger, 2006). The processing and simulation files are provided as MATLAB code. Additionally, the analysis files for the figures in Appendix A can be found in the data repository.

B1 Fig. 2 (SL Rabi nutations)

SL Rabi nutations obtained with the sequence in Fig. 2a had the following parameters: T=100 K, $t_{\mathit{\pi }/\mathrm{2}}=\mathrm{18}\mathrm{ns}$, t_π=36 ns, τ=1 µs, τ₀=1 µs, τ₁=2 µs, a_PM=0.4, shot repetition time (SRT) =40 ms, shots per point (SPP) = 4, and scans = 3. A four-step phase cycle of the π refocusing pulse was used:

$${\mathit{\varphi }}_{\mathit{\pi }}=[\mathrm{0}\mathit{°},\mathrm{90}\mathit{°},\mathrm{180}\mathit{°},\mathrm{270}\mathit{°}],{\mathit{\varphi }}_{\det }=[+,-,+,-].$$

Rabi traces were obtained by recording the integrated echo intensity as a function of t_PM. For the Fourier transformation, the traces were corrected for phase and baseline and apodized with a half-cosine bell function. This apodized data set was filled with zeros with twice the number of data points and Fourier transformed using the fft command in MATLAB. The simulation procedure is illustrated in Fig. 2c. The same overall trace length and FT procedure as in the experiment was used for the simulations. Experimental data sets: 999a–999f. Processing file: SL_Rabi_nutations_analysis.m.

B2 Fig. 3 (SL echoes and WOP SL echoes)

SL echoes and WOP SL echoes obtained with the sequences in Fig. 3a and b, respectively, used the following parameters: T=100 K, $t_{\mathit{\pi }/\mathrm{2}}=\mathrm{18}\mathrm{ns}$, t_π=36 ns, τ=0.3 µs, τ₀=2 µs, τ₁=3 µs, ${\mathit{\tau }}_{\mathrm{3}}=\mathrm{10.7}\mathrm{\mu }\mathrm{s}-{\mathit{\tau }}_{\mathrm{2}}$, a_PM=0.4, ν_PM=15 MHz, $t_{\mathrm{PM}\mathrm{1}}=t_{\mathrm{PM}\mathrm{3}}=\mathrm{80}\mathrm{ns}$, t_PM2=160 ns, SRT =40 ms, SPP = 1–4, and scans = 5. In addition to the PCs denoted in the figure, a two-step PC of PM3 was performed in all cases:

$${\mathit{\varphi }}_{\mathrm{PM}\mathrm{3}}=[\mathrm{0}\mathit{°},\mathrm{180}\mathit{°}],{\mathit{\varphi }}_{\det }=[+,-].$$

Echoes were obtained by recording the integrated echo intensity as a function of τ₂. The echo traces were phase corrected so that the dominant signal in the trace without PC of PM2 was maximized in the real part of the signal. The experiments with different PCs were subsequently corrected with the same phase factor. Experimental data sets: 967a–967f. Processing files: SL_echo_analysis.m and WOP_SL_echo_analysis.m.

B3 Fig. 4 (SL relaxation measurements)

SL echo decay traces were obtained using the sequence in Fig. 3a. The following parameters were used for the trace in Fig. 4a: T=100 K, $t_{\mathit{\pi }/\mathrm{2}}=\mathrm{20}\mathrm{ns}$, t_π=40 ns, τ=0.3 µs, τ₀=100 µs, τ₁=τ₂, ${\mathit{\tau }}_{\mathrm{3}}=\mathrm{32.2}\mathrm{\mu }\mathrm{s}-{\mathit{\tau }}_{\mathrm{1}}-{\mathit{\tau }}_{\mathrm{2}}$, a_PM=0.4, ν_PM=15 MHz, $t_{\mathrm{PM}\mathrm{1}}=t_{\mathrm{PM}\mathrm{3}}=\mathrm{80}\mathrm{ns}$, t_PM2=160 ns, SRT = 40 ms, SPP = 4, and scans = 3. The trace was obtained by recording the integrated echo intensity as a function of τ₁=τ₂. A two-step PC of PM3 was performed. The bare-state echo decay was obtained using the Hahn echo sequence ($\mathit{\pi }/\mathrm{2}-\mathit{\tau }-\mathit{\pi }-\mathit{\tau }-\mathrm{echo}$) with $t_{\mathit{\pi }/\mathrm{2}}=\mathrm{22}\mathrm{ns}$ and t_π=44 ns. The experiments in (b) were performed with similar parameters, but t_PM was, for each experiment, adjusted to the different SL Rabi frequencies. Hahn echoes with lower ν₁ were detected with a $\mathrm{50}/\mathrm{100}\mathrm{ns}$ and $\mathrm{200}/\mathrm{400}\mathrm{ns}$ echo. Carr–Purcell experiments were performed with $t_{\mathit{\pi }/\mathrm{2}}=\mathrm{22}\mathrm{ns}$, and t_π=44 ns. A 16-step PC of the last two π pulses was applied:

$$\begin{array}{c}{\displaystyle }{\mathit{\varphi }}_{\mathit{\pi },N}=[\mathrm{0}\mathit{°}{]}_{\mathrm{4}},[\mathrm{90}\mathit{°}{]}_{\mathrm{4}},[\mathrm{180}\mathit{°}{]}_{\mathrm{4}},[\mathrm{270}\mathit{°}{]}_{\mathrm{4}},\\ {\displaystyle }{\mathit{\varphi }}_{\mathit{\pi },N+\mathrm{1}}=[\mathrm{0}\mathit{°},\mathrm{90}\mathit{°},\mathrm{180}\mathit{°},\mathrm{270}\mathit{°}{]}_{\mathrm{4}},\\ {\displaystyle }{\mathit{\varphi }}_{\det }=[+,-,+,-,-,+,-,+{]}_{\mathrm{2}}.\end{array}$$

Experimental data sets: 1018, 1131, and 1134. Processing file: T2rho_analysis.m.

B4 Figs. 5 and 6 (CHEESY z^′ profiles)

CHEESY z^′ profiles were measured using the sequence in Fig. 5a and the following parameters: T=100 K, $t_{\mathit{\pi }/\mathrm{2},\mathrm{chirp}}=\mathrm{1}\mathrm{\mu }\mathrm{s}$, $t_{\mathit{\pi },\mathrm{chirp}}=\mathrm{0.5}\mathrm{\mu }\mathrm{s}$, τ=3 µs, SRT =40 ms, SPP = 150, and scans = 1. An eight-step PC of the chirp echo was used:

$$\begin{array}{c}{\displaystyle }{\mathit{\varphi }}_{\mathit{\pi }/\mathrm{2},\mathrm{chirp}}=[\mathrm{0}\mathit{°}{]}_{\mathrm{4}},[\mathrm{180}\mathit{°}{]}_{\mathrm{4}},\\ {\displaystyle }{\mathit{\varphi }}_{\mathit{\pi },\mathrm{chirp}}=[\mathrm{0}\mathit{°},\mathrm{90}\mathit{°},\mathrm{180}\mathit{°},\mathrm{270}\mathit{°}{]}_{\mathrm{2}},\\ {\displaystyle }{\mathit{\varphi }}_{\det }=[+,-{]}_{\mathrm{2}},[-,+{]}_{\mathrm{2}}.\end{array}$$

The chirp pulses had a sweep width of 200 MHz centered at the maximum of the resonator dip. The frequency sweep of the chirp pulses was adjusted to compensate for the resonator profile (FWHM≈140 MHz), yielding an approximately constant adiabaticity (Doll and Jeschke, 2014). A wideband, uniform rate, smooth truncation (WURST) amplitude modulation (Kupce and Freeman, 1995) was applied to yield the final chirp pulse.

Before the experiment, the magnetic field was adjusted so that the BDPA EPR line was centered at ${\mathit{\nu }}_{\mathrm{off}}=-\mathrm{31}\mathrm{MHz}$, and the analyzed pulse was applied with the same frequency offset. FT was performed using the fft command in MATLAB after choosing a FT window of 5 µs symmetrically around the echo, apodization with a cosine bell function and zero filling with the number of zeros equal to the length of the data set. Zeroth-order phase correction was performed so that the integral of the FT was maximized in the real part of the signal. Experimental data set: 1025. Processing files: CHEESY_z_profile_analysis.m and CHEESY_z_profile_SimSpinach.m.

B5 Fig. 7 (CHEESY $x^{\prime }{y}^{\prime }$ profiles)

The CHEESY $x^{\prime }{y}^{\prime }$ profiles in Fig. 7c and d were obtained using the pulse sequence in (a) and T=100 K. The chirp pulses were identical to the π chirp pulse in the previous section. The delays were set to τ₁=1 µs and τ₂=4 µs, together with SRT =40 ms, SPP = 250, scans = 1, and an offset ${\mathit{\nu }}_{\mathrm{off}}=-\mathrm{29}\mathrm{MHz}$. A 32-step PC was used:

$$\begin{array}{c}{\displaystyle }{\mathit{\varphi }}_{\mathrm{analyzed}}=[\mathrm{0}\mathit{°}{]}_{\mathrm{16}},[\mathrm{180}\mathit{°}{]}_{\mathrm{16}};\\ {\displaystyle }{\mathit{\varphi }}_{\mathit{\pi },\mathrm{chirp}\left(\mathrm{1}\right)}=\mathit{\left\{}\right[\mathrm{0}\mathit{°}{]}_{\mathrm{4}},[\mathrm{90}\mathit{°}{]}_{\mathrm{4}},[\mathrm{180}\mathit{°}{]}_{\mathrm{4}},[\mathrm{270}\mathit{°}{]}_{\mathrm{4}}{\mathit{\}}}_{\mathrm{2}};\\ {\displaystyle }{\mathit{\varphi }}_{\mathit{\pi },\mathrm{chirp}\left(\mathrm{2}\right)}=[\mathrm{0}\mathit{°},\mathrm{90}\mathit{°},\mathrm{180}\mathit{°},\mathrm{270}\mathit{°}{]}_{\mathrm{8}}:\\ {\displaystyle }{\mathit{\varphi }}_{\det }=\mathit{\left\{}\right[+,-{]}_{\mathrm{2}},[-,+{]}_{\mathrm{2}}{\mathit{\}}}_{\mathrm{2}},\mathit{\{}[-,+{]}_{\mathrm{2}},[+,-{]}_{\mathrm{2}}{\mathit{\}}}_{\mathrm{2}})\end{array}$$

After manual selection of an FT window with 4 µs length, FT was performed as in the previous section. Zeroth-order phase correction was performed so that the dispersive signal was maximized in the imaginary channel.

The measurement of T_1ρ (e) was performed with the sequence in (b), where $t_{\mathit{\pi }/\mathrm{2}}=\mathrm{14}\mathrm{ns}$, t_π=28 ns, τ=1 µs, and ${\mathit{\nu }}_{\mathrm{1},\mathrm{SL}}=\mathrm{0.8}\mathrm{MHz}$. The intensity was scaled to overlap with the data points from the R-CHEESY experiment. Experimental data sets: 875, 876, and 968. Processing file: CHEESY_xy_profile_analysis.m.

B6 Fig. 8 (Bloch equations)

The experimental data shown in (a) and (b) are the same as in Figs. 5c and 7c. Processing files: CHEESY_z_profile_analysis.m, CHEESY_xy_profile_analysis.m., and Bloch_simulations.m