Ultrafast spectroscopy is often collected in the mixed frequency/time domain, where pulse durations are similar to system dephasing times. In these experiments, expectations derived from the familiar driven and impulsive limits are not valid. This work simulates the mixed-domain four-wave mixing response of a model system to develop expectations for this more complex field-matter interaction. We explore frequency and delay axes. We show that these line shapes are exquisitely sensitive to excitation pulse widths and delays. Near pulse overlap, the excitation pulses induce correlations that resemble signatures of dynamic inhomogeneity. We describe these line shapes using an intuitive picture that connects to familiar field-matter expressions. We develop strategies for distinguishing pulse-induced correlations from true system inhomogeneity. These simulations provide a foundation for interpretation of ultrafast experiments in the mixed domain.

Ultrafast spectroscopy is based on using nonlinear interactions, created by multiple ultrashort (10−9–10−15 s) pulses, to resolve spectral information on time scales as short as the pulses themselves.1,2 The ultrafast spectra can be collected in the time domain or in the frequency domain.3 

Time-domain methods scan the pulse delays to resolve the free induction decay (FID).4 The Fourier transform of the FID gives the ultrafast spectrum. Ideally, these experiments are performed in the impulsive limit where FID dominates the measurement. FID occurs at the frequency of the transition that has been excited by a well-defined, time-ordered sequence of pulses. Time-domain methods are compromised when the dynamics occur on faster time scales than the ultrafast excitation pulses. As the pulses temporally overlap, FID from other pulse time-orderings and emission driven by the excitation pulses both become important. These factors are responsible for the complex “coherent artifacts” that are often ignored in pump-probe and related methods.5–8 Dynamics faster than the pulse envelopes are best measured using line shapes in frequency domain methods.

Frequency-domain methods scan pulse frequencies to resolve the ultrafast spectrum directly.9,10 Ideally, these experiments are performed in the driven limit where the steady state dominates the measurement. In the driven limit, all time-orderings of the pulse interactions are equally important and FID decay is negligible. The output signal is driven at the excitation pulse frequencies during the excitation pulse width. Frequency-domain methods are compromised when the spectral line shape is narrower than the frequency bandwidth of the excitation pulses. Dynamics that are slower than the pulse envelopes can be measured in the time domain by resolving the phase oscillations of the output signal during the entire FID decay.

There is also the hybrid mixed-time/frequency-domain approach, where pulse delays and pulse frequencies are both scanned to measure the system response. This approach is uniquely suited for experiments where the dephasing time is comparable to the pulse durations so that neither frequency-domain nor time-domain approaches excel on their own.10–12 In this regime, both FID and driven processes are important.13 Their relative importance depends on pulse frequencies and delays. Extracting the correct spectrum from the measurement then requires a more complex analysis that explicitly treats the excitation pulses and the different time-orderings.14–16 Despite these complications, mixed-domain methods have a practical advantage: the dual frequency- and delay-scanning capabilities allow these methods to address a wide variety of dephasing rates.

The relative importance of FID and driven processes and the changing importance of different coherence pathways are important factors for understanding spectral features in all ultrafast methods. These methods include partially coherent methods involving intermediate populations such as pump-probe spectroscopy,17 transient grating spectroscopy,18–20 transient absorption/reflection spectroscopy,21,22 photon echo spectroscopy,23–25 two dimensional-infrared (2D-IR) spectroscopy,26–28 2D-electronic spectroscopy (2D-ES),29,30 and three-pulse photon echo peak shift (3PEPS)23,31–34 spectroscopy. These methods also include fully coherent methods involving only coherences such as Stimulated Raman Spectroscopy (SRS),35,36 Doubly Vibrationally Enhanced (DOVE),37–43 Triply Resonant Sum Frequency (TRSF),44–46 Sum Frequency Generation (SFG),47 Coherent Anti-Stokes Raman Spectroscopy (CARS),48–50 and other coherent Raman methods.51 

This paper focuses on understanding the nature of the spectral changes that occur in Coherent Multidimensional Spectroscopy (CMDS) as experiments transition between the two limits of frequency- and time-domain methods. CMDS is a family of spectroscopies that use multiple delay and/or frequency axes to extract homogeneous and inhomogeneous broadening, as well as detailed information about spectral diffusion and chemical changes.52,53 For time-domain CMDS (2D-IR, 2D-ES), the complications that occur when the impulsive approximation does not strictly hold have only recently been addressed.54,55

Frequency-domain CMDS methods, referred to herein as multi-resonant CMDS (MR-CMDS), have similar capabilities for measuring homogeneous and inhomogeneous broadening. Although these experiments are typically described in the driven limit,4,19,20 many of the experiments involve pulse widths that are comparable to the widths of the system.15,37,40,41,56,57 MR-CMDS then becomes a mixed-domain experiment whereby resonances are characterized with marginal resolution in both frequency and time. For example, DOVE spectroscopy involves three different pathways58 whose relative importance depends on the relative importance of FID and driven responses.59 In the driven limit, the DOVE line shape depends on the difference between the first two pulse frequencies, so the line shape has a diagonal character that mimics the effects of inhomogeneous broadening. In the FID limit where the coherence frequencies are defined instead by the transition, the diagonal character is lost. Understanding these effects is crucial for interpreting experiments, yet these effects have not been characterized for MR-CMDS.

This work considers the third-order MR-CMDS response of a 3-level model system using three ultrafast excitation beams with the commonly used four-wave mixing (FWM) phase-matching condition, kout=k1k2+k2. Here, the subscripts represent the excitation pulse frequencies, ω1 and ω2=ω2. These experimental conditions were recently used to explore line shapes of excitonic systems15,57 and have been developed on vibrational states as well.60 Although MR-CMDS forms the context of this model, the treatment is quite general because the phase-matching condition can describe any of the spectroscopies mentioned above with the exception of SFG and TRSF, for which the model can be easily extended. We numerically simulate the MR-CMDS response with pulse durations at, above, and below the system coherence time. To highlight the role of pulse effects, we build an interpretation of the full MR-CMDS response by first showing how finite pulses affect the evolution of a coherence and then how finite pulses affect an isolated third-order pathway. When considering the full MR-CMDS response, we show that spectral features change dramatically as a function of delay, even for a homogeneous system with elementary dynamics. Importantly, the line shape can exhibit correlations that mimic inhomogeneity, and the temporal evolution of this line shape can mimic spectral diffusion. We identify key signatures that can help differentiate true inhomogeneity and spectral diffusion from these measurement artifacts.

We consider a simple three-level system (states n = 0, 1, 2) that highlights the multidimensional line shape changes resulting from choices of the relative dephasing and detuning of the system and the temporal and spectral widths of the excitation pulses. For simplicity, we will ignore population relaxation effects: Γ11=Γ00=0.

The electric field pulses, El, are given by

El(t;ωl,τl,klz)=12cl(tτl)eiklzeiωl(tτl)+c.c.,
(1)

where ωl is the field carrier frequency, kl is the wavevector, τl is the pulse delay, and cl is a slowly varying envelope. In this work, we assume normalized (real-valued) Gaussian envelopes

cl(t)=1Δt2ln22πexpln2tΔt2,
(2)

where Δt is the temporal FWHM of the envelope intensity. We neglect non-linear phase effects such as chirp, so the FWHM of the frequency bandwidth is transform limited: ΔωΔt=4ln22.77, where Δω is the spectral FWHM (intensity scale).

The Liouville-von Neumann equation propagates the density matrix, ρ,

dρdt=iH0+μl=1,2,2El(t),ρ+Γρ.
(3)

Here H0 is the time-independent Hamiltonian, μ is the dipole superoperator, and Γ contains the pure dephasing rate of the system. We perform the standard perturbative expansion of Eq. (3) to third order in the electric field interaction10,61–64 and restrict ourselves only to the terms that have the correct spatial wave vector kout=k1k2+k2. This approximation narrows the scope to sets of three interactions, one from each field, that result in the correct spatial dependence. The set of three interactions have 3! = 6 unique time-ordered sequences, and each time-ordering produces either two or three unique system-field interactions for our system, for a total of sixteen unique system-field interaction sequences, or Liouville pathways, to consider. Figure 1 shows these sixteen pathways as Wave Mixing Energy Level (WMEL) diagrams.65 

FIG. 1.

The sixteen triply resonant Liouville pathways for the third-order response of the system used here. Time flows from left to right. Each excitation is labeled by the pulse stimulating the transition; excitations with ω1 are yellow, excitations with ω2=ω2 are purple, and the final emission is gray.

FIG. 1.

The sixteen triply resonant Liouville pathways for the third-order response of the system used here. Time flows from left to right. Each excitation is labeled by the pulse stimulating the transition; excitations with ω1 are yellow, excitations with ω2=ω2 are purple, and the final emission is gray.

Close modal

We first focus on a single interaction in these sequences, where an excitation pulse, x, forms ρij from ρik or ρkj. For brevity, we use =1 and abbreviate the initial and final density matrix elements as ρi and ρf, respectively. Using the natural frequency rotating frame, ρ̃ij=ρijeiωijt, the formation of ρf using pulse x is written as

dρ̃fdt=Γfρ̃f+i2λfμfcx(tτx)eiκfkxz+ωxτxeiκfΩfxtρ̃i(t),
(4)

where Ωfx=κf1ωfωx(=|ωf|ωx) is the detuning, ωf is the transition frequency of the ith transition, μf is the transition dipole, and Γf is the dephasing/relaxation rate for ρf. The λf and κf parameters describe the phases of the interaction: λf=+1 for ket-side transitions and −1 for bra-side transitions, and κf depends on whether ρf is formed via absorption (κf=λf) or emission (κf=λf).66 In the following equations, we neglect spatial dependence (z = 0).

Equation (4) forms the basis for our simulations. It provides a general expression for arbitrary values of the dephasing rate and excitation pulse bandwidth. The integral solution is

ρ̃f(t)=i2λfμfeiκfωxτxeiκfΩfxtcx(tuτx)×ρ̃i(tu)Θ(u)eΓf+iκfΩfxudu,
(5)

where Θ is the Heaviside step function. Equation (5) becomes the steady state limit expression when Δt|Γf+iκfΩfx|>>1, and the impulsive limit expression results when Δt|Γf+iκfΩfx|<<1. Both limits are important for understanding the multidimensional line shape changes discussed in this paper. The steady state and impulsive limits of Eq. (5) are discussed in  Appendix A.

Figure 2 gives an overview of the simulations done in this work. Figure 2(a) shows an excitation pulse (gray-shaded) and examples of a coherent transient for three different dephasing rates. The color bindings to dephasing rates introduced in Fig. 2(a) will be used consistently throughout this work. Our simulations use systems with dephasing rates quantified relative to the pulse duration: Γ10Δt=0.5,1, or 2. The temporal axes are normalized to the pulse duration, Δt. The Γ10Δt=2 transient is mostly driven by the excitation pulse, while Γ10Δt=0.5 has a substantial free induction decay (FID) component at late times. Figure 2(b) shows a pulse sequence (pulses are shaded orange and pink) and the resulting system evolution of pathway Vγ (00201211110out00) with Γ10Δt=1. The final polarization (teal) is responsible for the emitted signal, which is then passed through a frequency bandpass filter to emulate monochromator detection [Fig. 2(c)]. The resulting signal is explored in 2D delay space [Fig. 2(d)], 2D frequency space [Fig. 2(f)], and hybrid delay-frequency space [Fig. 2(e)]. The detuning frequency axes are also normalized by the pulse bandwidth, Δω.

FIG. 2.

Overview of the MR-CMDS simulation. (a) The temporal profile of a coherence under pulsed excitation depends on how quickly the coherence dephases. In all subsequent panels, the relative dephasing rate is kept constant at Γ10Δt=1. (b) Simulated evolution of the density matrix elements of a third-order Liouville pathway Vγ under fully resonant excitation. Pulses can be labeled both by their time of arrival (x, y, z) and by the lab lasers used to stimulate the transitions (2, 2, 1). The final coherence (teal) creates the output electric field. (c) The frequency profile of the output electric field is filtered by a monochromator gating function, M(ω), and the passed components (shaded) are measured. [(d)–(f)] Signal is viewed against two laser parameters, either as 2D delay (d), mixed delay-frequency (e), or 2D frequency plots (f). The six time-orderings are labeled in (d) to help introduce our delay convention.

FIG. 2.

Overview of the MR-CMDS simulation. (a) The temporal profile of a coherence under pulsed excitation depends on how quickly the coherence dephases. In all subsequent panels, the relative dephasing rate is kept constant at Γ10Δt=1. (b) Simulated evolution of the density matrix elements of a third-order Liouville pathway Vγ under fully resonant excitation. Pulses can be labeled both by their time of arrival (x, y, z) and by the lab lasers used to stimulate the transitions (2, 2, 1). The final coherence (teal) creates the output electric field. (c) The frequency profile of the output electric field is filtered by a monochromator gating function, M(ω), and the passed components (shaded) are measured. [(d)–(f)] Signal is viewed against two laser parameters, either as 2D delay (d), mixed delay-frequency (e), or 2D frequency plots (f). The six time-orderings are labeled in (d) to help introduce our delay convention.

Close modal

We now consider the generalized Liouville pathway L:ρ0xρ1yρ2zρ3outρ4, where x, y, and z denote properties of the first, second, and third pulses, respectively, and indices 0, 1, 2, 3, and 4 define the properties of the ground state, first, second, third, and fourth density matrix elements, respectively. Figure 2(b) demonstrates the correspondence between x, y, z notation and 1, 2, 2 notation for the laser pulses with pathway Vγ.67 

The electric field emitted from a Liouville pathway is proportional to the polarization created by the third-order coherence,

EL(t)=iμ4ρ3(t).
(6)

Equation (6) assumes perfect phase-matching and no pulse distortions through propagation. Equation (5) shows that the output field for this Liouville pathway is

EL(t)=i8λ1λ2λ3μ1μ2μ3μ4eiκ1ωxτx+κ2ωyτy+κ3ωzτz×eiκ3ωz+κ2ωy+κ1ωxtcz(tuτz)×cy(tuvτy)cx(tuvwτx)×RL(u,v,w)dwdvdu,
(7)
RL(u,v,w)=Θ(w)eΓ1+iκ1Ω1xwΘ(v)eΓ2+iκ1Ω1x+κ2Ω2yv×Θ(u)eΓ3+iκ1Ω1x+κ2Ω2y+κ3Ω3zu,
(8)

where RL is the third-order response function for the Liouville pathway. The total electric field will be the superposition of all the Liouville pathways,

Etot=LEL(t).
(9)

For the superposition of Eq. (9) to be non-canceling, certain symmetries between the pathways must be broken. In general, this requires one or more of the following inequalities: Γ10Γ21, ω10ω21, and/or 2μ10μ21. Our simulations use the last inequality, which is important in two-level systems (μ21=0) and in systems where state-filling dominates the non-linear response, such as in semiconductor excitons. The exact ratio between μ10 and μ21 affects the absolute amplitude of the field but does not affect the multidimensional line shape. Importantly, the dipole inequality does not break the symmetry of double quantum coherence pathways (time-orderings II and IV), so such pathways are not present in our analysis.

In MR-CMDS, a monochromator resolves the driven output frequency from other nonlinear output frequencies, which in our case is ωm=ω1ω2+ω2=ω1. The monochromator can also enhance spectral resolution, as we show in Sec. IV A. In this simulation, the detection is emulated by transforming Etot(t) into the frequency domain, applying a narrow bandpass filter, M(ω), about ω1, and applying amplitude-scaled detection,

Stot(ω1,ω2,τ21,τ22)=|M(ωω1)Etot(ω)|2dω,
(10)

where Etot(ω) denotes the Fourier transform of Etot(t) [see Fig. 2(c)]. For M we used a rectangular function of width 0.408Δω. The arguments of Stot refer to the experimental degrees of freedom. The signal delay dependence is parameterized with the relative delays τ21 and τ22, where τnm=τnτm [see Fig. 2(b)]. Table S1 of the supplementary material summarizes the arguments for each Liouville pathway. Figure 2(f) shows the 2D (ω1,ω2)Stot spectrum resulting from the pulse delay times represented in Fig. 2(b).

Inhomogeneity is isolated in CMDS through both spectral signatures, such as line-narrowing,10,51,68–70 and temporal signatures, such as photon echoes.71,72 We simulate the effects of static inhomogeneous broadening by convolving the homogeneous response with a Gaussian distribution function. Further details of the convolution are in  Appendix B. Dynamic broadening effects such as spectral diffusion are beyond the scope of this work.

A matrix representation of differential equations of the type in Eq. (7) was numerically integrated for parallel computation of Liouville elements (see the supplementary material for details).73,74 The lower bound of integration was 2Δt before the first pulse, and the upper bound was 5Γ101 after the last pulse, with step sizes much shorter than the pulse durations. Integration was performed in the FID rotating frame; the time steps were chosen so that both the system-pulse difference frequencies and the pulse envelope were well-sampled.

The following simulations explore the four-dimensional (ω1,ω2,τ21,τ22) variable space. Both frequencies are scanned about the resonance, and both delays are scanned about pulse overlap. We explored the role of the sample dephasing rate by calculating the signal for systems with dephasing rates such that Γ10Δt=0.5,1, and 2. Inhomogeneous broadening used a spectral FWHM, Δinhom, that satisfied Δinhom/Δω=0,0.5,1, and 2 for the three dephasing rates. For all these dimensions, both ρ3(t;ω1,ω2,τ21,τ22) and Stot(ω1,ω2,τ21,τ22) are recorded for each unique Liouville pathway. Our simulations were done using the open-source SciPy library.75 

We now present portions of our simulated data that highlight the dependence of the spectral line shapes and transients on excitation pulse width, dephasing rate, detuning from resonance, pulse delay times, and inhomogeneous broadening.

It is illustrative to first consider the evolution of single coherences, ρ0xρ1, under various excitation conditions. Figure 3 shows the temporal evolution of ρ1 with various dephasing rates under Gaussian excitation. The value of ρ1 differs only by phase factors between various Liouville pathways [this can be verified by inspection of Eq. (5) under the various conditions in Table S1 of the supplementary material], so the profiles in Fig. 3 apply for the first interaction of any pathway. The pulse frequency was detuned from resonance so that frequency changes could be visualized by the color bar, but the detuning was kept slight so that it did not appreciably change the dimensionless product, ΔtΓf+iκfΩfxΓ10Δt. In this case, the evolution demonstrates the maximum impulsive character the transient can achieve. The instantaneous frequency, dφ/dt, is defined as

dφdt=ddttan1Im ρ1(t)Re ρ1(t).
(11)

The cases of Γ10Δt=0() agree with the impulsive (driven) expressions derived in  Appendix A. For Γ10Δt=0, the signal rises as the integral of the pulse and has instantaneous frequency close to that of the pulse [Eq. (A3)], but as the pulse vanishes, the signal adopts the natural system frequency and decay rate [Eq. (A2)]. For Γ10Δt=, the signal follows the amplitude and frequency of the pulse for all times [the driven limit, Eq. (A1)].

FIG. 3.

The relative importance of FID and driven responses for a single quantum coherence as a function of the relative dephasing rate (values of Γ10Δt are shown in insets). The black line shows the coherence amplitude profile, while the shaded color indicates the instantaneous frequency (see colorbar). For all cases, the pulsed excitation field (gray line, shown as electric field amplitude) is slightly detuned (relative detuning, Ωfx/Δω=0.1).

FIG. 3.

The relative importance of FID and driven responses for a single quantum coherence as a function of the relative dephasing rate (values of Γ10Δt are shown in insets). The black line shows the coherence amplitude profile, while the shaded color indicates the instantaneous frequency (see colorbar). For all cases, the pulsed excitation field (gray line, shown as electric field amplitude) is slightly detuned (relative detuning, Ωfx/Δω=0.1).

Close modal

The other three cases show a smooth interpolation between limits. As Γ10Δt increases from the impulsive limit, the coherence within the pulse region conforms less to a pulse integral profile and more to a pulse envelope profile. Accordingly, the FID component after the pulse becomes less prominent, and the instantaneous frequency pins to the driving frequency more strongly through the course of evolution. The trends can be understood by considering the differential form of evolution [Eq. (4)] and the time-dependent balance of optical coupling and system relaxation. We note that our choices of Γ10Δt=2.0,1.0, and 0.5 give coherences that have mainly driven, roughly equal driven and FID parts, and mainly FID components, respectively. The FID character is difficult to isolate when Γ10Δt=2.0.

Figure 4(a) shows the temporal evolution of ρ1 at several values of Ω1x/Δω with Γ10Δt=1 (supplementary material Fig. S3 shows the Fourier domain representation of Fig. 4(a)). As detuning increases, total amplitude decreases, FID character vanishes, and ρ1 assumes a more driven character, as expected. During the excitation, ρ1 maintains a phase relationship with the input field [as seen by the instantaneous frequency in Fig. 4(a)]. The coherence will persist beyond the pulse duration only if the pulse transfers energy into the system; FID evolution equates to absorption. The FID is therefore sensitive to the absorptive (imaginary) line shape of a transition, while the driven response is the composite of both absorptive and dispersive components. If the experiment isolates the latent FID response, there is consequently a narrower spectral response. This spectral narrowing can be seen in Fig. 4(a) by comparing the coherence amplitudes at t = 0 (driven) and at t/Δt=2 (FID); amplitudes for all Ωfx/Δω values shown are comparable at t = 0, but the lack of FID formation for Ωfx/Δω=±2 manifests as a visibly disproportionate amplitude decay (supplementary materials’ Fig. S4 shows explicit plots of ρ1(Ωfx/Δω) at discrete t/Δt values). Many ultrafast spectroscopies take advantage of the latent FID to suppress non-resonant background, improving signal to noise.42,47,59,76

FIG. 4.

Pulsed excitation of a single quantum coherence and its dependence on the pulse detuning. In all cases, relative dephasing is kept at Γ10Δt=1. (a) The relative importance of FID and driven responses for a single quantum coherence as a function of the detuning (values of relative detuning, Ωfx/Δω, are shown in insets). The color indicates the instantaneous frequency (scale bar on right), while the black line shows the amplitude profile. The gray line is the electric field amplitude. (b) The time-integrated coherence amplitude as a function of the detuning. The integrated amplitude is collected both with (teal) and without (magenta) a tracking monochromator that isolates the driven frequency components. For comparison, Green’s function of the single quantum coherence is also shown (amplitude is black, hashed; imaginary is black, solid). In all plots, the gray line is the electric field amplitude.

FIG. 4.

Pulsed excitation of a single quantum coherence and its dependence on the pulse detuning. In all cases, relative dephasing is kept at Γ10Δt=1. (a) The relative importance of FID and driven responses for a single quantum coherence as a function of the detuning (values of relative detuning, Ωfx/Δω, are shown in insets). The color indicates the instantaneous frequency (scale bar on right), while the black line shows the amplitude profile. The gray line is the electric field amplitude. (b) The time-integrated coherence amplitude as a function of the detuning. The integrated amplitude is collected both with (teal) and without (magenta) a tracking monochromator that isolates the driven frequency components. For comparison, Green’s function of the single quantum coherence is also shown (amplitude is black, hashed; imaginary is black, solid). In all plots, the gray line is the electric field amplitude.

Close modal

In driven experiments, the output frequency and line shape are fully constrained by the excitation beams. In such experiments, there is no additional information to be resolved in the output spectrum. The situation changes in the mixed domain, where Etot contains the FID signal that lasts longer than the pulse duration. Figure 4(a) provides insight on how frequency-resolved detection of coherent output can enhance resolution when pulses are spectrally broad. Without frequency-resolved detection, mixed-domain resonance enhancement occurs in two ways: (1) the peak amplitude increases, and (2) the coherence duration increases due to the FID transient. Frequency-resolved detection can further discriminate against detuning by requiring that the driving frequency agrees with latent FID. The implications of discrimination are most easily seen in Fig. 4(a) with Ω1x/Δω=±1, where the system frequency moves from the driving frequency to the FID frequency. When the excitation pulse frequency is scanned, the resonance will be more sensitive to detuning by isolating the driven frequency (tracking the monochromator with the excitation source).

The functional form of the measured line shape can be deduced by considering the frequency domain form of Eq. (5) (assume ρi=1 and τx=0),

ρ̃f(ω)=iλfμf22πFCxωκfΩfxΓf+iω,
(12)

where FCxω is the Fourier transform of cx,

FCxω=12πe(Δtω)24ln2.
(13)

For squared-law detection of ρf, the importance of the tracking monochromator is highlighted by two limits of Eq. (12) as follows:

  • When the transient is not frequency resolved, sig|ρ̃f(ω)|2dω and the measured line shape will be the convolution of the pulse envelope and the intrinsic (Green’s function) response [Fig. 4(b), magenta].

  • When the driven frequency is isolated, sig|ρ̃f(κfΩfx)|2 and the measured line shape will give the un-broadened Green’s function [Fig. 4(b), teal].

Monochromatic detection can remove broadening effects due to the pulse bandwidth. For large Γ10Δt values, FID evolution is negligible at all Ωfx/Δω values and the monochromator is not useful. Figure 4(b) shows the various detection methods for the relative dephasing rate of Γ10Δt=1.

We now consider the multidimensional response of a single Liouville pathway involving three pulse interactions. In a multi-pulse experiment, ρ1 acts as a source term for ρ2 (and subsequent excitations). The spectral and temporal features of ρ1 that are transferred to ρ2 depend on when the subsequent pulse arrives. Time-gating later in ρ1 evolution will produce responses with FID behavior, while time-gating ρ1 in the presence of the initial pulse will produce driven responses. An analogous relationship holds for ρ3 with its source term ρ2. As discussed above, the signal that time-gates FID evolution gives narrower spectra than the driven-gated signal. As a result, the spectra of even single Liouville pathways will change based on pulse delays.

The final coherence will also be frequency-gated by the monochromator. The monochromator isolates the signal at the fully driven frequency ωout=ω1. The monochromator will induce line-narrowing to the extent that FID takes place. It effectively enforces a frequency constraint that acts as an additional resonance condition, ωout=ω1. The driven frequency will be ω1 if E1 is the last pulse interaction (time-orderings V and VI), and the monochromator tracks the coherence frequency effectively. If E1 is not the last interaction, the output frequency may not be equal to the driven frequency, and the monochromator plays a more complex role.

We demonstrate this delay dependence using the multidimensional response of the Iγ Liouville pathway as an example (see Fig. 1). Figure 5 shows the resulting 2D delay profile of pathway Iγ signals for Γ10Δt=1 (left) and the corresponding ω1,ω2 2D spectra at several pulse delay values (right). The spectral changes result from changes in the relative importance of driven and FID components. The prominence of the FID signal can change the resonance conditions; Table I summarizes the changing resonance conditions for each of the four delay coordinates studied. Since E1 is not the last pulse in pathway Iγ, the tracking monochromator must also be considered. Supplementary materials Fig. S5 shows a simulation of Fig. 5 without monochromator frequency filtering (M(ωω1)=1 in Eq. 10).

FIG. 5.

Changes to the 2D frequency response of a single Liouville pathway (Iγ) at different delay values. The normalized dephasing rate is Γ10Δt=1. Left: The 2D delay response of pathway Iγ at triple resonance. Right: The 2D frequency response of pathway Iγ at different delay values. The delays at which the 2D frequency plots are collected are indicated on the delay plot; compare 2D spectrum frame color with dot color on the 2D delay plot.

FIG. 5.

Changes to the 2D frequency response of a single Liouville pathway (Iγ) at different delay values. The normalized dephasing rate is Γ10Δt=1. Left: The 2D delay response of pathway Iγ at triple resonance. Right: The 2D frequency response of pathway Iγ at different delay values. The delays at which the 2D frequency plots are collected are indicated on the delay plot; compare 2D spectrum frame color with dot color on the 2D delay plot.

Close modal
TABLE I.

Conditions for peak intensity at different pulse delays for pathway Iγ.

DelayApproximate resonance conditions
τ21/Δtτ22/Δtρ01ρ1ρ12ρ2ρ22ρ3ρ3 detection at ωm=ω1
ω1=ω10 ω1=ω2 ω1=ω10 … 
−2.4 ω1=ω10 ω1=ω2 ω2=ω10 ω1=ω2 
2.4 ω1=ω10 ω2=ω10 … ω1=ω10 
2.4 −2.4 ω1=ω10 ω2=ω10 ω2=ω10 ω1=ω2 
DelayApproximate resonance conditions
τ21/Δtτ22/Δtρ01ρ1ρ12ρ2ρ22ρ3ρ3 detection at ωm=ω1
ω1=ω10 ω1=ω2 ω1=ω10 … 
−2.4 ω1=ω10 ω1=ω2 ω2=ω10 ω1=ω2 
2.4 ω1=ω10 ω2=ω10 … ω1=ω10 
2.4 −2.4 ω1=ω10 ω2=ω10 ω2=ω10 ω1=ω2 

When the pulses are all overlapped (τ21=τ22=0, lower right, orange), all transitions in the Liouville pathway are simultaneously driven by the incident fields. This spectrum strongly resembles the driven limit spectrum. For this time-ordering, the first, second, and third density matrix elements have driven resonance conditions of ω1=ω10, ω1ω2=0, and ω1ω2+ω2=ω10, respectively. The second resonance condition causes elongation along the diagonal, and since ω2=ω2, the first and third resonance conditions are identical, effectively making ω1 doubly resonant at ω10 and resulting in the vertical elongation along ω1=ω10.

The other three spectra in Fig. 5 separate the pulse sequence over time so that not all interactions are driven. At τ21=0, τ22=2.4Δt (lower left, pink), the first two resonances remain the same as at pulse overlap (orange), but the last resonance is different. The final pulse, E2, is latent and probes ρ2 during its FID evolution after the memory of the driven frequency is lost. There are two important consequences. First, the third driven resonance condition is now approximated by ω2=ω10, which makes ω1 only singly resonant at ω1=ω10. Second, the driven portion of the signal frequency is determined only by the latent pulse: ωout=ω2. Since our monochromator gates ω1, we have the detection-induced correlation ω1=ω2. The net result is double resonance along ω1=ω2, and the vertical elongation of pulse overlap is strongly attenuated.

At τ21=2.4Δt,τ22=0 (upper right, purple), the first pulse E1 precedes the latter two, which makes the two resonance conditions for the input fields ω1=ω10 and ω2=ω10. The signal depends on the FID conversion of ρ1, which gives vertical elongation at ω1=ω10. Furthermore, ρ1 has no memory of ω1 when E2 interacts, which has two important implications. First, this means the second resonance condition is ω1=ω2 and the associated diagonal elongation is now absent. Second, the final output polarization frequency content is no longer functional of ω1. Coupled with the fact that E2 and E2 are coincident so that the final coherence can be approximated as driven by these two, we can approximate the final frequency as ωout=ω10ω2+ω2=ω10. Surprisingly, the frequency content of the output is strongly independent of all pulse frequencies. The monochromator narrows the ω1=ω10 resonance. The ω1=ω10 resonance condition now depends on the monochromator slit width, the FID propagation of ρ1, and the spectral bandwidth of ρ3; its spectral width is not easily related to material parameters. This resonance demonstrates the importance of the detection scheme for experiments and how the optimal detection can change depending on the pulse delay time.

Finally, when all pulses are well-separated (τ21=τ22=2.4Δt, upper left, cyan), each resonance condition is independent and both E1 and E2 require FID buildup to produce the final output. The resulting line shape is narrow in all directions. Again, the emitted frequency does not depend on ω1, yet the monochromator resolves the final coherence at frequency ω1. Since the driven part of the final interaction comes from E2, and since the monochromator track ω1, the output signal will increase when ω1=ω2. As a result, the line shape acquires a diagonal character.

The changes in the line shape seen in Fig. 5 have significant ramifications for the interpretations and strategies of MR-CMDS in the mixed domain. Time-gating has been used to isolate the 2D spectra of a certain time-ordering,13,41,60 but here we show that time-gating itself causes significant line shape changes to the isolated pathways. The phenomenon of time-gating can cause frequency and delay axes to become a functional of each other in unexpected ways.

In Sec. IV B, we showed how a single pathway’s spectra can evolve with delay due to pulse effects and time gating. In real experiments, evolution with delay is further complicated by the six time-orderings/sixteen pathways present in our three-beam experiment (see Fig. 1). Each time-ordering has different resonance conditions. When the signal is collected near pulse overlap, multiple time-orderings contribute. To identify these effects, we start by considering how strongly time-orderings are isolated at each delay coordinate.

While the general idea of using time delays to enhance certain time-ordered regions is widely applied, quantitation of this discrimination is rarely explored. Because the temporal profile of the signal is dependent on both the excitation pulse profile and the decay dynamics of the coherence itself, quantitation of pathway discrimination requires simulation.

Figure 6 shows the 2D delay space with all pathways present for ω1,ω2=ω10. It illustrates the interplay of pulse width and system decay rates on the isolation of time-ordered pathways. The color bar shows the signal amplitude. The signal is symmetric about the τ21=τ22 line because when ω1=ω2, E1 and E2 interactions are interchangeable: Stot(τ21,τ22)=Stot(τ22,τ21). The overlaid black contours represent signal “purity,” P, defined as the relative amount of the signal that comes from the dominant pathway at that delay value,

P(τ21,τ22)=maxSLτ21,τ22LSLτ21,τ22.
(14)

The dominant pathway (maxSLτ21,τ22) at given delays can be inferred by the time-ordered region defined in Fig. 2(d). The contours of purity generally run parallel to the time-ordering boundaries with the exception of time-ordered regions II and IV, which involve the double quantum coherences that have been neglected.

FIG. 6.

Comparison of the 2D delay response for different relative dephasing rates (labeled atop each column). All pulses are tuned to exact resonance. In each 2D delay plot, the signal amplitude is depicted by the colors. The black contour lines show signal purity, P [see Eq. (14)], with purity values denoted on each contour. The small plots above each 2D delay plot examine a τ22 slice of the delay response (τ21=0). The plot shows the total signal (black) as well as the component time-orderings VI (orange), V (purple), and III or I (teal).

FIG. 6.

Comparison of the 2D delay response for different relative dephasing rates (labeled atop each column). All pulses are tuned to exact resonance. In each 2D delay plot, the signal amplitude is depicted by the colors. The black contour lines show signal purity, P [see Eq. (14)], with purity values denoted on each contour. The small plots above each 2D delay plot examine a τ22 slice of the delay response (τ21=0). The plot shows the total signal (black) as well as the component time-orderings VI (orange), V (purple), and III or I (teal).

Close modal

A commonly employed metric for temporal selectivity is how definitively the pulses are ordered. This metric agrees with our simulations. The purity contours have a weak dependence on ΔtΓ10 for |τ22|/Δt<1 or |τ21|/Δt<1, where there is significant pathway overlap and a stronger dependence at larger values where the pathways are well-isolated. Because responses decay exponentially, while pulses decay as Gaussians, there always exist delays where temporal discrimination is possible. As ΔtΓ10, however, such discrimination is only achieved at vanishing signal intensities; the contour of P = 0.99 across our systems highlights this trend.

In Secs. IV A–IV C, we showed how pathway spectra and weights evolve with delay. This section ties the two concepts together by exploring the evolution of the spectral line shape over a span of τ21 delay times that include all pathways. It is a common practice to explore spectral evolution against τ21 because this delay axis shows population evolution in a manner analogous to pump-probe spectroscopies. The k2 and k2 interactions correspond to the pump, and the k1 interaction corresponds to the probe. Time-orderings V and VI are the normal pump-probe time-orderings, time-ordering III is a mixed pump-probe-pump ordering (so-called pump polarization coupling), and time-ordering I is the probe-pump ordering (so-called perturbed FID). Scanning τ21 through pulse overlap complicates interpretation of the line shape due to the changing nature and balance of the contributing time-orderings. At τ21>0, time-ordering I dominates; at τ21=0, all time-orderings contribute equally; at τ21<0, time-orderings V and VI dominate (Fig. 6). Conventional pump-probe techniques recognized these complications long ago,77,78 but the extension of these effects to MR-CMDS has not previously been done.

Figure 7 shows the MR-CMDS spectra as well as histograms of the pathway weights, while scanning τ21 through pulse overlap. The colored histogram bars and line shape contours correspond to different values of the relative dephasing rate, Γ10Δt. The contour is the half-maximum of the line shape (supplementary materials Fig. S6 shows fully colored contour plots of each 2D frequency spectrum). The dependence of the line shape amplitude on τ21 can be inferred from Fig. 6.

FIG. 7.

Evolution of the 2D frequency response as a function of τ21 (labeled inset) and the influence of the relative dephasing rate [Γ10Δt=0.5 (red), 1.0 (green), and 2.0 (blue)]. In all plots, τ22=0. To ease comparison between different dephasing rates, the colored line contours (showing the half-maximum) for all three relative dephasing rates are overlaid. The colored histograms below each 2D frequency plot show the relative weights of each time-ordering for each relative dephasing rate. Contributions from V and VI are grouped together because they have equal weights at τ22=0.

FIG. 7.

Evolution of the 2D frequency response as a function of τ21 (labeled inset) and the influence of the relative dephasing rate [Γ10Δt=0.5 (red), 1.0 (green), and 2.0 (blue)]. In all plots, τ22=0. To ease comparison between different dephasing rates, the colored line contours (showing the half-maximum) for all three relative dephasing rates are overlaid. The colored histograms below each 2D frequency plot show the relative weights of each time-ordering for each relative dephasing rate. Contributions from V and VI are grouped together because they have equal weights at τ22=0.

Close modal

The qualitative trend, as τ21 goes from positive to negative delays, is a change from diagonal/compressed line shapes to much broader resonances with no correlation (ω1 and ω2 interact with independent resonances). Such spectral changes could be misinterpreted as spectral diffusion, where the line shape changes from correlated to uncorrelated as population time increases due to system dynamics. The system dynamics included here, however, contain no structure that would allow for such diffusion. Rather, the spectral changes reflect the changes in the majority pathway contribution, starting with time-ordering I pathways, proceeding to an equal admixture of I, III, V, and VI, and finishing at an equal balance of V and VI when E1 arrives well after E2 and E2. Time-orderings I and III both exhibit a spectral correlation in ω1 and ω2 when driven, but time-orderings V and VI do not. Moreover, such a spectral correlation is forced near zero delay because the pulses time-gate the driven signals of the first two induced polarizations. The monochromator detection also plays a dynamic role because time-orderings V and VI always emit a signal at the monochromator frequency, while in time-orderings I and III, the emitted frequency is not defined by ω1, as discussed above.

When we isolate time-orderings V and VI, we can maintain the proper scaling of the FID bandwidth in the ω1 direction because our monochromator can gate the final coherence. This gating is not possible in time-orderings I and III because the final coherence frequency is determined by ω2, which is identical to ω2.

There are differences in the line shapes for the different values of the relative dephasing rate, Γ10Δt. The spectral correlation at τ21/Δt=2 decreases in strength as Γ10Δt decreases. As we illustrated in Fig. 5, this spectral correlation is a signature of the driven signal from the temporal overlap of E1 and E2; the loss of spectral correlation reflects the increased prominence of FID in the first coherence as the field-matter interactions become more impulsive. This increased prominence of FID also reflects an increase in signal strength, as shown by τ21 traces in Fig. 6. When all pulses are completely overlapped (τ21=0), each of the line shapes exhibit spectral correlation. At τ21/Δt=2, the line shape shrinks as Γ10Δt decreases, with the elongation direction changing from horizontal to vertical. The general shrinking reflects the narrowing homogeneous linewidth of the ω10 resonance. In all cases, the horizontal line shape corresponds to the homogeneous linewidth because the narrow bandpass monochromator resolves the final ω1 resonance. The change in the elongation direction is due to the resolving power of ω2. At Γ10Δt=0.5, the resonance is broader than our pulse bandwidth and is fully resolved vertically. It is narrower than the ω1 resonance because time-orderings V and VI interfere to isolate only the absorptive line shape along ω2. This narrowing, however, is unresolvable when the pulse bandwidth becomes broader than that of the resonance, which gives rise to a vertically elongated signal when Γ10Δt=0.5.

It is also common to represent data as “Wigner plots,” where one axis is delay and the other is frequency.14,15,21,57 In Fig. 8, we show five τ21,ω1 plots for varying ω2 with τ22=0. Supplementary materials’ Fig. S8 shows Wigner plots for other Γ10Δt values. The plots are the analog to the most common multidimensional experiment of transient absorption spectroscopy, where the non-linear probe spectrum is plotted as a function of the pump-probe delay. For each plot, the ω2 frequency is denoted by a vertical gray line. Each Wigner plot is scaled to its own dynamic range to emphasize the dependence on ω2. The dramatic line shape changes between positive and negative delays can be seen. This representation also highlights the asymmetric broadening of the ω1 line shape near pulse overlap when ω2 becomes non-resonant. Again, these features can resemble spectral diffusion even though our system is homogeneous.

FIG. 8.

Transient (ω1) line shapes and their dependence on ω2 frequency. The relative dephasing rate is Γ10Δt=1 and τ22=0. For each plot, the corresponding ω2 value is shown as a light gray vertical line.

FIG. 8.

Transient (ω1) line shapes and their dependence on ω2 frequency. The relative dephasing rate is Γ10Δt=1 and τ22=0. For each plot, the corresponding ω2 value is shown as a light gray vertical line.

Close modal

With the homogeneous system characterized, we can now consider the effect of inhomogeneity. For inhomogeneous systems, time-orderings III and V are enhanced because their final coherence will rephase to form a photon echo, whereas time-orderings I and VI will not. In delay space, this rephasing appears as a shift of the signal to time-ordered regions III and V that persists for all population times. Figure 9 shows the calculated spectra for relative dephasing rate Γ10Δt=1 with a frequency broadening function of width Δinhom=0.441Γ10. The inhomogeneity makes it easier to temporally isolate the rephasing pathways and harder to isolate the non-rephasing pathways, as shown by the purity contours.

FIG. 9.

2D delay response for Γ10Δt=1 with sample inhomogeneity. All pulses are tuned to exact resonance. The colors depict the signal amplitude. The black contour lines show signal purity, P [see Eq. (14)], with purity values denoted on each contour. The thick yellow line denotes the peak amplitude position that is used for 3PEPS analysis. The small plot above each 2D delay plot examines a τ22 slice of the delay response (τ21=0). The plot shows the total signal (black) as well as the component time-orderings VI (orange), V (purple), III (teal, dashed), and I (teal, solid).

FIG. 9.

2D delay response for Γ10Δt=1 with sample inhomogeneity. All pulses are tuned to exact resonance. The colors depict the signal amplitude. The black contour lines show signal purity, P [see Eq. (14)], with purity values denoted on each contour. The thick yellow line denotes the peak amplitude position that is used for 3PEPS analysis. The small plot above each 2D delay plot examines a τ22 slice of the delay response (τ21=0). The plot shows the total signal (black) as well as the component time-orderings VI (orange), V (purple), III (teal, dashed), and I (teal, solid).

Close modal

A common metric of rephasing in delay space is the 3PEPS measurement.71,79–81 In 3PEPS, one measures the signal as the first coherence time, τ, is scanned across both rephasing and non-rephasing pathways while keeping population time, T, constant. The position of the peak is measured; a peak shifted away from τ=0 reflects the rephasing ability of the system. An inhomogeneous system will emit a photon echo in the rephasing pathway, enhancing the signal in the rephasing time-ordering and creating the peak shift. In our 2D delay space, the τ trace can be defined if we assume that E2 and E2 create the population (time-orderings V and VI). The trace runs parallel to the III-V time-ordering boundary (diagonal) if τ22<0 and runs parallel to the IV-VI time-ordering boundary (horizontal) if τ22>0, and both intersect at τ22=0; the τ21 value at this intersection is T (supplementary materials Fig. S9 illustrates how 3PEPS shifts are measured from a 2D delay plot). In our 2D delay plots (Figs. 6 and 9), the peak shift is seen as the diagonal displacement of the signal peak from the τ21=0 vertical line (supplementary materials Fig. S10 shows the 3PEPS measurements of all 12 combinations of Γ10Δt and Δinhom for every population delay surveyed). Figure 9 highlights the peak shift profile as a function of population time with the yellow trace; it is easily verified that our static inhomogeneous system exhibits a non-zero peak shift value for all population times.

The unanticipated feature of the 3PEPS analysis is the dependence on T. Even though our inhomogeneity is static, the peak shift is maximal at T = 0 and dissipates as T increases, mimicking spectral diffusion. This dynamic arises from signal overlap with time-ordering III, which uses E2 and E1 as the first two interactions. At T = 0, the τ trace gives two ways to make a rephasing pathway (time-orderings III and V) and only one way to make a non-rephasing pathway (time-ordering VI). This pathway asymmetry shifts the signal away from τ=0 into the rephasing direction. At large T (large τ21), time-ordering III is not viable and pathway asymmetry disappears. Peak shifts imply inhomogeneity only when time-orderings V and III are minimally contaminated by each other, i.e., at population times that exceed pulse overlap. This fact is easily illustrated by the dynamics of the homogeneous system (Fig. 6); even though the homogeneous system cannot rephase, there is a non-zero peak shift near T = 0. The contamination of the 3PEPS measurement at pulse overlap is well known and is described in some studies,23,72 but the dependence of pulse and system properties on the distortion has not been investigated previously. Peak-shifting due to pulse overlap is less important when ω1ω2 because time-ordering III is decoupled by detuning.

In frequency space, spectral elongation along the diagonal is the signature of inhomogeneous broadening. Figure 10 shows the line shape changes of a Gaussian inhomogeneous distribution (supplementary materials Fig. S7 shows all contours for Fig. 10). All systems are broadened by a distribution proportional to their dephasing bandwidth. As expected, the sequence again shows a gradual broadening along the ω1 axis, with a strong spectral correlation at early delays (τ21>0) for the more driven signals. The anti-diagonal width at early delays (e.g., Fig. 10, τ21=2.0Δt) again depends on the pulse bandwidth and the monochromator slit width. At delay values that isolate time-orderings V and VI, however, the line shapes retain diagonal character, showing the characteristic balance of homogeneous and inhomogeneous widths.

FIG. 10.

The same as Fig. 7, but each system has inhomogeneity (Δinhom=0.441Γ10). Relative dephasing rates are Γ10Δt=0.5 (red), 1.0 (green), and 2.0 (blue). In all plots, τ22=0. To ease comparison between different dephasing rates, the colored line shapes of all three systems are overlaid. Each 2D plot shows a single representative contour (half-maximum) for each Γ10Δt value. The colored histograms below each 2D frequency plot show the relative weights of each time-ordering for each 2D frequency plot. In contrast to Fig. 7, inhomogeneity makes the relative contributions of time-orderings V and VI unequal.

FIG. 10.

The same as Fig. 7, but each system has inhomogeneity (Δinhom=0.441Γ10). Relative dephasing rates are Γ10Δt=0.5 (red), 1.0 (green), and 2.0 (blue). In all plots, τ22=0. To ease comparison between different dephasing rates, the colored line shapes of all three systems are overlaid. Each 2D plot shows a single representative contour (half-maximum) for each Γ10Δt value. The colored histograms below each 2D frequency plot show the relative weights of each time-ordering for each 2D frequency plot. In contrast to Fig. 7, inhomogeneity makes the relative contributions of time-orderings V and VI unequal.

Close modal

Our chosen values of the relative dephasing time, Γ10Δt, describe experiments where neither the impulsive nor the driven limit unilaterally applies. We have illustrated that in this intermediate regime, the multidimensional spectra contain attributes of both limits and that it is possible to judge when these attributes apply. In our three-pulse experiment, the second and third pulses time-gate coherences and populations produced by the previous pulse(s), and the monochromator frequency-gates the final coherence. Time-gating isolates different properties of the coherences and populations. Consequently, spectra evolve against delay. For any delay coordinate, one can develop qualitative line shape expectations by considering the following three principles:

  1. When time-gating during the pulse, the system pins to the driving frequency with a buildup efficiency determined by resonance.

  2. When time-gating after the pulse, the FID dominates the system response.

  3. The emitted signal field contains both FID and driven components; the ωout=ω1 component is isolated by the tracking monochromator.

Figure 3 illustrates principles 1 and 2 and Fig. 4 illustrates principle 2 and 3. Figure 5 provides a detailed example of the relationship between these principles and the multidimensional line shape changes for different delay times.

The principles presented above apply to a single pathway. For rapidly dephasing systems, it is difficult to achieve complete pathway discrimination, as shown in Fig. 6. In such situations, the interference between pathways must be considered to predict the line shape. The relative weight of each pathway to the interference can be approximated by the extent of pulse overlap. The pathway weights exchange when scanning across pulse overlap, which creates the dramatic line shape changes observed in Figs. 7 and 10.

We have shown that the driven limit misses details of the line shape if Γ10Δt1, but we have also reasoned that in certain conditions, the driven limit can approximate the response well (see principle 1). Here we examine the line shape at delay values that demonstrate this agreement. Figure 11 compares the results of our numerical simulation (third column) with the driven limit expressions for populations where Γ11Δt=0 (first column) or 1 (second column). The top and bottom rows compare the line shapes when τ22,τ21=(0,0) and (0,4Δt), respectively. The third column demonstrates the agreement between the driven limit approximations with the simulation by comparing the diagonal and anti-diagonal cross sections of the 2D spectra.

FIG. 11.

Comparing approximate expressions of the 2D frequency response with the directly integrated response. Γ10Δt=1. The top row compares the 2D response of all time-orderings (τ21=0) and the bottom row compares the response of time-orderings V and VI (τ21=4Δt). First column: The driven limit response. Note the narrow diagonal resonance for τ21=0. Second column: The same as the first column, but with ad hoc substitution Γ11=Δt. Third column: The directly integrated response.

FIG. 11.

Comparing approximate expressions of the 2D frequency response with the directly integrated response. Γ10Δt=1. The top row compares the 2D response of all time-orderings (τ21=0) and the bottom row compares the response of time-orderings V and VI (τ21=4Δt). First column: The driven limit response. Note the narrow diagonal resonance for τ21=0. Second column: The same as the first column, but with ad hoc substitution Γ11=Δt. Third column: The directly integrated response.

Close modal

Note the very sharp diagonal feature that appears for (τ21,τ22)=(0,0) and Γ11=0; this is due to population resonance in time-orderings I and III. This expression is inaccurate: the narrow resonance is only observed when pulse durations are much longer than the coherence time. A comparison of picosecond and femtosecond studies of quantum dot exciton line shapes (Yurs et al.82 and Kohler et al.,15 respectively) demonstrates this difference well. The driven equation fails to reproduce our numerical simulations here because resonant excitation of the population is impulsive; the experiment time-gates only the rise time of the population, yet driven theory predicts the resonance to be vanishingly narrow (Γ11=0). In light of this, one can approximate this time-gating effect by substituting population lifetime with the pulse duration (Γ11Δt=1), which gives good agreement with the numerical simulation (third column).

When τ22=0 and τ21<Δt, signals can also be approximated by the driven signal (Fig. 11, bottom row). Only time-orderings V and VI are relevant. The intermediate population resonance is still impulsive, but it depends on ω2ω2, which is not explored in our 2D frequency space.

We have shown that pulse effects mimic the qualitative signatures of inhomogeneity. Here we address how one can extract true system inhomogeneity in light of these effects. We focus on two ubiquitous metrics of inhomogeneity: 3PEPS for the time domain and ellipticity83 for the frequency domain.52,84 In the driven (impulsive) limit, ellipticity (3PEPS) corresponds to the frequency correlation function and uniquely extracts the inhomogeneity of the models presented here. In their respective limits, the metrics give values proportional to the inhomogeneity.

Figure 12 shows the results of this characterization for all Δinhom and Γ10Δt values explored in this work. We study how the correlations between the two metrics depend on the relative dephasing rate, Γ10Δt, the absolute inhomogeneity, Δinhom/Δω, the relative inhomogeneity Δinhom/Γ10, and the population time delay (supplementary materials Figs. S10–S12 show simulations for each value of the 3PEPS and ellipticity data in Fig. 12). The top row shows the correlations of the Δ3PEPS/Δt 3PEPS metric that represents the normalized coherence delay time required to reach the peak intensity. The upper right graph shows the correlations for a population time delay of T=4Δt that isolates the V and VI time-orderings. For this time delay, the Δ3PEPS/Δt metric works well for all dephasing times of Γ10Δt when the relative inhomogeneity is Δinhom/Δω<<1. It becomes independent of Δinhom/Δω when Δinhom/Δω>1. This saturation results because the frequency bandwidth of the excitation pulses becomes smaller than the inhomogeneous width and only a portion of the inhomogeneous ensemble contributes to the 3PEPS experiment.71 The corresponding graph for T = 0 shows that a large peak shift occurs, even without inhomogeneity. In this case, the peak shift depends on pathway overlap, as discussed in Sec. IV E.

FIG. 12.

Temporal (3PEPS) and spectral (ellipticity) metrics of correlation and their relation to the true system inhomogeneity. The left column plots the relationship at pulse overlap (T = 0) and the right column plots the relationship at a delay where driven correlations are removed (T=4Δt). For the ellipticity measurements, τ22=0. In each case, the two metrics are plotted directly against system inhomogeneity (top and middle row) and against each other (bottom row). Colored lines guide the eyes for systems with equal relative dephasing rates (Γ10Δt, see upper legend), while the area of the data point marker indicates the relative inhomogeneity (Δinhom/Γ10, see lower legend). Gray lines indicate contours of constant relative inhomogeneity (scatter points with the same area are connected).

FIG. 12.

Temporal (3PEPS) and spectral (ellipticity) metrics of correlation and their relation to the true system inhomogeneity. The left column plots the relationship at pulse overlap (T = 0) and the right column plots the relationship at a delay where driven correlations are removed (T=4Δt). For the ellipticity measurements, τ22=0. In each case, the two metrics are plotted directly against system inhomogeneity (top and middle row) and against each other (bottom row). Colored lines guide the eyes for systems with equal relative dephasing rates (Γ10Δt, see upper legend), while the area of the data point marker indicates the relative inhomogeneity (Δinhom/Γ10, see lower legend). Gray lines indicate contours of constant relative inhomogeneity (scatter points with the same area are connected).

Close modal

The middle row in Fig. 12 shows the ellipticity dependence on the relative dephasing rate and inhomogeneity assuming that the measurement is performed when the first two pulses are temporally overlapped (τ22=0). For a T=4Δt population time, the ellipticity is proportional to the inhomogeneity until Δinhom/Δω<<1, where the excitation bandwidth is wide compared with the inhomogeneity. Unlike 3PEPS, saturation is not observed because the pulse bandwidth does not limit the frequency range scanned. The 3PEPS and ellipticity metrics are therefore complementary since 3PEPS works well for Δinhom/Δω<<1 and ellipticity works well for Δinhom/Δω>>1. When all pulses are temporally overlapped at T = 0, the ellipticity is only weakly dependent on the inhomogeneity and dephasing rate. The ellipticity is instead dominated by the dependence on the excitation pulse frequency differences of time-orderings I and III that become important at pulse overlap.

It is clear from the previous discussion that both metrics depend on the dephasing and inhomogeneity. The dephasing can be measured independently in the frequency or time domain, depending upon whether the dephasing is very fast or slow, respectively. In the mixed frequency/time domain, the measurement of the dephasing becomes more difficult. One strategy to address this challenge is to use both the 3PEPS and ellipticity metrics. The bottom row in Fig. 12 plots 3PEPS against ellipticity to show how the relationship between the metrics changes for different amounts of dephasing and inhomogeneity. The anti-diagonal contours of constant relative inhomogeneity show that these metrics are complementary and can serve to extract the system correlation parameters.

Importantly, the metrics are uniquely mapped both in the presence and absence of pulse-induced effects (demonstrated by T = 0 and T=4Δt, respectively). The combined metrics can be used to determine correlation at T = 0, but the correlation-inducing pulse effects give a mapping significantly different than that at T=4Δt. At T = 0, 3PEPS is almost nonresponsive to inhomogeneity; instead, it is an almost independent characterization of the pure dephasing. In fact, the T = 0 trace is equivalent to the original photon echo traces used to resolve pure dephasing rates.85 Both metrics are offset due to the pulse overlap effects. Accordingly, the region to the left of the homogeneous contour is non-physical because it represents observed correlations that are less than those given by pulse overlap effects. If the metrics are measured as a function of T, the mapping gradually changes from the left figure to the right figure in accordance with the pulse overlap. Both metrics will show a decrease, even with static inhomogeneity. If a system has spectral diffusion, the mapping at late times will disagree with the mapping at early times; both ellipticity and 3PEPS will be smaller at later times than predicted by the change in mappings alone.

This study provides a framework to describe and disentangle the influence of the excitation pulses in mixed-domain ultrafast spectroscopy. We analyzed the features of mixed-domain spectroscopy through detailed simulations of MR-CMDS signals. When pulse durations are similar to coherence times, resolution is compromised by time-bandwidth uncertainty and the complex mixture of driven and FID responses. The dimensionless quantity ΔtΓf+iκfΩfx captures the balance of driven and FID character in a single field-matter interaction. In the nonlinear experiment, with multiple field-matter interactions, this balance is also controlled by pulse delays and frequency-resolved detection. Our analysis shows how these effects can be intuitive.

The dynamic nature of pulse effects can lead to misleading changes to spectra when delays are changed. When delays separate pulses, the spectral line shapes of individual pathways qualitatively change because the delays isolate FID contributions and de-emphasize the driven response. When delays are scanned across pulse overlap, the weights of individual pathways change, further changing the line shapes. In a real system, these changes would all be present in addition to actual dynamics and spectral changes of the material.

Finally, we find that in either the frequency or time domain, pulse effects mimic signatures of ultrafast inhomogeneity. Even homogeneous systems take on these signatures. Driven character gives rise to pathway overlap peak shifting in the 2D delay response, which artificially produces rephasing near pulse overlap. Driven character also produces resonances that depend on ω1ω2 near pulse overlap. Determination of the homogeneous and inhomogeneous broadening at ultrashort times is only possible by performing correlation analysis in both the frequency and time domains.

See supplementary material for the simulation results that are a sparse sampling of the total data simulated. The full dataset, as well as the software used to simulate, is available for download at http://dx.doi.org/10.17605/OSF.IO/EJ2XE. Additional derivations and figures can be found in the associated PDF. See Fig. S3 for a Fourier domain representation of Fig. 4(a). See Fig. S4 for explicit plots of ρ1(Ωfx/Δω) at discrete t/Δt values. See Fig. S5 for a representation of Fig. 5 simulated without monochromator frequency filtering [M(ωω1)=1 for Eq. (10)]. Figure S6 shows fully colored contour plots of each 2D frequency spectrum. See Fig. S8 for Wigner plots for all Γ10Δt values. See Fig. S9 for an illustration of how 3PEPS shifts are measured from a 2D delay plot. Figure S10 shows the 3PEPS measurements of all 12 combinations of Γ10Δt and Δinhom, for every population delay surveyed. As in Fig. 7, Fig. 10 shows only the contours at the half-maximum amplitude. See Fig. S7 for all contours. The simulations for each value of the 3PEPS and ellipticity data in Fig. 12 appear in Figs. S10–S12.

This work was supported by the U.S. Department of Energy, Office of Science, Basic Energy Sciences, under Award No. DE-FG02-09ER46664.

The changes in the spectral line shapes described in this work are best understood by examining the driven/continuous wave (CW) and impulsive limits of Eqs. (5) and (7). The driven limit is achieved when pulse durations are much longer than the response function dynamics: Δt|Γf+iκfΩfx|>>1. In this limit, the system will adopt a steady state over excitation: dρ/dt0. Neglecting phase factors, the driven solution to Eq. (5) will be

ρ̃f(t)=λfμf2cx(tτx)eiκfΩfxtκfΩfxρ̃i(t).
(A1)

The frequency and temporal envelopes of the excitation pulse control the coherence time evolution, and the relative amplitude and phase of the coherence are directly related to detuning from resonance.

The impulsive limit is achieved when the excitation pulses are much shorter than response function dynamics: Δt|Γf+iκfΩfx|<<1. The full description of the temporal evolution has two separate expressions: one for times when the pulse is interacting with the system and one for times after pulse interaction. Both expressions are important when describing CMDS experiments.

For times after the pulse interaction, tτx+Δt, the field-matter coupling is negligible. The evolution for these times, on resonance, is given by

ρ̃f(t)=iλfμf2ρ̃i(τx)cx(u)dueΓf(tτx).
(A2)

This is classic free induction decay (FID) evolution: the system evolves at its natural frequency and decays at rate Γf. It is important to note that while this expression is explicitly derived from the impulsive limit, FID behavior is not exclusive to impulsive excitation, as we have defined it. A latent FID will form if the pulse vanishes at a fast rate relative to the system dynamics.

For evaluating times near pulse excitation, tτx+Δt, we implement a Taylor expansion in the response function about zero: e(Γf+iκfΩfx)u=1(Γf+iκfΩfx)u+. Our impulsive criterion requires that a low order expansion suffices; it is instructive to consider the result of the first order expansion of Eq. (5),

ρ̃f(t)=iλfμf2eiκfωxτxeiκfΩfxtρ̃i(τx)×1(Γf+iκfΩfx)(tτx)tτxcx(u)du+(Γf+iκfΩfx)tτxcx(u)udu.
(A3)

During this time, ρ̃f builds up roughly according to the integration of the pulse envelope. The build-up is integrated because the pulse transfers energy before appreciable dephasing or detuning occurs. Contrary to the expectation of impulsive evolution, the evolution of ρ̃f is explicitly affected by the pulse frequency, and the temporal profile evolves according to the pulse.

It is important to recognize that the impulsive limit is defined not only by having slow relaxation relative to the pulse duration but also by small detuning relative to the pulse bandwidth (as is stated in the inequality). As detuning increases, the higher orders of the response function Taylor expansion will be needed to describe the rise time, and the driven limit of Eq. (A1) will become valid. The details of this build-up time can often be neglected in impulsive approximations because build-up contributions are often negligible in the analysis; the period over which the initial excitation occurs is small in comparison to the free evolution of the system. The build-up behavior can be emphasized by the measurement, which makes Eq. (A3) important.

We now consider full Liouville pathways in the impulsive and driven limits of Eq. (7). For the driven limit, Eq. (7) can be reduced to

EL(t)=18λ1λ2λ3μ1μ2μ3μ4ei(κ1ωxτx+κ2ωyτy+κ3ωzτz)×ei(κ3ωz+κ2ωy+κ1ωx)tcz(tτz)cy(tτy)cx(tτx)×1κ1Ω1xiΓ11κ1Ω1x+κ2Ω2yiΓ2×1κ1Ω1x+κ2Ω2y+κ3Ω3ziΓ3.
(A4)

It is important to note that the signal depends on the multiplication of all the fields; pathway discrimination based on pulse time-ordering is not achievable because polarizations exist only when all pulses are overlapped. This limit is the basis for frequency-domain techniques. Frequency axes, however, are not independent because the system is forced to the laser frequency and influences the resonance criterion for subsequent excitations. As an example, observe that the first two resonant terms in Eq. (A4) are maximized when ωx=|ω1| and ωy=|ω2|. If ωx is detuned by some value ε, however, the occurrence of the second resonance shifts to ωy=|ω2|+ε, effectively compensating for the ωx detuning. This shifting of the resonance results in 2D line shape correlations.

If the pulses do not temporally overlap (τx+Δtτy+Δtτz+Δtt), then the impulsive solution to the full Liouville pathway of Eq. (7) is

EL(t)=i8λ1λ2λ3μ1μ2μ3μ4ei(ω1+ω2+ω3)t×cx(w)dwcy(v)dvcz(u)du×eΓ1(τyτx)eΓ2(τzτy)eΓ(tτz).
(A5)

Pathway discrimination is demonstrated here because the signal is sensitive to the time-ordering of the pulses. This limit is suited for delay-scanning techniques. The emitted signal frequency is determined by the system and can be resolved by scanning a monochromator.

The driven and impulsive limits can qualitatively describe our simulated signals at certain frequency and delay combinations. Of the three expressions, the FID limit most resembles the signal when pulses are near resonance and well-separated in time (so that build-up behavior is negligible). The build-up limit is approximated well when pulses are near-resonant and arrive together (so that build-up behavior is emphasized). The driven limit holds for large detuning, regardless of delay.

Here we describe how to transform the data of a single reference oscillator signal to that of an inhomogeneous distribution. The oscillators in the distribution are allowed to have arbitrary energies for their states, which will cause frequency shifts in the resonances. To show this, we start with a modified, but equivalent, form of Eq. (4) as follows:

dρ̃fdt=Γfρ̃f+i2λfμfcx(tτx)×eiκfkz+ωxτxeiκfωx|ωf|tρ̃i(t).
(B1)

We consider two oscillators with transition frequencies ωf and ωf=ωf+δ. So long as |δ|ωf (so that |ωf+δ|=|ωf|+δ and thus the rotating wave approximation does not change), Eq. (B1) shows that the two are related by

dρ̃fdt(t;ωx)=dρ̃fdt(t;ωxδ)eiκfδτx.
(B2)

Because both coherences are assumed to have the same initial conditions [ρ0()=ρ0()=0], the equality also holds when both sides of the equation are integrated. The phase factor eiκfδτx in the substitution arises from Eq. (1), where the pulse carrier frequency maintains its phase within the pulse envelope for all delays.

The resonance translation can be extended to higher order signals as well. For a third-order signal, we compare systems with transition frequencies ω10=ω10+a and ω21=ω21+b. The extension of Eq. (B2) to pathway Vβ gives

ρ̃3(t;ω2,ω2,ω1)=ρ̃3(t;ω2a,ω2a,ω1b)×eiκ2aτ2eiκ2aτ2eiκ1bτ1.
(B3)

The translation of each laser coordinate depends on which transition is made (e.g., a for transitions between |0 and |1 or b for transitions between |1 and |2), so the exact translation relation differs between pathways. We can now compute the ensemble average of the signal for pathway Vβ as a convolution between the distribution function of the system, K(a, b), and the single oscillator response,

ρ̃3(t;ω2,ω2,ω1)=K(a,b)×ρ̃3(t;ω2+a,ω2+a,ω1+b)×eiκ2aτ2eiκ2aτ2eiκ1bτ1dadb.
(B4)

For this work, we restrict ourselves to a simpler ensemble where all oscillators have equally spaced levels (i.e., a = b). This makes the translation identical for all pathways and reduces the dimensionality of the convolution. Since pathways follow the same convolution, we may also perform the convolution on the total signal field,

Etot(t)=Lμ4,LK(a,a)×ρ̃3,L(t;ωxa,ωyaωza)×eiaκxτx+κyτy+κzτzda.
(B5)

Furthermore, since κ=1 for E1 and E2, while κ=1 for E2, we have eiaκxτx+κyτy+κzτz=eiaτ1τ2+τ2 for all pathways. Equivalently, if the electric field is parameterized in terms of laser coordinates ω1 and ω2, the ensemble field can be calculated as

Etot(t;ω1,ω2)=K(a,a)Etot(t;ω1a,ω2a)×eiaτ1τ2+τ2da,
(B6)

which is a 1D convolution along the diagonal axis in frequency space. Figure 13 demonstrates the use of Eq. (B6) on a homogeneous line shape.

FIG. 13.

Overview of the convolution. (a) The homogeneous line shape. (b) The distribution function, K, mapped onto laser coordinates. (c) The resulting ensemble line shape computed from the convolution. The thick black line represents the FWHM of the distribution function.

FIG. 13.

Overview of the convolution. (a) The homogeneous line shape. (b) The distribution function, K, mapped onto laser coordinates. (c) The resulting ensemble line shape computed from the convolution. The thick black line represents the FWHM of the distribution function.

Close modal
2.
3.
K.
Park
and
M.
Cho
,
J. Chem. Phys.
109
,
10559
(
1998
).
4.
S. M.
Gallagher
,
A. W.
Albrecht
,
J. D.
Hybl
,
B. L.
Landin
,
B.
Rajaram
, and
D. M.
Jonas
,
J. Opt. Soc. Am. B
15
,
2338
(
1998
).
5.
M. V.
Lebedev
,
O. V.
Misochko
,
T.
Dekorsy
, and
N.
Georgiev
,
J. Exp. Theor. Phys.
100
,
272
(
2005
).
6.
Z.
Vardeny
and
J.
Tauc
,
Opt. Commun.
39
,
396
(
1981
).
7.
M.
Joffre
,
C. B.
à la Guillaume
,
N.
Peyghambarian
,
M.
Lindberg
,
D.
Hulin
,
A.
Migus
,
S. W.
Koch
, and
A.
Antonetti
,
Opt. Lett.
13
,
276
(
1988
).
8.
9.
S.
Druet
,
J.-P.
Taran
, and
C. J.
Bordé
,
J. Phys.
40
,
819
(
1979
).
10.
J.-L.
Oudar
and
Y. R.
Shen
,
Phys. Rev. A
22
,
1141
(
1980
).
11.
J. C.
Wright
,
P. C.
Chen
,
J. P.
Hamilton
,
A.
Zilian
, and
M. J.
Labuda
,
Appl. Spectrosc.
51
,
949
(
1997
).
12.
J. C.
Wright
,
R. J.
Carlson
,
G. B.
Hurst
,
J. K.
Steehler
,
M. T.
Riebe
,
B. B.
Price
,
D. C.
Nguyen
, and
S. H.
Lee
,
Int. Rev. Phys. Chem.
10
,
349
(
1991
).
13.
A. V.
Pakoulev
,
M. A.
Rickard
,
K. A.
Meyer
,
K.
Kornau
,
N. A.
Mathew
,
D. E.
Thompson
, and
J. C.
Wright
,
J. Phys. Chem. A
110
,
3352
(
2006
).
14.
A. V.
Pakoulev
,
M. A.
Rickard
,
N. A.
Mathew
,
K. M.
Kornau
, and
J. C.
Wright
,
J. Phys. Chem. A
111
,
6999
(
2007
).
15.
D. D.
Kohler
,
S. B.
Block
,
S.
Kain
,
A. V.
Pakoulev
, and
J. C.
Wright
,
J. Phys. Chem. C
118
,
5020
(
2014
).
16.
M. F.
Gelin
,
D.
Egorova
, and
W.
Domcke
,
Acc. Chem. Res.
42
,
1290
(
2009
).
17.
P.
Hamm
,
M.
Lim
,
W. F.
DeGrado
, and
R. M.
Hochstrasser
,
J. Chem. Phys.
112
,
1907
(
2000
).
18.
J. R.
Salcedo
,
A. E.
Siegman
,
D. D.
Dlott
, and
M. D.
Fayer
,
Phys. Rev. Lett.
41
,
131
(
1978
).
19.
J. T.
Fourkas
,
R.
Trebino
, and
M. D.
Fayer
,
J. Chem. Phys.
97
,
69
(
1992
).
20.
J. T.
Fourkas
,
R.
Trebino
, and
M. D.
Fayer
,
J. Chem. Phys.
97
,
78
(
1992
).
21.
G.
Auböck
,
C.
Consani
,
F.
van Mourik
, and
M.
Chergui
,
Opt. Lett.
37
,
2337
(
2012
).
22.
H. J.
Bakker
,
H.-K.
Nienhuys
,
G.
Gallot
,
N.
Lascoux
,
G. M.
Gale
,
J.-C.
Leicknam
, and
S.
Bratos
,
J. Chem. Phys.
116
,
2592
(
2002
).
23.
W. P.
de Boeij
,
M. S.
Pshenichnikov
, and
D. A.
Wiersma
,
Chem. Phys. Lett.
253
,
53
(
1996
).
24.
F.
Patterson
,
H.
Lee
,
W. L.
Wilson
, and
M.
Fayer
,
Chem. Phys.
84
,
51
(
1984
).
25.
A.
Tokmakoff
and
M. D.
Fayer
,
Acc. Chem. Res.
28
,
437
(
1995
).
26.
P.
Hamm
,
M.
Lim
,
W. F.
DeGrado
, and
R. M.
Hochstrasser
,
Proc. Natl. Acad. Sci. U. S. A.
96
,
2036
(
1999
).
27.
M. C.
Asplund
,
M. T.
Zanni
, and
R. M.
Hochstrasser
,
Proc. Natl. Acad. Sci. U. S. A.
97
,
8219
(
2000
).
28.
M. T.
Zanni
,
M. C.
Asplund
, and
R. M.
Hochstrasser
,
J. Chem. Phys.
114
,
4579
(
2001
).
29.
J. D.
Hybl
,
Y.
Christophe
, and
D. M.
Jonas
,
Chem. Phys.
266
,
295
(
2001
).
30.
T.
Brixner
,
T.
Mančal
,
I. V.
Stiopkin
, and
G. R.
Fleming
,
J. Chem. Phys.
121
,
4221
(
2004
).
31.
M. F.
Emde
,
A.
Baltuska
,
A.
Kummrow
,
M. S.
Pshenichnikov
, and
D. A.
Wiersma
,
Phys. Rev. Lett.
80
,
4645
(
1998
).
32.
W. P.
de Boeij
,
M. S.
Pshenichnikov
, and
D. A.
Wiersma
,
Chem. Phys. Lett.
238
,
1
(
1995
).
33.
M.
Cho
,
N. F.
Scherer
,
G. R.
Fleming
, and
S.
Mukamel
,
J. Chem. Phys.
96
,
5618
(
1992
).
34.
S. A.
Passino
,
Y.
Nagasawa
,
T.
Joo
, and
G. R.
Fleming
,
J. Phys. Chem. A
101
,
725
(
1997
).
35.
S.
Yoon
,
D. W.
McCamant
,
P.
Kukura
,
R. A.
Mathies
,
D.
Zhang
, and
S.-Y.
Lee
,
J. Chem. Phys.
122
,
024505
(
2005
).
36.
D. W.
McCamant
,
P.
Kukura
, and
R. A.
Mathies
,
J. Phys. Chem. B
109
,
10449
(
2005
).
37.
W.
Zhao
and
J. C.
Wright
,
Phys. Rev. Lett.
83
,
1950
(
1999
).
38.
W.
Zhao
and
J. C.
Wright
,
J. Am. Chem. Soc.
121
,
10994
(
1999
).
39.
W.
Zhao
,
K. M. M.
Urdoch
,
D. M.
Besemann
,
N. J.
Condon
,
K. A.
Meyer
, and
J. C.
Wright
,
Appl. Spectrosc.
54
,
1000
(
2000
).
40.
K. A.
Meyer
and
J. C.
Wright
,
J. Phys. Chem. A
107
,
8388
(
2003
).
41.
P. M.
Donaldson
,
R.
Guo
,
F.
Fournier
,
E. M.
Gardner
,
L. M. C.
Barter
,
C. J.
Barnett
,
I. R.
Gould
,
D. R.
Klug
,
D. J.
Palmer
, and
K. R.
Willison
,
J. Chem. Phys.
127
,
114513
(
2007
).
42.
P. M.
Donaldson
,
R.
Guo
,
F.
Fournier
,
E. M.
Gardner
,
I. R.
Gould
, and
D. R.
Klug
,
Chem. Phys.
350
,
201
(
2008
).
43.
F.
Fournier
,
E. M.
Gardner
,
R.
Guo
,
P. M.
Donaldson
,
L. M.
Barter
,
D. J.
Palmer
,
C. J.
Barnett
,
K. R.
Willison
,
I. R.
Gould
, and
D. R.
Klug
,
Anal. Biochem.
374
,
358
(
2008
).
44.
E. S.
Boyle
,
N. A.
Neff-Mallon
, and
J. C.
Wright
,
J. Phys. Chem. A
117
,
12401
(
2013
).
45.
E. S.
Boyle
,
A. V.
Pakoulev
, and
J. C.
Wright
,
J. Phys. Chem. A
117
,
5578
(
2013
).
46.
E. S.
Boyle
,
N. A.
Neff-Mallon
,
J. D.
Handali
, and
J. C.
Wright
,
J. Phys. Chem. A
118
,
3112
(
2014
).
47.
A.
Lagutchev
,
S. A.
Hambir
, and
D. D.
Dlott
,
J. Phys. Chem. C
111
,
13645
(
2007
).
48.
R. J.
Carlson
,
D. C.
Nguyen
, and
J. C.
Wright
,
J. Chem. Phys.
92
,
1538
(
1990
).
49.
R. J.
Carlson
and
J. C.
Wright
,
J. Chem. Phys.
92
,
5186
(
1990
).
50.
R. J.
Carlson
and
J. C.
Wright
,
Anal. Chem.
63
,
1449
(
1991
).
51.
J. K.
Steehler
and
J. C.
Wright
,
J. Chem. Phys.
83
,
3200
(
1985
).
52.
K.
Kwac
and
M.
Cho
,
J. Phys. Chem. A
107
,
5903
(
2003
).
53.
J. C.
Wright
,
Chem. Phys. Lett.
662
,
1
(
2016
).
54.
V.
Perlík
,
J.
Hauer
, and
F.
Šanda
,
J. Opt. Soc. Am. B
34
,
430
(
2017
).
55.
C. L.
Smallwood
,
T. M.
Autry
, and
S. T.
Cundiff
,
J. Opt. Soc. Am. B
34
,
419
(
2017
).
56.
A. V.
Pakoulev
,
M. A.
Rickard
,
K. M.
Kornau
,
N. A.
Mathew
,
L. A.
Yurs
,
S. B.
Block
, and
J. C.
Wright
,
Acc. Chem. Res.
42
,
1310
(
2009
).
57.
K. J.
Czech
,
B. J.
Thompson
,
S.
Kain
,
Q.
Ding
,
M. J.
Shearer
,
R. J.
Hamers
,
S.
Jin
, and
J. C.
Wright
,
ACS Nano
9
,
12146
(
2015
).
58.
J. C.
Wright
,
N. J.
Condon
,
K. M.
Murdoch
,
D. M.
Besemann
, and
K. A.
Meyer
,
J. Phys. Chem. A
107
,
8166
(
2003
).
59.
P. M.
Donaldson
,
K. R.
Willison
, and
D. R.
Klug
,
J. Phys. Chem. B
114
,
12175
(
2010
).
60.
K. A.
Meyer
,
J. C.
Wright
, and
D. E.
Thompson
,
J. Phys. Chem. A
108
,
11485
(
2004
).
61.
S.
Mukamel
,
Principles of Nonlinear Optical Spectroscopy
, Oxford Series in Optical and Imaging Sciences (
Oxford University Press
,
1995
).
62.
T. K.
Yee
and
T. K.
Gustafson
,
Phys. Rev. A
18
,
1597
(
1978
).
63.
J. A.
Armstrong
,
N.
Bloembergen
,
J.
Ducuing
, and
P. S.
Pershan
,
Phys. Rev.
127
,
1918
(
1962
).
64.
I. V.
Schweigert
and
S.
Mukamel
,
Phys. Rev. A
77
,
033802
(
2008
).
65.
D.
Lee
and
A.
Albrecht Ferro
,
Advances in Infrared and Raman Spectroscopies
(
Wiley Heyden
,
1985
), Vol. 12, pp.
179
213
.
66.

κf also has a direct relationship to the phase-matching relationship: for transitions with E2, κf=1, and for E1 or E2, κf=1.

67.

For elucidation of the relationship between the generalized Liouville pathway notation and the specific parameters for each Liouville pathway, see Table S1 in the supplementary material.

68.
D. M.
Besemann
,
K. A.
Meyer
, and
J. C.
Wright
,
J. Phys. Chem. B
108
,
10493
(
2004
).
69.
R. J.
Carlson
and
J. C.
Wright
,
J. Mol. Spectrosc.
143
,
1
(
1990
).
70.
M. T.
Riebe
and
J. C.
Wright
,
J. Chem. Phys.
88
,
2981
(
1988
).
71.
A. M.
Weiner
,
S.
De Silvestri
, and
E. P.
Ippen
,
J. Opt. Soc. Am. B
2
,
654
(
1985
).
72.
R.
Agarwal
,
B. S.
Prall
,
A. H.
Rizvi
,
M.
Yang
, and
G. R.
Fleming
,
J. Chem. Phys.
116
,
6243
(
2002
).
73.
B.
Dick
and
R.
Hochstrasser
,
Chem. Phys.
75
,
133
(
1983
).
74.
M. F.
Gelin
,
D.
Egorova
, and
W.
Domcke
,
J. Chem. Phys.
123
,
164112
(
2005
).
75.
T. E.
Oliphant
,
Comput. Sci. Eng.
9
,
10
(
2007
).
76.
A.
Lagutchev
,
A.
Lozano
,
P.
Mukherjee
,
S. A.
Hambir
, and
D. D.
Dlott
,
Spectrochim. Acta, Part A
75
,
1289
(
2010
).
77.
C.
Brito Cruz
,
J.
Gordon
,
P.
Becker
,
R.
Fork
, and
C.
Shank
,
IEEE J. Quantum Electron.
24
,
261
(
1988
).
78.
S. L.
Palfrey
and
T. F.
Heinz
,
J. Opt. Soc. Am. B
2
,
674
(
1985
).
79.
G. R.
Fleming
,
S. A.
Passino
, and
Y.
Nagasawa
,
Philos. Trans. R. Soc., A
356
,
389
(
1998
).
80.
W. P.
de Boeij
,
M. S.
Pshenichnikov
, and
D. A.
Wiersma
,
Annu. Rev. Phys. Chem.
49
,
99
(
1998
).
81.
M. R.
Salvador
,
M. A.
Hines
, and
G. D.
Scholes
,
J. Chem. Phys.
118
,
9380
(
2003
).
82.
L. A.
Yurs
,
S. B.
Block
,
A. V.
Pakoulev
,
R. S.
Selinsky
,
S.
Jin
, and
J.
Wright
,
J. Phys. Chem. C
115
,
22833
(
2011
).
83.

There are many ways to characterize the ellipticity of a peak shape. We adopt the convention E=(a2b2)/(a2+b2), where a is the diagonal width and b is the anti-diagonal width.

84.
K.
Okumura
,
A.
Tokmakoff
, and
Y.
Tanimura
,
Chem. Phys. Lett.
314
,
488
(
1999
).
85.
T. J.
Aartsma
and
D. A.
Wiersma
,
Chem. Phys. Lett.
42
,
520
(
1976
).

Supplementary Material