Photoluminescence decomposition analysis: a technique to characterize NV creation in diamond
Scott Alsid John Barry Linh Pham Jennifer Schloss Michael O'Keefe Paola Cappellaro Danielle Braje


APS/123-QED




Lincoln Laboratory, Massachusetts Institute of Technology, Lexington, Massachusetts 02421, USA
Department of Nuclear Science and Engineering, Massachusetts Institute of Technology, Cambridge, Massachusetts 02139, USA

[Corresponding author.
]john.barry@ll.mit.edu
Lincoln Laboratory, Massachusetts Institute of Technology, Lexington, Massachusetts 02421, USA

Lincoln Laboratory, Massachusetts Institute of Technology, Lexington, Massachusetts 02421, USA

Department of Physics, Massachusetts Institute of Technology, Cambridge, Massachusetts 02139, USA
Center for Brain Science, Harvard University, Cambridge, Massachusetts 02138, USA

Lincoln Laboratory, Massachusetts Institute of Technology, Lexington, Massachusetts 02421, USA

Department of Nuclear Science and Engineering, Massachusetts Institute of Technology, Cambridge, Massachusetts 02139, USA

Lincoln Laboratory, Massachusetts Institute of Technology, Lexington, Massachusetts 02421, USA




Treatment of lab-grown diamond by electron irradiation and annealing has enabled quantum sensors based on negatively-charged nitrogen-vacancy (NV^-) centers to demonstrate record sensitivities. 
Here we investigate the irradiation and annealing process applied to 28 diamond samples using a new ambient-temperature, all-optical approach. As the presence of the neutrally-charged nitrogen-vacancy (NV^0) center is deleterious to sensor performance, this photoluminescence decomposition analysis (PDA) is first employed to determine the concentration ratio of NV^- to NV^0 in diamond samples from the measured photoluminescence spectrum. The analysis hinges on (i) isolating each NV charge state's emission spectrum and (ii) measuring the NV^- to NV^0 emission ratio, which is found to be 2.5±0.5 under low-intensity 532 nm illumination. Using the PDA method, we measure the effects of irradiation and annealing on conversion of substitutional nitrogen to NV centers. Combining these measurements with a phenomenological model for diamond irradiation and annealing, we extract an estimated monovacancy creation rate of 0.52± 0.26 cm^-1 for 1 MeV electron irradiation and an estimated monovacancy diffusion coefficient of 1.8 nm^2/s at 850 ^∘C. Finally we find that irradiation doses ≳ 10^18 e^-/cm^2 deteriorate the NV^- decoherence time T_2 whereas T_1 is unaffected up to the the maximum investigated dose of 5× 10^18 e^-/cm^2.



                             

                              
Photoluminescence decomposition analysis: a technique to characterize NV creation in diamond
    Danielle A. Braje
    December 30, 2023
============================================================================================






§ INTRODUCTION


The negatively-charged nitrogen vacancy (NV^-) center, with its rich spin dynamics even at room temperature, has captured the interest of the sensor community. Importantly, quantum sensors based on NV^- ensembles in diamond require efficient conversion of substitutional nitrogen to NV^- centers to achieve peak performance <cit.>. One common approach to creating NV^--rich diamond material involves electron-irradiating and subsequently annealing nitrogen-doped diamond <cit.>. However, maximizing the NV^- population while minimizing unwanted defects (e.g., neutral NV^0 centers, vacancy complexes, and other paramagnetic impurities) is critically dependent on employing proper irradiation and annealing processes.

Although high-quality diamond processing studies are available <cit.> to guide post-growth treatment, literature discrepancies (see Appendix Sec. <ref> and Refs. <cit.>) and a shortage of systematic studies impede full understanding of these post-growth treatments affect defect composition. In particular, variation in the base diamond material employed in the literature studies constitutes a prominent and continuing obstacle to optimizing treatment parameters. For example, diamond origin (natural or synthetic), growth method (CVD or HPHT), nitrogen concentration <cit.>, lattice defect density, impurity concentration <cit.>, and strain <cit.> have all been shown to affect the irradiation and annealing process. As diamond synthesis capabilities develop, producing material with a range of compositions, a comprehensive study of irradiation and annealing will enable targeted diamond engineering for diverse quantum sensing applications <cit.>.

To fine-tune the diamond treatment process, accurate measurements of the effects of irradiation and annealing on key diamond material properties are needed. In particular, optimizing NV-diamond quantum sensors necessitates post-treatment characterization of (i) the conversion efficiency from nitrogen to NV centers, (ii) the NV^- to NV^0 charge state ratio, and (iii) the relevant NV^- spin relaxation times. Such measurements can also provide physical insight into the vacancy creation and diffusion dynamics central to formation of NV^- centers. 
Here we present an all-optical, ambient-temperature technique for rapidly characterizing the NV^- to NV^0 charge state ratio and total NV concentration. This technique provides an appealing alternative to slower conventional methods, such as electron paramagnetic resonance (EPR) and ultraviolet-visible (UV-Vis) spectroscopy. Further, the technique is suitable for analyzing NV-doped layers or other small ensembles, regimes typically precluded from EPR or UV-Vis analysis by signal-to-noise requirements. Additionally, accurate charge state ratio measurements hinge on knowledge of the photoluminescence (PL) rate of NV^- relative to NV^0 for equal concentrations. Here we determine for the first time this PL ratio under low-intensity 532 nm excitation.


We demonstrate the presented technique's utility in a systematic study of irradiation and annealing on 28 diamonds. We focus on a diamond base material that is commercially available and relevant to many quantum sensing applications <cit.>. To support additional analysis of our measurements, we introduce a phenomenological model for the irradiation and annealing process and extract an estimated irradiation-induced vacancy creation rate of 0.52± 0.26 cm^-1 for 1 MeV electron irradiation. From the model, we also determine an estimated monovacancy diffusion coefficient of 1.8 nm^2/s at 850 ^∘C. Further, we examine the effects of irradiation and annealing on the NV^- coherence time T_2 and longitudinal relaxation time T_1 in the employed base material. 
Finally, we comment on applications of the methods and analysis detailed in this work to furthering NV physics research and diamond materials science. Example near-term studies include characterization of NV^0 optical dynamics and investigation into the effects of base material variation on vacancy migration.




§ DIAMOND MATERIAL OPTIMIZATION


The presence of various electronic spin species in diamond (e.g., NV^-, NV^0, substitutional nitrogen N_S, and other paramagnetic impurities) can limit the achievable sensitivity of an NV-diamond quantum sensor; this effect motivates the careful engineering of diamond material to contain optimal concentrations of each of these spin defects. In particular, spin defect concentrations affect three key parameters that can limit sensitivity: measurement time, measurement contrast, and number of photons collected per measurement. For example, measurement time is limited by the dephasing, coherence, or relaxation times of the NV^- spins, which in turn are typically limited by other spin species generating magnetic noise in the vicinity of the sensing NV^- centers. Measurement contrast can be degraded by the presence of non-information-containing background fluorescence emitters like NV^0.

The number of photons collected per measurement is dictated by the optical collection efficiency of the experimental apparatus and ultimately by the number of NV centers probed during the measurement. 

In essence, the sensitivity of an NV-diamond quantum sensor can generally be enhanced by maximizing the number of NV^- centers while minimizing the presence of all other spin defects. Since nitrogen is an integral constituent of the NV^- center, nitrogen-containing defects often dominate in diamond material employed for NV-ensemble-based sensors. Thus, optimal values of the NV conversion efficiency χ≡ [NV^T]/[N^T] are desired to achieve optimal sensitivity in these devices. Here [NV^T] indicates the total NV population, summed over all charge states, and [N^T] denotes the total nitrogen concentration in all forms including all charge states (e.g., N_S^0, N_S^+, NV^0, NV^-, etc.).

Note that there is an important distinction between the NV^- centers, which are the fundamental sensing qubits that comprise the NV-diamond quantum sensor, and NV^0 centers. The presence of NV^0 degrades sensor performance via several mechanisms. First, the photoluminescence from NV^0 spectrally overlaps with that of NV^-, contributing background signal and thereby reducing measurement contrast. This effect is particularly deleterious when [NV^0] ≳ [NV^-] <cit.>. Additionally, the NV^0 center's electronic S=1/2 paramagnetic ground state <cit.> may increase magnetic noise in the diamond and thereby degrade NV^- spin-coherence times. 
For these reasons, maximization of the charge state efficiency ζ≡ [NV^-]/[NV^T] is also desirable. 

A simple model suggests NV^- formation requires two nitrogen atoms and a vacancy. One nitrogen atom captures the vacancy to form the NV center, while the other nitrogen donates a valence electron to negatively charge the NV center <cit.>. For typical CVD diamonds, such as those employed in this work, the native grown-in vacancy concentration is too low to achieve optimal N^T-to-NV^T conversion efficiency χ∼ 1/2. One straightforward approach to improve χ in these diamonds is to perform electron irradiation to create vacancies in the crystal lattice. A subsequent annealing step mobilizes the vacancies and enables their capture by N_S defects to form NV centers. Achieving the optimal spin defect composition in the diamond following this two-step process is critically dependent on optimizing irradiation and annealing parameters. 
For example, too low an irradiation dose will introduce too few vacancies relative to N_S defects, limiting the attainable post-anneal NV conversion efficiency χ, which in turn limits the number of NV^- sensors available to perform measurements. Furthermore, the remaining under-utilized paramagnetic N_S defects can degrade NV^- spin-coherence times. On the other hand, too high an irradiation dose will introduce too many vacancies relative to N_S, resulting in overproduction of undesirable neutral NV^0 centers <cit.> and thereby degrading the charge state efficiency ζ.

Proper selection of irradiation energy and dose, annealing time and temperature, and other processing parameters is imperative for achieving NV-diamond optimized for high-sensitivity quantum sensing. 
Consequently, quantitative characterization of key material properties, such as the NV^- to NV^0 charge state ratio, is likewise important in determining these optimal processing parameters.
Conventional techniques for quantitatively characterizing relevant spin defect concentrations include ultraviolet-visible (UV-Vis) spectroscopy, which measures [NV^0] and [NV^-],
and electron paramagnetic resonance (EPR) spectroscopy, which measures [NV^-] and [N_S]. Unfortunately, UV-Vis spectroscopy is performed at cryogenic temperatures and thus requires additional time for thermal cycling. EPR spectroscopy of typically small diamond volume likewise requires significant averaging time. Since both of these conventional techniques typically take ≳1 hour to perform measurements on a single sample, characterizing the large number (≳10) of samples typically involved in optimizing processing parameters can be quite onerous.
Furthermore, both techniques employ specialized equipment, which is not readily accessible in many NV-diamond research labs.



§ EXPERIMENTAL DESIGN

This study employs 28 commercially available diamonds (Element Six, 145-500-0274-01) grown via chemical vapor deposition (CVD). The chip-sized samples are 2.6 mm × 2.6 mm × 0.3 mm, with <100> edges and {100} face orientations. The diamond material contains carbon in natural isotopic abundance (1.1% ^13C, 98.9% ^12C) and exhibits a mean total nitrogen concentration [N^T] ∼ 118 ppb (see Appendix Sec. <ref>). Prior to irradiation, the samples were pre-annealed at 850 ^∘C for one hour in a 760 Torr nitrogen atmosphere.

Diamond irradiation is performed at a commercial facility (Prism Gem LLC) using a 1 MeV electron beam with a 1 cm nominal diameter. Diamonds are placed on a water-cooled plate maintained at a temperature 66 ^∘C or below; the diamonds are estimated to remain below 120 ^∘C during irradiation, suggesting that interstitial-vacancy recombination during irradiation can be neglected <cit.>. Twenty-five diamonds received doses systematically varied between 1 × 10^15 e^-/cm^2 and 5 × 10^19 e^-/cm^2; the other three diamonds remained un-irradiated and were retained as controls. Nineteen of the irradiated diamonds were subsequently annealed at 850 ^∘C in a 760 Torr nitrogen atmosphere for one hour at the same commercial facility. The other six irradiated diamonds and a subset of the 850 ^∘C annealed diamonds were then annealed at 1250 ^∘C under 100 Torr hydrogen gas for one hour. Each diamond's specific irradiation and annealing parameters are given in the Appendix Sec. <ref>.

Diamond characterization is performed on NV-ensembles using a custom-built confocal microscope with attached spectrometer (See Appendix Sec. <ref>). While this confocal microscope follows standard designs <cit.> in most ways, a notable modification allows enhanced functionality for NV ensemble excitation power/intensity studies. Specifically, the microscope employs a detection volume centered within a much larger excitation volume, ensuring all NV centers within the detection volume receive similar optical illumination intensity (see Appendix Sec. <ref>, Fig. <ref>). Without this modification, the intensity gradient associated with the wings of a Gaussian excitation profile can complicate interpretation of intensity-induced dynamics, such as photoionization and charge conversion.

A spectrometer attached to the confocal microscope allows spectral characterization of the photoluminescence (PL) light from 560 nm to 818 nm. This range allows collection of nearly all emitted PL light from NV^- and NV^0; analysis of PL data in the literature from other groups <cit.> suggests this wavelength range covers ≳99% of NV^- and NV^0 photoluminescence.



§ CHARGE STATE PHOTOLUMINESCENCE STUDY



In this section, we present our method for quantitatively determining the charge state efficiency ζ using PL spectroscopy. 
Note that, historically, approximations of ζ have been extracted from PL spectra <cit.> using the Debye-Waller factor S̃ <cit.> (see Appendix Sec. <ref>).

However, this method is non-ideal because it results in large uncertainties in the extracted ζ and does not account for the difference in PL emission of the different charge states. We employ an alternative method for isolating the PL contributions due to NV^- and NV^0 for a given PL spectrum. We then use this alternative decomposition method to determine the ratio of PL light produced by a single NV^- to that of a single NV^0 per unit time. The charge-state photoluminescence ratio of NV^- relative to NV^0 is denoted κ_λ for optical excitation at wavelength λ (in this work, λ=532 nm) Using laser intensities where NV saturation effects can be neglected, the concentration of each charge state is then proportional to its respective PL when correcting for κ_532.  Therefore the charge state efficiency may be re-written as

    ζ≡[NV^-]/[NV^T]= [NV^-]/[NV^-]+[NV^0] = c_-/c_-+κ_532c_0,

where c_- and c_0 denote the fractional contributions to total PL from NV^- and NV^0 respectively. In the absence of PL from defects besides NV^0 and NV^-, c_-+c_0=1.  Post-analysis of our data shows little if any PL contributions in the 550-818 nm range from defects other than NV^- and NV^0 for irradiation doses ≲ 1 × 10^19 e^-/cm^2.



 §.§ Photoluminescence spectrum decomposition


To determine the charge state efficiency ζ using Eqn. (<ref>), the fractional contributions of each charge state in a given PL spectrum (c_0, c_-) are determined by fitting the acquired spectrum I(λ) 

to a linear combination of NV^- and NV^0 emission:  

    I(λ) = M[ c_-Î_-(λ)+c_0 Î_0(λ)].

Here, M is a proportionality constant, and Î_-(λ) and Î_0(λ) denote area-normalized PL emission spectra for NV^- and NV^0 respectively; the normalization satisfies ∫_λÎ_-(λ) dλ = 1 and ∫_λÎ_0 (λ) dλ=1, where the integral is taken over the wavelengths measured by the spectrometer. Note that, throughout this paper, the hat denotes dimensionless spectra which are normalized to unity area.

 

 

We now describe the process employed to determine Î_0(λ) and Î_-(λ), the accuracy of which are critical for the decomposition denoted in Eqn. (<ref>) to be effective.  
PL spectra are recorded for all 28 diamonds using 18 different 532 nm  optical excitation intensities between 1 kW/cm^2 and 4 MW/cm^2 (see Appendix Sec. <ref> for intensity calculations). 
The PL spectra are area-normalized and sorted by the fractional PL contribution in the 560-600 nm band, which includes the NV^0 zero-phonon line (ZPL) at 575 nm. 
The area-normalized PL spectrum with the greatest 560-600 nm content is denoted Î_0^pre(λ) and the area-normalized PL spectrum with the least 560-600 nm content is denoted Î_-^pre(λ). 
The superscripted pre label denotes the spectra are preliminary and require additional processing, as now  described. The spectra Î_0^pre(λ) and Î_-^pre(λ) are nominally dominated by PL emission from NV^0 and NV^- respectively. 
Post-analysis investigation shows Î_0^pre(λ) contains ∼94% NV^0 emission and ∼6% NV^- emission, while  Î_-^pre(λ) contains ∼84% NV^- emission and ∼16% NV^0 emission. 
We assume that PL from other defects (SiV, V^0, etc.,) can be ignored; no ZPLs from such defects were observed for irradiation doses ≲ 10^19 e^-/cm^2. 

The un-normalized NV^- I_-(λ) function is constructed using

    I_-(λ) = Î_-^pre(λ)- δ_0 Î_0^pre(λ),

where the value of δ_0 is adjusted so that I_-(λ) exhibits no evidence of the NV^0 ZPL feature at 575 nm. The I_-(λ) function is then area-normalized to yield Î_-(λ). The Î_-(λ) function is determined to be free of any NV^0 contribution to the 1% level or better (See Appendix Sec. <ref>). The un-normalized NV^0 Î_0(λ) function is then similarly constructed using

    I_0(λ) = Î_0^pre(λ)- δ_-Î_-(λ),

where the value of δ_- is adjusted so the resulting I_0(λ) exhibits no evidence of the NV^- ZPL feature at 637 nm. This I_0(λ) function is then area-normalized to yield Î_0(λ).  Ensuring the Î_0(λ) function is free of any NV^- contribution is non-trivial; the NV^- ZPL coincides with near-peak NV^0 vibronic sideband emission and is thus difficult to isolate and remove. Nevertheless, δ_- is determined to be 6.4%± 3.3%, suggesting the resultant Î_0(λ) function should exhibit at most a few percent contribution from NV^-. The extracted Î_-(λ) and Î_0(λ) functions for bulk NV-diamond are shown in Figs. <ref> and <ref>, respectively.

As an independent check of the extracted Î_-(λ) basis function, we adjust the low-temperature, one-phonon emission spectrum of the ^3A_2 ← ^3E transition <cit.> for higher temperatures and multi-phonon transitions, as outlined in Ref. <cit.>.

Specifically, the 300 K trace (Fig. 8) from Ref. <cit.> is digitized <cit.> and adjusted by the necessary ω^3 weighting <cit.> to yield an expected room temperature PL spectrum (less any ZPL contribution). The results are shown in Fig. <ref> and are in good agreement with our experimental determination of Î_-(λ). As the spectrum derived from Ref. <cit.> excludes the ZPL contribution, a comparison of the expected and measured data yields residuals whose main feature is the NV^- ZPL; this residual data is well-fit by a Lorentzian with a Debye-Waller factor of 4.3 and a full width at half maximum (FWHM) linewidth of 2.6 nm. These values are well within the range of values observed in the literature <cit.> (see Table <ref>). 

Using the determined Î_-(λ) and  Î_0(λ) basis functions, the PL spectra for all 28 diamonds at the 18 investigated optical excitation intensities are decomposed via Eqn. <ref>. An example decomposition is shown in Fig. <ref>. To verify our approach, we compare the PDA method presented here to the Debye-Waller decomposition and observe that the PDA method is more robust to noise and produces more consistent results (Fig. <ref>). 
This improvement is likely due to the PDA method harnessing the entire ∼200 nm PL spectrum, compared to the Debye-Waller decomposition, which utilizes only data in the ZPL vicinity (a few nm). Further information on the Debye Waller decomposition is detailed in Appendix Sec. <ref>.  As an additional benefit, the PDA method does not require a priori knowledge of the Debye-Waller factors, which are observed to vary in the literature (see Table <ref>). 

 





 §.§ Charge state photoluminescence ratio


Accurate determination of the charge state efficiency ζ using Eqn. <ref> requires first determining κ_532, the relative PL emission rate of NV^- to NV^0 under low-intensity (i.e. well below saturation for either species) 532 nm optical excitation. Equivalently, κ_λ can be formally defined as

    κ_λ = [∫_λ' I_-(λ')dλ'/∫_λ' I_0(λ')dλ']_[NV^0]=[NV^-],

where κ_λ is evaluated in the low optical excitation intensity limit. 

We estimate the value of κ_532 using three independent methods.  First, the total integrated PL is proportional to

    ∫_λ I(λ)  dλ ∝ ([NV^0]+[NV^-]κ_532) 2P/π w_0^2,

where P and w_0 are the optical excitation power at the confocal volume and the 1/e^2 intensity radius at the origin, respectively. The factor 2P/π w_0^2 corresponds to the peak intensity of the exciting laser beam (see Appendix Sec. <ref>). The concentrations [NV^0] and [NV^-] are constrained at any given time to satisfy 

    [NV^0]+[NV^-]=constant

for a given sample. Note that Equation <ref>  implicitly assumes [NV^+] can be neglected. Even well below saturation intensity, the ratio [NV^-]/[NV^0] varies as a function of the applied 532 nm optical excitation  intensity, as shown in Fig. <ref>. Therefore, by varying the 532 nm optical excitation intensity to vary the ratio [NV^-]/[NV^0], the value of κ_532 can be determined. This method was employed for seven diamonds where c_0 ≈ c_- for 532 nm laser intensities < 20 kW/cm^2 (i.e. well below saturation). Analyzing the combined data, we extract κ_532 = 2.5±0.5; this result represents the first reported estimation of the relative PL emission rate of NV^- to NV^0 and is one of the key results of this paper.

Second, the above-derived value for κ_532 is checked using third party data where the [NV^-]/[NV^0] ratio was varied via electrical control. The authors of Ref. <cit.> fabricate NV centers inside the insulating layer of a diamond-based PIN junction diode. By adjusting the bias voltage applied across the device, the authors directly controlled the [NV^-]/[NV^0] ratio (see Ref. <cit.> Fig. 4). Using PL spectra taken for a variety of bias voltages (obtained directly from the authors of Ref. <cit.>), the charge state PL ratio κ_532 is estimated (see Appendix Sec. <ref>). This secondary method finds κ_532 = 2.17. We note that this method implicitly assumes that the electric field used to vary the [NV^-]/[NV^0] ratio does not affect the dynamics of NV^- or NV^0 beyond their charge state ratio. Both methods implicitly assume that [NV^0] + [NV^-] is conserved. 

Third, as a final check on the determined value of κ_532, we evaluated the [NV^-]/[NV^0] ratio for a test diamond not part of the rest of this study. UV-Vis measurements on this diamond suggest [NV^-]/[NV^0]=1.33 ± 0.2. PL measurements at 1 kW/cm^2 and 2 kW/cm^2 intensity find the c_-/c_0 ratio to be = 2.6 and 2.1 respectively; a linear extrapolation (which we note is likely to underestimate c_-/c_0) to zero illumination suggests c_-/c_0= 3.1 in the dark. This analysis yields an estimate of κ_532 = 2.3 ± 0.2, although we note this linear extrapolation is likely to underestimate the actual value of κ_532. As all three methods are in good agreement, we estimate κ_532 = 2.5±0.5.

The  determination of the relative florescence ratios of NV^- to NV^0 is expected to aid in uncovering presently unknown dynamics of the NV^0, which is less well studied than NV^-. For example, this ratio may help determine the NV^0 ^2A state stimulated emission cross section or the radiative quantum efficiency. 




§ ANNEALING AND NITROGEN CONVERSION STUDY


For CVD-grown diamonds, the fraction of total nitrogen incorporated as NV centers during growth is typically quite low, with reported values ranging roughly from 0.0007 to 0.03 <cit.>. 
As discussed previously, this low as-grown NV conversion efficiency χ can be increased through electron irradiation and subsequent annealing. At the MeV-level energies typically employed in electron irradiation, NV creation via this mechanism depends significantly on irradiation dose, with behavior that can be simply described in the following three regimes:



  * For irradiation doses which introduce substantially fewer vacancies than the grown-in native NV population, the observed quantity  [NV^T] is predominantly independent of the electron irradiation dose D_e.

  * At very high irradiation doses, the quantity [NV^T] saturates as every nitrogen is paired with a vacancy (to an order unity correction factor <cit.>).

  * At certain intermediate irradiation doses, the value of [NV^T] is expected to vary linearly with the electron irradiation dose D_e <cit.>. In this regime, NV centers created by irradiation and annealing outnumber the native NV centers, yet [V] ≪  [N^T] so a given additional vacancy is expected to create an additional NV.




We propose a simple phenomenological model to approximate NV annealing dynamics consistent with the above three limiting cases,

    [NV^T]    = [NV^T]_gro
             +(F_c[N^T]-[NV^T]_gro)(1-e^- B D_e/F_c[N^T]-[NV^T]_gro)

where [NV^T]_gro is the total native NV concentration immediately following growth (i.e. before irradiation and annealing); the parameter B is the vacancy creation rate per unit length; and the constant F_c accounts for the empirical observation that [NV^T]/[N^T] saturates at 0.6 ± 0.07 <cit.> in the limit of high irradiation dose. The following analysis therefore employs the estimate F_c=0.6 with the caveat that the true value of F_c likely depends on the underlying diamond material and annealing conditions. For the diamonds in this study, [N^T] ≳ 100× [NV^T]_gro; as a result, Eqn. (<ref>) may be simplified to <cit.>:

    [NV^T] ≈ [NV^T]_gro+F_c[N^T](1-e^- B D_e/F_c[N^T]).


To test the validity of Eqn. (<ref>), a correlate of the quantity [NV^T] is determined for each diamond as described herein. For each diamond, the PL spectrum is recorded under non-saturating 532 nm illumination (≤ 20 kW/cm^2) and is decomposed into the NV^- and NV^0 components (see Sec. <ref>). Using the relative charge state PL emission ratio κ_532 at low excitation intensity and the fractional distribution between NV^- and NV^0 in the PL (c_0, c_-), we adjust the PL contributions into a correlate of [NV^T]. That is, we set [NV^T] =A(c_-+κ_532c_0), where A is proportional to the total PL. Such treatment is expected to be justified in the low intensity limit where the PL of NV^- and NV^0 are both proportional to the optical excitation intensity (see Appendix Sec. <ref>). An example plot depicting [NV^T] for all 28 diamonds  is shown in Fig. <ref> for a 532 nm intensity of 6 kW/cm^2. 

The annealing data shown in Fig. <ref> illustrate that for intermediate irradiation doses ∼ 10^17 e^-/cm^2, 
the 1250 ^∘C one-hour anneal was more effective at creating NV centers than the 850 ^∘C one-hour anneal. The results suggest that the 1 hour 850 ^∘C anneal was insufficient for most vacancies to diffuse to a substitutional nitrogen. In a simple model where non-nitrogen vacancy traps are ignored, we measure an approximate diffusion coefficient of D = 1.8 nm^2/s (see Appendix Sec. <ref> for additional details). 

Examination of the fit parameters in the phenomenological model given by Eqn. (<ref>) can give further information about the annealing and irradiation process. In particular, the data for samples annealed at 1250 ^∘C corresponds to a fitting parameter of B = 0.52 ± 0.26 cm^-1 for [N^T]=118 ppb. This value of B represents the vacancy creation rate for electron irradiation at 1 MeV when all recombination, including that which occurs during subsequent annealing, is accounted for. The value of B determined here is approximate; if there are substantial vacancy sinks, the obtained value of B will underestimate the true vacancy creation rate. However, conversely, the value of B can then be used as a probe of the prevalence of vacancy sinks in diamond material of interest. Finally, we note that the value of B obtained here falls within the range of values reported in the literature, 0.2 cm^-1 to 1.5 cm^-1 for electron energies between 1 and 5 MeV. Independent measurements using UV-Vis suggest a vacancy creation rate of 0.30±0.05 cm^-1,  consistent with the measured value of B. See Appendix Sec. <ref> for further discussion.

Additionally, Fig. <ref> shows that the fits for the 850 ^∘C and 1250 ^∘C anneals asymptote to different values of [NV^T], with the former resulting in ≈ 50% higher value than the latter. A comparison of means test yields a p-value of 0.05. Statistical significance aside, one plausible interpretation is that the increased mobility of one or more species at 1250 ^∘C is detrimental to the formation (or preservation) of NV centers. Another plausible hypothesis is that NV formation is less favorable at higher temperatures due to thermo-chromic charge transfer effects <cit.>.



§ T_1 AND T_2 RELAXATION STUDY


Following initialization into a given spin state, NV spins relax into the thermally mixed state from spin-lattice-induced interactions or other mechanisms <cit.>. The characteristic T_1 longitudinal relaxation time is particularly relevant for T_1 relaxometry measurements <cit.>. The value of T_1 is evaluated for a randomly selected diamond subset (14 of the 28 samples) after irradiation and annealing. All 14 diamond samples exhibit values of T_1 between 4.8 and 6.5 ms.  For doses ≤ 5× 10^18 e^-/cm^2,  the data suggest T_1 exhibits little if any dependence on irradiation dose, as shown in Fig. <ref>.





The T_2 coherence time is an important parameter for AC magnetometry <cit.>,  quantum information <cit.>, nanoscale NMR <cit.>, single protein sensing <cit.>, and other applications involving NV-diamond quantum sensors. The T_2 of a given species is determined in part by the surrounding paramagnetic spin bath <cit.>, and can be approximated as

    1/T_2≃∑_X1/T_2,NV^-/X+1/2 T_1,

where T_2,NV^-/X corresponds to the T_2 coherence time limited by NV^- dipolar interactions with species X, T_1 is the longitudinal relaxation time <cit.>, and the summation is evaluated over all paramagnetic species X in the diamond. Similar to ion implantation <cit.>, the electron irradiation process can introduce paramagnetic impurities into the diamond lattice <cit.>, which degrade the T_2 coherence time. While annealing can partially mitigate this damage, empirical observations suggest some lattice damage will always remain <cit.>. Therefore, the irradiation doses at which T_2 and T_1 are compromised are quantities of interest for many applications employing NV ensembles.

The value of  T_2 is measured for all 28 diamonds using a Hahn-Echo sequence (See Appendix Sec. <ref>) and is plotted versus irradiation dose in Fig. <ref>. To characterize the effects of the irradiation dose D_e and annealing on T_2, the data is first fit to a linear model of the form 1/T_2 = 1/W_1[1+D_e/W_2].  We determine W_1=436± 11 s, which can be interpreted as the T_2 value in the absence of irradiation. However the poor fit of this linear model at higher irradiation doses suggests the model does not capture the relevant dynamics. Instead, fitting the data to the quadratic model 1/T_2 = 1/W_1[1+D_e^2/W_3^2] produces a better fit.  This fit does not improve with addition of a term linear in D_e. The best-fit parameters for this quadratic model are W_1=420±10 s and W_3 =  6.8±0.8×10^18 cm^-2. The value of W_3 may be interpreted as the irradiation dose which decreases T_2 by 2× relative to the un-irradiated value T_2 = W_1. 

Given that W_3 ∼ 400 ×⟨ D_e^3dB⟩, where ⟨ D_e^3dB⟩=1.6± 0.8× 10^16 cm^-2 is the dose that converts half the initial nitrogen concentration to NV centers, the data suggest any desired NV conversion (i.e. [NV]/[N^T]) can be obtained with negligible decrease in T_2 from irradiation-induced defects. The failure of the linear model at the highest doses suggests the dominant source of decoherence observed at high doses is not the negatively charged vacancy. Further study is required to isolate the responsible decoherence mechanism.



§ SUMMARY AND OUTLOOK


In this work, we presented a new room-temperature, all-optical technique for rapidly and quantitatively characterizing NV^--to-NV^0 charge state ratio. In particular, we isolated the emission spectra of NV^- and NV^0, allowing arbitrary diamond spectra to be decomposed into a linear combination of the NV^- and NV^0 spectral contributions. This photoluminescence decomposition analysis exhibits better reproducibility and lower noise than the preexisting Debye-Waller-factor approach and does not require multiple excitation wavelengths <cit.>. Furthermore, to enable the quantitative extraction of the [NV^-]/[NV^0] charge state ratio from these PL measurements, we also empirically determined for the first time the relative PL emission rate of NV^- compared to that of NV^0, finding κ_532=2.5±0.5 under low-intensity 532-nm excitation. 

In addition to providing an expedient method to optimize irradiation and annealing of NV-diamond material targeting high-performance quantum sensing applications, the PDA technique and analyses presented are expected to help uncover unknown or poorly understood NV-diamond behavior. For example, the value κ_λ can likely be used to place bounds on the branching ratios and quantum efficiency of the NV^0 ^2E ↔ ^2A transition. Further, extensions of the PDA technique may be used to measure NV^- and NV^0 saturation intensities, which are not precisely known and may vary for NV ensembles in different strain or charge environments where reported branching ratios differ <cit.>. Finally, extending the findings of this work—through additional study of vacancy creation, diffusion, and annealing as well as through further investigation of the mechanism for the NV^- T_2 spin coherence time's quadratic dependence on irradiation dose—is expected to provide useful insight into the defect content of different diamond materials. Using these techniques and analyses to develop a more comprehensive understanding of diamond material composition and the effects of processing thereon is critical for advancing quantum sensing based on ensembles of NV centers in diamond.

We gratefully acknowledge the authors of Ref. <cit.>, Hiromitsu Kato, Marco Wolfer, Christoph Schreyvogel, Michael Kunzer, Wolfgang Müller-Sebert, Harald Obloh, Satoshi Yamasaki, and Christoph Nebel for providing raw data allowing us to confirm the measured value of κ_532. We thank Daniel Twitchen for insights into literature variation in irradiation outcomes.








 §.§ Diamond characteristics



  §.§.§ Diamond nitrogen concentration estimation



The primary contributors to NV^- T_2 coherence times are ^13C, residual neutral nitrogen, and spin-lattice interactions (i.e. T_1-related processes). Consequently, Eqn. (<ref>) can be used to estimate the value of [N^T].  The data in Fig. <ref> suggest T_2 = 420 ± 10 s for unirradiated diamonds in our study.  From Ref. <cit.> we expect T_2,NV^-/N = 165 ± 15 ppm×s, where [N^T] is evaluated in ppm (i.e., a diamond with [N^T] = 1 ppm will have T_2,NV^-/N = 165 s, while a diamond with [N^T] = 2 ppm will have T_2,NV^-/N = 165/2 s ≈ 83 s). For the ^13C decoherence contribution, we use values reported in the literature for samples exhibiting natural isotopic abundance: Hall et al. theoretically predicted 780 s <cit.>; Mizuochi et al. measured 650 s for a single NV <cit.>; Balasubranian et al. <cit.> measured 1.8 ms for a NV in 0.3% ^13C diamond, which corresponds to a 504 s T_2 for natural abundance ^13C; and Stanwix et al. reported 600 s for NV ensembles in high-purity diamonds with a natural abundance of ^13C <cit.>.   Taking the average of these four values yields T_2,NV^-/^13C=634 ± 114 s for diamonds with natural abundance ^13C. We evaluate T_1 from Fig. <ref> which yields T_1 = 5.5 ± 0.5 ms and is in good agreement with data presented in Ref. <cit.>. Employing Eqn. (<ref>) with the above estimates suggests a mean nitrogen concentration [N^T] =118 ± 48 ppb for the diamonds studied here.

This estimate is consistent with measurements in Ref. <cit.> for standard-grade diamonds from Element Six. Test results of twenty diamonds of the same part number exhibit a mean nitrogen concentration of 120 ppb with a standard deviation of 16 ppb and minimum and maximum concentrations of 95 and 162 ppb respectively <cit.>. Since the nitrogen concentration strongly affects CVD diamond growth rate, optical absorption, and other diamond properties, this parameter is carefully controlled during diamond growth, and it is perhaps unsurprising that the diamond samples employed in the work presented here exhibit similar nitrogen concentrations. 





  §.§.§ Irradiation parameters



Diamonds were electron irradiated by a commercial supplier (Prism Gem, LLC). The electron energy of the Van de Graaff accelerator was fixed at 1 MeV. Doses of 1 × 10^15 e^-/cm^2 employ 2.5 mA beam current for 36 seconds.  Doses of 1 × 10^16 e^-/cm^2 and 1 × 10^17 e^-/cm^2 employ the same 2.5 mA beam current for 6 minutes and one hour respectively.  Finally, doses of 1 × 10^18 e^-/cm^2 and 1 × 10^19 e^-/cm^2 employ 25 mA of beam current for one hour and 10 hours respectively.  Irradiation doses are calculated by the formula D_e = I t/A q_e where I is the beam current, A is the scan area over which the electron beam is rastered, q_e is the electron charge, and t is the irradiation time. Table <ref> details irradiation doses for the 28 diamonds used in this study. 



 §.§ Experimental details



  §.§.§ Confocal microscope



A home-built confocal microscope, shown diagrammatically in Fig. <ref>, allows PL spectroscopy and coherence measurements of the diamond samples.

Green 532 nm laser light (Coherent Verdi G5) is sent through an acousto-optic modulator (Crystal Technology/Gooch & Housego 3200-147). The zeroth order beam is discarded while the first order beam  passes through a half-wave plate, reflects off a dichroic mirror, and is focused onto the diamond with a commercial microscope objective (Mitutoyo M Plan Apo 20X 378-804-3) with 0.42 NA, a f_obj=10 mm focal length, and a D_ba=8.4 mm back aperture. A permanent magnet mounted on a translation stage creates a static bias magnetic field of moderate strength 40-60 gauss aligned along one of the four diamond crystallographic directions, allowing individual addressing of the magnetic resonances associated with the NV sub-ensemble oriented along the bias field. This bias field alignment prevents T_2 measurements from being skewed by the effects of ^13C.   

The diamond is mounted over a hole in a 170-m-thick glass coverslip attached to a stainless steel ring (Thorlabs, RS1M). Microwaves are applied to the diamond from above using a loop antenna. For T_2 and T_1  measurements, the interrogation volume is located at the top of diamond, closest to the microwave loop to achieve maximum MW drive strength. For PL measurements, the interrogation volume is located near the bottom diamond surface to minimize spherical aberration effects from the diamond's n_d≈2.41 refractive index.

Light emitted from the NV-ensemble is collected by the 20X objective and subsequently coupled using a 200 mm focal length lens into 10-m-core multi-mode fiber (Thorlabs, M64L02) that serves as a pinhole.  This fiber is connected to either a spectrometer (Princeton Instruments, IsoPlane SCT 320) for photoluminescence measurements or to a fiber connected to an avalanche photodiode single photon counter (Excelitas SPCM-AQRH-FC). The spectrometer is both wavelength and intensity calibrated (Princeton Instruments Intelical). The intensity calibration removes instrumental artifacts in the PL spectra that arise from the differing spectral responses from various optical elements in the spectrometer (mirrors, diffraction grating, CCD camera, etc.). Counts registered by the single photon counter are recorded by a DAQ (National Instruments, PCIe-6251), which is in turn gated by a pulse generator (SpinCore Pulse Blaster ESR-PRO-500). The pulse generator is used to gate the optical excitation, microwave radiation, and single photon counter measurement.





  §.§.§ Intensity calculations and saturation


To calculate the intensity associated with a given excitation laser power, we model the intensity of the green laser beam focused by the objective as a Gaussian profile:

    I(r,z) = 2P/π w(z)^2exp(-2r^2/w(z)^2),

where P is the incident 532 nm power corrected for any transmission losses, w(z) = w_0 √(1+(z/z_R)^2) is the 1/e^2 intensity radius, z_R = π n_d w_0^2/λ is the Rayleigh range (adjusted for the diamond index of refraction n_d), r is the radial distance in the plane transverse to the beam, λ is the wavelength in free space, and z is the axial distance from the origin, defined as the center of the collection volume <cit.>. The 1/e^2 intensity diameter at the origin is 2w_0. The peak intensity I_0 is

    I_0 = 2P/π w_0^2.

A beam profiler (Newport LBP2-HR-VIS) measures the spot size of the Gaussian laser beam right before the objective to be about 2.25 mm in diameter (1/e^2). For an input diameter D_in=2.25 mm, laser wavelength λ = 532 nm, and objective focal length f_obj = 10 mm, the focused excitation beam waist is w_0^ex = 2λ f_obj/π D_in= 1.5 m <cit.>. This yields a Rayleigh range z_R = 32 m for the excitation beam.

The saturation intensity of the NV^- center is on the order of 100 kW/cm^2 <cit.> with the saturation intensity of the NV^0 likely to be higher. For measurements relating collected PL to [NV^T], excitation intensities of 3.3 kW/cm^2 or less were employed.











  §.§.§ Excitation and collection volume


Interpretation of photoluminescence and related derived data is most straightforward if the interrogated NVs receive near-uniform optical excitation intensity. To ensure that the NVs in the collection volume are near-uniformly illuminated with 532 nm light, the collection volume is restricted to the central region of the excitation volume.  In contrast to the 532 nm excitation, the collected fluorescence makes full use of the objective 0.42 numerical aperture. As the objective back aperture D_ba is 8.4 mm, the collection beam waist is w_0^col = 2λ f_obj/π D_ba, which is w_0^col = 0.53m for 700 nm light. A radial slice of the excitation and collection profiles is shown in Fig. <ref>. 




 §.§ Charge state photoluminescence details



  §.§.§ Photoluminescence measurement protocol

PL spectroscopy measurements are performed with a spectrometer and an attached CCD camera (Princeton Instruments, PIXIS 100). A half wave plate tunes the incident light polarization so that all four NV orientations are excited equally. To avoid signal to noise degradation by modal noise (i.e. speckle) in the multi-mode fiber delivering light to the spectrometer, the fiber is taped to a speaker-head driven by a 200 Hz sinusoidal voltage <cit.>.  Multiple PL measurements are averaged together and normalized by the exposure time. The homogeneity of the commercial diamonds employed in this study was characterized by measuring PL spectra in 5 different spots in each of a subset of 8 diamonds; the PL intensity and extracted NV charge state ratios were found to vary by a maximum of 10% across each diamond sample, indicating reasonable homogeneity within each diamond.



  §.§.§ The Huang-Rhys and Debye-Waller factors


A commonly used method to approximate the charge state efficiency is the Debye-Waller (DW) factor S̃, which is defined in terms of the relative weight of the ZPL fluorescence:

    S̃_X = -ln(I_ZPL,X/I_total,X),

where I_ZPL,X is the integrated fluorescence of the zero-phonon line of NV charge state X's PL spectra and I_total,X is the total integrated fluorescence of that charge state <cit.>. In the limit where the ZPL intensities can be measured with sufficient accuracy, the DW factor would enable accurate determination of the relative fluorescence ratio of NV^- and NV^0 for a given sample. However, the relative fluorescence ratio alone is insufficient to accurately determine the charge state efficiency, due to different dynamics (lifetimes, absorption cross sections, etc.,) of the NV^- and NV^0 systems. The Debye-Waller factor can be derived from a theoretical treatment of the electron-phonon coupling for defects in the diamond lattice by approximating the coupling to be linear in the nuclear displacement <cit.>.  



Frequently in the literature the Debye-Waller factor is confused with the Huang-Rhys (HR) factor. In Ref. <cit.>, Huang and Rhys defined the Huang-Rhys factor to be the average number of photons of a k-phonon transition emitted in an optical transition <cit.>. Here the ZPL is weighted with an additional pre-factor of the cube of the optical emission frequency, which decreases the value of the DW factor compared to the HR factor.  To avoid confusion, we follow the naming convention suggested by Ref. <cit.>, denoting the Huang-Rhys factor by S and the Debye-Waller factor by S̃.

During this study, we found determination of the charge state efficiency using the Debye-Waller decomposition to be non-ideal, since PL features from other defects can obscure the ZPL, making accurate determination of the DW (or HR) factor difficult. This is most problematic for the NV^- ZPL, which lies within the NV^0 phonon sideband emission spectra. Additionally, the ZPL prominence relative to total PL is small at room temperature for both NV^- and NV^0.  Due in part to both non-idealities, there exists substantial variation in the DW factors reported in the literature.  Measured values of S̃ for NV^- in the literature are shown in Table <ref>.  The DW factor for NV^0 is less-studied and reported to be S̃ = 3.3 at low temperature <cit.>. 

To calculate the fractional contributions of NV^0 and NV^- to total PL using the Debye-Waller factor, we fit and subtract a local linear background to the ZPL of each charge state and fit a Gaussian to the modified ZPL.  The ZPL areas I_ZPL,0 and I_ZPL,- for NV^0 and NV^- respectively are weighted with the DW factors reported in Ref. <cit.>: S̃_0 = 3.3 for NV^0 and S̃_-=4.3 for NV^-.  Using the form of ζ in Eq. <ref> with the weighted ZPL areas (normalized by the charge state PL emission rate κ_532) gives the charge state efficiency as determined with the Debye-Waller decomposition, ζ^DW:

    ζ^DW = I_ZPL,-e^S̃_-/κ_532I_ZPL,0e^S̃_0+I_ZPL,-e^S̃_-.


The Debye-Waller decomposition provides fair performance for diamonds with low NV^0 fluorescence, which allow for easy extraction of the NV^- ZPL area <cit.>.  However, for diamonds with a high fraction of NV^0 fluorescence (in this study, diamonds with higher irradiation doses), determining the NV^- ZPL area becomes difficult. Figure <ref> shows the PL ratio (i.e. c_-/c_0) determined with the Debye-Waller decomposition for sample S for various laser excitation intensities. The increasingly large error bars and systematic errors at higher laser intensities reflects the difficulty of determining the NV^- ZPL area when the NV^0 photoluminescence dominates. Not only does the basis function decomposition lead to lower statistical errors when determining the charge state fluorescence ratio, but the method is particularly advantageous at the higher laser intensities where the Debye-Waller decomposition becomes unreliable. 




  §.§.§ Construction of basis functions



As the NV^0 ZPL lies outside the NV^- PL spectrum, the NV^- basis function can be created nearly free from NV^0. The little if any impurity is bounded in two ways and the basis function shape can be checked itself. First, by comparing peak deviations of the NV^- basis function around 575 nm, it is determined that the NV^- basis function contains <0.015 NV^0. Second, fitting a Lorentzian to the vicinity of 575 nm with FWHM bounded to lie in the range [2,5] nm suggests at most the NV^- basis function contains 0.003 NV^0. Third, the extracted NV^- basis function is in good agreement with theoretical predictions <cit.>, derived from experimentally determined single phonon sideband spectra from the ^3A_2↔ ^3E transition (see main text Fig. <ref>). Given the good fit with theory and the lack of any discernable NV^0 ZPL in the NV^- basis functions, we conclude that any contamination of the NV^- basis function with NV^0 is on the 1% level or below. Residual fits to the Î_0^pre function suggest 0.064 ± 0.033 is NV^-. As this component is removed, we expect the any residual component of NV^- in the NV^0 basis function is estimated to occur at the 3% level or less.

Further note that the process for constructing basis functions described in the main text is robust to the choice of initial spectra for Î_0^pre(λ) and Î_-^pre(λ). For example, repeating the basis function construction process with different yet comparable spectra to those employed in the main text yields a maximum 3% difference in the extracted basis functions and a maximum 3% difference in the PL decomposition results from a subset of 8 samples. However, due to the difficulty in fully isolating and removing NV^- residual contribution from the NV^0 basis function, selecting initial spectra with maximal NV^0 and NV^- contributions for the Î_0^pre(λ) and Î_-^pre(λ) functions, respectively, is recommended.

 






  §.§.§ Secondary determination of κ_532


To obtain the total NV content for the irradiation model (Eqn. (<ref>)) or to determine the charge state efficiency, we need to correct for the relative PL emission rate of NV^- to NV^0. In the main text, we find κ_532 = 2.5± 0.5 by varying the [NV^-]/[NV^0] ratio using different well-below-saturation 532 nm optical excitation intensities. To independently check the measured value of κ_532, we exploit 3rd party data, where the [NV^-]/[NV^0] ratio is varied electrically. Specifically, we employ the data from Ref. <cit.> Fig. 4a where a bias voltage is used to vary the [NV^-]/[NV^0] ratio. More positive bias voltages skew the favored charge state towards NV^0 while more negative bias voltages skew the favored charge state towards NV^-. 

We denote the PL spectra with 0 Volts applied as f'_NV^-(λ) and the PL spectra with 35 Volts applied as f'_NV^0(λ), which correspond to when the diamond under study in Kato et al. is primarily in the NV^- and NV^0 charge states respectively. Note that these spectra are not area-normalized. From the raw data for Ref. <cit.> Fig. 4a provided by the authors, we calculate a preliminary charge state PL ratio 

    κ'_532 = [ ∫_λ f'_NV^-(λ) dλ/∫_λ f'_NV^0(λ) dλ]_[NV^0]=[NV^-].

The above value of κ_532' will need to be corrected since the voltage only partially influences the charge state. Straightforward evaluation gives κ_532' = 1.687.


We now describe how κ_532' is corrected to obtain the true value of κ_532.
We desire to obtain the true charge state PL ratio κ_532 = ∫_λ I_-(λ) dλ /∫_λ I_0(λ) dλ.  We first express the normalized spectra from Ref. <cit.> Fig. 4a as a linear combination of the normalized basis spectra, Î_0(λ) and Î_-(λ) (see Figs. <ref> and <ref>): 

    f̂'_NV^-(λ)    = αÎ_0(λ) + (1-α) Î_-(λ) 
    f̂'_NV^0(λ)    = (1-β) Î_0(λ) + βÎ_-(λ).

We find α = 0.0853 and β = 0.273. 
The non-normalized PL spectra can be written as

    f'_NV^-(λ)    = C [a_0 Î_0(λ) + a_-κ_532Î_-(λ) ]
    
    f'_NV^0(λ)    = C [b_0 Î_0(λ) + b_-κ_532Î_-(λ) ]

where C is an overall scaling constant and the constants a_0, a_-, b_0, and b_- satisfy

    a_0+a_- = b_0+b_-,

which should be interpreted as the requirement that the quantity [NV^-]+[NV^0] be independent of the applied voltage. From Eqns. <ref>, <ref> we can rewrite Eqn. (<ref>) as 

    κ_532' = a_0 + a_-κ_532/b_0 + b_-κ_532.

Comparison of Eqns. <ref>, <ref> with Eqns. <ref>, <ref> allows us to note

    α   = a_0/a_0 +a_-κ_532
    β   = b_0/b_0 +b_-κ_532

Solving Eqns. <ref>, <ref>, <ref>, <ref> yields κ_532 = 2.17. This value of κ should be taken as a lower limit on the true value of κ_532 at low 532 nm optical excitation intensity since low optical excitation intensity was not employed in the work of Ref. <cit.>.






∘ 4.5 MeV from Ref. <cit.>

♢ 4.5 MeV from Ref. <cit.>

 4.5 MeV from Ref. <cit.>

♢ 3 MeV Iakoubovskii2005 <cit.>

 3 MeV Iakoubovskii2005 <cit.>

3 MeV Collins2003 <cit.>

△3 MeV Collins2003 <cit.>

1 MeV Prism Gem 2017 (unpublished)

∘ 2 MeV Twitchen 1999 Physica B <cit.>

▾ 1.9 MeV Twitchen 1999 DRM ] <cit.>

▽ 1.9 MeV Twitchen 1999 DRM <cit.>

 2 MeV Newton 2002 <cit.>



◂ 1.9 MeV Lawson 1998 (using Twitchen1999DRM fGR1+fND1) <cit.> using fGR1+fND1 from Ref. <cit.>

 2 MeV Allers 1998 (using Twitchen1999DRM fGR1+fND1) <cit.> using fGR1+fND1 from Ref. <cit.>






 §.§ Annealing and nitrogen conversion study details



  §.§.§ Irradiation studies in the literature



This section reviews prior literature studies reporting monovacancy creation rates by electron irradiation.  Figure <ref> displays vacancies measured in the diamond literature via UV-Vis absorption spectrophotometry for varying irradiation doses <cit.>, alongside simulated vacancy creation rates from Ref. <cit.>. In type IIa diamonds with low concentrations of electron donors, irradiation primarily creates neutral monovacancies (V^0, GR1), whereas in diamonds with a moderate to high concentration of electron donors, such as type Ib diamonds, negative monovacancies (V^-, ND1) are created at similar or greater rates than neutral monovacancies <cit.>. Fig. <ref> depicts only the V^0 concentration for the former sample type and the total vacancy concentration (V^0 + V^-) for the latter type. 

The GR1 center displays a zero-phonon absorption doublet at 740.9 nm and 744.4 nm <cit.>, while the ND1 center displays a zero-phonon line at 393.6 nm <cit.>. For the data depicted in Fig. <ref>, vacancy concentration is calculated by integrating under the absorption features to obtain A_GR1 and A_ND1, the absorption strengths in meV cm^-1, and then solving the equations A_GR1 = f_GR1[V^0] and A_ND1 = f_ND1[V^-], where f represents the transition’s oscillator strength <cit.>. The authors of Ref. <cit.> measured the ratio f_ND1/f_GR1 = 4.0 using detailed balance arguments applied to an annealing study. Later, the authors of Ref. <cit.> determined that f_ND1 = 4.8(2) × 10^-16 meV cm^2 by correlating ND1 absorption measurements with well-calibrated EPR measurements of [V^-]. Although V^0 is spinless and thus EPR inactive, f_GR1 is determined to be 1.2(3) × 10^-16 meV cm^2 by dividing f_ND1 by f_ND1/f_GR1 = 4.0. All vacancy concentrations reported in Fig. <ref> use these calibration constants. 

Large variation is apparent between data sets from different publications, and, more specifically, between different irradiation facilities. The diamonds from Refs. <cit.> were all irradiated at Reading University in the UK and all show significantly higher vacancy creation rates than diamonds irradiated at other facilities. Strikingly, Refs. <cit.> and <cit.>, published by the same group within the same year using approximately the same electron energy (≈2 MeV), report vacancy creation rates of 0.50 cm^-1 and 1.53 cm^-1 respectively, a difference of 3×. This discrepancy leaves open the possibility of a calibration error in the applied dosage at one or more irradiation facilities, or a variable which is not controlled for. 

Furthermore, although Ref. <cit.> observed the vacancy creation rate to be independent of temperature up to 600 K, such temperatures could be exceeded depending on the irradiation setup; inadequate heat sinking combined with high beam current could cause diamonds to be heated above the temperature at which vacancies become mobile and recombine with interstitials, thereby reducing the measured post-irradiation vacancy concentration for a given dose. The diamonds represented in Fig. <ref> are irradiated with beam currents ranging ranging from A <cit.> to tens of mA <cit.>.  In Ref. <cit.>, the vacancies created per dose is seen to vary by more than 2× between different diamonds. It is unclear whether this variation arises from differences in the samples <cit.>, or inconsistency in the irradiation conditions (such as electron beam inhomogeneity). It is possible that more than one of these surmised explanations combine to describe the variation depicted in Fig. <ref>. 



Moreover, as pointed out in Ref. <cit.>, there exists disagreement in the literature on the value of f_ND1/f_GR1. In particular, absorption measurements employing reversible charge interconversion between GR1 and ND1 find the ratio of oscillator strengths f_ND1/f_GR1 to range from 2 to 10 <cit.>. This variation suggests a possible uncertainty in the y-axis of Fig. <ref> of up to ∼ 5×. Such uncertainty complicates the determination of the required irradiation dose to create a fixed vacancy concentration, and thus generate a desired NV concentration after annealing.

We further note that all measured data points lie below the predicted relationship between dose and vacancy concentration from Ref. <cit.>, which uses SRIM (Stopping and Range of Ions in Matter) calculations <cit.>. These simulations ignore spontaneous recombination, which is crudely estimated to occur 30%-50% of the time for diamonds irradiated with ∼2 MeV electrons <cit.> and may occur at a higher rate for lower-energy electrons, where the vacancy-interstitial distance is reduced <cit.>. Uncertainty in the value of f_ND1/f_GR1 <cit.> could help account for the discrepancy between the simulations and much of reported data.

Indirect methods for determining the vacancy concentration produced by a given irradiation dose offer an independent alternative to UV-Vis absorption measurements. The indirect method presented in this paper finds a monovacancy creation rate of 0.52±0.26 cm^-1 for 1 MeV electron irradiation by Prism Gem. In comparison, direct vacancy concentration measurements by UV-Vis on five similar diamonds irradiated at the same facility (the  symbols in Fig. <ref>, excluding the point at 3× 10^18 e^-/cm^2), suggest a V^0 creation rate of 0.30± 0.05 cm^-1 and no production of V^- above the device's detection limit of [V^-]∼ 100 ppb. 
The indirect measurements detailed in the main text and these UV-Vis measurements are consistent to within the uncertainty in f_GR1, which demonstrates the need for additional detailed studies to reduce uncertainty and resolve the aforementioned discrepancies. 

The fit to Eqn. (<ref>) gives B = 0.52±0.26 cm^-1. 




  §.§.§ Diffusion


For one dimensional diffusion, the RMS distance σ from the starting point after N jumps of distance L is given by σ = L √(N). The diffusion coefficient is defined as D ≡σ^2/2t_ann where t_ann is the anneal time during which the diffusion occurs. Combining the two equations allows the diffusion coefficient to be rewritten as D = L^2Γ/2 where Γ≡ N/t_ann is the vacancy jump frequency. It can be shown that the distance between two adjacent carbon atoms in the diamond lattice is a√(3)/4 where a is the diamond lattice constant. Therefore, in one dimension L = a/4. Substituting in gives

    D = Γa^2/32,

which is consistent with alternative derivations <cit.>. We note that the diffusion constant is independent of the dimensionality so the above diffusion constant holds for three dimensional diffusion in the diamond lattice <cit.>. 

Comparison of samples with 1 hour post-irradiation anneal at 850 ^∘C with samples receiving a 1250 ^∘C post-irradiation anneal, suggests that 1250 ^∘C 1 hour anneal created approximately 4.5 times as many NVs as the 850 ^∘C anneal for irradiation doses ∼ 10^16 e^-/cm^2 (see Fig. <ref>). Under the likely assumption that the 1 hour 1250 ^∘C anneal proceeded to completion (no monovacancies are left), this finding suggests that employing only the 1 hour 850 ^∘C anneal limits NV^T to χ≈ 1/4.5 of that expected from the introduced vacancy concentration alone. From the observed value of χ we can estimate the diffusion constant D as we now describe. 

We assume vacancy diffusion can be modeled as a random walk on the diamond lattice and ignore sinks for vacancies other than substitutional nitrogen. We do not correct for any time the monovacancy spends in the negative charge state which is immobile at 850 ^∘ C <cit.>. The number of unique lattice sites a vacancy is expected to visit in order to create an NV is 

    M_need =[C]/4[N]F_c,

where [C]=1.76×10^17 cm^-3 is the density of carbon atoms in diamond and the factor 4 accounts for the four lattice sites adjacent to a substitutional nitrogen, each of which will result in NV formation if visited by a vacancy. The number of unique sites visited is 

    M_visited = Γ t_annη_uni = 32 D t_annη_uni/a^2

where η_uni characterizes the ratio of unique sites visited to discrete steps during the random walk. Ref. <cit.> finds the fraction of unique sites visited for a random walk on a diamond lattice to be n_uni≈ 0.56 n + 0.63√(n) where n is the number of steps. This result is similar to the number of distinct sites visited after n steps for other lattices where: s.c. n_unique≈ 0.66 + 0.57√(n), b.c.c. n_unique≈ 0.72 + 0.52√(n), and f.c.c. n_unique≈ 0.74 + 0.52√(n) <cit.>. 

We assume the conversion efficiency of the annealing follows

    χ = 1 - e^-M_visited/M_need.

Combining Eqns. <ref>, <ref>, <ref> allows the diffusion constant to be written as

    D = -[C] a^2 Log[1-χ]/32 t_ann 4 [N] F_c η_uni

where χ≈ 1/4.5 is the observed annealing efficiency for the 1 hour 850 ^∘C anneal. We obtain D = 1.8 nm^2/s for monovacancies at 850 ^∘C, subject to the assumptions of our diffusion model. This value should be compared to the vacancy diffusion constant determined in Ref. <cit.> which observes D = 1.1 nm^2/s for annealing at 750 ^∘C and Ref. <cit.> which observes the upper limit of D ≤ 40 nm^2/s for annealing at 1050 ^∘C. 




  §.§.§ Comparison with crude diffusion theory


A neutral monovacancy in diamond is surrounded by four carbon atoms. The average time for one of the carbon atoms to jump into the vacancy is given by <cit.>

    τ_carbon = e^E_a/k_B T/ν_0

where E_a = 2.3 ± 0.3 eV <cit.> is the activation energy, k_B is the Boltzmann constant, T is the temperature, and ν_0 is the attempt frequency for an individual carbon atom. As there are four carbons, the average time for the vacancy to move one lattice site is 

    τ_vacancy = e^E_a/k_B T/4ν_0

so that we can write the vacancy jumping frequency as

    Γ = 4ν_0 e^-E_a/k_B T.

Plugging the above into Eqn. (<ref>) gives

    D=1/8ν_0 a^2  e^-E_a/k_B T.

Although the above equation is approximately correct, accurate values of the attempt frequency ν_0 are difficult to determine a priori <cit.>. In particular, Ref. <cit.> notes, “For want of a suitable value, either the Debye frequency or Einstein frequency can be used but there is no fundamental reason why either is applicable”. Nevertheless, employing the diamond Debye frequency of 40 THz for ν_0 results in D = 30 nm^2/s, much higher than the experimental estimate of 1.8 nm^2/s. We note that the uncertainty on E_a = 2.3±0.3 results in a roughly 500× uncertainty in D at 850 ^∘C, so the discrepancy between the experimentally estimated D and that calculated from theory cannot be completely attributed to uncertainty in the value of ν_0.  




 §.§ T_2 and T_1 coherence study details




  §.§.§ T_2 measurement protocol

The T_2 coherence time for each diamond is measured as follows: A moderate magnetic field of 40-60 gauss is applied along the <111> diamond crystallographic axis, allowing NV centers of only a single orientation to be addressed. NV fluorescence collected from a π/2-π-π/2 spin-echo sequence and a subsequent π/2-π-3π/2 sequence is subtracted and normalized to eliminate common-mode laser intensity noise and other slow drifts.  To avoid anisotropic-hyperfine effects between the NV spin and the ^13C nuclei <cit.>, special care is taken to minimize misalignment of the static magnetic field with the NV symmetry axis.  To achieve this condition, the bias magnetic field is aligned along one of the four diamond crystallographic axes, and only NVs parallel to this direction are addressed with resonant microwave (MW) radiation. The envelope of the spin-echo signal is modeled as a stretched exponential e^-(τ/T_2)^p with stretched exponential parameter p. While the collapses and revivals in the fluorescence signal due to the ^13C nuclei are theoretically described by a fourth-power sinusoid <cit.> for a single NV spin, we instead used a linear combination of Lorentzian and Gaussian shapes for each revival to modulate the stretched exponential envelope.  The curve fits for our diamonds yielded 0.5 < p < 2. From this fit we extract T_2.



  §.§.§ T_1 measurement protocol

The longitudinal relaxation time T_1 is measured as follows: NVs are first optically initialized into the m=0 state.  Thereafter a near-resonant MW π pulse transfers the population to the m=1 state, after which the population relaxes for a variable time τ toward the thermal equilibrium mixed state before being read out with a second optical pulse.  This yields the fluorescence trace m_1. The NV ensemble is optically initialized   and again allowed to relax again to thermal equilibrium with no MW pulse, giving fluorescence trace m_0. The population difference m_0 - m_1 is fit to an exponential decay with decay time T_1.  The subtraction of the signals eliminates common-mode noise from laser fluctuations during the measurement and the common fluorescence from the NV orientations not addressed by the MW π pulse. 


