Local density of states and scanning tunneling currents in graphene

We present exact analytical calculations of scanning tunneling currents in locally disordered graphene using a multimode description of the microscope tip. Analytical expressions for the local density of states (LDOS) are given for energies beyond the Dirac cone approximation. We show that the LDOS at the $A$ and $B$ sublattices of graphene are out of phase by $\pi$ implying that the averaged LDOS, as one moves away from the impurity, shows no trace of the $2q_F$ (with $q_F$ the Fermi momentum) Friedel modulation. This means that a STM experiment lacking atomic resolution at the sublattice level will not be able of detecting the presence of the Friedel oscillations [this seems to be the case in the experiments reported in Phys. Rev. Lett. {\bf 101}, 206802 (2008)]. The momentum maps of the LDOS for different types of impurities are given. In the case of the vacancy, $2q_F$ features are seen in these maps. In all momentum space maps, $K$ and $K+K^\prime$ features are seen. The $K+K^\prime$ features are different from what is seen around zero momentum. An interpretation for these features is given. The calculations reported here are valid for chemical substitution impurities, such as boron and nitrogen atoms, as well as for vacancies. It is shown that the density of states close to the impurity is very sensitive to type of disorder: diagonal, non-diagonal, or vacancies. In the case of weakly coupled (to the carbon atoms) impurities, the local density of states presents strong resonances at finite energies, which leads to steps in the scanning tunneling currents and to suppression of the Fano factor.


Introduction
Graphene [1,2] consists of a monolayer of covalently bonded carbon atoms forming a twodimensional honeycomb lattice [3,4]. Low-energy electronic excitations in graphene are well described as massless Dirac fermions with an additional pseudospin degree of freedom. Because of the Dirac spectrum, impurities can have a strong effect on the local electronic structure of graphene when the Fermi energy is near the Dirac point [5,6,7,8,9,10,11]. Impurities in graphene can be in the substrate, in the form of adatoms, or as imperfections in the lattice itself. In one hand, there has been very significant progress in decreasing the amount of disorder introduced in graphene, for example by fabrication of suspended samples [12]. On the other hand, impurity effects have been explored to modify and tailor the electronic, thermal and chemical properties of graphene. Examples of the later include experiments with graphane [13], chemical substitution of some of graphene's carbon atoms by boro! n an! d nitrogen atoms [14,15], and doping of graphene with metals on top [16,17].
Since graphene is an atomically thin membrane, it can be easily accessed with Scanning Tunneling Microscopy (STM) measurements. Impurity effects can be studied with atomic resolution and the local spectrum can be obtained by STM spectroscopy. In addition, atomic manipulation can also be performed with STM. There has been several STM studies of graphene grown epitaxially on SiC [19,20,21], mechanically exfoliated graphene on SiO 2 [22,18,23,24,25], and graphene flakes on graphite [26]. In fact, STM experiments have proved instrumental in mapping the topography of corrugated graphene and determining the existence of charge puddles [27]. Additionally, STM experiments are also able to probe the chiral nature of the electrons in graphene when they scatter from impurities [28]. This experimental work showed that intravalley backscattering is virtually absent in graphene. In particular, the STM experiment showed the lack of the 2q F (q F is the Fermi momentum) Friedel modulation on the local density of states of graphene. As we show explicitly below, this lack of modulation can be traced to the fact that the local density of states at the A and B sublattices are out of phase by π, an aspect already noted in passing previously [29,30]. Therefore, the local density of states, when averaged over the unit cell, shows no trace of the 2q F oscillation. Additionally, as we show below, at distances d close to the impurity, d ≪ 1/q F , there is a strong departure from the 1/r 2 spatial dependence [31,8] of the local density of states. Moreover the form of the density of states close to the impurity is very sensitive to type of disorder: diagonal, non-diagonal, or vacancies.
In this work, we present calculations of STM currents in locally disordered graphene. We focus on the case of chemical substitution (by boron or nitrogen atoms, for example) and use a multimode description for the STM tip. We obtain exact analytical expressions for the local density of states, and also present results for energies beyond the Dirac cone approximation. We model the substitutional impurity by both an on-site impurity potential and local hopping disorder. We find that inclusion of the hopping disorder term leads to additional higher harmonics oscillations in real space for the density of states for the sub-lattice that does not contain the impurity. The main oscillations in the two sub-lattices are out-of-phase away from the impurity. For the regime in which the electronic hopping between the impurity and the nearest neighbor carbon atoms is decreased in relation to the hopping between carbon atoms in the clean system, the local density of states presents strong resonances. A vacancy is a extreme case of this regime. These resonances lead to the appearance of steps in the STM current, a signature that should be observable experimentally. These resonances also leads to open channels for tunneling between the STM tip and graphene, and lead to a decrease in the Fano factor.
Another issue addressed in this paper relates to the effect of the tip on the measured STM currents. In the usual analysis, the electrons in the tip are represented by jellium model with constant density of states. In this type of model neither the real part of the self-energy due to the tip-system coupling nor the variation of the density of states with energy is included (wide band limit). The tip, however, is not an infinite metal. In fact it has a structure where the number of atoms in the atomic planes reduces as we approach the tip. In a previous publication [32] we have modeled the tip as a one dimensional model (in that work we have also considered the simplification of zero on-site energy at the impurity), which corresponds essentially to the case of a constant density of state too a good approximation. In that case we found the STM current to be symmetric around zero energy. When we generalize to the case of a multimode tip this symmetry is lost, as we show in this work. Comparing the results of Ref. [32] with those given here it is possible to disentangle the effects due to graphene and to the impurities from those due to the tip. This work shows that some care has to be taken when interpreting the STM currents directly. The present manuscript is organized in the following way: In Section 2, the Green's function formalism is presented, with analytical results for the Green's function for graphene with a substitutional impurity, and for the STM tip modeled by a multimode system. Local density of states results are presented in Section 3, and results for the STM current are presented in Section 4. Section 5 contains final discussions and conclusions.

Graphene
The honeycomb lattice has two carbon atoms per unit cell, one from sublattice A and one from sublattice B, as depicted in Fig. 1. The unit cell vectors are a 1 and a 2 , with magnitudes |a 1 | = |a 2 | = a, where a = √ 3a 0 ≃ 2.461Å, and a 0 is the carbon-carbon distance. Any lattice vector r can be represented in this basis as r = na 1 + ma 2 , with n, m integers. In Cartesian coordinates, a 1 = a 0 (3, √ 3, 0)/2 and a 2 = a 0 (3, − √ 3, 0)/2, and the reciprocal lattice vectors are given by . The vectors connecting any A atom to its nearest neighbors are: We consider here the case where a substituting atom replaces a carbon atom in the A sublattice, say. When this happens two effects take place: (i) the on-site energy ǫ i at the impurity site is different from that at the carbon atoms; (ii) the hopping from and to the impurity atom, t i , changes relatively to that of pristine graphene. In the latter case, we model the change in the hopping by introducing an additional non-diagonal term to the Hamiltonian, such that t i = −t + t 0 (see below). Using these definitions the Hamiltonian can be written as: is the kinetic energy operator and a † , a (b † , b) are fermion creation and annihilation operators in the A (B) sites. The spin index is omitted for simplicity. We consider an isolated impurity located at r = (0, 0, 0), on sublattice A, so that its contribution to the Hamiltonian has two terms: and for hopping and potential disorder, respectively. In the case of zero chemical potential, when the Fermi level crosses the Dirac point, the system is most susceptible to the presence of impurities. In the particular case t 0 = t, hopping to the impurity site is completely suppressed, and the scattering term V t represents a vacancy. It is well known [34] that the formation of a vacancy will lead to some local distortion of the carbon-carbon bonds. This effect is not incorporated in our Hamiltonian, which in that case would not be exactly solvable. The substitution of carbon atoms by boron or nitrogen has the main consequence of changing the local hopping and the onsite energy, both effects included in our description. The particular choice for boron or nitrogen is due to size restrictions imposed by the unit cell of graphene. We note here that replacement of carbon atoms by boron and nitrogen was already experimentally achieved [15]. We first calculate the single particle Green's functions for the system comprised of graphene and a single impurity, described by the Hamiltonian H above. The single particle Green's functions carry sub-lattice indices, and are defined as: The equations of motion for the Green's functions are given by: where and N c is the total number of unit cells in the lattice, and ω n are fermionic Matsubara frequencies. Note that λ pq = λ qp . This is a consequence of the impurity hopping term V t , which breaks sub-lattice symmetry. The impurity potential term V i also breaks sub-lattice symmetry, and therefore ǫ i appears in an asymmetric way in the equations above. The set of equations of motions can be solved exactly. The presence of the scattering term V t leads to the appearance of the phases φ k and a more complex form for the T -matrix than usual. The exact solution for the Green's functions can be written, after a lengthy calculation, in the form [10]: where all the terms (G, G 0 , g, h and T ) also depend on ω n (omitted here for brevity). The terms g, h, and T correspond to sums over infinite series of Feynman diagrams for impurity scattering, and are given by: and where the denominator D(ω n ) is defined as with the diagonal component of the Green's function for the clean system given by ( = 1) which is translationally invariant. The expressions for the Green's functions, Eqs. (14)- (19), are exact analytic solutions for graphene with one isolated substitutional impurity, including contributions from both the on-site energy ǫ i and the hopping parameter t 0 . Inclusion of the off-diagonal disorder t 0 leads to additional terms, and additional ω-dependence of the graphene Green's function. Since single particle properties, such as local electronic spectra, can be obtained directly from the Green's functions, this ω-dependence has direct experimental consequences, such as for STM spectroscopy measurements. These results for the graphene Green's functions have been obtained in Ref. [10], where local density of states maps, electronic spectra and Friedel oscillations have been calculated for the cases of boron and nitrogen substitution. We include here a brief derivation of the analytical expressions for the Green's functions, Eqs. (14)- (21), for completeness. Upon closer inspection, the physical meaning of the extra terms in G aa become clear. The significance of the term g(ω n ) in (14) which only appears in G aa , is more easily interpreted if we do a double Fourier transform to real space. This term corresponds to the return amplitude to the impurity site for an electron starting at the impurity site. The factor 1/D(ω n ) contains a sum over an infinite series of intermediate scattering events, but the overall process is bounded and the t 2 0 factor denotes hopping from the impurity to the nearest neighbor B-sites and back to the impurity site. Likewise, an interpretation can be given to the other term which only appears in G aa , namely, h(ω n )G 0 k (ω n ). A double Fourier transform shows that this term contributes to G aa (r, 0) and describes the amplitude of propagation between the impurity site and another A site, again with an infinite series of intermediate scatterings. Similarly, the h(ω n )G 0 p (ω n ) term contributes to G aa (0, r). No such terms can, of course, appear in G bb when the inpurity is at a A site. And no such term can be present when there is only the impurity potential term term, which appears in both G aa and G bb , is the usual term also present in simple on-site impurity potential problems, but in this case the T -matrix contains contributions from both t 0 and ǫ i .

STM Tip
Let us consider a model for the STM tip represented by a multimode system. The bulk of the tip is modeled by a square lattice with two atoms in the transverse direction. The end of the tip is represented by a single atom. This choice renders the system multimode, with two transverse modes. It is as simple to include a truly three dimensional tip, but the current will not be much affected by it. The schematic atomic structure of the tip is represented in Fig. 2.
The Hamiltonian for the tip can be written as H t = H b + H 0 , where H b represents the bulk of the tip and H 0 the tip's last atom. These two parts of the Hamiltonian are defined as and where c † (c) are creation (annihilation) operators for fermions in the tip. Let us now consider the case of the bulk part of the Hamiltonian's tip, H b . In the case of a square lattice the wave function is separable and can be written as |ψ l,t = |φ l |φ t , with the longitudinal part of the wave function |φ l given by and the transverse part |φ t given by In Eqs. (25) and (26), the states |n, m = |n |m are position states and the numbers θ l and α t are given by and The resolvent for the Hamiltonian H b is defined asĜ where the + superscript denotes the retarded function. In the eigenstate basis, it has the form where E l,t are the eigenvalues of H b , with H b |ψ l,t = E l,t |ψ l,t , given by For the calcution of the STM current we will need the surface Green's functions defined as The calculation of (31) and (32) requires the evaluation of the integral which is easily done by contour integration methods [33]. The final results are In the case β 2 j < 1, the Green's functions are obtained from (34) and (35) by removing the factor sgn(β j ), and choosing the positive sign for the square root of the negative argument: Using Eq. (34), the local density of states at the sites n = −1, m = 1, 2 given as usual by Fig. 3. The multimode nature of the tip is clearly seen in the form of the density of states.

Graphene's local density of states
The properties of the STM current depend on the local density of states below the tip of the microscope. In this section we compute the Green's function of graphene in real space, from which the local density of states can be obtained. In Sec. 2.1, the position vector r denotes the position of the unit cell, which contains two atoms. The local density of states (per spin) at the sub-lattice A and sub-lattice B atoms of the unit cell localized in the position r is defined as where x = a, b, and G xx (r, r, ω) is obtained from after the usual analytical continuation ω n → ω + i0 + of the Matsubara Green's function. In the unit cell r = 0, which contains the impurity in its A site, we have simple expressions for the Green's functions for the A and B sites. These read whereG 0 (ω n ) is given bỹ andḠ 0 (ω n ) is defined in Eq. (20). As is clear from the above equations, the central quantity that needs to be calculated is the integrated Green's functionḠ 0 (ω n ). Since we want to compute the density of states for energies beyond the Dirac cone approximation we need to include in the density of states powers of the energy beyond the usual linear term. In a previous work [35], we have derived an expansion for the density of states (per unit cell, per spin) valid for energies up to ∼ 2.5 eV, reading (E = ω) Using this expression for the density of states, a close form for the retarded function G 0 (ω n → ω + i0 + ) can be derived. The imaginary part ofḠ 0 (ω) reads and the real part has the form where P 1 (x) and P 2 (x) are polynomial functions given by The energy D c is a cut-off energy chosen as D 2 c = √ 3πt 2 . It is possible to derive simple analytical expressions for the local density of states at the unit cell r = 0, where the impurity is located, using the Dirac cone approximation. If the full forms, Eqs. (44) and (45), of the Green's are used, an analytic close form is still possible, but is somewhat cumbersome.
For the calculation of the local density of states, the case of a vacancy and the case where ǫ i = 0 and t 0 = t have to be treated separately. For the vacancy, the density of states at the neigbouring B atom is with For the case of a general substituting atom, the local density of states at the impurity atom ρ a (0, ω) is obtained from the imaginary part of the Green's function which reads The imaginary part of the Green's function at the B site, next nearest neighbor to the impurity, reads  It is well-known that the density of states due to a vacancy has a strong departure from the pristine value close to the Dirac point. This is due to the logarithmic singularity seen in Eq. (48). In the case where there is an enhancement of the hopping amplitude between the impurity and the neighboring atoms (negative t 0 ) the local density of states retains its linear behavior close to the Dirac point, but its value at the A and B sublattices are different, as expected (Fig. 4, upper right panel). This case would mimic a boron impurity atom. Boron has a larger atomic radius (R ≃ 0.85Å) than carbon (R ≃ 0.7Å), and there should be an increase in the absolute value of the hopping amplitude when it substitutes a carbon in graphene. In the case of a decreasing of the electron hopping between the impurity and the carbon atoms (positive t 0 ) the behavior of the density of states is more interesting since resonances start to develop around the Dirac point (Fig. 4, two lower panels). Note that the density of states still goes to zero at the Dirac point. This behaviour is reminiscent of the vacancy, since we can picture the two resonances developed in both sides of the Dirac point as a splitting of the divergent peak for the vacancy due to the departure of t 0 from its vacancy value t 0 = t. This case would mimic a nitrogen atom, which has a smaller atomic radius (R ≃ 0.65Å) than carbon.
The calculation of the local density of states for finite r requires the calculation of the Fourier transform entering the definition in Eq. (39). For G bb (r, r, ω) we use an approximation for φ(k), which reads φ(k) ≃ 3a 0 (k y − ik x )/2. Carrying out the Fourier transform we obtain where F n (ω, r) (n = 0, 1) is defined as F n (ω, r) = (2n + 1)A c a n with J n (x) the Bessel function of integer order n; and k c = 2 √ π/( 3 √ 3a 0 ), α = |ω|r/v F , A c = 3 √ 3a 2 0 /2, and v F = 3ta 0 /2. The Cauchy principal value of the integral in Eq. (54) is computed using numerical methods.
In Fig. 5, we plot the local density of states at both sub-lattices A and B. The typical oscillations due to the presence of the impurity are present. Note that close to the impurity the B sublattice density of states presents higher harmonics as function of r; these are due to the cut-off momentum k c . On the other hand, the large wavelength oscillations are the 2q F Friedel oscillations: from Fig. 5 the wavelength is about λ ≃ 28a 0 ; on the other hand, for the energy ω = 0.5 eV the Fermi momentum is q F = 1/(9a 0 ), implying λ ≃ π/q F ≃ 28a 0 . At large values of r the two density of states are out of phase by a factor of π. Therefore, when we average over the unit cell the result is essentially the pristine density of states. This result has strong consequences for STM experiments. If the STM experiment lacks atomic resolution at the A and B sublattices level, the experimental data will show a very faint trace of the 2q F Friedel oscillations. This seems to be the case in the experiments reported in Ref. [28]. Closer to the impurity there is a strong departure from the asymptotic behavior [9] ρ(r) ∝ 1 and the A and B LDOS behave quite differently. STM measurements of a material surface are ideal for studying real-space local features with atomic resolution. In particular, real space modulations of the STM intensity can be observed when impurities are present at the surface of a given material. In general, the impurities lead to elastic scattering between the momentum q F and −q F , which is the most efficient process due to phase space restrictions, leading to 2q F Friedel oscillations. In the case of graphene, the chiral nature of its electronic spectrum changes this general behavior. The Fermi surface has disconnected pieces at different points of the Brillouin zone -the K and K ′ points. The Fermi surface consists of circumferences of radius q F around each K and K ′ points. The scattering process is then characterized by two channels: an intra-cone scattering (within the same K or K ′ points) of momentum change 2q F , and an inter-cone scattering (between the K and K ′ ) points.
A Fourier transform of the real-space STM-intensity currents, proportional to the LDOS, will produce bright spots at the momentum values seen in the real space modulations of the LDOS. When impurities are present, the momentum values characterizing the real space modulation are related to the momentum change associated with a given scattering process. In Ref. [28] it was found that the intra-cone scattering, which would give rise to a bright spot of radius 2q F , was absent in the momentum map of the density of states obtained by a Fourier transform of their STM data.
In what follows, we give Fourier transforms of the local density of states of graphene for the different types of impurities discussed previously in the text. We note that our derivation of the Fourier transform of LDOS uses the full Green's functions for the calculation of the LDOS, and therefore no approximation has been made in the calculation. The full real-space map of density of states can be obtained numerically from Eqs. (38) and (39), using the exact expressions Eqs. (14)-(21), as was done in Ref. [10]. We now perform a Fourier transform of a) b) c) d) the density of states, for each sublattice x = a, b. The results are shown in Fig. 6, where the three rows correspond to ρ a (k, ω), ρ b (k, ω), and the sum ρ a (k, ω) + ρ b (k, ω), from top to bottom; and the four columns correspond to four types of impurities discussed in Fig. 4. We recall that positive t 0 reduces the hopping from the impurity site to its nearest neighbor carbon atoms -the particular case of t 0 = t represents a vacancy -and negative t 0 increases the hopping of the electrons from the carbon atoms to the impurity site. This latter case would correspond to an atom with a a) b) c) d) radius larger than carbon, such as boron, leading to an increase of the hopping relatively to the hopping t between nearest neighbor carbons. The column (a) of Fig. 6 refers to a vacancy. In this case, it is clear that a 2q F circumference is seen around the k = (0, 0) point for the ρ a and ρ b plots, consistent with the modulations shown in Fig. 5. However, the intensity at the 2q F circle around k = (0, 0) is suppressed when looking at the ρ a +ρ b plot. Features at six spots at a distance |K| around k = (0, 0) are also seen, these correspond to the K and K ′ points at the corners of the first Brillouin zone, and represent inter-cone scattering. Additionally, six bright spots at distance |K +K ′ | are also present. These vectors with modulus |K +K ′ | are reciprocal lattice vectors G. For a pristine material the local density of states has the periodicity of the underlying lattice, that is, ρ(r) = ρ(r + R), and therefore, the Fourier transform of ρ(r) must show the same intensity at k = (0, 0) and k = G. While the presence of an impurity breaks the periodicity of the real-space lattice, the reason for k = (0, 0) and k = G being different is the fact that the impurity considered here is not on a whole unit cell, but on only one of the sites of the unit cell. If there were only one site per unit cell and a single short-range impurity, then the periodicity in k-space would be maintained. In the extreme case, if there was a line of impurities of the vacancy type, this would correspond to cutting the system into half, and k = (0, 0) and k = G would still be the same. Note that in all our cases ρ A is still periodic in k-space. The impurity on a site A that we consider here, introduces structure within the R = 0 unit cell, therefore, in k-space there is information going beyond the first Brillouin zone. Mathematically this appears because the B-site is located at R + δ 3 , where R denotes the position of the A sites, and δ 3 is not a lattice vector of triangular Bravais lattice (see Fig.  1). When we take the Fourier transform, the periodicity in k-space does not happen at k = G anymore, but at at larger k. All figures show this signature, which is rather clear in column (c) of Fig. 6, since the Fermi surface energy has been chosen as large as ω = −1.5 eV.
Inter-cone scatterings, represented by the region around the K and K ′ points, are highly angular-dependent [29,30], as can be seen particularly in the ρ b plots. While the scattering around the G vectors do not have such a strong angular dependence, some trigonal warping is observed in ρ B . The scatterings around k = (0, 0) are rotationally symmetric. Fig. 7 corresponds to the case t 0 = 2 eV and ǫ i = −1 eV (case (c) of Fig. 6), showing how the k-space LDOS map evolves with increasing energy.

STM current
In this Section we present calculations of the tunneling current between the STM tip and graphene, when the tip is close to an impurity atom. We model the tip by the multimode tight-binding model, as described in Sec. 2.2. This choice departs from the more simplified approach where the tip is modeled by a one dimensional system [36,37].
There is a number of ways one can use to describe the tunneling of the electrons between the STM tip and graphene. Here we assume that the coupling is made directly either to the impurity atom or to the next neighbor carbon atom. This choice corresponds to probing the local electronic properties at or around the impurity. More general types of coupling are easily included in the formalism. We write this coupling as where the operator d(0) can represent either an electron at the impurity atom in the A sublattice or at the carbon atom in the B sub-lattice.
Since the Hamiltonian of the problem is bilinear we can write it in matrix form as where the matrices V L and V R represent the coupling of the last atom in the tip of the STM microscope to the bulk of the tip and to graphene, respectively, and H b and H g stand for the bulk Hamiltonians of the tip and of graphene, respectively. H g also includes the impurity terms Eqs. (2) and (3). The matrix H is of infinite dimension due to H b and H g . The matrices where 0 represents an infinite dimensional null row vector. The tunneling is a local property, controlled by the coupling of the last atom of the tip to the bulk atoms and to graphene. Since we want to compute local quantities, this is best accomplished using Green's functions in real space. The full Green's function of the system is defined by where 1 is the identity matrix. The matrix form of the Green's function is The quantity of interest is G 00 , which can be shown to have the form where the matrices Σ + L and Σ + R are the self energies and have the form where G + xx is the surface Green's function of the Hamiltonian H g at the impurity unit cell (x = a, b), respectively. Note that the quantity G + xx is computed using Eq. (39) and setting r = 0.
The study of non-equilibrium transport is done using the non-equilibrium Green's function method, or Keldysh formalism. This method is particularly suited to study the regime where the system has a strong departure from equilibrium, such as when the bias potential on the STM tip, V b , is large. In this work, we consider, however, that the system is in the steady state. Since the seminal paper of Caroli et al. on non-equilibrium quantum transport [38], that the method of non-equilibrium Green's functions started to be generalized to the calculation of transport quantities of nanostructures. There are many places where one can find a description of the method [39,40], but a recent and elegant one was introduced in the context of transport through systems having bound states, showing that the problem can be reduced to the solution of an equation similar to a quantum Langevin equation [41]. The general idea of this method is that two perfect leads are coupled to the system, which is usually called the device. In our case the device is defined by the last atom of the tip of the microscope. The Green's function of the device has to be computed in the presence of the bulk of the tip and of graphene. This corresponds to our G + 00 Green's function. Besides the Green's function we need the effective coupling between the last atom of the tip and the bulk atoms as well as that to the graphene atoms, which are determined in terms of the self-energies Therefore the effective coupling Γ L/R depends on the surface Green's function of the tip and of graphene. According to the general theory, the two systems (bulk of the tip and graphene) are in thermal equilibrium at temperatures T L/R and chemical potential µ L/R and are connected to the system at some time t 0 . The total current through the device is then given by where the factor of 2 is due to the spin degrees of freedom, f (x) is the Fermi-Dirac distribution and the transmission T (E) is given by In Fig. 8 we depict T (E) in different cases. Note the asymmetry of the density of states which is exhibited even by the pristine case (Fig. 8, upper left panel). This asymmetry has a two fold nature: (i) it comes from the fact the bulk of the tip has two transverse atoms but the tip has only one; (ii) the fact that the atom at the tip has a different local energy from those in the bulk. This asymmetry carries on to the disordered cases. Additionaly, for the disordered cases the resonances seen in the local density of states has a strong impact on the transition probability T (E), leading to open transport channels with large values of T (E). This is specially true for the vacancy and for the weakly coupled impurity case, that is, when hopping between impurity and carbon atoms is suppressed in relation to the hopping between carbon atoms, as expected for nitrogen substitution.
Since we want to probe the properties of the STM current at zero doping we choose µ L = eV /2 and µ R = −eV /2. Also T L = T R = 0. This renders the calculation of the current   Fig. 4. The transmission function was computed at finite bias; the zero bias case is given in Fig. 8. The four panels here correspond to the same ones in that figure 8.
to a simple one-dimensional integral of T (E) over the energy. The form of the current will reflect the properties of T (E) as function of energy, and, as we have seen, these are markedly different for the different cases, depending strongly on the value and sign of t 0 . The presence of resonances in T (E) leads to steps in the STM current. This is seen in Fig. 9 for the cases of the vacancy and to the case of weak coupling (positive t 0 ) between the impurity and the neighboring carbon atoms.
To fully characterize the STM current, another important quantity is the shot noise [42]. For interacting systems, this quantity contains information on the nature of the quasi-particles, including for example the possible existence of quasi-particles with fractional charge. In disordered systems with no interactions, information on transport through open channels can be obtained. For non-interacting electrons at zero temperature the shot noise is defined as [43] The relevant quantity is not S directly but the Fano factor [42] defined as F = S eJ . When the transmission T (E) is strongly reduced we have F → 1, and the noise is said to be Poissonian. On the other hand, if the system has a finite density of open channels, T (E) → 1, we have F < 1 due to [1 − T (E)] ≪ 1. The resonances that we obtain in T (E) play a role in the resulting form of F . They lead to an enhancement of the current, and to a significant decrease of the Fano factor due to the opening of a transport channel. Note that the opening of the channels is a consequence of the local disorder induced by the impurity.

Conclusions
In this paper we have studied the STM currents through locally disordered graphene. We have considered a tip with transverse modes. Although the tip is strictly quasi-one-dimensional it  Fig. 4. The four panels here correspond to the four panels of the current given in Fig. 9.
still departs from the widely used model of a strictly one-dimensional model. Generalizing now the calculations to a truly three-dimensional tip is reasonably straightforward. The modifications would require introducing a three dimensional square lattice for describing the bulk of the tip, and a decreasing number of atoms for each transverse plane to describe the part of the tip in contact with graphene. This last part would lead to the most significant change in the calculation, since the device would not be a single atom as in our calculations here, but would be represented by a finite number of them, and therefore the Green's function for the device would be a matrix instead of a c-number. Nevertheless, as long as we take the dispersion of the electrons in the tip to have large bandwidth, the current should not depend much on the local density of states of the tip, since this would be essentially constant. This corresponds to the usual wide band limit.
We have also seen that tunneling through either impurity atoms, or their neighboring carbon atoms, depends on the local density of states of graphene. For certain circunstances -vacancy or weakly coupled impurities -there is a development of resonances at or close to the Dirac point. These resonances lead to a strong enhancement of the tunnelling probability which appear as steps in the tunneling current. It is conceivable that graphene could be locally modified in order to take advantage of these strong resonances developed close to the Dirac point. Clearly the substituting atoms would also locally distort the graphene lattice, an effect not included in our description. How much the STM current would depart from the values computed here would depend on the change in the values of the hopping parameter, and on additional features on the density of states due to the disorder. If future research will pursue the route of modifying graphene locally, our results will be important for the characterization of the surface. Even in the present state of affairs, our results could be used to interpret STM current due to local impurities.