Charge State Dynamics and Optically Detected Electron Spin Resonance Contrast of Shallow Nitrogen-Vacancy Centers in Diamond

Nitrogen-vacancy (NV) centers in diamond can be used for nanoscale sensing with atomic resolution and sensitivity; however, it has been observed that their properties degrade as they approach the diamond surface. Here we report that in addition to degraded spin coherence, NV centers within nanometers of the surface can also exhibit decreased fluorescence contrast for optically detected electron spin resonance (OD-ESR). We demonstrate that this decreased OD-ESR contrast arises from charge state dynamics of the NV center, and that it is strongly surface-dependent, indicating that surface engineering will be critical for nanoscale sensing applications based on color centers in diamond.

Nitrogen-vacancy (NV) centers in diamond are actively explored for a number of applications in quantum information processing and sensing because they exhibit long spin coherence times at room temperature, and their spin states can be optically initialized and read out with offresonant excitation [1][2][3][4][5][6].In order to achieve strong interactions with materials and molecules that are external to the diamond, NV centers must be placed close to the diamond surface.It has been well established that the diamond surface can host contaminants, magnetic defects, and electronic defects that give rise to noise, leading to short spin coherence times [7][8][9][10][11][12][13][14][15][16], and recent work has shown that careful preparation of the diamond surface can mitigate this noise, leading to extended spin coherence times [6,7,9].Most NV center applications are based on the negative charge state of NV centers, NV − , and it is also known that NV centers can exist in the neutral charge state, NV 0 .It has been recently shown that both the steady-state charge populations under illumination and the equilibrium charge state in the dark of shallow NV centers can vary significantly from bulk properties [17][18][19][20].Previous work has shown that charge state initialization can influence spin readout [21], however a detailed understanding of how these charge state dynamics can affect schemes for nanoscale sensing and its dependence on sample surface has not yet been established.
Here we demonstrate that the diamond surface can strongly affect both the steady-state charge state distribution and the charge state dynamics of NV centers within nanometers of the diamond surface, and that these charge state properties can significantly degrade optically detected electron spin resonance (OD-ESR) contrast.We focus on two diamond samples (samples A and F) that contain shallow NV centers introduced by ion implantation, which exhibit distinctly different OD-ESR contrast and charge state behavior.We find that sample A exhibits higher contrast, greater charge state stability in the dark, lower NV 0 population under green illumination, and slower ionization and recombination dynamics than sample F.More specifically, we find that sample F exhibits charge state conversion rates that are comparable to the internal spin-dependent dynamics of the NV center, thus leading to decreased OD-ESR contrast.Finally, we measure time-resolved, spin-dependent photoluminescence (PL) and compare to a coupled rate equation model that incorporates charge conversion rates to quantitatively understand the impact of charge dynamics on OD-ESR contrast.
A cartoon schematic of the interactions between shallow NV centers and charge traps is depicted in Fig. 1(a).The surface can host defects that act as charge traps, and tunneling from NV centers to these charge traps gives rise to charge state instability, leading to a decrease in OD-ESR contrast of shallow NV centers [17].
Following Ref. [22], we define OD-ESR contrast as where α i is the average number of photons collected in the readout pulse when the NV − spin state is initialized to the m s = i state.Figure 1(b) shows C ESR plotted for many NV centers as a function of depth across five samples, where depth is measured by detecting the proton NMR signal arising from the microscope immersion oil [23].NV centers within 10 nm of the surface display a wider distribution in C ESR .
In order to investigate the origin of the lower OD-ESR contrast, we identify two samples that have undergone nominally similar surface processing (see Supplemental Material Sec.I [24]), but exhibit markedly different distributions in C ESR , samples A and F [Fig. 1(c)].Most NV centers in sample F exhibit C ESR below 0.3, as well as short coherence times (T 2 < 4 µs), precluding many NV characterization measurements, such as single-shot charge state readout [25] and using proton NMR to measure depth.We use this comparison to quantify the impact of charge state conversion on C ESR .
Both samples were prepared using nitrogen ion implan- tation followed by thermal annealing.Sample A was implanted with a nitrogen ion energy of 3 keV while sample F was implanted with an energy of 1.5 keV.Despite undergoing nominally similar surface processing, sample F contains persistent boron contamination comprising approximately 4% of a surface monolayer as measured by X-ray photoelectron spectroscopy [24], likely arising from contamination in the furnace during vacuum thermal annealing.We note that although sample F is expected to have a shallower distribution of NV centers because of its lower implantation energy, the distribution of C ESR is much lower than another sample that was prepared with the same ion implantation energy that does not exhibit surface contamination, sample B [24].
To understand the origin of the lower OD-ESR contrast in sample F, we examine the steady-state charge distributions in both samples.First we measure the PL spectrum under green illumination.We observe a significantly higher fraction of emission between the zero phonon line (ZPL) of NV 0 (575 nm) and the ZPL of NV − (637 nm) in sample F, as well as a peak at 575 nm [Fig.2(a)], indicating a higher steady-state population of NV 0 [26].This higher NV 0 steady-state population contributes to background fluorescence, decreasing the OD-ESR contrast proportionally.NV centers are initialized into the negative charge state by a green pulse, and the charge state is read out using a chargestate-selective orange pulse with variable delay time.The PL is normalized to the value for the shortest dark time interval.Each NV center in sample F (orange triangles) exhibits a decay in the relative NV − population with a time constant in the range of 11 -300 ms, while NV centers in sample A (blue circles) exhibit much less decay out to 1 s.Solid lines are fits to the data.(c) Stability of NV 0 in the dark in sample A. NV centers are initialized into NV 0 using a green initialization pulse followed by an orange pulse, and then read out with another orange pulse after a variable delay time.The PL is normalized to the value after the green pulse.No change in the population is observed out to 1 s, indicating that recombination kinetics are also slow in sample A.
Recent work has shown that shallow NV centers can exhibit spontaneous conversion to NV 0 without illumination [17].Here, we perform a similar measurement on several NV centers in both samples.First, the NV center charge state is initialized using a green pulse of sufficient length (5 ms at 510 µW) to achieve a steady-state distribution.Then, after a variable time interval in the dark, we measure PL during an orange (590 nm) pulse, which preferentially excites NV − , allowing for charge state readout.Normalizing to the PL at the shortest dark time interval, we observe that NV centers in sample F exhibit a decay in the NV − population over timescales between 11 ms and 300 ms to less than half of the initial value, while a decay of less than 25% is observed in sample A for dark times up to 1 second [Fig.2(b)].In sample A, we also perform the inverse measurement by initializing NV centers into NV 0 using an orange pulse (5 ms at 18 µW, after a green pulse as Fig. 2(b) for consistency [19]) to look for evidence of spontaneous conversion to NV − .We observe no change in the NV − population out to 1 second [Fig.2(c)], indicating that spontaneous charge conversion is very slow in this sample, regardless of the initial charge state population.This is consistent with previous hypotheses that spontaneous charge conversion is mediated by the availability of nearby electron traps at or near the surface [Fig.1(a)] that strongly modify charge state kinetics [17].Now we turn to the charge state conversion dynamics under illumination for both samples.The rate of change of the charge populations can be expressed in the following simplified model [25]: where ρ − is the NV − population, ρ 0 is the NV 0 population, and r ion and r rec are ionization and recombination rates, respectively.Since the total population is conserved (ρ − +ρ 0 = 1), we can express the time dependence of the population in terms of a total charge conversion rate, ρ −,0 (t) ∝ e −rtott , where r tot = r ion + r rec .
We measure charge conversion rates under green illumination by using a fixed orange pulse to initialize the NV centers primarily into NV 0 .Then under green illumination, the population shifts to NV − , and we fit the time-resolved PL to an exponential to extract r tot .Similarly, we measure conversion rates under orange illumination by first initializing with a green pulse and measuring the overall ionization under orange, extracting r tot in the same manner.Both data sets are plotted in Fig. 3 for a set of NV centers in samples A and F. We observe a clear quadratic power dependence for NV centers in both samples, consistent with both ionization and recombination being two-photon processes [25] (see Supplemental Material Sec.VI [24] for a simplified derivation).NV centers in sample F have higher charge conversion rates under both green and orange illumination.Importantly, when the green power is comparable to the typical NV saturation power in our setup (≈ 300 µW), the charge state conversion rate can be faster than 1 × 10 6 s −1 , comparable to the spin polarization rate for NV − , suggesting that fast ionization and recombination processes can occur during the spin initialization and readout times.Since these processes are not spin conserving [27], this fast charge conversion can lead to degradation of the OD-ESR contrast.
To quantify how charge state conversion rates affect OD-ESR contrast for NV centers under distinctly different surface environments, we measure the time-resolved PL during the green readout pulse.Example PL traces with initial NV − spin states m s = 0 (upper traces) and m s = −1 (lower traces) for NV centers in both samples are shown in Fig. 4(a).The area between the two sets of curves is proportional to C ESR , and in a typical sensing experiment, the readout conditions would be optimized to maximize this area.We observe that the total area between the curves for sample A is ∼3.6 times the area between the curves for sample F. We fit the data using a rate equation model for NV − spin dynamics [27][28][29][30], modifying the model to additionally incorporate NV 0 states and charge state conversion rates Γ ion and Γ rec (see a complete description of the rate equation model in Supplemental Material Sec.III [24]).In order to constrain the model, we measure PL traces for both m s = 0 and m s = −1 initial states for eight different laser powers and fit all data sets simultaneously to a single model.We constrain the fit to incorporate the same intrinsic NV parameters (excited state lifetime, singlet state lifetime, etc.) across all powers, while allowing the photo-induced transition rates to be free parameters.Representative PL traces for a subset of these powers for an NV center from sample A are shown in Fig. 4(b), and corresponding datasets from sample F are shown in Fig. 4(c).Our fitted intrinsic NV parameters are consistent with previously reported measurements [27,29].We also find that extracted excited state and singlet lifetimes are consistent across NV centers.A summary of the extracted parameters for four NV centers in each sample can be found in Supplemental Material Sec.IV [24].
From the data it is clear that there are two different contributions of charge state dynamics to decreased C ESR : the increased background due to the steady-state NV 0 population, and charge conversion rates interfering with spin dynamics.Armed with a model that quantitatively accounts for both spin and charge dynamics of the NV center, we use this model to disentangle the relative We have shown that samples with different surface conditions can exhibit drastically different charge state dynamics for shallow NV centers, and that charge state stability can influence OD-ESR contrast in two ways: first, a high steady-state population in NV 0 will increase the PL background and second, when the charge state conversion time is comparable to the NV spin polarization and readout time, ionization results in a loss of spin information and initialization, and consequently a decrease of C ESR .These effects can be drastically different for different diamond samples, as well as among NV centers in a given sample.These results point to the importance of surface preparation and engineering in utilizing NV centers for nanoscale sensing and magnetometry.On-going and future work includes establishing microscopic mechanisms for charge state instabilities, such as surface trap states and contaminants, as well as exploring optimal surface terminations [31][32][33][34]  Most diamond samples in this work are prepared following the procedure outlined in Ref. [1].For samples A, C, D and F, we started with electronic grade diamonds from Element Six which were scaife polished to a RMS roughness of less than 1 nm.We then performed inductively-coupled plasma reactive ion etching (ICP-RIE) followed by high temperature annealing at 1200 • C in vacuum to remove surface and subsurface damage.After annealing, samples were cleaned in a refluxing mixture of 1:1:1 concentrated sulfuric, nitric, and perchloric acids (triacid clean) to remove surface amorphous carbon formed during annealing.Samples B and E did not go through this polish and pre-etch process and so their surfaces were 'as-grown' before implantation.Next, NV centers were formed by nitrogen ion implantation followed by annealing at 800 • C in vacuum.Further details about implantation parameters for different samples are listed in Table S1.Afterwards, samples were triacid cleaned again.This was the last step performed for sample B. For samples A, C, D, E and F we subsequently annealed them in an oxygen atmosphere at temperatures around 440 • C − 460 • C to create well-ordered oxygen-terminated diamond surfaces.After oxygen annealing, these samples were cleaned in Piranha solution (1:2 mixture of hydrogen peroxide and concentrated sulfuric acid).

Sample
Surface condition Implantation dose (cm −2 ) Implantation energy (keV) A Polished + Pre-etched, 12   Due to poor OD-ESR contrast in sample F, we are unable to use proton NMR to measure the depths of individual NV centers.Instead we provide a comparison sample that was implanted at the same energy, sample B. Histograms of the contrast and depth distributions for sample B are presented in Fig. S1.Although sample B was implanted with the same ion implantation energy as sample F, we observe significantly higher OD-ESR contrast for NV centers in this sample, which suggests the lower contrast in sample F mainly results from the diamond surface condition.X-ray photoelectron spectroscopy (XPS) characterization of the diamond surfaces was performed with a Thermo Fisher K-Alpha spectrometer at the Imaging and Analysis Center (IAC) at Princeton University.XPS spectra for samples A, B and F are shown in Fig. S2.A peak associated with boron 1s is evident in sample F suggesting boron contamination on the surface, which has an atomic percentage of 0.31%.We estimate that this atomic percentage corresponds to approximately 4% of a surface monolayer of boron using the method discussed in Ref. [1].

II. EXPERIMENTAL SETUP
NV photoluminescence (PL) measurements are performed in a home-built confocal microscope setup.A Nikon Plan Fluor 100×, N.A. = 1.30, oil immersion objective is used for focusing the excitation laser and collecting the PL.For the green illumination, we use a 532 nm optically pumped solid-state laser (Coherent Sapphire LP 532-300), and for the orange illumination we use a NKT SuperK laser (repetition rate 78 MHz, pulse width 5 ps) with two bandpass filters with transmission wavelength around 590 nm (Thorlabs FB590-10 and Semrock FF01-589/18-25).Both lasers are modulated by acousto-optic modulators (AOMs) (Isomet 1205C-1) and the beam is scanned with X-Y galvo mirrors (Thorlabs GVS012).A dichroic mirror (Semrock BLP01-647R-25) is used to separate the excitation and collection pathways, and the PL is measured using a fiber-coupled avalanche photodiode (Excelitas SPCM-AQRH-44-FC).
For the PL spectrum measurements shown in the main text Fig. 2, we use a Princeton Instruments Monochromator (Acton SP-2300i) with a CCD camera (Pixis 100).Each spectrum is obtained by first acquiring the background signal of a non-NV spot close to the target NV and subtracting this background spectrum from the acquired NV spectrum.This technique allows us to remove background fluorescence and Raman lines of the driving laser.
All the pulse sequences used in the experiments are programmed by a PulseBlaster (SpinCore ESR-PRO500 with a timing resolution of 2 ns).For time-resolved measurements of NV PL, the avalanche photodiode signal is sent to a PicoHarp (PicoQuant PicoHarp 300 with a highest resolution of 4 ps).

III. TIME-RESOLVED MEASUREMENTS: RATE EQUATION MODEL AND DATA FITTING
In our time-resolved measurements of NV PL, a green pulse is used to initialize the NV − spin state to m s = 0. Then we either leave the NV in m s = 0 or apply a calibrated microwave π pulse to initialize the NV into m s = −1.
A subsequent green pulse is used to measure the PL which is read out by a PicoHarp to achieve high time resolution.We vary the power of both initialization and readout green pulses together to get the power dependence of the time-resolved PL.It is worth noting that this potentially changes the initial state that our model uses to predict the time-resolved PL.We address this in our modeling.
In order to quantify the optically induced spin dynamics of our NV centers in the presence of NV ionization and recombination, we expand upon a model introduced in previous works [2][3][4][5].A representation of the states of our system and the corresponding transition rates can be found in Fig. S3.In our model, we use five levels to describe the NV − spin states, which include two states to describe the ground state triplet (GS 0 and GS 1 ), two states from the excited state triplet (ES 0 and ES 1 ), and a single state for the singlet manifold (A 1 ).We also include ground and excited states for NV 0 (GS N V 0 and ES N V 0 ).This allows us to describe the system as a state vector where the entries of the vector ρ ψ describe the probability of finding the system in each of the seven states in our model.
To model the dynamics of the system, we consider the set of transitions that couple all of these spin states, ignoring spin non-conserving transitions (aside from the intersystem crossing and ionization/recombination which is addressed later).This includes a transition rate Γ 532 between the ground triplet states and the excited triplet states, a transition rate Γ ES between the excited triplet states and the ground triplet states, and a transition rate Γ ES0,1→A1 between the triplet excited states and the singlet manifold.We also include a transition rate Γ A1 for the decay from the singlet state to the GS 0 and GS 1 states where the so-called "branching ratio" P A1→GS1 describes the probability of decay into GS 1 (with 1 − P A1→GS1 describing the probability of decay into GS 0 ).Finally, we include the transition rate Γ ion from the NV − excited states (ES 0 and ES 1 ) to the NV 0 ground state (GS N V 0 ) and the transition rate Γ rec from the NV 0 excited state (ES N V 0 ) to NV − ground states (GS 1 and GS 0 ).We assume that both ES 0 and ES 1 ionize at the same rate Γ ion and that recombination from the NV 0 excited state goes into the three spin states of NV − with equal probabilities [6].With these transition rates the dynamics of the system can be described by the following coupled rate equations: which can also be expressed in a matrix form For the readout, the NV center is initialized with the same green pulse.This means that the initial charge state can be estimated by solving for the steady state of the rate equation matrix ( ρψ (t) = 0).The solution can be found as the zero eigenvector of R M , ρ ψ,S.S. .We then sum the components of ρ ψ,S.S. from the NV 0 manifold and the NV − manifold, which we call P N V 0 and P N V − , respectively.Right before the readout, we assume that the NV is in the ground state because we have waited for a time longer than the excited state lifetimes.This means that if we initialize in m s = 0 then: and if we initialize in m s = −1 then: These initial states are then evolved in time according to the rate equation matrix to give The PL measured in the experiment is predominantly from transitions from the excited state triplet to the ground state triplet of NV − (which occurs at a rate Γ ES ) or from the excited state to ground state of NV 0 (which occurs at a rate Γ ES N V 0 ).This means that given the model predictions for the population in each of the states as a function of time, we can calculate the PL as To obtain the actual PL, we need to take into account the total collection efficiency of our measurement setup.This will introduce an overall scale factor to Eq. 6 and convert the proportionality (∝) to an equality (=).However, in order to reduce the number of fit parameters in our model, we normalize both the experimental data and the predicted PL from the model.

VI. RELATIONSHIP BETWEEN Γ ion,rec AND r ion,rec
To elucidate the relationship between Γ ion,rec and r ion,rec we consider a 4-level system, which consists of the ground state and excited state for NV − (GS N V − : level 1 and ES N V − : level 2) and for NV 0 (GS N V 0 : level 3 and ES N V 0 : level 4).We use ρ i to represent the population in level i and from the conservation of total population we have We assume the conversion rate from level i to level j is Γ ij .For each NV charge state we have photon-induced excitation rates Γ 12 , Γ 34 and emission rates Γ 21 , Γ 43 .Also, we include the ionization rate from the NV − excited state to the NV 0 ground state Γ 23 and the recombination rate from the NV 0 excited state to the NV − ground state Γ 41 .
Starting from this full picture, we want to consider the conversion rates between the two charge states: ρ N V − ≡ ρ 1 + ρ 2 and ρ N V 0 ≡ ρ 3 + ρ 4 .Thus, the time derivative of the two charge state populations will be In general if the ratios ρ 2 /ρ N V − and ρ 4 /ρ N V 0 change with time, we will not observe the exponential decay or increase in the populations.However, if the charge state conversion rates are much slower than the internal transition rates for each charge state, i.e. {Γ 23 , Γ 41 } ≪ {Γ 12 , Γ 21 , Γ 34 , Γ 43 }, then in a timescale shorter than the charge state conversion timescale, the populations will reach a steady state within each charge state: This leaves us with Meaning that we can write the rate equation for two charge states as   where r ion = Γ 23 Γ 12 /(Γ 12 + Γ 21 ) and r rec = Γ 41 Γ 34 /(Γ 43 + Γ 34 ).Now that we have derived effective rate equations for the populations in each of the charge states, we can explore how each of the individual transition rates affects the overall transition rate between charge states.The excitation from the ground state to the excited state, ionization, and recombination are each one-photon processes, which means that the rates Γ 12 , Γ 34 , Γ 23 , and Γ 41 should each have linear power dependence.Γ 21 and Γ 43 , which are the spontaneous emission rates analogous to Γ ES and Γ ES N V 0 , are decided by the intrinsic properties of NV centers and should not depend on power.If the laser power is much lower than the saturation power then Γ 12 ≪ Γ 21 and Γ 34 ≪ Γ 43 .This means that r ion and r rec will each have a quadratic power dependence.If the laser power is much higher than the saturation power, then Γ 12 ≫ Γ 21 and Γ 34 ≫ Γ 43 .This means that r ion and r rec should have linear power dependence, which agrees with the measurement results in Ref. [8].This analysis demonstrates an intuitive relationship between individual transition rates in our rate equation model and r ion and r rec in various parameter regimes.It also qualitatively explains the power dependence of r ion and r rec in the main text Fig. 3.It is worth noting, however, that the model here is a minimal model of the full rate equation model used in this paper but we believe the same intuition holds for the full model albeit with slightly more cumbersome expressions.

FIG. 1 .
FIG. 1. Variation in OD-ESR contrast for shallow NV centers.(a) Cartoon depicting the interaction between shallow NV centers and defects at the diamond surface.Defects can act as charge traps, leading to changes in ionization and recombination kinetics for nearby NV centers.(b) OD-ESR contrast, CESR, as a function of depth for NV centers across five diamond samples.A large spread and an average decrease in contrast is apparent within around 10 nm of the surface.(c) Histograms of CESR from the two selected samples, A (blue bars) and F (orange bars).99 NV centers in sample A and 108 NV centers in sample F were measured under the same experimental conditions including microwave pulse power and green laser power.An external magnetic field of approximately 27 G was aligned to the NV axis.

FIG. 2 .
FIG.2.Steady-state charge distributions.(a) PL spectrum under green illumination for representative NV centers in sample A (blue trace) and sample F (orange trace).The higher level of PL between 575 nm (dashed line) and 637 nm (dash-dotted line) in sample F indicates a higher steady-state population in NV 0 .(b) Charge state conversion in the dark.NV centers are initialized into the negative charge state by a green pulse, and the charge state is read out using a chargestate-selective orange pulse with variable delay time.The PL is normalized to the value for the shortest dark time interval.Each NV center in sample F (orange triangles) exhibits a decay in the relative NV − population with a time constant in the range of 11 -300 ms, while NV centers in sample A (blue circles) exhibit much less decay out to 1 s.Solid lines are fits to the data.(c) Stability of NV 0 in the dark in sample A. NV centers are initialized into NV 0 using a green initialization pulse followed by an orange pulse, and then read out with another orange pulse after a variable delay time.The PL is normalized to the value after the green pulse.No change in the population is observed out to 1 s, indicating that recombination kinetics are also slow in sample A.

)FIG. 3 .
FIG. 3. Charge state dynamics for NV centers in sample A (blue circles) and sample F (orange triangles).In both plots, points are experimental data and lines are fit using a quadratic function.(a) Power dependence of the charge state conversion rates under green illumination.(b) Power dependence of the charge state conversion rates under orange illumination.Laser powers are measured at the back of the objective.Under both illumination conditions, the charge state conversion rates are faster for each NV center in sample F compared to those in sample A.

FIG. 4 .
FIG. 4. Time-resolved PL traces and fitting to rate equation model.(a) Comparison of PL traces for representative NV centers from the two samples, A (blue) and F (orange), under 664 µW of green illumination.PL counts are normalized to steady-state values for each trace.The shaded area is proportional to the OD-ESR contrast.(b,c) Power dependence of NV PL traces from sample A (b) and sample F (c) after initialization into the NV − ms = 0 state and ms = −1 state.Each trace represents a different green illumination power (from bottom to top: 330 µW, 450 µW, 660 µW).Solid lines indicate fits to the data using the rate equation model.(d) Simulated contrast as a function of NV 0 population using extracted parameters from the fit to the same NV as in (c), where the upper bound (dark yellow points) represents the highest attainable contrast for a particular NV 0 population and the shaded region represents lower CESR caused by faster charge state conversion.The red cross indicates the extracted CESR and NV 0 population for this NV center, and the highest attainable contrast for this NV 0 population is indicated by the red circle.The total decrease in contrast compared to CESR with no charge dynamics is parameterized by the change due to the finite NV 0 population of 0.486, ∆CESR,static, and the change due to charge conversion rates, ∆C ESR,dynamic .(e) CESR as a function of Γion and Γrec for a fixed NV 0 population (±1%).The red cross indicates the extracted Γion and Γrec and the red circle indicates the point of highest CESR.
FIG. S1: Characteristics of NV centers in sample B, which was implanted at the same energy as sample F but has an as-grown surface.(a) Histogram of OD-ESR contrast, CESR.(b) Histogram of measured NV depths.
FIG. S2: XPS characterization of surfaces for samples A, B and F. (a) XPS survey scans showing dominant carbon and oxygen peaks, shifted vertically for clarity.Samples A and F were annealed in an oxygen atmospher and show higher oxygen 1s peaks.(b) High-resolution XPS boron 1s spectrum showing a small but clear peak in sample F. Intensity is normalized by setting the height of carbon 1s peak to 1.
FIG.S3: NV level diagram showing the states considered in our model, which includes the ground state triplet ms = 0 state (GS0), the ground state triplet ms = ±1 states (GS1), the excited state triplet ms = 0 state (ES0), the excited state triplet ms = ±1 states (ES1), the singlet state A1 and the NV 0 ground (GS NV 0 ) and excited (ES NV 0 ) states, along with all corresponding transition rates.
FIG. S4: (a), (b) Predicted OD-ESR contrast CESR as a function of Γion and Γrec using the extracted parameters from representative NV centers from main text Fig. 4 from sample A and sample F, respectively.(c), (d) Predicted SNR as a function of Γion and Γrec for the same NV centers.(e), (f) Steady-state NV 0 population based on Γion and Γrec for the same NV centers.(g), (h), Collected points from data in (a) and (e) or (b) and (f) that show CESR as a function of the NV 0 population.The finite lower bound of CESR is a result of sampling from a finite range of Γion and Γrec.Red crosses correspond to the fitted Γion and Γrec for the two NV centers.Based on the experimental contrast and the extracted parameter values we calculate a change in a change in contrast due to finite NV 0 population, ∆CESR,static, and due to the magnitude of the ionization and recombination rates ∆C ESR,dynamic .
for stabilizing the negative charge state of shallow NV centers.(Agreement No. D18AC00015).M. F. was supported by an appointment to the Intelligence Community Postdoctoral Research Fellowship Program by Oak Ridge Institute for Science and Education (ORISE) through an interagency agreement between the U.S. Department of Energy and the Office of the Director of National Intelligence (ODNI).L. V. H. R. acknowledges support from the Department of Defense through the National Defense Science and Engineering Graduate Fellowship Program.Material for "Charge State Dynamics and Electron Spin Resonance Contrast of Shallow Nitrogen-Vacancy Centers in Diamond" Zhiyang Yuan, Mattias Fitzpatrick, Lila V. H. Rodgers, Sorawis Sangtawesin, Srikanth Srinivasan, and Nathalie P. de Leon *

TABLE S2 :
Fitted parameters based on the rate equation model incorporating ionization and recombination.Listed CESR, P NV 0 , Γ532, Γion, and Γrec are the values at the maximum CESR across p532 (and resultant Γ532).We also include effective lifetimes ES0 τ = 1/(ΓES + ΓES 0 →A 1 ), ES1 τ = 1/(ΓES + ΓES 1 →A 1 ), and A1 τ = 1/ΓA 1 .It is worth noting that the actual excited state lifetime should also include ionization and recombination but we have selected these definitions to compare to models without ionization and recombination.ΓES is constrained to be < 75 MHz and we are saturating the bound in our fit procedure.

TABLE S3 :
Fitted parameters based on the rate equation model without ionization and recombination.