A density matrix treatment of plasmon-enhanced (PE) stimulated Raman spectroscopies is developed. Specifically, PE stimulated Raman Gain/Loss (PE-SRG/L) and coherent anti-Stokes Raman scattering (PE-CARS) due to monochromatic excitation and PE femtosecond stimulated Raman spectroscopy (PE-FSRS) are considered. A Lorentz oscillator model is used to explicitly describe the time dependence of plasmon-enhanced optical fields. These temporal characteristics are required for a density matrix based description of all plasmon-enhanced nonlinear molecular spectroscopies. Dispersive vibrational line shapes in PE-SRG/L and PE-FSRS spectra are shown to result primarily from terms proportional to the square of the complex optical field enhancement factor. The dependence on the plasmon resonance, picosecond and femtosecond pulse characteristics, and molecular vibrational properties are evident in the density matrix derived PE-FSRS intensity expression. The difference in signal detection mechanisms accounts for the lack of dispersive line shapes in PE spontaneous Raman spectroscopy. This density matrix treatment of PE-FSRS line shapes is compared with prior coupled wave results.

## I. INTRODUCTION

The ability of spontaneous and stimulated Raman spectroscopies to provide label-free, diffraction limit imaging capabilities, particularly for a range of biomedical applications, has been well established and offers advantages relative to traditional imaging techniques owing to its multiplexing molecular specificity.^{1–12} However, the inherent low cross section of non-electronically resonant Raman scattering often limits both spontaneous and stimulated Raman-based imaging techniques to probing distributions and dynamics of more highly concentrated analytes (∼mM and greater) and thus limits the impact of these techniques to provide near “real time” imaging that is important for many biological and biomedical applications. Electronic resonance^{13,14} or labeling/tagging with vibrators of large Raman cross sections^{1,15} are approaches that have been exploited to overcome this sensitivity limitation. Another approach that has been famously and widely used to boost effective spontaneous Raman signal intensities by as much as ∼10^{9}/molecule (ensemble average) is surface-enhanced Raman spectroscopy (SERS). As has been well established, on noble metals (e.g., Au or Ag), most of the observed enhancement is due to resonance of the excitation and scattered optical fields with the surface plasmon excitations of the nanostructured metal for molecules within a few nm of these surfaces.^{16–21} In addition, in favorable excitation conditions, physi-adsorbed analytes can have their SERS intensities boosted by another ∼10^{2} owing to chemical enhancement effects arising from resonant or near resonant weak metal–molecule charge transfer excitations.^{22–24} While SERS has been exploited in some microscopy/imaging studies,^{25–29} surface-enhanced stimulated Raman scattering methods, analogously deriving most of their signal amplification from plasmon enhanced (PE) optical fields, can provide higher sensitivities, thus enabling more rapid chemical imaging for “real time” applications and the ability to monitor lower concentration analytes. Furthermore, PE stimulated Raman techniques may further our understanding of molecule–surface interactions due to the optical radiation phase and amplitude control that can result from rational plasmonic design efforts. Thus, there is considerable current interest in characterizing and exploiting PE stimulated Raman techniques.

The enhancement of several stimulated Raman spectroscopies due to plasmonic effects has been experimentally demonstrated, albeit with varying degrees of enhancement efficiencies.^{30} These include surface-enhanced coherent anti-stokes Raman scattering (CARS),^{31–38} stimulated Raman gain (SRG) or loss (SRL),^{39–42} femtosecond stimulated Raman spectroscopy (FSRS),^{43–47} sum frequency generation (SFG),^{48–50} and impulsive stimulated pump–probe Raman scattering (ISRS).^{51} Both time and frequency domain, heterodyne and homodyne detected techniques are among this group of stimulated Raman methods. While the nomenclature of such experimental studies are rightfully characterized as surface-enhanced (SE) spectroscopies owing to all the potential enhancement mechanisms contributing on a given nanostructured surface, theoretical or modeling treatments explicitly considering local plasmonic effects are accurately describing plasmon-enhanced (PE) stimulated Raman spectroscopies, and that nomenclature will be employed here.

The appearance of dispersive vibrational line shapes are often observed in homodyne detected spectroscopies with Raman resonances, such as CARS and SFG. It is well established that these spectral features arise from the contributions of cross-terms between the complex Raman susceptibility and other vibrationally resonant and non-resonant electronic background contributions to the observed homodyne signal. However, more recently, dispersive line shapes have been observed in heterodyne detected SE-FSRS and SE-SRG/SE-SRL spectra.^{40,43,44,46,47} The origins of these complex vibrational line shapes has been discussed in some previous reports,^{45,52,53} and the dependence of these surface-enhanced coherent Raman signals on plasmonic enhancement factors has been described.^{34,36,45,54,55} No such dispersive vibrational line shapes are evident in observed plasmon-enhanced spontaneous Raman spectra, i.e., SERS.

Two mechanisms have been proposed as possible causes for these complex line shapes in heterodyne detected PE stimulated Raman spectra. It was shown that dispersive Raman line shapes could result from interference or heterodyne cross-terms between the plasmon-enhanced stimulated Raman signal field generated by molecules located near metal nanoparticles and a stimulated nonlinear emission/scattering field due to an optically resonant surface plasmon resonance (SPR) emission field acting as a local oscillator at the detector.^{52} A standard perturbative polarization approach was used to treat the FRSR and SPR emission signal fields and captured the qualitative features of the reported SE-FSRS spectra and their wavelength dependence.

Concurrently, a semi-classical coupled wave equation approach^{56,57} was employed to describe PE-FSRS. An analytical expression for PE-FSRS line shapes was derived for an interpulse delay (*τ*) of 0. In the coupled wave treatment dispersively shaped PE-FSRS, vibrational line shapes were shown to fundamentally result from the complex character of the frequency dependent field enhancement factor, *g*(*ω*).^{45} More precisely, this coupled wave analysis found that the dispersive contribution to the PE-FSRS, viewed along the SRG direction, is due to $Img(\Omega S)2$, where Ω_{S} is the Stokes Raman frequency. Qualitative agreement between experiment and theory was achieved by this treatment, although observed dispersion effects were not quantitatively captured.^{45,46}

More recently, Raman vibrational line shapes with dispersive character were also observed in heterodyne detected, surface enhanced SRG and SRL spectra.^{40,53} In a theoretical treatment of these PE-SFG/L spectra, a phasor based framework was used to describe the observed Raman line shapes for monochromatic pump and Stokes (anti-Stokes) radiation, and the surface plasmon resonance was modeled with an *ad hoc* Lorentzian shaped spectrum. The dispersive character of the PE-SRG/PE-SRL line shapes were attributed to Im *g*(Ω_{S})/Im *g*(Ω_{P}), respectively, in contrast to the dependence on the imaginary part of the *square* of the complex enhancement factor in the PE-FSRS coupled wave analysis.^{53} Finally, we note the enhancement factor dependence of PE-CARS has been described by several groups^{36,54} and generally identified as $g(\Omega P4g(\Omega S)2g(\Omega AS)2$.

The purpose here is to model the time dependence of the local plasmon enhanced optical fields and use these PE-fields in a density matrix approach to describe plasmon-enhanced SRG/SRL, CARS, and FSRS. While formally a time domain description of plasmonic fields is not necessary for stimulated PE Raman spectroscopies due to monochromatic incident pump and probe pulses, it is useful for developing a density matrix based consistent framework for descriptions of all PE spectroscopies and, in particular, for those experimental methods employing ultrafast pulses. A quantum mechanical, density matrix based description of the effects of a complex enhancement factor for surface enhanced stimulated Raman spectroscopies is applicable, in general, for both time and frequency domain, heterodyne and homodyne detected spectroscopies with Raman resonances. Although the general, density matrix framework outlined here is applicable to all types of plasmon-enhanced Raman spectrsocopies, this paper is limited to expressions for cw PE-SRG/L and PE-CARS and PE-FSRS where the complex character of the plasmonic enhancement is the sole potential cause of dispersive Raman line shapes for electronically non-resonant radiation. The formal Raman signal detection differences for spontaneous and stimulated Raman spectroscopies is also noted here.

## II. DESCRIPTIONS OF PLASMON-ENHANCED OPTICAL FIELDS

Formally, a quantum mechanical description of the time evolution of a system in response to incident radiation requires an explicit characterization of the time dependence of all the active optical electric fields. An incident electric field with carrier frequency, Ω_{j}, and pulse duration and amplitude given by $E\u0302j(t)$ is described classically by

The corresponding spectrum of this optical electric field is determined by Fourier transform,

### A. Relevant constitutive relations

In analogy to the well-established constitutive relations relating harmonic electric and magnetic fields to macroscopic responses, such as electric polarization, magnetic induction, and current density,^{58,59} we can describe a susceptibility or enhancement factor, *g*(*ω*), giving rise to the plasmon-enhanced complex electric field amplitude or induced polarization component at a given frequency by

where *g*(*ω*) is proportional to the polarizability of the nanostructured metal and its magnitude is also dependent on position and orientation factors. The focus in this discussion, however, is on the time and frequency dependence of plasmon-enhanced optical fields, and thus, explicit consideration of these additional factors will be neglected here. In the time domain, this constitutive relation [Eq. (2.3)] is equivalently described for harmonic fields by a convolution of the incident field time dependence with a materials specific plasmonic response function, *g*(*t*).^{58,59} Consequently, the time dependence of a plasmon-enhanced field [$E\u0302PE(t)$] relative to the incident field [$E\u0302(t)$] is given by

and the corresponding PE optical electric field spectrum is

### B. Enhancement response function derived from a Lorentz oscillator dielectric constant model

The dielectric response of a metal as a function of frequency due to a single electronic resonance modeled as a periodic force driven Lorentz oscillator has the well-known form^{59}

where *ω*_{p} and *ω*_{o} are the plasma resonance and the resonant frequency of the system, respectively, and 2Γ is the damping coefficient. Equation (2.6) also captures the free electron Drude contribution to a dielectric for *ω*_{o} = 0. Using the well-established enhancement factor resulting from the surface plasmon excitation of a spherical metal nanoparticle as a model plasmonic structure, the field enhancement factor is given by^{16,60–62}

where *ɛ*_{m} is the dielectric of the ambient medium and taken to be 1 here and *a* is a size and distance dependent factor. Other values for these constants can capture the effects of nanoparticle shape and solvent medium. For *ɛ*(*ω*) defined above [Eq. (2.6)], *g*(*ω*) resulting for this SPR is then

The resonant frequency of the surface plasmon is defined by $\Omega PL2=\omega o2+\omega p2/2$, and the damping constant (Γ) is more specifically identified as Γ_{PL}. The corresponding time dependent susceptibility is given by the Fourier transform,

where $A=\pi 21/2a\omega p2\Omega PL$ and we take $\Omega PL2\u226b\Gamma PL2$. Consequently, $G(t)=12g(t)+12g*(t)$ and

In the region of interest near the plasmon resonance for an underdamped plasmon response ($\Omega PL2\u226b\Gamma PL2$), Fourier transform of *g*(*t*) [Eq. (2.11)] yields a slightly simpler form for the enhancement factor spectrum, *g*(*ω*), than defined in Eq. (2.8),

and Δ_{ω} = Ω_{PL} − *ω*. [Alternatively, Eq. (2.12) results from Eq. (2.8) for the near resonance condition Ω_{PL} + *ω* ≈ 2*ω*.] This frequency dependent *g*(*ω*) has been previously used *ad hoc* to model the plasmon line shape in some prior PE nonlinear spectroscopy treatments.^{53,55,63} As shown here, the frequency dependent dielectric constants given by a single Lorentz driven oscillator model reveals the same Lorentzian enhancement factor. It also provides the corresponding response function, *g*(*t*), for describing the local time dependence of all plasmonic enhanced fields that interact with this metal nanoparticle. The magnitude of the enhanced optical field will be assumed throughout this treatment to be much greater than the incident field, i.e., 1 + *g*(*ω*) ≈ *g*(*ω*).

### C. Description of plasmon-enhanced radiation

#### 1. cw radiation

For an incident optical electric field as defined in Eq. (2.1) for cw monochromatic radiation, $E\u0302j(t)=E\u0302jo$, where the electric field amplitude $E\u0302jo$ is a constant. Thus, the plasmon-enhanced radiation at incident frequency Ω_{j} is simply [Eq. (2.4)]

The monochromatic field amplitude is enhanced by a Lorentzian factor dependent on the detuning between the surface plasmon resonance and the incident radiation frequencies (Δ_{j} = Ω_{PL} − Ω_{j}), and each Fourier component (±Ω_{j}) amplitude is now complex as a result of the plasmon enhancement.

#### 2. Picosecond pulses

For an incident Gaussian pulse with duration $2ln2\sigma j$, the time dependent amplitude [Eq. (2.1)] is $E\u0302j(t)=E\u0302joe\u2212t2/2\sigma j2$. The corresponding PE pulse is given by [Eq. (2.4)]

However, the picosecond scale pulse duration is much longer than characteristic surface plasmon relaxation times ($\Gamma PL\u22121$ ∼10 fs or less^{62,64}), and thus, the time dependence of the plasmon-enhanced picosecond pulse is given by the same expression as for cw excitation,

The PE picosecond scale pulse has the same temporal profile as the incident pulse but with an enhanced amplitude that is now complex whose phase is dependent on the detuning between the pulse and plasmon resonance frequencies described by *g*(Δ_{P}).

#### 3. Femtosecond pulse

For an incident ultrafast femtosecond pulse, the pulse duration and the surface plasmon damping time are of the same order, and thus, the required time domain convolution integral needs to be rigorously solved to derive an analytical expression for the amplitude and phase of the plasmonic enhanced ultrashort pulse. For a Gaussian incident pulse with a femtosecond scale pulse duration (FWHM) $2ln2\sigma j$, the plasmon-enhanced pulse is given by (see Sec. II of the supplementary material)

where

Again, Δ_{j} = Ω_{PL} − Ω_{j} is the detuning of the ultrafast pulse carrier frequency from the plasmon resonance and *s*′ is the complex argument in the error function (*erf*) defined by

As an illustration of the influence of a plasmonic resonance on the electric field envelop of an ultrafast pulse, the effects of a plasmon resonance with 1/Γ_{PL} = 10 fs, corresponding to a 1000 cm^{−1} FWHM, on the electric field of a 40 fs incident pulse are shown in Fig. 1. The incident, all real 40 fs pulse field and the corresponding normalized, complex plasmon-enhanced fields for carrier frequencies coincident with the plasmon resonant frequency (Δ_{F} = 0) and detuned 500 cm^{−1} from the plasmon resonance (Δ_{F} = 500 cm^{−1}) are shown in this figure. The changes in the temporal profile of this incident pulse, most notably the phase effects, are evident in this figure. A fuller description of the effects of plasmon responses on optical fields and the consequent ultrafast spectroscopic signal responses will be given in a subsequent report.

## III. PLASMON-ENHANCED STIMULATED RAMAN GAIN/LOSS (PE-SRG/L) AND CARS (PE-CARS)

### A. PE-SRG/PE-SRL

Expressions for PE-SRG/SRL spectra corresponding to an electronically non-resonant, three-level system (Fig. 2) resulting from two cw beams can be derived via a standard density matrix approach (supplementary material, Sec. I). The RWA (rotating wave approximation) density matrix histories contributing to these signals are represented by the WMEL (wave-mixing energy level) diagrams^{65} shown in Fig. 2. Following the discussion in Sec. II, the plasmon-enhanced optical electric fields at the monochromatic pump (Ω_{P}) and Stokes (Ω_{S}) frequencies are given by

where *g*(Δ_{P}) and *g*(Δ_{S}) are the complex electric field enhancement factors [Eq. (2.13)] for the Raman pump (P) and Stokes (S) radiation, where Δ_{P} = Ω_{PL} − Ω_{P} and Δ_{S} = Ω_{PL} − Ω_{S}, respectively, and $E\u0302P/So$ is the corresponding incident P/S field amplitude. With the carrier frequency dependence aside, a real incident pulse amplitude acquires a complex amplified amplitude with a phase dependent on the detuning from the plasmon resonant frequency given by *g*(Δ_{P/S}).

The third-order induced optical polarization at the signal frequency is defined by

The signal frequency, Ω_{sig}, is at Ω_{S}/Ω_{P} for SRG/SRL, and the complex third-order polarization amplitudes corresponding to the density matrix time evolution histories described in Fig. 2(a) for PE-SRG and in Fig. 2(b) PE-SRL are given by (see Sec. I of the supplementary material)

where $\rho \u0302gg(0)$ is the thermal population of the initial ground state level *g*, $\mu \u0302eg=e\mu \u0302g$ is the electronic transition moment between molecular states *e* and *g*, etc., and the electronic and vibrational detunings are defined by Δ_{1} = *ω*_{eg} − Ω_{P} and Δ_{v} = *ω*_{fg} − Ω_{P} + Ω_{S}, respectively. *ω*_{eg} is the energy gap between the ground state *g* and the representative virtual level *e* (Fig. 2). The RWA and Condon approximation have been used to describe the Raman susceptibilities [Eqs. (2.3) and (2.4)] for this non-electronic resonant, three-level system. Finally, *γ*_{fg} is the homogeneous (exponentially decaying) vibrational dephasing rate corresponding to the decay of the Raman mode coherence at frequency *ω*_{fg} = (*E*_{f} − *E*_{g})/*ℏ*.

The optically induced third-order polarization acts as a source term in the optical wave equation locally generating a new signal field with wavevector *k*_{sig} in the macroscopically phase matched direction. Thus, *k*_{sig} = *k*_{S} or *k*_{sig} = *k*_{P} for SRG or SRL, respectively, as demanded by phase-matching constraints, and thus, the signal field is perfectly aligned for heterodyne detection in typical pump–probe configurations. The locally generated coherent signal field is, however, also plasmon-enhanced, in addition to the incident pump and Stokes fields, and analogously described in this semiclassical framework by

where $E\u0302sigPE(t)=g(\Delta sig)E\u0302sig(t)$ at the signal frequency. By the optical wave equation solution, the signal field amplitude is *π*/2 out of phase with signal polarization amplitude,^{66}

where *L* is the length of the sample giving rise to this coherent signal. At the detector in the *k*_{sig} beam path, the gain/loss in the Stokes/pump beam intensity due to the stimulated Raman transition is given by

This mechanism of signal generation presumes that the Stokes or pump field ($E\u0302S/Po$) acting as the local oscillator in Eq. (3.7) is not concurrently amplified and consequently phase shifted by the plasmonic resonance as it passes through the sample medium. If it were, then additional interference terms could play a role in determining the observed intensity and line shape of the stimulated Raman spectrum.^{52}

Neglecting the weaker homodyne term, the plasmon-enhanced stimulated Raman signals are

Combining the results summarized in Eq. (3.3) with the signal field definition (3.6), the plasmon-enhanced stimulated Raman gain signal spectrum at the detector is

where $Reg(\Delta S)2=(Reg(\Delta S))2\u2212(Img(\Delta S))2$ and $Img(\Delta S)2=2Reg(\Delta S)Img(\Delta S)$ and

Thus, the line shape of PE-SRG features will be dependent on $g(\Delta P)2g(\Delta S)2$, the square of the complex enhancement factor at the Stokes frequency, not the complex conjugate square $g(\Delta P)2g(\Delta S)2$ as in spontaneous Raman (*vide infra*) and not the relative magnitudes of Re *g*(Δ_{S}) and Im *g*(Δ_{S}) as given in a recent SRG/SRL treatment.^{39} This density matrix derived result for plasmon-enhanced cw optical fields is consistent with the classical coupled wave result developed for the description of SRG within the context of the PE-FSRS treatment.^{45} The main conclusion here and there is that the dispersive Raman line shapes are due to the complex $g(\Delta S)2$ factor, which mixes the real and imaginary components of a Raman Lorentzian line shape function. As seen in Eq. (3.10), the SRG spectrum will recover its purely positive definite Lorentzian character for an enhancement factor *g*(Δ_{S}) that is all real.

Analogously, the plasmon-enhanced stimulated Raman loss signal is

Thus, a dispersive contribution to the PE-SRL vibrational line shape will be determined by $Img(\Delta P)2$ not Im *g*(Δ_{P}).^{53} Thus, regardless of additional complex heterodyning field effects,^{52} PE-SRG/SRL line shapes can be expected to exhibit dispersive character, in general, due to just the complex character of the plasmonic enhancement factor.

### B. PE-CARS

For completeness, this density matrix treatment can analogously be applied to the description of PE-CARS. A WMEL density matrix time evolution history for CARS is shown in Fig. 3. For this homodyne signal, the Raman response results from two pump fields and a Stokes field all acting on the bra or ket-side and a new signal anti-Stokes to the pump frequency results: Ω_{AS} = 2 Ω_{P} − Ω_{S}. The generated CARS third-order polarization when all the incident cw fields are enhanced by the plasmonic resonance is

where (see Sec. I of the supplementary material)

As described above for SRL/SRG, this third-order polarization acts as a source term in the wave equation generating a coherent signal field at the Ω_{AS} frequency [Eq. (3.6)] with an amplitude proportional to the third order polarization amplitude [Eq. (3.14)],

For PE-CARS resulting from cw fields, the locally generated anti-Stokes signal field, $E\u0302sigCARS$, is also enhanced in addition to the incident pump and Stokes laser fields,

Thus, the observed intensity of the homodyne detected PE-CARS is given by

The complex character of the field enhancement factors, *g*(Δ_{P,S,AS}), does not distort the observed PE-CARS line shape from its indicated Lorentzian character if the contributions of all other resonant and non-resonant optical contributions along the phase matched direction can be neglected. As stated previously^{34,36,54} and as found here, the plasmonic enhancement of the CARS signal is thus given by

Accordingly, whether a PE-CARS line shape shows dispersive character is dependent on the magnitude of the non-resonant background terms arising from electronic and nuclear degrees of freedom in the specific sample system and discussed in prior reports.^{54,55} The experimental challenge for being able to maximally exploit the potential amplification power of the ∼*g*^{8} dependence [Eq. (3.18)] for PE-CARS may be in finding a plasmonic surface with sufficiently rapid relaxation to provide enhancement that can encompass a frequency range on the order of 2 *ω*_{fg} (∼Δ_{AS} − Δ_{S}).

## IV. PLASMON-ENHANCED FEMTOSECOND STIMULATED RAMAN SPECTROSCOPY (PE-FSRS)

### A. Density matrix description of PE-FSRS spectra

A similar density matrix treatment can be followed for the description of PE-FSRS. The relevant time evolution history generating the FSRS optical polarization due to a Raman resonance is represented in Fig. 2(a). However, unlike the SRG treatment above, for FSRS, the pump and Stokes probe radiation are plasmon-enhanced picosecond and femtosecond pulses. Furthermore, there may be a delay, τ, between the arrival time of the temporal peak of these pulses. As discussed above, the temporal profile of the plasmon-enhanced optical fields are described by the convolution of the time dependence of the incident pulse and the plasmon response function here given by Eq. (2.11). For this PE-FSRS analysis, the incident pulses are assumed to be Gaussian shaped with pulse durations (FWHM) $2ln2\sigma j$. From above (Sec. II), the plasmon-enhanced picosecond pump (*P*) and femtosecond probe (*F*) pulses are

where Δ_{P/F} = Ω_{PL} − Ω_{P/F} is the detuning from the carrier frequency of the incident optical pulses from the plasmon resonance and *s*′ is the complex argument in the error function (*erf*) defined in Eq. (2.18).

Given these PE temporal field profiles, the density matrix treatment for PE-FSRS follows the steps essentially described previously for FSRS^{52,56} and shown in detail in the supplementary material (supplementary material, Sec. III). When the picosecond pump and femtosecond probe pulses are overlapped in the sample with a delay τ between the peaks of these pulses, the Raman resonance at *ω*_{fg} generates a third-order polarization response at frequency Ω_{sig} = Ω_{P} − *ω*_{fg},

The complex polarization amplitude for the signal frequency Ω_{P} − *ω*_{fg} is given by

where

Here, $\Delta 2v=\omega fg\u2212(\Omega P\u2212\Omega F)=\omega fg\u2212\Delta PF$ is a measure of the off-resonance between the pump and probe pulse carrier frequency difference (Δ_{PF} = Ω_{P} − Ω_{F}) and the vibrational frequency of the Raman active mode. When $\Delta 2v\u22600$, a femtosecond pulse effect will contribute to the observed dispersive PE-FSRS line shapes in addition to the plasmon enhancement factors [*g*(Δ_{P/F})] (*vide infra*) as evident from Eq. (4.6). In passing, we note that the *τ* dependence in $F(\Delta 2v,\tau )$ was not explicitly considered in prior FSRS density matrix treatments.^{45,52}

As for the other stimulated Raman spectroscopies described above, this optically induced polarization is a source term in Maxwell’s Equations for a generated signal field in the *k*_{F} phase-matched direction: $E\u0302sig(t,\tau )=2\pi ic\Omega sigLP\u0302(3)(t,\tau )$. Analogously, this signal field is also enhanced by the local plasmonic response of the nearby nanoparticles, and thus,

Taking advantage of the Fourier transform relationship for convolutions and substituting for the plasmon enhanced pump pulse in (4.5), the spectrum of the PE-FSRS signal field is

This signal field spectrum is peaked in the vicinity of Ω_{sig} = Ω_{P} − *ω*_{fg}. The probe pulse acts as the local oscillator for heterodyne detection of the FSRS signal field and is assumed in this treatment to be unaffected by plasmonic enhancement effects. The same assumption was made in the previous coupled wave description of PE-FSRS.^{45} The resulting local oscillator spectrum is thus given by

The corresponding heterodyne detected PE-FSRS spectrum obtained after both the signal and local oscillator (femtosecond probe) have been dispersed through the spectrometer and thus temporally and spatially overlapped on the detector is

The resulting PE-FSRS spectrum for a pump–probe delay of *τ* is given by

### B. PE-FSRS dependence on plasmonic enhancement factors

In order to compare this density matrix based expression for PE-FSRS with the corresponding results for cw PE-SRG and the prior coupled wave PE-FSRS result,^{45} it is useful to fully expose the PE-FRSR Raman intensity dependence on plasmonic enhancement factors in the derived density matrix intensity expression [Eq. (4.11)]. Two plasmonic enhancement factors are explicitly evident in Eq. (4.11). The other two plasmon enhancement factors are in $F*(\Delta 2v,\tau )$ [Eq. (4.6)]. For *τ* = 0, it is straightforward, taking advantage of the convolution theorem, to show that for a Gaussian shaped ultrafast probe pulse (supplementary material Sec. IV),

Thus, the value of this integral is determined by the product of the enhanced optical field values at the pump (Ω_{P}) and Raman signal (Ω_{sig} = Ω_{P} − *ω*_{fg}) frequencies and proportional to a complex Gaussian, which is dependent on the detuning of the signal frequency from the ultrafast probe pulse carrier frequency, $\Delta 2v=\omega fg\u2212\Delta PF=\Omega F\u2212\Omega sig$, and the vibrational dephasing rate (*γ*_{fg}). Combining this result with the heterodyne detected signal expression [Eq. (4.11)], the PE-FSRS spectrum for Gaussian pulses when the picosecond pump and femtosecond probe pulses are temporally overlapped (*τ* = 0) yields

For typical vibration dephasing times on the order of picoseconds (*γ*_{fg}) and for a pump and probe carrier frequency detuning nearly resonant with the Raman vibrational frequency Δ_{PF} ≈ *ω*_{fg}, $\Delta 2v\u22480$ and the complex Gaussian appearing in Eq. (4.13) is nearly all real (see Fig. S1). Defining a pump pulse dependent Raman susceptibility by

the *τ* = 0 PE-FSRS vibrational line shape (*ω* ≈ Ω_{sig} = Ω_{P} − *ω*_{fg}) when $\Delta 2v\u22480$ is thus

The real and imaginary parts of *χ*_{R}(*σ*_{P}) are Lorentzian (positive definite) and dispersive-like, respectively, and exhibit only a weak but observable dependence on the picosecond pulse width because the pump pulse duration is the same timescale as vibrational dephasing (see Fig. S2.) Thus, the density matrix description of PE-FSRS analogously finds, as shown for the coupled wave approach,^{45} that, in general, a dispersive line shape can be anticipated as found above for PE-SRG/SRL. It is the complex character of $g(\Delta sig)2$, where Δ_{sig} = Ω_{PL} − (Ω_{P} − *ω*_{fg}), that is primarily responsible for mixed real and imaginary *χ*_{R}(*σ*_{P}) contributions to the PE-FSRS spectrum. More specifically, it the imaginary part of $g(\Delta sig)2$ that accounts for this source of the dispersive PE-FSRS vibrational line shape component [$Im\chi R(\sigma P)$] in the absence of additional competing complex local oscillator fields at the detector.^{52} However, unlike PE-SRG/SRL, when the pump–probe frequency difference is sufficiently off-resonant with the vibrational frequency, the complex Gaussian quantity $e\u2212(\Delta 2v\u2212i\gamma fg)2\sigma F2/2$ offers a phase factor that also contributes to the mix of the real and imaginary *χ*_{R}(*σ*_{P}) components, as seen in Eq. (15), and thus affects the PE-FSRS line shapes, in addition to the complex enhancement factor $g(\Delta sig)2$ as shown below. Fundamentally, this additional phase contribution, evident as well in Eq. (4.11) when *τ* ≠ 0, is due to the signal polarization generated by the additional Fourier components of the stimulating femtosecond probe pulse.

Vibrational line shapes predicted by this PE-FSRS density matrix treatment [Eq. (4.11)] are first illustrated for the optical and materials parameters shown in Fig. 4. The decay time of the modeled plasmon resonance ($\Gamma PL\u22121$) is 10 fs corresponding to a plasmon spectral width (FWHM) of ∼1000 cm^{−1}, and dashed lines show the real and imaginary spectral components of this complex plasmon enhancement factor, *g*(*ω*), in Fig. 4. Other key parameters for the modeled PE-FSRS calculations pictured in Fig. 4 are an incident 40 fs Gaussian probe pulse with a carrier frequency (Ω_{F}) that is 780 cm^{−1} (Δ_{PF}) to the red of a 1 ps Gaussian pump pulse and 500 cm^{−1} to the red of the plasmon resonance peak (Δ_{F} = Ω_{PL} − Ω_{F} = +500 cm^{−1}). These excitation and plasmon parameters (Fig. 4) illustrate one set of parameters used for the calculated PE-FSRS line shapes shown in Fig. 5.

PE-FSRS spectra when pump and probe pulses are *τ* = 0 overlapped for a 1000 cm^{−1} (*ω*_{fg}) Raman active mode with a corresponding 1 ps vibration dephasing time (*γ*_{fg}) [Eq. (4.11)] are shown as a function of the femtosecond pulse detuning from the plasmon resonance (Δ_{F} = 0, ±500, ±1000, ±5000 cm^{−1}) in Fig. 5(a) (Δ_{PF} = Ω_{P} − Ω_{F} = 780 cm^{−1}) and Fig. 5(b) (Δ_{PF} = 1000 cm^{−1}). A 780 cm^{−1} pump–probe detuning (Δ_{PF}) was used in our prior PE-FSRS analysis of Raman line shapes resulting from a plasmon-enhanced local oscillator^{52} and in experiments.^{43} Δ_{PF} = 1000 cm^{−1} [Fig. 5(b)] was chosen because this pump–probe detuning matched the vibrational frequency in this model calculation. For these sets of typical molecular and material parameters, the calculated line shapes exhibit dispersive character resembling the experimental SE-FSRS spectra.^{43,44} As evident in Figs. 5(a) and 5(b), for each of pump–probe carrier detunings (Δ_{PF}), the modeled PE-FSRS line shapes of the 1000 cm^{−1} mode are strongly dependent on the detuning of the incident pulse carrier frequencies (Δ_{F/P}) from the plasmon resonance as anticipated. The sign of the dispersive component of the PE-FSRS line shape generally flips (Fig. 5) as the plasmon resonance (Ω_{PL}) is tuned from one side of the femtosecond probe pulse carrier frequency (±Δ_{F}) to the other (∓Δ_{F}) due to the sign carrying part of $Img(\Delta sig)2=2Reg(\Delta sig)Img(\Delta sig)$ [see Eq. (4.15)].

A comparison of Fig. 5(a) (Δ_{PF} = 780 cm^{−1}) with Fig. 5(b) (Δ_{PF} = 1000 cm^{−1}) shows that a relatively modest (220 cm^{−1}) change in detuning between the pump and probe carrier frequencies, while all other parameters remain the same, can typically result in observable line shape changes. When the pump–probe detuning matches the Raman frequency, i.e., Δ_{PF} = Ω_{P} − Ω_{F} = *ω*_{fg} = 1000 cm^{−1}, and the fs probe carrier frequency is coincident with the plasmon resonant frequency (Δ_{F} = Ω_{PL} − Ω_{F} = 0), the PE-FSRS vibrational feature is entirely negative going [Fig. 5(b), Δ_{F} = 0 panel] because $Img(\Delta sig)2$ = 0 and $Reg(\Delta sig)2<0$ at these resonant conditions for a Lorentzian shaped plasmon resonance. Alternatively, the all imaginary character of the plasmon enhanced femtosecond probe field (Fig. 1, panel on right) can be said to equivalently account for this π signal phase flip.

Furthermore, as indicated above, the complex Gaussian in Eq. (4.13) (*τ* = 0) also contributes a phase factor that will additionally mix real and imaginary *χ*_{R}(*σ*_{P}) when $\Delta 2v=\omega fg\u2212\Delta PF\u22600$ aside from the effects of Im $g(\Delta sig)2$. PE-FSRS line shapes calculated without inclusion of this complex Gaussian and just considering $g(\Delta sig)2$ effects [Eq. (4.15)] are shown in Fig. S3 for the same excitation and optical parameters used in Fig. 5(a) ($\Delta 2v=\omega fg\u2212\Delta PF=220cm\u22121$) and Fig. 5(b) ($\Delta 2v=\omega fg\u2212\Delta PF=0$). Figure S4 shows this same comparison when $\Delta 2v=\omega fg\u2212\Delta PF=500cm\u22121$ and all other parameters are the same. These calculated results (Figs. 5, S3, and S4) demonstrate that the influence of this additional Raman susceptibility phase mixing term is small, albeit evident as indicated above when $\Delta 2v=\omega fg\u2212\Delta PF\u22600$, and has a slightly larger effect on PE-FSRS line shapes as the detuning of the ultrafast probe carrier frequency from the Raman frequency increases. Finally, we note that even for an ultrafast probe pulse detuned from the plasmon resonance by five times the plasmon FWHM spectral width (Δ_{F} = ±5000 cm^{−1}), the effects of the imaginary parts of $g(\Delta sig)2$ can still be detected in the calculated PE-FSRS spectral line shapes [Figs. 5(a), 5(b), S3, and S4].

### C. Three-mode PE-FSRS spectrum

The PE-FSRS line shape dependence on vibrational frequency (*ω*_{fg}) when all other optical pulse and materials parameters are identical is illustrated in Fig. 6 for a molecule with three Raman active modes at 500, 800, and 1100 cm^{−1} that have 1:0.2:1 relative intensities, respectively, in the spontaneous Raman spectrum and 1 ps vibrational dephasing times. The three-mode spontaneous Raman spectrum (black line) and the real and imaginary components of $g(\Delta sig)2$ (dashed red and blue lines) are shown in Fig. 6(a). The calculated corresponding three-mode PE-FSRS spectrum is given in Fig. 6(b) for a pump–probe pulse detuning of 780 cm^{−1} (Δ_{PF}), a 1000 cm^{−1} plasmon spectral width, and a 40 fs probe pulse whose carrier frequency is 200 cm^{−1} to the blue of the plasmon resonance frequency (Δ_{F} = Ω_{PL} − Ω_{F} = −200 cm^{−1}). Thus, for this set of parameters, the plasmon resonance is centered 980 cm^{−1} to the red of the pump pulse (Ω_{PL} = Ω_{P} − (Δ_{PF} − Δ_{F}) = Ω_{P} − 980 [see Fig. 6(a)]. The three Raman features exhibit different line shapes resulting from this set of typical materials and excitation parameters. Such effects are observed in experimental PE-FSRS spectra.^{43,44,46} The 1100 cm^{−1} mode is closest to the plasmon resonance and is nearly all negative going due to the sign of Re $(g(\Delta sig)2)$, and the other two modes have the sign of their dispersive Raman contribution flipped relative to the 1100 cm^{−1} band because they are on opposite sides of the plasmon resonance [see Fig. 6(a)]. Furthermore, the dispersive-like contribution to the 500 cm^{−1} mode line shape is more dominant than for the other two modes because the $Re(g(\Delta sig)2)$ vanishes at 500 cm^{−1} from the peak for this Lorentzian shaped plasmon resonance (FWHH 1000 cm^{−1}). The relative intensities of the PE-FSRS three mode spectrum also differ from the spontaneous spectrum owing to the intensity profile of the femtosecond probe pulse centered in this model at Δ_{PF} = 780 cm^{−1}. While the PE-FSRS line shapes of the three-mode spectrum have been discussed with respect to effects of the complex enhancement factor, the calculated line shapes are also slightly affected by the complex Gaussian quantity in Eq. (4.13) as discussed above.

### D. Pump–probe delay (*τ*) dependence of PE-FSRS line shape

The effect of the pump–probe delay (*τ*) on the shape and intensity of the PE-FSRS has been previously reported^{44} and has been theoretically described for “normal” FSRS via the coupled wave approach.^{57} Calculated FSRS and PE-FSRS spectra are shown in Fig. 7 for a 1000 cm^{−1} mode with a 2 ps dephasing time, 1 ps probe pulse, and 780 cm^{−1} red shifted 40 fs probe pulse. τ > 0 indicates the maximum of picosecond pump arrives at the sample before the femtosecond Stokes probe. As seen in Fig. 7, the FSRS spectrum weakens and broadens as the ultrafast probe arrives increasing *after* the peak in the picosecond probe, i.e., for increasingly positive τ. In contrast, the FRSR spectra narrow and exhibit some negative going wings for increasingly negative τ where the probe *precedes* the temporal peak of the fs probe. Furthermore, the maximum FSRS signal intensity is found for slightly negative τ values, where the overlap of the vibrational coherence decay and the probe pulse are optimal as a function of time. These results are in agreement with the previous coupled wave results describing FSRS and can be understood in terms of the decreased effective decay time for the vibrational coherence when the probe acts after the peak of the ps pulse, whereas the natural vibrational coherence is more fully measured when the probe acts before the pump.^{44,57}

The corresponding calculated PE-FSRS spectra as a function of pump–probe delay (τ) are shown for probe pulse carrier frequencies detuned from the plasmon resonance (Δ_{F}) by +200 and −200 cm^{−1} in the middle and lower panels of Fig. 7, respectively. The calculated τ dependence of these PE-FSRS spectra is not very different from the FSRS effects noted above on the respective band shapes. Spectra get broader and weaker when the pump precedes the probe (*τ* > 0); they appear narrower for *τ* < 0 and peak when the probe pulse arrives slightly before the peak of the probe pulse. Plasmonic enhancement does not alter the *τ* dependent behavior of these spectra significantly, at least for the set of typical parameters considered here. This result follows because the PE ps and fs pulses remain disparate in pulse duration, any additional pulse delay due to the plasmonic interaction is ∼10 fs only or less, and the vibrational dephasing dynamics are not affected, at least in this treatment, by the plasmon resonance.

## V. SPONTANEOUS VS STIMULATED PLASMON-ENHANCED RAMAN SPECTROSCOPIES

A final point considered here is to underscore how the fundamental signal detection processes account for the vibrational line shape differences observed for plasmon-enhanced, heterodyne detected, stimulated Raman spectra and plasmon-enhanced spontaneous Raman (i.e., SERS) spectra in this semiclassical, polarization expansion framework. Formally, spontaneous emission is most properly captured by a quantized harmonic oscillator description of the radiation field.^{66,67} However, as shown previously,^{65} spontaneous emission spectroscopies can be incorporated into nonlinear density matrix treatments by including the incoherent spontaneous field with an amplitude given by the zero point electric field fluctuations, $E\u0302Ssp$ = ($\u210f\Omega S3/$4*π*^{2}*c*^{3})^{1/2}, and the *ad hoc* restriction that this spontaneous field only contributes to emission-like matter–radiation interactions. Thus, for the spontaneous Raman Stokes optical field,

When plasmon-enhanced, this field amplitude is correspondingly $g(\Delta S)E\u0302Ssp$. With this *ansatz*, the resulting third-order, PE-Raman resonant polarization amplitude for the three-level system shown in Fig. 2(a) is the same as given above for PE-SRG (electronically non-resonant, RWA and Condon approximation) with the downward acting spontaneous Stokes field described by

For stimulated Raman spectroscopies, this corresponding third-order polarization response of the Raman active medium [Eq. (5.2)] acts as a source term in Maxwell’s equations and generates a local signal field, which, in turn, is plasmon-enhanced as described above. In contrast, observed spontaneous Raman signals are generated by the incoherent responses of individual scatterers and can be given by the time averaged power, i.e., work done by the Stokes scattered field, Ω_{S}, in the sample per unit time per volume,^{60,65}

Cycle averaging is indicated by the angular brackets, and inelastic scattering corresponds to *W* < 0.^{65} This electromagnetic work term describing spontaneous emission may also be viewed as a consequence of the Poynting theorem.^{59} In contrast to PE-SRL/PE-SRG and PE-FSRS, a spontaneous signal field is not generated by a solution of optical wave equation. All the active fields in the spontaneous Raman process (applied and spontaneous) act locally to generate the observed incoherent emission, and thus, the SERS intensity for a homogeneous vibrational line shape function found by evaluating the work statement above [Eq. (5.3)] is given by

Unlike the heterodyne PE-stimulated Raman spectroscopies, the plasmon-enhanced spontaneous signal is positive definite and proportional to $g(\Delta S)2$ not $g(\Delta S)2$. Thus, due to the fundamentally different signal detection processes of the *same* third-order Raman resonant polarization response of a material, the complex plasmonic enhancement factors can lead to dispersive line shapes in PE-SRG, PE-SRL, and PE-FSRS but not for spontaneous Raman enhanced by plasmonic resonances (SERS).

## VI. CONCLUSIONS

Plasmonic enhancement can significantly improve optical imaging capabilities. A density matrix framework for treating plasmon-enhanced molecular spectroscopies is presented here. The time description of PE active optical electric fields is developed for a metal nanoparticle whose dielectric response is given by a Lorentz oscillator model. Formally, the temporal dependence of the incident fields, plasmon-enhanced or not, is required for a unified quantum mechanical description of nonlinear spectroscopies. A key element for describing plasmon-enhanced stimulated spectroscopies, including those that have Raman resonances, is that the locally generated signal field is also itself enhanced by a plasmonic nanoparticle. Convolution with the plasmonic response function is required to properly treat these stimulated spectroscopies, whether the signal is in the time or frequency domain. This density matrix framework has been applied to a description of PE-SRG/L and PE-CARS resulting from monochromatic radiation. The complex character of the field enhancement factor will contribute to the dispersive vibrational line shapes in agreement with prior reports. However, as first shown in the coupled wave approach treatment of PE-FSRS,^{45} it is the imaginary part of $g(\Delta S)2$/$g(\Delta P)2$, the square of the enhancement factors at the Stokes or pump frequency that determine the extent of the dispersive line shapes for PE-SRG/L respectively, not Im *g*(Δ_{S/P}).^{53} The density matrix approach outlined here for the PE-FSRS description readily incorporates the effects of interpulse (*τ*) delay as shown and will allow descriptions of incident radiation that is not only coincident with plasmon resonances but also resonant with molecular electronic resonances. These may be excitations inherent to the analyte or those resulting from new metal–molecule charge transfer excitations that contribute to chemical enhancement effects. Although both spontaneous and stimulated Raman spectroscopies can be described by a density matrix description of *P*^{(3)}, the different signal detection processes are shown to formally account for the different enhancement factor dependencies for plasmon-enhanced stimulated and spontaneous Raman scattering.

The locally generated coherent signal field is itself plasmon-enhanced as first pointed out in the classical coupled wave description of PE-FSRS,^{45} and the complex character of that enhancement factor primarily leads to dispersive PE-FSRS Raman line shapes as illustrated here as well. However, other mechanisms may additionally contribute to this observed effect for heterodyne detected stimulated spectroscopies when an optically phase distorted local oscillator field is temporally and spatially overlapped with the complex signal field at the detector. For example, such a field resulting from plasmon-enhanced nanoparticle emission as discussed previously^{52} may additionally contribute to the observed mixed real and imaginary Raman vibrational line shape character for PE stimulated Raman spectroscopies, further distorting vibrational line shapes described here.

The observed consequences of the complex character of plasmonic enhancement factors may be partially mitigated by inhomogeneous broadening effects for an ensemble of nanoparticles where size and relaxation rate distributions are an experimental reality. The effects of inhomogeneities, as well as more realistic models for dispersion of the dielectric constant of metal materials, will be considered in a subsequent report. Such models will include the effects of multiple contributions to the material dielectric, enhancement factors given by FDTD solutions to Maxwell’s equations incorporating experimental dielectric measurements, as well as considerations for enhanced fields of different nanostructures, i.e., nanoparticle dimers.

Finally, while the dispersion effects due to plasmonic enhancement on the amplitude and phase of the local optical fields made explicit here are required for the description of PE spectra via a density matrix approach, it is especially crucial for treatments of plasmon-enhanced ultrafast spectroscopies employing pulses capable of impulsive vibrational excitation and spectroscopies where the observed nonlinear signals are temporal responses. This will be explicitly treated in a subsequent publication.

## SUPPLEMENTARY MATERIAL

The supplementary material includes derivations of the third-order optical polarization responses corresponding to cw PE-SRG, PE-SRL, PE-CARS, and PE-FSRS, the analytical expression for the temporal dependence of plasmon-enhanced incident Gaussian ultrafast pulses, and an explicit form for the $F*(\Delta 2v,0)$ integral. Plots of the vibrational relaxation time dependence of the complex Gaussian in the PE-FSRS intensity expression, the picosecond pump pulse width dependence of the Raman susceptibility integral [*χ*_{R}(*σ*_{P})], and PE-FSRS spectral effects due to $g(\Delta sig)2$ alone are also shown in Figs. S1–S4.

## ACKNOWLEDGMENTS

The support of NSF Grant No. CHE 1609952 is gratefully acknowledged.

## AUTHOR DECLARATIONS

### Conflict of Interest

The authors have no conflicts of interest to disclose.

## DATA AVAILABILITY

The data that support the findings of this study are available within the article and its supplementary material.

## REFERENCES

^{5}–10

^{6}enhancement