Main

Topographic corrugations and charge puddles in graphene are two of the most significant types of disorder in this new material. Topographic corrugations2,3,4, for example, have been suggested as a cause for the suppression of anticipated antilocalization5. Electron and hole puddles6 have similarly been blamed for obscuring universal conductivity in graphene7. These issues are part of a puzzle regarding the factors that limit graphene’s mobility8,9,10,11,12. In order for graphene to fulfil its promise as a next-generation nanodevice substrate it is important to understand the origin of the disorder and the influence it has on Dirac fermions. We have made new progress in this direction by using the techniques of scanning tunnelling microscopy (STM) and spectroscopy to simultaneously probe topographic and electronic disorder in graphene with an electron-density spatial resolution two orders of magnitude higher than previous scanning single-electron transistor microscopy measurements6.

Figure 1a shows the STM topography of a typical 30×30 nm2 area of a graphene monolayer on SiO2. We observe random corrugations with lateral dimension of a few nanometres and a vertical dimension of ∼1.5 Å (r.m.s.), probably due to roughness in the underlying SiO2 surface and/or intrinsic ripples of the graphene sheet2,3,4,13. STM imaging at the atomic scale clearly resolves the graphene honeycomb lattice on top of the broader surface corrugation all over the sample surface (inset).

Figure 1: STM topography and charge-puddle profile of graphene.
figure 1

a, STM topograph (sample–tip bias Vb=−0.25 V, I=20 pA) of a 30×30 nm2 patch of graphene resting on a SiO2 substrate. Inset: Close-up topograph of the graphene honeycomb lattice. b, dI/dV spectra taken at two points separated by 17 nm on a graphene surface. Positional change in Dirac-point energy can be seen. Inset: Sketch showing how changes in the Dirac-point energy (ΔED) are proportional to changes in dI/dV signal intensity (ΔdI/dV) at a fixed sample–tip bias. c, Dirac-point energy (ED) map of a single charge puddle lying in the same area as shown in a (gate voltage Vg=15 V). This is converted to a local charge-density map of graphene (an average background charge density of 0.9×1012 cm−2 has been subtracted). d, Fixed-bias dI/dV map (Vb=−0.25 V, Vg=15 V) over the same area as a and c showing the same puddle profile.

Full size image

We explored the inhomogeneous graphene charge density by spatially mapping the Dirac point (that is, the charge neutral point in the density of states of undoped graphene). The graphene local density of states at the Dirac point shows a local minimum, which is reflected by a dip in the tunnelling spectra of graphene (Fig. 1b). The energy position of the dip, e VD, marks the Dirac-point energy, ED, offset by a constant ∼63 meV shift that arises from the loss of energy experienced by electrons as they inelastically tunnel into graphene by generating a phonon14. Spatial variation in the measured value of ED reflects the spatial profile of charge inhomogeneity in graphene (two spectra taken at points separated by 17 nm, for example, are shown in Fig. 1b). Charge puddles can thus be mapped by measuring the tunnel spectrum at every pixel over a given area and identifying ED at each point. A Dirac-point map, ED(x,y), can be converted into a charge-density map, n(x,y), through the relation n(x,y)=ED2(x,y)/π(ℏvF)2, where vF is the Fermi velocity in graphene. Figure 1c shows such a map of ED for the same area as shown in Fig. 1a at an applied gate voltage of Vg=15 V. We clearly resolve 30 meV fluctuations in the Dirac-point energy, corresponding to charge-density fluctuations of ∼4×1011 e cm−2, where e is the charge of an electron. A single puddle of electrons with a width of ∼20 nm can be seen over this area. Integration of n(x,y) over the puddle area yields a total charge inside this puddle of 0.3±0.2e (the average background charge density has been subtracted).

Charge puddles can also be probed by spatially mapping the tunnelling differential conductance, dI/dV, for a fixed sample–tip bias held slightly below VD. This technique reduces data acquisition time by an order of magnitude and is particularly suited for measuring large graphene areas containing multiple charge puddles. The basis for using this second technique to measure charge puddles is illustrated in Fig. 1b, inset. In the vicinity of the Dirac point (that is, VD) variations in dI/dV are proportional to variations in the electronic local density of states of graphene14. dI/dV maps taken at a fixed bias close to VD can thus produce a map of ED, up to a multiplicative factor. This is demonstrated by the fixed-bias (Vb=−0.25 V, Vg=15 V) dI/dV map in Fig. 1d, which shows the same charge puddle obtained from direct ED mapping (Fig. 1c). Applying this method to a larger area (topography shown in Fig. 2a), we are able to map the profile of multiple charge puddles as seen in Fig. 2b. Individual puddles with an average lateral dimension of 〈L〉∼20 nm are clearly resolved (the electron-rich puddle outlined by a dashed black box is the same as that shown in Fig. 1c). Such puddles are prevalent in graphene, and we have used this technique to explore 23 electron-rich charge puddles over an area of 23,000 nm2 for three different graphene samples.

Figure 2: Large-area image of graphene topography and charge puddles.
figure 2

a, 60×60 nm2 constant-current STM topograph of graphene (Vb=−0.225 V, I=20 pA). b, dI/dV map (Vb=−0.225 V, I=20 pA, Vg=15 V) taken simultaneously with a reveals electron puddles with a characteristic length of ∼20 nm. c, Curvature of surface obtained by calculating the Laplacian of the topographic image shown in a. Upper left dashed boxes indicate the same area as shown in Fig. 1.

Full size image

The same perturbations that create graphene charge puddles also act as scattering sites for the Dirac fermions in graphene, leading to quasiparticle interference (QPI) patterns15,16. This can be seen in Fig. 3a, which shows a dI/dV map taken with Vb=0.35 V (Vg=15 V) over the same area as shown in Fig. 2a. Standing wave patterns in electronic local density of states with a smaller feature size than the charge puddles are clearly resolved on top of the smooth background provided by the puddle profile shown in Fig. 2b. Dispersion in the QPI can be seen in Fig. 3b,c, which shows the interference wavelength decrease as sample–tip bias is increased to 0.6 and 0.85 V respectively (Vg=15 V). We emphasize that the charge puddles are a separate phenomenon from the QPI and that the puddle size scale, 〈L〉∼20 nm, is unrelated to the energy-dependent QPI wavelength.

Figure 3: Quasiparticle scattering on a graphene surface.
figure 3

a–c, dI/dV maps of the same area as shown in Fig. 2 obtained at Vb=0.35, 0.6 and 0.85 V respectively. The tunnel current is held at I=50, 60 and 70 pA respectively and the gate voltage is fixed at Vg=15 V for all three measurements. Lower right insets: two-dimensional (2D) Fourier transform of each image. Upper left dashed boxes indicate the same area as shown in Fig. 1. Red arrows in c point to localized scattering centres. d, Radial averaged intensity profiles of the 2D Fourier transforms shown in a–c plotted as a function of k. Red lines indicate Lorentzian fits. Curves are vertically displaced for clarity. e, Quasiparticle energy dispersion above and below the Dirac point (VD=−0.2 V). Each point is extracted from a Fourier analysis as in a–d. Solid red lines show fitted linear curves yielding vF=1.5±0.2 and 1.4±0.2×106 m s−1 for upper and lower branches. Blue dashed lines indicate theoretical dispersion for vF=1.1×106 m s−1 assuming 63 meV offsets due to phonon-assisted inelastic tunnelling14. f, Gate dependence of k as a function of Vg at a constant sample–tip bias of Vb=−0.75 V (each point is extracted from a Fourier-analysed dI/dV map). The solid red line shows the calculated dispersion obtained using equation (1) and vF as measured in e. The dashed blue line shows the theoretical dispersion arising when vF=1.1×106 m s−1. Error bars in e and f indicate the standard deviation of the mean obtained from fittings of the peak position from Fourier transforms such as those shown in d. g, Schematic diagram of the 2D Brillouin zone of graphene, with orange circles indicating constant-energy contours for states around the K and K′ points near the Fermi energy. The scattering wavevector for an intravalley backscattering process is given by q.

Full size image

These results raise two fundamental questions: (1) what specifically causes the charge puddles? and (2) how do the graphene Dirac fermions scatter from them, thus causing QPI? We now address these questions by first discussing electron scattering from the charge puddles and then by determining the actual origin of the charge puddles (we find it convenient to answer question (2) before answering question (1)).

The observed QPI patterns can be understood as the result of quasiparticle scattering from a disordered potential. This is schematically illustrated in the reciprocal-space sketch of Fig. 3g, where constant-energy contours cut through conical graphene bands produce circles with energy-dependent radius k around the Dirac points at K and K′. Intravalley scattering processes caused by long-range disorder scatter the electrons across the diameter of a single constant-energy circle through a scattering wavevector q (red arrow in Fig. 3g). This results in |q|=2k, that is, electrons are backscattered. 2D Fourier transforms of the dI/dV maps in Fig. 3a–c (shown in the insets) convert the observed spatial oscillations to reciprocal space and reveal constant-energy rings of radius 2k.

Probing QPI as a function of Vb enables us to map the 2D band structure of graphene17. Figure 3d shows a radial average of the Fourier transforms in Fig. 3a–c, and it is clear that the dominant wavevector, |q|=2k, of the observed QPI (the radius of the ring) varies significantly as a function of Vb. Figure 3e plots electron tunnel energy E=e Vb versus k=q/2 (red dots) from such analysis and reveals a linear dispersion relation for states above (e Vb>VD) and below (e Vb<VD) the Dirac point (Vg=15 V leads to a constant average VD=−0.2 V for this measurement). Fitting this data with the expected graphene dispersion relation, E=ℏvFk, we obtain vF=1.5±0.2×106 and 1.4±0.2×106 m s−1 for states above and below the Dirac point, respectively.

As the first and so far the only gate-tunable 2D system that is accessible to STM study, graphene provides a unique opportunity to probe the energy dependence of the QPI without changing the STM sample–tip bias. QPI patterns obtained in this way for fixed Vb=0.75 V and a changing Vg were Fourier analysed as above, resulting in a k versus Vg dispersion that is plotted in Fig. 3f. From the linear band structure of graphene (including the inelastic phonon offset, ℏω0=63 meV) we expect this gate-dependent dispersion to have the following form:

where n is the net charge carrier density induced by both the gate (Vg) and the environment (V0) assuming a simple parallel capacitor model. Here α=7.1×1010 cm−2 V−1 is estimated from the device geometry and V0≈0 V can be obtained from gate-dependent spectroscopic measurement14. Using the value for vF obtained from the data in Fig. 3e, we find that equation (1) fits our measured gate-dependent dispersion fairly well with no adjustable parameters (Fig. 3f, solid red line).

Our QPI-based electronic dispersion measurement differs from theoretical expectations and previous experimental measurements of epitaxial monolayer18 and bilayer17 graphene grown on SiC. There are three main points to notice. First, a gap exists in our experimental dispersion relation (Fig. 3e) at k=0, which we attribute to energy loss to the ℏω0=63 meV phonon modes during inelastic electron tunnelling14. Second, our extracted band slopes are ∼32±18% bigger than what we expect from the commonly accepted value of the graphene Fermi velocity, vF=1.1×106 m s−1 (this ‘theoretical’ value of vF leads to the poorly fitting dashed blue lines in Fig. 3e and f, with the inelastic phonon energy loss taken into account). We note that angle-resolved photoemission spectroscopy measurements of graphene also result in a similarly increased slope in the filled states, which has been attributed to band renormalization due to plasmons19.

The third intriguing aspect of our observed QPI is the fact that we see backscattering at all. Theoretical models that take Dirac fermion pseudospin into account suggest that intravalley backscattering processes are forbidden in monolayer graphene20 (in sharp contrast to bilayer graphene, where intravalley backscattering processes are allowed17). Intravalley backscattering was recently reported as absent from monolayer graphene epitaxially grown on SiC, and pseudospin-suppressed backscattering was provided as an explanation18. However, some intravalley backscattering is in general expected to occur as a second-order process, which is predicted to lead to a fast decay of standing wave patterns21. Such a second-order process, as well as other symmetry-breaking mechanisms, may explain our observed backscattering in graphene flake samples22.

We are now poised to explain the origin of the charge puddles, which is also the origin of the scattering-induced QPI that we observe. We first rule out the hypothesis that topographic corrugations in graphene are a primary cause of the charge puddles. A comparison between the geometry of the charge puddles we observe (Fig. 2b) and topographic corrugations over the same area (Fig. 2a) yields no apparent correlation, as the puddles are an order of magnitude larger than the size of the topographic corrugations. We have also computed the curvature of the graphene monolayer characterized by the Laplacian of the topography, . Figure 2c shows a map of the curvature over the same surface area as Fig. 2a. The average feature size in the curvature map is more than an order of magnitude smaller than that of the charge puddles, further ruling out surface corrugation as the cause of the puddles.

There is, however, a strong correlation between highly localized features seen in our large bias dI/dV maps and the charge puddles. These localized scattering centres show up as ‘dots’ in the QPI patterns and occur only in electron-rich charge puddles when the electron wavelength is reduced by large bias, as shown by the red arrows in Fig. 3c. We have observed such localized scattering centres in all of the electron-rich puddles that we have tested. For example, in Fig. 4a we show STM topography of a different region of the graphene surface that shows typical charge puddles in dI/dV maps obtained at sample–tip biases very close to the Dirac point (Fig. 4b). When the bias is moved away from the Dirac point, as shown in Fig. 4c, we clearly see local scattering centres in the electron-rich regions of these charge puddles. Because the scattering centres do not coincide with any clear topographical features, we believe that they arise from individual charged impurities located beneath the graphene. This interpretation is supported by recent experiments on suspended graphene sheets23,24, as well as theoretical predictions for impurity-based charge inhomogeneity in graphene25,26.

Figure 4: Impurities in graphene.
figure 4

a, STM topography of 50×50 nm2 area of graphene. b, dI/dV map at bias near Dirac point (Vb=−0.29 V, I=25 pA, Vg=15 V) showing electron puddles due to charge fluctuations over the same region of graphene as a. Red crosses indicate the location of quasiparticle scattering-centre impurities observed in c. c, dI/dV map of the same area at larger bias (Vb=−0.75 V, I=80 pA and Vg=60 V) revealing impurity scattering centres in electron-rich charge-density puddles (red crosses). d, Integrated charge per electron puddle plotted as a function of the number of observed impurities in each puddle (puddles are defined as the electron-rich regions left after subtracting the average background charge density). Error bars represent the total uncertainty arising from both the evaluation of VD (that is, standard deviation of mean fits) and possible extraneous tip-gating effects (that is, range of estimated effect for different possible tip geometrices). A linear fit to the data (black line) gives the charge associated with each impurity as 0.07±0.03e (e is the charge of an electron).

Full size image

To gain deeper insight into the origin of these subsurface impurities, we performed numerical integration of the charge in five different charge puddles and compared the total amount of charge per puddle to the number of impurities per puddle. In Fig. 4d we plot the total charge of the puddles as a function of the number of impurities they contain. This data falls roughly on a line, the slope of which enables us to estimate that the average charge fluctuation associated with an individual impurity is ∼0.07±0.03e. Interestingly, some calculations27,28 predict a charge transfer of this order when molecules from air (such as N2 and H2O) are physisorbed onto graphene (although the measured charge fluctuation may deviate from the actual amount of charge transfer owing to the specific criteria by which we define a charge puddle (see Fig. 4 caption)). This finding, combined with the fact that our samples are prepared in ambient conditions, suggests that molecules from air trapped between graphene and the SiO2 substrate are one likely origin of the charge puddles that we observe in graphene flake nanodevices.

We have imaged the nanometre-scale charge landscape that Dirac fermions experience as they move through graphene. We show directly that charge puddles with an average length scale of 20 nm arise from charge-donating impurities. Electronic scattering from these charge fluctuations leads to unexpected backscattering processes. These findings give us new insight into the microscopic processes that limit electron mobility in graphene flakes, and point toward new strategies for improving graphene nanodevice behaviour.

During preparation of our manuscript, we became aware of additional recent STM work on graphene flakes29.

Methods

Our graphene monolayer flakes were prepared on an oxidized Si wafer in a similar fashion as described in ref. 1. We made electrical contact to graphene by direct deposition of 12-nm-thick Ti (or Au) electrodes through a stencil mask to avoid photoresist contamination. Heavily doped Si under a 285 nm SiO2 layer was used as a back gate, enabling us to vary the carrier density in the graphene sample. As part of a cleaning procedure the samples were annealed at 180 ∘C in ultrahigh vacuum (background pressure<10−10 mbar) for ∼10 h. In situ electronic transport measurements have shown that the graphene samples prepared in this way have a typical mobility of ∼6,000 cm2 V−1 s−1 (see the Supplementary Information).

Experiments were conducted with a modified Omicron LT-STM at low temperature (T=4.8 K) and in an ultrahigh-vacuum environment with base pressure<10−11 mbar. We find that the preparation of STM tips is crucial for reliable spectroscopic measurement on graphene. To ensure that our STM tips were free of anomalies in their electronic structure, we calibrated the tips by performing tunnelling differential conductance (dI/dV) measurements on a clean Au(111) surface both before and after graphene measurement. dI/dV spectra were measured using lock-in detection of the a.c. tunnel current, I, after adding an 8 meV (r.m.s.) modulation at 517 Hz to the sample bias voltage Vb.