Coherent anti-Stokes Raman Spectroscopy of single and multi-layer graphene

TL;DR

This paper demonstrates the use of Coherent anti-Stokes Raman Spectroscopy (CARS) to analyze the vibrational properties of single and multi-layer graphene, highlighting methods to suppress background noise for clearer spectral measurements.

Contribution

It introduces a technique to measure the vibrationally resonant CARS peak in graphene by reducing pulse overlap, improving vibrational sensitivity in imaging.

Findings

CARS signal dominated by non-resonant background obscures Raman lineshape.

Reducing laser pulse overlap suppresses the non-resonant background.

CARS can be effectively used for vibrational imaging of graphene.

Abstract

We report stimulated Raman spectroscopy of the G phonon in both single and multi-layer graphene, through Coherent anti-Stokes Raman Scattering (CARS). The signal generated by the third order nonlinearity is dominated by a vibrationally non-resonant background (NVRB), which obscures the Raman lineshape. We demonstrate that the vibrationally resonant CARS peak can be measured by reducing the temporal overlap of the laser excitation pulses, suppressing the NVRB. We model the observed spectra, taking into account the electronically resonant nature of both CARS and NVRB. We show that CARS can be used for graphene imaging with vibrational sensitivity.

Figures

Equations

$$I({\omega }_{as})\propto |P_{CARS}^{(3)}({\omega }_{as})+P_{NVRB}^{(3)}({\omega }_{as}){|}^2.$$

$$\begin{gathered} P_{CARS}^{(3)}(\omega ,\Delta T)\propto -{\eta }_{CARS} \\ {\int }_{-\infty }^{\infty }d{\omega }_3{\int }_{-\infty }^{\infty }d{\omega }_2{\int }_{-\infty }^{\infty }d{\omega }_1\frac{{\hat{A}}_P({\omega }_3,\Delta T){\hat{A}}_P({\omega }_1,\Delta T){\hat{A}}_S^{\ast }({\omega }_2,0)\delta (\omega -2{\omega }_P+{\omega }_S-{\omega }_3-{\omega }_1+{\omega }_2)}{\left({\omega }_P+{\omega }_3-{\bar{\omega }}_{ba}\right)\left({\omega }_P-{\omega }_S+{\omega }_3-{\omega }_2-{\bar{\omega }}_{ca}\right)\left(2{\omega }_P-{\omega }_S+{\omega }_3-{\omega }_2+{\omega }_1-{\bar{\omega }}_{da}\right)} \end{gathered}$$

$$\begin{gathered} P_{CARS}^{(3)}(\omega ,\Delta T)\propto -{\eta }_{CARS}{\int }_{-\infty }^{\infty }d{\omega }_1{\int }_{-\infty }^{\infty }d{\omega }_3 \\ \frac{{\hat{A}}_P({\omega }_3,\Delta T){\hat{A}}_P({\omega }_1,\Delta T){\hat{A}}_S^{\ast }(2{\omega }_P-{\omega }_S-\omega +{\omega }_3+{\omega }_1,0)}{\left({\omega }_P+{\omega }_3-{\bar{\omega }}_{ba}\right)\left(\omega -{\omega }_P-{\omega }_1-{\bar{\omega }}_{ca}\right)\left(\omega -{\bar{\omega }}_{da}\right)} \end{gathered}$$

$$\begin{gathered} P_{NVRB}^{(3)}(\omega ,\Delta T)\propto -{\eta }_{NVRB}{\int }_{-\infty }^{\infty }d{\omega }_1{\int }_{-\infty }^{\infty }d{\omega }_2 \\ \frac{{\hat{A}}_P({\omega }_1,\Delta T){\hat{A}}_P({\omega }_2,\Delta T){\hat{A}}_S^{\ast }(2{\omega }_P-{\omega }_S-\omega +{\omega }_1+{\omega }_2,0)}{\left({\omega }_P+{\omega }_1-{\bar{\omega }}_{ea}\right)\left(2{\omega }_P+{\omega }_1+{\omega }_2-{\bar{\omega }}_{ea}\right)\left(\omega -{\bar{\omega }}_{ea}\right)} \end{gathered}$$

$$\begin{array}{rl}{P}_{CARS}^{(3)}(\omega )& \propto -\frac{{\eta }_{CARS}{E}_P^2{E}_S^{\ast }}{({\omega }_P-{\bar{\omega }}_{ba})(\omega -{\omega }_P-{\bar{\omega }}_{ca})(\omega -{\bar{\omega }}_{da})}=\\ & ={\chi }_{CARS}^{(3)}{E}_P^2{E}_S^{\ast }\end{array}$$

$$\begin{array}{rl}{P}_{NVRB}^{(3)}(\omega )& \propto -\frac{{\eta }_{NVRB}{E}_P^2{E}_S^{\ast }}{({\omega }_P-{\bar{\omega }}_{ea})(2{\omega }_P-{\bar{\omega }}_{ea})(\omega -{\bar{\omega }}_{ea})}=\\ & ={\chi }_{NVRB}^{(3)}{E}_P^2{E}_S^{\ast }\end{array}$$

$$\begin{array}{rl}I({\omega }_{as}{)}^{NR}& \sim |P_{NVRB}^{(3)}{|}^2+|P_{CARS}^{(3)}{|}^2+2\Re (P_{NVRB}^{(3)})\Re (P_{CARS}^{(3)})\\ & \propto |{\chi }_{NVRB}^{(3)}{|}^2+|{\chi }_{CARS}^{(3)}{|}^2+2\Re ({\chi }_{NVRB}^{(3)})\Re ({\chi }_{CARS}^{(3)})\end{array}$$

$$\begin{array}{rl}I({\omega }_{as}{)}^R& \sim |P_{NVRB}^{(3)}{|}^2+|P_{CARS}^{(3)}{|}^2+2\Im (P_{NVRB}^{(3)})\Im (P_{CARS}^{(3)})\\ & \propto |{\chi }_{NVRB}^{(3)}{|}^2+|{\chi }_{CARS}^{(3)}{|}^2+2\Im ({\chi }_{NVRB}^{(3)})\Im ({\chi }_{CARS}^{(3)})\end{array}$$

$$I=I_1-I_2+\frac{{\tilde{\nu }}_2-{\tilde{\nu }}_1}{{\tilde{\nu }}_3-{\tilde{\nu }}_1}(I_3-I_1)$$

$$\begin{gathered} P^{(3)}(t)\propto N{\int }_0^{\infty }d{\tau }_3{\int }_0^{\infty }d{\tau }_2{\int }_0^{\infty }d{\tau }_1\mathcal{E}(t-{\tau }_3) \\ \mathcal{E}(t-{\tau }_2-{\tau }_3)\mathcal{E}(t-{\tau }_1-{\tau }_2-{\tau }_3)S^{(3)}({\tau }_1,{\tau }_2,{\tau }_3) \end{gathered}$$

$$S^{(3)}(\tau_{1},\tau_{2},\tau_{3})\propto\left(i\right)^{3}Tr\Big{\{}\mu(\tau_{1}+\tau_{2}+\tau_{3})\\ \big{[}\mu(\tau_{1}+\tau_{2}),\big{[}\mu(\tau_{1}),\big{[}\mu(0),\rho(-\infty)\big{]}\big{]}\big{]}\Big{\}}$$

$$\mathcal{E}(t)=\sum_{i=P,S}\big{[}E_{i}(t,\Delta t_{i})+c.c.\big{]}=\sum_{i=P,S}\big{[}\hat{A}_{i}(t,\Delta_{i}t)e^{-i\omega_{i}t}+c.c.\big{]}$$

$$\begin{gathered} P_{CARS}^{(3)}(t)\propto {\left(i\right)}^3{n}_{CARS}{\mu }_{ba}{\mu }_{cb}{\mu }_{cd}{\mu }_{ad}{\int }_0^{\infty }d{\tau }_3{\int }_0^{\infty }d{\tau }_2{\int }_0^{\infty }d{\tau }_1{A}_P(t-{\tau }_1-{\tau }_2-{\tau }_3,\Delta t)A_S^{\ast }(t-{\tau }_2-{\tau }_3)A_P(t-{\tau }_3,\Delta t) \\ {e}^{-i{\omega }_P(t-{\tau }_1-{\tau }_2-{\tau }_3,\Delta t)}{e}^{+i{\omega }_S(t-{\tau }_2-{\tau }_3)}{e}^{-i{\omega }_P(t-{\tau }_3,\Delta t)}{e}^{-i{\bar{\omega }}_{ba}{\tau }_1}{e}^{-i{\bar{\omega }}_{ca}{\tau }_2}{e}^{-i{\bar{\omega }}_{da}{\tau }_3} \end{gathered}$$

$$\begin{gathered} P_{NVRB}^{(3)}(t)\propto {\left(i\right)}^3{n}_{NVRB}|{\mu }_{ea}{|}^4{\int }_0^{\infty }d{\tau }_3{\int }_0^{\infty }d{\tau }_2{\int }_0^{\infty }d{\tau }_1{A}_P(t-{\tau }_1-{\tau }_2-{\tau }_3)A_P(t-{\tau }_2-{\tau }_3)A_S^{\ast }(t-{\tau }_3) \\ {e}^{-i{\omega }_P(t-{\tau }_1-{\tau }_2-{\tau }_3,\Delta t)}{e}^{-i{\omega }_P(t-{\tau }_2-{\tau }_3,\Delta t)}{e}^{+i{\omega }_S(t-{\tau }_3)}{e}^{-i{\bar{\omega }}_{ea}{\tau }_1}{e}^{-i{\bar{\omega }}_{ea}{\tau }_2}{e}^{-i{\bar{\omega }}_{ea}{\tau }_3} \end{gathered}$$

$$P^{(3)}(\omega )={\int }_{-\infty }^{\infty }{P}^{(3)}(t)e^{i\omega t}dt.$$

$$\begin{gathered} P_{CARS}^{(3)}(\omega )\propto {\eta }_{CARS}{\left(i\right)}^3{\int }_{-\infty }^{\infty }dt{e}^{i\omega t}{\int }_0^{\infty }d{\tau }_3{\int }_0^{\infty }d{\tau }_2{\int }_0^{\infty } \\ d{\tau }_1{\int }_{-\infty }^{\infty }d{\omega }_1{\int }_{-\infty }^{\infty }d{\omega }_2{\int }_{-\infty }^{\infty }d{\omega }_3{\hat{A}}_P({\omega }_1,\Delta t){\hat{A}}_S^{\ast }({\omega }_2,0) \\ {\hat{A}}_P({\omega }_3,\Delta t)e^{-i({\omega }_P+{\omega }_1)(t-{\tau }_1-{\tau }_2-{\tau }_3)}{e}^{+i({\omega }_S+{\omega }_2)(t-{\tau }_2-{\tau }_3)} \\ {e}^{-i({\omega }_P+{\omega }_3)(t-{\tau }_3)}{e}^{-i{\bar{\omega }}_{ba}{\tau }_1}{e}^{-i{\bar{\omega }}_{ca}{\tau }_2}{e}^{-i{\bar{\omega }}_{da}{\tau }_3} \end{gathered}$$

$$\begin{gathered} P_{CARS}^{(3)}(\omega )\propto -{\eta }_{CARS}{\int }_{-\infty }^{\infty }d{\omega }_1{\int }_{-\infty }^{\infty }d{\omega }_2{\int }_{-\infty }^{\infty }d{\omega }_3\delta (\omega -2{\omega }_P+{\omega }_S+{\omega }_1-{\omega }_2+{\omega }_3) \\ \frac{{\hat{A}}_P({\omega }_1,\Delta t){\hat{A}}_S^{\ast }({\omega }_2,0){\hat{A}}_P({\omega }_3,\Delta t)}{({\omega }_P+{\omega }_1-{\bar{\omega }}_{ba})({\omega }_P-{\omega }_S+{\omega }_1-{\omega }_2-{\bar{\omega }}_{ca})(2{\omega }_P-{\omega }_S+{\omega }_1-{\omega }_2+{\omega }_3-{\bar{\omega }}_{da})} \end{gathered}$$

$$\delta (\omega -2{\omega }_P+{\omega }_S-{\omega }_1+{\omega }_2-{\omega }_3)={\int }_{-\infty }^{\infty }{e}^{i(\omega -2{\omega }_P+{\omega }_S-{\omega }_1+{\omega }_2-{\omega }_3)t},$$

$$\begin{gathered} P_{CARS}^{(3)}(\omega ,\Delta t)\propto -{\eta }_{CARS}{\int }_{-\infty }^{\infty }d{\omega }_1{\int }_{-\infty }^{\infty }d{\omega }_3 \\ \frac{{\hat{A}}_P({\omega }_3,\Delta t){\hat{A}}_P({\omega }_1,\Delta t){\hat{A}}_S^{\ast }(2{\omega }_P-{\omega }_S-\omega +{\omega }_3+{\omega }_1,0)}{\left({\omega }_P+{\omega }_3-{\bar{\omega }}_{ba}\right)\left(\omega -{\omega }_P-{\omega }_1-{\bar{\omega }}_{ca}\right)\left(\omega -{\bar{\omega }}_{da}\right)} \end{gathered}$$

$$\begin{gathered} P_{NVRB}^{(3)}(\omega ,\Delta t)\propto -{\eta }_{NVRB}{\int }_{-\infty }^{\infty }d{\omega }_1{\int }_{-\infty }^{\infty }d{\omega }_2 \\ \frac{{\hat{A}}_P({\omega }_1,\Delta t){\hat{A}}_P({\omega }_2,\Delta t){\hat{A}}_S^{\ast }(2{\omega }_P-{\omega }_S-\omega +{\omega }_1+{\omega }_2,0)}{\left({\omega }_P+{\omega }_1-{\bar{\omega }}_{ea}\right)\left(2{\omega }_P+{\omega }_1+{\omega }_2-{\bar{\omega }}_{ea}\right)\left(\omega -{\bar{\omega }}_{ea}\right)} \end{gathered}$$

Full text

Single-layer graphene (SLG) has a high nonlinear third-order susceptibility: $\lvert\chi^{(3)}\rvert$ $\sim 10^{-10}$ e.s.u for harmonic generationSoavi2018 and $\lvert\chi^{(3)}\rvert$ $\sim 10^{-7}$ e.s.u. for frequency mixingHendry2010 , where one electrostatic unit of charge (1 e.s.u), in standard units (SI) isBoydBook : $\chi^{(3)}$(SI)/$\chi^{(3)}$(e.s.u.) = $4\pi/(3\times 10^{4})^{2}$. This is up to seven orders of magnitude greater than those of dielectric materials such as silica ($\chi^{(3)}=1.4\times 10^{-14}$ e.s.uPotmaMukamel2013 ). This property is due to optical resonance with interband electronic transitionsKhurgin2014 and has led to the observation of gate-tunable third-harmonic generationSoavi2018 and nonlinear four-wave mixingHendry2010 ; Wu2011 ; Gu2012 (FWM, i.e. the third-order processes whereby an electromagnetic field is emitted by the nonlinear polarization induced by three field-matter interactionsBoydBook ). FWM can be exploited for for graphene imaging, with an image contrast of up to seven orders of magnitudeHendry2010 higher than that of optical reflection microscopyduong_probing_2012 . However, FWM-based imaging reported to date in grapheneHendry2010 lacks chemical selectivity and does not provide the same wealth of information brought about by the vibrational sensitivity of Raman spectroscopyFerrariPRL2006 ; FerrariReview2013 .

Coherent anti-Stokes Raman scattering (CARS)Begley1974 ; Duncan1982 ; Zumbusch1999 ; Cheng2002 is a FWM process that exploits the nonlinear interaction of two laser beams, the pump field $E_{P}$ at frequency $\omega_{P}$ and the Stokes field $E_{S}$ at frequency $\omega_{S}<\omega_{P}$, to access the vibrational properties of a material. As shown in Fig.1a, when the energy difference between the two photons matches a phonon energy ($\hbar\omega_{P}-\hbar\omega_{S}=\hbar\omega_{v}$), the interaction of the laser pulses and the sample results in the generation of vibrational coherencesPotmaMukamel2013 . While spontaneous Raman (SR) scattering is an incoherent signalMukamel1999 , since the phases of the electromagnetic fields emitted by individual scatterers are uncorrelatedMukamel1999 , in CARS, atomic vibrations are coherently stimulated, i.e. atoms oscillate with the same phasePotmaMukamel2013 , potentially leading to a signal enhancement of several orders of magnitude depending on incident power and scatterer densityCui2009 ; Pestov265 .

The same combination of optical fields used for CARS can generate another FWM signal, a non-vibrationally resonant background (NVRB)Hendry2010 , Fig.1b. In both processes, the optical response consists of a field emitted at the anti-Stokes frequency $\omega_{as}=2\omega_{P}-\omega_{s}$PotmaMukamel2013 . However, the interference of the two effects usually generates an an additional contribution which is dispersive with respect to the emitted optical frequency, i.e. shaped as the first derivative of a peaked function (resembling the real part of the refractive index around a resonance), which introduces an asymmetric distortion of the Raman peak profile in the region $\omega_{as}=\omega_{P}+\omega_{v}$ evans_coherent_2008 .

In the biological fieldCui2009 ; Evans16807 , a wealth of studies has demonstrated the potential of CARS for fast imagingCui2009 ; Pestov265 ; CiceroneNat2015 , with pixel acquisition times as low as$\sim 0.16\mu s$Evans16807 , thus allowing for video-rate microscopyEvans16807 . By contrast, there are only a few reports to date of CARS imaging of micro-structured materials (such as polyethylene blendLee2012 , multicomponent polymersLim2006 , cholesterol micro-crystalsLim2011 ) and nano-structured ones (patterned gold surfacesSteuwe2011 , single wall nanotubesIkeda2009 ; Lee2010 , highly oriented pyrolytic graphiteDovbeshko2010 ). Such studies, performed with pixel acquisition times down to$\sim 2\mu s$Baldacchini2010 , have shown the ability of CARS to identify chemical heterogeneities on sub-micrometer scales and characterize single particles that are part of a larger domain, enabling e.g. to visualize microscopic domains (polystyrene, poly-methyl methacrylate, and poly-ethylene terephthalate) in the case of the above mentioned polymer mixturesPotma2013 , or to provide detailed maps of microcrystal orientation in organic matrices (e.g. cholesterols in atherosclerotic plaquesLim2011 ).

In graphene, despite the large $\chi^{(3)}$Hendry2010 ; Soavi2018 , no CARS peak profiles, equivalent to those measured in SR, have been observed to date, to the best of our knowledge. We reported SR with single-color pulsed excitationferrante_nc , using the same picosecond lasers usually adopted for CARSEvans16807 . However, in order to measure CARS, a combination of pulses with different colors must be adoptedTRCARS2017 . By scanning the pulse frequency detuning in a two-color experiment, a dip has been observed in the third order nonlinear spectral response of SLG at the G phonon frequency. This was interpreted as an anomalous antiresonance and phenomenologically described in terms of a Fano lineshapeLafetaArx2017 .

Here we use two 1ps pulses (see inset of Fig.2) to explore FWM in SLG and few-layer graphene (FLG). We experimentally demonstrate and theoretically describe how the inter-pulse delay, $\Delta T$ (Figs.1c-d) can be used to modify the relative weight of CARS and NVRB components that simultaneously contribute to the FWM, thus recovering the G-band Raman peak profile. We show that the dip in the nonlinear optical response around the vibrational resonance is due to the interplay of CARS and NVRB under electronically resonant conditions, which allows vibrational imaging with signal levels as large as those of the third-order nonlinear response.

SLG is grown on a 35$\mu$m copper (Cu) foil, following the process described in Ref.33. The substrate is heated up to 1000*∘C and annealed in hydrogen atmosphere (H2, 20 sccm) for 30 minutes at $\sim$200 mTorr. Then, 5 sccm of methane (CH4) are let into the chamber for the following 30 minutes to enable growthBaeNN2010 ; LiS2009 . The sample is then cooled back to room temperature in vacuum ($\sim$1 mTorr) and unloaded from the chamber. SLG is subsequently transferred on a glass substrate by a wet method. Poly-methyl methacrylate (PMMA) is spin coated on the SLG/Cu and floated on a solution of ammonium persulfate (APS) and deionized water. When Cu is etchedBaeNN2010 ; BonaMT2012 , the PMMA membrane with attached SLG is then moved to a beaker with deionized water for cleaning APS residuals. The membrane is subsequently lifted with the target substrate. After drying, PMMA is removed in acetone leaving SLG on glass. SLG is characterized by SR after transfer using a Renishaw InVia spectrometer at 514nm. The position of the G peak, Pos(G), is$\sim$1582cm-1*, while FWHM(G)$\sim$14cm*-1*. The 2D to G peak area ratio is$\sim$5.3, indicating a p-doping after transfer$\sim$250meVBaskPRB2009 ; DasNN2008 , which corresponds to a carrier concentration$\sim$5$\cdot$1012cm*-2*. FLG flakes are produced by micromechanical cleavage from bulk graphiteNovoselov2005PNAS . The bulk crystal is exfoliated on Nitto Denko tape. The FLG G peak is$\sim$1580cm*-1*. The D peak is negligible. The 2D peak shape indicates this is Bernal stacked FLGFerrariPRL2006 ; FerrariReview2013 .

For CARS experiments, we use a two-modules Toptica FemtoFiber Pro source, with two Erbium-Doped fiber amplifiers (EDFA) at$\sim 1550$nm generating 90fs pulses at 40MHz, seeded by a common mode-locked Er:fiber oscillatorKrauss:09 , Fig.2. In the first branch (FemtoFiber pro NIR), 1ps pulses at 784nm (pump pulse, PP) are produced by second-harmonic spectral compressionMarangoni:09 in a 1cm Periodically Poled Lithium Niobate (PPLN) crystal. In the second branch (FemtoFiber pro TNIR), the amplified laser passes through a nonlinear fiber (NLF), wherein a supercontinuum (SC) output is generated. The SC spectral intensity can be tuned with a motorized Si-prism-pair compressor. A PPLN crystal with a fan-out grating (a poling period changing along the transverse direction) is exploited to produce broadly tunable (from 840 to 1100nm) narrowband 1ps Stokes pulses (SP), with a power$<10$mW. A dichroic mirror (DM) is used to combine the two beams, whose relative temporal delay is tuned with an optical delay line (DL). A long-working distance 20x objective (O, numerical aperture NA=0.4) focuses the pulses onto the sample (S). The generated and transmitted light is collected by a condenser (C) and the PP and SP are filtered out by a short-pass filter (F). The total FWM signal is collected with an optical multichannel analyser (OMA, Photon Control SPM-002-E). A dichroic mirror reflects the SP in order to measure its intensity ($I_{s}$) with a powermeter (P). FWM spectra are obtained by scanning $\omega_{S}$ around $\omega_{P}-\omega_{v}$ (at fixed $\omega_{P}$) to probe the G-band phonon frequency, $\omega_{v}=\omega_{G}$. Fig.3 displays the FWM intensity, normalized to $I_{s}$, for different $\Delta T$.

In both SLG and FLG measurements, for $\Delta T$ shorter than the vibrational dephasing time $\tau\sim$1psBonini2007 , i.e. the characteristic time of coherence lossMukamel1999 , a Lorentzian dip at $\omega_{as}=\omega_{P}+\omega_{G}$ appears on top of a backgroundLafetaArx2017 . For $\Delta T>2$ps, while the total FWM signal decreases by nearly two orders of magnitude, the dip observed in FLG at $\Delta T\sim 0$ps evolves into a Raman peak shape at the G-phonon energy. No dispersive features are observed at any $\Delta T$, unlike what normally expected for the interference between NVRB and CARSevans_coherent_2008 . Here we use pulses with duration $\delta t\sim 1$ps since this allows us to scan the inter-pulse delay across the vibrational dephasing time $\tau$ to suppress the NVRB cross section more than the vibrational contribution, while minimizing the spectral broadening due to the finite pulse duration $1/\delta t\simeq$15cm*-1* Mukamel1999 .

Since both CARS and NVRB signals depend quadratically on the number of scatterersPotma2013 , the SLG signal intensity is significantly reduced with respect to FLG (Fig.3), with a lower signal-to-noise ratio, hampering the detection of peak-shaped vibrational resonances expected for $\Delta T>1.4$ps.

The data in Fig.3 can be qualitatively understood as follows. The anti-Stokes signal, $I(\omega_{as})$, generated by CARS and NVRB, is proportional to the square modulus of the electric field emitted by the third-order polarization, $P^{(3)}$PotmaMukamel2013 , as:

[TABLE]

CARS and NVRB signals are simultaneously generated by two light-matter interactions with the PP and a single interaction with the SP, with different time ordering, Figs.1a-b. Consequently, $P^{(3)}\propto E_{P}^{2}E_{S}^{*}$, where * indicates the complex conjugate. However, the temporal constraints for such interactions are significantly different for the two cases. As shown in Figs.1c-d, in the case of NVRB the three interactions must take place within the dephasing time of the involved electronic excitation, which in SLG is$\sim$10fsTomadin2013 ; Brida2013 , i.e. much shorter than the pulses duration ($\delta t\sim$1ps). Hence, $P^{(3)}_{NVRB}(\omega_{as})$ is only generated during the temporal overlap between the two pulses $P^{(3)}_{NVRB}\propto E_{P}^{2}(t-\Delta T)E_{S}^{*}(t)$ (the three field interactions, in a representative NVRB event, are indicated by three nearly coincident dots in Fig.1d). In CARS, the electronic dephasing time only constrains the lag between the first two interactions that generate the vibrational coherence (the two stimulating-field interactions are represented by the two nearly coincident dots in Fig.1c). This can be read out by the third field interaction within the phonon dephasing time, $\tau\sim 1$psBonini2007 (indicated, for a representative CARS event, by the third dot in Fig.1d). Thus, $P^{(3)}_{CARS}\propto E_{P}(t-\Delta T)E_{S}^{*}(t)\int_{-\infty}^{t}E_{P}(t^{\prime}-\Delta T)e^{-t^{\prime}/\tau}dt^{\prime}$Volkmer2002 . Therefore, $\Delta T$ can be used to control the relative weights of $P^{(3)}_{CARS}$ and $P^{(3)}_{NVRB}$Kamga:80 ; Volkmer2002 ; Sidorov-Biryukov:06 ; Pestov265 ; Selm:10 ; Kumar:11 . For positive time delays within a few $\tau$, $P^{(3)}_{CARS}/P^{(3)}_{NVRB}$ is progressively enhanced, as shown in Fig.4a,b,c.

The system response can be evaluated through a density-matrix description of $P^{(3)}(\omega,\Delta T)$Mukamel1999 . In the presence of a temporal delay between PP and SP, their electric fields can be written asBoydBook : $E_{P}(t,\Delta T)={A}_{P}(t,\Delta T)e^{-i\omega_{P}t}$ and $E_{S}(t,0)={A}_{S}(t,0)e^{-i\omega_{S}t}$, where ${A}_{P/S}(t,\Delta T)$ indicates the PP/SP temporal envelope. By Fourier transform, the fields can be expressed in the frequency domain as: $\hat{E}_{P}(\omega,\Delta T)=\int_{-\infty}^{+\infty}{E}_{P}(t,\Delta T)e^{i\omega t}dt$ and $\hat{E}_{S}(\omega,0)=\int_{-\infty}^{+\infty}{E}_{S}(t,0)e^{i\omega t}dt$, which can be used to calculate $P^{(3)}_{CARS}(\omega,\Delta T)$ asMukamel1999 ; Batignani2016 :

[TABLE]

where $\eta_{\hbox{\footnotesize{CARS}}}=n_{\hbox{\footnotesize{CARS}}}\mu_{ba}\mu_{cb}\mu_{cd}\mu_{ad}$, $\mu_{ij}$ is the transition dipole moment between the $i$ and $j$ states, $n_{\hbox{\footnotesize{CARS}}}$ is the number of scatterers involved in the CARS process, $\bar{\omega}_{ij}=\omega_{ij}-i\gamma_{ij}=\omega_{i}-\omega_{j}-i\gamma_{ij}$, $\omega_{ij}=\omega_{i}-\omega_{j}$ is the energy difference between the levels i and j, and $\gamma_{ij}=\tau_{ij}^{-1}$ is the dephasing rate of the $\ket{i}\bra{j}$ coherenceMukamel1999 ; $a$ and $c$ denote the vibrational ground state $|g,0\rangle$, and the first vibrational excited level, $|g,1\rangle$, with respect to the electronic ground state $|g\rangle$ ($\pi$ band). In our experiments, $c$ corresponds to the G phonon at $q\sim 0$, $b$ and $d$ indicate the vibrational ground state, $|e,0\rangle$, and the first vibrational excited level, $|e,1\rangle$, with respect to the excited electronic state $|e\rangle$ ($\pi^{*}$ band).

Using the conservation of energy represented by the $\delta$-distribution in Eq.2 and integrating over $\omega_{2}$, we get:

[TABLE]

Defining $\tilde{\nu}=\omega/(2\pi c)$, the third-order nonlinear polarization can be expressed as a function of the Raman shift ($\tilde{\nu}-\tilde{\nu}_{P}$) as $P^{(3)}(\omega,\Delta T)=P^{(3)}(2\pi c\tilde{\nu},\Delta T)$.

The ${\omega}_{ca}$ level in the denominator of Eq.3 is the frequency of the Raman mode coherently stimulated in CARS, while ${\omega}_{ba}$ and ${\omega}_{da}$ are frequency differences between the electronic levels. In the case of real levels, resonance enhancement occursMukamel1999 . In view of the optical nature of the involved phonons ($q\sim 0$), and due to momentum conservation, only one electronic level must be included in the calculation and, consequently, the nonlinear response can be derived for $\omega_{ba}=\omega_{dc}=\omega_{P}$. In a similar manner, $P^{(3)}_{NVRB}$ can be expressed asMukamel1999 :

[TABLE]

where $\eta_{\hbox{\footnotesize{NVRB}}}=n_{\hbox{\footnotesize{NVRB}}}|\mu_{ea}|^{4}$, $n_{\hbox{\footnotesize{NVRB}}}$ is the number of scatterers involved in the NVRB process, and $\omega_{ea}$ is the energy of the electronic excited level involved in the NVRB process. Since the cross section of third-order nonlinear processes in graphene is enhanced by increasing the photon energyWang2013146 ; concado2007 , we consider only the dominant case, i.e., $\tilde{\nu}_{ea}=2\tilde{\nu}_{P}$.

We describe the spectral FWM response assuming monochromatic fields with no inter-pulse delay: $\hat{E}_{P}(\omega)=E_{P}\cdot\delta(\omega-\omega_{P})$, $\hat{E}_{S}(\omega)=E_{S}\cdot\delta(\omega-\omega_{S})$. From Eqs.3,4, the CARS and NVRB nonlinear polarizations can be expressed asMukamel1999 ; PotmaMukamel2013 :

[TABLE]

[TABLE]

which can be used to calculate the total FWM spectrum according to Eq.1. Fig.4a plots the electronically non-resonant case. The CARS polarization, defined by Eq.5, is a complex quantity: the real part has a dispersive lineshape, while the imaginary part peaks at the phonon frequency $\omega_{ca}$. The NVRB polarization, defined by Eq.6, is a positive real quantity. Accordingly, the FWM spectrum in the electronically non-resonant condition, $I(\omega_{as})^{NR}$, can be written asMukamel1999 :

[TABLE]

and it can be significantly distorted by the third term in Eq.7 depending on the relative weight of the two corresponding susceptibilities. The maximum of the signal, when the dispersive contribution is dominant, can be frequency shifted from the phonon frequency. This is the most common scenario, in which the dispersive lineshapes hampering a direct access to the vibrational characterization of the sample in terms of phonon frequencies and lifetimes. Such limitation is particularly severe when $\chi^{(3)}_{NRVB}$ is comparable to $\chi^{(3)}_{CARS}$ and the NVRB and CARS contributions have the same order of magnitude. This condition is common in the case of a weak vibrational resonant contribution ($\frac{\mu_{ba}\mu_{cb}\mu_{cd}\mu_{ad}}{|\mu_{ea}|^{4}}<<1$), as in the case of low concentrations of oscillators ($\frac{n_{\hbox{\footnotesize{CARS}}}}{n_{\hbox{\footnotesize{NVRB}}}}<<1$). Hence, this produces an intense NVRB signal and reduces the vibrational contrast, hindering the imaging of electronically non-resonant samples. This is the case for cells and tissues which need to be excited in the near infrared to avoid radiation damageMuller2007 .

For SLG, the linear dispersion of the massless Dirac Fermions makes the response always electronically resonant. In the case of FLG, absorption has a complex dependence on wavelength, as well as on the number of layers and their relative orientation, exhibiting, for instance, a tunable band gap in twisted bilayer grapheneZhang2009 . This is also reflected in the resonant nature of SRWu2014 ; Wu2015 . However, approaching visible wavelengths, the absorption spectrum flattens above$\sim 0.8$eV and it is $\sim(1-\pi e^{2}/2h)^{N}$ for Bernal-stacked $N$-layer grapheneMak2010 . Our exfoliated FLG are Bernal stacked, as also confirmed by the measured SR 2D peak shape in SRFerrariPRL2006 ; FerrariReview2013 . Accordingly, at the typical CARS wavelengths used here (784 and 894nm), SLG and Bernal FLG are electronically resonant, unlike the situation for most biological samplesMuller2007 . This results in an opposite sign in the CARS response, i.e. a spectral dip in $\Im(\chi^{(3)}_{CARS})$, related to two additional imaginary unit contributions in the denominator of Eq.5, $(\omega_{P}-\bar{\omega}_{ba})$ and $(\omega-\bar{\omega}_{da})$, wherein the $i\gamma_{ba}$, $i\gamma_{da}$ components dominate. Further, the $-i$ contribution from $(2\omega_{P}-\bar{\omega}_{ea})$ results in a NVRB dominated by the imaginary part, as illustrated in Fig.4c,d,e.

Thus, the third term in Eq.7 must be replaced with the contribution from the interference of the spectral dip $\Im(\chi^{(3)}_{CARS})$ with the imaginary part $\Im(\chi^{(3)}_{NVRB})$. This leads to a signal that, under the electronically resonant regime, becomesMukamel1999 :

[TABLE]

which indicates that the total FWM, at the phonon frequency, can be either a negative dip or a positive peak depending on the ratio between the vibrationally resonant and the non-resonant susceptibilities ($\chi^{(3)}_{CARS}/\chi^{(3)}_{NVRB}$), which depends only on the material under examination and not on the pulses used in the experiment. Such a qualitatively different interplay between NVRB and CARS, compared with the experimental lineshapes for $\Delta T=0$ in Fig.3, unambiguously indicates the presence of electronic resonance in SLG and Bernal FLG. For a given material, the relative weight of the two FWM contributions can be in general modified by using pulsed excitation and tuning the temporal overlap between the PP and SP fieldsVolkmer2002 , i.e. changing $\Delta T$. The experimentally observed evolution of the FWM signal in FLG as a function of PP-SP delay in Fig.3 has a trend similar to that shown in Fig.4e,f as function of $\eta_{NRVB}/\eta_{CARS}$, validating the resonance-dominated scenario.

A more quantitative picture can be derived from Eqs.3,4, where the PP and SP temporal profiles are taken into account, matching those retrieved from the experimentally measured autocorrelation (Fig.2).

As model parameters for the FLG we use the experimental SR value $\tilde{\nu}_{ca}$=1580cm*-1*, with fitted $\tau_{ca}=1.1\pm 0.1$psFerrariPRL2006 ; Bonini2007 (corresponding to FWHM(G)$=10$cm*-1*) and $\tau_{da}=\tau_{ba}=\tau_{ea}=10\pm 2$fs in agreement with the value measured for SLGBrida2013 . The ratio between NVRB and CARS contributions $\frac{\eta_{CARS}}{\eta_{NVRB}}=(3.0\pm 0.7)\times 10^{-5}$ is obtained by fitting to the experimental data in Fig.3 with Eqs.1,19,20. The resulting spectra (colored dashed lines in Fig.3), evaluated by tuning only the PP-SP delays, are in good agreement with the experimental data, with $\frac{\eta_{CARS}}{\eta_{NVRB}}$ as the only adjustable parameter. This ratio, combined with Eqs.5,6, allows us to extract the ratio between the third-order nonlinear susceptibilities for CARS and NVRB: $\frac{|\chi^{(3)}_{CARS}|}{|\chi^{(3)}_{NVRB}|}\sim 1.3$ at the G-phonon resonance.

The dependence of our spectral profiles on the inter-pulse delay, $\Delta T$, indicates that the peculiar FWM lineshapes for SLG and FLG originates from the interference between two electronically resonant radiation-matter interactions (NVRB and CARS) rather than from a matter-only Hamiltonian coupling the electronic continuum and a discrete phonon state (implying a resonance between the corresponding energies), resulting in the Fano resonancefano1961 suggested in Ref.LafetaArx2017 .

In the electronically non-resonant case, CARS provides access to the real part of $\chi^{(3)}$evans_coherent_2008 . However, due to the dispersive nature of the $\chi^{(3)}$ real partevans_coherent_2008 , it distorts the phonon lineshapesBoydBook , unlike SR. In SLG and FLG the FWM signal arises from the imaginary (non dispersive) CARS susceptibility, and is amplified by its NVRB (third term in Eq.8). Thus, the signal can be used for vibrational imaging, unlike the non-resonant caseevans_coherent_2008 .

The vibrationally resonant contribution $I$ can be isolated by subtracting from the $I_{2}$ FWM signal at $\tilde{\nu}_{2}-\tilde{\nu}_{P}\sim\tilde{\nu}_{G}$, the NRVB obtained by linear interpolation of the spectral intensities measured at the two frequencies at the opposite sides of vibrational resonance:

[TABLE]

where the indexes i=1,2,3 refer to data at $\tilde{\nu}_{1}=\tilde{\nu}_{P}+1545$cm*-1*, $\tilde{\nu}_{2}=\tilde{\nu}_{P}+1576$cm*-1*, $\tilde{\nu}_{3}=\tilde{\nu}_{P}+1607$cm*-1* (i.e. with $\tilde{\nu}_{2}$ near to the G phonon frequency and $|\tilde{\nu}_{1,3}-\tilde{\nu}_{G}|$ greater than two half-widths at half maximum of the measured profiles, as shown in Fig.3).

This combination of electronically resonant NVRB and CARS non-linear responses gives CARS images (i.e. retaining vibrational sensitivity) with signal intensities comparable to those of NVRB, for which sub-ms pixel dwell times have been demonstrated with the use of a point detector, e.g. photomultiplierHendry2010 . In our case, the images in Fig.5 are obtained with a pixel dwell time$\sim 200$ms using a Si-Charge-Coupled Device (CCD) array, aiming at a complete spectral characterization, and scanning the sample at fixed $\Delta T$ with stepper-motor translation stages.

Figs.5(a–c) display nonlinear optical images measured at two different $\omega_{S}$, corresponding to vibrationally non-resonant and resonant conditions. Extracting for each image pixel $I_{1}$ (Fig.5a), $I_{2}$ (Fig.5c), and $I_{3}$, required for Eq.9, we obtain an image with suppression of the signal not generated by FLG, as in Fig.5e.

To obtain a quantitative comparison of the different images, we plot the pixel intensity histogram in Figs.5b,d,f. This gives a bimodal distribution: one peak corresponds to the most intense pixels, associated with FLG ($I_{g}$) and the other is related to the weaker substrate signal ($I_{s}$). The ability to discriminate sample from substrate can be quantified in terms of 1) $I_{g}$ compared to $I_{s}$, evaluated as the difference $I_{g}-I_{s}$, and 2) the proximity of $I_{s}$ to $I=0$ in the histogram (dashed black line in Figs.5b,d,f). These two features can be quantified by the contrast $C$ in order to compare the images: $C=\frac{\overline{I}_{g}-\overline{I}_{s}}{\overline{I}_{s}}$, where $\overline{I}_{g}$ and $\overline{I}_{s}$ are the mean FLG and substrate intensities, corresponding to the local maxima in the histograms in Figs.5b,d,f. In Figs.5g,h we plot the intensity profiles along two scanning paths, one inside (dashed) and the other adjacent to (full line) the FLG flake.

Comparing the three histograms (Figs.5b,d,f), the vibrationally off resonant FWM image (NVRB only, Fig.5b) has the highest $\overline{I}_{g}$. The visibility of the flakes is limited by the noise of the detector and by $\chi^{(3)}$ nonlinearity of the substrate. NVRB, lacking vibrational specificity, can also originate from the glass substrate outside the FLG flake ($\overline{I}_{s}\gg 0$), as indicated by the scanning profile in Fig.5h (red line). This may become a critical limitation in those substrates with $\chi^{(3)}$ much larger than Si ($\lvert\chi^{(3)}\rvert$ $\sim 2.5\times 10^{-10}$ e.s.u.Potma2013 ), such as Au ($\chi^{(3)}=4\times 10^{-9}$e.s.uBOYD201474 ; XENOGIANNOPOULOU2007217 ). Similarly, the vibrationally resonant FWM, $I_{2}$, originating from concurrent CARS and NVRB processes (Fig.5d), has a $\overline{I}_{s}\gg 0$ related to NRVB. The depth of the FWM dip (Fig.5f) is related to the CARS signal intensity, and its vibrational sensitivity brings about a substantial contrast increase, as demonstrated by the close-to-zero average value of the (green) scanning profile in Fig.5h.

In summary, by using an experimental time-delayed FWM scheme, CARS peaks equivalent to those seen in spontaneous Raman were obtained from graphene. By explaining the physical mechanism responsible for the FWM signal, we demonstrated that the spectral response can be described in terms of joint CARS and NVRB contributions concurring to the overall signal. Unlike non-resonant FWM, where dispersive lineshapes hamper vibrational imaging of biological systems, the resonant nature of FWM in graphene, which can be traced back to its peculiar electronic properties, mixes CARS and NVRB, resulting in Lorentzian profiles which are either peaks or dips depending on their relative strength. We also demonstrated that CARS can be used for vibrational imaging with contrast equivalent to spontaneous Raman microscopy and signal levels as large as those of the third order nonlinear response.

I acknowledgements

We thank M. Polini for useful discussions. We acknowledge funding from the EU Graphene Flagship, ERC grant Hetero2D, EPSRC grants EP/L016087/1, EP/K01711X/1, EP/K017144/1. This project has received funding from the European Union’s Horizon 2020 research and innovation programme under grant agreement No. 785219 - GrapheneCore2.

II Methods

The third -order response for the SLG and FLG samples can be obtained from the third-order polarizationMukamel1999 :

[TABLE]

where $N$ is the number of scatterers, $S^{(3)}(\tau_{1},\tau_{2},\tau_{3})$ may be expressed asMukamel1999 :

[TABLE]

and $\mathcal{E}(t)$ is the total electric field on the sample

[TABLE]

Consider a SP at $\Delta t_{S}=0$ with $\Delta t=\Delta t_{P}$. The energy level diagrams in Fig.1schematically illustrate the CARS and NVRB processesMukamel1999 :

[TABLE]

[TABLE]

where $\bar{\omega}_{ij}=\omega_{i}-\omega_{j}-i\gamma_{ij}$.

By Fourier transform, the frequency dispersed signal can be expressed as:

[TABLE]

In order to reduce the computational effort to calculate Eqs.13,14, we also write the pulse fields in terms of their Fourier transforms, obtaining:

[TABLE]

where $\eta_{\hbox{\footnotesize{CARS}}}=n_{\hbox{\footnotesize{CARS}}}\mu_{ba}\mu_{cb}\mu_{cd}\mu_{ad}$. In this way all the temporal integrals can be solved analytically:

[TABLE]

using the energy conservation, represented by the delta distribution:

[TABLE]

the $\omega_{2}$ integral can be simplified:

[TABLE]

In a similar way, using $\eta_{\hbox{\footnotesize{NVRB}}}=n_{\hbox{\footnotesize{NVRB}}}|\mu_{ea}|^{4}$, Eq.14 can be written as:

[TABLE]